ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitHzz4l/LeptonSelection/src/MuonSelection.cc
Revision: 1.35
Committed: Tue Jul 17 10:17:22 2012 UTC (12 years, 10 months ago) by dkralph
Content type: text/plain
Branch: MAIN
Changes since 1.34: +4 -2 lines
Log Message:
1) implemented a way to store the bitset of all matched triggers for the leptons in the FO branches (using TBits, since root doesn't seem to yet support bitsets). See some commented lines in applyZPlusX.cc and ZPlusX.cc (more or less moot -- see below)

2) These bitsets are almost entirely zeroes, but root didn't compress them as much as I was hoping -- adding them tripled the file size, so I had to implement another way to store trigger info in the ntuples. I decided to store info for only a few  triggers in an unsigned value in SimpleLepton, with the def of the bits at the bottom of TriggerUtils.cc in the function initAnalysisTriggers (http://cmssw.cvs.cern.ch/cgi-bin/cmssw.cgi/UserCode/MitHzz4l/Util/src/TriggerUtils.cc?r1=1.15&r2=1.16)

3) added the same (info for a few triggers) for HLT Info to an unsigned value in InfoStruct

4) added all versions of the emu cross triggers to TriggerUtils.cc

File Contents

# Content
1 #include <math.h>
2 #include <iostream>
3
4 #include "MuonSelection.h"
5 #include "ParseArgs.h"
6 #include "SelectionStatus.h"
7
8 #include "Track.h"
9 #include "MuonTools.h"
10 #include "MuonQuality.h"
11 #include "MuonIDMVA.h"
12
13 mithep::MuonIDMVA * muIDMVA;
14 mithep::MuonTools muTools;
15
16 extern vector<bool> PFnoPUflag;
17 extern Float_t computePFMuonIso(const mithep::Muon *mu,
18 const mithep::Vertex * vtx,
19 const mithep::Array<mithep::PFCandidate> * pfCandidates,
20 const Double_t dRMax);
21
22 bool muonPOG2012CutBasedIDTight(const mithep::Muon *mu, const mithep::Vertex * vtx, const mithep::Array<mithep::PFCandidate> * pfCandidates, const mithep::Array<mithep::PileupEnergyDensity> * puEnergyDensity)
23 {
24
25 Float_t charged_iso = 0;
26
27 Float_t sumChargedHadronPt = 0;
28 Float_t sumNeutralHadronPt = 0;
29 Float_t sumPhotonPt = 0;
30 Float_t sumPUPt = 0;
31
32 for(unsigned k=0; k<pfCandidates->GetEntries(); ++k) {
33 const mithep::PFCandidate *pf = (mithep::PFCandidate*)((*pfCandidates)[k]);
34 Double_t dr = mithep::MathUtils::DeltaR(mu->Phi(),mu->Eta(), pf->Phi(), pf->Eta());
35 if(dr > 0.4) continue;
36 if(pf->TrackerTrk() && mu->TrackerTrk() && pf->TrackerTrk() == mu->TrackerTrk()) continue;
37
38 if(PFnoPUflag[k]){ //pf no pileup
39 if(abs(pf->PFType()) == mithep::PFCandidate::eHadron)
40 sumChargedHadronPt += pf->Pt();
41 else if (pf->Pt() > 0.5 && abs(pf->PFType()) == mithep::PFCandidate::eNeutralHadron)
42 sumNeutralHadronPt += pf->Pt();
43 else if (pf->Pt() > 0.5 && abs(pf->PFType()) == mithep::PFCandidate::eGamma)
44 sumPhotonPt += pf->Pt();
45 }
46 else {//pf pile up
47 if(pf->Charge() != 0)
48 sumPUPt += pf->Pt();
49 }
50
51 }
52 double iso_sum = sumChargedHadronPt+ max(0.,sumNeutralHadronPt+sumPhotonPt-0.5*sumPUPt);
53
54 bool pass = kTRUE;
55 if(!mu->IsGlobalMuon()) pass = kFALSE;
56
57 //there is not necessarily a global track or a tracker track if is global muon is false
58 if(!pass) return pass;
59
60 if(!mu->IsPFMuon()) pass = kFALSE;
61 if(mu->GlobalTrk()->RChi2() > 10) pass= kFALSE;
62 if(mu->NValidHits() == 0 ) pass= kFALSE;
63 if(mu->NMatches() <= 1) pass = kFALSE;
64 if(fabs(mu->TrackerTrk()->D0Corrected(*vtx)) > 0.2) pass = kFALSE;
65 if(fabs(mu->TrackerTrk()->DzCorrected(*vtx)) > 0.5) pass = kFALSE;
66 if(mu->BestTrk()->NPixelHits() == 0) pass = kFALSE;
67 if(mu->NTrkLayersHit() <= 5) pass = kFALSE;
68 if(iso_sum/mu->Pt() > 0.12 ) pass = kFALSE;
69
70 return pass;
71
72
73 }
74 //----------------------------------------------------------------------------------------
75 bool muon2012CutBasedIDTightVersionWithOldIsoThatWorksOn2011(ControlFlags &ctrl, const mithep::Muon *mu, const mithep::Vertex * vtx, const mithep::Array<mithep::PFCandidate> * pfCandidates,
76 const mithep::Array<mithep::PileupEnergyDensity> * puEnergyDensity, mithep::MuonTools::EMuonEffectiveAreaTarget EffectiveAreaVersion
77 )
78 {
79
80 mithep::MuonTools muT;
81
82 Float_t charged_iso = 0;
83 Float_t neutral_iso = 0;
84 Float_t gamma_iso = 0;
85
86 for(int k=0; k<pfCandidates->GetEntries(); ++k) {
87
88 if( !(PFnoPUflag[k]) ) continue;
89 const mithep::PFCandidate *pf = (mithep::PFCandidate*)((*pfCandidates)[k]);
90
91 Double_t dr = mithep::MathUtils::DeltaR(mu->Phi(),mu->Eta(), pf->Phi(), pf->Eta());
92 if (dr > 0.4) continue;
93
94 if (pf->Charge() != 0 && (pf->HasTrackerTrk()||pf->HasGsfTrk()) ) {
95
96 if (abs(pf->PFType()) == mithep::PFCandidate::eElectron || abs(pf->PFType()) == mithep::PFCandidate::eMuon) continue;
97
98 charged_iso += pf->Pt();
99 }
100
101 else if (abs(pf->PFType()) == mithep::PFCandidate::eGamma) {
102
103 gamma_iso += pf->Pt();
104 }
105
106 else {
107 neutral_iso += pf->Pt();
108 }
109 }
110
111 double iso_sum=999;
112 double rho;
113 if(ctrl.era == 2012) {
114 EffectiveAreaVersion = mithep::MuonTools::kMuEAData2012;
115 rho = puEnergyDensity->At(0)->RhoKt6PFJetsCentralNeutral();
116
117 iso_sum = charged_iso + fmax(0.0,(gamma_iso + neutral_iso
118 -rho*muT.MuonEffectiveArea(muT.kMuGammaAndNeutralHadronIso04,
119 mu->Eta(),EffectiveAreaVersion)));
120 } else if(ctrl.era == 2011) {
121 setMuonRhoAndEaVersion(rho,EffectiveAreaVersion,puEnergyDensity);
122 iso_sum = charged_iso + fmax(0.0,(gamma_iso + neutral_iso
123 -rho*muT.MuonEffectiveArea(muT.kMuGammaAndNeutralHadronIso04,
124 mu->Eta(),EffectiveAreaVersion)));
125 } else assert(0);
126
127 bool pass = kTRUE;
128 if(!mu->IsGlobalMuon()) pass = kFALSE;
129
130 //there is not necessarily a global track or a tracker track if is global muon is false
131 if(!pass) return pass;
132
133 if(!mu->IsPFMuon() && ctrl.era==2012) pass = kFALSE;
134 if(mu->GlobalTrk()->RChi2() > 10) pass = kFALSE;
135 if(mu->NValidHits() == 0 ) pass = kFALSE;
136 if(mu->NMatches() <= 1) pass = kFALSE;
137 if(fabs(mu->TrackerTrk()->D0Corrected(*vtx)) > 0.2) pass = kFALSE;
138 if(fabs(mu->TrackerTrk()->DzCorrected(*vtx)) > 0.5) pass = kFALSE;
139 if(mu->BestTrk()->NPixelHits() == 0) pass = kFALSE;
140 if(mu->NTrkLayersHit() <= 5) pass = kFALSE;
141 if(iso_sum/mu->Pt() > 0.12 ) pass = kFALSE;
142
143 return pass;
144
145
146 }
147
148 bool trackToPFMuonSelection( const mithep::Track *track, const mithep::MuonCol * muoncol )
149 {
150 bool pass = false;
151 for(UInt_t i=0; i<muoncol->GetEntries(); ++i) {
152 const mithep::Muon *mu = muoncol->At(i);
153 if(mu->TrackerTrk() == track && (mu->IsGlobalMuon() || mu->IsTrackerMuon()) && mu->IsPFMuon())
154 pass = true;
155 }
156
157 return pass;
158
159 }
160
161 //--------------------------------------------------------------------------------------------------
162 SelectionStatus muonDummyVeto(ControlFlags &ctrl,
163 const mithep::Muon *muon,
164 const mithep::Vertex * vtx)
165 //--------------------------------------------------------------------------------------------------
166 {
167 SelectionStatus status;
168 status.setStatus(SelectionStatus::PRESELECTION);
169 return status;
170 }
171
172 //--------------------------------------------------------------------------------------------------
173 SelectionStatus muonCutBasedVeto(ControlFlags &ctrl,
174 const mithep::Muon *muon,
175 const mithep::Vertex *vtx)
176 //--------------------------------------------------------------------------------------------------
177 {
178 //
179 // Loose cut-based ID for isolation veto
180 //
181 bool ret = true;
182
183 if(!(muon->IsGlobalMuon() || muon->IsTrackerMuon())) ret = false;
184 if( muon->NValidHits() < 11 ) ret = false;
185 if( fabs(muon->Ip3dPVSignificance()) >= 4 ) ret = false;
186
187 SelectionStatus status;
188 if( ret ) status.setStatus(SelectionStatus::PRESELECTION);
189 return status;
190 }
191
192
193 //--------------------------------------------------------------------------------------------------
194 SelectionStatus noPreSelection( ControlFlags &ctrl, const mithep::Muon * mu )
195
196 //--------------------------------------------------------------------------------------------------
197 {
198 SelectionStatus status;
199 status.setStatus(SelectionStatus::PRESELECTION);
200 if(ctrl.debug) cout << "muon presel returning status : " << status.getStatus() << endl;
201 return status;
202 }
203
204 //--------------------------------------------------------------------------------------------------
205 SelectionStatus muonPreSelection( ControlFlags &ctrl,
206 const mithep::Muon * mu,
207 const mithep::Vertex * vtx,
208 const mithep::Array<mithep::PFCandidate> * pfCandidates )
209 //--------------------------------------------------------------------------------------------------
210 {
211 bool ret = true;
212
213 ret &= ( fabs(mu->Ip3dPVSignificance()) < 100 );
214 ret &= ( mu->Pt() >= 5 );
215 ret &= ( fabs(mu->Eta()) <= 2.4 );
216 ret &= ( mu->IsTrackerMuon() || mu->IsGlobalMuon() );
217 // if( mu->IsTrackerMuon() && mu->HasTrackerTrk() && !mu->IsGlobalMuon() )
218 // ret &= ( mu->Quality().QualityMask().Mask() & mithep::MuonQuality::AllArbitrated);
219 ret &= ( mu->IsoR03SumPt()/mu->Pt() < 0.7 );
220
221 SelectionStatus status;
222 if( ret ) status.setStatus(SelectionStatus::PRESELECTION);
223 if(ctrl.debug) cout << "muon presel returning status : " << status.getStatus() << endl;
224 return status;
225 }
226
227
228
229 //--------------------------------------------------------------------------------------------------
230 SelectionStatus muonReferencePreSelection( ControlFlags &ctrl,
231 const mithep::Muon * mu,
232 const mithep::Vertex * vtx,
233 const mithep::Array<mithep::PFCandidate> * pfCandidates )
234 //--------------------------------------------------------------------------------------------------
235 {
236 bool ret = true;
237
238 if(ctrl.debug) cout << "inside muonpresel ... " << endl;
239 ret &= ( fabs(mu->Ip3dPVSignificance()) < 100 );
240 if(ctrl.debug) cout << "and pass IP (" << mu->Ip3dPVSignificance() << ") ? ... " << ret << endl;
241 ret &= ( mu->Pt() >= 5 );
242 if(ctrl.debug) cout << "and >5 GeV (" << mu->Pt() << ") ? ... " << ret << endl;
243 ret &= ( fabs(mu->Eta()) <= 2.4 );
244 if(ctrl.debug) cout << "and < 2.4 eta (" << mu->Eta() << ")? ... " << ret << endl;
245 ret &= ( mu->IsTrackerMuon() || mu->IsGlobalMuon() );
246 if(ctrl.debug) cout << "is Tracker or Global? ... " << ret << endl;
247 ret &= (mu->HasTrackerTrk() &&
248 fabs(mu->TrackerTrk()->D0Corrected(*vtx)) < 0.5 &&
249 fabs(mu->TrackerTrk()->DzCorrected(*vtx)) < 1.0);
250 if(ctrl.debug) cout << "d0 and dz? ... " << ret << endl;
251
252 // if( mu->IsTrackerMuon() && mu->HasTrackerTrk() && !mu->IsGlobalMuon() )
253 // ret &= ( mu->Quality().QualityMask().Mask() & mithep::MuonQuality::AllArbitrated);
254 // if(ctrl.debug) cout << "if isTrackerMuon, arbitrated ? ... " << ret << endl;
255 // ret &= ( mu->BestTrk()->DzCorrected(vtx) < 0.1);
256 // if( ctrl.debug ) cout << "elepresel : ret after dZcorr : " << ret << endl;
257
258
259 SelectionStatus status;
260 if( ret ) status.setStatus(SelectionStatus::PRESELECTION);
261 if(ctrl.debug) cout << "muon presel returning status : " << status.getStatus() << endl;
262 return status;
263 }
264 //--------------------------------------------------------------------------------------------------
265 SelectionStatus muonPreSelectionNoD0DzIP( ControlFlags &ctrl,
266 const mithep::Muon * mu,
267 const mithep::Vertex * vtx,
268 const mithep::Array<mithep::PFCandidate> * pfCandidates )
269 //--------------------------------------------------------------------------------------------------
270 {
271 bool ret = true;
272
273 if(ctrl.debug) cout << "inside muonpresel ... " << endl;
274 ret &= ( mu->Pt() >= 5 );
275 if(ctrl.debug) cout << "and >5 GeV (" << mu->Pt() << ") ? ... " << ret << endl;
276 ret &= ( fabs(mu->Eta()) <= 2.4 );
277 if(ctrl.debug) cout << "and < 2.4 eta (" << mu->Eta() << ")? ... " << ret << endl;
278 ret &= ( mu->IsTrackerMuon() || mu->IsGlobalMuon() );
279 if(ctrl.debug) cout << "is Tracker or Global? ... " << ret << endl;
280
281 // if( mu->IsTrackerMuon() && mu->HasTrackerTrk() && !mu->IsGlobalMuon() )
282 // ret &= ( mu->Quality().QualityMask().Mask() & mithep::MuonQuality::AllArbitrated);
283 // if(ctrl.debug) cout << "if isTrackerMuon, arbitrated ? ... " << ret << endl;
284 // ret &= ( mu->BestTrk()->DzCorrected(vtx) < 0.1);
285 // if( ctrl.debug ) cout << "elepresel : ret after dZcorr : " << ret << endl;
286
287
288 SelectionStatus status;
289 if( ret ) status.setStatus(SelectionStatus::PRESELECTION);
290 if(ctrl.debug) cout << "muon presel returning status : " << status.getStatus() << endl;
291 return status;
292 }
293 //--------------------------------------------------------------------------------------------------
294 SelectionStatus muonPreSelectionNoD0IP( ControlFlags &ctrl,
295 const mithep::Muon * mu,
296 const mithep::Vertex * vtx,
297 const mithep::Array<mithep::PFCandidate> * pfCandidates )
298 //--------------------------------------------------------------------------------------------------
299 {
300 bool ret = true;
301
302 ret &= ( mu->Pt() >= 5 );
303 ret &= ( fabs(mu->Eta()) <= 2.4 );
304 ret &= ( mu->IsTrackerMuon() || mu->IsGlobalMuon() );
305 // ret &= (mu->HasTrackerTrk() &&
306 // fabs(mu->TrackerTrk()->D0Corrected(*vtx)) < 0.5 &&
307 // fabs(mu->TrackerTrk()->DzCorrected(*vtx)) < 1.0);
308 ret &= (mu->HasTrackerTrk() &&
309 fabs(mu->TrackerTrk()->DzCorrected(*vtx)) < 1.0);
310
311 SelectionStatus status;
312 if( ret ) status.setStatus(SelectionStatus::PRESELECTION);
313 if(ctrl.debug) cout << "muon presel returning status : " << status.getStatus() << endl;
314 return status;
315 }
316 //--------------------------------------------------------------------------------------------------
317 SelectionStatus muonIDMVASelection(ControlFlags &ctrl,
318 const mithep::Muon *mu,
319 const mithep::Vertex * vtx )
320 //--------------------------------------------------------------------------------------------------
321 {
322 double global_rchi2 = (mu->IsGlobalMuon()) ? mu->GlobalTrk()->RChi2() : 0.;
323
324 //
325 // WARNING !!! KH hacked out to sync
326 //
327 /*
328 double mvaval = muIDMVA->MVAValue_ID(mu->Pt(),
329 mu->Eta(),
330 mu->IsGlobalMuon(),
331 mu->IsTrackerMuon(),
332 mu->TrackerTrk()->RChi2(),
333 global_rchi2,
334 (Double_t)(mu->NValidHits()),
335 (Double_t)(mu->TrackerTrk()->NHits()),
336 (Double_t)(mu->TrackerTrk()->NPixelHits()),
337 (Double_t)(mu->NMatches()),
338 mu->TrkKink(),
339 muTools.GetSegmentCompatability(mu),
340 muTools.GetCaloCompatability(mu,true,true),
341 mu->HadEnergy(),
342 mu->EmEnergy(),
343 mu->HadS9Energy(),
344 mu->EmS9Energy(),
345 (Bool_t)(ctrl.debug) );
346 */
347 double mvaval = 0.0;
348
349 SelectionStatus status;
350 bool pass = false;
351
352 if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
353 && fabs(mu->Eta()) <= 1.5 && mu->Pt() <= 10 && mvaval >= MUON_IDMVA_CUT_BIN0) pass = true;
354 else if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
355 && fabs(mu->Eta()) <= 1.5 && mu->Pt() > 10 && mvaval >= MUON_IDMVA_CUT_BIN1) pass = true;
356 else if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
357 && fabs(mu->Eta()) > 1.5 && mu->Pt() <= 10 && mvaval >= MUON_IDMVA_CUT_BIN2) pass = true;
358 else if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
359 && fabs(mu->Eta()) > 1.5 && mu->Pt() > 10 && mvaval >= MUON_IDMVA_CUT_BIN3) pass = true;
360 else if( !(mu->IsGlobalMuon()) && mu->IsTrackerMuon() && mvaval >= MUON_IDMVA_CUT_BIN4 ) pass = true;
361 else if( mu->IsGlobalMuon() && !(mu->IsTrackerMuon()) && mvaval >= MUON_IDMVA_CUT_BIN5 ) pass = true;
362
363
364
365 if( pass ) {
366 status.orStatus(SelectionStatus::LOOSEID);
367 status.orStatus(SelectionStatus::TIGHTID);
368 }
369
370 return status;
371 }
372
373 //--------------------------------------------------------------------------------------------------
374 void initMuonIDMVA()
375 //--------------------------------------------------------------------------------------------------
376 {
377 muIDMVA = new mithep::MuonIDMVA();
378 vector<string> weightFiles;
379 weightFiles.push_back("./data/MuonIDMVAWeights/MuonIDMVA_BDTG_barrel_lowpt_V2.weights.xml");
380 weightFiles.push_back("./data/MuonIDMVAWeights/MuonIDMVA_BDTG_barrel_highpt_V2.weights.xml");
381 weightFiles.push_back("./data/MuonIDMVAWeights/MuonIDMVA_BDTG_endcap_lowpt_V2.weights.xml");
382 weightFiles.push_back("./data/MuonIDMVAWeights/MuonIDMVA_BDTG_endcap_highpt_V2.weights.xml");
383 weightFiles.push_back("./data/MuonIDMVAWeights/MuonIDMVA_BDTG_tracker_V2.weights.xml");
384 weightFiles.push_back("./data/MuonIDMVAWeights/MuonIDMVA_BDTG_global_V2.weights.xml");
385 muIDMVA->Initialize( "MuonIDMVA",
386 mithep::MuonIDMVA::kIDV0,
387 kTRUE, weightFiles);
388 }
389
390 //--------------------------------------------------------------------------------------------------
391 SelectionStatus muonIDPFSelection(ControlFlags &ctrl,
392 const mithep::Muon *mu,
393 const mithep::Vertex * vtx,
394 const mithep::Array<mithep::PFCandidate> * pfCandidates )
395 //--------------------------------------------------------------------------------------------------
396 {
397 bool pass = false;
398 SelectionStatus status; // init'ed to FAIL
399
400 // check that it matches to a PF muon
401 for( int i=0; i<pfCandidates->GetEntries(); i++ ) {
402 // tmp
403 if( !(PFnoPUflag[i]) ) continue; // my PF no PU hack
404
405 const mithep::PFCandidate *pf = (mithep::PFCandidate*)((*pfCandidates)[i]);
406 if( (pf->TrackerTrk() == mu->TrackerTrk()) && abs(pf->PFType()) == mithep::PFCandidate::eMuon ) {
407 pass = true;
408 break;
409 }
410 }
411 //
412 // if(mu->IsPFMuon()) pass = true;
413 // else pass = false;
414
415 if( pass ) {
416 status.orStatus(SelectionStatus::LOOSEID);
417 status.orStatus(SelectionStatus::TIGHTID);
418 }
419
420 return status;
421 }
422 //----------------------------------------------------------------------------------------
423 SelectionStatus PassWwMuonSel(const mithep::Muon *mu, const mithep::Vertex *vtx,
424 const mithep::Array<mithep::PFCandidate> *pfArr
425 )
426 {
427 bool passID=true;
428 bool passIso=true;
429
430 if (mu->Pt() < 5)
431 passID = false;
432 if (fabs(mu->Eta()) > 2.4)
433 passID = false;
434 if (mu->BestTrk()->PtErr()/mu->Pt() > 0.1)
435 passID = false;
436 if (fabs(mu->BestTrk()->DzCorrected(*vtx)) > 0.1)
437 passID = false;
438 Bool_t isGlobal = (mu->IsGlobalMuon()) && (mu->BestTrk()->RChi2() < 10) &&
439 (mu->NMatches() > 1) && (mu->NValidHits() > 0);
440 Bool_t isTracker = (mu->IsTrackerMuon()) &&
441 (mu->Quality().Quality(mithep::MuonQuality::TMLastStationTight));
442 if (!isGlobal && !isTracker)
443 passID = false;
444 int ntrkhits = (mu->HasTrackerTrk()) ? mu->TrackerTrk()->NHits() : 0;
445 if (ntrkhits < 10)
446 passID = false;
447 if (mu->BestTrk()->NPixelHits() < 1)
448 passID = false;
449 if (fabs(mu->BestTrk()->D0Corrected(*vtx)) > 0.02)
450 passID = false;
451
452 // note: this isn't really ww muon isolation
453 float pfiso = computePFMuonIso(mu,vtx,pfArr,0.3);
454
455 if (pfiso > 0.2*mu->Pt())
456 passIso = false;
457
458 SelectionStatus status;
459 if(passID) status |= SelectionStatus::TIGHTIDANDPRE;
460 if(passIso) status |= SelectionStatus::TIGHTISOANDPRE;
461
462 return status;
463
464 }