ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/FastOpenGlDisplayer/src/OpenGLDisplayer.cc
Revision: 1.9
Committed: Fri Mar 21 15:11:29 2008 UTC (17 years, 1 month ago) by querten
Content type: text/plain
Branch: MAIN
Changes since 1.8: +67 -9 lines
Log Message:
Add Calo Geometry

File Contents

# User Rev Content
1 querten 1.1 // -*- C++ -*-
2     //
3     // Package: OpenGLDisplayer
4     // Class: OpenGLDisplayer
5     //
6     /**\class OpenGLDisplayer OpenGLDisplayer.cc Visualisation/OpenGLDisplayer/src/OpenGLDisplayer.cc
7    
8     Description: <one line class summary>
9    
10     Implementation:
11     <Notes on implementation>
12     */
13     //
14     // Original Author: Loic QUERTENMONT
15     // Created: Fri Oct 26 07:22:12 CEST 2007
16 querten 1.3 // $Id: OpenGLDisplayer.cc,v 1.2 2008/03/15 17:02:20 querten Exp $
17 querten 1.1 //
18     //
19    
20    
21     // system include files
22     #include <memory>
23    
24     // user include files
25     #include "FWCore/Framework/interface/Frameworkfwd.h"
26     #include "FWCore/Framework/interface/EDAnalyzer.h"
27    
28     #include "FWCore/Framework/interface/Event.h"
29     #include "FWCore/Framework/interface/MakerMacros.h"
30    
31     #include "FWCore/ParameterSet/interface/ParameterSet.h"
32     #include "FWCore/ServiceRegistry/interface/Service.h"
33     #include "FWCore/Framework/interface/ESHandle.h"
34    
35     #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
36     #include "Geometry/CommonDetUnit/interface/GeomDetUnit.h"
37     #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
38     #include "Geometry/TrackerNumberingBuilder/interface/GeometricDet.h"
39     #include "Geometry/CommonTopologies/interface/PixelTopology.h"
40     #include "Geometry/CommonTopologies/interface/StripTopology.h"
41     #include "Geometry/TrackerGeometryBuilder/interface/PixelGeomDetType.h"
42     #include "Geometry/TrackerGeometryBuilder/interface/StripGeomDetType.h"
43     #include "Geometry/TrackerGeometryBuilder/interface/PixelGeomDetUnit.h"
44     #include "Geometry/TrackerGeometryBuilder/interface/StripGeomDetUnit.h"
45     #include "DataFormats/GeometrySurface/interface/BoundSurface.h"
46     #include "DataFormats/DetId/interface/DetId.h"
47    
48 querten 1.9 #include "Geometry/CaloEventSetup/interface/CaloTopologyRecord.h"
49     #include "Geometry/Records/interface/IdealGeometryRecord.h"
50     #include "Geometry/CaloGeometry/interface/CaloSubdetectorGeometry.h"
51     #include "Geometry/CaloGeometry/interface/CaloGeometry.h"
52     #include "Geometry/EcalBarrelAlgo/interface/EcalBarrelGeometry.h"
53     #include "Geometry/EcalEndcapAlgo/interface/EcalEndcapGeometry.h"
54     #include "Geometry/CaloGeometry/interface/CaloCellGeometry.h"
55     #include "DataFormats/EcalDetId/interface/EcalSubdetector.h"
56    
57    
58    
59 querten 1.1 #include "SimDataFormats/TrackingHit/interface/PSimHitContainer.h"
60     #include "SimDataFormats/TrackingHit/interface/PSimHit.h"
61     #include "SimDataFormats/CrossingFrame/interface/CrossingFrame.h"
62     #include "SimDataFormats/CrossingFrame/interface/MixCollection.h"
63     #include "SimDataFormats/CaloHit/interface/PCaloHit.h"
64     #include "SimDataFormats/Track/interface/SimTrack.h"
65     #include "SimDataFormats/Vertex/interface/SimVertex.h"
66    
67     #include "DataFormats/Provenance/interface/BranchDescription.h"
68     #include "DataFormats/Provenance/interface/Provenance.h"
69     #include "DataFormats/Candidate/interface/Candidate.h"
70     #include "DataFormats/HepMCCandidate/interface/GenParticleCandidate.h"
71     #include "DataFormats/RecoCandidate/interface/RecoChargedCandidate.h"
72     #include "DataFormats/TrackReco/interface/Track.h"
73    
74     #include "DataFormats/Math/interface/Point3D.h"
75 querten 1.7 #include "DataFormats/GeometrySurface/interface/TrapezoidalPlaneBounds.h"
76     #include "DataFormats/GeometrySurface/interface/RectangularPlaneBounds.h"
77 querten 1.1
78     #include "Visualisation/OpenGLDisplayer/interface/SimEvent.h"
79 querten 1.3 #include "Visualisation/OpenGLDisplayer/interface/Geometry.h"
80 querten 1.1
81     using namespace edm;
82     using namespace std;
83    
84     //
85     // class decleration
86     //
87    
88     class OpenGLDisplayer : public edm::EDAnalyzer {
89     public:
90     explicit OpenGLDisplayer(const edm::ParameterSet&);
91     ~OpenGLDisplayer();
92    
93    
94     private:
95 querten 1.3 virtual void beginJob(const edm::EventSetup& iSetup);
96     virtual void analyze (const edm::Event&, const edm::EventSetup&);
97     virtual void endJob ();
98 querten 1.1
99    
100     MySimEvents* MyEvents;
101 querten 1.3 Geometry* Geom;
102 querten 1.1
103    
104     std::vector<std::string> SimHitSubdetectors;
105    
106     std::string OutputFile;
107    
108    
109     // ----------member data ---------------------------
110     };
111    
112     //
113     // constructors and destructor
114     //
115     OpenGLDisplayer::OpenGLDisplayer(const edm::ParameterSet& iConfig)
116     {
117     SimHitSubdetectors = iConfig.getParameter<std::vector<std::string> >("SimHitSubdetectors");
118    
119     OutputFile = iConfig.getParameter<std::string >("OutputFile");
120     }
121    
122    
123     OpenGLDisplayer::~OpenGLDisplayer()
124     {
125     }
126    
127     // ------------ method called once each job just before starting event loop ------------
128     void
129 querten 1.3 OpenGLDisplayer::beginJob(const edm::EventSetup& iSetup)
130 querten 1.1 {
131 querten 1.9 DetId Detid;
132     int SubDet;
133 querten 1.1
134 querten 1.3 MyEvents = new MySimEvents();
135     Geom = new Geometry();
136    
137     edm::ESHandle<TrackerGeometry> tkGeom;
138     iSetup.get<TrackerDigiGeometryRecord>().get( tkGeom );
139 querten 1.9 vector<GeomDet*> TrackerDets = tkGeom->dets();
140 querten 1.3
141 querten 1.9 for(unsigned int i=0;i<TrackerDets.size();i++){
142     Detid = TrackerDets[i]->geographicalId();
143     // SubDet = Detid.subdetId();
144 querten 1.4
145 querten 1.9 GeomDet* DetUnit = TrackerDets[i];
146 querten 1.3 if(!DetUnit)continue;
147 querten 1.7 const BoundPlane plane = DetUnit->surface();
148     const TrapezoidalPlaneBounds* trapezoidalBounds( dynamic_cast<const TrapezoidalPlaneBounds*>(&(plane.bounds())));
149     const RectangularPlaneBounds* rectangularBounds( dynamic_cast<const RectangularPlaneBounds*>(&(plane.bounds())));
150 querten 1.3
151 querten 1.9 float width = 0;
152     float length = 0;
153     float thickness = 0;
154     float TrapezoidalParam = 0;
155 querten 1.7
156     if(trapezoidalBounds)
157     {
158     std::vector<float> const & parameters = (*trapezoidalBounds).parameters();
159 querten 1.8 width = parameters[0]*2;
160     length = parameters[3]*2;
161 querten 1.7 thickness = (*trapezoidalBounds).thickness();
162     TrapezoidalParam = parameters[1]/parameters[0];
163     }else if(rectangularBounds){
164     width = DetUnit->surface().bounds().width();
165     length = DetUnit->surface().bounds().length();
166     thickness = DetUnit->surface().bounds().thickness();
167     TrapezoidalParam = 1;
168     }
169 querten 1.4
170     Surface::GlobalPoint WidthVector = plane.toGlobal( LocalPoint(width/2, 0, 0) );
171     Surface::GlobalPoint LengthVector = plane.toGlobal( LocalPoint(0, length/2, 0) );
172     Surface::GlobalPoint ThickVector = plane.toGlobal( LocalPoint(0, 0, thickness/2) );
173    
174 querten 1.7 GlobalVector Pos = GlobalVector(DetUnit->position().basicVector());
175 querten 1.4
176 querten 1.7 Geom->Add_TrackerDet(Detid.rawId(), TrapezoidalParam,
177     Pos.x(), Pos.y(), Pos.z(),
178     WidthVector.x() -Pos.x(), WidthVector.y() -Pos.y(), WidthVector.z() -Pos.z(),
179     LengthVector.x()-Pos.x(), LengthVector.y()-Pos.y(), LengthVector.z()-Pos.z(),
180     ThickVector.x() -Pos.x(), ThickVector.y() -Pos.y(), ThickVector.z() -Pos.z());
181    
182 querten 1.3
183     }
184 querten 1.9
185    
186     // ### CALO GEOMETRY ###
187    
188     edm::ESHandle<CaloGeometry> CaloGeom;
189     iSetup.get<IdealGeometryRecord>().get( CaloGeom );
190     const vector<DetId> CaloDets = CaloGeom->getValidDetIds();
191    
192     for(unsigned int i=0;i<CaloDets.size();i++)
193     {
194     Detid = CaloDets[i];
195     SubDet = Detid.subdetId();
196    
197     if(! (Detid.det()==DetId::Ecal||Detid.det()==DetId::Hcal)) continue;
198    
199     //EcalBarrel=1, EcalEndcap=2, EcalPreshower=3, EcalTriggerTower=4, EcalLaserPnDiode=5
200     // if(SubDet!=EcalBarrel)continue;
201    
202     const CaloCellGeometry* CellGeom = CaloGeom->getGeometry(Detid);
203     GlobalPoint CellPos = CellGeom->getPosition();
204    
205     const CaloCellGeometry::CornersVec CellCorners = CellGeom->getCorners();
206    
207     float Coord[24];
208     for(unsigned int i=0;i<8;i++){
209     Coord[i*3+0] = CellCorners[i].x();
210     Coord[i*3+1] = CellCorners[i].y();
211     Coord[i*3+2] = CellCorners[i].z();
212     }
213    
214     Geom->Add_CaloDet(Detid.rawId(),
215     CellPos.x(), CellPos.y(), CellPos.z(),
216     Coord[0 ] , Coord[1 ] , Coord[2 ] ,
217     Coord[3 ] , Coord[4 ] , Coord[5 ] ,
218     Coord[6 ] , Coord[7 ] , Coord[8 ] ,
219     Coord[9 ] , Coord[10] , Coord[11] ,
220     Coord[12] , Coord[13] , Coord[14] ,
221     Coord[15] , Coord[16] , Coord[17] ,
222     Coord[18] , Coord[19] , Coord[20] ,
223     Coord[21] , Coord[22] , Coord[23] );
224    
225     printf("Det = %i Subdet = %i --> %6.2f %6.2f %6.2f - %i\n", Detid.det(), SubDet, CellPos.x(), CellPos.y(), CellPos.z(), CellCorners.size());
226     }
227    
228    
229 querten 1.3 Geom->Save("Tracker.geom");
230 querten 1.1
231    
232     }
233    
234     // ------------ method called once each job just after ending the event loop ------------
235     void
236     OpenGLDisplayer::endJob() {
237    
238     MyEvents->Save((char*) OutputFile.c_str());
239     MyEvents->Load((char*) OutputFile.c_str());
240     }
241    
242    
243    
244     //
245     // member functions
246     //
247    
248     // ------------ method called to for each event ------------
249     void
250     OpenGLDisplayer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup)
251     {
252    
253     edm::Handle<std::vector< SimTrack > > h_SimTracks;
254     iEvent.getByLabel("g4SimHits", h_SimTracks);
255     std::vector< SimTrack > SimTrackColl = *h_SimTracks.product();
256    
257     edm::Handle<std::vector< SimVertex > > h_Vertex;
258     iEvent.getByLabel("g4SimHits", h_Vertex);
259     std::vector< SimVertex > VertexColl = *h_Vertex.product();
260    
261    
262    
263     // access the tracker
264     edm::ESHandle<TrackerGeometry> theTrackerGeometry;
265     iSetup.get<TrackerDigiGeometryRecord>().get(theTrackerGeometry);
266     const TrackerGeometry& theTracker(*theTrackerGeometry);
267    
268     // the DetUnits
269     TrackingGeometry::DetContainer theDetUnits = theTracker.dets();
270    
271    
272    
273    
274     // Save the data
275    
276     MySimEvent* MyEvent = new MySimEvent;
277    
278     for ( unsigned int a = 0; a < SimTrackColl.size(); ++a ) {
279     MySimTrack MysimTrack;
280     SimTrack simTrack =SimTrackColl[a];
281    
282     MysimTrack.track_id =simTrack.trackId();
283     MysimTrack.Type =simTrack.type();
284     MysimTrack.parent_vertex=simTrack.vertIndex();
285     MysimTrack.Px =simTrack.momentum().x();
286     MysimTrack.Py =simTrack.momentum().y();
287     MysimTrack.Pz =simTrack.momentum().z();
288     MysimTrack.E =simTrack.momentum().e();
289     MysimTrack.charge =simTrack.charge();
290    
291     MyEvent->MySimTrackCollection.push_back(MysimTrack);
292     }
293    
294    
295     for (unsigned int b = 0; b < VertexColl.size(); ++b ) {
296     MySimVertex MyVertex;
297     SimVertex Vertex = VertexColl[b];
298    
299     MyVertex.parentTrack_id =Vertex.parentIndex ();
300     MyVertex.x =Vertex.position().x();
301     MyVertex.y =Vertex.position().y();
302     MyVertex.z =Vertex.position().z();
303    
304     MyEvent->MySimVertexCollection.push_back(MyVertex);
305     }
306    
307     for(unsigned int i=0; i<SimHitSubdetectors.size(); i++)
308     {
309     edm::Handle<std::vector< PSimHit > > h_Hits;
310     iEvent.getByLabel("g4SimHits",SimHitSubdetectors[i].c_str(), h_Hits);
311     std::vector< PSimHit > Hits = *h_Hits.product();
312    
313     for(unsigned int h=0; h<Hits.size(); h++)
314     {
315     DetId theDetUnitId(Hits[h].detUnitId());
316     const GeomDet * theDet = theTracker.idToDet(theDetUnitId);
317    
318     MyPSimHit Hit;
319     Hit.x = theDet->surface().toGlobal(Hits[h].localPosition()).x();
320     Hit.y = theDet->surface().toGlobal(Hits[h].localPosition()).y();
321     Hit.z = theDet->surface().toGlobal(Hits[h].localPosition()).z();
322     Hit.ProcessType = Hits[h].processType();
323     Hit.dEdX = Hits[h].energyLoss();
324    
325     MyEvent->MyPSimHitCollection.push_back(Hit);
326    
327     }
328     }
329    
330     edm::Handle<std::vector< reco::Track > > h_Tracks;
331     iEvent.getByLabel("ctfWithMaterialTracks", h_Tracks);
332     std::vector< reco::Track > TrackColl = *h_Tracks.product();
333    
334     for ( unsigned int t = 0; t < TrackColl.size(); ++t ) {
335     MyRecoTrack MyrecoTrack;
336     reco::Track recoTrack =TrackColl[t];
337    
338 querten 1.2 // MyrecoTrack.Hits = new MyRecoHit[MyrecoTrack.N];
339 querten 1.1
340 querten 1.2 for(unsigned int h=0;h<recoTrack.recHitsSize();h++){
341 querten 1.1 TrackingRecHitRef h_it = recoTrack.recHit(h);
342     if(!h_it->isValid() )continue;
343     DetId detId = h_it->geographicalId();
344     const GeomDet * theDet = theTracker.idToDet(detId);
345     LocalPoint localPos = h_it->localPosition();
346 querten 1.2
347     MyRecoHit hit;
348     hit.x = theDet->surface().toGlobal(localPos).x();
349     hit.y = theDet->surface().toGlobal(localPos).y();
350     hit.z = theDet->surface().toGlobal(localPos).z();
351     hit.DetId = detId.rawId();
352     hit.Charge = -1;
353     MyrecoTrack.Hits.push_back(hit);
354    
355    
356     // printf("%8.2f %8.2f %8.2f\n",theDet->surface().toGlobal(localPos).x(),theDet->surface().toGlobal(localPos).y(),theDet->surface().toGlobal(localPos).z());
357 querten 1.1 }
358     MyEvent->MyRecoTrackCollection.push_back(MyrecoTrack);
359     }
360    
361     MyEvents->Events.push_back(MyEvent);
362    
363     }
364    
365    
366     //define this as a plug-in
367     DEFINE_FWK_MODULE(OpenGLDisplayer);
368    
369