CMS 3D CMS Logo

CosmicSplitterValidation.cc
Go to the documentation of this file.
1 // -*- C++ -*-
2 //
3 // Package: CosmicSplitterValidation
4 // Class: CosmicSplitterValidation
5 //
13 //
14 // Original Author: Nhan Tran
15 // Created: Mon Jul 16m 16:56:34 CDT 2007
16 // $Id: CosmicSplitterValidation.cc,v 1.11 2010/03/29 13:18:43 mussgill Exp $
17 //
18 //
19 
20 // system include files
21 #include <algorithm>
22 #include <TTree.h>
23 
25 
34 
47 
49 
51 
57 
58 //
59 // class decleration
60 //
61 
63 public:
66 
67 
68 private:
69  virtual void beginJob() override;
70  virtual void analyze(const edm::Event &iEvent, const edm::EventSetup &iSetup) override;
71  virtual void endJob() override ;
72 
73  bool is_gold_muon(const edm::Event& e);
74 
79 
83  int goldenCtr;
87 
89  // ----------member data ---------------------------
90  //std::vector<AlignTransform> m_align;
91  // tree
92  TTree* splitterTree_;
93  // tree vars
95  // split track variables
110  // split track errors
118 
119  // original track variables
121  double dxy_org_, dz_org_;
124  double norchi2_org_;
127  // original track errors
131 
132  // split sta variables
143  // split sta_ errors
151 
152  // split glb_ variables
164  // split glb_ errors
172  // original glb muon variables
174  double dxy_orm_, dz_orm_;
177  double norchi2_orm_;
178  // original glb muon errors
182 
183  //consumes tokens
189 
190 };
191 
192 //
193 // constants, enums and typedefs
194 //
195 
196 //
197 // static data member definitions
198 //
199 
200 //
201 // constructors and destructor
202 //
204  splitTracks_(iConfig.getParameter<edm::InputTag>("splitTracks")),
205  splitGlobalMuons_(iConfig.getParameter<edm::InputTag>("splitGlobalMuons")),
206  originalTracks_(iConfig.getParameter<edm::InputTag>("originalTracks")),
207  originalGlobalMuons_(iConfig.getParameter<edm::InputTag>("originalGlobalMuons")),
208  checkIfGolden_(iConfig.getParameter<bool>("checkIfGolden")),
209  splitMuons_(iConfig.getParameter<bool> ("ifSplitMuons")),
211  goldenCtr(0),
212  twoTracksCtr(0),
215  splitterTree_(0),
217  dcaX1_spl_(0), dcaY1_spl_(0), dcaZ1_spl_(0),
218  dcaX2_spl_(0), dcaY2_spl_(0), dcaZ2_spl_(0),
219  dxy1_spl_(0), dxy2_spl_(0), ddxy_spl_(0),
220  dz1_spl_(0), dz2_spl_(0), ddz_spl_(0),
222  eta1_spl_(0), eta2_spl_(0), deta_spl_(0),
223  phi1_spl_(0), phi2_spl_(0), dphi_spl_(0),
224  pt1_spl_(0), pt2_spl_(0), dpt_spl_(0),
225  p1_spl_(0), p2_spl_(0), dp_spl_(0),
231 
232  dxy1Err_spl_(0), dxy2Err_spl_(0),
233  dz1Err_spl_(0), dz2Err_spl_(0),
235  eta1Err_spl_(0), eta2Err_spl_(0),
236  phi1Err_spl_(0), phi2Err_spl_(0),
237  pt1Err_spl_(0), pt2Err_spl_(0),
239 
240  dcaX_org_(0), dcaY_org_(0), dcaZ_org_(0),
241  dxy_org_(0), dz_org_(0),
242  theta_org_(0), eta_org_(0), phi_org_(0),
243  pt_org_(0), p_org_(0), qoverpt_org_(0),
244  norchi2_org_(0),
247 
248  dxyErr_org_(0), dzErr_org_(0),
251 
252  dcaX1_sta_(0), dcaY1_sta_(0), dcaZ1_sta_(0),
253  dcaX2_sta_(0), dcaY2_sta_(0), dcaZ2_sta_(0),
254  dxy1_sta_(0), dxy2_sta_(0), ddxy_sta_(0),
255  dz1_sta_(0), dz2_sta_(0), ddz_sta_(0),
257  eta1_sta_(0), eta2_sta_(0), deta_sta_(0),
258  phi1_sta_(0), phi2_sta_(0), dphi_sta_(0),
259  pt1_sta_(0), pt2_sta_(0), dpt_sta_(0),
260  p1_sta_(0), p2_sta_(0), dp_sta_(0),
262 
263  dxy1Err_sta_(0), dxy2Err_sta_(0),
264  dz1Err_sta_(0), dz2Err_sta_(0),
266  eta1Err_sta_(0), eta2Err_sta_(0),
267  phi1Err_sta_(0), phi2Err_sta_(0),
268  pt1Err_sta_(0), pt2Err_sta_(0),
270 
271  dcaX1_glb_(0), dcaY1_glb_(0), dcaZ1_glb_(0),
272  dcaX2_glb_(0), dcaY2_glb_(0), dcaZ2_glb_(0),
273  dxy1_glb_(0), dxy2_glb_(0), ddxy_glb_(0),
274  dz1_glb_(0), dz2_glb_(0), ddz_glb_(0),
276  eta1_glb_(0), eta2_glb_(0), deta_glb_(0),
277  phi1_glb_(0), phi2_glb_(0), dphi_glb_(0),
278  pt1_glb_(0), pt2_glb_(0), dpt_glb_(0),
279  p1_glb_(0), p2_glb_(0), dp_glb_(0),
281  norchi1_glb_(0), norchi2_glb_(0),
282 
283  dxy1Err_glb_(0), dxy2Err_glb_(0),
284  dz1Err_glb_(0), dz2Err_glb_(0),
286  eta1Err_glb_(0), eta2Err_glb_(0),
287  phi1Err_glb_(0), phi2Err_glb_(0),
288  pt1Err_glb_(0), pt2Err_glb_(0),
290 
291  dcaX_orm_(0), dcaY_orm_(0), dcaZ_orm_(0),
292  dxy_orm_(0), dz_orm_(0),
293  theta_orm_(0), eta_orm_(0), phi_orm_(0),
294  pt_orm_(0), p_orm_(0), qoverpt_orm_(0),
295  norchi2_orm_(0),
296  dxyErr_orm_(0), dzErr_orm_(0),
299 {
301  splitTracksToken_ = iC.consumes<std::vector<reco::Track>>(splitTracks_);
302  originalTracksToken_ = iC.consumes<std::vector<reco::Track>>(originalTracks_);
303  if (splitMuons_){
306  }
307  if (checkIfGolden_)
309 }
310 
311 
313 {}
314 
315 
316 //
317 // member functions
318 //
319 
320 // ------------ method called to for each event ------------
322 {
323 
324  // check if golden muon
325  bool isGolden = true;
326  if (checkIfGolden_) isGolden = is_gold_muon( iEvent );
327 
328  // grab collections
332  edm::Handle<reco::MuonCollection> originalGlobalMuons;
333  iEvent.getByToken(splitTracksToken_, tracks);
334  iEvent.getByToken(originalTracksToken_, originalTracks);
335  if (splitMuons_){
336  iEvent.getByToken(splitGlobalMuonsToken_, globalMuons);
337  iEvent.getByToken(originalGlobalMuonsToken_, originalGlobalMuons);
338  }
339 
342 
343 
345  if (isGolden) goldenCtr++;
346  if (tracks->size() == 2) twoTracksCtr++;
347  if (tracks->size() == 2 && originalTracks->size() == 1 && isGolden){
349 
350  int gmCtr = 0;
351  bool topGlobalMuonFlag = false;
352  bool bottomGlobalMuonFlag = false;
353  int topGlobalMuon = -1;
354  int bottomGlobalMuon = -1;
355  double topGlobalMuonNorchi2 = 1e10;
356  double bottomGlobalMuonNorchi2 = 1e10;
357 
358  if (splitMuons_){
359  // check if split global muons are good
360  for (std::vector<reco::Muon>::const_iterator gmI = globalMuons->begin(); gmI != globalMuons->end(); gmI++){
361 
362  if ( gmI->isTrackerMuon() && gmI->isStandAloneMuon() && gmI->isGlobalMuon() ){
363 
364  reco::TrackRef trackerTrackRef1( tracks, 0 );
365  reco::TrackRef trackerTrackRef2( tracks, 1 );
366 
367  if (gmI->innerTrack() == trackerTrackRef1){
368  if (gmI->globalTrack()->normalizedChi2() < topGlobalMuonNorchi2){
369  topGlobalMuonFlag = true;
370  topGlobalMuonNorchi2 = gmI->globalTrack()->normalizedChi2();
371  topGlobalMuon = gmCtr;
372  }
373  }
374  if (gmI->innerTrack() == trackerTrackRef2){
375  if (gmI->globalTrack()->normalizedChi2() < bottomGlobalMuonNorchi2){
376  bottomGlobalMuonFlag = true;
377  bottomGlobalMuonNorchi2 = gmI->globalTrack()->normalizedChi2();
378  bottomGlobalMuon = gmCtr;
379  }
380  }
381  }
382  gmCtr++;
383  }
384  }
385 
386  if ( (!splitMuons_) || (splitMuons_ && bottomGlobalMuonFlag && topGlobalMuonFlag) ){
387 
389 
390  // split tracks calculations
391  reco::Track track1 = tracks->at(0);
392  reco::Track track2 = tracks->at(1);
393 
394  math::XYZPoint dca1 = track1.referencePoint();
395  math::XYZPoint dca2 = track2.referencePoint();
396 
397  // looping through the hits for track 1
398  int Nrechits1 =0;
399  int nhitinTIB1 =0;
400  int nhitinTOB1 =0;
401  int nhitinTID1 =0;
402  int nhitinTEC1 =0;
403  int nhitinBPIX1 =0;
404  int nhitinFPIX1 =0;
405  for (trackingRecHit_iterator iHit = track1.recHitsBegin(); iHit != track1.recHitsEnd(); ++iHit) {
406 
407  if((*iHit)->isValid()) {
408 
409  Nrechits1++;
410 
411  int type =(*iHit)->geographicalId().subdetId();
412 
413  if(type==int(StripSubdetector::TIB)){++nhitinTIB1;}
414  if(type==int(StripSubdetector::TOB)){++nhitinTOB1;}
415  if(type==int(StripSubdetector::TID)){++nhitinTID1;}
416  if(type==int(StripSubdetector::TEC)){++nhitinTEC1;}
417  if(type==int( kBPIX)){++nhitinBPIX1;}
418  if(type==int( kFPIX)){++nhitinFPIX1;}
419 
420  }
421  }
422  nHits1_spl_ = Nrechits1;
423  nHitsTIB1_spl_ = nhitinTIB1;
424  nHitsTOB1_spl_ = nhitinTOB1;
425  nHitsTID1_spl_ = nhitinTID1;
426  nHitsTEC1_spl_ = nhitinTEC1;
427  nHitsPXB1_spl_ = nhitinBPIX1;
428  nHitsPXF1_spl_ = nhitinFPIX1;
429 
430  // looping through the hits for track 2
431  int Nrechits2 =0;
432  int nhitinTIB2 =0;
433  int nhitinTOB2 =0;
434  int nhitinTID2 =0;
435  int nhitinTEC2 =0;
436  int nhitinBPIX2 =0;
437  int nhitinFPIX2 =0;
438  for (trackingRecHit_iterator iHit = track2.recHitsBegin(); iHit != track2.recHitsEnd(); ++iHit) {
439 
440  if((*iHit)->isValid()) {
441 
442  Nrechits2++;
443 
444  int type =(*iHit)->geographicalId().subdetId();
445 
446  if(type==int(StripSubdetector::TIB)){++nhitinTIB2;}
447  if(type==int(StripSubdetector::TOB)){++nhitinTOB2;}
448  if(type==int(StripSubdetector::TID)){++nhitinTID2;}
449  if(type==int(StripSubdetector::TEC)){++nhitinTEC2;}
450  if(type==int( kBPIX)){++nhitinBPIX2;}
451  if(type==int( kFPIX)){++nhitinFPIX2;}
452 
453  }
454  }
455  nHits2_spl_ = Nrechits2;
456  nHitsTIB2_spl_ = nhitinTIB2;
457  nHitsTOB2_spl_ = nhitinTOB2;
458  nHitsTID2_spl_ = nhitinTID2;
459  nHitsTEC2_spl_ = nhitinTEC2;
460  nHitsPXB2_spl_ = nhitinBPIX2;
461  nHitsPXF2_spl_ = nhitinFPIX2;
462 
463 
464  double ddxy_Val = track1.d0() - track2.d0();
465  double ddz_Val = track1.dz() - track2.dz();
466  double dtheta_Val = track1.theta() - track2.theta();
467  double deta_Val = track1.eta() - track2.eta();
468  double dphi_Val = track1.phi() - track2.phi();
469  double dpt_Val = track1.pt() - track2.pt();
470  double dp_Val = track1.p() - track2.p();
471  double dqoverpt_Val = track1.charge() / track1.pt() - track2.charge() / track2.pt();
472 
473  // original tracks calculations
474  reco::Track origTrack = originalTracks->at(0);
475  math::XYZPoint dca_org = origTrack.referencePoint();
476 
477  // looping through the hits for the original track
478  int Nrechitsorg =0;
479  int nhitinTIBorg =0;
480  int nhitinTOBorg =0;
481  int nhitinTIDorg =0;
482  int nhitinTECorg =0;
483  int nhitinBPIXorg =0;
484  int nhitinFPIXorg =0;
485  for (trackingRecHit_iterator iHit = origTrack.recHitsBegin(); iHit != origTrack.recHitsEnd(); ++iHit) {
486 
487  if((*iHit)->isValid()) {
488 
489  Nrechitsorg++;
490 
491  int type =(*iHit)->geographicalId().subdetId();
492 
493  if(type==int(StripSubdetector::TIB)){++nhitinTIBorg;}
494  if(type==int(StripSubdetector::TOB)){++nhitinTOBorg;}
495  if(type==int(StripSubdetector::TID)){++nhitinTIDorg;}
496  if(type==int(StripSubdetector::TEC)){++nhitinTECorg;}
497  if(type==int( kBPIX)){++nhitinBPIXorg;}
498  if(type==int( kFPIX)){++nhitinFPIXorg;}
499 
500  }
501  }
502  nHits_org_ = Nrechitsorg;
503  nHitsTIB_org_ = nhitinTIBorg;
504  nHitsTOB_org_ = nhitinTOBorg;
505  nHitsTID_org_ = nhitinTIDorg;
506  nHitsTEC_org_ = nhitinTECorg;
507  nHitsPXB_org_ = nhitinBPIXorg;
508  nHitsPXF_org_ = nhitinFPIXorg;
509 
510  // fill tree
511  runNumber_ = iEvent.id().run();
512  luminosityBlock_ = iEvent.id().luminosityBlock();
513  eventNumber_ = iEvent.id().event();
514  // split tracks
515  dcaX1_spl_ = dca1.x();
516  dcaY1_spl_ = dca1.y();
517  dcaZ1_spl_ = dca1.z();
518  dcaX2_spl_ = dca2.x();
519  dcaY2_spl_ = dca2.y();
520  dcaZ2_spl_ = dca2.z();
521  dxy1_spl_ = track1.d0();
522  dxy2_spl_ = track2.d0();
523  ddxy_spl_ = ddxy_Val;
524  dz1_spl_ = track1.dz();
525  dz2_spl_ = track2.dz();
526  ddz_spl_ = ddz_Val;
527  theta1_spl_ = track1.theta();
528  theta2_spl_ = track2.theta();
529  dtheta_spl_ = dtheta_Val;
530  eta1_spl_ = track1.eta();
531  eta2_spl_ = track2.eta();
532  deta_spl_ = deta_Val;
533  phi1_spl_ = track1.phi();
534  phi2_spl_ = track2.phi();
535  dphi_spl_ = dphi_Val;
536  pt1_spl_ = track1.pt();
537  pt2_spl_ = track2.pt();
538  dpt_spl_ = dpt_Val;
539  p1_spl_ = track1.p();
540  p2_spl_ = track2.p();
541  dp_spl_ = dp_Val;
542  qoverpt1_spl_ = track1.charge() / track1.pt();
543  qoverpt2_spl_ = track2.charge() / track2.pt();
544  dqoverpt_spl_ = dqoverpt_Val;
545  dxy1Err_spl_ = track1.d0Error();
546  dxy2Err_spl_ = track2.d0Error();
547  dz1Err_spl_ = track1.dzError();
548  dz2Err_spl_ = track2.dzError();
549  theta1Err_spl_ = track1.thetaError();
550  theta2Err_spl_ = track2.thetaError();
551  eta1Err_spl_ = track1.etaError();
552  eta2Err_spl_ = track2.etaError();
553  phi1Err_spl_ = track1.phiError();
554  phi2Err_spl_ = track2.phiError();
555  pt1Err_spl_ = track1.ptError();
556  pt2Err_spl_ = track2.ptError();
559 
560  // original tracks
561  dcaX_org_ = dca_org.x();
562  dcaY_org_ = dca_org.y();
563  dcaZ_org_ = dca_org.z();
564  dxy_org_ = origTrack.d0();
565  dz_org_ = origTrack.dz();
566  theta_org_ = origTrack.theta();
567  eta_org_ = origTrack.eta();
568  phi_org_ = origTrack.phi();
569  pt_org_ = origTrack.pt();
570  p_org_ = origTrack.p();
571  qoverpt_org_ = origTrack.charge() / origTrack.pt();
572  norchi2_org_ = origTrack.normalizedChi2();
573 
574  dxyErr_org_ = origTrack.d0Error();
575  dzErr_org_ = origTrack.dzError();
576  thetaErr_org_ = origTrack.thetaError();
577  etaErr_org_ = origTrack.etaError();
578  phiErr_org_ = origTrack.phiError();
579  ptErr_org_ = origTrack.ptError();
581 
582  // split muons calculations
583  if (splitMuons_){
584 
585  reco::Muon muonTop = globalMuons->at( topGlobalMuon );
586  reco::Muon muonBottom = globalMuons->at( bottomGlobalMuon );
587 
588  reco::TrackRef glb1 = muonTop.globalTrack();
589  reco::TrackRef glb2 = muonBottom.globalTrack();
590  reco::TrackRef sta1 = muonTop.outerTrack();
591  reco::TrackRef sta2 = muonBottom.outerTrack();
592 
593  // standalone muon variables
594  dcaX1_sta_ = sta1->referencePoint().x();
595  dcaY1_sta_ = sta1->referencePoint().y();
596  dcaZ1_sta_ = sta1->referencePoint().z();
597  dcaX2_sta_ = sta2->referencePoint().x();
598  dcaY2_sta_ = sta2->referencePoint().y();
599  dcaZ2_sta_ = sta2->referencePoint().z();
600  dxy1_sta_ = sta1->d0();
601  dxy2_sta_ = sta2->d0();
602  ddxy_sta_ = sta1->d0() - sta2->d0();
603  dz1_sta_ = sta1->dz();
604  dz2_sta_ = sta2->dz();
605  ddz_sta_ = sta1->dz() - sta2->dz();
606  theta1_sta_ = sta1->theta();
607  theta2_sta_ = sta2->theta();
608  dtheta_sta_ = sta1->theta() - sta2->theta();
609  eta1_sta_ = sta1->eta();
610  eta2_sta_ = sta2->eta();
612  phi1_sta_ = sta1->phi();
613  phi2_sta_ = sta2->phi();
614  dphi_sta_ = sta1->phi() - sta2->phi();
615  pt1_sta_ = sta1->pt();
616  pt2_sta_ = sta2->pt();
617  dpt_sta_ = sta1->pt() - sta2->pt();
618  p1_sta_ = sta1->p();
619  p2_sta_ = sta2->p();
620  dp_sta_ = sta1->p() - sta2->p();
621  qoverpt1_sta_ = sta1->charge() / sta1->pt();
622  qoverpt2_sta_ = sta2->charge() / sta2->pt();
624  dxy1Err_sta_ = sta1->dxyError();
625  dxy2Err_sta_ = sta2->dxyError();
626  dz1Err_sta_ = sta1->dzError();
627  dz2Err_sta_ = sta2->dzError();
628  theta1Err_sta_ = sta1->thetaError();
629  theta2Err_sta_ = sta2->thetaError();
630  eta1Err_sta_ = sta1->etaError();
631  eta2Err_sta_ = sta2->etaError();
632  phi1Err_sta_ = sta1->phiError();
633  phi2Err_sta_ = sta2->phiError();
634  pt1Err_sta_ = sta1->ptError();
635  pt2Err_sta_ = sta2->ptError();
637  qoverpt2Err_sta_ = qoverpt2_sta_ * pt2Err_sta_ / pt2_sta_;
638 
639  // global muon variables
640  dcaX1_glb_ = glb1->referencePoint().x();
641  dcaY1_glb_ = glb1->referencePoint().y();
642  dcaZ1_glb_ = glb1->referencePoint().z();
643  dcaX2_glb_ = glb2->referencePoint().x();
644  dcaY2_glb_ = glb2->referencePoint().y();
645  dcaZ2_glb_ = glb2->referencePoint().z();
646  dxy1_glb_ = glb1->d0();
647  dxy2_glb_ = glb2->d0();
648  ddxy_glb_ = glb1->d0() - glb2->d0();
649  dz1_glb_ = glb1->dz();
650  dz2_glb_ = glb2->dz();
651  ddz_glb_ = glb1->dz() - glb2->dz();
652  theta1_glb_ = glb1->theta();
653  theta2_glb_ = glb2->theta();
654  dtheta_glb_ = glb1->theta() - glb2->theta();
655  eta1_glb_ = glb1->eta();
656  eta2_glb_ = glb2->eta();
658  phi1_glb_ = glb1->phi();
659  phi2_glb_ = glb2->phi();
660  dphi_glb_ = glb1->phi() - glb2->phi();
661  pt1_glb_ = glb1->pt();
662  pt2_glb_ = glb2->pt();
663  dpt_glb_ = glb1->pt() - glb2->pt();
664  p1_glb_ = glb1->p();
665  p2_glb_ = glb2->p();
666  dp_glb_ = glb1->p() - glb2->p();
667  qoverpt1_glb_ = glb1->charge() / glb1->pt();
668  qoverpt2_glb_ = glb2->charge() / glb2->pt();
670  norchi1_glb_ = glb1->normalizedChi2();
671  norchi2_glb_ = glb2->normalizedChi2();
672 
673  dxy1Err_glb_ = glb1->d0Error();
674  dxy2Err_glb_ = glb2->d0Error();
675  dz1Err_glb_ = glb1->dzError();
676  dz2Err_glb_ = glb2->dzError();
677  theta1Err_glb_ = glb1->thetaError();
678  theta2Err_glb_ = glb2->thetaError();
679  eta1Err_glb_ = glb1->etaError();
680  eta2Err_glb_ = glb2->etaError();
681  phi1Err_glb_ = glb1->phiError();
682  phi2Err_glb_ = glb2->phiError();
683  pt1Err_glb_ = glb1->ptError();
684  pt2Err_glb_ = glb2->ptError();
686  qoverpt2Err_glb_ = qoverpt2_glb_ * pt2Err_glb_ / pt2_glb_;
687 
688  }
689 
690 
691  splitterTree_->Fill();
692  }
693  }
694 
695 
696 }
697 
698 
699 // ------------ method called once each job just before starting event loop ------------
701 {
702  edm::LogInfo("beginJob") << "Begin Job" << std::endl;
703 
704  splitterTree_ = tfile->make<TTree>("splitterTree","splitterTree");
705 
706  splitterTree_->Branch("runNumber", &runNumber_, "runNumber/I");
707  splitterTree_->Branch("eventNumber", &eventNumber_, "eventNumber/I");
708  splitterTree_->Branch("luminosityBlock", &luminosityBlock_, "luminosityBlock/I");
709 
710  // split track variables
711  splitterTree_->Branch("dcaX1_spl", &dcaX1_spl_, "dcaX1_spl/D");
712  splitterTree_->Branch("dcaY1_spl", &dcaY1_spl_, "dcaY1_spl/D");
713  splitterTree_->Branch("dcaZ1_spl", &dcaZ1_spl_, "dcaZ1_spl/D");
714  splitterTree_->Branch("dcaX2_spl", &dcaX2_spl_, "dcaX2_spl/D");
715  splitterTree_->Branch("dcaY2_spl", &dcaY2_spl_, "dcaY2_spl/D");
716  splitterTree_->Branch("dcaZ2_spl", &dcaZ2_spl_, "dcaZ2_spl/D");
717  splitterTree_->Branch("dxy1_spl", &dxy1_spl_, "dxy1_spl/D");
718  splitterTree_->Branch("dxy2_spl", &dxy2_spl_, "dxy2_spl/D");
719  splitterTree_->Branch("dz1_spl", &dz1_spl_, "dz1_spl/D");
720  splitterTree_->Branch("dz2_spl", &dz2_spl_, "dz2_spl/D");
721  splitterTree_->Branch("theta1_spl", &theta1_spl_, "theta1_spl/D");
722  splitterTree_->Branch("theta2_spl", &theta2_spl_, "theta2_spl/D");
723  splitterTree_->Branch("eta1_spl", &eta1_spl_, "eta1_spl/D");
724  splitterTree_->Branch("eta2_spl", &eta2_spl_, "eta2_spl/D");
725  splitterTree_->Branch("phi1_spl", &phi1_spl_, "phi1_spl/D");
726  splitterTree_->Branch("phi2_spl", &phi2_spl_, "phi2_spl/D");
727  splitterTree_->Branch("pt1_spl", &pt1_spl_, "pt1_spl/D");
728  splitterTree_->Branch("pt2_spl", &pt2_spl_, "pt2_spl/D");
729  splitterTree_->Branch("p1_spl", &p1_spl_, "p1_spl/D");
730  splitterTree_->Branch("p2_spl", &p2_spl_, "p2_spl/D");
731  splitterTree_->Branch("qoverpt1_spl", &qoverpt1_spl_, "qoverpt1_spl/D");
732  splitterTree_->Branch("qoverpt2_spl", &qoverpt2_spl_, "qoverpt2_spl/D");
733  splitterTree_->Branch("nHits1_spl", &nHits1_spl_, "nHits1_spl/I");
734  splitterTree_->Branch("nHitsPXB1_spl", &nHitsPXB1_spl_, "nHitsPXB1_spl/I");
735  splitterTree_->Branch("nHitsPXF1_spl", &nHitsPXF1_spl_, "nHitsPXF1_spl/I");
736  splitterTree_->Branch("nHitsTIB1_spl", &nHitsTIB1_spl_, "nHitsTIB1_spl/I");
737  splitterTree_->Branch("nHitsTOB1_spl", &nHitsTOB1_spl_, "nHitsTOB1_spl/I");
738  splitterTree_->Branch("nHitsTID1_spl", &nHitsTID1_spl_, "nHitsTID1_spl/I");
739  splitterTree_->Branch("nHitsTEC1_spl", &nHitsTEC1_spl_, "nHitsTEC1_spl/I");
740  splitterTree_->Branch("nHits2_spl", &nHits2_spl_, "nHits2_spl/I");
741  splitterTree_->Branch("nHitsPXB2_spl", &nHitsPXB2_spl_, "nHitsPXB2_spl/I");
742  splitterTree_->Branch("nHitsPXF2_spl", &nHitsPXF2_spl_, "nHitsPXF2_spl/I");
743  splitterTree_->Branch("nHitsTIB2_spl", &nHitsTIB2_spl_, "nHitsTIB2_spl/I");
744  splitterTree_->Branch("nHitsTOB2_spl", &nHitsTOB2_spl_, "nHitsTOB2_spl/I");
745  splitterTree_->Branch("nHitsTID2_spl", &nHitsTID2_spl_, "nHitsTID2_spl/I");
746  splitterTree_->Branch("nHitsTEC2_spl", &nHitsTEC2_spl_, "nHitsTEC2_spl/I");
747 
748 
749  splitterTree_->Branch("dxy1Err_spl", &dxy1Err_spl_, "dxy1Err_spl/D");
750  splitterTree_->Branch("dxy2Err_spl", &dxy2Err_spl_, "dxy2Err_spl/D");
751  splitterTree_->Branch("dz1Err_spl", &dz1Err_spl_, "dz1Err_spl/D");
752  splitterTree_->Branch("dz2Err_spl", &dz2Err_spl_, "dz2Err_spl/D");
753  splitterTree_->Branch("theta1Err_spl", &theta1Err_spl_, "theta1Err_spl/D");
754  splitterTree_->Branch("theta2Err_spl", &theta2Err_spl_, "theta2Err_spl/D");
755  splitterTree_->Branch("eta1Err_spl", &eta1Err_spl_, "eta1Err_spl/D");
756  splitterTree_->Branch("eta2Err_spl", &eta2Err_spl_, "eta2Err_spl/D");
757  splitterTree_->Branch("phi1Err_spl", &phi1Err_spl_, "phi1Err_spl/D");
758  splitterTree_->Branch("phi2Err_spl", &phi2Err_spl_, "phi2Err_spl/D");
759  splitterTree_->Branch("pt1Err_spl", &pt1Err_spl_, "pt1Err_spl/D");
760  splitterTree_->Branch("pt2Err_spl", &pt2Err_spl_, "pt2Err_spl/D");
761  splitterTree_->Branch("qoverpt1Err_spl", &qoverpt1Err_spl_, "qoverpt1Err_spl/D");
762  splitterTree_->Branch("qoverpt2Err_spl", &qoverpt2Err_spl_, "qoverpt2Err_spl/D");
763 
764  splitterTree_->Branch("dcaX_org", &dcaX_org_, "dcaX_org/D");
765  splitterTree_->Branch("dcaY_org", &dcaY_org_, "dcaY_org/D");
766  splitterTree_->Branch("dcaZ_org", &dcaZ_org_, "dcaZ_org/D");
767  splitterTree_->Branch("dxy_org", &dxy_org_, "dxy_org/D");
768  splitterTree_->Branch("dz_org", &dz_org_, "dz_org/D");
769  splitterTree_->Branch("theta_org", &theta_org_, "theta_org/D");
770  splitterTree_->Branch("eta_org", &eta_org_, "eta_org/D");
771  splitterTree_->Branch("phi_org", &phi_org_, "phi_org/D");
772  splitterTree_->Branch("pt_org", &pt_org_, "pt_org/D");
773  splitterTree_->Branch("p_org", &p_org_, "p_org/D");
774  splitterTree_->Branch("qoverpt_org", &qoverpt_org_, "qoverpt_org/D");
775  splitterTree_->Branch("norchi2_org", &norchi2_org_, "norchi2_org/D");
776 
777  splitterTree_->Branch("nHits_org", &nHits_org_, "nHits_org/I");
778  splitterTree_->Branch("nHitsPXB_org", &nHitsPXB_org_, "nHitsPXB_org/I");
779  splitterTree_->Branch("nHitsPXF_org", &nHitsPXF_org_, "nHitsPXF_org/I");
780  splitterTree_->Branch("nHitsTIB_org", &nHitsTIB_org_, "nHitsTIB_org/I");
781  splitterTree_->Branch("nHitsTOB_org", &nHitsTOB_org_, "nHitsTOB_org/I");
782  splitterTree_->Branch("nHitsTID_org", &nHitsTID_org_, "nHitsTID_org/I");
783  splitterTree_->Branch("nHitsTEC_org", &nHitsTEC_org_, "nHitsTEC_org/I");
784 
785  splitterTree_->Branch("dxyErr_org", &dxyErr_org_, "dxyErr_org/D");
786  splitterTree_->Branch("dzErr_org", &dzErr_org_, "dzErr_org/D");
787  splitterTree_->Branch("thetaErr_org", &thetaErr_org_, "thetaErr_org/D");
788  splitterTree_->Branch("etaErr_org", &etaErr_org_, "etaErr_org/D");
789  splitterTree_->Branch("phiErr_org", &phiErr_org_, "phiErr_org/D");
790  splitterTree_->Branch("ptErr_org", &ptErr_org_, "ptErr_org/D");
791  splitterTree_->Branch("qoverptErr_org", &qoverptErr_org_, "qoverptErr_org/D");
792 
793  //put the Deltas at the end of the tree, since they're the focus of the validation.
794  //this way, if you use splitterTree->Show(), they are immediately visible
795  splitterTree_->Branch("Delta_dxy", &ddxy_spl_, "Delta_dxy/D");
796  splitterTree_->Branch("Delta_dz", &ddz_spl_, "Delta_dz/D");
797  splitterTree_->Branch("Delta_theta", &dtheta_spl_, "Delta_theta/D");
798  splitterTree_->Branch("Delta_eta", &deta_spl_, "Delta_eta/D");
799  splitterTree_->Branch("Delta_phi", &dphi_spl_, "Delta_phi/D");
800  splitterTree_->Branch("Delta_pt", &dpt_spl_, "Delta_pt/D");
801  splitterTree_->Branch("Delta_p", &dp_spl_, "Delta_p/D");
802  splitterTree_->Branch("Delta_qoverpt", &dqoverpt_spl_, "Delta_qoverpt/D");
803 
804  if (splitMuons_){
805 
806  // standalone split
807  splitterTree_->Branch("dcaX1_sta", &dcaX1_sta_, "dcaX1_sta/D");
808  splitterTree_->Branch("dcaY1_sta", &dcaY1_sta_, "dcaY1_sta/D");
809  splitterTree_->Branch("dcaZ1_sta", &dcaZ1_sta_, "dcaZ1_sta/D");
810  splitterTree_->Branch("dcaX2_sta", &dcaX2_sta_, "dcaX2_sta/D");
811  splitterTree_->Branch("dcaY2_sta", &dcaY2_sta_, "dcaY2_sta/D");
812  splitterTree_->Branch("dcaZ2_sta", &dcaZ2_sta_, "dcaZ2_sta/D");
813  splitterTree_->Branch("dxy1_sta", &dxy1_sta_, "dxy1_sta/D");
814  splitterTree_->Branch("dxy2_sta", &dxy2_sta_, "dxy2_sta/D");
815  splitterTree_->Branch("ddxy_sta", &ddxy_sta_, "ddxy_sta/D");
816  splitterTree_->Branch("dz1_sta", &dz1_sta_, "dz1_sta/D");
817  splitterTree_->Branch("dz2_sta", &dz2_sta_, "dz2_sta/D");
818  splitterTree_->Branch("ddz_sta", &ddz_sta_, "ddz_sta/D");
819  splitterTree_->Branch("theta1_sta", &theta1_sta_, "theta1_sta/D");
820  splitterTree_->Branch("theta2_sta", &theta2_sta_, "theta2_sta/D");
821  splitterTree_->Branch("dtheta_sta", &dtheta_sta_, "dtheta_sta/D");
822  splitterTree_->Branch("eta1_sta", &eta1_sta_, "eta1_sta/D");
823  splitterTree_->Branch("eta2_sta", &eta2_sta_, "eta2_sta/D");
824  splitterTree_->Branch("deta_sta", &deta_sta_, "deta_sta/D");
825  splitterTree_->Branch("phi1_sta", &phi1_sta_, "phi1_sta/D");
826  splitterTree_->Branch("phi2_sta", &phi2_sta_, "phi2_sta/D");
827  splitterTree_->Branch("dphi_sta", &dphi_sta_, "dphi_sta/D");
828  splitterTree_->Branch("pt1_sta", &pt1_sta_, "pt1_sta/D");
829  splitterTree_->Branch("pt2_sta", &pt2_sta_, "pt2_sta/D");
830  splitterTree_->Branch("dpt_sta", &dpt_sta_, "dpt_sta/D");
831  splitterTree_->Branch("p1_sta", &p1_sta_, "p1_sta/D");
832  splitterTree_->Branch("p2_sta", &p2_sta_, "p2_sta/D");
833  splitterTree_->Branch("dp_sta", &dp_sta_, "dp_sta/D");
834 
835  splitterTree_->Branch("dxy1Err_sta", &dxy1Err_sta_, "dxy1Err_sta/D");
836  splitterTree_->Branch("dxy2Err_sta", &dxy2Err_sta_, "dxy2Err_sta/D");
837  splitterTree_->Branch("dz1Err_sta", &dz1Err_sta_, "dz1Err_sta/D");
838  splitterTree_->Branch("dz2Err_sta", &dz2Err_sta_, "dz2Err_sta/D");
839  splitterTree_->Branch("theta1Err_sta", &theta1Err_sta_, "theta1Err_sta/D");
840  splitterTree_->Branch("theta2Err_sta", &theta2Err_sta_, "theta2Err_sta/D");
841  splitterTree_->Branch("eta1Err_sta", &eta1Err_sta_, "eta1Err_sta/D");
842  splitterTree_->Branch("eta2Err_sta", &eta2Err_sta_, "eta2Err_sta/D");
843  splitterTree_->Branch("phi1Err_sta", &phi1Err_sta_, "phi1Err_sta/D");
844  splitterTree_->Branch("phi2Err_sta", &phi2Err_sta_, "phi2Err_sta/D");
845  splitterTree_->Branch("pt1Err_sta", &pt1Err_sta_, "pt1Err_sta/D");
846  splitterTree_->Branch("pt2Err_sta", &pt2Err_sta_, "pt2Err_sta/D");
847  splitterTree_->Branch("qoverpt1Err_sta", &qoverpt1Err_sta_, "qoverpt1Err_sta/D");
848  splitterTree_->Branch("qoverpt2Err_sta", &qoverpt2Err_sta_, "qoverpt2Err_sta/D");
849 
850  // global split
851  splitterTree_->Branch("dcaX1_glb", &dcaX1_glb_, "dcaX1_glb/D");
852  splitterTree_->Branch("dcaY1_glb", &dcaY1_glb_, "dcaY1_glb/D");
853  splitterTree_->Branch("dcaZ1_glb", &dcaZ1_glb_, "dcaZ1_glb/D");
854  splitterTree_->Branch("dcaX2_glb", &dcaX2_glb_, "dcaX2_glb/D");
855  splitterTree_->Branch("dcaY2_glb", &dcaY2_glb_, "dcaY2_glb/D");
856  splitterTree_->Branch("dcaZ2_glb", &dcaZ2_glb_, "dcaZ2_glb/D");
857  splitterTree_->Branch("dxy1_glb", &dxy1_glb_, "dxy1_glb/D");
858  splitterTree_->Branch("dxy2_glb", &dxy2_glb_, "dxy2_glb/D");
859  splitterTree_->Branch("ddxy_glb", &ddxy_glb_, "ddxy_glb/D");
860  splitterTree_->Branch("dz1_glb", &dz1_glb_, "dz1_glb/D");
861  splitterTree_->Branch("dz2_glb", &dz2_glb_, "dz2_glb/D");
862  splitterTree_->Branch("ddz_glb", &ddz_glb_, "ddz_glb/D");
863  splitterTree_->Branch("theta1_glb", &theta1_glb_, "theta1_glb/D");
864  splitterTree_->Branch("theta2_glb", &theta2_glb_, "theta2_glb/D");
865  splitterTree_->Branch("dtheta_glb", &dtheta_glb_, "dtheta_glb/D");
866  splitterTree_->Branch("eta1_glb", &eta1_glb_, "eta1_glb/D");
867  splitterTree_->Branch("eta2_glb", &eta2_glb_, "eta2_glb/D");
868  splitterTree_->Branch("deta_glb", &deta_glb_, "deta_glb/D");
869  splitterTree_->Branch("phi1_glb", &phi1_glb_, "phi1_glb/D");
870  splitterTree_->Branch("phi2_glb", &phi2_glb_, "phi2_glb/D");
871  splitterTree_->Branch("dphi_glb", &dphi_glb_, "dphi_glb/D");
872  splitterTree_->Branch("pt1_glb", &pt1_glb_, "pt1_glb/D");
873  splitterTree_->Branch("pt2_glb", &pt2_glb_, "pt2_glb/D");
874  splitterTree_->Branch("dpt_glb", &dpt_glb_, "dpt_glb/D");
875  splitterTree_->Branch("p1_glb", &p1_glb_, "p1_glb/D");
876  splitterTree_->Branch("p2_glb", &p2_glb_, "p2_glb/D");
877  splitterTree_->Branch("dp_glb", &dp_glb_, "dp_glb/D");
878  splitterTree_->Branch("qoverpt1_glb", &qoverpt1_glb_, "qoverpt1_glb/D");
879  splitterTree_->Branch("qoverpt2_glb", &qoverpt2_glb_, "qoverpt2_glb/D");
880  splitterTree_->Branch("dqoverpt_glb", &dqoverpt_glb_, "dqoverpt_glb/D");
881  splitterTree_->Branch("norchi1_glb", &norchi1_glb_, "norchi1_glb/D");
882  splitterTree_->Branch("norchi2_glb", &norchi2_glb_, "norchi2_glb/D");
883 
884  splitterTree_->Branch("dxy1Err_glb", &dxy1Err_glb_, "dxy1Err_glb/D");
885  splitterTree_->Branch("dxy2Err_glb", &dxy2Err_glb_, "dxy2Err_glb/D");
886  splitterTree_->Branch("dz1Err_glb", &dz1Err_glb_, "dz1Err_glb/D");
887  splitterTree_->Branch("dz2Err_glb", &dz2Err_glb_, "dz2Err_glb/D");
888  splitterTree_->Branch("theta1Err_glb", &theta1Err_glb_, "theta1Err_glb/D");
889  splitterTree_->Branch("theta2Err_glb", &theta2Err_glb_, "theta2Err_glb/D");
890  splitterTree_->Branch("eta1Err_glb", &eta1Err_glb_, "eta1Err_glb/D");
891  splitterTree_->Branch("eta2Err_glb", &eta2Err_glb_, "eta2Err_glb/D");
892  splitterTree_->Branch("phi1Err_glb", &phi1Err_glb_, "phi1Err_glb/D");
893  splitterTree_->Branch("phi2Err_glb", &phi2Err_glb_, "phi2Err_glb/D");
894  splitterTree_->Branch("pt1Err_glb", &pt1Err_glb_, "pt1Err_glb/D");
895  splitterTree_->Branch("pt2Err_glb", &pt2Err_glb_, "pt2Err_glb/D");
896  splitterTree_->Branch("qoverpt1Err_glb", &qoverpt1Err_glb_, "qoverpt1Err_glb/D");
897  splitterTree_->Branch("qoverpt2Err_glb", &qoverpt2Err_glb_, "qoverpt2Err_glb/D");
898 
899  }
900 
902  goldenCtr = 0;
903  twoTracksCtr = 0;
906 }
907 
908 
909 // ------------ method called once each job just after ending the event loop ------------
911 
912  //std::cout << "totalTracksToAnalyzer: " << totalTracksToAnalyzer_ << std::endl;
913  std::cout << "golden: " << goldenCtr << ", two tracks: " << twoTracksCtr << ", "
914  << "golden+twotracks: " << goldenPlusTwoTracksCtr << ", "
915  << "tracks+muons cuts: " << _passesTracksPlusMuonsCuts << std::endl;
916 }
917 
919 
921  e.getByToken(STAMuonsToken_, muHandle);
922  const reco::MuonCollection & muons = *(muHandle.product());
923  // make sure there are 2 muons
924  if ( 2 != muons.size() ) return false;
925 
926  double mudd0=0., mudphi=0., muddsz=0., mudeta=0.;
927  for ( unsigned int bindex = 0; bindex < muons.size(); ++bindex ) {
928  reco::Muon mymuon = muons[bindex];
929  // deprecated in 21x (now outerTrack)
930  //reco::TrackRef mutrackref = mymuon.standAloneMuon();
931  reco::TrackRef mutrackref = mymuon.outerTrack();
932  const reco::Track* mutrack = mutrackref.get();
933  if (0 == bindex){
934  mudd0+=mutrack->d0();
935  mudphi+=mutrack->phi();
936  muddsz+=mutrack->dsz();
937  mudeta+=mymuon.eta();
938  }
939  if (1 == bindex){
940  mudd0-=mutrack->d0();
941  mudphi-=mutrack->phi();
942  muddsz-=mutrack->dsz();
943  mudeta-=mymuon.eta();
944  }
945  }
946  if ((fabs(mudd0)<15.0)&&(fabs(mudphi)<0.045)&&(fabs(muddsz)<20.0)&&(fabs(mudeta)<0.060)) return true;
947  return false;
948 }
949 
950 
951 //define this as a plug-in
RunNumber_t run() const
Definition: EventID.h:39
double p() const
momentum vector magnitude
Definition: TrackBase.h:610
EDGetTokenT< ProductType > consumes(edm::InputTag const &tag)
type
Definition: HCALResponse.h:21
edm::Service< TFileService > tfile
EventNumber_t event() const
Definition: EventID.h:41
const Point & referencePoint() const
Reference point on the track.
Definition: TrackBase.h:676
double d0Error() const
error on d0
Definition: TrackBase.h:797
edm::EDGetTokenT< reco::MuonCollection > STAMuonsToken_
edm::EDGetTokenT< reco::MuonCollection > splitGlobalMuonsToken_
double d0() const
dxy parameter in perigee convention (d0 = -dxy)
Definition: TrackBase.h:592
double normalizedChi2() const
chi-squared divided by n.d.o.f. (or chi-squared * 1e6 if n.d.o.f. is zero)
Definition: TrackBase.h:556
CosmicSplitterValidation(const edm::ParameterSet &)
virtual double eta() const final
momentum pseudorapidity
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:460
virtual void analyze(const edm::Event &iEvent, const edm::EventSetup &iSetup) override
double theta() const
polar angle
Definition: TrackBase.h:574
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:17
double etaError() const
error on eta
Definition: TrackBase.h:779
double phi() const
azimuthal angle of momentum vector
Definition: TrackBase.h:640
LuminosityBlockNumber_t luminosityBlock() const
Definition: EventID.h:40
double dsz() const
dsz parameter (THIS IS NOT the SZ impact parameter to (0,0,0) if refPoint is far from (0...
Definition: TrackBase.h:598
std::vector< Muon > MuonCollection
collection of Muon objects
Definition: MuonFwd.h:9
edm::EDGetTokenT< std::vector< reco::Track > > splitTracksToken_
const int kFPIX
int iEvent
Definition: GenABIO.cc:230
double eta() const
pseudorapidity of momentum vector
Definition: TrackBase.h:646
virtual void beginJob() override
ConsumesCollector consumesCollector()
Use a ConsumesCollector to gather consumes information from helper functions.
double pt() const
track transverse momentum
Definition: TrackBase.h:616
Definition: tfile.py:1
double ptError() const
error on Pt (set to 1000 TeV if charge==0 for safety)
Definition: TrackBase.h:758
double phiError() const
error on phi
Definition: TrackBase.h:785
edm::EDGetTokenT< reco::MuonCollection > originalGlobalMuonsToken_
T const * get() const
Returns C++ pointer to the item.
Definition: Ref.h:244
trackingRecHit_iterator recHitsBegin() const
Iterator to first hit on the track.
Definition: Track.h:104
virtual TrackRef outerTrack() const
reference to Track reconstructed in the muon detector only
Definition: Muon.h:51
double dz() const
dz parameter (= dsz/cos(lambda)). This is the track z0 w.r.t (0,0,0) only if the refPoint is close to...
Definition: TrackBase.h:604
double dzError() const
error on dz
Definition: TrackBase.h:809
T const * product() const
Definition: Handle.h:81
XYZPointD XYZPoint
point in space with cartesian internal representation
Definition: Point3D.h:12
edm::EventID id() const
Definition: EventBase.h:60
HLT enums.
virtual void endJob() override
bool is_gold_muon(const edm::Event &e)
int charge() const
track electric charge
Definition: TrackBase.h:562
edm::EDGetTokenT< std::vector< reco::Track > > originalTracksToken_
TrackingRecHitCollection::base::const_iterator trackingRecHit_iterator
iterator over a vector of reference to TrackingRecHit in the same collection
double thetaError() const
error on theta
Definition: TrackBase.h:767
virtual TrackRef globalTrack() const
reference to Track reconstructed in both tracked and muon detector
Definition: Muon.h:54
const int kBPIX
trackingRecHit_iterator recHitsEnd() const
Iterator to last hit on the track.
Definition: Track.h:109