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