ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitHzz4l/LeptonSelection/src/IsolationSelection.cc
Revision: 1.28
Committed: Mon May 28 18:27:38 2012 UTC (12 years, 11 months ago) by khahn
Content type: text/plain
Branch: MAIN
CVS Tags: synched2
Changes since 1.27: +2 -2 lines
Log Message: 
2012 data/mc should use 2012 EA's

File Contents

# Content
1 #include <math.h>
2
3 #include "IsolationSelection.h"
4 #include "IsolationSelectionDefs.h"
5
6 #include "MathUtils.h"
7 #include "MuonTools.h"
8 #include "MuonIDMVA.h"
9 #include "ElectronTools.h"
10 #include "ElectronIDMVA.h"
11
12 using namespace mithep;
13
14 mithep::MuonIDMVA * muIsoMVA;
15 mithep::MuonTools muT;
16 mithep::ElectronIDMVA * eleIsoMVA;
17 mithep::ElectronTools eleT;
18
19 // global hack to sync
20 double gChargedIso;
21 double gGammaIso;
22 double gNeutralIso;
23
24 extern vector<bool> PFnoPUflag;
25
26 //--------------------------------------------------------------------------------------------------
27 Float_t computePFMuonIso(const mithep::Muon *muon,
28 const mithep::Vertex * vtx,
29 const mithep::Array<mithep::PFCandidate> * fPFCandidates,
30 const Double_t dRMax)
31 //--------------------------------------------------------------------------------------------------
32 {
33 const Double_t dRMin = 0;
34 const Double_t neuPtMin = 1.0;
35 const Double_t dzMax = 0.1;
36
37 Double_t zLepton = (muon->BestTrk()) ? muon->BestTrk()->DzCorrected(*vtx) : 0.0;
38
39 Float_t iso=0;
40 for(UInt_t ipf=0; ipf<fPFCandidates->GetEntries(); ipf++) {
41 const PFCandidate *pfcand = fPFCandidates->At(ipf);
42
43 if(!pfcand->HasTrk() && (pfcand->Pt()<=neuPtMin)) continue; // pT cut on neutral particles
44
45 // exclude THE muon
46 if(pfcand->TrackerTrk() && muon->TrackerTrk() && (pfcand->TrackerTrk()==muon->TrackerTrk())) continue;
47
48 // dz cut
49 Double_t dz = (pfcand->BestTrk()) ? fabs(pfcand->BestTrk()->DzCorrected(*vtx) - zLepton) : 0;
50 if(dz >= dzMax) continue;
51
52 // check iso cone
53 Double_t dr = MathUtils::DeltaR(muon->Mom(), pfcand->Mom());
54 if(dr<dRMax && dr>=dRMin)
55 iso += pfcand->Pt();
56 }
57
58 return iso;
59 }
60
61 //--------------------------------------------------------------------------------------------------
62 Float_t computePFEleIso(const mithep::Electron *electron,
63 const mithep::Vertex * fVertex,
64 const mithep::Array<mithep::PFCandidate> * fPFCandidates,
65 const Double_t dRMax)
66 //--------------------------------------------------------------------------------------------------
67 {
68 const Double_t dRMin = 0;
69 const Double_t neuPtMin = 1.0;
70 const Double_t dzMax = 0.1;
71
72 Double_t zLepton = (electron->BestTrk()) ? electron->BestTrk()->DzCorrected(*fVertex) : 0.0;
73
74 Float_t iso=0;
75 for(UInt_t ipf=0; ipf<fPFCandidates->GetEntries(); ipf++) {
76 const PFCandidate *pfcand = (PFCandidate*)(fPFCandidates->At(ipf));
77
78 if(!pfcand->HasTrk() && (pfcand->Pt()<=neuPtMin)) continue; // pT cut on neutral particles
79
80 // dz cut
81 Double_t dz = (pfcand->BestTrk()) ? fabs(pfcand->BestTrk()->DzCorrected(*fVertex) - zLepton) : 0;
82 if(dz >= dzMax) continue;
83
84 // remove THE electron
85 if(pfcand->TrackerTrk() && electron->TrackerTrk() && (pfcand->TrackerTrk()==electron->TrackerTrk())) continue;
86 if(pfcand->GsfTrk() && electron->GsfTrk() && (pfcand->GsfTrk()==electron->GsfTrk())) continue;
87
88 // check iso cone
89 Double_t dr = MathUtils::DeltaR(electron->Mom(), pfcand->Mom());
90 if(dr<dRMax && dr>=dRMin) {
91 // eta-strip veto for photons
92 if((pfcand->PFType() == PFCandidate::eGamma) && fabs(electron->Eta() - pfcand->Eta()) < 0.025) continue;
93
94 // Inner cone (one tower = dR < 0.07) veto for non-photon neutrals
95 if(!pfcand->HasTrk() && (pfcand->PFType() == PFCandidate::eNeutralHadron) &&
96 (MathUtils::DeltaR(electron->Mom(), pfcand->Mom()) < 0.07)) continue;
97
98 iso += pfcand->Pt();
99 }
100 }
101
102 return iso;
103 };
104
105 //--------------------------------------------------------------------------------------------------
106 bool pairwiseIsoSelection( ControlFlags &ctrl,
107 vector<SimpleLepton> &lepvec,
108 float rho )
109 //--------------------------------------------------------------------------------------------------
110 {
111
112 bool passiso=true;
113
114 for( int i=0; i<lepvec.size(); i++ )
115 {
116
117 if( !(lepvec[i].is4l) ) continue;
118
119 float effArea_ecal_i, effArea_hcal_i;
120 if( lepvec[i].isEB ) {
121 if( lepvec[i].type == 11 ) {
122 effArea_ecal_i = 0.101;
123 effArea_hcal_i = 0.021;
124 } else {
125 effArea_ecal_i = 0.074;
126 effArea_hcal_i = 0.022;
127 }
128 } else {
129 if( lepvec[i].type == 11 ) {
130 effArea_ecal_i = 0.046;
131 effArea_hcal_i = 0.040;
132 } else {
133 effArea_ecal_i = 0.045;
134 effArea_hcal_i = 0.030;
135 }
136 }
137
138 float isoEcal_corr_i = lepvec[i].isoEcal - (effArea_ecal_i*rho);
139 float isoHcal_corr_i = lepvec[i].isoHcal - (effArea_hcal_i*rho);
140
141 for( int j=i+1; j<lepvec.size(); j++ )
142 {
143
144 if( !(lepvec[j].is4l) ) continue;
145
146 float effArea_ecal_j, effArea_hcal_j;
147 if( lepvec[j].isEB ) {
148 if( lepvec[j].type == 11 ) {
149 effArea_ecal_j = 0.101;
150 effArea_hcal_j = 0.021;
151 } else {
152 effArea_ecal_j = 0.074;
153 effArea_hcal_j = 0.022;
154 }
155 } else {
156 if( lepvec[j].type == 11 ) {
157 effArea_ecal_j = 0.046;
158 effArea_hcal_j = 0.040;
159 } else {
160 effArea_ecal_j = 0.045;
161 effArea_hcal_j = 0.030;
162 }
163 }
164
165 float isoEcal_corr_j = lepvec[j].isoEcal - (effArea_ecal_j*rho);
166 float isoHcal_corr_j = lepvec[j].isoHcal - (effArea_hcal_j*rho);
167 float RIso_i = (lepvec[i].isoTrk+isoEcal_corr_i+isoHcal_corr_i)/lepvec[i].vec.Pt();
168 float RIso_j = (lepvec[j].isoTrk+isoEcal_corr_j+isoHcal_corr_j)/lepvec[j].vec.Pt();
169 float comboIso = RIso_i + RIso_j;
170
171 if( comboIso > 0.35 ) {
172 if( ctrl.debug ) cout << "combo failing for indices: " << i << "," << j << endl;
173 passiso = false;
174 return passiso;
175 }
176 }
177 }
178
179 return passiso;
180 }
181
182 //--------------------------------------------------------------------------------------------------
183 SelectionStatus muonIsoSelection(ControlFlags &ctrl,
184 const mithep::Muon * mu,
185 const mithep::Vertex * vtx,
186 const mithep::Array<mithep::PFCandidate> * fPFCandidateCol )
187 //--------------------------------------------------------------------------------------------------
188 {
189 float reliso = computePFMuonIso(mu,vtx,fPFCandidateCol,0.3)/mu->Pt();
190 bool isEB = (fabs(mu->Eta()) < 1.479 ? 1 : 0 );
191 bool failiso = false;
192 if( isEB && mu->Pt() > 20 && reliso > PFISO_MU_LOOSE_EB_HIGHPT ) {
193 failiso = true;
194 }
195 if( isEB && mu->Pt() < 20 && reliso > PFISO_MU_LOOSE_EB_LOWPT ) {
196 failiso = true;
197 }
198 if( !(isEB) && mu->Pt() > 20 && reliso > PFISO_MU_LOOSE_EE_HIGHPT ) {
199 failiso = true;
200 }
201 if( !(isEB) && mu->Pt() < 20 && reliso > PFISO_MU_LOOSE_EE_LOWPT ) {
202 failiso = true;
203 }
204
205 SelectionStatus status;
206 if( !failiso ) status.setStatus(SelectionStatus::LOOSEISO);
207 if( !failiso ) status.setStatus(SelectionStatus::TIGHTISO);
208 return status;
209
210 };
211
212 //--------------------------------------------------------------------------------------------------
213 SelectionStatus electronIsoSelection(ControlFlags &ctrl,
214 const mithep::Electron * ele,
215 const mithep::Vertex *fVertex,
216 const mithep::Array<mithep::PFCandidate> * fPFCandidates)
217 //--------------------------------------------------------------------------------------------------
218 {
219
220 bool failiso=false;
221
222 float reliso = computePFEleIso(ele,fVertex,fPFCandidates,0.4)/ele->Pt();
223
224 if( ele->IsEB() && ele->Pt() > 20 && reliso > PFISO_ELE_LOOSE_EB_HIGHPT ) {
225 failiso = true;
226 }
227 if( ele->IsEB() && ele->Pt() < 20 && reliso > PFISO_ELE_LOOSE_EB_LOWPT ) {
228 failiso = true;
229 }
230 if( !(ele->IsEB()) && ele->Pt() > 20 && reliso > PFISO_ELE_LOOSE_EE_HIGHPT ) {
231 failiso = true;
232 }
233 if( !(ele->IsEB()) && ele->Pt() < 20 && reliso > PFISO_ELE_LOOSE_EE_LOWPT ) {
234 failiso = true;
235 }
236
237 SelectionStatus status;
238 if( !failiso ) {
239 status.orStatus(SelectionStatus::LOOSEISO);
240 status.orStatus(SelectionStatus::TIGHTISO);
241 }
242 if(ctrl.debug) cout << "returning status : " << hex << status.getStatus() << dec << endl;
243 return status;
244
245 }
246
247
248 bool noIso(ControlFlags &, vector<SimpleLepton> &, float rho) {
249
250 return true;
251 }
252
253
254 //--------------------------------------------------------------------------------------------------
255 SelectionStatus muonIsoMVASelection(ControlFlags &ctrl,
256 const mithep::Muon * mu,
257 const mithep::Vertex * vtx,
258 const mithep::Array<mithep::PFCandidate> * fPFCandidates,
259 const mithep::Array<mithep::PileupEnergyDensity> * fPUEnergyDensity,
260 mithep::MuonTools::EMuonEffectiveAreaTarget EffectiveAreaVersion,
261 vector<const mithep::Muon*> muonsToVeto,
262 vector<const mithep::Electron*> electronsToVeto)
263 //--------------------------------------------------------------------------------------------------
264 {
265
266 if( ctrl.debug ) {
267 cout << "muonIsoMVASelection :: muons to veto " << endl;
268 for( int i=0; i<muonsToVeto.size(); i++ ) {
269 const mithep::Muon * vmu = muonsToVeto[i];
270 cout << "\tpt: " << vmu->Pt()
271 << "\teta: " << vmu->Eta()
272 << "\tphi: " << vmu->Phi()
273 << endl;
274 }
275 cout << "muonIsoMVASelection :: electrson to veto " << endl;
276 for( int i=0; i<electronsToVeto.size(); i++ ) {
277 const mithep::Electron * vel = electronsToVeto[i];
278 cout << "\tpt: " << vel->Pt()
279 << "\teta: " << vel->Eta()
280 << "\tphi: " << vel->Phi()
281 << endl;
282 }
283 }
284 bool failiso=false;
285
286 //
287 // tmp iso rings
288 //
289 Double_t tmpChargedIso_DR0p0To0p1 = 0;
290 Double_t tmpChargedIso_DR0p1To0p2 = 0;
291 Double_t tmpChargedIso_DR0p2To0p3 = 0;
292 Double_t tmpChargedIso_DR0p3To0p4 = 0;
293 Double_t tmpChargedIso_DR0p4To0p5 = 0;
294 Double_t tmpChargedIso_DR0p5To0p7 = 0;
295
296 Double_t tmpGammaIso_DR0p0To0p1 = 0;
297 Double_t tmpGammaIso_DR0p1To0p2 = 0;
298 Double_t tmpGammaIso_DR0p2To0p3 = 0;
299 Double_t tmpGammaIso_DR0p3To0p4 = 0;
300 Double_t tmpGammaIso_DR0p4To0p5 = 0;
301 Double_t tmpGammaIso_DR0p5To0p7 = 0;
302
303 Double_t tmpNeutralHadronIso_DR0p0To0p1 = 0;
304 Double_t tmpNeutralHadronIso_DR0p1To0p2 = 0;
305 Double_t tmpNeutralHadronIso_DR0p2To0p3 = 0;
306 Double_t tmpNeutralHadronIso_DR0p3To0p4 = 0;
307 Double_t tmpNeutralHadronIso_DR0p4To0p5 = 0;
308 Double_t tmpNeutralHadronIso_DR0p5To0p7 = 0;
309
310
311
312 //
313 // final rings for the MVA
314 //
315 Double_t fChargedIso_DR0p0To0p1;
316 Double_t fChargedIso_DR0p1To0p2;
317 Double_t fChargedIso_DR0p2To0p3;
318 Double_t fChargedIso_DR0p3To0p4;
319 Double_t fChargedIso_DR0p4To0p5;
320 Double_t fChargedIso_DR0p5To0p7;
321
322 Double_t fGammaIso_DR0p0To0p1;
323 Double_t fGammaIso_DR0p1To0p2;
324 Double_t fGammaIso_DR0p2To0p3;
325 Double_t fGammaIso_DR0p3To0p4;
326 Double_t fGammaIso_DR0p4To0p5;
327 Double_t fGammaIso_DR0p5To0p7;
328
329 Double_t fNeutralHadronIso_DR0p0To0p1;
330 Double_t fNeutralHadronIso_DR0p1To0p2;
331 Double_t fNeutralHadronIso_DR0p2To0p3;
332 Double_t fNeutralHadronIso_DR0p3To0p4;
333 Double_t fNeutralHadronIso_DR0p4To0p5;
334 Double_t fNeutralHadronIso_DR0p5To0p7;
335
336
337 //
338 //Loop over PF Candidates
339 //
340 for(int k=0; k<fPFCandidates->GetEntries(); ++k) {
341
342 if( !(PFnoPUflag[k]) ) continue; // my PF no PU hack
343
344 const mithep::PFCandidate *pf = (mithep::PFCandidate*)((*fPFCandidates)[k]);
345
346 Double_t deta = (mu->Eta() - pf->Eta());
347 Double_t dphi = mithep::MathUtils::DeltaPhi(Double_t(mu->Phi()),Double_t(pf->Phi()));
348 Double_t dr = mithep::MathUtils::DeltaR(mu->Phi(),mu->Eta(), pf->Phi(), pf->Eta());
349 if (dr > 1.0) continue;
350
351 if (pf->HasTrackerTrk() && (pf->TrackerTrk() == mu->TrackerTrk()) ) continue;
352
353 //
354 // Lepton Footprint Removal
355 //
356 Bool_t IsLeptonFootprint = kFALSE;
357 if (dr < 1.0) {
358
359 //
360 // Check for electrons
361 //
362 for (Int_t q=0; q < electronsToVeto.size(); ++q) {
363 const mithep::Electron *tmpele = electronsToVeto[q];
364 // 4l electron
365 if( pf->HasTrackerTrk() ) {
366 if( pf->TrackerTrk() == tmpele->TrackerTrk() )
367 IsLeptonFootprint = kTRUE;
368 }
369 if( pf->HasGsfTrk() ) {
370 if( pf->GsfTrk() == tmpele->GsfTrk() )
371 IsLeptonFootprint = kTRUE;
372 }
373 // PF charged
374 if (pf->Charge() != 0 && fabs(tmpele->SCluster()->Eta()) >= 1.479
375 && mithep::MathUtils::DeltaR(tmpele->Phi(),tmpele->Eta(), pf->Phi(), pf->Eta()) < 0.015)
376 IsLeptonFootprint = kTRUE;
377 // PF gamma
378 if (abs(pf->PFType()) == PFCandidate::eGamma && fabs(tmpele->SCluster()->Eta()) >= 1.479
379 && mithep::MathUtils::DeltaR(tmpele->Phi(),tmpele->Eta(), pf->Phi(), pf->Eta()) < 0.08)
380 IsLeptonFootprint = kTRUE;
381 } // loop over electrons
382
383 /* KH - commented for sync
384 //
385 // Check for muons
386 //
387 for (Int_t q=0; q < muonsToVeto.size(); ++q) {
388 const mithep::Muon *tmpmu = muonsToVeto[q];
389 // 4l muon
390 if( pf->HasTrackerTrk() ) {
391 if( pf->TrackerTrk() == tmpmu->TrackerTrk() )
392 IsLeptonFootprint = kTRUE;
393 }
394 // PF charged
395 if (pf->Charge() != 0 && mithep::MathUtils::DeltaR(tmpmu->Phi(),tmpmu->Eta(), pf->Phi(), pf->Eta()) < 0.01)
396 IsLeptonFootprint = kTRUE;
397 } // loop over muons
398 */
399
400 if (IsLeptonFootprint)
401 continue;
402
403 //
404 // Charged Iso Rings
405 //
406 if (pf->Charge() != 0 && (pf->HasTrackerTrk()||pf->HasGsfTrk()) ) {
407
408 if( dr < 0.01 ) continue; // only for muon iso mva?
409 if (abs(pf->PFType()) == PFCandidate::eElectron || abs(pf->PFType()) == PFCandidate::eMuon) continue;
410
411 // if( pf->HasTrackerTrk() ) {
412 // if (abs(pf->TrackerTrk()->DzCorrected(vtx)) > 0.2) continue;
413 // if( ctrl.debug ) cout << "charged:: " << pf->PFType() << " " << pf->Pt() << " "
414 // << abs(pf->TrackerTrk()->DzCorrected(vtx)) << " "
415 // << dr << endl;
416 // }
417 // if( pf->HasGsfTrk() ) {
418 // if (abs(pf->GsfTrk()->DzCorrected(vtx)) > 0.2) continue;
419 // if( ctrl.debug ) cout << "charged:: " << pf->PFType() << " " << pf->Pt() << " "
420 // << abs(pf->GsfTrk()->DzCorrected(vtx)) << " "
421 // << dr << endl;
422 // }
423
424 // Footprint Veto
425 if (dr < 0.1) tmpChargedIso_DR0p0To0p1 += pf->Pt();
426 if (dr >= 0.1 && dr < 0.2) tmpChargedIso_DR0p1To0p2 += pf->Pt();
427 if (dr >= 0.2 && dr < 0.3) tmpChargedIso_DR0p2To0p3 += pf->Pt();
428 if (dr >= 0.3 && dr < 0.4) tmpChargedIso_DR0p3To0p4 += pf->Pt();
429 if (dr >= 0.4 && dr < 0.5) tmpChargedIso_DR0p4To0p5 += pf->Pt();
430 if (dr >= 0.5 && dr < 0.7) tmpChargedIso_DR0p5To0p7 += pf->Pt();
431 }
432
433 //
434 // Gamma Iso Rings
435 //
436 else if (abs(pf->PFType()) == PFCandidate::eGamma) {
437 if (dr < 0.1) tmpGammaIso_DR0p0To0p1 += pf->Pt();
438 if (dr >= 0.1 && dr < 0.2) tmpGammaIso_DR0p1To0p2 += pf->Pt();
439 if (dr >= 0.2 && dr < 0.3) tmpGammaIso_DR0p2To0p3 += pf->Pt();
440 if (dr >= 0.3 && dr < 0.4) tmpGammaIso_DR0p3To0p4 += pf->Pt();
441 if (dr >= 0.4 && dr < 0.5) tmpGammaIso_DR0p4To0p5 += pf->Pt();
442 if (dr >= 0.5 && dr < 0.7) tmpGammaIso_DR0p5To0p7 += pf->Pt();
443 }
444
445 //
446 // Other Neutral Iso Rings
447 //
448 else {
449 if (dr < 0.1) tmpNeutralHadronIso_DR0p0To0p1 += pf->Pt();
450 if (dr >= 0.1 && dr < 0.2) tmpNeutralHadronIso_DR0p1To0p2 += pf->Pt();
451 if (dr >= 0.2 && dr < 0.3) tmpNeutralHadronIso_DR0p2To0p3 += pf->Pt();
452 if (dr >= 0.3 && dr < 0.4) tmpNeutralHadronIso_DR0p3To0p4 += pf->Pt();
453 if (dr >= 0.4 && dr < 0.5) tmpNeutralHadronIso_DR0p4To0p5 += pf->Pt();
454 if (dr >= 0.5 && dr < 0.7) tmpNeutralHadronIso_DR0p5To0p7 += pf->Pt();
455 }
456
457 }
458
459 }
460
461 fChargedIso_DR0p0To0p1 = fmin((tmpChargedIso_DR0p0To0p1)/mu->Pt(), 2.5);
462 fChargedIso_DR0p1To0p2 = fmin((tmpChargedIso_DR0p1To0p2)/mu->Pt(), 2.5);
463 fChargedIso_DR0p2To0p3 = fmin((tmpChargedIso_DR0p2To0p3)/mu->Pt(), 2.5);
464 fChargedIso_DR0p3To0p4 = fmin((tmpChargedIso_DR0p3To0p4)/mu->Pt(), 2.5);
465 fChargedIso_DR0p4To0p5 = fmin((tmpChargedIso_DR0p4To0p5)/mu->Pt(), 2.5);
466
467
468 double rho = 0;
469 if (!(isnan(fPUEnergyDensity->At(0)->Rho()) || isinf(fPUEnergyDensity->At(0)->Rho())))
470 rho = fPUEnergyDensity->At(0)->Rho();
471
472 // if (!(isnan(fPUEnergyDensity->At(0)->RhoLowEta()) || isinf(fPUEnergyDensity->At(0)->RhoLowEta())))
473 // rho = fPUEnergyDensity->At(0)->RhoLowEta();
474
475 // WARNING!!!!
476 // hardcode for sync ...
477 EffectiveAreaVersion = muT.kMuEAData2011;
478 // WARNING!!!!
479
480
481 fGammaIso_DR0p0To0p1 = fmax(fmin((tmpGammaIso_DR0p0To0p1
482 -rho*muT.MuonEffectiveArea(muT.kMuGammaIsoDR0p0To0p1,mu->Eta(),EffectiveAreaVersion))/mu->Pt()
483 ,2.5)
484 ,0.0);
485 fGammaIso_DR0p1To0p2 = fmax(fmin((tmpGammaIso_DR0p1To0p2
486 -rho*muT.MuonEffectiveArea(muT.kMuGammaIsoDR0p1To0p2,mu->Eta(),EffectiveAreaVersion))/mu->Pt()
487 ,2.5)
488 ,0.0);
489 fGammaIso_DR0p2To0p3 = fmax(fmin((tmpGammaIso_DR0p2To0p3
490 -rho*muT.MuonEffectiveArea(muT.kMuGammaIsoDR0p2To0p3,mu->Eta(),EffectiveAreaVersion))/mu->Pt()
491 ,2.5)
492 ,0.0);
493 fGammaIso_DR0p3To0p4 = fmax(fmin((tmpGammaIso_DR0p3To0p4
494 -rho*muT.MuonEffectiveArea(muT.kMuGammaIsoDR0p3To0p4,mu->Eta(),EffectiveAreaVersion))/mu->Pt()
495 ,2.5)
496 ,0.0);
497 fGammaIso_DR0p4To0p5 = fmax(fmin((tmpGammaIso_DR0p4To0p5
498 -rho*muT.MuonEffectiveArea(muT.kMuGammaIsoDR0p4To0p5,mu->Eta(),EffectiveAreaVersion))/mu->Pt()
499 ,2.5)
500 ,0.0);
501
502
503
504 fNeutralHadronIso_DR0p0To0p1 = fmax(fmin((tmpNeutralHadronIso_DR0p0To0p1
505 -rho*muT.MuonEffectiveArea(muT.kMuNeutralHadronIsoDR0p0To0p1,
506 mu->Eta(),EffectiveAreaVersion))/mu->Pt()
507 , 2.5)
508 , 0.0);
509 fNeutralHadronIso_DR0p1To0p2 = fmax(fmin((tmpNeutralHadronIso_DR0p1To0p2
510 -rho*muT.MuonEffectiveArea(muT.kMuNeutralHadronIsoDR0p1To0p2,
511 mu->Eta(),EffectiveAreaVersion))/mu->Pt()
512 , 2.5)
513 , 0.0);
514 fNeutralHadronIso_DR0p2To0p3 = fmax(fmin((tmpNeutralHadronIso_DR0p2To0p3
515 -rho*muT.MuonEffectiveArea(muT.kMuNeutralHadronIsoDR0p2To0p3,
516 mu->Eta(),EffectiveAreaVersion))/mu->Pt()
517 , 2.5)
518 , 0.0);
519 fNeutralHadronIso_DR0p3To0p4 = fmax(fmin((tmpNeutralHadronIso_DR0p3To0p4
520 -rho*muT.MuonEffectiveArea(muT.kMuNeutralHadronIsoDR0p3To0p4,
521 mu->Eta(), EffectiveAreaVersion))/mu->Pt()
522 , 2.5)
523 , 0.0);
524 fNeutralHadronIso_DR0p4To0p5 = fmax(fmin((tmpNeutralHadronIso_DR0p4To0p5
525 -rho*muT.MuonEffectiveArea(muT.kMuNeutralHadronIsoDR0p4To0p5,
526 mu->Eta(), EffectiveAreaVersion))/mu->Pt()
527 , 2.5)
528 , 0.0);
529
530
531 double mvaval = muIsoMVA->MVAValue_IsoRings( mu->Pt(),
532 mu->Eta(),
533 mu->IsGlobalMuon(),
534 mu->IsTrackerMuon(),
535 fChargedIso_DR0p0To0p1,
536 fChargedIso_DR0p1To0p2,
537 fChargedIso_DR0p2To0p3,
538 fChargedIso_DR0p3To0p4,
539 fChargedIso_DR0p4To0p5,
540 fGammaIso_DR0p0To0p1,
541 fGammaIso_DR0p1To0p2,
542 fGammaIso_DR0p2To0p3,
543 fGammaIso_DR0p3To0p4,
544 fGammaIso_DR0p4To0p5,
545 fNeutralHadronIso_DR0p0To0p1,
546 fNeutralHadronIso_DR0p1To0p2,
547 fNeutralHadronIso_DR0p2To0p3,
548 fNeutralHadronIso_DR0p3To0p4,
549 fNeutralHadronIso_DR0p4To0p5,
550 ctrl.debug);
551
552 SelectionStatus status;
553 bool pass;
554
555 pass = false;
556 if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
557 && fabs(mu->Eta()) <= 1.5 && mu->Pt() <= 10 && mvaval >= MUON_ISOMVA_LOOSE_FORPFID_CUT_OPT_BIN0) pass = true;
558 else if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
559 && fabs(mu->Eta()) <= 1.5 && mu->Pt() > 10 && mvaval >= MUON_ISOMVA_LOOSE_FORPFID_CUT_OPT_BIN1) pass = true;
560 else if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
561 && fabs(mu->Eta()) > 1.5 && mu->Pt() <= 10 && mvaval >= MUON_ISOMVA_LOOSE_FORPFID_CUT_OPT_BIN2) pass = true;
562 else if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
563 && fabs(mu->Eta()) > 1.5 && mu->Pt() > 10 && mvaval >= MUON_ISOMVA_LOOSE_FORPFID_CUT_OPT_BIN3) pass = true;
564 else if( !(mu->IsGlobalMuon()) && mu->IsTrackerMuon() && mvaval >= MUON_ISOMVA_LOOSE_FORPFID_CUT_OPT_BIN4) pass = true;
565 else if( mu->IsGlobalMuon() && !(mu->IsTrackerMuon()) && mvaval >= MUON_ISOMVA_LOOSE_FORPFID_CUT_OPT_BIN5) pass = true;
566 if( pass ) status.orStatus(SelectionStatus::LOOSEISO);
567
568 /*
569 pass = false;
570 if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
571 && fabs(mu->Eta()) <= 1.5 && mu->Pt() <= 10 && mvaval >= MUON_ISOMVA_TIGHT_FORPFID_CUT_BIN0) pass = true;
572 else if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
573 && fabs(mu->Eta()) <= 1.5 && mu->Pt() > 10 && mvaval >= MUON_ISOMVA_TIGHT_FORPFID_CUT_BIN1) pass = true;
574 else if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
575 && fabs(mu->Eta()) > 1.5 && mu->Pt() <= 10 && mvaval >= MUON_ISOMVA_TIGHT_FORPFID_CUT_BIN2) pass = true;
576 else if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
577 && fabs(mu->Eta()) > 1.5 && mu->Pt() > 10 && mvaval >= MUON_ISOMVA_TIGHT_FORPFID_CUT_BIN3) pass = true;
578 else if( !(mu->IsGlobalMuon()) && mu->IsTrackerMuon() && mvaval >= MUON_ISOMVA_TIGHT_FORPFID_CUT_BIN4) pass = true;
579 else if( mu->IsGlobalMuon() && !(mu->IsTrackerMuon()) && mvaval >= MUON_ISOMVA_TIGHT_FORPFID_CUT_BIN5) pass = true;
580 if( pass ) status.orStatus(SelectionStatus::TIGHTISO);
581 */
582
583 // pass &= (fChargedIso_DR0p0To0p1 + fChargedIso_DR0p1To0p2 + fChargedIso_DR0p2To0p3 < 0.7);
584
585 status.isoMVA = mvaval;
586
587 if(ctrl.debug) {
588 cout << "returning status : " << hex << status.getStatus() << dec << endl;
589 cout << "MVAVAL : " << status.isoMVA << endl;
590 }
591 return status;
592
593 }
594
595
596 //--------------------------------------------------------------------------------------------------
597 SelectionStatus muonIsoMVASelection(ControlFlags &ctrl,
598 const mithep::Muon * mu,
599 const mithep::Vertex * vtx,
600 const mithep::Array<mithep::PFCandidate> * fPFCandidates,
601 float rho,
602 //const mithep::Array<mithep::PileupEnergyDensity> * fPUEnergyDensity,
603 mithep::MuonTools::EMuonEffectiveAreaTarget EffectiveAreaVersion,
604 vector<const mithep::Muon*> muonsToVeto,
605 vector<const mithep::Electron*> electronsToVeto)
606 //--------------------------------------------------------------------------------------------------
607 // hacked version
608 {
609
610 if( ctrl.debug ) {
611 cout << "muonIsoMVASelection :: muons to veto " << endl;
612 for( int i=0; i<muonsToVeto.size(); i++ ) {
613 const mithep::Muon * vmu = muonsToVeto[i];
614 cout << "\tpt: " << vmu->Pt()
615 << "\teta: " << vmu->Eta()
616 << "\tphi: " << vmu->Phi()
617 << endl;
618 }
619 cout << "muonIsoMVASelection :: electrson to veto " << endl;
620 for( int i=0; i<electronsToVeto.size(); i++ ) {
621 const mithep::Electron * vel = electronsToVeto[i];
622 cout << "\tpt: " << vel->Pt()
623 << "\teta: " << vel->Eta()
624 << "\tphi: " << vel->Phi()
625 << endl;
626 }
627 }
628 bool failiso=false;
629
630 //
631 // tmp iso rings
632 //
633 Double_t tmpChargedIso_DR0p0To0p1 = 0;
634 Double_t tmpChargedIso_DR0p1To0p2 = 0;
635 Double_t tmpChargedIso_DR0p2To0p3 = 0;
636 Double_t tmpChargedIso_DR0p3To0p4 = 0;
637 Double_t tmpChargedIso_DR0p4To0p5 = 0;
638 Double_t tmpChargedIso_DR0p5To0p7 = 0;
639
640 Double_t tmpGammaIso_DR0p0To0p1 = 0;
641 Double_t tmpGammaIso_DR0p1To0p2 = 0;
642 Double_t tmpGammaIso_DR0p2To0p3 = 0;
643 Double_t tmpGammaIso_DR0p3To0p4 = 0;
644 Double_t tmpGammaIso_DR0p4To0p5 = 0;
645 Double_t tmpGammaIso_DR0p5To0p7 = 0;
646
647 Double_t tmpNeutralHadronIso_DR0p0To0p1 = 0;
648 Double_t tmpNeutralHadronIso_DR0p1To0p2 = 0;
649 Double_t tmpNeutralHadronIso_DR0p2To0p3 = 0;
650 Double_t tmpNeutralHadronIso_DR0p3To0p4 = 0;
651 Double_t tmpNeutralHadronIso_DR0p4To0p5 = 0;
652 Double_t tmpNeutralHadronIso_DR0p5To0p7 = 0;
653
654
655
656 //
657 // final rings for the MVA
658 //
659 Double_t fChargedIso_DR0p0To0p1;
660 Double_t fChargedIso_DR0p1To0p2;
661 Double_t fChargedIso_DR0p2To0p3;
662 Double_t fChargedIso_DR0p3To0p4;
663 Double_t fChargedIso_DR0p4To0p5;
664 Double_t fChargedIso_DR0p5To0p7;
665
666 Double_t fGammaIso_DR0p0To0p1;
667 Double_t fGammaIso_DR0p1To0p2;
668 Double_t fGammaIso_DR0p2To0p3;
669 Double_t fGammaIso_DR0p3To0p4;
670 Double_t fGammaIso_DR0p4To0p5;
671 Double_t fGammaIso_DR0p5To0p7;
672
673 Double_t fNeutralHadronIso_DR0p0To0p1;
674 Double_t fNeutralHadronIso_DR0p1To0p2;
675 Double_t fNeutralHadronIso_DR0p2To0p3;
676 Double_t fNeutralHadronIso_DR0p3To0p4;
677 Double_t fNeutralHadronIso_DR0p4To0p5;
678 Double_t fNeutralHadronIso_DR0p5To0p7;
679
680
681 //
682 //Loop over PF Candidates
683 //
684 for(int k=0; k<fPFCandidates->GetEntries(); ++k) {
685
686 if( !(PFnoPUflag[k]) ) continue; // my PF no PU hack
687
688 const mithep::PFCandidate *pf = (mithep::PFCandidate*)((*fPFCandidates)[k]);
689
690 Double_t deta = (mu->Eta() - pf->Eta());
691 Double_t dphi = mithep::MathUtils::DeltaPhi(Double_t(mu->Phi()),Double_t(pf->Phi()));
692 Double_t dr = mithep::MathUtils::DeltaR(mu->Phi(),mu->Eta(), pf->Phi(), pf->Eta());
693 if (dr > 1.0) continue;
694
695 if (pf->HasTrackerTrk() && (pf->TrackerTrk() == mu->TrackerTrk()) ) continue;
696
697 //
698 // Lepton Footprint Removal
699 //
700 Bool_t IsLeptonFootprint = kFALSE;
701 if (dr < 1.0) {
702
703 //
704 // Check for electrons
705 //
706 for (Int_t q=0; q < electronsToVeto.size(); ++q) {
707 const mithep::Electron *tmpele = electronsToVeto[q];
708 // 4l electron
709 if( pf->HasTrackerTrk() ) {
710 if( pf->TrackerTrk() == tmpele->TrackerTrk() )
711 IsLeptonFootprint = kTRUE;
712 }
713 if( pf->HasGsfTrk() ) {
714 if( pf->GsfTrk() == tmpele->GsfTrk() )
715 IsLeptonFootprint = kTRUE;
716 }
717 // PF charged
718 if (pf->Charge() != 0 && fabs(tmpele->SCluster()->Eta()) >= 1.479
719 && mithep::MathUtils::DeltaR(tmpele->Phi(),tmpele->Eta(), pf->Phi(), pf->Eta()) < 0.015)
720 IsLeptonFootprint = kTRUE;
721 // PF gamma
722 if (abs(pf->PFType()) == PFCandidate::eGamma && fabs(tmpele->SCluster()->Eta()) >= 1.479
723 && mithep::MathUtils::DeltaR(tmpele->Phi(),tmpele->Eta(), pf->Phi(), pf->Eta()) < 0.08)
724 IsLeptonFootprint = kTRUE;
725 } // loop over electrons
726
727 /* KH - commented for sync
728 //
729 // Check for muons
730 //
731 for (Int_t q=0; q < muonsToVeto.size(); ++q) {
732 const mithep::Muon *tmpmu = muonsToVeto[q];
733 // 4l muon
734 if( pf->HasTrackerTrk() ) {
735 if( pf->TrackerTrk() == tmpmu->TrackerTrk() )
736 IsLeptonFootprint = kTRUE;
737 }
738 // PF charged
739 if (pf->Charge() != 0 && mithep::MathUtils::DeltaR(tmpmu->Phi(),tmpmu->Eta(), pf->Phi(), pf->Eta()) < 0.01)
740 IsLeptonFootprint = kTRUE;
741 } // loop over muons
742 */
743
744 if (IsLeptonFootprint)
745 continue;
746
747 //
748 // Charged Iso Rings
749 //
750 if (pf->Charge() != 0 && (pf->HasTrackerTrk()||pf->HasGsfTrk()) ) {
751
752 if( dr < 0.01 ) continue; // only for muon iso mva?
753 if (abs(pf->PFType()) == PFCandidate::eElectron || abs(pf->PFType()) == PFCandidate::eMuon) continue;
754
755 // if( pf->HasTrackerTrk() ) {
756 // if (abs(pf->TrackerTrk()->DzCorrected(vtx)) > 0.2) continue;
757 // if( ctrl.debug ) cout << "charged:: " << pf->PFType() << " " << pf->Pt() << " "
758 // << abs(pf->TrackerTrk()->DzCorrected(vtx)) << " "
759 // << dr << endl;
760 // }
761 // if( pf->HasGsfTrk() ) {
762 // if (abs(pf->GsfTrk()->DzCorrected(vtx)) > 0.2) continue;
763 // if( ctrl.debug ) cout << "charged:: " << pf->PFType() << " " << pf->Pt() << " "
764 // << abs(pf->GsfTrk()->DzCorrected(vtx)) << " "
765 // << dr << endl;
766 // }
767
768 // Footprint Veto
769 if (dr < 0.1) tmpChargedIso_DR0p0To0p1 += pf->Pt();
770 if (dr >= 0.1 && dr < 0.2) tmpChargedIso_DR0p1To0p2 += pf->Pt();
771 if (dr >= 0.2 && dr < 0.3) tmpChargedIso_DR0p2To0p3 += pf->Pt();
772 if (dr >= 0.3 && dr < 0.4) tmpChargedIso_DR0p3To0p4 += pf->Pt();
773 if (dr >= 0.4 && dr < 0.5) tmpChargedIso_DR0p4To0p5 += pf->Pt();
774 if (dr >= 0.5 && dr < 0.7) tmpChargedIso_DR0p5To0p7 += pf->Pt();
775 }
776
777 //
778 // Gamma Iso Rings
779 //
780 else if (abs(pf->PFType()) == PFCandidate::eGamma) {
781 if (dr < 0.1) tmpGammaIso_DR0p0To0p1 += pf->Pt();
782 if (dr >= 0.1 && dr < 0.2) tmpGammaIso_DR0p1To0p2 += pf->Pt();
783 if (dr >= 0.2 && dr < 0.3) tmpGammaIso_DR0p2To0p3 += pf->Pt();
784 if (dr >= 0.3 && dr < 0.4) tmpGammaIso_DR0p3To0p4 += pf->Pt();
785 if (dr >= 0.4 && dr < 0.5) tmpGammaIso_DR0p4To0p5 += pf->Pt();
786 if (dr >= 0.5 && dr < 0.7) tmpGammaIso_DR0p5To0p7 += pf->Pt();
787 }
788
789 //
790 // Other Neutral Iso Rings
791 //
792 else {
793 if (dr < 0.1) tmpNeutralHadronIso_DR0p0To0p1 += pf->Pt();
794 if (dr >= 0.1 && dr < 0.2) tmpNeutralHadronIso_DR0p1To0p2 += pf->Pt();
795 if (dr >= 0.2 && dr < 0.3) tmpNeutralHadronIso_DR0p2To0p3 += pf->Pt();
796 if (dr >= 0.3 && dr < 0.4) tmpNeutralHadronIso_DR0p3To0p4 += pf->Pt();
797 if (dr >= 0.4 && dr < 0.5) tmpNeutralHadronIso_DR0p4To0p5 += pf->Pt();
798 if (dr >= 0.5 && dr < 0.7) tmpNeutralHadronIso_DR0p5To0p7 += pf->Pt();
799 }
800
801 }
802
803 }
804
805 fChargedIso_DR0p0To0p1 = fmin((tmpChargedIso_DR0p0To0p1)/mu->Pt(), 2.5);
806 fChargedIso_DR0p1To0p2 = fmin((tmpChargedIso_DR0p1To0p2)/mu->Pt(), 2.5);
807 fChargedIso_DR0p2To0p3 = fmin((tmpChargedIso_DR0p2To0p3)/mu->Pt(), 2.5);
808 fChargedIso_DR0p3To0p4 = fmin((tmpChargedIso_DR0p3To0p4)/mu->Pt(), 2.5);
809 fChargedIso_DR0p4To0p5 = fmin((tmpChargedIso_DR0p4To0p5)/mu->Pt(), 2.5);
810
811
812 // double rho = 0;
813 // if (!(isnan(fPUEnergyDensity->At(0)->Rho()) || isinf(fPUEnergyDensity->At(0)->Rho())))
814 // rho = fPUEnergyDensity->At(0)->Rho();
815 // if (!(isnan(fPUEnergyDensity->At(0)->RhoLowEta()) || isinf(fPUEnergyDensity->At(0)->RhoLowEta())))
816 // rho = fPUEnergyDensity->At(0)->RhoLowEta();
817
818 // WARNING!!!!
819 // hardcode for sync ...
820 EffectiveAreaVersion = muT.kMuEAData2011;
821 // WARNING!!!!
822
823
824 fGammaIso_DR0p0To0p1 = fmax(fmin((tmpGammaIso_DR0p0To0p1
825 -rho*muT.MuonEffectiveArea(muT.kMuGammaIsoDR0p0To0p1,mu->Eta(),EffectiveAreaVersion))/mu->Pt()
826 ,2.5)
827 ,0.0);
828 fGammaIso_DR0p1To0p2 = fmax(fmin((tmpGammaIso_DR0p1To0p2
829 -rho*muT.MuonEffectiveArea(muT.kMuGammaIsoDR0p1To0p2,mu->Eta(),EffectiveAreaVersion))/mu->Pt()
830 ,2.5)
831 ,0.0);
832 fGammaIso_DR0p2To0p3 = fmax(fmin((tmpGammaIso_DR0p2To0p3
833 -rho*muT.MuonEffectiveArea(muT.kMuGammaIsoDR0p2To0p3,mu->Eta(),EffectiveAreaVersion))/mu->Pt()
834 ,2.5)
835 ,0.0);
836 fGammaIso_DR0p3To0p4 = fmax(fmin((tmpGammaIso_DR0p3To0p4
837 -rho*muT.MuonEffectiveArea(muT.kMuGammaIsoDR0p3To0p4,mu->Eta(),EffectiveAreaVersion))/mu->Pt()
838 ,2.5)
839 ,0.0);
840 fGammaIso_DR0p4To0p5 = fmax(fmin((tmpGammaIso_DR0p4To0p5
841 -rho*muT.MuonEffectiveArea(muT.kMuGammaIsoDR0p4To0p5,mu->Eta(),EffectiveAreaVersion))/mu->Pt()
842 ,2.5)
843 ,0.0);
844
845
846
847 fNeutralHadronIso_DR0p0To0p1 = fmax(fmin((tmpNeutralHadronIso_DR0p0To0p1
848 -rho*muT.MuonEffectiveArea(muT.kMuNeutralHadronIsoDR0p0To0p1,
849 mu->Eta(),EffectiveAreaVersion))/mu->Pt()
850 , 2.5)
851 , 0.0);
852 fNeutralHadronIso_DR0p1To0p2 = fmax(fmin((tmpNeutralHadronIso_DR0p1To0p2
853 -rho*muT.MuonEffectiveArea(muT.kMuNeutralHadronIsoDR0p1To0p2,
854 mu->Eta(),EffectiveAreaVersion))/mu->Pt()
855 , 2.5)
856 , 0.0);
857 fNeutralHadronIso_DR0p2To0p3 = fmax(fmin((tmpNeutralHadronIso_DR0p2To0p3
858 -rho*muT.MuonEffectiveArea(muT.kMuNeutralHadronIsoDR0p2To0p3,
859 mu->Eta(),EffectiveAreaVersion))/mu->Pt()
860 , 2.5)
861 , 0.0);
862 fNeutralHadronIso_DR0p3To0p4 = fmax(fmin((tmpNeutralHadronIso_DR0p3To0p4
863 -rho*muT.MuonEffectiveArea(muT.kMuNeutralHadronIsoDR0p3To0p4,
864 mu->Eta(), EffectiveAreaVersion))/mu->Pt()
865 , 2.5)
866 , 0.0);
867 fNeutralHadronIso_DR0p4To0p5 = fmax(fmin((tmpNeutralHadronIso_DR0p4To0p5
868 -rho*muT.MuonEffectiveArea(muT.kMuNeutralHadronIsoDR0p4To0p5,
869 mu->Eta(), EffectiveAreaVersion))/mu->Pt()
870 , 2.5)
871 , 0.0);
872
873
874 double mvaval = muIsoMVA->MVAValue_IsoRings( mu->Pt(),
875 mu->Eta(),
876 mu->IsGlobalMuon(),
877 mu->IsTrackerMuon(),
878 fChargedIso_DR0p0To0p1,
879 fChargedIso_DR0p1To0p2,
880 fChargedIso_DR0p2To0p3,
881 fChargedIso_DR0p3To0p4,
882 fChargedIso_DR0p4To0p5,
883 fGammaIso_DR0p0To0p1,
884 fGammaIso_DR0p1To0p2,
885 fGammaIso_DR0p2To0p3,
886 fGammaIso_DR0p3To0p4,
887 fGammaIso_DR0p4To0p5,
888 fNeutralHadronIso_DR0p0To0p1,
889 fNeutralHadronIso_DR0p1To0p2,
890 fNeutralHadronIso_DR0p2To0p3,
891 fNeutralHadronIso_DR0p3To0p4,
892 fNeutralHadronIso_DR0p4To0p5,
893 ctrl.debug);
894
895 SelectionStatus status;
896 bool pass;
897
898 pass = false;
899 if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
900 && fabs(mu->Eta()) <= 1.5 && mu->Pt() <= 10 && mvaval >= MUON_ISOMVA_LOOSE_FORPFID_CUT_BIN0) pass = true;
901 else if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
902 && fabs(mu->Eta()) <= 1.5 && mu->Pt() > 10 && mvaval >= MUON_ISOMVA_LOOSE_FORPFID_CUT_BIN1) pass = true;
903 else if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
904 && fabs(mu->Eta()) > 1.5 && mu->Pt() <= 10 && mvaval >= MUON_ISOMVA_LOOSE_FORPFID_CUT_BIN2) pass = true;
905 else if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
906 && fabs(mu->Eta()) > 1.5 && mu->Pt() > 10 && mvaval >= MUON_ISOMVA_LOOSE_FORPFID_CUT_BIN3) pass = true;
907 else if( !(mu->IsGlobalMuon()) && mu->IsTrackerMuon() && mvaval >= MUON_ISOMVA_LOOSE_FORPFID_CUT_BIN4) pass = true;
908 else if( mu->IsGlobalMuon() && !(mu->IsTrackerMuon()) && mvaval >= MUON_ISOMVA_LOOSE_FORPFID_CUT_BIN5) pass = true;
909 if( pass ) status.orStatus(SelectionStatus::LOOSEISO);
910
911 /*
912 pass = false;
913 if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
914 && fabs(mu->Eta()) <= 1.5 && mu->Pt() <= 10 && mvaval >= MUON_ISOMVA_TIGHT_FORPFID_CUT_BIN0) pass = true;
915 else if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
916 && fabs(mu->Eta()) <= 1.5 && mu->Pt() > 10 && mvaval >= MUON_ISOMVA_TIGHT_FORPFID_CUT_BIN1) pass = true;
917 else if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
918 && fabs(mu->Eta()) > 1.5 && mu->Pt() <= 10 && mvaval >= MUON_ISOMVA_TIGHT_FORPFID_CUT_BIN2) pass = true;
919 else if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
920 && fabs(mu->Eta()) > 1.5 && mu->Pt() > 10 && mvaval >= MUON_ISOMVA_TIGHT_FORPFID_CUT_BIN3) pass = true;
921 else if( !(mu->IsGlobalMuon()) && mu->IsTrackerMuon() && mvaval >= MUON_ISOMVA_TIGHT_FORPFID_CUT_BIN4) pass = true;
922 else if( mu->IsGlobalMuon() && !(mu->IsTrackerMuon()) && mvaval >= MUON_ISOMVA_TIGHT_FORPFID_CUT_BIN5) pass = true;
923 if( pass ) status.orStatus(SelectionStatus::TIGHTISO);
924 */
925
926 // pass &= (fChargedIso_DR0p0To0p1 + fChargedIso_DR0p1To0p2 + fChargedIso_DR0p2To0p3 < 0.7);
927
928 status.isoMVA = mvaval;
929
930 if(ctrl.debug) {
931 cout << "returning status : " << hex << status.getStatus() << dec << endl;
932 cout << "MVAVAL : " << status.isoMVA << endl;
933 }
934 return status;
935
936 }
937
938
939 //--------------------------------------------------------------------------------------------------
940 void initMuonIsoMVA() {
941 //--------------------------------------------------------------------------------------------------
942 muIsoMVA = new mithep::MuonIDMVA();
943 vector<string> weightFiles;
944 weightFiles.push_back("./data/MuonIsoMVAWeights/MuonIsoMVA_BDTG_V0_barrel_lowpt.weights.xml");
945 weightFiles.push_back("./data/MuonIsoMVAWeights/MuonIsoMVA_BDTG_V0_barrel_highpt.weights.xml");
946 weightFiles.push_back("./data/MuonIsoMVAWeights/MuonIsoMVA_BDTG_V0_endcap_lowpt.weights.xml");
947 weightFiles.push_back("./data/MuonIsoMVAWeights/MuonIsoMVA_BDTG_V0_endcap_highpt.weights.xml");
948 weightFiles.push_back("./data/MuonIsoMVAWeights/MuonIsoMVA_BDTG_V0_tracker.weights.xml");
949 weightFiles.push_back("./data/MuonIsoMVAWeights/MuonIsoMVA_BDTG_V0_global.weights.xml");
950 muIsoMVA->Initialize( "MuonIsoMVA",
951 mithep::MuonIDMVA::kIsoRingsV0,
952 kTRUE, weightFiles);
953 }
954
955
956
957
958 //--------------------------------------------------------------------------------------------------
959 double muonPFIso04(ControlFlags &ctrl,
960 const mithep::Muon * mu,
961 const mithep::Vertex * vtx,
962 const mithep::Array<mithep::PFCandidate> * fPFCandidates,
963 const mithep::Array<mithep::PileupEnergyDensity> * fPUEnergyDensity,
964 mithep::MuonTools::EMuonEffectiveAreaTarget EffectiveAreaVersion,
965 vector<const mithep::Muon*> muonsToVeto,
966 vector<const mithep::Electron*> electronsToVeto)
967 //--------------------------------------------------------------------------------------------------
968 {
969
970 extern double gChargedIso;
971 extern double gGammaIso;
972 extern double gNeutralIso;
973
974 if( ctrl.debug ) {
975 cout << "muonIsoMVASelection :: muons to veto " << endl;
976 for( int i=0; i<muonsToVeto.size(); i++ ) {
977 const mithep::Muon * vmu = muonsToVeto[i];
978 cout << "\tpt: " << vmu->Pt()
979 << "\teta: " << vmu->Eta()
980 << "\tphi: " << vmu->Phi()
981 << endl;
982 }
983 cout << "muonIsoMVASelection :: electrson to veto " << endl;
984 for( int i=0; i<electronsToVeto.size(); i++ ) {
985 const mithep::Electron * vel = electronsToVeto[i];
986 cout << "\tpt: " << vel->Pt()
987 << "\teta: " << vel->Eta()
988 << "\tphi: " << vel->Phi()
989 << endl;
990 }
991 }
992
993 //
994 // final iso
995 //
996 Double_t fChargedIso = 0.0;
997 Double_t fGammaIso = 0.0;
998 Double_t fNeutralHadronIso = 0.0;
999
1000 //
1001 //Loop over PF Candidates
1002 //
1003 for(int k=0; k<fPFCandidates->GetEntries(); ++k) {
1004
1005 if( !(PFnoPUflag[k]) ) continue; // my PF no PU hack
1006 const mithep::PFCandidate *pf = (mithep::PFCandidate*)((*fPFCandidates)[k]);
1007
1008 Double_t deta = (mu->Eta() - pf->Eta());
1009 Double_t dphi = mithep::MathUtils::DeltaPhi(Double_t(mu->Phi()),Double_t(pf->Phi()));
1010 Double_t dr = mithep::MathUtils::DeltaR(mu->Phi(),mu->Eta(), pf->Phi(), pf->Eta());
1011 if (dr > 0.4) continue;
1012
1013 if (pf->HasTrackerTrk() && (pf->TrackerTrk() == mu->TrackerTrk()) ) continue;
1014
1015 //
1016 // Lepton Footprint Removal
1017 //
1018 Bool_t IsLeptonFootprint = kFALSE;
1019 if (dr < 1.0) {
1020
1021 //
1022 // Check for electrons
1023 //
1024 for (Int_t q=0; q < electronsToVeto.size(); ++q) {
1025 const mithep::Electron *tmpele = electronsToVeto[q];
1026 // 4l electron
1027 if( pf->HasTrackerTrk() ) {
1028 if( pf->TrackerTrk() == tmpele->TrackerTrk() )
1029 IsLeptonFootprint = kTRUE;
1030 }
1031 if( pf->HasGsfTrk() ) {
1032 if( pf->GsfTrk() == tmpele->GsfTrk() )
1033 IsLeptonFootprint = kTRUE;
1034 }
1035 // PF charged
1036 if (pf->Charge() != 0 && fabs(tmpele->SCluster()->Eta()) > 1.479
1037 && mithep::MathUtils::DeltaR(tmpele->Phi(),tmpele->Eta(), pf->Phi(), pf->Eta()) < 0.015) {
1038 if( ctrl.debug) cout << "\tcharged trk, dR ("
1039 << mithep::MathUtils::DeltaR(tmpele->Phi(),tmpele->Eta(), pf->Phi(), pf->Eta())
1040 << " matches 4L ele ..." << endl;
1041 IsLeptonFootprint = kTRUE;
1042 }
1043 // PF gamma
1044 if (abs(pf->PFType()) == PFCandidate::eGamma && fabs(tmpele->SCluster()->Eta()) > 1.479
1045 && mithep::MathUtils::DeltaR(tmpele->Phi(),tmpele->Eta(), pf->Phi(), pf->Eta()) < 0.08)
1046 IsLeptonFootprint = kTRUE;
1047 } // loop over electrons
1048
1049 /* KH - comment for sync
1050 //
1051 // Check for muons
1052 //
1053 for (Int_t q=0; q < muonsToVeto.size(); ++q) {
1054 const mithep::Muon *tmpmu = muonsToVeto[q];
1055 // 4l muon
1056 if( pf->HasTrackerTrk() ) {
1057 if( pf->TrackerTrk() == tmpmu->TrackerTrk() )
1058 IsLeptonFootprint = kTRUE;
1059 }
1060 // PF charged
1061 if (pf->Charge() != 0 && mithep::MathUtils::DeltaR(tmpmu->Phi(),tmpmu->Eta(), pf->Phi(), pf->Eta()) < 0.01)
1062 IsLeptonFootprint = kTRUE;
1063 } // loop over muons
1064 */
1065
1066 if (IsLeptonFootprint)
1067 continue;
1068
1069 //
1070 // Charged Iso
1071 //
1072 if (pf->Charge() != 0 && (pf->HasTrackerTrk()||pf->HasGsfTrk()) ) {
1073
1074 //if( dr < 0.01 ) continue; // only for muon iso mva?
1075 if (abs(pf->PFType()) == PFCandidate::eElectron || abs(pf->PFType()) == PFCandidate::eMuon) continue;
1076
1077
1078 // if( pf->HasTrackerTrk() ) {
1079 // if (abs(pf->TrackerTrk()->DzCorrected(vtx)) > 0.2) continue;
1080 // if( ctrl.debug ) cout << "charged:: " << pf->PFType() << " " << pf->Pt() << " "
1081 // << abs(pf->TrackerTrk()->DzCorrected(vtx)) << " "
1082 // << dr << endl;
1083 // }
1084 // if( pf->HasGsfTrk() ) {
1085 // if (abs(pf->GsfTrk()->DzCorrected(vtx)) > 0.2) continue;
1086 // if( ctrl.debug ) cout << "charged:: " << pf->PFType() << " " << pf->Pt() << " "
1087 // << abs(pf->GsfTrk()->DzCorrected(vtx)) << " "
1088 // << dr << endl;
1089 // }
1090
1091
1092 fChargedIso += pf->Pt();
1093 }
1094
1095 //
1096 // Gamma Iso
1097 //
1098 else if (abs(pf->PFType()) == PFCandidate::eGamma) {
1099 // KH, add to sync
1100 if( pf->Pt() > 0.5 )
1101 fGammaIso += pf->Pt();
1102 }
1103
1104 //
1105 // Other Neutrals
1106 //
1107 else {
1108 // KH, add to sync
1109 if( pf->Pt() > 0.5 )
1110 fNeutralHadronIso += pf->Pt();
1111 }
1112
1113 }
1114
1115 }
1116
1117 double rho=0;
1118 if( (EffectiveAreaVersion == mithep::MuonTools::kMuEAFall11MC) ||
1119 (EffectiveAreaVersion == mithep::MuonTools::kMuEAData2011) ) {
1120 if (!(isnan(fPUEnergyDensity->At(0)->RhoKt6PFJetsForIso25()) ||
1121 isinf(fPUEnergyDensity->At(0)->RhoKt6PFJetsForIso25())))
1122 rho = fPUEnergyDensity->At(0)->RhoKt6PFJetsForIso25();
1123 // !!!!!!!!!!!!! TMP HACK FOR SYNC !!!!!!!!!!!!!!!!!!!!!
1124 EffectiveAreaVersion = mithep::MuonTools::kMuEAData2011;
1125 // !!!!!!!!!!!!! TMP HACK FOR SYNC !!!!!!!!!!!!!!!!!!!!!
1126 } else {
1127 if (!(isnan(fPUEnergyDensity->At(0)->RhoKt6PFJetsCentralNeutral()) ||
1128 isinf(fPUEnergyDensity->At(0)->RhoKt6PFJetsCentralNeutral())))
1129 rho = fPUEnergyDensity->At(0)->RhoKt6PFJetsCentralNeutral();
1130 // !!!!!!!!!!!!! TMP HACK FOR SYNC !!!!!!!!!!!!!!!!!!!!!
1131 EffectiveAreaVersion = mithep::MuonTools::kMuEAData2012;
1132 // !!!!!!!!!!!!! TMP HACK FOR SYNC !!!!!!!!!!!!!!!!!!!!!
1133 }
1134 if(ctrl.debug) cout << "rho: " << rho << endl;
1135
1136 double pfIso = fChargedIso + fmax(0.0,(fGammaIso + fNeutralHadronIso
1137 -rho*muT.MuonEffectiveArea(muT.kMuGammaAndNeutralHadronIso04,
1138 mu->Eta(),EffectiveAreaVersion)));
1139 gChargedIso = fChargedIso;
1140 gGammaIso = fGammaIso;
1141 gNeutralIso = fNeutralHadronIso;
1142
1143 return pfIso;
1144 }
1145
1146
1147
1148
1149 //--------------------------------------------------------------------------------------------------
1150 // hacked version
1151 double muonPFIso04(ControlFlags &ctrl,
1152 const mithep::Muon * mu,
1153 const mithep::Vertex * vtx,
1154 const mithep::Array<mithep::PFCandidate> * fPFCandidates,
1155 float rho,
1156 mithep::MuonTools::EMuonEffectiveAreaTarget EffectiveAreaVersion,
1157 vector<const mithep::Muon*> muonsToVeto,
1158 vector<const mithep::Electron*> electronsToVeto)
1159 //--------------------------------------------------------------------------------------------------
1160 {
1161
1162 extern double gChargedIso;
1163 extern double gGammaIso;
1164 extern double gNeutralIso;
1165
1166 if( ctrl.debug ) {
1167 cout << "muonIsoMVASelection :: muons to veto " << endl;
1168 for( int i=0; i<muonsToVeto.size(); i++ ) {
1169 const mithep::Muon * vmu = muonsToVeto[i];
1170 cout << "\tpt: " << vmu->Pt()
1171 << "\teta: " << vmu->Eta()
1172 << "\tphi: " << vmu->Phi()
1173 << endl;
1174 }
1175 cout << "muonIsoMVASelection :: electrson to veto " << endl;
1176 for( int i=0; i<electronsToVeto.size(); i++ ) {
1177 const mithep::Electron * vel = electronsToVeto[i];
1178 cout << "\tpt: " << vel->Pt()
1179 << "\teta: " << vel->Eta()
1180 << "\tphi: " << vel->Phi()
1181 << endl;
1182 }
1183 }
1184
1185 //
1186 // final iso
1187 //
1188 Double_t fChargedIso = 0.0;
1189 Double_t fGammaIso = 0.0;
1190 Double_t fNeutralHadronIso = 0.0;
1191
1192 //
1193 //Loop over PF Candidates
1194 //
1195 for(int k=0; k<fPFCandidates->GetEntries(); ++k) {
1196
1197 if( !(PFnoPUflag[k]) ) continue; // my PF no PU hack
1198 const mithep::PFCandidate *pf = (mithep::PFCandidate*)((*fPFCandidates)[k]);
1199
1200 Double_t deta = (mu->Eta() - pf->Eta());
1201 Double_t dphi = mithep::MathUtils::DeltaPhi(Double_t(mu->Phi()),Double_t(pf->Phi()));
1202 Double_t dr = mithep::MathUtils::DeltaR(mu->Phi(),mu->Eta(), pf->Phi(), pf->Eta());
1203 if (dr > 0.4) continue;
1204
1205 if (pf->HasTrackerTrk() && (pf->TrackerTrk() == mu->TrackerTrk()) ) continue;
1206
1207 //
1208 // Lepton Footprint Removal
1209 //
1210 Bool_t IsLeptonFootprint = kFALSE;
1211 if (dr < 1.0) {
1212
1213 //
1214 // Check for electrons
1215 //
1216 for (Int_t q=0; q < electronsToVeto.size(); ++q) {
1217 const mithep::Electron *tmpele = electronsToVeto[q];
1218 // 4l electron
1219 if( pf->HasTrackerTrk() ) {
1220 if( pf->TrackerTrk() == tmpele->TrackerTrk() )
1221 IsLeptonFootprint = kTRUE;
1222 }
1223 if( pf->HasGsfTrk() ) {
1224 if( pf->GsfTrk() == tmpele->GsfTrk() )
1225 IsLeptonFootprint = kTRUE;
1226 }
1227 // PF charged
1228 if (pf->Charge() != 0 && fabs(tmpele->SCluster()->Eta()) > 1.479
1229 && mithep::MathUtils::DeltaR(tmpele->Phi(),tmpele->Eta(), pf->Phi(), pf->Eta()) < 0.015)
1230 IsLeptonFootprint = kTRUE;
1231 // PF gamma
1232 if (abs(pf->PFType()) == PFCandidate::eGamma && fabs(tmpele->SCluster()->Eta()) > 1.479
1233 && mithep::MathUtils::DeltaR(tmpele->Phi(),tmpele->Eta(), pf->Phi(), pf->Eta()) < 0.08)
1234 IsLeptonFootprint = kTRUE;
1235 } // loop over electrons
1236
1237 /* KH - comment for sync
1238 //
1239 // Check for muons
1240 //
1241 for (Int_t q=0; q < muonsToVeto.size(); ++q) {
1242 const mithep::Muon *tmpmu = muonsToVeto[q];
1243 // 4l muon
1244 if( pf->HasTrackerTrk() ) {
1245 if( pf->TrackerTrk() == tmpmu->TrackerTrk() )
1246 IsLeptonFootprint = kTRUE;
1247 }
1248 // PF charged
1249 if (pf->Charge() != 0 && mithep::MathUtils::DeltaR(tmpmu->Phi(),tmpmu->Eta(), pf->Phi(), pf->Eta()) < 0.01)
1250 IsLeptonFootprint = kTRUE;
1251 } // loop over muons
1252 */
1253
1254 if (IsLeptonFootprint)
1255 continue;
1256
1257 //
1258 // Charged Iso
1259 //
1260 if (pf->Charge() != 0 && (pf->HasTrackerTrk()||pf->HasGsfTrk()) ) {
1261
1262 //if( dr < 0.01 ) continue; // only for muon iso mva?
1263 if (abs(pf->PFType()) == PFCandidate::eElectron || abs(pf->PFType()) == PFCandidate::eMuon) continue;
1264
1265
1266 // if( pf->HasTrackerTrk() ) {
1267 // if (abs(pf->TrackerTrk()->DzCorrected(vtx)) > 0.2) continue;
1268 // if( ctrl.debug ) cout << "charged:: " << pf->PFType() << " " << pf->Pt() << " "
1269 // << abs(pf->TrackerTrk()->DzCorrected(vtx)) << " "
1270 // << dr << endl;
1271 // }
1272 // if( pf->HasGsfTrk() ) {
1273 // if (abs(pf->GsfTrk()->DzCorrected(vtx)) > 0.2) continue;
1274 // if( ctrl.debug ) cout << "charged:: " << pf->PFType() << " " << pf->Pt() << " "
1275 // << abs(pf->GsfTrk()->DzCorrected(vtx)) << " "
1276 // << dr << endl;
1277 // }
1278
1279
1280 fChargedIso += pf->Pt();
1281 }
1282
1283 //
1284 // Gamma Iso
1285 //
1286 else if (abs(pf->PFType()) == PFCandidate::eGamma) {
1287 // KH, add to sync
1288 if( pf->Pt() > 0.5 )
1289 fGammaIso += pf->Pt();
1290 }
1291
1292 //
1293 // Other Neutrals
1294 //
1295 else {
1296 // KH, add to sync
1297 if( pf->Pt() > 0.5 )
1298 fNeutralHadronIso += pf->Pt();
1299 }
1300
1301 }
1302
1303 }
1304
1305 // double rho = 0;
1306 // if (!(isnan(fPUEnergyDensity->At(0)->Rho()) || isinf(fPUEnergyDensity->At(0)->Rho())))
1307 // rho = fPUEnergyDensity->At(0)->Rho();
1308
1309 // WARNING!!!!
1310 // hardcode for sync ...
1311 EffectiveAreaVersion = muT.kMuEAData2011;
1312 // WARNING!!!!
1313
1314
1315 double pfIso = fChargedIso + fmax(0.0,(fGammaIso + fNeutralHadronIso
1316 -rho*muT.MuonEffectiveArea(muT.kMuGammaAndNeutralHadronIso04,
1317 mu->Eta(),EffectiveAreaVersion)));
1318 gChargedIso = fChargedIso;
1319 gGammaIso = fGammaIso;
1320 gNeutralIso = fNeutralHadronIso;
1321
1322 return pfIso;
1323 }
1324
1325
1326 //--------------------------------------------------------------------------------------------------
1327 SelectionStatus muonReferenceIsoSelection(ControlFlags &ctrl,
1328 const mithep::Muon * mu,
1329 const mithep::Vertex * vtx,
1330 const mithep::Array<mithep::PFCandidate> * fPFCandidates,
1331 const mithep::Array<mithep::PileupEnergyDensity> * fPUEnergyDensity,
1332 mithep::MuonTools::EMuonEffectiveAreaTarget EffectiveAreaVersion,
1333 vector<const mithep::Muon*> muonsToVeto,
1334 vector<const mithep::Electron*> electronsToVeto)
1335 //--------------------------------------------------------------------------------------------------
1336 {
1337
1338 SelectionStatus status;
1339
1340 double pfIso = muonPFIso04( ctrl, mu, vtx, fPFCandidates, fPUEnergyDensity,
1341 EffectiveAreaVersion, muonsToVeto ,electronsToVeto );
1342 // cout << "--------------> setting muon isoPF04 to" << pfIso << endl;
1343 status.isoPF04 = pfIso;
1344 status.chisoPF04 = gChargedIso;
1345 status.gaisoPF04 = gGammaIso;
1346 status.neisoPF04 = gNeutralIso;
1347
1348 bool pass = false;
1349 if( (pfIso/mu->Pt()) < MUON_REFERENCE_PFISO_CUT ) pass = true;
1350
1351 if( pass ) {
1352 status.orStatus(SelectionStatus::LOOSEISO);
1353 status.orStatus(SelectionStatus::TIGHTISO);
1354 }
1355 if(ctrl.debug) cout << "returning status : " << hex << status.getStatus() << dec << endl;
1356 return status;
1357
1358 }
1359
1360
1361 //--------------------------------------------------------------------------------------------------
1362 // hacked version
1363 SelectionStatus muonReferenceIsoSelection(ControlFlags &ctrl,
1364 const mithep::Muon * mu,
1365 const mithep::Vertex * vtx,
1366 const mithep::Array<mithep::PFCandidate> * fPFCandidates,
1367 float rho,
1368 mithep::MuonTools::EMuonEffectiveAreaTarget EffectiveAreaVersion,
1369 vector<const mithep::Muon*> muonsToVeto,
1370 vector<const mithep::Electron*> electronsToVeto)
1371 //--------------------------------------------------------------------------------------------------
1372 {
1373
1374 SelectionStatus status;
1375
1376 double pfIso = muonPFIso04( ctrl, mu, vtx, fPFCandidates, rho,
1377 EffectiveAreaVersion, muonsToVeto ,electronsToVeto );
1378
1379 status.isoPF04 = pfIso;
1380 status.chisoPF04 = gChargedIso;
1381 status.gaisoPF04 = gGammaIso;
1382 status.neisoPF04 = gNeutralIso;
1383
1384 bool pass = false;
1385 if( (pfIso/mu->Pt()) < MUON_REFERENCE_PFISO_CUT ) pass = true;
1386
1387 if( pass ) {
1388 status.orStatus(SelectionStatus::LOOSEISO);
1389 status.orStatus(SelectionStatus::TIGHTISO);
1390 }
1391 if(ctrl.debug) cout << "returning status : " << hex << status.getStatus() << dec << endl;
1392 return status;
1393
1394 }
1395
1396
1397
1398
1399 //--------------------------------------------------------------------------------------------------
1400 SelectionStatus electronIsoMVASelection(ControlFlags &ctrl,
1401 const mithep::Electron * ele,
1402 const mithep::Vertex * vtx,
1403 const mithep::Array<mithep::PFCandidate> * fPFCandidates,
1404 const mithep::Array<mithep::PileupEnergyDensity> * fPUEnergyDensity,
1405 mithep::ElectronTools::EElectronEffectiveAreaTarget EffectiveAreaVersion,
1406 vector<const mithep::Muon*> muonsToVeto,
1407 vector<const mithep::Electron*> electronsToVeto)
1408 //--------------------------------------------------------------------------------------------------
1409 {
1410
1411 if( ctrl.debug ) {
1412 cout << "electronIsoMVASelection :: muons to veto " << endl;
1413 for( int i=0; i<muonsToVeto.size(); i++ ) {
1414 const mithep::Muon * vmu = muonsToVeto[i];
1415 cout << "\tpt: " << vmu->Pt()
1416 << "\teta: " << vmu->Eta()
1417 << "\tphi: " << vmu->Phi()
1418 << endl;
1419 }
1420 cout << "electronIsoMVASelection :: electrson to veto " << endl;
1421 for( int i=0; i<electronsToVeto.size(); i++ ) {
1422 const mithep::Electron * vel = electronsToVeto[i];
1423 cout << "\tpt: " << vel->Pt()
1424 << "\teta: " << vel->Eta()
1425 << "\tphi: " << vel->Phi()
1426 << "\ttrk: " << vel->TrackerTrk()
1427 << endl;
1428 }
1429 }
1430
1431 bool failiso=false;
1432
1433 //
1434 // tmp iso rings
1435 //
1436 Double_t tmpChargedIso_DR0p0To0p1 = 0;
1437 Double_t tmpChargedIso_DR0p1To0p2 = 0;
1438 Double_t tmpChargedIso_DR0p2To0p3 = 0;
1439 Double_t tmpChargedIso_DR0p3To0p4 = 0;
1440 Double_t tmpChargedIso_DR0p4To0p5 = 0;
1441
1442 Double_t tmpGammaIso_DR0p0To0p1 = 0;
1443 Double_t tmpGammaIso_DR0p1To0p2 = 0;
1444 Double_t tmpGammaIso_DR0p2To0p3 = 0;
1445 Double_t tmpGammaIso_DR0p3To0p4 = 0;
1446 Double_t tmpGammaIso_DR0p4To0p5 = 0;
1447
1448
1449 Double_t tmpNeutralHadronIso_DR0p0To0p1 = 0;
1450 Double_t tmpNeutralHadronIso_DR0p1To0p2 = 0;
1451 Double_t tmpNeutralHadronIso_DR0p2To0p3 = 0;
1452 Double_t tmpNeutralHadronIso_DR0p3To0p4 = 0;
1453 Double_t tmpNeutralHadronIso_DR0p4To0p5 = 0;
1454
1455
1456
1457 //
1458 // final rings for the MVA
1459 //
1460 Double_t fChargedIso_DR0p0To0p1;
1461 Double_t fChargedIso_DR0p1To0p2;
1462 Double_t fChargedIso_DR0p2To0p3;
1463 Double_t fChargedIso_DR0p3To0p4;
1464 Double_t fChargedIso_DR0p4To0p5;
1465
1466 Double_t fGammaIso_DR0p0To0p1;
1467 Double_t fGammaIso_DR0p1To0p2;
1468 Double_t fGammaIso_DR0p2To0p3;
1469 Double_t fGammaIso_DR0p3To0p4;
1470 Double_t fGammaIso_DR0p4To0p5;
1471
1472 Double_t fNeutralHadronIso_DR0p0To0p1;
1473 Double_t fNeutralHadronIso_DR0p1To0p2;
1474 Double_t fNeutralHadronIso_DR0p2To0p3;
1475 Double_t fNeutralHadronIso_DR0p3To0p4;
1476 Double_t fNeutralHadronIso_DR0p4To0p5;
1477
1478
1479 //
1480 //Loop over PF Candidates
1481 //
1482 for(int k=0; k<fPFCandidates->GetEntries(); ++k) {
1483
1484 if( !(PFnoPUflag[k]) ) continue; // my PF no PU hack
1485
1486 const mithep::PFCandidate *pf = (mithep::PFCandidate*)((*fPFCandidates)[k]);
1487 Double_t deta = (ele->Eta() - pf->Eta());
1488 Double_t dphi = mithep::MathUtils::DeltaPhi(Double_t(ele->Phi()),Double_t(pf->Phi()));
1489 Double_t dr = mithep::MathUtils::DeltaR(ele->Phi(),ele->Eta(), pf->Phi(), pf->Eta());
1490 if (dr > 1.0) continue;
1491
1492 if(ctrl.debug) {
1493 cout << "pf :: type: " << pf->PFType() << "\tpt: " << pf->Pt();
1494 if( pf->HasTrackerTrk() ) cout << "\tdZ: " << pf->TrackerTrk()->DzCorrected(*vtx);
1495 cout << endl;
1496 }
1497
1498
1499 if ( (pf->HasTrackerTrk() && (pf->TrackerTrk() == ele->TrackerTrk())) ||
1500 (pf->HasGsfTrk() && (pf->GsfTrk() == ele->GsfTrk()))) continue;
1501
1502
1503 //
1504 // Lepton Footprint Removal
1505 //
1506 Bool_t IsLeptonFootprint = kFALSE;
1507 if (dr < 1.0) {
1508
1509
1510 //
1511 // Check for electrons
1512 //
1513
1514 for (Int_t q=0; q < electronsToVeto.size(); ++q) {
1515 const mithep::Electron *tmpele = electronsToVeto[q];
1516 double tmpdr = mithep::MathUtils::DeltaR(tmpele->Phi(),tmpele->Eta(), pf->Phi(), pf->Eta());
1517
1518 // 4l electron
1519 if( pf->HasTrackerTrk() ) {
1520 if( pf->TrackerTrk() == tmpele->TrackerTrk() ) {
1521 if( ctrl.debug) cout << "\tcharged tktrk, matches 4L ele ..." << endl;
1522 IsLeptonFootprint = kTRUE;
1523 }
1524 }
1525 if( pf->HasGsfTrk() ) {
1526 if( pf->GsfTrk() == tmpele->GsfTrk() ) {
1527 if( ctrl.debug) cout << "\tcharged gsftrk, matches 4L ele ..." << endl;
1528 IsLeptonFootprint = kTRUE;
1529 }
1530 }
1531 // PF charged
1532 if (pf->Charge() != 0 && fabs(tmpele->SCluster()->Eta()) >= 1.479 && tmpdr < 0.015) {
1533 if( ctrl.debug) cout << "\tcharged trk, dR matches 4L ele ..." << endl;
1534 IsLeptonFootprint = kTRUE;
1535 }
1536 // PF gamma
1537 if (abs(pf->PFType()) == PFCandidate::eGamma && fabs(tmpele->SCluster()->Eta()) >= 1.479
1538 && tmpdr < 0.08) {
1539 if( ctrl.debug) cout << "\tPF gamma, matches 4L ele ..." << endl;
1540 IsLeptonFootprint = kTRUE;
1541 }
1542 } // loop over electrons
1543
1544
1545 /* KH - comment for sync
1546 //
1547 // Check for muons
1548 //
1549 for (Int_t q=0; q < muonsToVeto.size(); ++q) {
1550 const mithep::Muon *tmpmu = muonsToVeto[q];
1551 // 4l muon
1552 if( pf->HasTrackerTrk() ) {
1553 if (pf->TrackerTrk() == tmpmu->TrackerTrk() ){
1554 if( ctrl.debug) cout << "\tmatches 4L mu ..." << endl;
1555 IsLeptonFootprint = kTRUE;
1556 }
1557 }
1558 // PF charged
1559 if (pf->Charge() != 0 && mithep::MathUtils::DeltaR(tmpmu->Phi(),tmpmu->Eta(), pf->Phi(), pf->Eta()) < 0.01) {
1560 if( ctrl.debug) cout << "\tcharged trk, dR matches 4L mu ..." << endl;
1561 IsLeptonFootprint = kTRUE;
1562 }
1563 } // loop over muons
1564 */
1565
1566 if (IsLeptonFootprint)
1567 continue;
1568
1569 //
1570 // Charged Iso Rings
1571 //
1572 if (pf->Charge() != 0 && (pf->HasTrackerTrk()||pf->HasGsfTrk()) ) {
1573
1574 // if( pf->HasGsfTrk() ) {
1575 // if (abs(pf->GsfTrk()->DzCorrected(vtx)) > 0.2) continue;
1576 // } else if( pf->HasTrackerTrk() ){
1577 // if (abs(pf->TrackerTrk()->DzCorrected(vtx)) > 0.2) continue;
1578 // }
1579
1580 // Veto any PFmuon, or PFEle
1581 if (abs(pf->PFType()) == PFCandidate::eElectron || abs(pf->PFType()) == PFCandidate::eMuon) continue;
1582
1583 // Footprint Veto
1584 if (fabs(ele->SCluster()->Eta()) > 1.479 && dr < 0.015) continue;
1585
1586 if( ctrl.debug) cout << "charged:: pt: " << pf->Pt()
1587 << "\ttype: " << pf->PFType()
1588 << "\ttrk: " << pf->TrackerTrk() << endl;
1589
1590 if (dr < 0.1) tmpChargedIso_DR0p0To0p1 += pf->Pt();
1591 if (dr >= 0.1 && dr < 0.2) tmpChargedIso_DR0p1To0p2 += pf->Pt();
1592 if (dr >= 0.2 && dr < 0.3) tmpChargedIso_DR0p2To0p3 += pf->Pt();
1593 if (dr >= 0.3 && dr < 0.4) tmpChargedIso_DR0p3To0p4 += pf->Pt();
1594 if (dr >= 0.4 && dr < 0.5) tmpChargedIso_DR0p4To0p5 += pf->Pt();
1595
1596 }
1597
1598 //
1599 // Gamma Iso Rings
1600 //
1601 else if (abs(pf->PFType()) == PFCandidate::eGamma) {
1602
1603 if (fabs(ele->SCluster()->Eta()) > 1.479 && dr < 0.08) continue;
1604
1605 if( ctrl.debug) cout << "gamma:: " << pf->Pt() << " "
1606 << dr << endl;
1607
1608 if (dr < 0.1) tmpGammaIso_DR0p0To0p1 += pf->Pt();
1609 if (dr >= 0.1 && dr < 0.2) tmpGammaIso_DR0p1To0p2 += pf->Pt();
1610 if (dr >= 0.2 && dr < 0.3) tmpGammaIso_DR0p2To0p3 += pf->Pt();
1611 if (dr >= 0.3 && dr < 0.4) tmpGammaIso_DR0p3To0p4 += pf->Pt();
1612 if (dr >= 0.4 && dr < 0.5) tmpGammaIso_DR0p4To0p5 += pf->Pt();
1613 }
1614
1615 //
1616 // Other Neutral Iso Rings
1617 //
1618 else {
1619 if( ctrl.debug) cout << "neutral:: " << pf->Pt() << " "
1620 << dr << endl;
1621 if (dr < 0.1) tmpNeutralHadronIso_DR0p0To0p1 += pf->Pt();
1622 if (dr >= 0.1 && dr < 0.2) tmpNeutralHadronIso_DR0p1To0p2 += pf->Pt();
1623 if (dr >= 0.2 && dr < 0.3) tmpNeutralHadronIso_DR0p2To0p3 += pf->Pt();
1624 if (dr >= 0.3 && dr < 0.4) tmpNeutralHadronIso_DR0p3To0p4 += pf->Pt();
1625 if (dr >= 0.4 && dr < 0.5) tmpNeutralHadronIso_DR0p4To0p5 += pf->Pt();
1626 }
1627
1628 }
1629
1630 }
1631
1632 fChargedIso_DR0p0To0p1 = fmin((tmpChargedIso_DR0p0To0p1)/ele->Pt(), 2.5);
1633 fChargedIso_DR0p1To0p2 = fmin((tmpChargedIso_DR0p1To0p2)/ele->Pt(), 2.5);
1634 fChargedIso_DR0p2To0p3 = fmin((tmpChargedIso_DR0p2To0p3)/ele->Pt(), 2.5);
1635 fChargedIso_DR0p3To0p4 = fmin((tmpChargedIso_DR0p3To0p4)/ele->Pt(), 2.5);
1636 fChargedIso_DR0p4To0p5 = fmin((tmpChargedIso_DR0p4To0p5)/ele->Pt(), 2.5);
1637
1638 if(ctrl.debug) {
1639 cout << "fChargedIso_DR0p0To0p1 : " << fChargedIso_DR0p0To0p1 << endl;
1640 cout << "fChargedIso_DR0p1To0p2 : " << fChargedIso_DR0p1To0p2 << endl;
1641 cout << "fChargedIso_DR0p2To0p3 : " << fChargedIso_DR0p2To0p3 << endl;
1642 cout << "fChargedIso_DR0p3To0p4 : " << fChargedIso_DR0p3To0p4 << endl;
1643 cout << "fChargedIso_DR0p4To0p5 : " << fChargedIso_DR0p4To0p5 << endl;
1644 }
1645
1646
1647 double rho = 0;
1648 if (!(isnan(fPUEnergyDensity->At(0)->Rho()) || isinf(fPUEnergyDensity->At(0)->Rho())))
1649 rho = fPUEnergyDensity->At(0)->Rho();
1650 // if (!(isnan(fPUEnergyDensity->At(0)->RhoLowEta()) || isinf(fPUEnergyDensity->At(0)->RhoLowEta())))
1651 // rho = fPUEnergyDensity->At(0)->RhoLowEta();
1652
1653 // WARNING!!!!
1654 // hardcode for sync ...
1655 EffectiveAreaVersion = eleT.kEleEAData2011;
1656 // WARNING!!!!
1657
1658 if( ctrl.debug) {
1659 cout << "RHO: " << rho << endl;
1660 cout << "eta: " << ele->SCluster()->Eta() << endl;
1661 cout << "target: " << EffectiveAreaVersion << endl;
1662 cout << "effA 0-1: " << eleT.ElectronEffectiveArea(eleT.kEleNeutralHadronIsoDR0p0To0p1,
1663 ele->SCluster()->Eta(),
1664 EffectiveAreaVersion)
1665 << endl;
1666 cout << "effA 1-2: " << eleT.ElectronEffectiveArea(eleT.kEleNeutralHadronIsoDR0p1To0p2,
1667 ele->SCluster()->Eta(),
1668 EffectiveAreaVersion)
1669 << endl;
1670 cout << "effA 2-3: " << eleT.ElectronEffectiveArea(eleT.kEleNeutralHadronIsoDR0p2To0p3,
1671 ele->SCluster()->Eta(),
1672 EffectiveAreaVersion)
1673 << endl;
1674 cout << "effA 3-4: " << eleT.ElectronEffectiveArea(eleT.kEleNeutralHadronIsoDR0p3To0p4,
1675 ele->SCluster()->Eta(),
1676 EffectiveAreaVersion)
1677 << endl;
1678 }
1679
1680 fGammaIso_DR0p0To0p1 = fmax(fmin((tmpGammaIso_DR0p0To0p1
1681 -rho*eleT.ElectronEffectiveArea(eleT.kEleGammaIsoDR0p0To0p1,
1682 ele->SCluster()->Eta(),
1683 EffectiveAreaVersion))/ele->Pt()
1684 ,2.5)
1685 ,0.0);
1686 fGammaIso_DR0p1To0p2 = fmax(fmin((tmpGammaIso_DR0p1To0p2
1687 -rho*eleT.ElectronEffectiveArea(eleT.kEleGammaIsoDR0p1To0p2,
1688 ele->SCluster()->Eta(),
1689 EffectiveAreaVersion))/ele->Pt()
1690 ,2.5)
1691 ,0.0);
1692 fGammaIso_DR0p2To0p3 = fmax(fmin((tmpGammaIso_DR0p2To0p3
1693 -rho*eleT.ElectronEffectiveArea(eleT.kEleGammaIsoDR0p2To0p3,
1694 ele->SCluster()->Eta()
1695 ,EffectiveAreaVersion))/ele->Pt()
1696 ,2.5)
1697 ,0.0);
1698 fGammaIso_DR0p3To0p4 = fmax(fmin((tmpGammaIso_DR0p3To0p4
1699 -rho*eleT.ElectronEffectiveArea(eleT.kEleGammaIsoDR0p3To0p4,
1700 ele->SCluster()->Eta(),
1701 EffectiveAreaVersion))/ele->Pt()
1702 ,2.5)
1703 ,0.0);
1704 fGammaIso_DR0p4To0p5 = fmax(fmin((tmpGammaIso_DR0p4To0p5
1705 -rho*eleT.ElectronEffectiveArea(eleT.kEleGammaIsoDR0p4To0p5,
1706 ele->SCluster()->Eta(),
1707 EffectiveAreaVersion))/ele->Pt()
1708 ,2.5)
1709 ,0.0);
1710
1711
1712 if( ctrl.debug) {
1713 cout << "fGammaIso_DR0p0To0p1: " << fGammaIso_DR0p0To0p1 << endl;
1714 cout << "fGammaIso_DR0p1To0p2: " << fGammaIso_DR0p1To0p2 << endl;
1715 cout << "fGammaIso_DR0p2To0p3: " << fGammaIso_DR0p2To0p3 << endl;
1716 cout << "fGammaIso_DR0p3To0p4: " << fGammaIso_DR0p3To0p4 << endl;
1717 cout << "fGammaIso_DR0p4To0p5: " << fGammaIso_DR0p4To0p5 << endl;
1718 }
1719
1720 fNeutralHadronIso_DR0p0To0p1 = fmax(fmin((tmpNeutralHadronIso_DR0p0To0p1
1721 -rho*eleT.ElectronEffectiveArea(eleT.kEleNeutralHadronIsoDR0p0To0p1,
1722 ele->SCluster()->Eta(),EffectiveAreaVersion))/ele->Pt()
1723 , 2.5)
1724 , 0.0);
1725 fNeutralHadronIso_DR0p1To0p2 = fmax(fmin((tmpNeutralHadronIso_DR0p1To0p2
1726 -rho*eleT.ElectronEffectiveArea(eleT.kEleNeutralHadronIsoDR0p1To0p2,
1727 ele->SCluster()->Eta(),EffectiveAreaVersion))/ele->Pt()
1728 , 2.5)
1729 , 0.0);
1730 fNeutralHadronIso_DR0p2To0p3 = fmax(fmin((tmpNeutralHadronIso_DR0p2To0p3
1731 -rho*eleT.ElectronEffectiveArea(eleT.kEleNeutralHadronIsoDR0p2To0p3,
1732 ele->SCluster()->Eta(),EffectiveAreaVersion))/ele->Pt()
1733 , 2.5)
1734 , 0.0);
1735 fNeutralHadronIso_DR0p3To0p4 = fmax(fmin((tmpNeutralHadronIso_DR0p3To0p4
1736 -rho*eleT.ElectronEffectiveArea(eleT.kEleNeutralHadronIsoDR0p3To0p4,
1737 ele->SCluster()->Eta(), EffectiveAreaVersion))/ele->Pt()
1738 , 2.5)
1739 , 0.0);
1740 fNeutralHadronIso_DR0p4To0p5 = fmax(fmin((tmpNeutralHadronIso_DR0p4To0p5
1741 -rho*eleT.ElectronEffectiveArea(eleT.kEleNeutralHadronIsoDR0p4To0p5,
1742 ele->SCluster()->Eta(), EffectiveAreaVersion))/ele->Pt()
1743 , 2.5)
1744 , 0.0);
1745
1746 if( ctrl.debug) {
1747 cout << "fNeutralHadronIso_DR0p0To0p1: " << fNeutralHadronIso_DR0p0To0p1 << endl;
1748 cout << "fNeutralHadronIso_DR0p1To0p2: " << fNeutralHadronIso_DR0p1To0p2 << endl;
1749 cout << "fNeutralHadronIso_DR0p2To0p3: " << fNeutralHadronIso_DR0p2To0p3 << endl;
1750 cout << "fNeutralHadronIso_DR0p3To0p4: " << fNeutralHadronIso_DR0p3To0p4 << endl;
1751 cout << "fNeutralHadronIso_DR0p4To0p5: " << fNeutralHadronIso_DR0p4To0p5 << endl;
1752 }
1753
1754 double mvaval = eleIsoMVA->MVAValue_IsoRings( ele->Pt(),
1755 ele->SCluster()->Eta(),
1756 fChargedIso_DR0p0To0p1,
1757 fChargedIso_DR0p1To0p2,
1758 fChargedIso_DR0p2To0p3,
1759 fChargedIso_DR0p3To0p4,
1760 fChargedIso_DR0p4To0p5,
1761 fGammaIso_DR0p0To0p1,
1762 fGammaIso_DR0p1To0p2,
1763 fGammaIso_DR0p2To0p3,
1764 fGammaIso_DR0p3To0p4,
1765 fGammaIso_DR0p4To0p5,
1766 fNeutralHadronIso_DR0p0To0p1,
1767 fNeutralHadronIso_DR0p1To0p2,
1768 fNeutralHadronIso_DR0p2To0p3,
1769 fNeutralHadronIso_DR0p3To0p4,
1770 fNeutralHadronIso_DR0p4To0p5,
1771 ctrl.debug);
1772
1773 SelectionStatus status;
1774 status.isoMVA = mvaval;
1775 bool pass = false;
1776
1777 Int_t subdet = 0;
1778 if (fabs(ele->SCluster()->Eta()) < 0.8) subdet = 0;
1779 else if (fabs(ele->SCluster()->Eta()) < 1.479) subdet = 1;
1780 else subdet = 2;
1781
1782 Int_t ptBin = 0;
1783 if (ele->Pt() >= 10.0) ptBin = 1;
1784
1785 Int_t MVABin = -1;
1786 if (subdet == 0 && ptBin == 0) MVABin = 0;
1787 if (subdet == 1 && ptBin == 0) MVABin = 1;
1788 if (subdet == 2 && ptBin == 0) MVABin = 2;
1789 if (subdet == 0 && ptBin == 1) MVABin = 3;
1790 if (subdet == 1 && ptBin == 1) MVABin = 4;
1791 if (subdet == 2 && ptBin == 1) MVABin = 5;
1792
1793 pass = false;
1794 if( MVABin == 0 && mvaval > ELECTRON_LOOSE_ISOMVA_CUT_OPT_BIN0 ) pass = true;
1795 if( MVABin == 1 && mvaval > ELECTRON_LOOSE_ISOMVA_CUT_OPT_BIN1 ) pass = true;
1796 if( MVABin == 2 && mvaval > ELECTRON_LOOSE_ISOMVA_CUT_OPT_BIN2 ) pass = true;
1797 if( MVABin == 3 && mvaval > ELECTRON_LOOSE_ISOMVA_CUT_OPT_BIN3 ) pass = true;
1798 if( MVABin == 4 && mvaval > ELECTRON_LOOSE_ISOMVA_CUT_OPT_BIN4 ) pass = true;
1799 if( MVABin == 5 && mvaval > ELECTRON_LOOSE_ISOMVA_CUT_OPT_BIN5 ) pass = true;
1800 // pass &= (fChargedIso_DR0p0To0p1 + fChargedIso_DR0p1To0p2 + fChargedIso_DR0p2To0p3 < 0.7);
1801 if( pass ) status.orStatus(SelectionStatus::LOOSEISO);
1802
1803 // pass = false;
1804 // if( MVABin == 0 && mvaval > ELECTRON_TIGHT_ISOMVA_CUT_BIN0 ) pass = true;
1805 // if( MVABin == 1 && mvaval > ELECTRON_TIGHT_ISOMVA_CUT_BIN1 ) pass = true;
1806 // if( MVABin == 2 && mvaval > ELECTRON_TIGHT_ISOMVA_CUT_BIN2 ) pass = true;
1807 // if( MVABin == 3 && mvaval > ELECTRON_TIGHT_ISOMVA_CUT_BIN3 ) pass = true;
1808 // if( MVABin == 4 && mvaval > ELECTRON_TIGHT_ISOMVA_CUT_BIN4 ) pass = true;
1809 // if( MVABin == 5 && mvaval > ELECTRON_TIGHT_ISOMVA_CUT_BIN5 ) pass = true;
1810 // if( pass ) status.orStatus(SelectionStatus::TIGHTISO);
1811
1812 if(ctrl.debug) cout << "returning status : " << hex << status.getStatus() << dec << endl;
1813 return status;
1814
1815 }
1816
1817
1818 //--------------------------------------------------------------------------------------------------
1819 SelectionStatus electronIsoMVASelection(ControlFlags &ctrl,
1820 const mithep::Electron * ele,
1821 const mithep::Vertex * vtx,
1822 const mithep::Array<mithep::PFCandidate> * fPFCandidates,
1823 float rho,
1824 //const mithep::Array<mithep::PileupEnergyDensity> * fPUEnergyDensity,
1825 mithep::ElectronTools::EElectronEffectiveAreaTarget EffectiveAreaVersion,
1826 vector<const mithep::Muon*> muonsToVeto,
1827 vector<const mithep::Electron*> electronsToVeto)
1828 //--------------------------------------------------------------------------------------------------
1829 // hacked version
1830 {
1831 if( ctrl.debug ) {
1832 cout << "================> hacked ele Iso MVA <======================" << endl;
1833 }
1834
1835 if( ctrl.debug ) {
1836 cout << "electronIsoMVASelection :: muons to veto " << endl;
1837 for( int i=0; i<muonsToVeto.size(); i++ ) {
1838 const mithep::Muon * vmu = muonsToVeto[i];
1839 cout << "\tpt: " << vmu->Pt()
1840 << "\teta: " << vmu->Eta()
1841 << "\tphi: " << vmu->Phi()
1842 << endl;
1843 }
1844 cout << "electronIsoMVASelection :: electrson to veto " << endl;
1845 for( int i=0; i<electronsToVeto.size(); i++ ) {
1846 const mithep::Electron * vel = electronsToVeto[i];
1847 cout << "\tpt: " << vel->Pt()
1848 << "\teta: " << vel->Eta()
1849 << "\tphi: " << vel->Phi()
1850 << "\ttrk: " << vel->TrackerTrk()
1851 << endl;
1852 }
1853 }
1854
1855 bool failiso=false;
1856
1857 //
1858 // tmp iso rings
1859 //
1860 Double_t tmpChargedIso_DR0p0To0p1 = 0;
1861 Double_t tmpChargedIso_DR0p1To0p2 = 0;
1862 Double_t tmpChargedIso_DR0p2To0p3 = 0;
1863 Double_t tmpChargedIso_DR0p3To0p4 = 0;
1864 Double_t tmpChargedIso_DR0p4To0p5 = 0;
1865
1866 Double_t tmpGammaIso_DR0p0To0p1 = 0;
1867 Double_t tmpGammaIso_DR0p1To0p2 = 0;
1868 Double_t tmpGammaIso_DR0p2To0p3 = 0;
1869 Double_t tmpGammaIso_DR0p3To0p4 = 0;
1870 Double_t tmpGammaIso_DR0p4To0p5 = 0;
1871
1872
1873 Double_t tmpNeutralHadronIso_DR0p0To0p1 = 0;
1874 Double_t tmpNeutralHadronIso_DR0p1To0p2 = 0;
1875 Double_t tmpNeutralHadronIso_DR0p2To0p3 = 0;
1876 Double_t tmpNeutralHadronIso_DR0p3To0p4 = 0;
1877 Double_t tmpNeutralHadronIso_DR0p4To0p5 = 0;
1878
1879
1880
1881 //
1882 // final rings for the MVA
1883 //
1884 Double_t fChargedIso_DR0p0To0p1;
1885 Double_t fChargedIso_DR0p1To0p2;
1886 Double_t fChargedIso_DR0p2To0p3;
1887 Double_t fChargedIso_DR0p3To0p4;
1888 Double_t fChargedIso_DR0p4To0p5;
1889
1890 Double_t fGammaIso_DR0p0To0p1;
1891 Double_t fGammaIso_DR0p1To0p2;
1892 Double_t fGammaIso_DR0p2To0p3;
1893 Double_t fGammaIso_DR0p3To0p4;
1894 Double_t fGammaIso_DR0p4To0p5;
1895
1896 Double_t fNeutralHadronIso_DR0p0To0p1;
1897 Double_t fNeutralHadronIso_DR0p1To0p2;
1898 Double_t fNeutralHadronIso_DR0p2To0p3;
1899 Double_t fNeutralHadronIso_DR0p3To0p4;
1900 Double_t fNeutralHadronIso_DR0p4To0p5;
1901
1902
1903 //
1904 //Loop over PF Candidates
1905 //
1906 for(int k=0; k<fPFCandidates->GetEntries(); ++k) {
1907
1908 if( !(PFnoPUflag[k]) ) continue; // my PF no PU hack
1909
1910 const mithep::PFCandidate *pf = (mithep::PFCandidate*)((*fPFCandidates)[k]);
1911 Double_t deta = (ele->Eta() - pf->Eta());
1912 Double_t dphi = mithep::MathUtils::DeltaPhi(Double_t(ele->Phi()),Double_t(pf->Phi()));
1913 Double_t dr = mithep::MathUtils::DeltaR(ele->Phi(),ele->Eta(), pf->Phi(), pf->Eta());
1914 if (dr > 1.0) continue;
1915
1916 if(ctrl.debug) {
1917 cout << "pf :: type: " << pf->PFType() << "\tpt: " << pf->Pt();
1918 if( pf->HasTrackerTrk() ) cout << "\tdZ: " << pf->TrackerTrk()->DzCorrected(*vtx);
1919 cout << endl;
1920 }
1921
1922
1923 if ( (pf->HasTrackerTrk() && (pf->TrackerTrk() == ele->TrackerTrk())) ||
1924 (pf->HasGsfTrk() && (pf->GsfTrk() == ele->GsfTrk()))) continue;
1925
1926
1927 //
1928 // Lepton Footprint Removal
1929 //
1930 Bool_t IsLeptonFootprint = kFALSE;
1931 if (dr < 1.0) {
1932
1933
1934 //
1935 // Check for electrons
1936 //
1937
1938 for (Int_t q=0; q < electronsToVeto.size(); ++q) {
1939 const mithep::Electron *tmpele = electronsToVeto[q];
1940 double tmpdr = mithep::MathUtils::DeltaR(tmpele->Phi(),tmpele->Eta(), pf->Phi(), pf->Eta());
1941
1942 // 4l electron
1943 if( pf->HasTrackerTrk() ) {
1944 if( pf->TrackerTrk() == tmpele->TrackerTrk() ) {
1945 if( ctrl.debug) cout << "\tcharged tktrk, matches 4L ele ..." << endl;
1946 IsLeptonFootprint = kTRUE;
1947 }
1948 }
1949 if( pf->HasGsfTrk() ) {
1950 if( pf->GsfTrk() == tmpele->GsfTrk() ) {
1951 if( ctrl.debug) cout << "\tcharged gsftrk, matches 4L ele ..." << endl;
1952 IsLeptonFootprint = kTRUE;
1953 }
1954 }
1955 // PF charged
1956 if (pf->Charge() != 0 && fabs(tmpele->SCluster()->Eta()) >= 1.479 && tmpdr < 0.015) {
1957 if( ctrl.debug) cout << "\tcharged trk, dR matches 4L ele ..." << endl;
1958 IsLeptonFootprint = kTRUE;
1959 }
1960 // PF gamma
1961 if (abs(pf->PFType()) == PFCandidate::eGamma && fabs(tmpele->SCluster()->Eta()) >= 1.479
1962 && tmpdr < 0.08) {
1963 if( ctrl.debug) cout << "\tPF gamma, matches 4L ele ..." << endl;
1964 IsLeptonFootprint = kTRUE;
1965 }
1966 } // loop over electrons
1967
1968
1969 /* KH - comment for sync
1970 //
1971 // Check for muons
1972 //
1973 for (Int_t q=0; q < muonsToVeto.size(); ++q) {
1974 const mithep::Muon *tmpmu = muonsToVeto[q];
1975 // 4l muon
1976 if( pf->HasTrackerTrk() ) {
1977 if (pf->TrackerTrk() == tmpmu->TrackerTrk() ){
1978 if( ctrl.debug) cout << "\tmatches 4L mu ..." << endl;
1979 IsLeptonFootprint = kTRUE;
1980 }
1981 }
1982 // PF charged
1983 if (pf->Charge() != 0 && mithep::MathUtils::DeltaR(tmpmu->Phi(),tmpmu->Eta(), pf->Phi(), pf->Eta()) < 0.01) {
1984 if( ctrl.debug) cout << "\tcharged trk, dR matches 4L mu ..." << endl;
1985 IsLeptonFootprint = kTRUE;
1986 }
1987 } // loop over muons
1988 */
1989
1990 if (IsLeptonFootprint)
1991 continue;
1992
1993 //
1994 // Charged Iso Rings
1995 //
1996 if (pf->Charge() != 0 && (pf->HasTrackerTrk()||pf->HasGsfTrk()) ) {
1997
1998 // if( pf->HasGsfTrk() ) {
1999 // if (abs(pf->GsfTrk()->DzCorrected(vtx)) > 0.2) continue;
2000 // } else if( pf->HasTrackerTrk() ){
2001 // if (abs(pf->TrackerTrk()->DzCorrected(vtx)) > 0.2) continue;
2002 // }
2003
2004 // Veto any PFmuon, or PFEle
2005 if (abs(pf->PFType()) == PFCandidate::eElectron || abs(pf->PFType()) == PFCandidate::eMuon) continue;
2006
2007 // Footprint Veto
2008 if (fabs(ele->SCluster()->Eta()) > 1.479 && dr < 0.015) continue;
2009
2010 if( ctrl.debug) cout << "charged:: pt: " << pf->Pt()
2011 << "\ttype: " << pf->PFType()
2012 << "\ttrk: " << pf->TrackerTrk() << endl;
2013
2014 if (dr < 0.1) tmpChargedIso_DR0p0To0p1 += pf->Pt();
2015 if (dr >= 0.1 && dr < 0.2) tmpChargedIso_DR0p1To0p2 += pf->Pt();
2016 if (dr >= 0.2 && dr < 0.3) tmpChargedIso_DR0p2To0p3 += pf->Pt();
2017 if (dr >= 0.3 && dr < 0.4) tmpChargedIso_DR0p3To0p4 += pf->Pt();
2018 if (dr >= 0.4 && dr < 0.5) tmpChargedIso_DR0p4To0p5 += pf->Pt();
2019
2020 }
2021
2022 //
2023 // Gamma Iso Rings
2024 //
2025 else if (abs(pf->PFType()) == PFCandidate::eGamma) {
2026
2027 if (fabs(ele->SCluster()->Eta()) > 1.479 && dr < 0.08) continue;
2028
2029 if( ctrl.debug) cout << "gamma:: " << pf->Pt() << " "
2030 << dr << endl;
2031
2032 if (dr < 0.1) tmpGammaIso_DR0p0To0p1 += pf->Pt();
2033 if (dr >= 0.1 && dr < 0.2) tmpGammaIso_DR0p1To0p2 += pf->Pt();
2034 if (dr >= 0.2 && dr < 0.3) tmpGammaIso_DR0p2To0p3 += pf->Pt();
2035 if (dr >= 0.3 && dr < 0.4) tmpGammaIso_DR0p3To0p4 += pf->Pt();
2036 if (dr >= 0.4 && dr < 0.5) tmpGammaIso_DR0p4To0p5 += pf->Pt();
2037 }
2038
2039 //
2040 // Other Neutral Iso Rings
2041 //
2042 else {
2043 if( ctrl.debug) cout << "neutral:: " << pf->Pt() << " "
2044 << dr << endl;
2045 if (dr < 0.1) tmpNeutralHadronIso_DR0p0To0p1 += pf->Pt();
2046 if (dr >= 0.1 && dr < 0.2) tmpNeutralHadronIso_DR0p1To0p2 += pf->Pt();
2047 if (dr >= 0.2 && dr < 0.3) tmpNeutralHadronIso_DR0p2To0p3 += pf->Pt();
2048 if (dr >= 0.3 && dr < 0.4) tmpNeutralHadronIso_DR0p3To0p4 += pf->Pt();
2049 if (dr >= 0.4 && dr < 0.5) tmpNeutralHadronIso_DR0p4To0p5 += pf->Pt();
2050 }
2051
2052 }
2053
2054 }
2055
2056 fChargedIso_DR0p0To0p1 = fmin((tmpChargedIso_DR0p0To0p1)/ele->Pt(), 2.5);
2057 fChargedIso_DR0p1To0p2 = fmin((tmpChargedIso_DR0p1To0p2)/ele->Pt(), 2.5);
2058 fChargedIso_DR0p2To0p3 = fmin((tmpChargedIso_DR0p2To0p3)/ele->Pt(), 2.5);
2059 fChargedIso_DR0p3To0p4 = fmin((tmpChargedIso_DR0p3To0p4)/ele->Pt(), 2.5);
2060 fChargedIso_DR0p4To0p5 = fmin((tmpChargedIso_DR0p4To0p5)/ele->Pt(), 2.5);
2061
2062 if(ctrl.debug) {
2063 cout << "fChargedIso_DR0p0To0p1 : " << fChargedIso_DR0p0To0p1 << endl;
2064 cout << "fChargedIso_DR0p1To0p2 : " << fChargedIso_DR0p1To0p2 << endl;
2065 cout << "fChargedIso_DR0p2To0p3 : " << fChargedIso_DR0p2To0p3 << endl;
2066 cout << "fChargedIso_DR0p3To0p4 : " << fChargedIso_DR0p3To0p4 << endl;
2067 cout << "fChargedIso_DR0p4To0p5 : " << fChargedIso_DR0p4To0p5 << endl;
2068 }
2069
2070
2071 // rho=0;
2072 // double rho = 0;
2073 // if (!(isnan(fPUEnergyDensity->At(0)->Rho()) || isinf(fPUEnergyDensity->At(0)->Rho())))
2074 // rho = fPUEnergyDensity->At(0)->Rho();
2075 // if (!(isnan(fPUEnergyDensity->At(0)->RhoLowEta()) || isinf(fPUEnergyDensity->At(0)->RhoLowEta())))
2076 // rho = fPUEnergyDensity->At(0)->RhoLowEta();
2077
2078 // WARNING!!!!
2079 // hardcode for sync ...
2080 EffectiveAreaVersion = eleT.kEleEAData2011;
2081 // WARNING!!!!
2082
2083 if( ctrl.debug) {
2084 cout << "RHO: " << rho << endl;
2085 cout << "eta: " << ele->SCluster()->Eta() << endl;
2086 cout << "target: " << EffectiveAreaVersion << endl;
2087 cout << "effA 0-1: " << eleT.ElectronEffectiveArea(eleT.kEleNeutralHadronIsoDR0p0To0p1,
2088 ele->SCluster()->Eta(),
2089 EffectiveAreaVersion)
2090 << endl;
2091 cout << "effA 1-2: " << eleT.ElectronEffectiveArea(eleT.kEleNeutralHadronIsoDR0p1To0p2,
2092 ele->SCluster()->Eta(),
2093 EffectiveAreaVersion)
2094 << endl;
2095 cout << "effA 2-3: " << eleT.ElectronEffectiveArea(eleT.kEleNeutralHadronIsoDR0p2To0p3,
2096 ele->SCluster()->Eta(),
2097 EffectiveAreaVersion)
2098 << endl;
2099 cout << "effA 3-4: " << eleT.ElectronEffectiveArea(eleT.kEleNeutralHadronIsoDR0p3To0p4,
2100 ele->SCluster()->Eta(),
2101 EffectiveAreaVersion)
2102 << endl;
2103 }
2104
2105 fGammaIso_DR0p0To0p1 = fmax(fmin((tmpGammaIso_DR0p0To0p1
2106 -rho*eleT.ElectronEffectiveArea(eleT.kEleGammaIsoDR0p0To0p1,
2107 ele->SCluster()->Eta(),
2108 EffectiveAreaVersion))/ele->Pt()
2109 ,2.5)
2110 ,0.0);
2111 fGammaIso_DR0p1To0p2 = fmax(fmin((tmpGammaIso_DR0p1To0p2
2112 -rho*eleT.ElectronEffectiveArea(eleT.kEleGammaIsoDR0p1To0p2,
2113 ele->SCluster()->Eta(),
2114 EffectiveAreaVersion))/ele->Pt()
2115 ,2.5)
2116 ,0.0);
2117 fGammaIso_DR0p2To0p3 = fmax(fmin((tmpGammaIso_DR0p2To0p3
2118 -rho*eleT.ElectronEffectiveArea(eleT.kEleGammaIsoDR0p2To0p3,
2119 ele->SCluster()->Eta()
2120 ,EffectiveAreaVersion))/ele->Pt()
2121 ,2.5)
2122 ,0.0);
2123 fGammaIso_DR0p3To0p4 = fmax(fmin((tmpGammaIso_DR0p3To0p4
2124 -rho*eleT.ElectronEffectiveArea(eleT.kEleGammaIsoDR0p3To0p4,
2125 ele->SCluster()->Eta(),
2126 EffectiveAreaVersion))/ele->Pt()
2127 ,2.5)
2128 ,0.0);
2129 fGammaIso_DR0p4To0p5 = fmax(fmin((tmpGammaIso_DR0p4To0p5
2130 -rho*eleT.ElectronEffectiveArea(eleT.kEleGammaIsoDR0p4To0p5,
2131 ele->SCluster()->Eta(),
2132 EffectiveAreaVersion))/ele->Pt()
2133 ,2.5)
2134 ,0.0);
2135
2136
2137 if( ctrl.debug) {
2138 cout << "fGammaIso_DR0p0To0p1: " << fGammaIso_DR0p0To0p1 << endl;
2139 cout << "fGammaIso_DR0p1To0p2: " << fGammaIso_DR0p1To0p2 << endl;
2140 cout << "fGammaIso_DR0p2To0p3: " << fGammaIso_DR0p2To0p3 << endl;
2141 cout << "fGammaIso_DR0p3To0p4: " << fGammaIso_DR0p3To0p4 << endl;
2142 cout << "fGammaIso_DR0p4To0p5: " << fGammaIso_DR0p4To0p5 << endl;
2143 }
2144
2145 fNeutralHadronIso_DR0p0To0p1 = fmax(fmin((tmpNeutralHadronIso_DR0p0To0p1
2146 -rho*eleT.ElectronEffectiveArea(eleT.kEleNeutralHadronIsoDR0p0To0p1,
2147 ele->SCluster()->Eta(),EffectiveAreaVersion))/ele->Pt()
2148 , 2.5)
2149 , 0.0);
2150 fNeutralHadronIso_DR0p1To0p2 = fmax(fmin((tmpNeutralHadronIso_DR0p1To0p2
2151 -rho*eleT.ElectronEffectiveArea(eleT.kEleNeutralHadronIsoDR0p1To0p2,
2152 ele->SCluster()->Eta(),EffectiveAreaVersion))/ele->Pt()
2153 , 2.5)
2154 , 0.0);
2155 fNeutralHadronIso_DR0p2To0p3 = fmax(fmin((tmpNeutralHadronIso_DR0p2To0p3
2156 -rho*eleT.ElectronEffectiveArea(eleT.kEleNeutralHadronIsoDR0p2To0p3,
2157 ele->SCluster()->Eta(),EffectiveAreaVersion))/ele->Pt()
2158 , 2.5)
2159 , 0.0);
2160 fNeutralHadronIso_DR0p3To0p4 = fmax(fmin((tmpNeutralHadronIso_DR0p3To0p4
2161 -rho*eleT.ElectronEffectiveArea(eleT.kEleNeutralHadronIsoDR0p3To0p4,
2162 ele->SCluster()->Eta(), EffectiveAreaVersion))/ele->Pt()
2163 , 2.5)
2164 , 0.0);
2165 fNeutralHadronIso_DR0p4To0p5 = fmax(fmin((tmpNeutralHadronIso_DR0p4To0p5
2166 -rho*eleT.ElectronEffectiveArea(eleT.kEleNeutralHadronIsoDR0p4To0p5,
2167 ele->SCluster()->Eta(), EffectiveAreaVersion))/ele->Pt()
2168 , 2.5)
2169 , 0.0);
2170
2171 if( ctrl.debug) {
2172 cout << "fNeutralHadronIso_DR0p0To0p1: " << fNeutralHadronIso_DR0p0To0p1 << endl;
2173 cout << "fNeutralHadronIso_DR0p1To0p2: " << fNeutralHadronIso_DR0p1To0p2 << endl;
2174 cout << "fNeutralHadronIso_DR0p2To0p3: " << fNeutralHadronIso_DR0p2To0p3 << endl;
2175 cout << "fNeutralHadronIso_DR0p3To0p4: " << fNeutralHadronIso_DR0p3To0p4 << endl;
2176 cout << "fNeutralHadronIso_DR0p4To0p5: " << fNeutralHadronIso_DR0p4To0p5 << endl;
2177 }
2178
2179 double mvaval = eleIsoMVA->MVAValue_IsoRings( ele->Pt(),
2180 ele->SCluster()->Eta(),
2181 fChargedIso_DR0p0To0p1,
2182 fChargedIso_DR0p1To0p2,
2183 fChargedIso_DR0p2To0p3,
2184 fChargedIso_DR0p3To0p4,
2185 fChargedIso_DR0p4To0p5,
2186 fGammaIso_DR0p0To0p1,
2187 fGammaIso_DR0p1To0p2,
2188 fGammaIso_DR0p2To0p3,
2189 fGammaIso_DR0p3To0p4,
2190 fGammaIso_DR0p4To0p5,
2191 fNeutralHadronIso_DR0p0To0p1,
2192 fNeutralHadronIso_DR0p1To0p2,
2193 fNeutralHadronIso_DR0p2To0p3,
2194 fNeutralHadronIso_DR0p3To0p4,
2195 fNeutralHadronIso_DR0p4To0p5,
2196 ctrl.debug);
2197
2198 SelectionStatus status;
2199 status.isoMVA = mvaval;
2200 bool pass = false;
2201
2202 Int_t subdet = 0;
2203 if (fabs(ele->SCluster()->Eta()) < 0.8) subdet = 0;
2204 else if (fabs(ele->SCluster()->Eta()) < 1.479) subdet = 1;
2205 else subdet = 2;
2206
2207 Int_t ptBin = 0;
2208 if (ele->Pt() >= 10.0) ptBin = 1;
2209
2210 Int_t MVABin = -1;
2211 if (subdet == 0 && ptBin == 0) MVABin = 0;
2212 if (subdet == 1 && ptBin == 0) MVABin = 1;
2213 if (subdet == 2 && ptBin == 0) MVABin = 2;
2214 if (subdet == 0 && ptBin == 1) MVABin = 3;
2215 if (subdet == 1 && ptBin == 1) MVABin = 4;
2216 if (subdet == 2 && ptBin == 1) MVABin = 5;
2217
2218 pass = false;
2219 if( MVABin == 0 && mvaval > ELECTRON_LOOSE_ISOMVA_CUT_BIN0 ) pass = true;
2220 if( MVABin == 1 && mvaval > ELECTRON_LOOSE_ISOMVA_CUT_BIN1 ) pass = true;
2221 if( MVABin == 2 && mvaval > ELECTRON_LOOSE_ISOMVA_CUT_BIN2 ) pass = true;
2222 if( MVABin == 3 && mvaval > ELECTRON_LOOSE_ISOMVA_CUT_BIN3 ) pass = true;
2223 if( MVABin == 4 && mvaval > ELECTRON_LOOSE_ISOMVA_CUT_BIN4 ) pass = true;
2224 if( MVABin == 5 && mvaval > ELECTRON_LOOSE_ISOMVA_CUT_BIN5 ) pass = true;
2225 if( pass ) status.orStatus(SelectionStatus::LOOSEISO);
2226
2227 // pass = false;
2228 // if( MVABin == 0 && mvaval > ELECTRON_TIGHT_ISOMVA_CUT_BIN0 ) pass = true;
2229 // if( MVABin == 1 && mvaval > ELECTRON_TIGHT_ISOMVA_CUT_BIN1 ) pass = true;
2230 // if( MVABin == 2 && mvaval > ELECTRON_TIGHT_ISOMVA_CUT_BIN2 ) pass = true;
2231 // if( MVABin == 3 && mvaval > ELECTRON_TIGHT_ISOMVA_CUT_BIN3 ) pass = true;
2232 // if( MVABin == 4 && mvaval > ELECTRON_TIGHT_ISOMVA_CUT_BIN4 ) pass = true;
2233 // if( MVABin == 5 && mvaval > ELECTRON_TIGHT_ISOMVA_CUT_BIN5 ) pass = true;
2234 // if( pass ) status.orStatus(SelectionStatus::TIGHTISO);
2235
2236 if(ctrl.debug) cout << "returning status : " << hex << status.getStatus() << dec << endl;
2237 return status;
2238
2239 }
2240
2241
2242 //--------------------------------------------------------------------------------------------------
2243 void initElectronIsoMVA() {
2244 //--------------------------------------------------------------------------------------------------
2245 eleIsoMVA = new mithep::ElectronIDMVA();
2246 vector<string> weightFiles;
2247 weightFiles.push_back("../MitPhysics/data/ElectronMVAWeights/ElectronIso_BDTG_V0_BarrelPt5To10.weights.xml");
2248 weightFiles.push_back("../MitPhysics/data/ElectronMVAWeights/ElectronIso_BDTG_V0_EndcapPt5To10.weights.xml");
2249 weightFiles.push_back("../MitPhysics/data/ElectronMVAWeights/ElectronIso_BDTG_V0_BarrelPt10ToInf.weights.xml");
2250 weightFiles.push_back("../MitPhysics/data/ElectronMVAWeights/ElectronIso_BDTG_V0_EndcapPt10ToInf.weights.xml");
2251 eleIsoMVA->Initialize( "ElectronIsoMVA",
2252 mithep::ElectronIDMVA::kIsoRingsV0,
2253 kTRUE, weightFiles);
2254 }
2255
2256
2257
2258
2259 //--------------------------------------------------------------------------------------------------
2260 float electronPFIso04(ControlFlags &ctrl,
2261 const mithep::Electron * ele,
2262 const mithep::Vertex * vtx,
2263 const mithep::Array<mithep::PFCandidate> * fPFCandidates,
2264 const mithep::Array<mithep::PileupEnergyDensity> * fPUEnergyDensity,
2265 mithep::ElectronTools::EElectronEffectiveAreaTarget EffectiveAreaVersion,
2266 vector<const mithep::Muon*> muonsToVeto,
2267 vector<const mithep::Electron*> electronsToVeto)
2268 //--------------------------------------------------------------------------------------------------
2269 {
2270 /*
2271 if( ctrl.debug ) {
2272 cout << "electronIsoMVASelection :: muons to veto " << endl;
2273 for( int i=0; i<muonsToVeto.size(); i++ ) {
2274 const mithep::Muon * vmu = muonsToVeto[i];
2275 cout << "\tpt: " << vmu->Pt()
2276 << "\teta: " << vmu->Eta()
2277 << "\tphi: " << vmu->Phi()
2278 << endl;
2279 }
2280 cout << "electronIsoMVASelection :: electrons to veto " << endl;
2281 for( int i=0; i<electronsToVeto.size(); i++ ) {
2282 const mithep::Electron * vel = electronsToVeto[i];
2283 cout << "\tpt: " << vel->Pt()
2284 << "\teta: " << vel->Eta()
2285 << "\tphi: " << vel->Phi()
2286 << "\ttrk: " << vel->TrackerTrk()
2287 << endl;
2288 }
2289 }
2290 */
2291
2292 //
2293 // final iso
2294 //
2295 Double_t fChargedIso = 0.0;
2296 Double_t fGammaIso = 0.0;
2297 Double_t fNeutralHadronIso = 0.0;
2298
2299
2300 //
2301 //Loop over PF Candidates
2302 //
2303 for(int k=0; k<fPFCandidates->GetEntries(); ++k) {
2304
2305
2306 const mithep::PFCandidate *pf = (mithep::PFCandidate*)((*fPFCandidates)[k]);
2307 Double_t deta = (ele->Eta() - pf->Eta());
2308 Double_t dphi = mithep::MathUtils::DeltaPhi(Double_t(ele->Phi()),Double_t(pf->Phi()));
2309 Double_t dr = mithep::MathUtils::DeltaR(ele->Phi(),ele->Eta(), pf->Phi(), pf->Eta());
2310
2311 if (dr > 0.4) continue;
2312 if( !(PFnoPUflag[k]) ) continue; // my PF no PU hack
2313
2314 if(ctrl.debug) {
2315 cout << "pf :: type: " << pf->PFType() << "\tpt: " << pf->Pt() << "\tdR: " << dr;
2316 if( pf->HasTrackerTrk() ) cout << "\tdZ: " << pf->TrackerTrk()->DzCorrected(*vtx)
2317 << "\ttrk: " << pf->HasTrackerTrk()
2318 << "\tgsf: " << pf->HasGsfTrk();
2319
2320 cout << endl;
2321 }
2322
2323
2324 //
2325 // sync : I don't think theyre doing this ...
2326 //
2327 // if ( (pf->HasTrackerTrk() && (pf->TrackerTrk() == ele->TrackerTrk())) ||
2328 // (pf->HasGsfTrk() && (pf->GsfTrk() == ele->GsfTrk()))) {
2329 // if( ctrl.debug ) cout << "\tskipping, matches to the electron ..." << endl;
2330 // continue;
2331 // }
2332
2333
2334 //
2335 // Lepton Footprint Removal
2336 //
2337 Bool_t IsLeptonFootprint = kFALSE;
2338 if (dr < 1.0) {
2339
2340 //
2341 // Check for electrons
2342 //
2343 for (Int_t q=0; q < electronsToVeto.size(); ++q) {
2344 const mithep::Electron *tmpele = electronsToVeto[q];
2345 /*
2346 // 4l electron
2347 if( pf->HasTrackerTrk() ) {
2348 if( pf->TrackerTrk() == tmpele->TrackerTrk() ) {
2349 if( ctrl.debug) cout << "\tcharged tktrk, matches 4L ele ..." << endl;
2350 IsLeptonFootprint = kTRUE;
2351 }
2352 }
2353 if( pf->HasGsfTrk() ) {
2354 if( pf->GsfTrk() == tmpele->GsfTrk() ) {
2355 if( ctrl.debug) cout << "\tcharged gsftrk, matches 4L ele ..." << endl;
2356 IsLeptonFootprint = kTRUE;
2357 }
2358 }
2359 */
2360 // PF charged
2361 if (pf->Charge() != 0 && fabs(tmpele->SCluster()->Eta()) > 1.479
2362 && mithep::MathUtils::DeltaR(tmpele->Phi(),tmpele->Eta(), pf->Phi(), pf->Eta()) < 0.015) {
2363 if( ctrl.debug) cout << "\tcharged trk, dR matches 4L ele ..." << endl;
2364 IsLeptonFootprint = kTRUE;
2365 }
2366 // PF gamma
2367 if (abs(pf->PFType()) == PFCandidate::eGamma && fabs(tmpele->SCluster()->Eta()) > 1.479
2368 && mithep::MathUtils::DeltaR(tmpele->Phi(),tmpele->Eta(), pf->Phi(), pf->Eta()) < 0.08) {
2369 if( ctrl.debug) cout << "\tPF gamma, matches 4L ele ..." << endl;
2370 IsLeptonFootprint = kTRUE;
2371 }
2372 } // loop over electrons
2373
2374 /* KH - comment for sync
2375 //
2376 // Check for muons
2377 //
2378 for (Int_t q=0; q < muonsToVeto.size(); ++q) {
2379 const mithep::Muon *tmpmu = muonsToVeto[q];
2380 // 4l muon
2381 if( pf->HasTrackerTrk() ) {
2382 if (pf->TrackerTrk() == tmpmu->TrackerTrk() ){
2383 if( ctrl.debug) cout << "\tmatches 4L mu ..." << endl;
2384 IsLeptonFootprint = kTRUE;
2385 }
2386 }
2387 // PF charged
2388 if (pf->Charge() != 0 && mithep::MathUtils::DeltaR(tmpmu->Phi(),tmpmu->Eta(), pf->Phi(), pf->Eta()) < 0.01) {
2389 if( ctrl.debug) cout << "\tcharged trk, dR matches 4L mu ..." << endl;
2390 IsLeptonFootprint = kTRUE;
2391 }
2392 } // loop over muons
2393 */
2394
2395 if (IsLeptonFootprint)
2396 continue;
2397
2398 //
2399 // Charged Iso
2400 //
2401 if (pf->Charge() != 0 && (pf->HasTrackerTrk()||pf->HasGsfTrk()) ) {
2402
2403 // if( pf->HasTrackerTrk() )
2404 // if (abs(pf->TrackerTrk()->DzCorrected(vtx)) > 0.2) continue;
2405 // if( pf->HasGsfTrk() )
2406 // if (abs(pf->GsfTrk()->DzCorrected(vtx)) > 0.2) continue;
2407
2408 // Veto any PFmuon, or PFEle
2409 if (abs(pf->PFType()) == PFCandidate::eElectron || abs(pf->PFType()) == PFCandidate::eMuon) {
2410 if( ctrl.debug ) cout << "\t skipping, pf is and ele or mu .." <<endl;
2411 continue;
2412 }
2413
2414 // Footprint Veto
2415 if (fabs(ele->SCluster()->Eta()) > 1.479 && dr < 0.015) continue;
2416
2417 if( ctrl.debug) cout << "charged:: pt: " << pf->Pt()
2418 << "\ttype: " << pf->PFType()
2419 << "\ttrk: " << pf->TrackerTrk() << endl;
2420
2421 fChargedIso += pf->Pt();
2422 }
2423
2424 //
2425 // Gamma Iso
2426 //
2427 else if (abs(pf->PFType()) == PFCandidate::eGamma) {
2428
2429 if (fabs(ele->SCluster()->Eta()) > 1.479) {
2430 if (mithep::MathUtils::DeltaR(ele->Phi(),ele->Eta(), pf->Phi(), pf->Eta()) < 0.08) continue;
2431 }
2432 if( ctrl.debug) cout << "gamma:: " << pf->Pt() << " "
2433 << dr << endl;
2434 // KH, add to sync
2435 // if( pf->Pt() > 0.5 )
2436 fGammaIso += pf->Pt();
2437 }
2438
2439 //
2440 // Neutral Iso
2441 //
2442 else {
2443 if( ctrl.debug) cout << "neutral:: " << pf->Pt() << " "
2444 << dr << endl;
2445 // KH, add to sync
2446 // if( pf->Pt() > 0.5 )
2447 fNeutralHadronIso += pf->Pt();
2448 }
2449
2450 }
2451
2452 }
2453
2454
2455 double rho=0;
2456 if( (EffectiveAreaVersion == mithep::ElectronTools::kEleEAFall11MC) ||
2457 (EffectiveAreaVersion == mithep::ElectronTools::kEleEAData2011) ) {
2458 if (!(isnan(fPUEnergyDensity->At(0)->RhoKt6PFJetsForIso25()) ||
2459 isinf(fPUEnergyDensity->At(0)->RhoKt6PFJetsForIso25())))
2460 rho = fPUEnergyDensity->At(0)->RhoKt6PFJetsForIso25();
2461 // !!!!!!!!!!!!! TMP HACK FOR SYNC !!!!!!!!!!!!!!!!!!!!!
2462 EffectiveAreaVersion = mithep::ElectronTools::kEleEAData2011;
2463 // !!!!!!!!!!!!! TMP HACK FOR SYNC !!!!!!!!!!!!!!!!!!!!!
2464 } else {
2465 if (!(isnan(fPUEnergyDensity->At(0)->RhoKt6PFJets()) ||
2466 isinf(fPUEnergyDensity->At(0)->RhoKt6PFJets())))
2467 rho = fPUEnergyDensity->At(0)->RhoKt6PFJets();
2468 // !!!!!!!!!!!!! TMP HACK FOR SYNC !!!!!!!!!!!!!!!!!!!!!
2469 EffectiveAreaVersion = mithep::ElectronTools::kEleEAData2012;
2470 // !!!!!!!!!!!!! TMP HACK FOR SYNC !!!!!!!!!!!!!!!!!!!!!
2471 }
2472 if(ctrl.debug) cout << "rho: " << rho << endl;
2473
2474 double pfIso = fChargedIso + fmax(0.0,(fGammaIso + fNeutralHadronIso
2475 -rho*eleT.ElectronEffectiveArea(eleT.kEleGammaAndNeutralHadronIso04,
2476 ele->Eta(),EffectiveAreaVersion)));
2477
2478
2479 gChargedIso = fChargedIso;
2480 gGammaIso = fGammaIso;
2481 gNeutralIso = fNeutralHadronIso;
2482 return pfIso;
2483 }
2484
2485
2486
2487 //--------------------------------------------------------------------------------------------------
2488 // hacked version
2489 float electronPFIso04(ControlFlags &ctrl,
2490 const mithep::Electron * ele,
2491 const mithep::Vertex * vtx,
2492 const mithep::Array<mithep::PFCandidate> * fPFCandidates,
2493 float rho,
2494 mithep::ElectronTools::EElectronEffectiveAreaTarget EffectiveAreaVersion,
2495 vector<const mithep::Muon*> muonsToVeto,
2496 vector<const mithep::Electron*> electronsToVeto)
2497 //--------------------------------------------------------------------------------------------------
2498 {
2499
2500 if( ctrl.debug ) {
2501 cout << "electronIsoMVASelection :: muons to veto " << endl;
2502 for( int i=0; i<muonsToVeto.size(); i++ ) {
2503 const mithep::Muon * vmu = muonsToVeto[i];
2504 cout << "\tpt: " << vmu->Pt()
2505 << "\teta: " << vmu->Eta()
2506 << "\tphi: " << vmu->Phi()
2507 << endl;
2508 }
2509 cout << "electronIsoMVASelection :: electrons to veto " << endl;
2510 for( int i=0; i<electronsToVeto.size(); i++ ) {
2511 const mithep::Electron * vel = electronsToVeto[i];
2512 cout << "\tpt: " << vel->Pt()
2513 << "\teta: " << vel->Eta()
2514 << "\tphi: " << vel->Phi()
2515 << "\ttrk: " << vel->TrackerTrk()
2516 << endl;
2517 }
2518 }
2519
2520
2521 //
2522 // final iso
2523 //
2524 Double_t fChargedIso = 0.0;
2525 Double_t fGammaIso = 0.0;
2526 Double_t fNeutralHadronIso = 0.0;
2527
2528
2529 //
2530 //Loop over PF Candidates
2531 //
2532 for(int k=0; k<fPFCandidates->GetEntries(); ++k) {
2533
2534
2535 const mithep::PFCandidate *pf = (mithep::PFCandidate*)((*fPFCandidates)[k]);
2536 Double_t deta = (ele->Eta() - pf->Eta());
2537 Double_t dphi = mithep::MathUtils::DeltaPhi(Double_t(ele->Phi()),Double_t(pf->Phi()));
2538 Double_t dr = mithep::MathUtils::DeltaR(ele->Phi(),ele->Eta(), pf->Phi(), pf->Eta());
2539
2540 if (dr > 0.4) continue;
2541 if( !(PFnoPUflag[k]) ) continue; // my PF no PU hack
2542
2543 if(ctrl.debug) {
2544 cout << "pf :: type: " << pf->PFType() << "\tpt: " << pf->Pt() << "\tdR: " << dr;
2545 if( pf->HasTrackerTrk() ) cout << "\tdZ: " << pf->TrackerTrk()->DzCorrected(*vtx)
2546 << "\ttrk: " << pf->HasTrackerTrk()
2547 << "\tgsf: " << pf->HasGsfTrk();
2548
2549 cout << endl;
2550 }
2551
2552
2553 //
2554 // sync : I don't think theyre doing this ...
2555 //
2556 // if ( (pf->HasTrackerTrk() && (pf->TrackerTrk() == ele->TrackerTrk())) ||
2557 // (pf->HasGsfTrk() && (pf->GsfTrk() == ele->GsfTrk()))) {
2558 // if( ctrl.debug ) cout << "\tskipping, matches to the electron ..." << endl;
2559 // continue;
2560 // }
2561
2562
2563 //
2564 // Lepton Footprint Removal
2565 //
2566 Bool_t IsLeptonFootprint = kFALSE;
2567 if (dr < 1.0) {
2568
2569 //
2570 // Check for electrons
2571 //
2572 for (Int_t q=0; q < electronsToVeto.size(); ++q) {
2573 const mithep::Electron *tmpele = electronsToVeto[q];
2574 /*
2575 // 4l electron
2576 if( pf->HasTrackerTrk() ) {
2577 if( pf->TrackerTrk() == tmpele->TrackerTrk() ) {
2578 if( ctrl.debug) cout << "\tcharged tktrk, matches 4L ele ..." << endl;
2579 IsLeptonFootprint = kTRUE;
2580 }
2581 }
2582 if( pf->HasGsfTrk() ) {
2583 if( pf->GsfTrk() == tmpele->GsfTrk() ) {
2584 if( ctrl.debug) cout << "\tcharged gsftrk, matches 4L ele ..." << endl;
2585 IsLeptonFootprint = kTRUE;
2586 }
2587 }
2588 */
2589 // PF charged
2590 if (pf->Charge() != 0 && fabs(tmpele->SCluster()->Eta()) > 1.479
2591 && mithep::MathUtils::DeltaR(tmpele->Phi(),tmpele->Eta(), pf->Phi(), pf->Eta()) < 0.015) {
2592 if( ctrl.debug) cout << "\tcharged trk, dR matches 4L ele ..." << endl;
2593 IsLeptonFootprint = kTRUE;
2594 }
2595 // PF gamma
2596 if (abs(pf->PFType()) == PFCandidate::eGamma && fabs(tmpele->SCluster()->Eta()) > 1.479
2597 && mithep::MathUtils::DeltaR(tmpele->Phi(),tmpele->Eta(), pf->Phi(), pf->Eta()) < 0.08) {
2598 if( ctrl.debug) cout << "\tPF gamma, matches 4L ele ..." << endl;
2599 IsLeptonFootprint = kTRUE;
2600 }
2601 } // loop over electrons
2602
2603 /* KH - comment for sync
2604 //
2605 // Check for muons
2606 //
2607 for (Int_t q=0; q < muonsToVeto.size(); ++q) {
2608 const mithep::Muon *tmpmu = muonsToVeto[q];
2609 // 4l muon
2610 if( pf->HasTrackerTrk() ) {
2611 if (pf->TrackerTrk() == tmpmu->TrackerTrk() ){
2612 if( ctrl.debug) cout << "\tmatches 4L mu ..." << endl;
2613 IsLeptonFootprint = kTRUE;
2614 }
2615 }
2616 // PF charged
2617 if (pf->Charge() != 0 && mithep::MathUtils::DeltaR(tmpmu->Phi(),tmpmu->Eta(), pf->Phi(), pf->Eta()) < 0.01) {
2618 if( ctrl.debug) cout << "\tcharged trk, dR matches 4L mu ..." << endl;
2619 IsLeptonFootprint = kTRUE;
2620 }
2621 } // loop over muons
2622 */
2623
2624 if (IsLeptonFootprint)
2625 continue;
2626
2627 //
2628 // Charged Iso
2629 //
2630 if (pf->Charge() != 0 && (pf->HasTrackerTrk()||pf->HasGsfTrk()) ) {
2631
2632 // if( pf->HasTrackerTrk() )
2633 // if (abs(pf->TrackerTrk()->DzCorrected(vtx)) > 0.2) continue;
2634 // if( pf->HasGsfTrk() )
2635 // if (abs(pf->GsfTrk()->DzCorrected(vtx)) > 0.2) continue;
2636
2637 // Veto any PFmuon, or PFEle
2638 if (abs(pf->PFType()) == PFCandidate::eElectron || abs(pf->PFType()) == PFCandidate::eMuon) {
2639 if( ctrl.debug ) cout << "\t skipping, pf is and ele or mu .." <<endl;
2640 continue;
2641 }
2642
2643 // Footprint Veto
2644 if (fabs(ele->SCluster()->Eta()) > 1.479 && dr < 0.015) continue;
2645
2646 if( ctrl.debug) cout << "charged:: pt: " << pf->Pt()
2647 << "\ttype: " << pf->PFType()
2648 << "\ttrk: " << pf->TrackerTrk() << endl;
2649
2650 fChargedIso += pf->Pt();
2651 }
2652
2653 //
2654 // Gamma Iso
2655 //
2656 else if (abs(pf->PFType()) == PFCandidate::eGamma) {
2657
2658 if (fabs(ele->SCluster()->Eta()) > 1.479) {
2659 if (mithep::MathUtils::DeltaR(ele->Phi(),ele->Eta(), pf->Phi(), pf->Eta()) < 0.08) continue;
2660 }
2661 if( ctrl.debug) cout << "gamma:: " << pf->Pt() << " "
2662 << dr << endl;
2663 // KH, add to sync
2664 // if( pf->Pt() > 0.5 )
2665 fGammaIso += pf->Pt();
2666 }
2667
2668 //
2669 // Neutral Iso
2670 //
2671 else {
2672 if( ctrl.debug) cout << "neutral:: " << pf->Pt() << " "
2673 << dr << endl;
2674 // KH, add to sync
2675 // if( pf->Pt() > 0.5 )
2676 fNeutralHadronIso += pf->Pt();
2677 }
2678
2679 }
2680
2681 }
2682
2683 // double rho = 0;
2684 // if (!(isnan(fPUEnergyDensity->At(0)->Rho()) || isinf(fPUEnergyDensity->At(0)->Rho())))
2685 // rho = fPUEnergyDensity->At(0)->Rho();
2686
2687 // WARNING!!!!
2688 // hardcode for sync ...
2689 EffectiveAreaVersion = eleT.kEleEAData2011;
2690 // WARNING!!!!
2691
2692
2693 double pfIso = fChargedIso + fmax(0.0,(fGammaIso + fNeutralHadronIso
2694 -rho*eleT.ElectronEffectiveArea(eleT.kEleGammaAndNeutralHadronIso04,
2695 ele->Eta(),EffectiveAreaVersion)));
2696
2697
2698 gChargedIso = fChargedIso;
2699 gGammaIso = fGammaIso;
2700 gNeutralIso = fNeutralHadronIso;
2701 return pfIso;
2702 }
2703
2704
2705 //--------------------------------------------------------------------------------------------------
2706 SelectionStatus electronReferenceIsoSelection(ControlFlags &ctrl,
2707 const mithep::Electron * ele,
2708 const mithep::Vertex * vtx,
2709 const mithep::Array<mithep::PFCandidate> * fPFCandidates,
2710 const mithep::Array<mithep::PileupEnergyDensity> * fPUEnergyDensity,
2711 mithep::ElectronTools::EElectronEffectiveAreaTarget EffectiveAreaVersion,
2712 vector<const mithep::Muon*> muonsToVeto,
2713 vector<const mithep::Electron*> electronsToVeto)
2714 //--------------------------------------------------------------------------------------------------
2715 {
2716
2717 SelectionStatus status;
2718
2719 double pfIso = electronPFIso04( ctrl, ele, vtx, fPFCandidates, fPUEnergyDensity,
2720 EffectiveAreaVersion, muonsToVeto ,electronsToVeto );
2721 // cout << "--------------> setting electron isoPF04 to " << pfIso << endl;
2722 status.isoPF04 = pfIso;
2723 status.chisoPF04 = gChargedIso;
2724 status.gaisoPF04 = gGammaIso;
2725 status.neisoPF04 = gNeutralIso;
2726
2727 bool pass = false;
2728 if( (pfIso/ele->Pt()) < ELECTRON_REFERENCE_PFISO_CUT ) pass = true;
2729
2730 if( pass ) {
2731 status.orStatus(SelectionStatus::LOOSEISO);
2732 status.orStatus(SelectionStatus::TIGHTISO);
2733 }
2734 if(ctrl.debug) cout << "returning status : " << hex << status.getStatus() << dec << endl;
2735 return status;
2736
2737 }
2738
2739
2740 //--------------------------------------------------------------------------------------------------
2741 // hacked version
2742 SelectionStatus electronReferenceIsoSelection(ControlFlags &ctrl,
2743 const mithep::Electron * ele,
2744 const mithep::Vertex * vtx,
2745 const mithep::Array<mithep::PFCandidate> * fPFCandidates,
2746 float rho,
2747 mithep::ElectronTools::EElectronEffectiveAreaTarget EffectiveAreaVersion,
2748 vector<const mithep::Muon*> muonsToVeto,
2749 vector<const mithep::Electron*> electronsToVeto)
2750 //--------------------------------------------------------------------------------------------------
2751 {
2752
2753 SelectionStatus status;
2754
2755 double pfIso = electronPFIso04( ctrl, ele, vtx, fPFCandidates, rho,
2756 EffectiveAreaVersion, muonsToVeto ,electronsToVeto );
2757 status.isoPF04 = pfIso;
2758 status.chisoPF04 = gChargedIso;
2759 status.gaisoPF04 = gGammaIso;
2760 status.neisoPF04 = gNeutralIso;
2761 bool pass = false;
2762 if( (pfIso/ele->Pt()) < ELECTRON_REFERENCE_PFISO_CUT ) pass = true;
2763
2764 if( pass ) {
2765 status.orStatus(SelectionStatus::LOOSEISO);
2766 status.orStatus(SelectionStatus::TIGHTISO);
2767 }
2768 if(ctrl.debug) cout << "returning status : " << hex << status.getStatus() << dec << endl;
2769 return status;
2770
2771 }