CMS 3D CMS Logo
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
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 "FWCore/Framework/interface/EDAnalyzer.h"
22 #include "FWCore/ParameterSet/interface/ParameterSet.h"
23 #include "FWCore/Framework/interface/MakerMacros.h"
24 
25 #include <algorithm>
26 
27 #include "FWCore/Framework/interface/ESHandle.h"
28 #include "FWCore/Framework/interface/EventSetup.h"
29 #include "FWCore/MessageLogger/interface/MessageLogger.h"
30 
31 #include "FWCore/ParameterSet/interface/ParameterSet.h"
32 #include "FWCore/Utilities/interface/InputTag.h"
33 
34 
35 #include "FWCore/Framework/interface/EDProducer.h"
36 #include "FWCore/Framework/interface/Event.h"
37 #include "FWCore/Framework/interface/EventSetup.h"
38 #include "FWCore/Framework/interface/ESHandle.h"
39 #include "CommonTools/UtilAlgos/interface/TFileService.h"
40 #include "FWCore/ServiceRegistry/interface/Service.h"
41 
42 
43 #include "DataFormats/Common/interface/View.h"
44 #include "DataFormats/TrackReco/interface/TrackFwd.h"
45 #include "DataFormats/TrackReco/interface/Track.h"
46 #include "DataFormats/TrackCandidate/interface/TrackCandidateCollection.h"
47 #include "DataFormats/TrackingRecHit/interface/InvalidTrackingRecHit.h"
48 #include "DataFormats/MuonReco/interface/Muon.h"
49 #include "DataFormats/MuonReco/interface/MuonFwd.h"
50 
51 #include "TrackingTools/TransientTrack/interface/BasicTransientTrack.h"
52 #include "TrackingTools/TrajectoryState/interface/FreeTrajectoryState.h"
53 #include "DataFormats/TrackReco/interface/Track.h"
54 #include "DataFormats/TrackReco/interface/TrackFwd.h"
55 #include "DataFormats/Common/interface/RefToBase.h"
56 #include "TrackingTools/TransientTrackingRecHit/interface/TransientTrackingRecHitBuilder.h"
57 
58 #include "TrackingTools/TrajectoryState/interface/FreeTrajectoryState.h"
59 #include "TrackingTools/TransientTrack/interface/TransientTrack.h"
60 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
61 #include "TrackingTools/TrajectoryState/interface/CopyUsingClone.h"
62 #include "RecoVertex/VertexTools/interface/PerigeeLinearizedTrackState.h"
63 
64 #include <TTree.h>
65 
66 //
67 // class decleration
68 //
69 
70 class CosmicSplitterValidation : public edm::EDAnalyzer {
71 public:
72  explicit CosmicSplitterValidation(const edm::ParameterSet&);
73  ~CosmicSplitterValidation();
74 
75 
76 private:
77  virtual void beginJob() override;
78  virtual void analyze(const edm::Event &iEvent, const edm::EventSetup &iSetup) override;
79  virtual void endJob() override ;
80 
81  bool is_gold_muon(const edm::Event& e);
82 
83  edm::InputTag splitTracks_;
84  edm::InputTag splitGlobalMuons_;
85  edm::InputTag originalTracks_;
86  edm::InputTag originalGlobalMuons_;
87 
88  bool checkIfGolden_;
89  bool splitMuons_;
90  int totalTracksToAnalyzer_;
91  int goldenCtr;
92  int twoTracksCtr;
93  int goldenPlusTwoTracksCtr;
94  int _passesTracksPlusMuonsCuts;
95 
96  edm::Service<TFileService> tfile;
97  // ----------member data ---------------------------
98  //std::vector<AlignTransform> m_align;
99  // tree
100  TTree* splitterTree_;
101  // tree vars
102  // split track variables
103  double dcaX1_spl_, dcaY1_spl_, dcaZ1_spl_;
104  double dcaX2_spl_, dcaY2_spl_, dcaZ2_spl_;
105  double dxy1_spl_, dxy2_spl_, dz1_spl_, dz2_spl_;
106  double theta1_spl_, theta2_spl_, phi1_spl_, phi2_spl_;
107  double ddxy_spl_, ddz_spl_, dtheta_spl_, dphi_spl_;
108  double pt1_spl_, pt2_spl_, dpt_spl_, p1_spl_, p2_spl_;
109  double eta1_spl_, eta2_spl_, deta_spl_;
110  double nHits1_spl_, nHitsPXB1_spl_, nHitsPXF1_spl_, nHitsTIB1_spl_;
111  double nHitsTOB1_spl_, nHitsTID1_spl_, nHitsTEC1_spl_;
112  double nHits2_spl_, nHitsPXB2_spl_, nHitsPXF2_spl_, nHitsTIB2_spl_;
113  double nHitsTOB2_spl_, nHitsTID2_spl_, nHitsTEC2_spl_;
114  // spl_it track errors
115  double pt1Err_spl_, pt2Err_spl_;
116  double theta1Err_spl_, theta2Err_spl_;
117  double phi1Err_spl_, phi2Err_spl_;
118  double d01Err_spl_, d02Err_spl_;
119  double dz1Err_spl_, dz2Err_spl_;
120  // original track variables
121  double dcaX_org_, dcaY_org_, dcaZ_org_;
122  double dxy_org_, dz_org_;
123  double theta_org_, phi_org_, eta_org_, pt_org_, p_org_;
124  double norchi2_org_;
125 
126  // split sta variables
127  double dcaX1_sta_, dcaY1_sta_, dcaZ1_sta_;
128  double dcaX2_sta_, dcaY2_sta_, dcaZ2_sta_;
129  double dxy1_sta_, dxy2_sta_, dz1_sta_, dz2_sta_;
130  double theta1_sta_, theta2_sta_, phi1_sta_, phi2_sta_;
131  double ddxy_sta_, ddz_sta_, dtheta_sta_, dphi_sta_;
132  double pt1_sta_, pt2_sta_, dpt_sta_, p1_sta_, p2_sta_;
133  double eta1_sta_, eta2_sta_, deta_sta_;
134  // split sta_ errors
135  double pt1Err_sta_, pt2Err_sta_;
136  double theta1Err_sta_, theta2Err_sta_;
137  double phi1Err_sta_, phi2Err_sta_;
138  double d01Err_sta_, d02Err_sta_;
139  double dz1Err_sta_, dz2Err_sta_;
140 
141  // split glb_ variables
142  double dcaX1_glb_, dcaY1_glb_, dcaZ1_glb_;
143  double dcaX2_glb_, dcaY2_glb_, dcaZ2_glb_;
144  double dxy1_glb_, dxy2_glb_, dz1_glb_, dz2_glb_;
145  double theta1_glb_, theta2_glb_, phi1_glb_, phi2_glb_;
146  double ddxy_glb_, ddz_glb_, dtheta_glb_, dphi_glb_;
147  double pt1_glb_, pt2_glb_, dpt_glb_, p1_glb_, p2_glb_;
148  double eta1_glb_, eta2_glb_, deta_glb_;
149  double norchi1_glb_, norchi2_glb_;
150  // split glb_ errors
151  double pt1Err_glb_, pt2Err_glb_;
152  double theta1Err_glb_, theta2Err_glb_;
153  double phi1Err_glb_, phi2Err_glb_;
154  double d01Err_glb_, d02Err_glb_;
155  double dz1Err_glb_, dz2Err_glb_;
156  // original glb muon variables
157  double dcaX_orm_, dcaY_orm_, dcaZ_orm_;
158  double dxy_orm_, dz_orm_;
159  double theta_orm_, phi_orm_, eta_orm_, pt_orm_, p_orm_;
160  double norchi2_orm_;
161 
162 };
163 
164 //
165 // constants, enums and typedefs
166 //
167 
168 //
169 // static data member definitions
170 //
171 
172 //
173 // constructors and destructor
174 //
175 CosmicSplitterValidation::CosmicSplitterValidation(const edm::ParameterSet& iConfig):
176  splitTracks_(iConfig.getParameter<edm::InputTag>("splitTracks")),
177  splitGlobalMuons_(iConfig.getParameter<edm::InputTag>("splitGlobalMuons")),
178  originalTracks_(iConfig.getParameter<edm::InputTag>("originalTracks")),
179  originalGlobalMuons_(iConfig.getParameter<edm::InputTag>("originalGlobalMuons")),
180  checkIfGolden_(iConfig.getParameter<bool>("checkIfGolden")),
181  splitMuons_(iConfig.getParameter<bool> ("ifSplitMuons")),
182  totalTracksToAnalyzer_(0),
183  goldenCtr(0),
184  twoTracksCtr(0),
185  goldenPlusTwoTracksCtr(0),
186  _passesTracksPlusMuonsCuts(0),
187  splitterTree_(0),
188  dcaX1_spl_(0), dcaY1_spl_(0), dcaZ1_spl_(0),
189  dcaX2_spl_(0), dcaY2_spl_(0), dcaZ2_spl_(0),
190  dxy1_spl_(0), dxy2_spl_(0), dz1_spl_(0), dz2_spl_(0),
191  theta1_spl_(0), theta2_spl_(0), phi1_spl_(0), phi2_spl_(0),
192  ddxy_spl_(0), ddz_spl_(0), dtheta_spl_(0), dphi_spl_(0),
193  pt1_spl_(0), pt2_spl_(0), dpt_spl_(0), p1_spl_(0), p2_spl_(0),
194  eta1_spl_(0), eta2_spl_(0), deta_spl_(0),
195  nHits1_spl_(0), nHitsPXB1_spl_(0), nHitsPXF1_spl_(0), nHitsTIB1_spl_(0),
196  nHitsTOB1_spl_(0), nHitsTID1_spl_(0), nHitsTEC1_spl_(0),
197  nHits2_spl_(0), nHitsPXB2_spl_(0), nHitsPXF2_spl_(0), nHitsTIB2_spl_(0),
198  nHitsTOB2_spl_(0), nHitsTID2_spl_(0), nHitsTEC2_spl_(0),
199  pt1Err_spl_(0), pt2Err_spl_(0),
200  theta1Err_spl_(0), theta2Err_spl_(0),
201  phi1Err_spl_(0), phi2Err_spl_(0),
202  d01Err_spl_(0), d02Err_spl_(0),
203  dz1Err_spl_(0), dz2Err_spl_(0),
204  dcaX_org_(0), dcaY_org_(0), dcaZ_org_(0),
205  dxy_org_(0), dz_org_(0),
206  theta_org_(0), phi_org_(0), eta_org_(0), pt_org_(0), p_org_(0),
207  norchi2_org_(0),
208  dcaX1_sta_(0), dcaY1_sta_(0), dcaZ1_sta_(0),
209  dcaX2_sta_(0), dcaY2_sta_(0), dcaZ2_sta_(0),
210  dxy1_sta_(0), dxy2_sta_(0), dz1_sta_(0), dz2_sta_(0),
211  theta1_sta_(0), theta2_sta_(0), phi1_sta_(0), phi2_sta_(0),
212  ddxy_sta_(0), ddz_sta_(0), dtheta_sta_(0), dphi_sta_(0),
213  pt1_sta_(0), pt2_sta_(0), dpt_sta_(0), p1_sta_(0), p2_sta_(0),
214  eta1_sta_(0), eta2_sta_(0), deta_sta_(0),
215  pt1Err_sta_(0), pt2Err_sta_(0),
216  theta1Err_sta_(0), theta2Err_sta_(0),
217  phi1Err_sta_(0), phi2Err_sta_(0),
218  d01Err_sta_(0), d02Err_sta_(0),
219  dz1Err_sta_(0), dz2Err_sta_(0),
220  dcaX1_glb_(0), dcaY1_glb_(0), dcaZ1_glb_(0),
221  dcaX2_glb_(0), dcaY2_glb_(0), dcaZ2_glb_(0),
222  dxy1_glb_(0), dxy2_glb_(0), dz1_glb_(0), dz2_glb_(0),
223  theta1_glb_(0), theta2_glb_(0), phi1_glb_(0), phi2_glb_(0),
224  ddxy_glb_(0), ddz_glb_(0), dtheta_glb_(0), dphi_glb_(0),
225  pt1_glb_(0), pt2_glb_(0), dpt_glb_(0), p1_glb_(0), p2_glb_(0),
226  eta1_glb_(0), eta2_glb_(0), deta_glb_(0),
227  norchi1_glb_(0), norchi2_glb_(0),
228  pt1Err_glb_(0), pt2Err_glb_(0),
229  theta1Err_glb_(0), theta2Err_glb_(0),
230  phi1Err_glb_(0), phi2Err_glb_(0),
231  d01Err_glb_(0), d02Err_glb_(0),
232  dz1Err_glb_(0), dz2Err_glb_(0),
233  dcaX_orm_(0), dcaY_orm_(0), dcaZ_orm_(0),
234  dxy_orm_(0), dz_orm_(0),
235  theta_orm_(0), phi_orm_(0), eta_orm_(0), pt_orm_(0), p_orm_(0),
236  norchi2_orm_(0)
237 {
238 
239 }
240 
241 
242 CosmicSplitterValidation::~CosmicSplitterValidation()
243 {}
244 
245 
246 //
247 // member functions
248 //
249 
250 // ------------ method called to for each event ------------
251 void CosmicSplitterValidation::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup)
252 {
253 
254  // check if golden muon
255  bool isGolden = true;
256  if (checkIfGolden_) isGolden = is_gold_muon( iEvent );
257 
258  // grab collections
259  edm::Handle<std::vector<reco::Track> > tracks;
260  edm::Handle<reco::MuonCollection> globalMuons;
261  edm::Handle<reco::MuonCollection> originalGlobalMuons;
262  edm::Handle<std::vector<reco::Track> > originalTracks;
263  iEvent.getByLabel(splitTracks_, tracks);
264  iEvent.getByLabel(originalTracks_, originalTracks);
265  if (splitMuons_){
266  iEvent.getByLabel(splitGlobalMuons_, globalMuons);
267  iEvent.getByLabel(originalGlobalMuons_, originalGlobalMuons);
268  }
269 
270  const int kBPIX = PixelSubdetector::PixelBarrel;
271  const int kFPIX = PixelSubdetector::PixelEndcap;
272 
273 
274  totalTracksToAnalyzer_ = totalTracksToAnalyzer_ + tracks->size();
275  if (isGolden) goldenCtr++;
276  if (tracks->size() == 2) twoTracksCtr++;
277  if (tracks->size() == 2 && originalTracks->size() == 1 && isGolden){
278  goldenPlusTwoTracksCtr++;
279 
280  int gmCtr = 0;
281  bool topGlobalMuonFlag = false;
282  bool bottomGlobalMuonFlag = false;
283  int topGlobalMuon = -1;
284  int bottomGlobalMuon = -1;
285  double topGlobalMuonNorchi2 = 1e10;
286  double bottomGlobalMuonNorchi2 = 1e10;
287 
288  if (splitMuons_){
289  // check if split global muons are good
290  for (std::vector<reco::Muon>::const_iterator gmI = globalMuons->begin(); gmI != globalMuons->end(); gmI++){
291 
292  if ( gmI->isTrackerMuon() && gmI->isStandAloneMuon() && gmI->isGlobalMuon() ){
293 
294  reco::TrackRef trackerTrackRef1( tracks, 0 );
295  reco::TrackRef trackerTrackRef2( tracks, 1 );
296 
297  if (gmI->innerTrack() == trackerTrackRef1){
298  if (gmI->globalTrack()->normalizedChi2() < topGlobalMuonNorchi2){
299  topGlobalMuonFlag = true;
300  topGlobalMuonNorchi2 = gmI->globalTrack()->normalizedChi2();
301  topGlobalMuon = gmCtr;
302  }
303  }
304  if (gmI->innerTrack() == trackerTrackRef2){
305  if (gmI->globalTrack()->normalizedChi2() < bottomGlobalMuonNorchi2){
306  bottomGlobalMuonFlag = true;
307  bottomGlobalMuonNorchi2 = gmI->globalTrack()->normalizedChi2();
308  bottomGlobalMuon = gmCtr;
309  }
310  }
311  }
312  gmCtr++;
313  }
314  }
315 
316  if ( (!splitMuons_) || (splitMuons_ && bottomGlobalMuonFlag && topGlobalMuonFlag) ){
317 
318  _passesTracksPlusMuonsCuts++;
319 
320  // split tracks calculations
321  reco::Track track1 = tracks->at(0);
322  reco::Track track2 = tracks->at(1);
323 
324  math::XYZPoint dca1 = track1.referencePoint();
325  math::XYZPoint dca2 = track2.referencePoint();
326 
327  // looping through the hits for track 1
328  double Nrechits1 =0;
329  double nhitinTIB1 =0;
330  double nhitinTOB1 =0;
331  double nhitinTID1 =0;
332  double nhitinTEC1 =0;
333  double nhitinBPIX1 =0;
334  double nhitinFPIX1 =0;
335  for (trackingRecHit_iterator iHit = track1.recHitsBegin(); iHit != track1.recHitsEnd(); ++iHit) {
336 
337  if((*iHit)->isValid()) {
338 
339  Nrechits1++;
340 
341  int type =(*iHit)->geographicalId().subdetId();
342 
343  if(type==int(StripSubdetector::TIB)){++nhitinTIB1;}
344  if(type==int(StripSubdetector::TOB)){++nhitinTOB1;}
345  if(type==int(StripSubdetector::TID)){++nhitinTID1;}
346  if(type==int(StripSubdetector::TEC)){++nhitinTEC1;}
347  if(type==int( kBPIX)){++nhitinBPIX1;}
348  if(type==int( kFPIX)){++nhitinFPIX1;}
349 
350  }
351  }
352  nHits1_spl_ = Nrechits1;
353  nHitsTIB1_spl_ = nhitinTIB1;
354  nHitsTOB1_spl_ = nhitinTOB1;
355  nHitsTID1_spl_ = nhitinTID1;
356  nHitsTEC1_spl_ = nhitinTEC1;
357  nHitsPXB1_spl_ = nhitinBPIX1;
358  nHitsPXF1_spl_ = nhitinFPIX1;
359 
360  // looping through the hits for track 2
361  double Nrechits2 =0;
362  double nhitinTIB2 =0;
363  double nhitinTOB2 =0;
364  double nhitinTID2 =0;
365  double nhitinTEC2 =0;
366  double nhitinBPIX2 =0;
367  double nhitinFPIX2 =0;
368  for (trackingRecHit_iterator iHit = track2.recHitsBegin(); iHit != track2.recHitsEnd(); ++iHit) {
369 
370  if((*iHit)->isValid()) {
371 
372  Nrechits2++;
373 
374  int type =(*iHit)->geographicalId().subdetId();
375 
376  if(type==int(StripSubdetector::TIB)){++nhitinTIB2;}
377  if(type==int(StripSubdetector::TOB)){++nhitinTOB2;}
378  if(type==int(StripSubdetector::TID)){++nhitinTID2;}
379  //\\if(type==int(StripSubdetector::TEC)){++nhitinTEC2;}
380  if(type==int( kBPIX)){++nhitinBPIX2;}
381  if(type==int( kFPIX)){++nhitinFPIX2;}
382 
383  }
384  }
385  nHits2_spl_ = Nrechits2;
386  nHitsTIB2_spl_ = nhitinTIB2;
387  nHitsTOB2_spl_ = nhitinTOB2;
388  nHitsTID2_spl_ = nhitinTID2;
389  nHitsTEC2_spl_ = nhitinTEC2;
390  nHitsPXB2_spl_ = nhitinBPIX2;
391  nHitsPXF2_spl_ = nhitinFPIX2;
392 
393 
394  double dtheta_Val = track1.theta() - track2.theta();
395  double dphi_Val = track1.phi() - track2.phi();
396  double ddxy_Val = track1.d0() - track2.d0();
397  double ddz_Val = track1.dz() - track2.dz();
398  double dpt_Val = track1.pt() - track2.pt();
399 
400  // original tracks calculations
401  reco::Track origTrack = originalTracks->at(0);
402  math::XYZPoint dca_org = origTrack.referencePoint();
403 
404  // fill tree
405  // split tracks
406  dcaX1_spl_ = dca1.x();
407  dcaY1_spl_ = dca1.y();
408  dcaZ1_spl_ = dca1.z();
409  dcaX2_spl_ = dca2.x();
410  dcaY2_spl_ = dca2.y();
411  dcaZ2_spl_ = dca2.z();
412  dxy1_spl_ = track1.d0();
413  dxy2_spl_ = track2.d0();
414  dz1_spl_ = track1.dz();
415  dz2_spl_ = track2.dz();
416  d01Err_spl_ = track1.d0Error();
417  d02Err_spl_ = track2.d0Error();
418  dz1Err_spl_ = track1.dzError();
419  dz2Err_spl_ = track2.dzError();
420  theta1_spl_ = track1.theta();
421  theta2_spl_ = track2.theta();
422  theta1Err_spl_ = track1.thetaError();
423  theta2Err_spl_ = track2.thetaError();
424  phi1_spl_ = track1.phi();
425  phi2_spl_ = track2.phi();
426  phi1Err_spl_ = track1.phiError();
427  phi2Err_spl_ = track2.phiError();
428  ddxy_spl_ = ddxy_Val;
429  ddz_spl_ = ddz_Val;
430  dtheta_spl_ = dtheta_Val;
431  dphi_spl_ = dphi_Val;
432  pt1_spl_ = track1.pt();
433  pt2_spl_ = track2.pt();
434  p1_spl_ = track1.p();
435  p2_spl_ = track2.p();
436  pt1Err_spl_ = track1.ptError();
437  pt2Err_spl_ = track2.ptError();
438  dpt_spl_ = dpt_Val;
439  eta1_spl_ = track1.eta();
440  eta2_spl_ = track2.eta();
441  deta_spl_ = eta1_spl_ - eta2_spl_;
442 
443  // original tracks
444  dcaX_org_ = dca_org.x();
445  dcaY_org_ = dca_org.y();
446  dcaZ_org_ = dca_org.z();
447  dxy_org_ = origTrack.d0();
448  dz_org_ = origTrack.dz();
449  theta_org_ = origTrack.theta();
450  phi_org_ = origTrack.phi();
451  eta_org_ = origTrack.eta();
452  pt_org_ = origTrack.pt();
453  p_org_ = origTrack.p();
454  norchi2_org_ = origTrack.normalizedChi2();
455 
456  // split muons calculations
457  if (splitMuons_){
458 
459  reco::Muon muonTop = globalMuons->at( topGlobalMuon );
460  reco::Muon muonBottom = globalMuons->at( bottomGlobalMuon );
461 
462  reco::TrackRef glb1 = muonTop.globalTrack();
463  reco::TrackRef glb2 = muonBottom.globalTrack();
464  reco::TrackRef sta1 = muonTop.outerTrack();
465  reco::TrackRef sta2 = muonBottom.outerTrack();
466 
467  // standalone muon variables
468  dcaX1_sta_ = sta1->referencePoint().x();
469  dcaY1_sta_ = sta1->referencePoint().y();
470  dcaZ1_sta_ = sta1->referencePoint().z();
471  dcaX2_sta_ = sta2->referencePoint().x();
472  dcaY2_sta_ = sta2->referencePoint().y();
473  dcaZ2_sta_ = sta2->referencePoint().z();
474  dxy1_sta_ = sta1->d0();
475  dxy2_sta_ = sta2->d0();
476  dz1_sta_ = sta1->dz();
477  dz2_sta_ = sta2->dz();
478  d01Err_sta_ = sta1->d0Error();
479  d02Err_sta_ = sta2->d0Error();
480  dz1Err_sta_ = sta1->dzError();
481  dz2Err_sta_ = sta2->dzError();
482  theta1_sta_ = sta1->theta();
483  theta2_sta_ = sta2->theta();
484  theta1Err_sta_ = sta1->thetaError();
485  theta2Err_sta_ = sta2->thetaError();
486  phi1_sta_ = sta1->phi();
487  phi2_sta_ = sta2->phi();
488  phi1Err_sta_ = sta1->phiError();
489  phi2Err_sta_ = sta2->phiError();
490  ddxy_sta_ = sta1->d0() - sta2->d0();
491  ddz_sta_ = sta1->dz() - sta2->dz();
492  dtheta_sta_ = sta1->theta() - sta2->theta();
493  dphi_sta_ = sta1->phi() - sta2->phi();
494  pt1_sta_ = sta1->pt();
495  pt2_sta_ = sta2->pt();
496  p1_sta_ = sta1->p();
497  p2_sta_ = sta2->p();
498  pt1Err_sta_ = sta1->ptError();
499  pt2Err_sta_ = sta2->ptError();
500  dpt_sta_ = sta1->pt() - sta2->pt();
501  eta1_sta_ = sta1->eta();
502  eta2_sta_ = sta2->eta();
503  deta_sta_ = eta1_sta_ - eta2_sta_;
504 
505  // global muon variables
506  dcaX1_glb_ = glb1->referencePoint().x();
507  dcaY1_glb_ = glb1->referencePoint().y();
508  dcaZ1_glb_ = glb1->referencePoint().z();
509  dcaX2_glb_ = glb2->referencePoint().x();
510  dcaY2_glb_ = glb2->referencePoint().y();
511  dcaZ2_glb_ = glb2->referencePoint().z();
512  dxy1_glb_ = glb1->d0();
513  dxy2_glb_ = glb2->d0();
514  dz1_glb_ = glb1->dz();
515  dz2_glb_ = glb2->dz();
516  d01Err_glb_ = glb1->d0Error();
517  d02Err_glb_ = glb2->d0Error();
518  dz1Err_glb_ = glb1->dzError();
519  dz2Err_glb_ = glb2->dzError();
520  theta1_glb_ = glb1->theta();
521  theta2_glb_ = glb2->theta();
522  theta1Err_glb_ = glb1->thetaError();
523  theta2Err_glb_ = glb2->thetaError();
524  phi1_glb_ = glb1->phi();
525  phi2_glb_ = glb2->phi();
526  phi1Err_glb_ = glb1->phiError();
527  phi2Err_glb_ = glb2->phiError();
528  ddxy_glb_ = glb1->d0() - glb2->d0();
529  ddz_glb_ = glb1->dz() - glb2->dz();
530  dtheta_glb_ = glb1->theta() - glb2->theta();
531  dphi_glb_ = glb1->phi() - glb2->phi();
532  pt1_glb_ = glb1->pt();
533  pt2_glb_ = glb2->pt();
534  p1_glb_ = glb1->p();
535  p2_glb_ = glb2->p();
536  pt1Err_glb_ = glb1->ptError();
537  pt2Err_glb_ = glb2->ptError();
538  dpt_glb_ = glb1->pt() - glb2->pt();
539  eta1_glb_ = glb1->eta();
540  eta2_glb_ = glb2->eta();
541  deta_glb_ = eta1_glb_ - eta2_glb_;
542  norchi1_glb_ = glb1->normalizedChi2();
543  norchi2_glb_ = glb2->normalizedChi2();
544 
545  }
546 
547 
548  splitterTree_->Fill();
549  }
550  }
551 
552 
553 }
554 
555 
556 // ------------ method called once each job just before starting event loop ------------
557 void CosmicSplitterValidation::beginJob()
558 {
559  edm::LogInfo("beginJob") << "Begin Job" << std::endl;
560 
561  splitterTree_ = tfile->make<TTree>("splitterTree","splitterTree");
562 
563  // split track variables
564  splitterTree_->Branch("dcaX1_spl", &dcaX1_spl_, "dcaX1_spl/D");
565  splitterTree_->Branch("dcaY1_spl", &dcaY1_spl_, "dcaY1_spl/D");
566  splitterTree_->Branch("dcaZ1_spl", &dcaZ1_spl_, "dcaZ1_spl/D");
567  splitterTree_->Branch("dcaX2_spl", &dcaX2_spl_, "dcaX2_spl/D");
568  splitterTree_->Branch("dcaY2_spl", &dcaY2_spl_, "dcaY2_spl/D");
569  splitterTree_->Branch("dcaZ2_spl", &dcaZ2_spl_, "dcaZ2_spl/D");
570  splitterTree_->Branch("dxy1_spl", &dxy1_spl_, "dxy1_spl/D");
571  splitterTree_->Branch("dz1_spl", &dz1_spl_, "dz1_spl/D");
572  splitterTree_->Branch("dxy2_spl", &dxy2_spl_, "dxy2_spl/D");
573  splitterTree_->Branch("dz2_spl", &dz2_spl_, "dz2_spl/D");
574  splitterTree_->Branch("theta1_spl", &theta1_spl_, "theta1_spl/D");
575  splitterTree_->Branch("theta2_spl", &theta2_spl_, "theta2_spl/D");
576  splitterTree_->Branch("phi1_spl", &phi1_spl_, "phi1_spl/D");
577  splitterTree_->Branch("phi2_spl", &phi2_spl_, "phi2_spl/D");
578  splitterTree_->Branch("ddxy_spl", &ddxy_spl_, "ddxy_spl/D");
579  splitterTree_->Branch("ddz_spl", &ddz_spl_, "ddz_spl/D");
580  splitterTree_->Branch("dphi_spl", &dphi_spl_, "dphi_spl/D");
581  splitterTree_->Branch("dtheta_spl", &dtheta_spl_, "dtheta_spl/D");
582  splitterTree_->Branch("pt1_spl", &pt1_spl_, "pt1_spl/D");
583  splitterTree_->Branch("pt2_spl", &pt2_spl_, "pt2_spl/D");
584  splitterTree_->Branch("dpt_spl", &dpt_spl_, "dpt_spl/D");
585  splitterTree_->Branch("p1_spl", &p1_spl_, "p1_spl/D");
586  splitterTree_->Branch("p2_spl", &p2_spl_, "p2_spl/D");
587  splitterTree_->Branch("eta1_spl", &eta1_spl_, "eta1_spl/D");
588  splitterTree_->Branch("eta2_spl", &eta2_spl_, "eta2_spl/D");
589  splitterTree_->Branch("deta_spl", &deta_spl_, "deta_spl/D");
590  splitterTree_->Branch("nHits1_spl", &nHits1_spl_, "nHits1_spl/D");
591  splitterTree_->Branch("nHitsPXB1_spl", &nHitsPXB1_spl_, "nHitsPXB1_spl/D");
592  splitterTree_->Branch("nHitsPXF1_spl", &nHitsPXF1_spl_, "nHitsPXF1_spl/D");
593  splitterTree_->Branch("nHitsTIB1_spl", &nHitsTIB1_spl_, "nHitsTIB1_spl/D");
594  splitterTree_->Branch("nHitsTOB1_spl", &nHitsTOB1_spl_, "nHitsTOB1_spl/D");
595  splitterTree_->Branch("nHitsTID1_spl", &nHitsTID1_spl_, "nHitsTID1_spl/D");
596  splitterTree_->Branch("nHitsTEC1_spl", &nHitsTEC1_spl_, "nHitsTEC1_spl/D");
597  splitterTree_->Branch("nHits2_spl", &nHits2_spl_, "nHits2_spl/D");
598  splitterTree_->Branch("nHitsPXB2_spl", &nHitsPXB2_spl_, "nHitsPXB2_spl/D");
599  splitterTree_->Branch("nHitsPXF2_spl", &nHitsPXF2_spl_, "nHitsPXF2_spl/D");
600  splitterTree_->Branch("nHitsTIB2_spl", &nHitsTIB2_spl_, "nHitsTIB2_spl/D");
601  splitterTree_->Branch("nHitsTOB2_spl", &nHitsTOB2_spl_, "nHitsTOB2_spl/D");
602  splitterTree_->Branch("nHitsTID2_spl", &nHitsTID2_spl_, "nHitsTID2_spl/D");
603  splitterTree_->Branch("nHitsTEC2_spl", &nHitsTEC2_spl_, "nHitsTEC2_spl/D");
604 
605 
606  splitterTree_->Branch("d01Err_spl", &d01Err_spl_, "d01Err_spl/D");
607  splitterTree_->Branch("d02Err_spl", &d02Err_spl_, "d02Err_spl/D");
608  splitterTree_->Branch("dz1Err_spl", &dz1Err_spl_, "dz1Err_spl/D");
609  splitterTree_->Branch("dz2Err_spl", &dz2Err_spl_, "dz2Err_spl/D");
610  splitterTree_->Branch("phi1Err_spl", &phi1Err_spl_, "phi1Err_spl/D");
611  splitterTree_->Branch("phi2Err_spl", &phi2Err_spl_, "phi2Err_spl/D");
612  splitterTree_->Branch("theta1Err_spl", &theta1Err_spl_, "theta1Err_spl/D");
613  splitterTree_->Branch("theta2Err_spl", &theta2Err_spl_, "theta2Err_spl/D");
614  splitterTree_->Branch("pt1Err_spl", &pt1Err_spl_, "pt1Err_spl/D");
615  splitterTree_->Branch("pt2Err_spl", &pt2Err_spl_, "pt2Err_spl/D");
616 
617  splitterTree_->Branch("dcaX_org", &dcaX_org_, "dcaX_org/D");
618  splitterTree_->Branch("dcaY_org", &dcaY_org_, "dcaY_org/D");
619  splitterTree_->Branch("dcaZ_org", &dcaZ_org_, "dcaZ_org/D");
620  splitterTree_->Branch("dxy_org", &dxy_org_, "dxy_org/D");
621  splitterTree_->Branch("dz_org", &dz_org_, "dz_org/D");
622  splitterTree_->Branch("theta_org", &theta_org_, "theta_org/D");
623  splitterTree_->Branch("phi_org", &phi_org_, "phi_org/D");
624  splitterTree_->Branch("eta_org", &eta_org_, "eta_org/D");
625  splitterTree_->Branch("pt_org", &pt_org_, "pt_org/D");
626  splitterTree_->Branch("p_org", &p_org_, "p_org/D");
627  splitterTree_->Branch("norchi2_org", &norchi2_org_, "norchi2_org/D");
628 
629  if (splitMuons_){
630 
631  // standalone split
632  splitterTree_->Branch("dcaX1_sta", &dcaX1_sta_, "dcaX1_sta/D");
633  splitterTree_->Branch("dcaY1_sta", &dcaY1_sta_, "dcaY1_sta/D");
634  splitterTree_->Branch("dcaZ1_sta", &dcaZ1_sta_, "dcaZ1_sta/D");
635  splitterTree_->Branch("dcaX2_sta", &dcaX2_sta_, "dcaX2_sta/D");
636  splitterTree_->Branch("dcaY2_sta", &dcaY2_sta_, "dcaY2_sta/D");
637  splitterTree_->Branch("dcaZ2_sta", &dcaZ2_sta_, "dcaZ2_sta/D");
638  splitterTree_->Branch("dxy1_sta", &dxy1_sta_, "dxy1_sta/D");
639  splitterTree_->Branch("dz1_sta", &dz1_sta_, "dz1_sta/D");
640  splitterTree_->Branch("dxy2_sta", &dxy2_sta_, "dxy2_sta/D");
641  splitterTree_->Branch("dz2_sta", &dz2_sta_, "dz2_sta/D");
642  splitterTree_->Branch("theta1_sta", &theta1_sta_, "theta1_sta/D");
643  splitterTree_->Branch("theta2_sta", &theta2_sta_, "theta2_sta/D");
644  splitterTree_->Branch("phi1_sta", &phi1_sta_, "phi1_sta/D");
645  splitterTree_->Branch("phi2_sta", &phi2_sta_, "phi2_sta/D");
646  splitterTree_->Branch("ddxy_sta", &ddxy_sta_, "ddxy_sta/D");
647  splitterTree_->Branch("ddz_sta", &ddz_sta_, "ddz_sta/D");
648  splitterTree_->Branch("dphi_sta", &dphi_sta_, "dphi_sta/D");
649  splitterTree_->Branch("dtheta_sta", &dtheta_sta_, "dtheta_sta/D");
650  splitterTree_->Branch("pt1_sta", &pt1_sta_, "pt1_sta/D");
651  splitterTree_->Branch("pt2_sta", &pt2_sta_, "pt2_sta/D");
652  splitterTree_->Branch("dpt_sta", &dpt_sta_, "dpt_sta/D");
653  splitterTree_->Branch("p1_sta", &p1_sta_, "p1_sta/D");
654  splitterTree_->Branch("p2_sta", &p2_sta_, "p2_sta/D");
655  splitterTree_->Branch("eta1_sta", &eta1_sta_, "eta1_sta/D");
656  splitterTree_->Branch("eta2_sta", &eta2_sta_, "eta2_sta/D");
657  splitterTree_->Branch("deta_sta", &deta_sta_, "deta_sta/D");
658 
659  splitterTree_->Branch("d01Err_sta", &d01Err_sta_, "d01Err_sta/D");
660  splitterTree_->Branch("d02Err_sta", &d02Err_sta_, "d02Err_sta/D");
661  splitterTree_->Branch("dz1Err_sta", &dz1Err_sta_, "dz1Err_sta/D");
662  splitterTree_->Branch("dz2Err_sta", &dz2Err_sta_, "dz2Err_sta/D");
663  splitterTree_->Branch("phi1Err_sta", &phi1Err_sta_, "phi1Err_sta/D");
664  splitterTree_->Branch("phi2Err_sta", &phi2Err_sta_, "phi2Err_sta/D");
665  splitterTree_->Branch("theta1Err_sta", &theta1Err_sta_, "theta1Err_sta/D");
666  splitterTree_->Branch("theta2Err_sta", &theta2Err_sta_, "theta2Err_sta/D");
667  splitterTree_->Branch("pt1Err_sta", &pt1Err_sta_, "pt1Err_sta/D");
668  splitterTree_->Branch("pt2Err_sta", &pt2Err_sta_, "pt2Err_sta/D");
669 
670  // global split
671  splitterTree_->Branch("dcaX1_glb", &dcaX1_glb_, "dcaX1_glb/D");
672  splitterTree_->Branch("dcaY1_glb", &dcaY1_glb_, "dcaY1_glb/D");
673  splitterTree_->Branch("dcaZ1_glb", &dcaZ1_glb_, "dcaZ1_glb/D");
674  splitterTree_->Branch("dcaX2_glb", &dcaX2_glb_, "dcaX2_glb/D");
675  splitterTree_->Branch("dcaY2_glb", &dcaY2_glb_, "dcaY2_glb/D");
676  splitterTree_->Branch("dcaZ2_glb", &dcaZ2_glb_, "dcaZ2_glb/D");
677  splitterTree_->Branch("dxy1_glb", &dxy1_glb_, "dxy1_glb/D");
678  splitterTree_->Branch("dz1_glb", &dz1_glb_, "dz1_glb/D");
679  splitterTree_->Branch("dxy2_glb", &dxy2_glb_, "dxy2_glb/D");
680  splitterTree_->Branch("dz2_glb", &dz2_glb_, "dz2_glb/D");
681  splitterTree_->Branch("theta1_glb", &theta1_glb_, "theta1_glb/D");
682  splitterTree_->Branch("theta2_glb", &theta2_glb_, "theta2_glb/D");
683  splitterTree_->Branch("phi1_glb", &phi1_glb_, "phi1_glb/D");
684  splitterTree_->Branch("phi2_glb", &phi2_glb_, "phi2_glb/D");
685  splitterTree_->Branch("ddxy_glb", &ddxy_glb_, "ddxy_glb/D");
686  splitterTree_->Branch("ddz_glb", &ddz_glb_, "ddz_glb/D");
687  splitterTree_->Branch("dphi_glb", &dphi_glb_, "dphi_glb/D");
688  splitterTree_->Branch("dtheta_glb", &dtheta_glb_, "dtheta_glb/D");
689  splitterTree_->Branch("pt1_glb", &pt1_glb_, "pt1_glb/D");
690  splitterTree_->Branch("pt2_glb", &pt2_glb_, "pt2_glb/D");
691  splitterTree_->Branch("dpt_glb", &dpt_glb_, "dpt_glb/D");
692  splitterTree_->Branch("p1_glb", &p1_glb_, "p1_glb/D");
693  splitterTree_->Branch("p2_glb", &p2_glb_, "p2_glb/D");
694  splitterTree_->Branch("eta1_glb", &eta1_glb_, "eta1_glb/D");
695  splitterTree_->Branch("eta2_glb", &eta2_glb_, "eta2_glb/D");
696  splitterTree_->Branch("deta_glb", &deta_glb_, "deta_glb/D");
697  splitterTree_->Branch("norchi1_glb", &norchi1_glb_, "norchi1_glb/D");
698  splitterTree_->Branch("norchi2_glb", &norchi2_glb_, "norchi2_glb/D");
699 
700  splitterTree_->Branch("d01Err_glb", &d01Err_glb_, "d01Err_glb/D");
701  splitterTree_->Branch("d02Err_glb", &d02Err_glb_, "d02Err_glb/D");
702  splitterTree_->Branch("dz1Err_glb", &dz1Err_glb_, "dz1Err_glb/D");
703  splitterTree_->Branch("dz2Err_glb", &dz2Err_glb_, "dz2Err_glb/D");
704  splitterTree_->Branch("phi1Err_glb", &phi1Err_glb_, "phi1Err_glb/D");
705  splitterTree_->Branch("phi2Err_glb", &phi2Err_glb_, "phi2Err_glb/D");
706  splitterTree_->Branch("theta1Err_glb", &theta1Err_glb_, "theta1Err_glb/D");
707  splitterTree_->Branch("theta2Err_glb", &theta2Err_glb_, "theta2Err_glb/D");
708  splitterTree_->Branch("pt1Err_glb", &pt1Err_glb_, "pt1Err_glb/D");
709  splitterTree_->Branch("pt2Err_glb", &pt2Err_glb_, "pt2Err_glb/D");
710 
711  }
712 
713  totalTracksToAnalyzer_ = 0;
714  goldenCtr = 0;
715  twoTracksCtr = 0;
716  goldenPlusTwoTracksCtr = 0;
717  _passesTracksPlusMuonsCuts = 0;
718 }
719 
720 
721 // ------------ method called once each job just after ending the event loop ------------
722 void CosmicSplitterValidation::endJob() {
723 
724  //std::cout << "totalTracksToAnalyzer: " << totalTracksToAnalyzer_ << std::endl;
725  std::cout << "golden: " << goldenCtr << ", two tracks: " << twoTracksCtr << ", golden+twotracks: " << goldenPlusTwoTracksCtr << ", tracks+muons cuts: " << _passesTracksPlusMuonsCuts << std::endl;
726 }
727 
728 bool CosmicSplitterValidation::is_gold_muon(const edm::Event& e){
729 
730  edm::Handle<reco::MuonCollection> muHandle;
731  e.getByLabel("STAMuons", muHandle);
732  const reco::MuonCollection & muons = *(muHandle.product());
733  // make sure there are 2 muons
734  if ( 2 != muons.size() ) return false;
735 
736  double mudd0=0., mudphi=0., muddsz=0., mudeta=0.;
737  for ( unsigned int bindex = 0; bindex < muons.size(); ++bindex ) {
738  reco::Muon mymuon = muons[bindex];
739  // deprecated in 21x (now outerTrack)
740  //reco::TrackRef mutrackref = mymuon.standAloneMuon();
741  reco::TrackRef mutrackref = mymuon.outerTrack();
742  const reco::Track* mutrack = mutrackref.get();
743  if (0 == bindex){
744  mudd0+=mutrack->d0();
745  mudphi+=mutrack->phi();
746  muddsz+=mutrack->dsz();
747  mudeta+=mymuon.eta();
748  }
749  if (1 == bindex){
750  mudd0-=mutrack->d0();
751  mudphi-=mutrack->phi();
752  muddsz-=mutrack->dsz();
753  mudeta-=mymuon.eta();
754  }
755  }
756  if ((fabs(mudd0)<15.0)&&(fabs(mudphi)<0.045)&&(fabs(muddsz)<20.0)&&(fabs(mudeta)<0.060)) return true;
757  return false;
758 }
759 
760 
761 //define this as a plug-in
762 DEFINE_FWK_MODULE(CosmicSplitterValidation);
reco::TrackBase::p
double p() const
momentum vector magnitude
Definition: TrackBase.h:127
type
type
Definition: HCALResponse.h:21
CosmicSplitterValidation::tfile
edm::Service< TFileService > tfile
Definition: CosmicSplitterValidation.cc:96
reco::TrackBase::d0Error
double d0Error() const
error on d0
Definition: TrackBase.h:209
reco::TrackBase::d0
double d0() const
dxy parameter in perigee convention (d0 = - dxy)
Definition: TrackBase.h:121
reco::TrackBase::normalizedChi2
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:109
CosmicSplitterValidation::CosmicSplitterValidation
CosmicSplitterValidation(const edm::ParameterSet &)
Definition: CosmicSplitterValidation.cc:175
CosmicSplitterValidation::analyze
virtual void analyze(const edm::Event &iEvent, const edm::EventSetup &iSetup) override
Definition: CosmicSplitterValidation.cc:251
reco::TrackBase::theta
double theta() const
polar angle
Definition: TrackBase.h:115
DEFINE_FWK_MODULE
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:17
reco::TrackBase::phi
double phi() const
azimuthal angle of momentum vector
Definition: TrackBase.h:137
TFileService::make
T * make(const Args &...args) const
make new ROOT object
Definition: TFileService.h:64
reco::TrackBase::referencePoint
const Point & referencePoint() const
Reference point on the track.
Definition: TrackBase.h:151
reco::TrackBase::dsz
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:123
reco::MuonCollection
std::vector< Muon > MuonCollection
collection of Muon objects
Definition: MuonFwd.h:9
kFPIX
const int kFPIX
Definition: AlignmentTrackSelector.cc:28
iEvent
int iEvent
Definition: GenABIO.cc:243
reco::TrackBase::eta
double eta() const
pseudorapidity of momentum vector
Definition: TrackBase.h:139
CosmicSplitterValidation::beginJob
virtual void beginJob() override
Definition: CosmicSplitterValidation.cc:557
reco::TrackBase::pt
double pt() const
track transverse momentum
Definition: TrackBase.h:129
reco::TrackBase::ptError
double ptError() const
error on Pt (set to 1000 TeV if charge==0 for safety)
Definition: TrackBase.h:192
reco::TrackBase::phiError
double phiError() const
error on phi
Definition: TrackBase.h:205
reco::Track::recHitsBegin
trackingRecHit_iterator recHitsBegin() const
Iterator to first hit on the track.
Definition: Track.h:62
edm::Event::getByLabel
bool getByLabel(InputTag const &tag, Handle< PROD > &result) const
Definition: Event.h:390
reco::Muon::outerTrack
virtual TrackRef outerTrack() const
reference to Track reconstructed in the muon detector only
Definition: Muon.h:51
reco::LeafCandidate::eta
virtual float eta() const GCC11_FINAL
momentum pseudorapidity
Definition: LeafCandidate.h:159
reco::TrackBase::dz
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:125
reco::TrackBase::dzError
double dzError() const
error on dz
Definition: TrackBase.h:213
testEve_cfg.tracks
tuple tracks
Definition: testEve_cfg.py:39
math::XYZPoint
XYZPointD XYZPoint
point in space with cartesian internal representation
Definition: Point3D.h:12
globalMuons_cfi.globalMuons
tuple globalMuons
Definition: globalMuons_cfi.py:6
edm::Handle::product
T const * product() const
Definition: Handle.h:81
patZpeak.muons
tuple muons
Definition: patZpeak.py:38
CosmicSplitterValidation::endJob
virtual void endJob() override
Definition: CosmicSplitterValidation.cc:722
CosmicSplitterValidation::is_gold_muon
bool is_gold_muon(const edm::Event &e)
Definition: CosmicSplitterValidation.cc:728
gather_cfg.cout
tuple cout
Definition: gather_cfg.py:121
edm::Ref::get
T const * get() const
Returns C++ pointer to the item.
Definition: Ref.h:242
reco::TrackBase::thetaError
double thetaError() const
error on theta
Definition: TrackBase.h:199
reco::Muon::globalTrack
virtual TrackRef globalTrack() const
reference to Track reconstructed in both tracked and muon detector
Definition: Muon.h:54
kBPIX
const int kBPIX
Definition: AlignmentTrackSelector.cc:27
reco::Track::recHitsEnd
trackingRecHit_iterator recHitsEnd() const
Iterator to last hit on the track.
Definition: Track.h:64
Generated for CMSSW Reference Manual by doxygen 1.8.5