test
CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
ZMuMu_efficiencyAnalyzer.cc
Go to the documentation of this file.
1 /* \class ZMuMu_efficieyAnalyzer
2  *
3  * author: Davide Piccolo
4  *
5  * ZMuMu efficiency analyzer:
6  * check muon reco efficiencies from MC truth,
7  *
8  */
9 
34 
35 #include "TH1.h"
36 #include "TH2.h"
37 #include "TH3.h"
38 #include <vector>
39 
40 using namespace edm;
41 using namespace std;
42 using namespace reco;
43 
45 
47 public:
49 private:
50  virtual void analyze(const edm::Event& event, const edm::EventSetup& setup) override;
51  bool check_ifZmumu(const Candidate * dauGen0, const Candidate * dauGen1, const Candidate * dauGen2);
52  float getParticlePt(const int ipart, const Candidate * dauGen0, const Candidate * dauGen1, const Candidate * dauGen2);
53  float getParticleEta(const int ipart, const Candidate * dauGen0, const Candidate * dauGen1, const Candidate * dauGen2);
54  float getParticlePhi(const int ipart, const Candidate * dauGen0, const Candidate * dauGen1, const Candidate * dauGen2);
55  Particle::LorentzVector getParticleP4(const int ipart, const Candidate * dauGen0, const Candidate * dauGen1, const Candidate * dauGen2);
56  virtual void endJob() override;
57 
68 
69  bool bothMuons_;
70 
71  double etamax_, ptmin_, massMin_, massMax_, isoMax_;
72 
73  // binning of entries array (at moment defined by hand and not in cfg file)
74  unsigned int etaBins;
75  unsigned int ptBins;
76  double etaRange[7];
77  double ptRange[5];
78 
79  reco::CandidateBaseRef globalMuonCandRef_, trackMuonCandRef_, standAloneMuonCandRef_;
81 
82  // general histograms
83  TH1D *h_zmm_mass, *h_zmm2HLT_mass;
84  TH1D *h_zmm1HLTplus_mass, *h_zmmNotIsoplus_mass, *h_zmsplus_mass, *h_zmtplus_mass;
85  TH1D *h_zmm1HLTminus_mass, *h_zmmNotIsominus_mass, *h_zmsminus_mass, *h_zmtminus_mass;
86 
87  // global counters
88  int nGlobalMuonsMatched_passed; // total number of global muons MC matched and passing cuts (and triggered)
89 
90  vector<TH1D *> hmumu2HLTplus_eta, hmumu1HLTplus_eta, hmustaplus_eta, hmutrackplus_eta, hmumuNotIsoplus_eta;
91  vector<TH1D *> hmumu2HLTplus_pt, hmumu1HLTplus_pt, hmustaplus_pt, hmutrackplus_pt, hmumuNotIsoplus_pt;
92  vector<TH1D *> hmumu2HLTminus_eta, hmumu1HLTminus_eta, hmustaminus_eta, hmutrackminus_eta, hmumuNotIsominus_eta;
93  vector<TH1D *> hmumu2HLTminus_pt, hmumu1HLTminus_pt, hmustaminus_pt, hmutrackminus_pt, hmumuNotIsominus_pt;
94 };
95 
101 #include <iostream>
102 #include <iterator>
103 #include <cmath>
104 
106  zMuMuToken_(consumes<CandidateView>(pset.getParameter<InputTag>("zMuMu"))),
107  zMuMuMatchMapToken_(mayConsume<GenParticleMatch>(pset.getParameter<InputTag>("zMuMuMatchMap"))),
108  zMuStandAloneToken_(consumes<CandidateView>(pset.getParameter<InputTag>("zMuStandAlone"))),
109  zMuStandAloneMatchMapToken_(mayConsume<GenParticleMatch>(pset.getParameter<InputTag>("zMuStandAloneMatchMap"))),
110  zMuTrackToken_(consumes<CandidateView>(pset.getParameter<InputTag>("zMuTrack"))),
111  zMuTrackMatchMapToken_(mayConsume<GenParticleMatch>(pset.getParameter<InputTag>("zMuTrackMatchMap"))),
112  muonsToken_(consumes<CandidateView>(pset.getParameter<InputTag>("muons"))),
113  tracksToken_(consumes<CandidateView>(pset.getParameter<InputTag>("tracks"))),
114  genParticlesToken_(consumes<GenParticleCollection>(pset.getParameter<InputTag>("genParticles"))),
115  primaryVerticesToken_(consumes<VertexCollection>(pset.getParameter<InputTag>("primaryVertices"))),
116 
117  bothMuons_(pset.getParameter<bool>("bothMuons")),
118 
119  etamax_(pset.getUntrackedParameter<double>("etamax")),
120  ptmin_(pset.getUntrackedParameter<double>("ptmin")),
121  massMin_(pset.getUntrackedParameter<double>("zMassMin")),
122  massMax_(pset.getUntrackedParameter<double>("zMassMax")),
123  isoMax_(pset.getUntrackedParameter<double>("isomax")) {
125 
126  // general histograms
127  h_zmm_mass = fs->make<TH1D>("zmm_mass","zmumu mass",100,0.,200.);
128  h_zmm2HLT_mass = fs->make<TH1D>("zmm2HLT_mass","zmumu 2HLT mass",100,0.,200.);
129  h_zmm1HLTplus_mass = fs->make<TH1D>("zmm1HLTplus_mass","zmumu 1HLT plus mass",100,0.,200.);
130  h_zmmNotIsoplus_mass = fs->make<TH1D>("zmmNotIsoplus_mass","zmumu a least One Not Iso plus mass",100,0.,200.);
131  h_zmsplus_mass = fs->make<TH1D>("zmsplus_mass","zmusta plus mass",100,0.,200.);
132  h_zmtplus_mass = fs->make<TH1D>("zmtplus_mass","zmutrack plus mass",100,0.,200.);
133  h_zmm1HLTminus_mass = fs->make<TH1D>("zmm1HLTminus_mass","zmumu 1HLT minus mass",100,0.,200.);
134  h_zmmNotIsominus_mass = fs->make<TH1D>("zmmNotIsominus_mass","zmumu a least One Not Iso minus mass",100,0.,200.);
135  h_zmsminus_mass = fs->make<TH1D>("zmsminus_mass","zmusta minus mass",100,0.,200.);
136  h_zmtminus_mass = fs->make<TH1D>("zmtminus_mass","zmutrack minus mass",100,0.,200.);
137 
138  cout << "primo" << endl;
139  // creating histograms for each Pt, eta interval
140 
141  TFileDirectory etaDirectory = fs->mkdir("etaIntervals"); // in this directory will be saved all the histos of different eta intervals
142  TFileDirectory ptDirectory = fs->mkdir("ptIntervals"); // in this directory will be saved all the histos of different pt intervals
143 
144  // binning of entries array (at moment defined by hand and not in cfg file)
145  etaBins = 6;
146  ptBins = 4;
147  double etaRangeTmp[7] = {-2.,-1.2,-0.8,0.,0.8,1.2,2.};
148  double ptRangeTmp[5] = {20.,40.,60.,80.,100.};
149  for (unsigned int i=0;i<=etaBins;i++) etaRange[i] = etaRangeTmp[i];
150  for (unsigned int i=0;i<=ptBins;i++) ptRange[i] = ptRangeTmp[i];
151 
152  // eta histograms creation
153  cout << "eta istograms creation " << endl;
154 
155  for (unsigned int i=0;i<etaBins;i++) {
156  cout << " bin eta plus " << i << endl;
157  // muon plus
158  double range0 = etaRange[i];
159  double range1= etaRange[i+1];
160  char ap[30], bp[50];
161  sprintf(ap,"zmumu2HLTplus_etaRange%d",i);
162  sprintf(bp,"zmumu2HLT plus mass eta Range %f to %f",range0,range1);
163  cout << ap << " " << bp << endl;
164  hmumu2HLTplus_eta.push_back(etaDirectory.make<TH1D>(ap,bp,200,0.,200.));
165  sprintf(ap,"zmumu1HLTplus_etaRange%d",i);
166  sprintf(bp,"zmumu1HLT plus mass eta Range %f to %f",range0,range1);
167  cout << ap << " " << bp << endl;
168  hmumu1HLTplus_eta.push_back(etaDirectory.make<TH1D>(ap,bp,200,0.,200.));
169  sprintf(ap,"zmustaplus_etaRange%d",i);
170  sprintf(bp,"zmusta plus mass eta Range %f to %f",range0,range1);
171  cout << ap << " " << bp << endl;
172  hmustaplus_eta.push_back(etaDirectory.make<TH1D>(ap,bp,50,0.,200.));
173  sprintf(ap,"zmutrackplus_etaRange%d",i);
174  sprintf(bp,"zmutrack plus mass eta Range %f to %f",range0,range1);
175  cout << ap << " " << bp << endl;
176  hmutrackplus_eta.push_back(etaDirectory.make<TH1D>(ap,bp,100,0.,200.));
177  sprintf(ap,"zmumuNotIsoplus_etaRange%d",i);
178  sprintf(bp,"zmumuNotIso plus mass eta Range %f to %f",range0,range1);
179  cout << ap << " " << bp << endl;
180  hmumuNotIsoplus_eta.push_back(etaDirectory.make<TH1D>(ap,bp,100,0.,200.));
181  // muon minus
182  cout << " bin eta minus " << i << endl;
183  char am[30], bm[50];
184  sprintf(am,"zmumu2HLTminus_etaRange%d",i);
185  sprintf(bm,"zmumu2HLT minus mass eta Range %f to %f",range0,range1);
186  cout << am << " " << bm << endl;
187  hmumu2HLTminus_eta.push_back(etaDirectory.make<TH1D>(am,bm,200,0.,200.));
188  sprintf(am,"zmumu1HLTminus_etaRange%d",i);
189  sprintf(bm,"zmumu1HLT minus mass eta Range %f to %f",range0,range1);
190  cout << am << " " << bm << endl;
191  hmumu1HLTminus_eta.push_back(etaDirectory.make<TH1D>(am,bm,200,0.,200.));
192  sprintf(am,"zmustaminus_etaRange%d",i);
193  sprintf(bm,"zmusta minus mass eta Range %f to %f",range0,range1);
194  cout << am << " " << bm << endl;
195  hmustaminus_eta.push_back(etaDirectory.make<TH1D>(am,bm,50,0.,200.));
196  sprintf(am,"zmutrackminus_etaRange%d",i);
197  sprintf(bm,"zmutrack minus mass eta Range %f to %f",range0,range1);
198  cout << am << " " << bm << endl;
199  hmutrackminus_eta.push_back(etaDirectory.make<TH1D>(am,bm,100,0.,200.));
200  sprintf(am,"zmumuNotIsominus_etaRange%d",i);
201  sprintf(bm,"zmumuNotIso minus mass eta Range %f to %f",range0,range1);
202  cout << am << " " << bm << endl;
203  hmumuNotIsominus_eta.push_back(etaDirectory.make<TH1D>(am,bm,100,0.,200.));
204  }
205 
206  // pt histograms creation
207  cout << "pt istograms creation " << endl;
208 
209  for (unsigned int i=0;i<ptBins;i++) {
210  double range0 = ptRange[i];
211  double range1= ptRange[i+1];
212  // muon plus
213  cout << " bin pt plus " << i << endl;
214  char ap1[30], bp1[50];
215  sprintf(ap1,"zmumu2HLTplus_ptRange%d",i);
216  sprintf(bp1,"zmumu2HLT plus mass pt Range %f to %f",range0,range1);
217  cout << ap1 << " " << bp1 << endl;
218  hmumu2HLTplus_pt.push_back(ptDirectory.make<TH1D>(ap1,bp1,200,0.,200.));
219  sprintf(ap1,"zmumu1HLTplus_ptRange%d",i);
220  sprintf(bp1,"zmumu1HLT plus mass pt Range %f to %f",range0,range1);
221  cout << ap1 << " " << bp1 << endl;
222  hmumu1HLTplus_pt.push_back(ptDirectory.make<TH1D>(ap1,bp1,200,0.,200.));
223  sprintf(ap1,"zmustaplus_ptRange%d",i);
224  sprintf(bp1,"zmusta plus mass pt Range %f to %f",range0,range1);
225  cout << ap1 << " " << bp1 << endl;
226  hmustaplus_pt.push_back(ptDirectory.make<TH1D>(ap1,bp1,50,0.,200.));
227  sprintf(ap1,"zmutrackplus_ptRange%d",i);
228  sprintf(bp1,"zmutrack plus mass pt Range %f to %f",range0,range1);
229  cout << ap1 << " " << bp1 << endl;
230  hmutrackplus_pt.push_back(ptDirectory.make<TH1D>(ap1,bp1,100,0.,200.));
231  sprintf(ap1,"zmumuNotIsoplus_ptRange%d",i);
232  sprintf(bp1,"zmumuNotIso plus mass pt Range %f to %f",range0,range1);
233  cout << ap1 << " " << bp1 << endl;
234  hmumuNotIsoplus_pt.push_back(ptDirectory.make<TH1D>(ap1,bp1,100,0.,200.));
235  // muon minus
236  cout << " bin pt minus " << i << endl;
237  char am1[30], bm1[50];
238  sprintf(am1,"zmumu2HLTminus_ptRange%d",i);
239  sprintf(bm1,"zmumu2HLT minus mass pt Range %f to %f",range0,range1);
240  cout << am1 << " " << bm1 << endl;
241  hmumu2HLTminus_pt.push_back(ptDirectory.make<TH1D>(am1,bm1,200,0.,200.));
242  sprintf(am1,"zmumu1HLTminus_ptRange%d",i);
243  sprintf(bm1,"zmumu1HLT minus mass pt Range %f to %f",range0,range1);
244  cout << am1 << " " << bm1 << endl;
245  hmumu1HLTminus_pt.push_back(ptDirectory.make<TH1D>(am1,bm1,200,0.,200.));
246  sprintf(am1,"zmustaminus_ptRange%d",i);
247  sprintf(bm1,"zmusta minus mass pt Range %f to %f",range0,range1);
248  cout << am1 << " " << bm1 << endl;
249  hmustaminus_pt.push_back(ptDirectory.make<TH1D>(am1,bm1,50,0.,200.));
250  sprintf(am1,"zmutrackminus_ptRange%d",i);
251  sprintf(bm1,"zmutrack minus mass pt Range %f to %f",range0,range1);
252  cout << am1 << " " << bm1 << endl;
253  hmutrackminus_pt.push_back(ptDirectory.make<TH1D>(am1,bm1,100,0.,200.));
254  sprintf(am1,"zmumuNotIsominus_ptRange%d",i);
255  sprintf(bm1,"zmumuNotIso minus mass pt Range %f to %f",range0,range1);
256  cout << am1 << " " << bm1 << endl;
257  hmumuNotIsominus_pt.push_back(ptDirectory.make<TH1D>(am1,bm1,100,0.,200.));
258  }
259 
260  // clear global counters
262 }
263 
266  Handle<GenParticleMatch> zMuMuMatchMap; //Map of Z made by Mu global + Mu global
267  Handle<CandidateView> zMuStandAlone;
268  Handle<GenParticleMatch> zMuStandAloneMatchMap; //Map of Z made by Mu + StandAlone
269  Handle<CandidateView> zMuTrack;
270  Handle<GenParticleMatch> zMuTrackMatchMap; //Map of Z made by Mu + Track
271  Handle<CandidateView> muons; //Collection of Muons
272  Handle<CandidateView> tracks; //Collection of Tracks
273 
274  Handle<GenParticleCollection> genParticles; // Collection of Generatd Particles
275  Handle<VertexCollection> primaryVertices; // Collection of primary Vertices
276 
277  event.getByToken(zMuMuToken_, zMuMu);
278  event.getByToken(zMuStandAloneToken_, zMuStandAlone);
279  event.getByToken(zMuTrackToken_, zMuTrack);
280  event.getByToken(genParticlesToken_, genParticles);
281  event.getByToken(primaryVerticesToken_, primaryVertices);
282  event.getByToken(muonsToken_, muons);
283  event.getByToken(tracksToken_, tracks);
284 
285  /*
286  cout << "********* zMuMu size : " << zMuMu->size() << endl;
287  cout << "********* zMuStandAlone size : " << zMuStandAlone->size() << endl;
288  cout << "********* zMuTrack size : " << zMuTrack->size() << endl;
289  cout << "********* muons size : " << muons->size() << endl;
290  cout << "********* tracks size : " << tracks->size() << endl;
291  cout << "********* vertices size : " << primaryVertices->size() << endl;
292  */
293 
294  // std::cout<<"Run-> "<<event.id().run()<<std::endl;
295  // std::cout<<"Event-> "<<event.id().event()<<std::endl;
296 
297 
298 
299  bool zMuMu_found = false;
300  // loop on ZMuMu
301  if (zMuMu->size() > 0 ) {
302  for(unsigned int i = 0; i < zMuMu->size(); ++i) { //loop on candidates
303  const Candidate & zMuMuCand = (*zMuMu)[i]; //the candidate
304  CandidateBaseRef zMuMuCandRef = zMuMu->refAt(i);
305 
306  const Candidate * lep0 = zMuMuCand.daughter( 0 );
307  const Candidate * lep1 = zMuMuCand.daughter( 1 );
308  const pat::Muon & muonDau0 = dynamic_cast<const pat::Muon &>(*lep0->masterClone());
309  double trkiso0 = muonDau0.trackIso();
310  const pat::Muon & muonDau1 = dynamic_cast<const pat::Muon &>(*lep1->masterClone());
311  double trkiso1 = muonDau1.trackIso();
312 
313  // kinemtic variables
314  double pt0 = zMuMuCand.daughter(0)->pt();
315  double pt1 = zMuMuCand.daughter(1)->pt();
316  double eta0 = zMuMuCand.daughter(0)->eta();
317  double eta1 = zMuMuCand.daughter(1)->eta();
318  double charge0 = zMuMuCand.daughter(0)->charge();
319  double charge1 = zMuMuCand.daughter(1)->charge();
320  double mass = zMuMuCand.mass();
321 
322  // HLT match
323  const pat::TriggerObjectStandAloneCollection mu0HLTMatches =
324  muonDau0.triggerObjectMatchesByPath( "HLT_Mu9" );
325  const pat::TriggerObjectStandAloneCollection mu1HLTMatches =
326  muonDau1.triggerObjectMatchesByPath( "HLT_Mu9" );
327 
328  bool trig0found = false;
329  bool trig1found = false;
330  if( mu0HLTMatches.size()>0 )
331  trig0found = true;
332  if( mu1HLTMatches.size()>0 )
333  trig1found = true;
334 
335  // kinematic selection
336 
337  bool checkOppositeCharge = false;
338  if (charge0 != charge1) checkOppositeCharge = true;
339  if (pt0>ptmin_ && pt1>ptmin_ && abs(eta0)<etamax_ && abs(eta1)<etamax_ && mass>massMin_ && mass<massMax_ && checkOppositeCharge) {
340  if (trig0found || trig1found) { // at least one muon match HLT
341  zMuMu_found = true; // Z found as global-global (so don't check Zms and Zmt)
342  if (trkiso0 < isoMax_ && trkiso1 < isoMax_) { // both muons are isolated
343  if (trig0found && trig1found) {
344 
345  // ******************** category zmm 2 HLT ****************
346 
347  h_zmm2HLT_mass->Fill(mass);
348  h_zmm_mass->Fill(mass);
349 
350  // check the cynematics to fill correct histograms
351 
352  for (unsigned int j=0;j<etaBins;j++) { // eta Bins loop
353  double range0 = etaRange[j];
354  double range1= etaRange[j+1];
355 
356  // eta histograms
357 
358  if (eta0>=range0 && eta0<range1)
359  {
360  if (charge0<0) hmumu2HLTminus_eta[j]->Fill(mass); // mu- in bin eta
361  if (charge0>0) hmumu2HLTplus_eta[j]->Fill(mass); // mu+ in bin eta
362  }
363  if (eta1>=range0 && eta1<range1)
364  {
365  if (charge1<0) hmumu2HLTminus_eta[j]->Fill(mass); // mu- in bin eta
366  if (charge1>0) hmumu2HLTplus_eta[j]->Fill(mass); // mu+ in bin eta
367  }
368  } // end loop etaBins
369 
370  for (unsigned int j=0;j<ptBins;j++) { // pt Bins loop
371  double range0pt = ptRange[j];
372  double range1pt = ptRange[j+1];
373  // pt histograms
374  if (pt0>=range0pt && pt0<range1pt)
375  {
376  if (charge0<0) hmumu2HLTminus_pt[j]->Fill(mass); // mu- in bin eta
377  if (charge0>0) hmumu2HLTplus_pt[j]->Fill(mass); // mu+ in bin eta
378  }
379  if (pt1>=range0pt && pt1<range1pt)
380  {
381  if (charge1<0) hmumu2HLTminus_pt[j]->Fill(mass); // mu- in bin eta
382  if (charge1>0) hmumu2HLTplus_pt[j]->Fill(mass); // mu+ in bin eta
383  }
384  } // end loop ptBins
385 
386  } // ******************* end category zmm 2 HLT ****************
387 
388  if (!trig0found || !trig1found) {
389  // ****************** category zmm 1 HLT ******************
390  h_zmm_mass->Fill(mass);
391  double eta = 9999;
392  double pt = 9999;
393  double charge = 0;
394  if (trig0found) {
395  eta = eta1; // check muon not HLT matched
396  pt = pt1;
397  charge = charge1;
398  } else {
399  eta = eta0;
400  pt =pt0;
401  charge = charge0;
402  }
403  if (charge<0) h_zmm1HLTminus_mass->Fill(mass);
404  if (charge>0) h_zmm1HLTplus_mass->Fill(mass);
405 
406  for (unsigned int j=0;j<etaBins;j++) { // eta Bins loop
407  double range0 = etaRange[j];
408  double range1= etaRange[j+1];
409  // eta histograms fill the bin of the muon not HLT matched
410  if (eta>=range0 && eta<range1)
411  {
412  if (charge<0) hmumu1HLTminus_eta[j]->Fill(mass);
413  if (charge>0) hmumu1HLTplus_eta[j]->Fill(mass);
414  }
415  } // end loop etaBins
416  for (unsigned int j=0;j<ptBins;j++) { // pt Bins loop
417  double range0 = ptRange[j];
418  double range1= ptRange[j+1];
419  // pt histograms
420  if (pt>=range0 && pt<range1)
421  {
422  if (charge<0) hmumu1HLTminus_pt[j]->Fill(mass);
423  if (charge>0) hmumu1HLTplus_pt[j]->Fill(mass);
424  }
425  } // end loop ptBins
426 
427  } // ****************** end category zmm 1 HLT ***************
428 
429  } else { // one or both muons are not isolated
430  // ***************** category zmumuNotIso **************** (per ora non studio iso vs eta e pt da capire meglio)
431 
432  } // end if both muons isolated
433 
434  } // end if at least 1 HLT trigger found
435  } // end if kinematic selection
436 
437 
438  } // end loop on ZMuMu cand
439  } // end if ZMuMu size > 0
440 
441  // loop on ZMuSta
442  bool zMuSta_found = false;
443  if (!zMuMu_found && zMuStandAlone->size() > 0 ) {
444  event.getByToken(zMuStandAloneMatchMapToken_, zMuStandAloneMatchMap);
445  for(unsigned int i = 0; i < zMuStandAlone->size(); ++i) { //loop on candidates
446  const Candidate & zMuStandAloneCand = (*zMuStandAlone)[i]; //the candidate
447  CandidateBaseRef zMuStandAloneCandRef = zMuStandAlone->refAt(i);
448  GenParticleRef zMuStandAloneMatch = (*zMuStandAloneMatchMap)[zMuStandAloneCandRef];
449 
450  const Candidate * lep0 = zMuStandAloneCand.daughter( 0 );
451  const Candidate * lep1 = zMuStandAloneCand.daughter( 1 );
452  const pat::Muon & muonDau0 = dynamic_cast<const pat::Muon &>(*lep0->masterClone());
453  double trkiso0 = muonDau0.trackIso();
454  const pat::Muon & muonDau1 = dynamic_cast<const pat::Muon &>(*lep1->masterClone());
455  double trkiso1 = muonDau1.trackIso();
456  double pt0 = zMuStandAloneCand.daughter(0)->pt();
457  double pt1 = zMuStandAloneCand.daughter(1)->pt();
458  double eta0 = zMuStandAloneCand.daughter(0)->eta();
459  double eta1 = zMuStandAloneCand.daughter(1)->eta();
460  double charge0 = zMuStandAloneCand.daughter(0)->charge();
461  double charge1 = zMuStandAloneCand.daughter(1)->charge();
462  double mass = zMuStandAloneCand.mass();
463 
464  // HLT match
465  const pat::TriggerObjectStandAloneCollection mu0HLTMatches =
466  muonDau0.triggerObjectMatchesByPath( "HLT_Mu9" );
467  const pat::TriggerObjectStandAloneCollection mu1HLTMatches =
468  muonDau1.triggerObjectMatchesByPath( "HLT_Mu9" );
469 
470  bool trig0found = false;
471  bool trig1found = false;
472  if( mu0HLTMatches.size()>0 )
473  trig0found = true;
474  if( mu1HLTMatches.size()>0 )
475  trig1found = true;
476 
477  // check HLT match of Global muon and save eta, pt of second muon (standAlone)
478  bool trigGlbfound = false;
479  double pt =999.;
480  double eta = 999.;
481  double charge = 0;
482  if (muonDau0.isGlobalMuon()) {
483  trigGlbfound = trig0found;
484  pt = pt1;
485  eta = eta1;
486  charge = charge1;
487  }
488  if (muonDau1.isGlobalMuon()) {
489  trigGlbfound = trig1found;
490  pt = pt0;
491  eta = eta0;
492  charge = charge0;
493  }
494 
495  bool checkOppositeCharge = false;
496  if (charge0 != charge1) checkOppositeCharge = true;
497 
498  if (checkOppositeCharge && trigGlbfound && pt0>ptmin_ && pt1>ptmin_ && abs(eta0)<etamax_ && abs(eta1)<etamax_ && mass>massMin_ && mass<massMax_ && trkiso0<isoMax_ && trkiso1<isoMax_ ) { // global mu match HLT + kinematic cuts + opposite charge
499 
500  if (charge<0) h_zmsminus_mass->Fill(mass);
501  if (charge>0) h_zmsplus_mass->Fill(mass);
502 
503  for (unsigned int j=0;j<etaBins;j++) { // eta Bins loop
504  double range0 = etaRange[j];
505  double range1= etaRange[j+1];
506  // eta histograms
507  if (eta>=range0 && eta<range1) {
508  if (charge<0) hmustaminus_eta[j]->Fill(mass);
509  if (charge>0) hmustaplus_eta[j]->Fill(mass);
510  }
511  } // end loop etaBins
512  for (unsigned int j=0;j<ptBins;j++) { // pt Bins loop
513  double range0 = ptRange[j];
514  double range1= ptRange[j+1];
515  // pt histograms
516  if (pt>=range0 && pt<range1) {
517  if (charge<0) hmustaminus_pt[j]->Fill(mass);
518  if (charge>0) hmustaplus_pt[j]->Fill(mass);
519  }
520  } // end loop ptBins
521 
522  } // end if trigGlbfound + kinecuts + OppostieCharge
523  } // end loop on ZMuStandAlone cand
524  } // end if ZMuStandAlone size > 0
525 
526 
527  // loop on ZMuTrack
528  // bool zMuTrack_found = false;
529  if (!zMuMu_found && !zMuSta_found && zMuTrack->size() > 0 ) {
530  event.getByToken(zMuTrackMatchMapToken_, zMuTrackMatchMap);
531  for(unsigned int i = 0; i < zMuTrack->size(); ++i) { //loop on candidates
532  const Candidate & zMuTrackCand = (*zMuTrack)[i]; //the candidate
533  CandidateBaseRef zMuTrackCandRef = zMuTrack->refAt(i);
534  const Candidate * lep0 = zMuTrackCand.daughter( 0 );
535  const Candidate * lep1 = zMuTrackCand.daughter( 1 );
536  const pat::Muon & muonDau0 = dynamic_cast<const pat::Muon &>(*lep0->masterClone());
537  double trkiso0 = muonDau0.trackIso();
538  const pat::GenericParticle & trackDau1 = dynamic_cast<const pat::GenericParticle &>(*lep1->masterClone());
539  double trkiso1 = trackDau1.trackIso();
540  double pt0 = zMuTrackCand.daughter(0)->pt();
541  double pt1 = zMuTrackCand.daughter(1)->pt();
542  double eta0 = zMuTrackCand.daughter(0)->eta();
543  double eta1 = zMuTrackCand.daughter(1)->eta();
544  double charge0 = zMuTrackCand.daughter(0)->charge();
545  double charge1 = zMuTrackCand.daughter(1)->charge();
546  double mass = zMuTrackCand.mass();
547 
548  // HLT match (check just dau0 the global)
549  const pat::TriggerObjectStandAloneCollection mu0HLTMatches =
550  muonDau0.triggerObjectMatchesByPath( "HLT_Mu9" );
551 
552  bool trig0found = false;
553  if( mu0HLTMatches.size()>0 )
554  trig0found = true;
555 
556  bool checkOppositeCharge = false;
557  if (charge0 != charge1) checkOppositeCharge = true;
558 
559  if (checkOppositeCharge && trig0found && pt0>ptmin_ && pt1>ptmin_ && abs(eta0)<etamax_ && abs(eta1)<etamax_ && mass>massMin_ && mass<massMax_ && trkiso0<isoMax_ && trkiso1<isoMax_ ) { // global mu match HLT + kinematic cuts + opposite charge
560 
561  if (charge1<0) h_zmtminus_mass->Fill(mass);
562  if (charge1>0) h_zmtplus_mass->Fill(mass);
563 
564  for (unsigned int j=0;j<etaBins;j++) { // eta Bins loop
565  double range0 = etaRange[j];
566  double range1= etaRange[j+1];
567  // eta histograms
568  if (eta1>=range0 && eta1<range1) {
569  if (charge1<0) hmutrackminus_eta[j]->Fill(mass); // just check muon1 (mu0 is global by definition)
570  if (charge1>0) hmutrackplus_eta[j]->Fill(mass); // just check muon1 (mu0 is global by definition)
571  }
572  } // end loop etaBins
573  for (unsigned int j=0;j<ptBins;j++) { // pt Bins loop
574  double range0 = ptRange[j];
575  double range1= ptRange[j+1];
576  // pt histograms
577  if (pt1>=range0 && pt1<range1) {
578  if (charge1<0) hmutrackminus_pt[j]->Fill(mass); // just check muon1 (mu0 is global by definition)
579  if (charge1>0) hmutrackplus_pt[j]->Fill(mass); // just check muon1 (mu0 is global by definition)
580  }
581  } // end loop ptBins
582 
583  } // end if trig0found
584 
585 
586  } // end loop on ZMuTrack cand
587  } // end if ZMuTrack size > 0
588 
589 } // end analyze
590 
591 bool ZMuMu_efficiencyAnalyzer::check_ifZmumu(const Candidate * dauGen0, const Candidate * dauGen1, const Candidate * dauGen2)
592 {
593  int partId0 = dauGen0->pdgId();
594  int partId1 = dauGen1->pdgId();
595  int partId2 = dauGen2->pdgId();
596  bool muplusFound=false;
597  bool muminusFound=false;
598  bool ZFound=false;
599  if (partId0==13 || partId1==13 || partId2==13) muminusFound=true;
600  if (partId0==-13 || partId1==-13 || partId2==-13) muplusFound=true;
601  if (partId0==23 || partId1==23 || partId2==23) ZFound=true;
602  return muplusFound*muminusFound*ZFound;
603 }
604 
605 float ZMuMu_efficiencyAnalyzer::getParticlePt(const int ipart, const Candidate * dauGen0, const Candidate * dauGen1, const Candidate * dauGen2)
606 {
607  int partId0 = dauGen0->pdgId();
608  int partId1 = dauGen1->pdgId();
609  int partId2 = dauGen2->pdgId();
610  float ptpart=0.;
611  if (partId0 == ipart) {
612  for(unsigned int k = 0; k < dauGen0->numberOfDaughters(); ++k) {
613  const Candidate * dauMuGen = dauGen0->daughter(k);
614  if(dauMuGen->pdgId() == ipart && dauMuGen->status() ==1) {
615  ptpart = dauMuGen->pt();
616  }
617  }
618  }
619  if (partId1 == ipart) {
620  for(unsigned int k = 0; k < dauGen1->numberOfDaughters(); ++k) {
621  const Candidate * dauMuGen = dauGen1->daughter(k);
622  if(dauMuGen->pdgId() == ipart && dauMuGen->status() ==1) {
623  ptpart = dauMuGen->pt();
624  }
625  }
626  }
627  if (partId2 == ipart) {
628  for(unsigned int k = 0; k < dauGen2->numberOfDaughters(); ++k) {
629  const Candidate * dauMuGen = dauGen2->daughter(k);
630  if(abs(dauMuGen->pdgId()) == ipart && dauMuGen->status() ==1) {
631  ptpart = dauMuGen->pt();
632  }
633  }
634  }
635  return ptpart;
636 }
637 
638 float ZMuMu_efficiencyAnalyzer::getParticleEta(const int ipart, const Candidate * dauGen0, const Candidate * dauGen1, const Candidate * dauGen2)
639 {
640  int partId0 = dauGen0->pdgId();
641  int partId1 = dauGen1->pdgId();
642  int partId2 = dauGen2->pdgId();
643  float etapart=0.;
644  if (partId0 == ipart) {
645  for(unsigned int k = 0; k < dauGen0->numberOfDaughters(); ++k) {
646  const Candidate * dauMuGen = dauGen0->daughter(k);
647  if(dauMuGen->pdgId() == ipart && dauMuGen->status() ==1) {
648  etapart = dauMuGen->eta();
649  }
650  }
651  }
652  if (partId1 == ipart) {
653  for(unsigned int k = 0; k < dauGen1->numberOfDaughters(); ++k) {
654  const Candidate * dauMuGen = dauGen1->daughter(k);
655  if(dauMuGen->pdgId() == ipart && dauMuGen->status() ==1) {
656  etapart = dauMuGen->eta();
657  }
658  }
659  }
660  if (partId2 == ipart) {
661  for(unsigned int k = 0; k < dauGen2->numberOfDaughters(); ++k) {
662  const Candidate * dauMuGen = dauGen2->daughter(k);
663  if(abs(dauMuGen->pdgId()) == ipart && dauMuGen->status() ==1) {
664  etapart = dauMuGen->eta();
665  }
666  }
667  }
668  return etapart;
669 }
670 
671 float ZMuMu_efficiencyAnalyzer::getParticlePhi(const int ipart, const Candidate * dauGen0, const Candidate * dauGen1, const Candidate * dauGen2)
672 {
673  int partId0 = dauGen0->pdgId();
674  int partId1 = dauGen1->pdgId();
675  int partId2 = dauGen2->pdgId();
676  float phipart=0.;
677  if (partId0 == ipart) {
678  for(unsigned int k = 0; k < dauGen0->numberOfDaughters(); ++k) {
679  const Candidate * dauMuGen = dauGen0->daughter(k);
680  if(dauMuGen->pdgId() == ipart && dauMuGen->status() ==1) {
681  phipart = dauMuGen->phi();
682  }
683  }
684  }
685  if (partId1 == ipart) {
686  for(unsigned int k = 0; k < dauGen1->numberOfDaughters(); ++k) {
687  const Candidate * dauMuGen = dauGen1->daughter(k);
688  if(dauMuGen->pdgId() == ipart && dauMuGen->status() ==1) {
689  phipart = dauMuGen->phi();
690  }
691  }
692  }
693  if (partId2 == ipart) {
694  for(unsigned int k = 0; k < dauGen2->numberOfDaughters(); ++k) {
695  const Candidate * dauMuGen = dauGen2->daughter(k);
696  if(abs(dauMuGen->pdgId()) == ipart && dauMuGen->status() ==1) {
697  phipart = dauMuGen->phi();
698  }
699  }
700  }
701  return phipart;
702 }
703 
704 Particle::LorentzVector ZMuMu_efficiencyAnalyzer::getParticleP4(const int ipart, const Candidate * dauGen0, const Candidate * dauGen1, const Candidate * dauGen2)
705 {
706  int partId0 = dauGen0->pdgId();
707  int partId1 = dauGen1->pdgId();
708  int partId2 = dauGen2->pdgId();
709  Particle::LorentzVector p4part(0.,0.,0.,0.);
710  if (partId0 == ipart) {
711  for(unsigned int k = 0; k < dauGen0->numberOfDaughters(); ++k) {
712  const Candidate * dauMuGen = dauGen0->daughter(k);
713  if(dauMuGen->pdgId() == ipart && dauMuGen->status() ==1) {
714  p4part = dauMuGen->p4();
715  }
716  }
717  }
718  if (partId1 == ipart) {
719  for(unsigned int k = 0; k < dauGen1->numberOfDaughters(); ++k) {
720  const Candidate * dauMuGen = dauGen1->daughter(k);
721  if(dauMuGen->pdgId() == ipart && dauMuGen->status() ==1) {
722  p4part = dauMuGen->p4();
723  }
724  }
725  }
726  if (partId2 == ipart) {
727  for(unsigned int k = 0; k < dauGen2->numberOfDaughters(); ++k) {
728  const Candidate * dauMuGen = dauGen2->daughter(k);
729  if(abs(dauMuGen->pdgId()) == ipart && dauMuGen->status() ==1) {
730  p4part = dauMuGen->p4();
731  }
732  }
733  }
734  return p4part;
735 }
736 
737 
738 
740 
741 
742 
743 }
744 
746 
748 
std::vector< GenParticle > GenParticleCollection
collection of GenParticles
int i
Definition: DBlmapReader.cc:9
virtual const Candidate * daughter(size_type i) const =0
return daughter at a given position, i = 0, ... numberOfDaughters() - 1 (read only mode) ...
EDGetTokenT< GenParticleMatch > zMuMuMatchMapToken_
reco::CandidateBaseRef trackMuonCandRef_
EDGetTokenT< CandidateView > muonsToken_
virtual double pt() const =0
transverse momentum
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:17
virtual double mass() const =0
mass
EDGetTokenT< GenParticleMatch > zMuStandAloneMatchMapToken_
virtual int status() const =0
status word
EDGetTokenT< GenParticleCollection > genParticlesToken_
bool isGlobalMuon() const
Definition: Muon.h:222
std::vector< Vertex > VertexCollection
collection of Vertex objects
Definition: VertexFwd.h:9
T * make(const Args &...args) const
make new ROOT object
Definition: TFileService.h:64
std::vector< TriggerObjectStandAlone > TriggerObjectStandAloneCollection
Collection of TriggerObjectStandAlone.
EDGetTokenT< CandidateView > tracksToken_
virtual void analyze(const edm::Event &event, const edm::EventSetup &setup) override
float trackIso() const
float trackIso() const
Definition: Muon.h:170
EDGetTokenT< GenParticleMatch > zMuTrackMatchMapToken_
bool check_ifZmumu(const Candidate *dauGen0, const Candidate *dauGen1, const Candidate *dauGen2)
EDGetTokenT< CandidateView > zMuTrackToken_
float getParticlePhi(const int ipart, const Candidate *dauGen0, const Candidate *dauGen1, const Candidate *dauGen2)
Analysis-level Generic Particle class (e.g. for hadron or muon not fully reconstructed) ...
virtual size_type numberOfDaughters() const =0
number of daughters
edm::ValueMap< float > IsolationCollection
tuple lep1
print &#39;MRbb(1b)&#39;,event.mr_bb
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
int j
Definition: DBlmapReader.cc:9
virtual int charge() const =0
electric charge
T * make(const Args &...args) const
make new ROOT object
const TriggerObjectStandAloneCollection triggerObjectMatchesByPath(const std::string &namePath, const bool pathLastFilterAccepted=false, const bool pathL3FilterAccepted=true) const
Definition: PATObject.h:612
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< CandidateView > zMuStandAloneToken_
virtual int pdgId() const =0
PDG identifier.
float getParticleEta(const int ipart, const Candidate *dauGen0, const Candidate *dauGen1, const Candidate *dauGen2)
TFileDirectory mkdir(const std::string &dir, const std::string &descr="")
create a new subdirectory
Definition: TFileService.h:69
tuple tracks
Definition: testEve_cfg.py:39
Particle::LorentzVector getParticleP4(const int ipart, const Candidate *dauGen0, const Candidate *dauGen1, const Candidate *dauGen2)
EDGetTokenT< VertexCollection > primaryVerticesToken_
virtual void endJob() override
tuple muons
Definition: patZpeak.py:38
tuple cout
Definition: gather_cfg.py:145
EDGetTokenT< CandidateView > zMuMuToken_
float getParticlePt(const int ipart, const Candidate *dauGen0, const Candidate *dauGen1, const Candidate *dauGen2)
Analysis-level muon class.
Definition: Muon.h:49
ZMuMu_efficiencyAnalyzer(const edm::ParameterSet &pset)
tuple zMuMu
zMuMu vector of PSet is common to all categories except zMuTrk category
math::PtEtaPhiELorentzVectorF LorentzVector
virtual double phi() const =0
momentum azimuthal angle
virtual double eta() const =0
momentum pseudorapidity
virtual const LorentzVector & p4() const =0
four-momentum Lorentz vector
virtual const CandidateBaseRef & masterClone() const =0