CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
ZMuMu_MCanalyzer.cc
Go to the documentation of this file.
1 /* \class ZMuMu_MCanalyzer
2  *
3  * author: Davide Piccolo
4  *
5  * ZMuMu MC analyzer:
6  * check muon reco efficiencies from MC truth,
7  *
8  */
9 
36 #include "TH1.h"
37 #include "TH2.h"
38 #include "TH3.h"
39 #include <vector>
40 
41 using namespace edm;
42 using namespace std;
43 using namespace reco;
44 using namespace reco;
45 using namespace isodeposit;
46 
48 
50 public:
52 private:
53  virtual void analyze(const edm::Event& event, const edm::EventSetup& setup) override;
54  bool check_ifZmumu(const Candidate * dauGen0, const Candidate * dauGen1, const Candidate * dauGen2);
55  float getParticlePt(const int ipart, const Candidate * dauGen0, const Candidate * dauGen1, const Candidate * dauGen2);
56  float getParticleEta(const int ipart, const Candidate * dauGen0, const Candidate * dauGen1, const Candidate * dauGen2);
57  float getParticlePhi(const int ipart, const Candidate * dauGen0, const Candidate * dauGen1, const Candidate * dauGen2);
58  Particle::LorentzVector getParticleP4(const int ipart, const Candidate * dauGen0, const Candidate * dauGen1, const Candidate * dauGen2);
59  virtual void endJob() override;
60 
70 
71  bool bothMuons_;
72 
73  double etamin_, etamax_, ptmin_, massMin_, massMax_, isoMax_;
74 
75  double ptThreshold_, etEcalThreshold_, etHcalThreshold_ ,dRVetoTrk_, dRTrk_, dREcal_ , dRHcal_, alpha_, beta_;
76 
78  string hltPath_;
79  reco::CandidateBaseRef globalMuonCandRef_, trackMuonCandRef_, standAloneMuonCandRef_;
81 
82  // general histograms
83  TH1D *h_trackProbe_eta, *h_trackProbe_pt, *h_staProbe_eta, *h_staProbe_pt, *h_ProbeOk_eta, *h_ProbeOk_pt;
84 
85  // global counters
86  int nGlobalMuonsMatched_passed; // total number of global muons MC matched and passing cuts (and triggered)
87  int nGlobalMuonsMatched_passedIso; // total number of global muons MC matched and passing cuts including Iso
88  int n2GlobalMuonsMatched_passedIso; // total number of Z->2 global muons MC matched and passing cuts including Iso
89  int nStaMuonsMatched_passedIso; // total number of sta only muons MC matched and passing cuts including Iso
90  int nTracksMuonsMatched_passedIso; // total number of tracks only muons MC matched and passing cuts including Iso
91  int n2GlobalMuonsMatched_passedIso2Trg; // total number of Z->2 global muons MC matched and passing cuts including Iso and both triggered
92  int nMu0onlyTriggered; // n. of events zMuMu with mu0 only triggered
93  int nMu1onlyTriggered; // n. of events zMuMu with mu1 only triggered
94 
95  int nZMuMu_matched; // n. of events zMuMu macthed
96  int nZMuSta_matched; // n of events zMuSta macthed
97  int nZMuTrk_matched; // n. of events zMuTrk mathced
98 };
99 
100 
101 template<typename T>
102 double isolation(const T * t, double ptThreshold, double etEcalThreshold, double etHcalThreshold , double dRVetoTrk, double dRTrk, double dREcal , double dRHcal, double alpha, double beta, bool relativeIsolation) {
103  // on 34X:
104  const pat::IsoDeposit * trkIso = t->isoDeposit(pat::TrackIso);
105  // const pat::IsoDeposit * trkIso = t->trackerIsoDeposit();
106  // on 34X
107  const pat::IsoDeposit * ecalIso = t->isoDeposit(pat::EcalIso);
108  // const pat::IsoDeposit * ecalIso = t->ecalIsoDeposit();
109  // on 34X
110  const pat::IsoDeposit * hcalIso = t->isoDeposit(pat::HcalIso);
111  // const pat::IsoDeposit * hcalIso = t->hcalIsoDeposit();
112 
113  Direction dir = Direction(t->eta(), t->phi());
114 
115  pat::IsoDeposit::AbsVetos vetosTrk;
116  vetosTrk.push_back(new ConeVeto( dir, dRVetoTrk ));
117  vetosTrk.push_back(new ThresholdVeto( ptThreshold ));
118 
119  pat::IsoDeposit::AbsVetos vetosEcal;
120  vetosEcal.push_back(new ConeVeto( dir, 0.));
121  vetosEcal.push_back(new ThresholdVeto( etEcalThreshold ));
122 
123  pat::IsoDeposit::AbsVetos vetosHcal;
124  vetosHcal.push_back(new ConeVeto( dir, 0. ));
125  vetosHcal.push_back(new ThresholdVeto( etHcalThreshold ));
126 
127  double isovalueTrk = (trkIso->sumWithin(dRTrk,vetosTrk));
128  double isovalueEcal = (ecalIso->sumWithin(dREcal,vetosEcal));
129  double isovalueHcal = (hcalIso->sumWithin(dRHcal,vetosHcal));
130 
131 
132  double iso = alpha*( ((1+beta)/2*isovalueEcal) + ((1-beta)/2*isovalueHcal) ) + ((1-alpha)*isovalueTrk) ;
133  if(relativeIsolation) iso /= t->pt();
134  return iso;
135 }
136 
137 
138 double candidateIsolation( const reco::Candidate* c, double ptThreshold, double etEcalThreshold, double etHcalThreshold , double dRVetoTrk, double dRTrk, double dREcal , double dRHcal, double alpha, double beta, bool relativeIsolation) {
139  const pat::Muon * mu = dynamic_cast<const pat::Muon *>(&*c->masterClone());
140  if(mu != 0) return isolation(mu, ptThreshold, etEcalThreshold, etHcalThreshold ,dRVetoTrk, dRTrk, dREcal , dRHcal, alpha, beta, relativeIsolation);
141  const pat::GenericParticle * trk = dynamic_cast<const pat::GenericParticle*>(&*c->masterClone());
142  if(trk != 0) return isolation(trk, ptThreshold, etEcalThreshold, etHcalThreshold ,dRVetoTrk, dRTrk, dREcal ,
143  dRHcal, alpha, beta, relativeIsolation);
145  << "Candidate daughter #0 is neither pat::Muons nor pat::GenericParticle\n";
146  return -1;
147 }
148 
149 
150 
151 
162 #include <iostream>
163 #include <iterator>
164 #include <cmath>
165 
167  zMuMuToken_(consumes<CandidateView>(pset.getParameter<InputTag>("zMuMu"))),
168  zMuMuMatchMapToken_(mayConsume<GenParticleMatch>(pset.getParameter<InputTag>("zMuMuMatchMap"))),
169  zMuStandAloneToken_(consumes<CandidateView>(pset.getParameter<InputTag>("zMuStandAlone"))),
170  zMuStandAloneMatchMapToken_(mayConsume<GenParticleMatch>(pset.getParameter<InputTag>("zMuStandAloneMatchMap"))),
171  zMuTrackToken_(consumes<CandidateView>(pset.getParameter<InputTag>("zMuTrack"))),
172  zMuTrackMatchMapToken_(mayConsume<GenParticleMatch>(pset.getParameter<InputTag>("zMuTrackMatchMap"))),
173  muonsToken_(consumes<CandidateView>(pset.getParameter<InputTag>("muons"))),
174  tracksToken_(consumes<CandidateView>(pset.getParameter<InputTag>("tracks"))),
175  genParticlesToken_(consumes<GenParticleCollection>(pset.getParameter<InputTag>("genParticles"))),
176 
177  bothMuons_(pset.getParameter<bool>("bothMuons")),
178  etamin_(pset.getUntrackedParameter<double>("etamin")),
179  etamax_(pset.getUntrackedParameter<double>("etamax")),
180  ptmin_(pset.getUntrackedParameter<double>("ptmin")),
181  massMin_(pset.getUntrackedParameter<double>("zMassMin")),
182  massMax_(pset.getUntrackedParameter<double>("zMassMax")),
183  isoMax_(pset.getUntrackedParameter<double>("isomax")),
184  ptThreshold_(pset.getUntrackedParameter<double>("ptThreshold")),
185  etEcalThreshold_(pset.getUntrackedParameter<double>("etEcalThreshold")),
186  etHcalThreshold_(pset.getUntrackedParameter<double>("etHcalThreshold")),
187  dRVetoTrk_(pset.getUntrackedParameter<double>("deltaRVetoTrk")),
188  dRTrk_(pset.getUntrackedParameter<double>("deltaRTrk")),
189  dREcal_(pset.getUntrackedParameter<double>("deltaREcal")),
190  dRHcal_(pset.getUntrackedParameter<double>("deltaRHcal")),
191  alpha_(pset.getUntrackedParameter<double>("alpha")),
192  beta_(pset.getUntrackedParameter<double>("beta")),
193  relativeIsolation_(pset.getUntrackedParameter<bool>("relativeIsolation")),
194  hltPath_(pset.getUntrackedParameter<std::string >("hltPath")) {
196 
197  // binning of entries array (at moment defined by hand and not in cfg file)
198  double etaRange[8] = {-2.5,-2.,-1.2,-0.8,0.8,1.2,2.,2.5};
199  double ptRange[4] = {20.,40.,60.,100.};
200 
201  // general histograms
202  h_trackProbe_eta = fs->make<TH1D>("trackProbeEta","Eta of tracks",7,etaRange);
203  h_trackProbe_pt = fs->make<TH1D>("trackProbePt","Pt of tracks",3,ptRange);
204  h_staProbe_eta = fs->make<TH1D>("standAloneProbeEta","Eta of standAlone",7,etaRange);
205  h_staProbe_pt = fs->make<TH1D>("standAloneProbePt","Pt of standAlone",3,ptRange);
206  h_ProbeOk_eta = fs->make<TH1D>("probeOkEta","Eta of probe Ok",7,etaRange);
207  h_ProbeOk_pt = fs->make<TH1D>("probeOkPt","Pt of probe ok",3,ptRange);
208 
209  // clear global counters
216  nMu0onlyTriggered = 0;
217  nMu1onlyTriggered = 0;
218  nZMuMu_matched = 0;
219  nZMuSta_matched = 0;
220  nZMuTrk_matched = 0;
221 }
222 
225  Handle<GenParticleMatch> zMuMuMatchMap; //Map of Z made by Mu global + Mu global
226  Handle<CandidateView> zMuStandAlone;
227  Handle<GenParticleMatch> zMuStandAloneMatchMap; //Map of Z made by Mu + StandAlone
228  Handle<CandidateView> zMuTrack;
229  Handle<GenParticleMatch> zMuTrackMatchMap; //Map of Z made by Mu + Track
230  Handle<CandidateView> muons; //Collection of Muons
231  Handle<CandidateView> tracks; //Collection of Tracks
232 
233  Handle<GenParticleCollection> genParticles; // Collection of Generatd Particles
234 
235  event.getByToken(zMuMuToken_, zMuMu);
236  event.getByToken(zMuStandAloneToken_, zMuStandAlone);
237  event.getByToken(zMuTrackToken_, zMuTrack);
238  event.getByToken(genParticlesToken_, genParticles);
239  event.getByToken(muonsToken_, muons);
240  event.getByToken(tracksToken_, tracks);
241 
242  /*
243  cout << "********* zMuMu size : " << zMuMu->size() << endl;
244  cout << "********* zMuStandAlone size : " << zMuStandAlone->size() << endl;
245  cout << "********* zMuTrack size : " << zMuTrack->size() << endl;
246  cout << "********* muons size : " << muons->size() << endl;
247  cout << "********* tracks size : " << tracks->size() << endl;
248  */
249  // std::cout<<"Run-> "<<event.id().run()<<std::endl;
250  // std::cout<<"Event-> "<<event.id().event()<<std::endl;
251 
252 
253  bool zMuMu_found = false;
254 
255  // loop on ZMuMu
256  if (zMuMu->size() > 0 ) {
257  event.getByToken(zMuMuMatchMapToken_, zMuMuMatchMap);
258  for(unsigned int i = 0; i < zMuMu->size(); ++i) { //loop on candidates
259  const Candidate & zMuMuCand = (*zMuMu)[i]; //the candidate
260  CandidateBaseRef zMuMuCandRef = zMuMu->refAt(i);
261 
262  const Candidate * lep0 = zMuMuCand.daughter( 0 );
263  const Candidate * lep1 = zMuMuCand.daughter( 1 );
264  const pat::Muon & muonDau0 = dynamic_cast<const pat::Muon &>(*lep0->masterClone());
265  // double trkiso0 = muonDau0.trackIso();
266  const pat::Muon & muonDau1 = dynamic_cast<const pat::Muon &>(*lep1->masterClone());
267  //double trkiso1 = muonDau1.trackIso();
268 
271 
272  double pt0 = zMuMuCand.daughter(0)->pt();
273  double pt1 = zMuMuCand.daughter(1)->pt();
274  double eta0 = zMuMuCand.daughter(0)->eta();
275  double eta1 = zMuMuCand.daughter(1)->eta();
276  double mass = zMuMuCand.mass();
277 
278  // HLT match
279  const pat::TriggerObjectStandAloneCollection mu0HLTMatches =
281  const pat::TriggerObjectStandAloneCollection mu1HLTMatches =
283 
284  bool trig0found = false;
285  bool trig1found = false;
286  if( mu0HLTMatches.size()>0 )
287  trig0found = true;
288  if( mu1HLTMatches.size()>0 )
289  trig1found = true;
290 
291  GenParticleRef zMuMuMatch = (*zMuMuMatchMap)[zMuMuCandRef];
292  if(zMuMuMatch.isNonnull()) { // ZMuMu matched
293  zMuMu_found = true;
294  nZMuMu_matched++;
295  if (pt0>ptmin_ && pt1>ptmin_ && abs(eta0)>etamin_ && abs(eta1) >etamin_ && abs(eta0)<etamax_ && abs(eta1) <etamax_ && mass >massMin_ && mass < massMax_ && (trig0found || trig1found)) { // kinematic and trigger cuts passed
296  nGlobalMuonsMatched_passed++; // first global Muon passed kine cuts
297  nGlobalMuonsMatched_passed++; // second global muon passsed kine cuts
298  if (iso0<isoMax_) nGlobalMuonsMatched_passedIso++; // first global muon passed the iso cut
299  if (iso1<isoMax_) nGlobalMuonsMatched_passedIso++; // second global muon passed the iso cut
300  if (iso0<isoMax_ && iso1<isoMax_) {
301  n2GlobalMuonsMatched_passedIso++; // both muons passed iso cut
302  if (trig0found && trig1found) n2GlobalMuonsMatched_passedIso2Trg++; // both muons have HLT
303  if (trig0found && !trig1found) nMu0onlyTriggered++;
304  if (trig1found && !trig0found) nMu1onlyTriggered++;
305  // histograms vs eta and pt
306  if (trig1found) { // check efficiency of muon0 not imposing the trigger on it
307  h_trackProbe_eta->Fill(eta0);
308  h_trackProbe_pt->Fill(pt0);
309  h_staProbe_eta->Fill(eta0);
310  h_staProbe_pt->Fill(pt0);
311  h_ProbeOk_eta->Fill(eta0);
312  h_ProbeOk_pt->Fill(pt0);
313  }
314  if (trig0found) { // check efficiency of muon1 not imposing the trigger on it
315  h_trackProbe_eta->Fill(eta1);
316  h_staProbe_eta->Fill(eta1);
317  h_trackProbe_pt->Fill(pt1);
318  h_staProbe_pt->Fill(pt1);
319  h_ProbeOk_eta->Fill(eta1);
320  h_ProbeOk_pt->Fill(pt1);
321  }
322  }
323  }
324  } // end MC match
325 
326  } // end loop on ZMuMu cand
327  } // end if ZMuMu size > 0
328 
329  // loop on ZMuSta
330  bool zMuSta_found = false;
331  if (!zMuMu_found && zMuStandAlone->size() > 0 ) {
332  event.getByToken(zMuStandAloneMatchMapToken_, zMuStandAloneMatchMap);
333  for(unsigned int i = 0; i < zMuStandAlone->size(); ++i) { //loop on candidates
334  const Candidate & zMuStandAloneCand = (*zMuStandAlone)[i]; //the candidate
335  CandidateBaseRef zMuStandAloneCandRef = zMuStandAlone->refAt(i);
336  GenParticleRef zMuStandAloneMatch = (*zMuStandAloneMatchMap)[zMuStandAloneCandRef];
337 
338  const Candidate * lep0 = zMuStandAloneCand.daughter( 0 );
339  const Candidate * lep1 = zMuStandAloneCand.daughter( 1 );
340  const pat::Muon & muonDau0 = dynamic_cast<const pat::Muon &>(*lep0->masterClone());
341  //double trkiso0 = muonDau0.trackIso();
342  // const pat::Muon & muonDau1 = dynamic_cast<const pat::Muon &>(*lep1->masterClone());
343  //double trkiso1 = muonDau1.trackIso();
344 
347 
348  double pt0 = zMuStandAloneCand.daughter(0)->pt();
349  double pt1 = zMuStandAloneCand.daughter(1)->pt();
350  double eta0 = zMuStandAloneCand.daughter(0)->eta();
351  double eta1 = zMuStandAloneCand.daughter(1)->eta();
352  double mass = zMuStandAloneCand.mass();
353 
354  // HLT match (check just dau0 the global)
355  const pat::TriggerObjectStandAloneCollection mu0HLTMatches =
357 
358  bool trig0found = false;
359  if( mu0HLTMatches.size()>0 )
360  trig0found = true;
361 
362  if(zMuStandAloneMatch.isNonnull()) { // ZMuStandAlone matched
363  zMuSta_found = true;
364  nZMuSta_matched++;
365  if (pt0>ptmin_ && pt1>ptmin_ && abs(eta0)>etamin_ && abs(eta1)>etamin_ && abs(eta0)<etamax_ && abs(eta1) <etamax_ && mass >massMin_ &&
366  mass < massMax_ && iso0<isoMax_ && iso1 < isoMax_ && trig0found) { // all cuts and trigger passed
368  // histograms vs eta and pt
369  h_staProbe_eta->Fill(eta1);
370  h_staProbe_pt->Fill(pt1);
371  }
372  } // end MC match
373  } // end loop on ZMuStandAlone cand
374  } // end if ZMuStandAlone size > 0
375 
376 
377  // loop on ZMuTrack
378  if (!zMuMu_found && !zMuSta_found && zMuTrack->size() > 0 ) {
379  event.getByToken(zMuTrackMatchMapToken_, zMuTrackMatchMap);
380  for(unsigned int i = 0; i < zMuTrack->size(); ++i) { //loop on candidates
381  const Candidate & zMuTrackCand = (*zMuTrack)[i]; //the candidate
382  CandidateBaseRef zMuTrackCandRef = zMuTrack->refAt(i);
383  const Candidate * lep0 = zMuTrackCand.daughter( 0 );
384  const Candidate * lep1 = zMuTrackCand.daughter( 1 );
385  const pat::Muon & muonDau0 = dynamic_cast<const pat::Muon &>(*lep0->masterClone());
386  //double trkiso0 = muonDau0.trackIso();
387  //const pat::GenericParticle & trackDau1 = dynamic_cast<const pat::GenericParticle &>(*lep1->masterClone());
388  //double trkiso1 = trackDau1.trackIso();
391 
392 
393  double pt0 = zMuTrackCand.daughter(0)->pt();
394  double pt1 = zMuTrackCand.daughter(1)->pt();
395  double eta0 = zMuTrackCand.daughter(0)->eta();
396  double eta1 = zMuTrackCand.daughter(1)->eta();
397  double mass = zMuTrackCand.mass();
398 
399  // HLT match (check just dau0 the global)
400  const pat::TriggerObjectStandAloneCollection mu0HLTMatches =
402 
403  bool trig0found = false;
404  if( mu0HLTMatches.size()>0 )
405  trig0found = true;
406 
407  GenParticleRef zMuTrackMatch = (*zMuTrackMatchMap)[zMuTrackCandRef];
408  if(zMuTrackMatch.isNonnull()) { // ZMuTrack matched
409  nZMuTrk_matched++;
410  if (pt0>ptmin_ && pt1>ptmin_ && abs(eta0)>etamin_ && abs(eta1)>etamin_ && abs(eta0)<etamax_ && abs(eta1) <etamax_ && mass >massMin_ &&
411  mass < massMax_ && iso0<isoMax_ && iso1 < isoMax_ && trig0found) { // all cuts and trigger passed
413  // histograms vs eta and pt
414  h_trackProbe_eta->Fill(eta1);
415  h_trackProbe_pt->Fill(pt1);
416  }
417  } // end MC match
418  } // end loop on ZMuTrack cand
419  } // end if ZMuTrack size > 0
420 
421 } // end analyze
422 
423 bool ZMuMu_MCanalyzer::check_ifZmumu(const Candidate * dauGen0, const Candidate * dauGen1, const Candidate * dauGen2)
424 {
425  int partId0 = dauGen0->pdgId();
426  int partId1 = dauGen1->pdgId();
427  int partId2 = dauGen2->pdgId();
428  bool muplusFound=false;
429  bool muminusFound=false;
430  bool ZFound=false;
431  if (partId0==13 || partId1==13 || partId2==13) muminusFound=true;
432  if (partId0==-13 || partId1==-13 || partId2==-13) muplusFound=true;
433  if (partId0==23 || partId1==23 || partId2==23) ZFound=true;
434  return muplusFound*muminusFound*ZFound;
435 }
436 
437 float ZMuMu_MCanalyzer::getParticlePt(const int ipart, const Candidate * dauGen0, const Candidate * dauGen1, const Candidate * dauGen2)
438 {
439  int partId0 = dauGen0->pdgId();
440  int partId1 = dauGen1->pdgId();
441  int partId2 = dauGen2->pdgId();
442  float ptpart=0.;
443  if (partId0 == ipart) {
444  for(unsigned int k = 0; k < dauGen0->numberOfDaughters(); ++k) {
445  const Candidate * dauMuGen = dauGen0->daughter(k);
446  if(dauMuGen->pdgId() == ipart && dauMuGen->status() ==1) {
447  ptpart = dauMuGen->pt();
448  }
449  }
450  }
451  if (partId1 == ipart) {
452  for(unsigned int k = 0; k < dauGen1->numberOfDaughters(); ++k) {
453  const Candidate * dauMuGen = dauGen1->daughter(k);
454  if(dauMuGen->pdgId() == ipart && dauMuGen->status() ==1) {
455  ptpart = dauMuGen->pt();
456  }
457  }
458  }
459  if (partId2 == ipart) {
460  for(unsigned int k = 0; k < dauGen2->numberOfDaughters(); ++k) {
461  const Candidate * dauMuGen = dauGen2->daughter(k);
462  if(abs(dauMuGen->pdgId()) == ipart && dauMuGen->status() ==1) {
463  ptpart = dauMuGen->pt();
464  }
465  }
466  }
467  return ptpart;
468 }
469 
470 float ZMuMu_MCanalyzer::getParticleEta(const int ipart, const Candidate * dauGen0, const Candidate * dauGen1, const Candidate * dauGen2)
471 {
472  int partId0 = dauGen0->pdgId();
473  int partId1 = dauGen1->pdgId();
474  int partId2 = dauGen2->pdgId();
475  float etapart=0.;
476  if (partId0 == ipart) {
477  for(unsigned int k = 0; k < dauGen0->numberOfDaughters(); ++k) {
478  const Candidate * dauMuGen = dauGen0->daughter(k);
479  if(dauMuGen->pdgId() == ipart && dauMuGen->status() ==1) {
480  etapart = dauMuGen->eta();
481  }
482  }
483  }
484  if (partId1 == ipart) {
485  for(unsigned int k = 0; k < dauGen1->numberOfDaughters(); ++k) {
486  const Candidate * dauMuGen = dauGen1->daughter(k);
487  if(dauMuGen->pdgId() == ipart && dauMuGen->status() ==1) {
488  etapart = dauMuGen->eta();
489  }
490  }
491  }
492  if (partId2 == ipart) {
493  for(unsigned int k = 0; k < dauGen2->numberOfDaughters(); ++k) {
494  const Candidate * dauMuGen = dauGen2->daughter(k);
495  if(abs(dauMuGen->pdgId()) == ipart && dauMuGen->status() ==1) {
496  etapart = dauMuGen->eta();
497  }
498  }
499  }
500  return etapart;
501 }
502 
503 float ZMuMu_MCanalyzer::getParticlePhi(const int ipart, const Candidate * dauGen0, const Candidate * dauGen1, const Candidate * dauGen2)
504 {
505  int partId0 = dauGen0->pdgId();
506  int partId1 = dauGen1->pdgId();
507  int partId2 = dauGen2->pdgId();
508  float phipart=0.;
509  if (partId0 == ipart) {
510  for(unsigned int k = 0; k < dauGen0->numberOfDaughters(); ++k) {
511  const Candidate * dauMuGen = dauGen0->daughter(k);
512  if(dauMuGen->pdgId() == ipart && dauMuGen->status() ==1) {
513  phipart = dauMuGen->phi();
514  }
515  }
516  }
517  if (partId1 == ipart) {
518  for(unsigned int k = 0; k < dauGen1->numberOfDaughters(); ++k) {
519  const Candidate * dauMuGen = dauGen1->daughter(k);
520  if(dauMuGen->pdgId() == ipart && dauMuGen->status() ==1) {
521  phipart = dauMuGen->phi();
522  }
523  }
524  }
525  if (partId2 == ipart) {
526  for(unsigned int k = 0; k < dauGen2->numberOfDaughters(); ++k) {
527  const Candidate * dauMuGen = dauGen2->daughter(k);
528  if(abs(dauMuGen->pdgId()) == ipart && dauMuGen->status() ==1) {
529  phipart = dauMuGen->phi();
530  }
531  }
532  }
533  return phipart;
534 }
535 
536 Particle::LorentzVector ZMuMu_MCanalyzer::getParticleP4(const int ipart, const Candidate * dauGen0, const Candidate * dauGen1, const Candidate * dauGen2)
537 {
538  int partId0 = dauGen0->pdgId();
539  int partId1 = dauGen1->pdgId();
540  int partId2 = dauGen2->pdgId();
541  Particle::LorentzVector p4part(0.,0.,0.,0.);
542  if (partId0 == ipart) {
543  for(unsigned int k = 0; k < dauGen0->numberOfDaughters(); ++k) {
544  const Candidate * dauMuGen = dauGen0->daughter(k);
545  if(dauMuGen->pdgId() == ipart && dauMuGen->status() ==1) {
546  p4part = dauMuGen->p4();
547  }
548  }
549  }
550  if (partId1 == ipart) {
551  for(unsigned int k = 0; k < dauGen1->numberOfDaughters(); ++k) {
552  const Candidate * dauMuGen = dauGen1->daughter(k);
553  if(dauMuGen->pdgId() == ipart && dauMuGen->status() ==1) {
554  p4part = dauMuGen->p4();
555  }
556  }
557  }
558  if (partId2 == ipart) {
559  for(unsigned int k = 0; k < dauGen2->numberOfDaughters(); ++k) {
560  const Candidate * dauMuGen = dauGen2->daughter(k);
561  if(abs(dauMuGen->pdgId()) == ipart && dauMuGen->status() ==1) {
562  p4part = dauMuGen->p4();
563  }
564  }
565  }
566  return p4part;
567 }
568 
569 
570 
572 
573 
575  double err_effIso = sqrt(eff_Iso*(1-eff_Iso)/nGlobalMuonsMatched_passed);
576 
580  double effSta_afterIso = (2*n2GlobalMuonsMatched_passedIso2Trg+nMu0onlyTriggered+nMu1onlyTriggered)/n1_afterIso;
581  double effTrk_afterIso = (2*n2GlobalMuonsMatched_passedIso2Trg+nMu0onlyTriggered+nMu1onlyTriggered)/n2_afterIso;
582  double effHLT_afterIso = (2.* n2GlobalMuonsMatched_passedIso2Trg)/(2.* n2GlobalMuonsMatched_passedIso2Trg + nMu0onlyTriggered + nMu1onlyTriggered);
583  double err_effHLT_afterIso= sqrt( effHLT_afterIso * (1 - effHLT_afterIso)/nGLB_afterIso);
584  double err_effsta_afterIso = sqrt(effSta_afterIso*(1-effSta_afterIso)/n1_afterIso);
585  double err_efftrk_afterIso = sqrt(effTrk_afterIso*(1-effTrk_afterIso)/n2_afterIso);
586 
587 
588  cout << "------------------------------------ Counters --------------------------------" << endl;
589 
590  cout << "number of events zMuMu matched " << nZMuMu_matched << endl;
591  cout << "number of events zMuSta matched " << nZMuSta_matched << endl;
592  cout << "number of events zMuTk matched " << nZMuTrk_matched << endl;
593  cout << "number of events zMuMu with mu0 only triggered " << nMu0onlyTriggered << endl;
594  cout << "number of events zMuMu with mu1 only triggered " << nMu1onlyTriggered << endl;
595  cout << "=========================================" << endl;
596  cout << "n. of global muons MC matched and passing cuts: " << nGlobalMuonsMatched_passed << endl;
597  cout << "n. of global muons MC matched and passing also Iso cut: " << nGlobalMuonsMatched_passedIso << endl;
598  cout << "n. of Z -> 2 global muons MC matched and passing ALL cuts: " << n2GlobalMuonsMatched_passedIso << endl;
599  cout << "n. of ZMuSta MC matched and passing ALL cuts: " << nStaMuonsMatched_passedIso << endl;
600  cout << "n. of ZmuTrck MC matched and passing ALL cuts: " << nTracksMuonsMatched_passedIso << endl;
601  cout << "n. of Z -> 2 global muons MC matched and passing ALL cuts and both triggered: " << n2GlobalMuonsMatched_passedIso2Trg << endl;
602  cout << "=================================================================================" << endl;
603  cout << "Iso efficiency: " << eff_Iso << " +/- " << err_effIso << endl;
604  cout << "HLT efficiency: " << effHLT_afterIso << " +/- " << err_effHLT_afterIso << endl;
605  cout << "eff StandAlone (after Isocut) : " << effSta_afterIso << "+/-" << err_effsta_afterIso << endl;
606  cout << "eff Tracker (after Isocut) : " << effTrk_afterIso << "+/-" << err_efftrk_afterIso << endl;
607 
608 }
609 
611 
613 
const double beta
std::vector< GenParticle > GenParticleCollection
collection of GenParticles
int i
Definition: DBlmapReader.cc:9
ZMuMu_MCanalyzer(const edm::ParameterSet &pset)
float alpha
Definition: AMPTWrapper.h:95
virtual const Candidate * daughter(size_type i) const =0
return daughter at a given position, i = 0, ... numberOfDaughters() - 1 (read only mode) ...
double candidateIsolation(const reco::Candidate *c, double ptThreshold, double etEcalThreshold, double etHcalThreshold, double dRVetoTrk, double dRTrk, double dREcal, double dRHcal, double alpha, double beta, bool relativeIsolation)
virtual float mass() const =0
mass
virtual float eta() const =0
momentum pseudorapidity
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:17
EDGetTokenT< CandidateView > zMuStandAloneToken_
EDGetTokenT< CandidateView > zMuMuToken_
virtual int status() const =0
status word
virtual float phi() const =0
momentum azimuthal angle
float getParticlePhi(const int ipart, const Candidate *dauGen0, const Candidate *dauGen1, const Candidate *dauGen2)
T * make(const Args &...args) const
make new ROOT object
Definition: TFileService.h:64
std::vector< TriggerObjectStandAlone > TriggerObjectStandAloneCollection
Collection of TriggerObjectStandAlone.
virtual void analyze(const edm::Event &event, const edm::EventSetup &setup) override
double sumWithin(double coneSize, const AbsVetos &vetos=AbsVetos(), bool skipDepositVeto=false) const
Definition: IsoDeposit.cc:138
bool check_ifZmumu(const Candidate *dauGen0, const Candidate *dauGen1, const Candidate *dauGen2)
bool isNonnull() const
Checks for non-null.
Definition: Ref.h:250
Analysis-level Generic Particle class (e.g. for hadron or muon not fully reconstructed) ...
virtual size_type numberOfDaughters() const =0
number of daughters
virtual float pt() const =0
transverse momentum
EDGetTokenT< CandidateView > muonsToken_
edm::ValueMap< float > IsolationCollection
EDGetTokenT< GenParticleMatch > zMuTrackMatchMapToken_
T sqrt(T t)
Definition: SSEVec.h:48
reco::CandidateBaseRef trackMuonCandRef_
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
const int mu
Definition: Constants.h:22
const TriggerObjectStandAloneCollection triggerObjectMatchesByPath(const std::string &namePath, const bool pathLastFilterAccepted=false, const bool pathL3FilterAccepted=true) const
Definition: PATObject.h:596
How EventSelector::AcceptEvent() decides whether to accept an event for output otherwise it is excluding the probing of A single or multiple positive and the trigger will pass if any such matching triggers are PASS or EXCEPTION[A criterion thatmatches no triggers at all is detected and causes a throw.] A single negative with an expectation of appropriate bit checking in the decision and the trigger will pass if any such matching triggers are FAIL or EXCEPTION A wildcarded negative criterion that matches more than one trigger in the trigger but the state exists so we define the behavior If all triggers are the negative crieriion will lead to accepting the event(this again matches the behavior of"!*"before the partial wildcard feature was incorporated).The per-event"cost"of each negative criterion with multiple relevant triggers is about the same as!*was in the past
EDGetTokenT< GenParticleCollection > genParticlesToken_
int k[5][pyjets_maxn]
virtual int pdgId() const =0
PDG identifier.
int n2GlobalMuonsMatched_passedIso2Trg
tuple tracks
Definition: testEve_cfg.py:39
EDGetTokenT< GenParticleMatch > zMuMuMatchMapToken_
OverlapChecker overlap_
tuple muons
Definition: patZpeak.py:38
float getParticleEta(const int ipart, const Candidate *dauGen0, const Candidate *dauGen1, const Candidate *dauGen2)
tuple cout
Definition: gather_cfg.py:121
EDGetTokenT< CandidateView > zMuTrackToken_
dbl *** dir
Definition: mlp_gen.cc:35
EDGetTokenT< CandidateView > tracksToken_
isodeposit::AbsVetos AbsVetos
Definition: IsoDeposit.h:51
long double T
math::XYZTLorentzVector LorentzVector
Lorentz vector.
Definition: Particle.h:27
Analysis-level muon class.
Definition: Muon.h:50
void setup(std::vector< TH2F > &depth, std::string name, std::string units="")
EDGetTokenT< GenParticleMatch > zMuStandAloneMatchMapToken_
tuple zMuMu
zMuMu vector of PSet is common to all categories except zMuTrk category
float getParticlePt(const int ipart, const Candidate *dauGen0, const Candidate *dauGen1, const Candidate *dauGen2)
virtual void endJob() override
Particle::LorentzVector getParticleP4(const int ipart, const Candidate *dauGen0, const Candidate *dauGen1, const Candidate *dauGen2)
virtual const LorentzVector & p4() const =0
four-momentum Lorentz vector
virtual const CandidateBaseRef & masterClone() const =0