ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/FastOpenGlDisplayer/src/OpenGLDisplayer.cc
Revision: 1.10
Committed: Fri Mar 21 19:40:46 2008 UTC (17 years, 1 month ago) by querten
Content type: text/plain
Branch: MAIN
Changes since 1.9: +3 -21 lines
Log Message:
Include Geometry of Calorimter

File Contents

# Content
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 // $Id: OpenGLDisplayer.cc,v 1.2 2008/03/15 17:02:20 querten Exp $
17 //
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 #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 #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 #include "DataFormats/GeometrySurface/interface/TrapezoidalPlaneBounds.h"
76 #include "DataFormats/GeometrySurface/interface/RectangularPlaneBounds.h"
77
78 #include "Visualisation/OpenGLDisplayer/interface/SimEvent.h"
79 #include "Visualisation/OpenGLDisplayer/interface/Geometry.h"
80
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 virtual void beginJob(const edm::EventSetup& iSetup);
96 virtual void analyze (const edm::Event&, const edm::EventSetup&);
97 virtual void endJob ();
98
99
100 MySimEvents* MyEvents;
101 Geometry* Geom;
102
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 OpenGLDisplayer::beginJob(const edm::EventSetup& iSetup)
130 {
131 DetId Detid;
132 int SubDet;
133
134 MyEvents = new MySimEvents();
135 Geom = new Geometry();
136
137 edm::ESHandle<TrackerGeometry> tkGeom;
138 iSetup.get<TrackerDigiGeometryRecord>().get( tkGeom );
139 vector<GeomDet*> TrackerDets = tkGeom->dets();
140
141 for(unsigned int i=0;i<TrackerDets.size();i++){
142 Detid = TrackerDets[i]->geographicalId();
143 // SubDet = Detid.subdetId();
144
145 GeomDet* DetUnit = TrackerDets[i];
146 if(!DetUnit)continue;
147 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
151 float width = 0;
152 float length = 0;
153 float thickness = 0;
154 float TrapezoidalParam = 0;
155
156 if(trapezoidalBounds)
157 {
158 std::vector<float> const & parameters = (*trapezoidalBounds).parameters();
159 width = parameters[0]*2;
160 length = parameters[3]*2;
161 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
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 GlobalVector Pos = GlobalVector(DetUnit->position().basicVector());
175
176 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
183 }
184
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) 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 printf("Det = %i Subdet = %i --> %6.2f %6.2f %6.2f - %i\n", Detid.det(), SubDet, CellPos.x(), CellPos.y(), CellPos.z(), CellCorners.size());
208 }
209
210
211 Geom->Save("Tracker.geom");
212
213
214 }
215
216 // ------------ method called once each job just after ending the event loop ------------
217 void
218 OpenGLDisplayer::endJob() {
219
220 MyEvents->Save((char*) OutputFile.c_str());
221 MyEvents->Load((char*) OutputFile.c_str());
222 }
223
224
225
226 //
227 // member functions
228 //
229
230 // ------------ method called to for each event ------------
231 void
232 OpenGLDisplayer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup)
233 {
234
235 edm::Handle<std::vector< SimTrack > > h_SimTracks;
236 iEvent.getByLabel("g4SimHits", h_SimTracks);
237 std::vector< SimTrack > SimTrackColl = *h_SimTracks.product();
238
239 edm::Handle<std::vector< SimVertex > > h_Vertex;
240 iEvent.getByLabel("g4SimHits", h_Vertex);
241 std::vector< SimVertex > VertexColl = *h_Vertex.product();
242
243
244
245 // access the tracker
246 edm::ESHandle<TrackerGeometry> theTrackerGeometry;
247 iSetup.get<TrackerDigiGeometryRecord>().get(theTrackerGeometry);
248 const TrackerGeometry& theTracker(*theTrackerGeometry);
249
250 // the DetUnits
251 TrackingGeometry::DetContainer theDetUnits = theTracker.dets();
252
253
254
255
256 // Save the data
257
258 MySimEvent* MyEvent = new MySimEvent;
259
260 for ( unsigned int a = 0; a < SimTrackColl.size(); ++a ) {
261 MySimTrack MysimTrack;
262 SimTrack simTrack =SimTrackColl[a];
263
264 MysimTrack.track_id =simTrack.trackId();
265 MysimTrack.Type =simTrack.type();
266 MysimTrack.parent_vertex=simTrack.vertIndex();
267 MysimTrack.Px =simTrack.momentum().x();
268 MysimTrack.Py =simTrack.momentum().y();
269 MysimTrack.Pz =simTrack.momentum().z();
270 MysimTrack.E =simTrack.momentum().e();
271 MysimTrack.charge =simTrack.charge();
272
273 MyEvent->MySimTrackCollection.push_back(MysimTrack);
274 }
275
276
277 for (unsigned int b = 0; b < VertexColl.size(); ++b ) {
278 MySimVertex MyVertex;
279 SimVertex Vertex = VertexColl[b];
280
281 MyVertex.parentTrack_id =Vertex.parentIndex ();
282 MyVertex.x =Vertex.position().x();
283 MyVertex.y =Vertex.position().y();
284 MyVertex.z =Vertex.position().z();
285
286 MyEvent->MySimVertexCollection.push_back(MyVertex);
287 }
288
289 for(unsigned int i=0; i<SimHitSubdetectors.size(); i++)
290 {
291 edm::Handle<std::vector< PSimHit > > h_Hits;
292 iEvent.getByLabel("g4SimHits",SimHitSubdetectors[i].c_str(), h_Hits);
293 std::vector< PSimHit > Hits = *h_Hits.product();
294
295 for(unsigned int h=0; h<Hits.size(); h++)
296 {
297 DetId theDetUnitId(Hits[h].detUnitId());
298 const GeomDet * theDet = theTracker.idToDet(theDetUnitId);
299
300 MyPSimHit Hit;
301 Hit.x = theDet->surface().toGlobal(Hits[h].localPosition()).x();
302 Hit.y = theDet->surface().toGlobal(Hits[h].localPosition()).y();
303 Hit.z = theDet->surface().toGlobal(Hits[h].localPosition()).z();
304 Hit.ProcessType = Hits[h].processType();
305 Hit.dEdX = Hits[h].energyLoss();
306
307 MyEvent->MyPSimHitCollection.push_back(Hit);
308
309 }
310 }
311
312 edm::Handle<std::vector< reco::Track > > h_Tracks;
313 iEvent.getByLabel("ctfWithMaterialTracks", h_Tracks);
314 std::vector< reco::Track > TrackColl = *h_Tracks.product();
315
316 for ( unsigned int t = 0; t < TrackColl.size(); ++t ) {
317 MyRecoTrack MyrecoTrack;
318 reco::Track recoTrack =TrackColl[t];
319
320 // MyrecoTrack.Hits = new MyRecoHit[MyrecoTrack.N];
321
322 for(unsigned int h=0;h<recoTrack.recHitsSize();h++){
323 TrackingRecHitRef h_it = recoTrack.recHit(h);
324 if(!h_it->isValid() )continue;
325 DetId detId = h_it->geographicalId();
326 const GeomDet * theDet = theTracker.idToDet(detId);
327 LocalPoint localPos = h_it->localPosition();
328
329 MyRecoHit hit;
330 hit.x = theDet->surface().toGlobal(localPos).x();
331 hit.y = theDet->surface().toGlobal(localPos).y();
332 hit.z = theDet->surface().toGlobal(localPos).z();
333 hit.DetId = detId.rawId();
334 hit.Charge = -1;
335 MyrecoTrack.Hits.push_back(hit);
336
337
338 // printf("%8.2f %8.2f %8.2f\n",theDet->surface().toGlobal(localPos).x(),theDet->surface().toGlobal(localPos).y(),theDet->surface().toGlobal(localPos).z());
339 }
340 MyEvent->MyRecoTrackCollection.push_back(MyrecoTrack);
341 }
342
343 MyEvents->Events.push_back(MyEvent);
344
345 }
346
347
348 //define this as a plug-in
349 DEFINE_FWK_MODULE(OpenGLDisplayer);
350
351