CMS 3D CMS Logo

TrackBuildingAnalyzer.cc
Go to the documentation of this file.
11 
16 #include <string>
17 #include "TMath.h"
18 
19 #include <iostream>
20 
22  : SeedPt(NULL)
23  , SeedEta(NULL)
24  , SeedPhi(NULL)
25  , SeedPhiVsEta(NULL)
26  , SeedTheta(NULL)
27  , SeedQ(NULL)
28  , SeedDxy(NULL)
29  , SeedDz(NULL)
30  , NumberOfRecHitsPerSeed(NULL)
31  , NumberOfRecHitsPerSeedVsPhiProfile(NULL)
32  , NumberOfRecHitsPerSeedVsEtaProfile(NULL)
33  , stoppingSource(NULL)
34  , stoppingSourceVSeta(NULL)
35  , stoppingSourceVSphi(NULL)
36 {
37 }
38 
40 {
41 }
42 
44 {
45 
46  // parameters from the configuration
47  std::string AlgoName = iConfig.getParameter<std::string>("AlgoName");
48  std::string MEFolderName = iConfig.getParameter<std::string>("FolderName");
49 
50  // std::cout << "[TrackBuildingAnalyzer::beginRun] AlgoName: " << AlgoName << std::endl;
51 
52  // use the AlgoName and Quality Name
53  std::string CatagoryName = AlgoName;
54 
55  // get binning from the configuration
56  int TrackPtBin = iConfig.getParameter<int>( "TrackPtBin");
57  double TrackPtMin = iConfig.getParameter<double>("TrackPtMin");
58  double TrackPtMax = iConfig.getParameter<double>("TrackPtMax");
59 
60  int PhiBin = iConfig.getParameter<int>( "PhiBin");
61  double PhiMin = iConfig.getParameter<double>("PhiMin");
62  double PhiMax = iConfig.getParameter<double>("PhiMax");
63 
64  int EtaBin = iConfig.getParameter<int>( "EtaBin");
65  double EtaMin = iConfig.getParameter<double>("EtaMin");
66  double EtaMax = iConfig.getParameter<double>("EtaMax");
67 
68  int ThetaBin = iConfig.getParameter<int>( "ThetaBin");
69  double ThetaMin = iConfig.getParameter<double>("ThetaMin");
70  double ThetaMax = iConfig.getParameter<double>("ThetaMax");
71 
72  int TrackQBin = iConfig.getParameter<int>( "TrackQBin");
73  double TrackQMin = iConfig.getParameter<double>("TrackQMin");
74  double TrackQMax = iConfig.getParameter<double>("TrackQMax");
75 
76  int SeedDxyBin = iConfig.getParameter<int>( "SeedDxyBin");
77  double SeedDxyMin = iConfig.getParameter<double>("SeedDxyMin");
78  double SeedDxyMax = iConfig.getParameter<double>("SeedDxyMax");
79 
80  int SeedDzBin = iConfig.getParameter<int>( "SeedDzBin");
81  double SeedDzMin = iConfig.getParameter<double>("SeedDzMin");
82  double SeedDzMax = iConfig.getParameter<double>("SeedDzMax");
83 
84  int SeedHitBin = iConfig.getParameter<int>( "SeedHitBin");
85  double SeedHitMin = iConfig.getParameter<double>("SeedHitMin");
86  double SeedHitMax = iConfig.getParameter<double>("SeedHitMax");
87 
88  int TCDxyBin = iConfig.getParameter<int>( "TCDxyBin");
89  double TCDxyMin = iConfig.getParameter<double>("TCDxyMin");
90  double TCDxyMax = iConfig.getParameter<double>("TCDxyMax");
91 
92  int TCDzBin = iConfig.getParameter<int>( "TCDzBin");
93  double TCDzMin = iConfig.getParameter<double>("TCDzMin");
94  double TCDzMax = iConfig.getParameter<double>("TCDzMax");
95 
96  int TCHitBin = iConfig.getParameter<int>( "TCHitBin");
97  double TCHitMin = iConfig.getParameter<double>("TCHitMin");
98  double TCHitMax = iConfig.getParameter<double>("TCHitMax");
99 
100 
101  edm::InputTag seedProducer = iConfig.getParameter<edm::InputTag>("SeedProducer");
102  edm::InputTag tcProducer = iConfig.getParameter<edm::InputTag>("TCProducer");
103 
104  doAllPlots = iConfig.getParameter<bool>("doAllPlots");
105  doAllSeedPlots = iConfig.getParameter<bool>("doSeedParameterHistos");
106  doTCPlots = iConfig.getParameter<bool>("doTrackCandHistos");
107  doAllTCPlots = iConfig.getParameter<bool>("doAllTrackCandHistos");
108  doPT = iConfig.getParameter<bool>("doSeedPTHisto");
109  doETA = iConfig.getParameter<bool>("doSeedETAHisto");
110  doPHI = iConfig.getParameter<bool>("doSeedPHIHisto");
111  doPHIVsETA = iConfig.getParameter<bool>("doSeedPHIVsETAHisto");
112  doTheta = iConfig.getParameter<bool>("doSeedThetaHisto");
113  doQ = iConfig.getParameter<bool>("doSeedQHisto");
114  doDxy = iConfig.getParameter<bool>("doSeedDxyHisto");
115  doDz = iConfig.getParameter<bool>("doSeedDzHisto");
116  doNRecHits = iConfig.getParameter<bool>("doSeedNRecHitsHisto");
117  doProfPHI = iConfig.getParameter<bool>("doSeedNVsPhiProf");
118  doProfETA = iConfig.getParameter<bool>("doSeedNVsEtaProf");
119  doStopSource = iConfig.getParameter<bool>("doStopSource");
120 
121  // if (doAllPlots){doAllSeedPlots=true; doTCPlots=true;}
122 
123  ibooker.setCurrentFolder(MEFolderName);
124 
125  // book the Seed histograms
126  // ---------------------------------------------------------------------------------//
127  // std::cout << "[TrackBuildingAnalyzer::beginRun] MEFolderName: " << MEFolderName << std::endl;
128  ibooker.setCurrentFolder(MEFolderName+"/TrackBuilding");
129 
130  if (doAllSeedPlots || doPT) {
131  histname = "SeedPt_"+seedProducer.label() + "_";
132  SeedPt = ibooker.book1D(histname+CatagoryName, histname+CatagoryName, TrackPtBin, TrackPtMin, TrackPtMax);
133  SeedPt->setAxisTitle("Seed p_{T} (GeV/c)", 1);
134  SeedPt->setAxisTitle("Number of Seeds", 2);
135  }
136 
137  if (doAllSeedPlots || doETA) {
138  histname = "SeedEta_"+seedProducer.label() + "_";
139  SeedEta = ibooker.book1D(histname+CatagoryName, histname+CatagoryName, EtaBin, EtaMin, EtaMax);
140  SeedEta->setAxisTitle("Seed #eta", 1);
141  SeedEta->setAxisTitle("Number of Seeds", 2);
142  }
143 
144  if (doAllSeedPlots || doPHI) {
145  histname = "SeedPhi_"+seedProducer.label() + "_";
146  SeedPhi = ibooker.book1D(histname+CatagoryName, histname+CatagoryName, PhiBin, PhiMin, PhiMax);
147  SeedPhi->setAxisTitle("Seed #phi", 1);
148  SeedPhi->setAxisTitle("Number of Seed", 2);
149  }
150 
151  if (doAllSeedPlots || doPHIVsETA) {
152  histname = "SeedPhiVsEta_"+seedProducer.label() + "_";
153  SeedPhiVsEta = ibooker.book2D(histname+CatagoryName, histname+CatagoryName, EtaBin, EtaMin, EtaMax, PhiBin, PhiMin, PhiMax);
154  SeedPhiVsEta->setAxisTitle("Seed #eta", 1);
155  SeedPhiVsEta->setAxisTitle("Seed #phi", 2);
156  }
157 
158  if (doAllSeedPlots || doTheta){
159  histname = "SeedTheta_"+seedProducer.label() + "_";
160  SeedTheta = ibooker.book1D(histname+CatagoryName, histname+CatagoryName, ThetaBin, ThetaMin, ThetaMax);
161  SeedTheta->setAxisTitle("Seed #theta", 1);
162  SeedTheta->setAxisTitle("Number of Seeds", 2);
163  }
164 
165  if (doAllSeedPlots || doQ){
166  histname = "SeedQ_"+seedProducer.label() + "_";
167  SeedQ = ibooker.book1D(histname+CatagoryName, histname+CatagoryName, TrackQBin, TrackQMin, TrackQMax);
168  SeedQ->setAxisTitle("Seed Charge", 1);
169  SeedQ->setAxisTitle("Number of Seeds",2);
170  }
171 
172  if (doAllSeedPlots || doDxy){
173  histname = "SeedDxy_"+seedProducer.label() + "_";
174  SeedDxy = ibooker.book1D(histname+CatagoryName, histname+CatagoryName, SeedDxyBin, SeedDxyMin, SeedDxyMax);
175  SeedDxy->setAxisTitle("Seed d_{xy} (cm)", 1);
176  SeedDxy->setAxisTitle("Number of Seeds",2);
177  }
178 
179  if (doAllSeedPlots || doDz){
180  histname = "SeedDz_"+seedProducer.label() + "_";
181  SeedDz = ibooker.book1D(histname+CatagoryName, histname+CatagoryName, SeedDzBin, SeedDzMin, SeedDzMax);
182  SeedDz->setAxisTitle("Seed d_{z} (cm)", 1);
183  SeedDz->setAxisTitle("Number of Seeds",2);
184  }
185 
186  if (doAllSeedPlots || doNRecHits){
187  histname = "NumberOfRecHitsPerSeed_"+seedProducer.label() + "_";
188  NumberOfRecHitsPerSeed = ibooker.book1D(histname+CatagoryName, histname+CatagoryName, SeedHitBin, SeedHitMin, SeedHitMax);
189  NumberOfRecHitsPerSeed->setAxisTitle("Number of RecHits per Seed", 1);
190  NumberOfRecHitsPerSeed->setAxisTitle("Number of Seeds",2);
191  }
192 
193  if (doAllSeedPlots || doProfPHI){
194  histname = "NumberOfRecHitsPerSeedVsPhiProfile_"+seedProducer.label() + "_";
195  NumberOfRecHitsPerSeedVsPhiProfile = ibooker.bookProfile(histname+CatagoryName, histname+CatagoryName, PhiBin, PhiMin, PhiMax, SeedHitBin, SeedHitMin, SeedHitMax,"s");
197  NumberOfRecHitsPerSeedVsPhiProfile->setAxisTitle("Number of RecHits of each Seed",2);
198  }
199 
200  if (doAllSeedPlots || doProfETA){
201  histname = "NumberOfRecHitsPerSeedVsEtaProfile_"+seedProducer.label() + "_";
202  NumberOfRecHitsPerSeedVsEtaProfile = ibooker.bookProfile(histname+CatagoryName, histname+CatagoryName, EtaBin, EtaMin, EtaMax, SeedHitBin, SeedHitMin, SeedHitMax,"s");
204  NumberOfRecHitsPerSeedVsEtaProfile->setAxisTitle("Number of RecHits of each Seed",2);
205  }
206 
207  if (doAllTCPlots || doStopSource) {
208  // DataFormats/TrackReco/interface/TrajectoryStopReasons.h
209  size_t StopReasonNameSize = sizeof(StopReasonName::StopReasonName)/sizeof(std::string);
210  if(StopReasonNameSize != static_cast<unsigned int>(StopReason::SIZE)) {
211  throw cms::Exception("Assert") << "StopReason::SIZE is " << static_cast<unsigned int>(StopReason::SIZE)
212  << " but StopReasonName's only for "
213  << StopReasonNameSize
214  << ". Please update DataFormats/TrackReco/interface/TrajectoryStopReasons.h.";
215  }
216 
217 
218  histname = "StoppingSource_"+seedProducer.label() + "_";
219  stoppingSource = ibooker.book1D(histname+CatagoryName,
220  histname+CatagoryName,
221  StopReasonNameSize,
222  0., double(StopReasonNameSize));
223  stoppingSource->setAxisTitle("stopping reason",1);
224  stoppingSource->setAxisTitle("Number of Tracks",2);
225 
226  histname = "StoppingSourceVSeta_"+seedProducer.label() + "_";
227  stoppingSourceVSeta = ibooker.book2D(histname+CatagoryName,
228  histname+CatagoryName,
229  EtaBin,
230  EtaMin,
231  EtaMax,
232  StopReasonNameSize,
233  0., double(StopReasonNameSize));
234  stoppingSourceVSeta->setAxisTitle("track #eta",1);
235  stoppingSourceVSeta->setAxisTitle("stopping reason",2);
236 
237  histname = "StoppingSourceVSphi_"+seedProducer.label() + "_";
238  stoppingSourceVSphi = ibooker.book2D(histname+CatagoryName,
239  histname+CatagoryName,
240  PhiBin,
241  PhiMin,
242  PhiMax,
243  StopReasonNameSize,
244  0., double(StopReasonNameSize));
245  stoppingSourceVSphi->setAxisTitle("track #phi",1);
246  stoppingSourceVSphi->setAxisTitle("stopping reason",2);
247 
248  for (size_t ibin=0; ibin<StopReasonNameSize; ibin++) {
252  }
253  }
254 
255 
256 
257  // book the TrackCandidate histograms
258  // ---------------------------------------------------------------------------------//
259 
260  if (doTCPlots){
261 
262  ibooker.setCurrentFolder(MEFolderName+"/TrackBuilding");
263 
264  histname = "TrackCandPt_"+tcProducer.label() + "_";
265  TrackCandPt = ibooker.book1D(histname+CatagoryName, histname+CatagoryName, TrackPtBin, TrackPtMin, TrackPtMax);
266  TrackCandPt->setAxisTitle("Track Candidate p_{T} (GeV/c)", 1);
267  TrackCandPt->setAxisTitle("Number of Track Candidates", 2);
268 
269  histname = "TrackCandEta_"+tcProducer.label() + "_";
270  TrackCandEta = ibooker.book1D(histname+CatagoryName, histname+CatagoryName, EtaBin, EtaMin, EtaMax);
271  TrackCandEta->setAxisTitle("Track Candidate #eta", 1);
272  TrackCandEta->setAxisTitle("Number of Track Candidates", 2);
273 
274  histname = "TrackCandPhi_"+tcProducer.label() + "_";
275  TrackCandPhi = ibooker.book1D(histname+CatagoryName, histname+CatagoryName, PhiBin, PhiMin, PhiMax);
276  TrackCandPhi->setAxisTitle("Track Candidate #phi", 1);
277  TrackCandPhi->setAxisTitle("Number of Track Candidates", 2);
278 
279  if (doTheta) {
280  histname = "TrackCandTheta_"+tcProducer.label() + "_";
281  TrackCandTheta = ibooker.book1D(histname+CatagoryName, histname+CatagoryName, ThetaBin, ThetaMin, ThetaMax);
282  TrackCandTheta->setAxisTitle("Track Candidate #theta", 1);
283  TrackCandTheta->setAxisTitle("Number of Track Candidates", 2);
284  }
285 
286  if (doAllTCPlots) {
287  histname = "TrackCandQ_"+tcProducer.label() + "_";
288  TrackCandQ = ibooker.book1D(histname+CatagoryName, histname+CatagoryName, TrackQBin, TrackQMin, TrackQMax);
289  TrackCandQ->setAxisTitle("Track Candidate Charge", 1);
290  TrackCandQ->setAxisTitle("Number of Track Candidates",2);
291 
292  histname = "TrackCandDxy_"+tcProducer.label() + "_";
293  TrackCandDxy = ibooker.book1D(histname+CatagoryName, histname+CatagoryName, TCDxyBin, TCDxyMin, TCDxyMax);
294  TrackCandDxy->setAxisTitle("Track Candidate d_{xy} (cm)", 1);
295  TrackCandDxy->setAxisTitle("Number of Track Candidates",2);
296 
297  histname = "TrackCandDz_"+tcProducer.label() + "_";
298  TrackCandDz = ibooker.book1D(histname+CatagoryName, histname+CatagoryName, TCDzBin, TCDzMin, TCDzMax);
299  TrackCandDz->setAxisTitle("Track Candidate d_{z} (cm)", 1);
300  TrackCandDz->setAxisTitle("Number of Track Candidates",2);
301 
302  histname = "NumberOfRecHitsPerTrackCand_"+tcProducer.label() + "_";
303  NumberOfRecHitsPerTrackCand = ibooker.book1D(histname+CatagoryName, histname+CatagoryName, TCHitBin, TCHitMin, TCHitMax);
304  NumberOfRecHitsPerTrackCand->setAxisTitle("Number of RecHits per Track Candidate", 1);
305  NumberOfRecHitsPerTrackCand->setAxisTitle("Number of Track Candidates",2);
306 
307  histname = "NumberOfRecHitsPerTrackCandVsPhiProfile_"+tcProducer.label() + "_";
308  NumberOfRecHitsPerTrackCandVsPhiProfile = ibooker.bookProfile(histname+CatagoryName, histname+CatagoryName, PhiBin, PhiMin, PhiMax, TCHitBin, TCHitMin, TCHitMax,"s");
309  NumberOfRecHitsPerTrackCandVsPhiProfile->setAxisTitle("Track Candidate #phi",1);
310  NumberOfRecHitsPerTrackCandVsPhiProfile->setAxisTitle("Number of RecHits of each Track Candidate",2);
311 
312  histname = "NumberOfRecHitsPerTrackCandVsEtaProfile_"+tcProducer.label() + "_";
313  NumberOfRecHitsPerTrackCandVsEtaProfile = ibooker.bookProfile(histname+CatagoryName, histname+CatagoryName, EtaBin, EtaMin, EtaMax, TCHitBin, TCHitMin, TCHitMax,"s");
314  NumberOfRecHitsPerTrackCandVsEtaProfile->setAxisTitle("Track Candidate #eta",1);
315  NumberOfRecHitsPerTrackCandVsEtaProfile->setAxisTitle("Number of RecHits of each Track Candidate",2);
316  }
317 
318  histname = "TrackCandPhiVsEta_"+tcProducer.label() + "_";
319  TrackCandPhiVsEta = ibooker.book2D(histname+CatagoryName, histname+CatagoryName, EtaBin, EtaMin, EtaMax, PhiBin, PhiMin, PhiMax);
320  TrackCandPhiVsEta->setAxisTitle("Track Candidate #eta", 1);
321  TrackCandPhiVsEta->setAxisTitle("Track Candidate #phi", 2);
322  }
323 
324 }
325 
326 // -- Analyse
327 // ---------------------------------------------------------------------------------//
329 (
330  const edm::Event& iEvent,
331  const edm::EventSetup& iSetup,
332  const TrajectorySeed& candidate,
333  const reco::BeamSpot& bs,
334  const edm::ESHandle<MagneticField>& theMF,
336 )
337 {
338  TSCBLBuilderNoMaterial tscblBuilder;
339 
340  //get parameters and errors from the candidate state
341  auto const & theG = ((TkTransientTrackingRecHitBuilder const *)(theTTRHBuilder.product()))->geometry();
342  auto const & candSS = candidate.startingState();
343  TrajectoryStateOnSurface state = trajectoryStateTransform::transientState( candSS, &(theG->idToDet(candSS.detId())->surface()), theMF.product());
344  TrajectoryStateClosestToBeamLine tsAtClosestApproachSeed = tscblBuilder(*state.freeState(),bs);//as in TrackProducerAlgorithm
345  if(!(tsAtClosestApproachSeed.isValid())) {
346  edm::LogVerbatim("TrackBuilding") << "TrajectoryStateClosestToBeamLine not valid";
347  return;
348  }
349  GlobalPoint v0 = tsAtClosestApproachSeed.trackStateAtPCA().position();
350  GlobalVector p = tsAtClosestApproachSeed.trackStateAtPCA().momentum();
351  GlobalPoint v(v0.x()-bs.x0(),v0.y()-bs.y0(),v0.z()-bs.z0());
352 
353  double pt = sqrt(state.globalMomentum().perp2());
354  double eta = state.globalPosition().eta();
355  double phi = state.globalPosition().phi();
356  double theta = state.globalPosition().theta();
357  //double pm = sqrt(state.globalMomentum().mag2());
358  //double pz = state.globalMomentum().z();
359  //double qoverp = tsAtClosestApproachSeed.trackStateAtPCA().charge()/p.mag();
360  //double theta = p.theta();
361  //double lambda = M_PI/2-p.theta();
362  double numberOfHits = candidate.recHits().second-candidate.recHits().first;
363  double dxy = (-v.x()*sin(p.phi())+v.y()*cos(p.phi()));
364  double dz = v.z() - (v.x()*p.x()+v.y()*p.y())/p.perp() * p.z()/p.perp();
365 
366  // fill the ME's
367  if (doAllSeedPlots || doQ)SeedQ->Fill( state.charge() );
368  if (doAllSeedPlots || doPT) SeedPt->Fill( pt );
369  if (doAllSeedPlots || doETA) SeedEta->Fill( eta );
370  if (doAllSeedPlots || doPHI) SeedPhi->Fill( phi );
371  if (doAllSeedPlots || doPHIVsETA) SeedPhiVsEta->Fill( eta, phi);
372  if (doAllSeedPlots || doTheta) SeedTheta->Fill( theta );
373  if (doAllSeedPlots || doDxy) SeedDxy->Fill( dxy );
374  if (doAllSeedPlots || doDz) SeedDz->Fill( dz );
375  if (doAllSeedPlots || doNRecHits) NumberOfRecHitsPerSeed->Fill( numberOfHits );
376  if (doAllSeedPlots || doProfETA) NumberOfRecHitsPerSeedVsEtaProfile->Fill( eta, numberOfHits );
377  if (doAllSeedPlots || doProfPHI) NumberOfRecHitsPerSeedVsPhiProfile->Fill( phi, numberOfHits );
378 
379 }
380 
381 // -- Analyse
382 // ---------------------------------------------------------------------------------//
384 (
385  const edm::Event& iEvent,
386  const edm::EventSetup& iSetup,
387  const TrackCandidate& candidate,
388  const reco::BeamSpot& bs,
389  const edm::ESHandle<MagneticField>& theMF,
391 )
392 {
393  TSCBLBuilderNoMaterial tscblBuilder;
394 
395  //get parameters and errors from the candidate state
396  auto const & theG = ((TkTransientTrackingRecHitBuilder const *)(theTTRHBuilder.product()))->geometry();
397  auto const & candSS = candidate.trajectoryStateOnDet();
398  TrajectoryStateOnSurface state = trajectoryStateTransform::transientState( candSS, &(theG->idToDet(candSS.detId())->surface()), theMF.product());
399  TrajectoryStateClosestToBeamLine tsAtClosestApproachTrackCand = tscblBuilder(*state.freeState(),bs);//as in TrackProducerAlgorithm
400  if(!(tsAtClosestApproachTrackCand.isValid())) {
401  edm::LogVerbatim("TrackBuilding") << "TrajectoryStateClosestToBeamLine not valid";
402  return;
403  }
404  GlobalPoint v0 = tsAtClosestApproachTrackCand.trackStateAtPCA().position();
405  GlobalVector p = tsAtClosestApproachTrackCand.trackStateAtPCA().momentum();
406  GlobalPoint v(v0.x()-bs.x0(),v0.y()-bs.y0(),v0.z()-bs.z0());
407 
408  double pt = sqrt(state.globalMomentum().perp2());
409  double eta = state.globalPosition().eta();
410  double phi = state.globalPosition().phi();
411  double theta = state.globalPosition().theta();
412  //double pm = sqrt(state.globalMomentum().mag2());
413  //double pz = state.globalMomentum().z();
414  //double qoverp = tsAtClosestApproachTrackCand.trackStateAtPCA().charge()/p.mag();
415  //double theta = p.theta();
416  //double lambda = M_PI/2-p.theta();
417  double numberOfHits = candidate.recHits().second-candidate.recHits().first;
418  double dxy = (-v.x()*sin(p.phi())+v.y()*cos(p.phi()));
419 
420  double dz = v.z() - (v.x()*p.x()+v.y()*p.y())/p.perp() * p.z()/p.perp();
421 
422  if (doAllTCPlots || doStopSource) {
423  // stopping source
424  int max = stoppingSource->getNbinsX();
425  double stop = candidate.stopReason() > max ? double(max-1) : static_cast<double>(candidate.stopReason());
426  stoppingSource ->Fill(stop);
427  stoppingSourceVSeta ->Fill(eta,stop);
428  stoppingSourceVSphi ->Fill(phi,stop);
429  }
430 
431  if (doTCPlots){
432  // fill the ME's
433  if (doAllTCPlots) TrackCandQ->Fill( state.charge() );
434  TrackCandPt->Fill( pt );
435  TrackCandEta->Fill( eta );
436  TrackCandPhi->Fill( phi );
437  TrackCandPhiVsEta->Fill( eta, phi );
438  if (doTheta) TrackCandTheta->Fill( theta );
439  if (doAllTCPlots) TrackCandDxy->Fill( dxy );
440  if (doAllTCPlots) TrackCandDz->Fill( dz );
441  if (doAllTCPlots) NumberOfRecHitsPerTrackCand->Fill( numberOfHits );
442  if (doAllTCPlots) NumberOfRecHitsPerTrackCandVsEtaProfile->Fill( eta, numberOfHits );
443  if (doAllTCPlots) NumberOfRecHitsPerTrackCandVsPhiProfile->Fill( phi, numberOfHits );
444  }
445 }
MonitorElement * NumberOfRecHitsPerTrackCandVsPhiProfile
MonitorElement * TrackCandPt
T getParameter(std::string const &) const
double z0() const
z coordinate
Definition: BeamSpot.h:68
T perp() const
Definition: PV3DBase.h:72
range recHits() const
MonitorElement * bookProfile(Args &&...args)
Definition: DQMStore.h:157
MonitorElement * TrackCandDz
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
Geom::Phi< T > phi() const
Definition: PV3DBase.h:69
Geom::Theta< T > theta() const
T y() const
Definition: PV3DBase.h:63
void setBinLabel(int bin, const std::string &label, int axis=1)
set bin label for x, y or z axis (axis=1, 2, 3 respectively)
MonitorElement * stoppingSourceVSphi
GlobalPoint globalPosition() const
#define NULL
Definition: scimark2.h:8
T perp2() const
Definition: PV3DBase.h:71
const double EtaMax[kNumberCalorimeter]
MonitorElement * TrackCandPhi
static const std::string StopReasonName[]
void initHisto(DQMStore::IBooker &ibooker, const edm::ParameterSet &)
MonitorElement * TrackCandPhiVsEta
void Fill(long long x)
Geom::Theta< T > theta() const
Definition: PV3DBase.h:75
MonitorElement * NumberOfRecHitsPerTrackCandVsEtaProfile
MonitorElement * stoppingSource
PTrajectoryStateOnDet const & trajectoryStateOnDet() const
MonitorElement * NumberOfRecHitsPerSeed
void analyze(const edm::Event &iEvent, const edm::EventSetup &iSetup, const TrajectorySeed &seed, const reco::BeamSpot &bs, const edm::ESHandle< MagneticField > &theMF, const edm::ESHandle< TransientTrackingRecHitBuilder > &theTTRHBuilder)
int iEvent
Definition: GenABIO.cc:230
MonitorElement * NumberOfRecHitsPerSeedVsPhiProfile
MonitorElement * TrackCandDxy
T sqrt(T t)
Definition: SSEVec.h:18
FreeTrajectoryState const * freeState(bool withErrors=true) const
T z() const
Definition: PV3DBase.h:64
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
MonitorElement * TrackCandEta
MonitorElement * book1D(Args &&...args)
Definition: DQMStore.h:115
const double EtaMin[kNumberCalorimeter]
TrackBuildingAnalyzer(const edm::ParameterSet &)
PTrajectoryStateOnDet const & startingState() const
MonitorElement * stoppingSourceVSeta
void setCurrentFolder(const std::string &fullpath)
Definition: DQMStore.cc:277
TrajectoryStateOnSurface transientState(const PTrajectoryStateOnDet &ts, const Surface *surface, const MagneticField *field)
MonitorElement * book2D(Args &&...args)
Definition: DQMStore.h:133
range recHits() const
MonitorElement * TrackCandTheta
MonitorElement * SeedPhiVsEta
T eta() const
Definition: PV3DBase.h:76
GlobalVector globalMomentum() const
int getNbinsX(void) const
get # of bins in X-axis
double y0() const
y coordinate
Definition: BeamSpot.h:66
MonitorElement * NumberOfRecHitsPerSeedVsEtaProfile
T x() const
Definition: PV3DBase.h:62
void setAxisTitle(const std::string &title, int axis=1)
set x-, y- or z-axis title (axis=1, 2, 3 respectively)
T const * product() const
Definition: ESHandle.h:86
MonitorElement * NumberOfRecHitsPerTrackCand
uint8_t stopReason() const
double x0() const
x coordinate
Definition: BeamSpot.h:64