ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitHzz4l/LeptonSelection/src/IsolationSelection.cc
Revision: 1.10
Committed: Sat May 5 13:09:05 2012 UTC (13 years ago) by khahn
Content type: text/plain
Branch: MAIN
Changes since 1.9: +497 -26 lines
Log Message:
stuff for reference selection

File Contents

# User Rev Content
1 khahn 1.1 #include <math.h>
2    
3     #include "IsolationSelection.h"
4     #include "IsolationSelectionDefs.h"
5    
6 khahn 1.4 #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     //--------------------------------------------------------------------------------------------------
20     Float_t computePFMuonIso(const mithep::Muon *muon,
21     const mithep::Vertex & vtx,
22     const mithep::Array<mithep::PFCandidate> * fPFCandidates,
23     const Double_t dRMax)
24     //--------------------------------------------------------------------------------------------------
25     {
26     const Double_t dRMin = 0;
27     const Double_t neuPtMin = 1.0;
28     const Double_t dzMax = 0.1;
29    
30     Double_t zLepton = (muon->BestTrk()) ? muon->BestTrk()->DzCorrected(vtx) : 0.0;
31    
32     Float_t iso=0;
33     for(UInt_t ipf=0; ipf<fPFCandidates->GetEntries(); ipf++) {
34     const PFCandidate *pfcand = fPFCandidates->At(ipf);
35    
36     if(!pfcand->HasTrk() && (pfcand->Pt()<=neuPtMin)) continue; // pT cut on neutral particles
37    
38     // exclude THE muon
39     if(pfcand->TrackerTrk() && muon->TrackerTrk() && (pfcand->TrackerTrk()==muon->TrackerTrk())) continue;
40    
41     // dz cut
42     Double_t dz = (pfcand->BestTrk()) ? fabs(pfcand->BestTrk()->DzCorrected(vtx) - zLepton) : 0;
43     if(dz >= dzMax) continue;
44    
45     // check iso cone
46     Double_t dr = MathUtils::DeltaR(muon->Mom(), pfcand->Mom());
47     if(dr<dRMax && dr>=dRMin)
48     iso += pfcand->Pt();
49     }
50    
51     return iso;
52     }
53    
54     //--------------------------------------------------------------------------------------------------
55     Float_t computePFEleIso(const mithep::Electron *electron,
56     const mithep::Vertex & fVertex,
57     const mithep::Array<mithep::PFCandidate> * fPFCandidates,
58     const Double_t dRMax)
59     //--------------------------------------------------------------------------------------------------
60     {
61     const Double_t dRMin = 0;
62     const Double_t neuPtMin = 1.0;
63     const Double_t dzMax = 0.1;
64    
65     Double_t zLepton = (electron->BestTrk()) ? electron->BestTrk()->DzCorrected(fVertex) : 0.0;
66    
67     Float_t iso=0;
68     for(UInt_t ipf=0; ipf<fPFCandidates->GetEntries(); ipf++) {
69     const PFCandidate *pfcand = (PFCandidate*)(fPFCandidates->At(ipf));
70    
71     if(!pfcand->HasTrk() && (pfcand->Pt()<=neuPtMin)) continue; // pT cut on neutral particles
72    
73     // dz cut
74     Double_t dz = (pfcand->BestTrk()) ? fabs(pfcand->BestTrk()->DzCorrected(fVertex) - zLepton) : 0;
75     if(dz >= dzMax) continue;
76    
77     // remove THE electron
78     if(pfcand->TrackerTrk() && electron->TrackerTrk() && (pfcand->TrackerTrk()==electron->TrackerTrk())) continue;
79     if(pfcand->GsfTrk() && electron->GsfTrk() && (pfcand->GsfTrk()==electron->GsfTrk())) continue;
80    
81     // check iso cone
82     Double_t dr = MathUtils::DeltaR(electron->Mom(), pfcand->Mom());
83     if(dr<dRMax && dr>=dRMin) {
84     // eta-strip veto for photons
85     if((pfcand->PFType() == PFCandidate::eGamma) && fabs(electron->Eta() - pfcand->Eta()) < 0.025) continue;
86    
87     // Inner cone (one tower = dR < 0.07) veto for non-photon neutrals
88     if(!pfcand->HasTrk() && (pfcand->PFType() == PFCandidate::eNeutralHadron) &&
89     (MathUtils::DeltaR(electron->Mom(), pfcand->Mom()) < 0.07)) continue;
90    
91     iso += pfcand->Pt();
92     }
93     }
94    
95     return iso;
96     };
97    
98     //--------------------------------------------------------------------------------------------------
99     bool pairwiseIsoSelection( ControlFlags &ctrl,
100     vector<SimpleLepton> &lepvec,
101     float rho )
102     //--------------------------------------------------------------------------------------------------
103     {
104 khahn 1.2
105     bool passiso=true;
106    
107     for( int i=0; i<lepvec.size(); i++ )
108     {
109    
110     if( !(lepvec[i].is4l) ) continue;
111    
112     float effArea_ecal_i, effArea_hcal_i;
113     if( lepvec[i].isEB ) {
114     if( lepvec[i].type == 11 ) {
115     effArea_ecal_i = 0.101;
116     effArea_hcal_i = 0.021;
117     } else {
118     effArea_ecal_i = 0.074;
119     effArea_hcal_i = 0.022;
120     }
121     } else {
122     if( lepvec[i].type == 11 ) {
123     effArea_ecal_i = 0.046;
124     effArea_hcal_i = 0.040;
125     } else {
126     effArea_ecal_i = 0.045;
127     effArea_hcal_i = 0.030;
128     }
129     }
130    
131     float isoEcal_corr_i = lepvec[i].isoEcal - (effArea_ecal_i*rho);
132     float isoHcal_corr_i = lepvec[i].isoHcal - (effArea_hcal_i*rho);
133    
134     for( int j=i+1; j<lepvec.size(); j++ )
135     {
136    
137     if( !(lepvec[j].is4l) ) continue;
138    
139     float effArea_ecal_j, effArea_hcal_j;
140     if( lepvec[j].isEB ) {
141     if( lepvec[j].type == 11 ) {
142     effArea_ecal_j = 0.101;
143     effArea_hcal_j = 0.021;
144     } else {
145     effArea_ecal_j = 0.074;
146     effArea_hcal_j = 0.022;
147     }
148     } else {
149     if( lepvec[j].type == 11 ) {
150     effArea_ecal_j = 0.046;
151     effArea_hcal_j = 0.040;
152     } else {
153     effArea_ecal_j = 0.045;
154     effArea_hcal_j = 0.030;
155     }
156     }
157    
158     float isoEcal_corr_j = lepvec[j].isoEcal - (effArea_ecal_j*rho);
159     float isoHcal_corr_j = lepvec[j].isoHcal - (effArea_hcal_j*rho);
160     float RIso_i = (lepvec[i].isoTrk+isoEcal_corr_i+isoHcal_corr_i)/lepvec[i].vec->Pt();
161     float RIso_j = (lepvec[j].isoTrk+isoEcal_corr_j+isoHcal_corr_j)/lepvec[j].vec->Pt();
162     float comboIso = RIso_i + RIso_j;
163    
164     if( comboIso > 0.35 ) {
165     if( ctrl.debug ) cout << "combo failing for indices: " << i << "," << j << endl;
166     passiso = false;
167     return passiso;
168     }
169     }
170     }
171    
172     return passiso;
173     }
174    
175 khahn 1.4 //--------------------------------------------------------------------------------------------------
176     SelectionStatus muonIsoSelection(ControlFlags &ctrl,
177     const mithep::Muon * mu,
178     const mithep::Vertex & vtx,
179     const mithep::Array<mithep::PFCandidate> * fPFCandidateCol )
180     //--------------------------------------------------------------------------------------------------
181     {
182     float reliso = computePFMuonIso(mu,vtx,fPFCandidateCol,0.3)/mu->Pt();
183     bool isEB = (fabs(mu->Eta()) < 1.479 ? 1 : 0 );
184 khahn 1.1 bool failiso = false;
185 khahn 1.4 if( isEB && mu->Pt() > 20 && reliso > PFISO_MU_LOOSE_EB_HIGHPT ) {
186 khahn 1.1 failiso = true;
187     }
188 khahn 1.4 if( isEB && mu->Pt() < 20 && reliso > PFISO_MU_LOOSE_EB_LOWPT ) {
189 khahn 1.1 failiso = true;
190     }
191 khahn 1.4 if( !(isEB) && mu->Pt() > 20 && reliso > PFISO_MU_LOOSE_EE_HIGHPT ) {
192 khahn 1.1 failiso = true;
193     }
194 khahn 1.4 if( !(isEB) && mu->Pt() < 20 && reliso > PFISO_MU_LOOSE_EE_LOWPT ) {
195 khahn 1.1 failiso = true;
196     }
197    
198     SelectionStatus status;
199     if( !failiso ) status.setStatus(SelectionStatus::LOOSEISO);
200     if( !failiso ) status.setStatus(SelectionStatus::TIGHTISO);
201     return status;
202    
203     };
204    
205 khahn 1.4 //--------------------------------------------------------------------------------------------------
206     SelectionStatus electronIsoSelection(ControlFlags &ctrl,
207     const mithep::Electron * ele,
208     const mithep::Vertex &fVertex,
209     const mithep::Array<mithep::PFCandidate> * fPFCandidates)
210     //--------------------------------------------------------------------------------------------------
211     {
212 khahn 1.1
213 khahn 1.4 bool failiso=false;
214 khahn 1.1
215 khahn 1.4 float reliso = computePFEleIso(ele,fVertex,fPFCandidates,0.4)/ele->Pt();
216 khahn 1.1
217 khahn 1.4 if( ele->IsEB() && ele->Pt() > 20 && reliso > PFISO_ELE_LOOSE_EB_HIGHPT ) {
218 khahn 1.1 failiso = true;
219     }
220 khahn 1.4 if( ele->IsEB() && ele->Pt() < 20 && reliso > PFISO_ELE_LOOSE_EB_LOWPT ) {
221 khahn 1.1 failiso = true;
222     }
223     if(ctrl.debug) cout << "before iso check ..." << endl;
224 khahn 1.4 if( !(ele->IsEB()) && ele->Pt() > 20 && reliso > PFISO_ELE_LOOSE_EE_HIGHPT ) {
225 khahn 1.1 if(ctrl.debug) cout << "\tit fails ..." << endl;
226     failiso = true;
227     }
228 khahn 1.4 if( !(ele->IsEB()) && ele->Pt() < 20 && reliso > PFISO_ELE_LOOSE_EE_LOWPT ) {
229 khahn 1.1 failiso = true;
230     }
231    
232     SelectionStatus status;
233     if( !failiso ) {
234 khahn 1.4 status.orStatus(SelectionStatus::LOOSEISO);
235     status.orStatus(SelectionStatus::TIGHTISO);
236 khahn 1.1 }
237     if(ctrl.debug) cout << "returning status : " << hex << status.getStatus() << dec << endl;
238     return status;
239    
240     }
241 anlevin 1.3
242 khahn 1.4
243 anlevin 1.3 bool noIso(ControlFlags &, vector<SimpleLepton> &, float rho) {
244    
245     return true;
246     }
247 khahn 1.4
248     //--------------------------------------------------------------------------------------------------
249     SelectionStatus muonIsoMVASelection(ControlFlags &ctrl,
250     const mithep::Muon * mu,
251     const mithep::Vertex & vtx,
252     const mithep::Array<mithep::PFCandidate> * fPFCandidates,
253     const mithep::Array<mithep::PileupEnergyDensity> * fPUEnergyDensity,
254     mithep::MuonTools::EMuonEffectiveAreaTarget EffectiveAreaVersion,
255     vector<const mithep::Muon*> muonsToVeto,
256     vector<const mithep::Electron*> electronsToVeto)
257     //--------------------------------------------------------------------------------------------------
258     {
259    
260 khahn 1.7 if( ctrl.debug ) {
261     cout << "muonIsoMVASelection :: muons to veto " << endl;
262     for( int i=0; i<muonsToVeto.size(); i++ ) {
263     const mithep::Muon * vmu = muonsToVeto[i];
264     cout << "\tpt: " << vmu->Pt()
265     << "\teta: " << vmu->Eta()
266     << "\tphi: " << vmu->Phi()
267     << endl;
268     }
269     cout << "muonIsoMVASelection :: electrson to veto " << endl;
270     for( int i=0; i<electronsToVeto.size(); i++ ) {
271     const mithep::Electron * vel = electronsToVeto[i];
272     cout << "\tpt: " << vel->Pt()
273     << "\teta: " << vel->Eta()
274     << "\tphi: " << vel->Phi()
275     << endl;
276     }
277     }
278 khahn 1.4 bool failiso=false;
279    
280     //
281     // tmp iso rings
282     //
283 khahn 1.7 Double_t tmpChargedIso_DR0p0To0p1 = 0;
284     Double_t tmpChargedIso_DR0p1To0p2 = 0;
285     Double_t tmpChargedIso_DR0p2To0p3 = 0;
286 khahn 1.4 Double_t tmpChargedIso_DR0p3To0p4 = 0;
287     Double_t tmpChargedIso_DR0p4To0p5 = 0;
288 khahn 1.7 Double_t tmpChargedIso_DR0p5To0p7 = 0;
289 khahn 1.4
290 khahn 1.7 Double_t tmpGammaIso_DR0p0To0p1 = 0;
291     Double_t tmpGammaIso_DR0p1To0p2 = 0;
292     Double_t tmpGammaIso_DR0p2To0p3 = 0;
293 khahn 1.4 Double_t tmpGammaIso_DR0p3To0p4 = 0;
294     Double_t tmpGammaIso_DR0p4To0p5 = 0;
295 khahn 1.7 Double_t tmpGammaIso_DR0p5To0p7 = 0;
296 khahn 1.4
297 khahn 1.7 Double_t tmpNeutralHadronIso_DR0p0To0p1 = 0;
298     Double_t tmpNeutralHadronIso_DR0p1To0p2 = 0;
299     Double_t tmpNeutralHadronIso_DR0p2To0p3 = 0;
300 khahn 1.4 Double_t tmpNeutralHadronIso_DR0p3To0p4 = 0;
301     Double_t tmpNeutralHadronIso_DR0p4To0p5 = 0;
302 khahn 1.7 Double_t tmpNeutralHadronIso_DR0p5To0p7 = 0;
303    
304 khahn 1.4
305    
306     //
307     // final rings for the MVA
308     //
309     Double_t fChargedIso_DR0p0To0p1;
310     Double_t fChargedIso_DR0p1To0p2;
311     Double_t fChargedIso_DR0p2To0p3;
312     Double_t fChargedIso_DR0p3To0p4;
313     Double_t fChargedIso_DR0p4To0p5;
314 khahn 1.7 Double_t fChargedIso_DR0p5To0p7;
315 khahn 1.4
316     Double_t fGammaIso_DR0p0To0p1;
317     Double_t fGammaIso_DR0p1To0p2;
318     Double_t fGammaIso_DR0p2To0p3;
319     Double_t fGammaIso_DR0p3To0p4;
320     Double_t fGammaIso_DR0p4To0p5;
321 khahn 1.7 Double_t fGammaIso_DR0p5To0p7;
322 khahn 1.4
323     Double_t fNeutralHadronIso_DR0p0To0p1;
324     Double_t fNeutralHadronIso_DR0p1To0p2;
325     Double_t fNeutralHadronIso_DR0p2To0p3;
326     Double_t fNeutralHadronIso_DR0p3To0p4;
327     Double_t fNeutralHadronIso_DR0p4To0p5;
328 khahn 1.7 Double_t fNeutralHadronIso_DR0p5To0p7;
329 khahn 1.4
330    
331     //
332     //Loop over PF Candidates
333     //
334     for(int k=0; k<fPFCandidates->GetEntries(); ++k) {
335     const mithep::PFCandidate *pf = (mithep::PFCandidate*)((*fPFCandidates)[k]);
336    
337     Double_t deta = (mu->Eta() - pf->Eta());
338     Double_t dphi = mithep::MathUtils::DeltaPhi(Double_t(mu->Phi()),Double_t(pf->Phi()));
339     Double_t dr = mithep::MathUtils::DeltaR(mu->Phi(),mu->Eta(), pf->Phi(), pf->Eta());
340     if (dr > 1.0) continue;
341    
342 khahn 1.7 if (pf->HasTrackerTrk() && (pf->TrackerTrk() == mu->TrackerTrk()) ) continue;
343 khahn 1.4
344     //
345     // Lepton Footprint Removal
346     //
347     Bool_t IsLeptonFootprint = kFALSE;
348     if (dr < 1.0) {
349    
350     //
351     // Check for electrons
352     //
353     for (Int_t q=0; q < electronsToVeto.size(); ++q) {
354     const mithep::Electron *tmpele = electronsToVeto[q];
355 khahn 1.7 // 4l electron
356     if( pf->HasTrackerTrk() ) {
357     if( pf->TrackerTrk() == tmpele->TrackerTrk() )
358     IsLeptonFootprint = kTRUE;
359     }
360     if( pf->HasGsfTrk() ) {
361     if( pf->GsfTrk() == tmpele->GsfTrk() )
362     IsLeptonFootprint = kTRUE;
363     }
364 khahn 1.4 // PF charged
365     if (pf->Charge() != 0 && fabs(tmpele->SCluster()->Eta()) > 1.479
366     && mithep::MathUtils::DeltaR(tmpele->Phi(),tmpele->Eta(), pf->Phi(), pf->Eta()) < 0.015)
367     IsLeptonFootprint = kTRUE;
368     // PF gamma
369 khahn 1.7 if (abs(pf->PFType()) == PFCandidate::eGamma && fabs(tmpele->SCluster()->Eta()) > 1.479
370 khahn 1.4 && mithep::MathUtils::DeltaR(tmpele->Phi(),tmpele->Eta(), pf->Phi(), pf->Eta()) < 0.08)
371     IsLeptonFootprint = kTRUE;
372     } // loop over electrons
373    
374     //
375     // Check for muons
376     //
377     for (Int_t q=0; q < muonsToVeto.size(); ++q) {
378     const mithep::Muon *tmpmu = muonsToVeto[q];
379 khahn 1.7 // 4l muon
380     if( pf->HasTrackerTrk() ) {
381     if( pf->TrackerTrk() == tmpmu->TrackerTrk() )
382     IsLeptonFootprint = kTRUE;
383     }
384 khahn 1.4 // PF charged
385     if (pf->Charge() != 0 && mithep::MathUtils::DeltaR(tmpmu->Phi(),tmpmu->Eta(), pf->Phi(), pf->Eta()) < 0.01)
386     IsLeptonFootprint = kTRUE;
387     } // loop over muons
388    
389    
390     if (IsLeptonFootprint)
391     continue;
392    
393     //
394     // Charged Iso Rings
395     //
396 khahn 1.7 if (pf->Charge() != 0 && (pf->HasTrackerTrk()||pf->HasGsfTrk()) ) {
397    
398     if( dr < 0.01 ) continue; // only for muon iso mva?
399     if (abs(pf->PFType()) == PFCandidate::eElectron || abs(pf->PFType()) == PFCandidate::eMuon) continue;
400 khahn 1.4
401 khahn 1.7 if( pf->HasTrackerTrk() ) {
402     if (abs(pf->TrackerTrk()->DzCorrected(vtx)) > 0.2) continue;
403     if( ctrl.debug ) cout << "charged:: " << pf->PFType() << " " << pf->Pt() << " "
404     << abs(pf->TrackerTrk()->DzCorrected(vtx)) << " "
405     << dr << endl;
406     }
407     if( pf->HasGsfTrk() ) {
408     if (abs(pf->GsfTrk()->DzCorrected(vtx)) > 0.2) continue;
409     if( ctrl.debug ) cout << "charged:: " << pf->PFType() << " " << pf->Pt() << " "
410     << abs(pf->GsfTrk()->DzCorrected(vtx)) << " "
411     << dr << endl;
412     }
413 khahn 1.4
414     // Footprint Veto
415 khahn 1.7 if (dr < 0.1) tmpChargedIso_DR0p0To0p1 += pf->Pt();
416     if (dr >= 0.1 && dr < 0.2) tmpChargedIso_DR0p1To0p2 += pf->Pt();
417     if (dr >= 0.2 && dr < 0.3) tmpChargedIso_DR0p2To0p3 += pf->Pt();
418     if (dr >= 0.3 && dr < 0.4) tmpChargedIso_DR0p3To0p4 += pf->Pt();
419     if (dr >= 0.4 && dr < 0.5) tmpChargedIso_DR0p4To0p5 += pf->Pt();
420     if (dr >= 0.5 && dr < 0.7) tmpChargedIso_DR0p5To0p7 += pf->Pt();
421 khahn 1.4 }
422    
423     //
424     // Gamma Iso Rings
425     //
426 khahn 1.7 else if (abs(pf->PFType()) == PFCandidate::eGamma) {
427     if (dr < 0.1) tmpGammaIso_DR0p0To0p1 += pf->Pt();
428     if (dr >= 0.1 && dr < 0.2) tmpGammaIso_DR0p1To0p2 += pf->Pt();
429     if (dr >= 0.2 && dr < 0.3) tmpGammaIso_DR0p2To0p3 += pf->Pt();
430     if (dr >= 0.3 && dr < 0.4) tmpGammaIso_DR0p3To0p4 += pf->Pt();
431     if (dr >= 0.4 && dr < 0.5) tmpGammaIso_DR0p4To0p5 += pf->Pt();
432     if (dr >= 0.5 && dr < 0.7) tmpGammaIso_DR0p5To0p7 += pf->Pt();
433 khahn 1.4 }
434    
435     //
436     // Other Neutral Iso Rings
437     //
438     else {
439 khahn 1.7 if (dr < 0.1) tmpNeutralHadronIso_DR0p0To0p1 += pf->Pt();
440     if (dr >= 0.1 && dr < 0.2) tmpNeutralHadronIso_DR0p1To0p2 += pf->Pt();
441     if (dr >= 0.2 && dr < 0.3) tmpNeutralHadronIso_DR0p2To0p3 += pf->Pt();
442     if (dr >= 0.3 && dr < 0.4) tmpNeutralHadronIso_DR0p3To0p4 += pf->Pt();
443     if (dr >= 0.4 && dr < 0.5) tmpNeutralHadronIso_DR0p4To0p5 += pf->Pt();
444     if (dr >= 0.5 && dr < 0.7) tmpNeutralHadronIso_DR0p5To0p7 += pf->Pt();
445 khahn 1.4 }
446    
447     }
448    
449     }
450    
451 khahn 1.7 fChargedIso_DR0p0To0p1 = min((tmpChargedIso_DR0p0To0p1)/mu->Pt(), 2.5);
452     fChargedIso_DR0p1To0p2 = min((tmpChargedIso_DR0p1To0p2)/mu->Pt(), 2.5);
453     fChargedIso_DR0p2To0p3 = min((tmpChargedIso_DR0p2To0p3)/mu->Pt(), 2.5);
454 khahn 1.4 fChargedIso_DR0p3To0p4 = min((tmpChargedIso_DR0p3To0p4)/mu->Pt(), 2.5);
455     fChargedIso_DR0p4To0p5 = min((tmpChargedIso_DR0p4To0p5)/mu->Pt(), 2.5);
456    
457 khahn 1.7
458 khahn 1.4 double rho = 0;
459     if (!(isnan(fPUEnergyDensity->At(0)->Rho()) || isinf(fPUEnergyDensity->At(0)->Rho())))
460     rho = fPUEnergyDensity->At(0)->Rho();
461    
462    
463 khahn 1.7 fGammaIso_DR0p0To0p1 = max(min((tmpGammaIso_DR0p0To0p1
464 khahn 1.4 -rho*muT.MuonEffectiveArea(muT.kMuGammaIsoDR0p0To0p1,mu->Eta(),EffectiveAreaVersion))/mu->Pt()
465     ,2.5)
466     ,0.0);
467 khahn 1.7 fGammaIso_DR0p1To0p2 = max(min((tmpGammaIso_DR0p1To0p2
468 khahn 1.4 -rho*muT.MuonEffectiveArea(muT.kMuGammaIsoDR0p1To0p2,mu->Eta(),EffectiveAreaVersion))/mu->Pt()
469     ,2.5)
470     ,0.0);
471 khahn 1.7 fGammaIso_DR0p2To0p3 = max(min((tmpGammaIso_DR0p2To0p3
472 khahn 1.4 -rho*muT.MuonEffectiveArea(muT.kMuGammaIsoDR0p2To0p3,mu->Eta(),EffectiveAreaVersion))/mu->Pt()
473     ,2.5)
474     ,0.0);
475     fGammaIso_DR0p3To0p4 = max(min((tmpGammaIso_DR0p3To0p4
476     -rho*muT.MuonEffectiveArea(muT.kMuGammaIsoDR0p3To0p4,mu->Eta(),EffectiveAreaVersion))/mu->Pt()
477     ,2.5)
478     ,0.0);
479     fGammaIso_DR0p4To0p5 = max(min((tmpGammaIso_DR0p4To0p5
480     -rho*muT.MuonEffectiveArea(muT.kMuGammaIsoDR0p4To0p5,mu->Eta(),EffectiveAreaVersion))/mu->Pt()
481     ,2.5)
482     ,0.0);
483    
484    
485 khahn 1.7
486     fNeutralHadronIso_DR0p0To0p1 = max(min((tmpNeutralHadronIso_DR0p0To0p1
487 khahn 1.4 -rho*muT.MuonEffectiveArea(muT.kMuNeutralHadronIsoDR0p0To0p1,
488     mu->Eta(),EffectiveAreaVersion))/mu->Pt()
489     , 2.5)
490     , 0.0);
491 khahn 1.7 fNeutralHadronIso_DR0p1To0p2 = max(min((tmpNeutralHadronIso_DR0p1To0p2
492 khahn 1.4 -rho*muT.MuonEffectiveArea(muT.kMuNeutralHadronIsoDR0p1To0p2,
493     mu->Eta(),EffectiveAreaVersion))/mu->Pt()
494     , 2.5)
495     , 0.0);
496 khahn 1.7 fNeutralHadronIso_DR0p2To0p3 = max(min((tmpNeutralHadronIso_DR0p2To0p3
497 khahn 1.4 -rho*muT.MuonEffectiveArea(muT.kMuNeutralHadronIsoDR0p2To0p3,
498     mu->Eta(),EffectiveAreaVersion))/mu->Pt()
499     , 2.5)
500     , 0.0);
501     fNeutralHadronIso_DR0p3To0p4 = max(min((tmpNeutralHadronIso_DR0p3To0p4
502     -rho*muT.MuonEffectiveArea(muT.kMuNeutralHadronIsoDR0p3To0p4,
503     mu->Eta(), EffectiveAreaVersion))/mu->Pt()
504     , 2.5)
505     , 0.0);
506     fNeutralHadronIso_DR0p4To0p5 = max(min((tmpNeutralHadronIso_DR0p4To0p5
507     -rho*muT.MuonEffectiveArea(muT.kMuNeutralHadronIsoDR0p4To0p5,
508     mu->Eta(), EffectiveAreaVersion))/mu->Pt()
509     , 2.5)
510     , 0.0);
511    
512 khahn 1.7
513 khahn 1.4 double mvaval = muIsoMVA->MVAValue_IsoRings( mu->Pt(),
514     mu->Eta(),
515     fChargedIso_DR0p0To0p1,
516     fChargedIso_DR0p1To0p2,
517     fChargedIso_DR0p2To0p3,
518     fChargedIso_DR0p3To0p4,
519     fChargedIso_DR0p4To0p5,
520     fGammaIso_DR0p0To0p1,
521     fGammaIso_DR0p1To0p2,
522     fGammaIso_DR0p2To0p3,
523     fGammaIso_DR0p3To0p4,
524     fGammaIso_DR0p4To0p5,
525     fNeutralHadronIso_DR0p0To0p1,
526     fNeutralHadronIso_DR0p1To0p2,
527     fNeutralHadronIso_DR0p2To0p3,
528     fNeutralHadronIso_DR0p3To0p4,
529     fNeutralHadronIso_DR0p4To0p5,
530     ctrl.debug);
531    
532     SelectionStatus status;
533 khahn 1.10 bool pass;
534 khahn 1.4
535 khahn 1.10 pass = false;
536 khahn 1.4 if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
537 khahn 1.10 && fabs(mu->Eta()) <= 1.5 && mu->Pt() <= 10 && mvaval >= MUON_ISOMVA_LOOSE_FORPFID_CUT_BIN0) pass = true;
538 khahn 1.4 else if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
539 khahn 1.10 && fabs(mu->Eta()) <= 1.5 && mu->Pt() > 10 && mvaval >= MUON_ISOMVA_LOOSE_FORPFID_CUT_BIN1) pass = true;
540 khahn 1.4 else if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
541 khahn 1.10 && fabs(mu->Eta()) > 1.5 && mu->Pt() <= 10 && mvaval >= MUON_ISOMVA_LOOSE_FORPFID_CUT_BIN2) pass = true;
542 khahn 1.4 else if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
543 khahn 1.10 && fabs(mu->Eta()) > 1.5 && mu->Pt() > 10 && mvaval >= MUON_ISOMVA_LOOSE_FORPFID_CUT_BIN3) pass = true;
544     else if( !(mu->IsGlobalMuon()) && mu->IsTrackerMuon() && mvaval >= MUON_ISOMVA_LOOSE_FORPFID_CUT_BIN4) pass = true;
545     else if( mu->IsGlobalMuon() && !(mu->IsTrackerMuon()) && mvaval >= MUON_ISOMVA_LOOSE_FORPFID_CUT_BIN5) pass = true;
546     if( pass ) status.orStatus(SelectionStatus::LOOSEISO);
547 khahn 1.4
548 khahn 1.10 pass = false;
549     if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
550     && fabs(mu->Eta()) <= 1.5 && mu->Pt() <= 10 && mvaval >= MUON_ISOMVA_TIGHT_FORPFID_CUT_BIN0) pass = true;
551     else if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
552     && fabs(mu->Eta()) <= 1.5 && mu->Pt() > 10 && mvaval >= MUON_ISOMVA_TIGHT_FORPFID_CUT_BIN1) pass = true;
553     else if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
554     && fabs(mu->Eta()) > 1.5 && mu->Pt() <= 10 && mvaval >= MUON_ISOMVA_TIGHT_FORPFID_CUT_BIN2) pass = true;
555     else if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
556     && fabs(mu->Eta()) > 1.5 && mu->Pt() > 10 && mvaval >= MUON_ISOMVA_TIGHT_FORPFID_CUT_BIN3) pass = true;
557     else if( !(mu->IsGlobalMuon()) && mu->IsTrackerMuon() && mvaval >= MUON_ISOMVA_TIGHT_FORPFID_CUT_BIN4) pass = true;
558     else if( mu->IsGlobalMuon() && !(mu->IsTrackerMuon()) && mvaval >= MUON_ISOMVA_TIGHT_FORPFID_CUT_BIN5) pass = true;
559     if( pass ) status.orStatus(SelectionStatus::TIGHTISO);
560 khahn 1.4
561 khahn 1.10 // pass &= (fChargedIso_DR0p0To0p1 + fChargedIso_DR0p1To0p2 + fChargedIso_DR0p2To0p3 < 0.7);
562 khahn 1.7
563 khahn 1.4 if(ctrl.debug) cout << "returning status : " << hex << status.getStatus() << dec << endl;
564     return status;
565    
566     }
567    
568     //--------------------------------------------------------------------------------------------------
569     void initMuonIsoMVA() {
570     //--------------------------------------------------------------------------------------------------
571     muIsoMVA = new mithep::MuonIDMVA();
572     vector<string> weightFiles;
573 khahn 1.6 weightFiles.push_back("./data/MuonIsoMVAWeights/MuonIsoMVA_BDTG_V0_barrel_lowpt.weights.xml");
574     weightFiles.push_back("./data/MuonIsoMVAWeights/MuonIsoMVA_BDTG_V0_barrel_highpt.weights.xml");
575     weightFiles.push_back("./data/MuonIsoMVAWeights/MuonIsoMVA_BDTG_V0_endcap_lowpt.weights.xml");
576     weightFiles.push_back("./data/MuonIsoMVAWeights/MuonIsoMVA_BDTG_V0_endcap_highpt.weights.xml");
577     weightFiles.push_back("./data/MuonIsoMVAWeights/MuonIsoMVA_BDTG_V0_tracker.weights.xml");
578     weightFiles.push_back("./data/MuonIsoMVAWeights/MuonIsoMVA_BDTG_V0_global.weights.xml");
579 khahn 1.4 muIsoMVA->Initialize( "MuonIsoMVA",
580     mithep::MuonIDMVA::kIsoRingsV0,
581     kTRUE, weightFiles);
582     }
583    
584    
585 khahn 1.10 //--------------------------------------------------------------------------------------------------
586     double muonPFIso04(ControlFlags &ctrl,
587     const mithep::Muon * mu,
588     const mithep::Vertex & vtx,
589     const mithep::Array<mithep::PFCandidate> * fPFCandidates,
590     const mithep::Array<mithep::PileupEnergyDensity> * fPUEnergyDensity,
591     mithep::MuonTools::EMuonEffectiveAreaTarget EffectiveAreaVersion,
592     vector<const mithep::Muon*> muonsToVeto,
593     vector<const mithep::Electron*> electronsToVeto)
594     //--------------------------------------------------------------------------------------------------
595     {
596    
597     if( ctrl.debug ) {
598     cout << "muonIsoMVASelection :: muons to veto " << endl;
599     for( int i=0; i<muonsToVeto.size(); i++ ) {
600     const mithep::Muon * vmu = muonsToVeto[i];
601     cout << "\tpt: " << vmu->Pt()
602     << "\teta: " << vmu->Eta()
603     << "\tphi: " << vmu->Phi()
604     << endl;
605     }
606     cout << "muonIsoMVASelection :: electrson to veto " << endl;
607     for( int i=0; i<electronsToVeto.size(); i++ ) {
608     const mithep::Electron * vel = electronsToVeto[i];
609     cout << "\tpt: " << vel->Pt()
610     << "\teta: " << vel->Eta()
611     << "\tphi: " << vel->Phi()
612     << endl;
613     }
614     }
615     bool failiso=false;
616    
617     //
618     // tmp iso
619     //
620     Double_t tmpChargedIso = 0;
621     Double_t tmpGammaIso = 0;
622     Double_t tmpNeutralHadronIso = 0;
623    
624     //
625     // final iso
626     //
627     Double_t fChargedIso;
628     Double_t fGammaIso;
629     Double_t fNeutralHadronIso;
630    
631     //
632     //Loop over PF Candidates
633     //
634     for(int k=0; k<fPFCandidates->GetEntries(); ++k) {
635     const mithep::PFCandidate *pf = (mithep::PFCandidate*)((*fPFCandidates)[k]);
636    
637     Double_t deta = (mu->Eta() - pf->Eta());
638     Double_t dphi = mithep::MathUtils::DeltaPhi(Double_t(mu->Phi()),Double_t(pf->Phi()));
639     Double_t dr = mithep::MathUtils::DeltaR(mu->Phi(),mu->Eta(), pf->Phi(), pf->Eta());
640     if (dr > 0.4) continue;
641    
642     if (pf->HasTrackerTrk() && (pf->TrackerTrk() == mu->TrackerTrk()) ) continue;
643    
644     //
645     // Lepton Footprint Removal
646     //
647     Bool_t IsLeptonFootprint = kFALSE;
648     if (dr < 1.0) {
649    
650     //
651     // Check for electrons
652     //
653     for (Int_t q=0; q < electronsToVeto.size(); ++q) {
654     const mithep::Electron *tmpele = electronsToVeto[q];
655     // 4l electron
656     if( pf->HasTrackerTrk() ) {
657     if( pf->TrackerTrk() == tmpele->TrackerTrk() )
658     IsLeptonFootprint = kTRUE;
659     }
660     if( pf->HasGsfTrk() ) {
661     if( pf->GsfTrk() == tmpele->GsfTrk() )
662     IsLeptonFootprint = kTRUE;
663     }
664     // PF charged
665     if (pf->Charge() != 0 && fabs(tmpele->SCluster()->Eta()) > 1.479
666     && mithep::MathUtils::DeltaR(tmpele->Phi(),tmpele->Eta(), pf->Phi(), pf->Eta()) < 0.015)
667     IsLeptonFootprint = kTRUE;
668     // PF gamma
669     if (abs(pf->PFType()) == PFCandidate::eGamma && fabs(tmpele->SCluster()->Eta()) > 1.479
670     && mithep::MathUtils::DeltaR(tmpele->Phi(),tmpele->Eta(), pf->Phi(), pf->Eta()) < 0.08)
671     IsLeptonFootprint = kTRUE;
672     } // loop over electrons
673    
674     //
675     // Check for muons
676     //
677     for (Int_t q=0; q < muonsToVeto.size(); ++q) {
678     const mithep::Muon *tmpmu = muonsToVeto[q];
679     // 4l muon
680     if( pf->HasTrackerTrk() ) {
681     if( pf->TrackerTrk() == tmpmu->TrackerTrk() )
682     IsLeptonFootprint = kTRUE;
683     }
684     // PF charged
685     if (pf->Charge() != 0 && mithep::MathUtils::DeltaR(tmpmu->Phi(),tmpmu->Eta(), pf->Phi(), pf->Eta()) < 0.01)
686     IsLeptonFootprint = kTRUE;
687     } // loop over muons
688    
689    
690     if (IsLeptonFootprint)
691     continue;
692    
693     //
694     // Charged Iso Rings
695     //
696     if (pf->Charge() != 0 && (pf->HasTrackerTrk()||pf->HasGsfTrk()) ) {
697    
698     if( dr < 0.01 ) continue; // only for muon iso mva?
699     if (abs(pf->PFType()) == PFCandidate::eElectron || abs(pf->PFType()) == PFCandidate::eMuon) continue;
700    
701     if( pf->HasTrackerTrk() ) {
702     if (abs(pf->TrackerTrk()->DzCorrected(vtx)) > 0.2) continue;
703     if( ctrl.debug ) cout << "charged:: " << pf->PFType() << " " << pf->Pt() << " "
704     << abs(pf->TrackerTrk()->DzCorrected(vtx)) << " "
705     << dr << endl;
706     }
707     if( pf->HasGsfTrk() ) {
708     if (abs(pf->GsfTrk()->DzCorrected(vtx)) > 0.2) continue;
709     if( ctrl.debug ) cout << "charged:: " << pf->PFType() << " " << pf->Pt() << " "
710     << abs(pf->GsfTrk()->DzCorrected(vtx)) << " "
711     << dr << endl;
712     }
713    
714    
715     tmpChargedIso += pf->Pt();
716     }
717    
718     //
719     // Gamma Iso
720     //
721     else if (abs(pf->PFType()) == PFCandidate::eGamma) {
722     tmpGammaIso += pf->Pt();
723     }
724    
725     //
726     // Other Neutral Iso Rings
727     //
728     else {
729     tmpNeutralHadronIso += pf->Pt();
730     }
731    
732     }
733    
734     }
735    
736     fChargedIso = mu->Pt() * min((tmpChargedIso)/mu->Pt(), 2.5);
737    
738    
739     double rho = 0;
740     if (!(isnan(fPUEnergyDensity->At(0)->Rho()) || isinf(fPUEnergyDensity->At(0)->Rho())))
741     rho = fPUEnergyDensity->At(0)->Rho();
742    
743    
744     fGammaIso = mu->Pt()*max(min((tmpGammaIso
745     -rho*muT.MuonEffectiveArea(muT.kMuGammaIso04,mu->Eta(),EffectiveAreaVersion))/mu->Pt()
746     ,2.5)
747     ,0.0);
748     fNeutralHadronIso = mu->Pt()*max(min((tmpNeutralHadronIso
749     -rho*muT.MuonEffectiveArea(muT.kMuNeutralIso04,
750     mu->Eta(),EffectiveAreaVersion))/mu->Pt()
751     , 2.5)
752     , 0.0);
753    
754     double pfIso = fChargedIso + fGammaIso + fNeutralHadronIso;
755    
756     return pfIso;
757     }
758    
759     //--------------------------------------------------------------------------------------------------
760     SelectionStatus muonIsoReferenceSelection(ControlFlags &ctrl,
761     const mithep::Muon * mu,
762     const mithep::Vertex & vtx,
763     const mithep::Array<mithep::PFCandidate> * fPFCandidates,
764     const mithep::Array<mithep::PileupEnergyDensity> * fPUEnergyDensity,
765     mithep::MuonTools::EMuonEffectiveAreaTarget EffectiveAreaVersion,
766     vector<const mithep::Muon*> muonsToVeto,
767     vector<const mithep::Electron*> electronsToVeto)
768     //--------------------------------------------------------------------------------------------------
769     {
770    
771     SelectionStatus status;
772    
773     double pfIso = muonPFIso04( ctrl, mu, vtx, fPFCandidates, fPUEnergyDensity,
774     EffectiveAreaVersion, muonsToVeto ,electronsToVeto );
775     bool pass = false;
776     if( (pfIso/mu->Pt()) < MUON_REFERENCE_PFISO_CUT ) pass = true;
777    
778     if( pass ) {
779     status.orStatus(SelectionStatus::LOOSEISO);
780     status.orStatus(SelectionStatus::TIGHTISO);
781     }
782     if(ctrl.debug) cout << "returning status : " << hex << status.getStatus() << dec << endl;
783     return status;
784    
785     }
786    
787    
788 khahn 1.4
789     //--------------------------------------------------------------------------------------------------
790     SelectionStatus electronIsoMVASelection(ControlFlags &ctrl,
791     const mithep::Electron * ele,
792     const mithep::Vertex & vtx,
793     const mithep::Array<mithep::PFCandidate> * fPFCandidates,
794     const mithep::Array<mithep::PileupEnergyDensity> * fPUEnergyDensity,
795     mithep::ElectronTools::EElectronEffectiveAreaTarget EffectiveAreaVersion,
796     vector<const mithep::Muon*> muonsToVeto,
797     vector<const mithep::Electron*> electronsToVeto)
798     //--------------------------------------------------------------------------------------------------
799     {
800    
801 khahn 1.7 if( ctrl.debug ) {
802     cout << "electronIsoMVASelection :: muons to veto " << endl;
803     for( int i=0; i<muonsToVeto.size(); i++ ) {
804     const mithep::Muon * vmu = muonsToVeto[i];
805     cout << "\tpt: " << vmu->Pt()
806     << "\teta: " << vmu->Eta()
807     << "\tphi: " << vmu->Phi()
808     << endl;
809     }
810     cout << "electronIsoMVASelection :: electrson to veto " << endl;
811     for( int i=0; i<electronsToVeto.size(); i++ ) {
812     const mithep::Electron * vel = electronsToVeto[i];
813     cout << "\tpt: " << vel->Pt()
814     << "\teta: " << vel->Eta()
815     << "\tphi: " << vel->Phi()
816     << "\ttrk: " << vel->TrackerTrk()
817     << endl;
818     }
819     }
820    
821 khahn 1.4 bool failiso=false;
822    
823     //
824     // tmp iso rings
825     //
826 khahn 1.7 Double_t tmpChargedIso_DR0p0To0p1 = 0;
827     Double_t tmpChargedIso_DR0p1To0p2 = 0;
828     Double_t tmpChargedIso_DR0p2To0p3 = 0;
829 khahn 1.4 Double_t tmpChargedIso_DR0p3To0p4 = 0;
830     Double_t tmpChargedIso_DR0p4To0p5 = 0;
831 khahn 1.7 Double_t tmpChargedIso_DR0p5To0p7 = 0;
832 khahn 1.4
833 khahn 1.7 Double_t tmpGammaIso_DR0p0To0p1 = 0;
834     Double_t tmpGammaIso_DR0p1To0p2 = 0;
835     Double_t tmpGammaIso_DR0p2To0p3 = 0;
836 khahn 1.4 Double_t tmpGammaIso_DR0p3To0p4 = 0;
837     Double_t tmpGammaIso_DR0p4To0p5 = 0;
838 khahn 1.7 Double_t tmpGammaIso_DR0p5To0p7 = 0;
839 khahn 1.4
840 khahn 1.7 Double_t tmpNeutralHadronIso_DR0p0To0p1 = 0;
841     Double_t tmpNeutralHadronIso_DR0p1To0p2 = 0;
842     Double_t tmpNeutralHadronIso_DR0p2To0p3 = 0;
843 khahn 1.4 Double_t tmpNeutralHadronIso_DR0p3To0p4 = 0;
844     Double_t tmpNeutralHadronIso_DR0p4To0p5 = 0;
845 khahn 1.7 Double_t tmpNeutralHadronIso_DR0p5To0p7 = 0;
846    
847 khahn 1.4
848    
849     //
850     // final rings for the MVA
851     //
852     Double_t fChargedIso_DR0p0To0p1;
853     Double_t fChargedIso_DR0p1To0p2;
854     Double_t fChargedIso_DR0p2To0p3;
855     Double_t fChargedIso_DR0p3To0p4;
856     Double_t fChargedIso_DR0p4To0p5;
857    
858     Double_t fGammaIso_DR0p0To0p1;
859     Double_t fGammaIso_DR0p1To0p2;
860     Double_t fGammaIso_DR0p2To0p3;
861     Double_t fGammaIso_DR0p3To0p4;
862     Double_t fGammaIso_DR0p4To0p5;
863    
864     Double_t fNeutralHadronIso_DR0p0To0p1;
865     Double_t fNeutralHadronIso_DR0p1To0p2;
866     Double_t fNeutralHadronIso_DR0p2To0p3;
867     Double_t fNeutralHadronIso_DR0p3To0p4;
868     Double_t fNeutralHadronIso_DR0p4To0p5;
869    
870    
871     //
872     //Loop over PF Candidates
873     //
874     for(int k=0; k<fPFCandidates->GetEntries(); ++k) {
875     const mithep::PFCandidate *pf = (mithep::PFCandidate*)((*fPFCandidates)[k]);
876     Double_t deta = (ele->Eta() - pf->Eta());
877     Double_t dphi = mithep::MathUtils::DeltaPhi(Double_t(ele->Phi()),Double_t(pf->Phi()));
878     Double_t dr = mithep::MathUtils::DeltaR(ele->Phi(),ele->Eta(), pf->Phi(), pf->Eta());
879 khahn 1.7 if (dr >= 0.5) continue;
880     if(ctrl.debug) {
881     cout << "pf :: type: " << pf->PFType() << "\tpt: " << pf->Pt();
882     if( pf->HasTrackerTrk() ) cout << "\tdZ: " << pf->TrackerTrk()->DzCorrected(vtx);
883     cout << endl;
884     }
885    
886 khahn 1.4
887 khahn 1.7 if ( (pf->HasTrackerTrk() && (pf->TrackerTrk() == ele->TrackerTrk())) ||
888     (pf->HasGsfTrk() && (pf->GsfTrk() == ele->GsfTrk()))) continue;
889    
890 khahn 1.4
891     //
892     // Lepton Footprint Removal
893     //
894     Bool_t IsLeptonFootprint = kFALSE;
895     if (dr < 1.0) {
896    
897     //
898     // Check for electrons
899     //
900     for (Int_t q=0; q < electronsToVeto.size(); ++q) {
901     const mithep::Electron *tmpele = electronsToVeto[q];
902 khahn 1.7 // 4l electron
903     if( pf->HasTrackerTrk() ) {
904     if( pf->TrackerTrk() == tmpele->TrackerTrk() ) {
905     if( ctrl.debug) cout << "\tcharged tktrk, matches 4L ele ..." << endl;
906     IsLeptonFootprint = kTRUE;
907     }
908     }
909     if( pf->HasGsfTrk() ) {
910     if( pf->GsfTrk() == tmpele->GsfTrk() ) {
911     if( ctrl.debug) cout << "\tcharged gsftrk, matches 4L ele ..." << endl;
912     IsLeptonFootprint = kTRUE;
913     }
914     }
915 khahn 1.4 // PF charged
916     if (pf->Charge() != 0 && fabs(tmpele->SCluster()->Eta()) > 1.479
917 khahn 1.7 && mithep::MathUtils::DeltaR(tmpele->Phi(),tmpele->Eta(), pf->Phi(), pf->Eta()) < 0.015) {
918     if( ctrl.debug) cout << "\tcharged trk, dR matches 4L ele ..." << endl;
919 khahn 1.4 IsLeptonFootprint = kTRUE;
920 khahn 1.7 }
921 khahn 1.4 // PF gamma
922 khahn 1.7 if (abs(pf->PFType()) == PFCandidate::eGamma && fabs(tmpele->SCluster()->Eta()) > 1.479
923     && mithep::MathUtils::DeltaR(tmpele->Phi(),tmpele->Eta(), pf->Phi(), pf->Eta()) < 0.08) {
924     if( ctrl.debug) cout << "\tPF gamma, matches 4L ele ..." << endl;
925 khahn 1.4 IsLeptonFootprint = kTRUE;
926 khahn 1.7 }
927 khahn 1.4 } // loop over electrons
928    
929     //
930     // Check for muons
931     //
932     for (Int_t q=0; q < muonsToVeto.size(); ++q) {
933     const mithep::Muon *tmpmu = muonsToVeto[q];
934 khahn 1.7 // 4l muon
935     if( pf->HasTrackerTrk() ) {
936     if (pf->TrackerTrk() == tmpmu->TrackerTrk() ){
937     if( ctrl.debug) cout << "\tmatches 4L mu ..." << endl;
938     IsLeptonFootprint = kTRUE;
939     }
940     }
941 khahn 1.4 // PF charged
942 khahn 1.7 if (pf->Charge() != 0 && mithep::MathUtils::DeltaR(tmpmu->Phi(),tmpmu->Eta(), pf->Phi(), pf->Eta()) < 0.01) {
943     if( ctrl.debug) cout << "\tcharged trk, dR matches 4L mu ..." << endl;
944 khahn 1.4 IsLeptonFootprint = kTRUE;
945 khahn 1.7 }
946 khahn 1.4 } // loop over muons
947    
948    
949     if (IsLeptonFootprint)
950     continue;
951    
952     //
953     // Charged Iso Rings
954     //
955 khahn 1.7 if (pf->Charge() != 0 && (pf->HasTrackerTrk()||pf->HasGsfTrk()) ) {
956    
957     if( pf->HasTrackerTrk() )
958     if (abs(pf->TrackerTrk()->DzCorrected(vtx)) > 0.2) continue;
959     if( pf->HasGsfTrk() )
960     if (abs(pf->GsfTrk()->DzCorrected(vtx)) > 0.2) continue;
961 khahn 1.4
962     // Veto any PFmuon, or PFEle
963 khahn 1.7 if (abs(pf->PFType()) == PFCandidate::eElectron || abs(pf->PFType()) == PFCandidate::eMuon) continue;
964 khahn 1.4
965     // Footprint Veto
966 khahn 1.7 if (fabs(ele->SCluster()->Eta()) > 1.479 && dr < 0.015) continue;
967    
968     if( ctrl.debug) cout << "charged:: pt: " << pf->Pt()
969     << "\ttype: " << pf->PFType()
970     << "\ttrk: " << pf->TrackerTrk() << endl;
971    
972     if (dr < 0.1) tmpChargedIso_DR0p0To0p1 += pf->Pt();
973     if (dr >= 0.1 && dr < 0.2) tmpChargedIso_DR0p1To0p2 += pf->Pt();
974     if (dr >= 0.2 && dr < 0.3) tmpChargedIso_DR0p2To0p3 += pf->Pt();
975     if (dr >= 0.3 && dr < 0.4) tmpChargedIso_DR0p3To0p4 += pf->Pt();
976     if (dr >= 0.4 && dr < 0.5) tmpChargedIso_DR0p4To0p5 += pf->Pt();
977     if (dr >= 0.5 && dr < 0.7) tmpChargedIso_DR0p5To0p7 += pf->Pt();
978 khahn 1.4
979     }
980    
981     //
982     // Gamma Iso Rings
983     //
984 khahn 1.7 else if (abs(pf->PFType()) == PFCandidate::eGamma) {
985    
986     if (fabs(ele->SCluster()->Eta()) > 1.479) {
987     if (mithep::MathUtils::DeltaR(ele->Phi(),ele->Eta(), pf->Phi(), pf->Eta()) < 0.08) continue;
988     }
989    
990     if( ctrl.debug) cout << "gamma:: " << pf->Pt() << " "
991     << dr << endl;
992    
993     if (dr < 0.1) tmpGammaIso_DR0p0To0p1 += pf->Pt();
994     if (dr >= 0.1 && dr < 0.2) tmpGammaIso_DR0p1To0p2 += pf->Pt();
995     if (dr >= 0.2 && dr < 0.3) tmpGammaIso_DR0p2To0p3 += pf->Pt();
996     if (dr >= 0.3 && dr < 0.4) tmpGammaIso_DR0p3To0p4 += pf->Pt();
997     if (dr >= 0.4 && dr < 0.5) tmpGammaIso_DR0p4To0p5 += pf->Pt();
998     if (dr >= 0.5 && dr < 0.7) tmpGammaIso_DR0p5To0p7 += pf->Pt();
999 khahn 1.4
1000     }
1001    
1002     //
1003     // Other Neutral Iso Rings
1004     //
1005     else {
1006 khahn 1.7 if( ctrl.debug) cout << "neutral:: " << pf->Pt() << " "
1007     << dr << endl;
1008     if (dr < 0.1) tmpNeutralHadronIso_DR0p0To0p1 += pf->Pt();
1009     if (dr >= 0.1 && dr < 0.2) tmpNeutralHadronIso_DR0p1To0p2 += pf->Pt();
1010     if (dr >= 0.2 && dr < 0.3) tmpNeutralHadronIso_DR0p2To0p3 += pf->Pt();
1011     if (dr >= 0.3 && dr < 0.4) tmpNeutralHadronIso_DR0p3To0p4 += pf->Pt();
1012     if (dr >= 0.4 && dr < 0.5) tmpNeutralHadronIso_DR0p4To0p5 += pf->Pt();
1013     if (dr >= 0.5 && dr < 0.7) tmpNeutralHadronIso_DR0p5To0p7 += pf->Pt();
1014 khahn 1.4 }
1015    
1016     }
1017    
1018     }
1019    
1020 khahn 1.7 fChargedIso_DR0p0To0p1 = min((tmpChargedIso_DR0p0To0p1)/ele->Pt(), 2.5);
1021     fChargedIso_DR0p1To0p2 = min((tmpChargedIso_DR0p1To0p2)/ele->Pt(), 2.5);
1022     fChargedIso_DR0p2To0p3 = min((tmpChargedIso_DR0p2To0p3)/ele->Pt(), 2.5);
1023 khahn 1.4 fChargedIso_DR0p3To0p4 = min((tmpChargedIso_DR0p3To0p4)/ele->Pt(), 2.5);
1024     fChargedIso_DR0p4To0p5 = min((tmpChargedIso_DR0p4To0p5)/ele->Pt(), 2.5);
1025    
1026     double rho = 0;
1027     if (!(isnan(fPUEnergyDensity->At(0)->Rho()) || isinf(fPUEnergyDensity->At(0)->Rho())))
1028     rho = fPUEnergyDensity->At(0)->Rho();
1029 khahn 1.7
1030     if( ctrl.debug) {
1031     cout << "RHO: " << rho << endl;
1032     cout << "eta: " << ele->SCluster()->Eta() << endl;
1033     cout << "target: " << EffectiveAreaVersion << endl;
1034     cout << "effA 0-1: " << eleT.ElectronEffectiveArea(eleT.kEleNeutralHadronIsoDR0p0To0p1,
1035     ele->SCluster()->Eta(),
1036     EffectiveAreaVersion)
1037     << endl;
1038     cout << "effA 1-2: " << eleT.ElectronEffectiveArea(eleT.kEleNeutralHadronIsoDR0p1To0p2,
1039     ele->SCluster()->Eta(),
1040     EffectiveAreaVersion)
1041     << endl;
1042     cout << "effA 2-3: " << eleT.ElectronEffectiveArea(eleT.kEleNeutralHadronIsoDR0p2To0p3,
1043     ele->SCluster()->Eta(),
1044     EffectiveAreaVersion)
1045     << endl;
1046     cout << "effA 3-4: " << eleT.ElectronEffectiveArea(eleT.kEleNeutralHadronIsoDR0p3To0p4,
1047     ele->SCluster()->Eta(),
1048     EffectiveAreaVersion)
1049     << endl;
1050     }
1051    
1052     fGammaIso_DR0p0To0p1 = max(min((tmpGammaIso_DR0p0To0p1
1053 khahn 1.4 -rho*eleT.ElectronEffectiveArea(eleT.kEleGammaIsoDR0p0To0p1,
1054 khahn 1.7 ele->SCluster()->Eta(),
1055 khahn 1.4 EffectiveAreaVersion))/ele->Pt()
1056     ,2.5)
1057     ,0.0);
1058 khahn 1.7 fGammaIso_DR0p1To0p2 = max(min((tmpGammaIso_DR0p1To0p2
1059 khahn 1.4 -rho*eleT.ElectronEffectiveArea(eleT.kEleGammaIsoDR0p1To0p2,
1060 khahn 1.7 ele->SCluster()->Eta(),
1061 khahn 1.4 EffectiveAreaVersion))/ele->Pt()
1062     ,2.5)
1063     ,0.0);
1064 khahn 1.7 fGammaIso_DR0p2To0p3 = max(min((tmpGammaIso_DR0p2To0p3
1065 khahn 1.4 -rho*eleT.ElectronEffectiveArea(eleT.kEleGammaIsoDR0p2To0p3,
1066 khahn 1.7 ele->SCluster()->Eta()
1067 khahn 1.4 ,EffectiveAreaVersion))/ele->Pt()
1068     ,2.5)
1069     ,0.0);
1070     fGammaIso_DR0p3To0p4 = max(min((tmpGammaIso_DR0p3To0p4
1071     -rho*eleT.ElectronEffectiveArea(eleT.kEleGammaIsoDR0p3To0p4,
1072 khahn 1.7 ele->SCluster()->Eta(),
1073 khahn 1.4 EffectiveAreaVersion))/ele->Pt()
1074     ,2.5)
1075     ,0.0);
1076     fGammaIso_DR0p4To0p5 = max(min((tmpGammaIso_DR0p4To0p5
1077     -rho*eleT.ElectronEffectiveArea(eleT.kEleGammaIsoDR0p4To0p5,
1078 khahn 1.7 ele->SCluster()->Eta(),
1079 khahn 1.4 EffectiveAreaVersion))/ele->Pt()
1080     ,2.5)
1081     ,0.0);
1082    
1083    
1084 khahn 1.7 fNeutralHadronIso_DR0p0To0p1 = max(min((tmpNeutralHadronIso_DR0p0To0p1
1085 khahn 1.4 -rho*eleT.ElectronEffectiveArea(eleT.kEleNeutralHadronIsoDR0p0To0p1,
1086 khahn 1.7 ele->SCluster()->Eta(),EffectiveAreaVersion))/ele->Pt()
1087 khahn 1.4 , 2.5)
1088     , 0.0);
1089 khahn 1.7 fNeutralHadronIso_DR0p1To0p2 = max(min((tmpNeutralHadronIso_DR0p1To0p2
1090 khahn 1.4 -rho*eleT.ElectronEffectiveArea(eleT.kEleNeutralHadronIsoDR0p1To0p2,
1091 khahn 1.7 ele->SCluster()->Eta(),EffectiveAreaVersion))/ele->Pt()
1092 khahn 1.4 , 2.5)
1093     , 0.0);
1094 khahn 1.7 fNeutralHadronIso_DR0p2To0p3 = max(min((tmpNeutralHadronIso_DR0p2To0p3
1095 khahn 1.4 -rho*eleT.ElectronEffectiveArea(eleT.kEleNeutralHadronIsoDR0p2To0p3,
1096 khahn 1.7 ele->SCluster()->Eta(),EffectiveAreaVersion))/ele->Pt()
1097 khahn 1.4 , 2.5)
1098     , 0.0);
1099     fNeutralHadronIso_DR0p3To0p4 = max(min((tmpNeutralHadronIso_DR0p3To0p4
1100     -rho*eleT.ElectronEffectiveArea(eleT.kEleNeutralHadronIsoDR0p3To0p4,
1101 khahn 1.7 ele->SCluster()->Eta(), EffectiveAreaVersion))/ele->Pt()
1102 khahn 1.4 , 2.5)
1103     , 0.0);
1104     fNeutralHadronIso_DR0p4To0p5 = max(min((tmpNeutralHadronIso_DR0p4To0p5
1105     -rho*eleT.ElectronEffectiveArea(eleT.kEleNeutralHadronIsoDR0p4To0p5,
1106 khahn 1.7 ele->SCluster()->Eta(), EffectiveAreaVersion))/ele->Pt()
1107 khahn 1.4 , 2.5)
1108     , 0.0);
1109    
1110     double mvaval = eleIsoMVA->MVAValue_IsoRings( ele->Pt(),
1111 khahn 1.7 ele->SCluster()->Eta(),
1112     fChargedIso_DR0p0To0p1,
1113     fChargedIso_DR0p1To0p2,
1114     fChargedIso_DR0p2To0p3,
1115     fChargedIso_DR0p3To0p4,
1116     fChargedIso_DR0p4To0p5,
1117     fGammaIso_DR0p0To0p1,
1118     fGammaIso_DR0p1To0p2,
1119     fGammaIso_DR0p2To0p3,
1120     fGammaIso_DR0p3To0p4,
1121     fGammaIso_DR0p4To0p5,
1122     fNeutralHadronIso_DR0p0To0p1,
1123     fNeutralHadronIso_DR0p1To0p2,
1124     fNeutralHadronIso_DR0p2To0p3,
1125     fNeutralHadronIso_DR0p3To0p4,
1126     fNeutralHadronIso_DR0p4To0p5,
1127     ctrl.debug);
1128 khahn 1.4
1129     SelectionStatus status;
1130     bool pass = false;
1131    
1132     Int_t subdet = 0;
1133     if (fabs(ele->SCluster()->Eta()) < 0.8) subdet = 0;
1134     else if (fabs(ele->SCluster()->Eta()) < 1.479) subdet = 1;
1135     else subdet = 2;
1136     Int_t ptBin = 0;
1137     if (ele->Pt() > 10.0) ptBin = 1;
1138    
1139     Int_t MVABin = -1;
1140     if (subdet == 0 && ptBin == 0) MVABin = 0;
1141     if (subdet == 1 && ptBin == 0) MVABin = 1;
1142     if (subdet == 2 && ptBin == 0) MVABin = 2;
1143     if (subdet == 0 && ptBin == 1) MVABin = 3;
1144     if (subdet == 1 && ptBin == 1) MVABin = 4;
1145     if (subdet == 2 && ptBin == 1) MVABin = 5;
1146    
1147 khahn 1.10 pass = false;
1148     if( MVABin == 0 && mvaval > ELECTRON_LOOSE_ISOMVA_CUT_BIN0 ) pass = true;
1149     if( MVABin == 1 && mvaval > ELECTRON_LOOSE_ISOMVA_CUT_BIN1 ) pass = true;
1150     if( MVABin == 2 && mvaval > ELECTRON_LOOSE_ISOMVA_CUT_BIN2 ) pass = true;
1151     if( MVABin == 3 && mvaval > ELECTRON_LOOSE_ISOMVA_CUT_BIN3 ) pass = true;
1152     if( MVABin == 4 && mvaval > ELECTRON_LOOSE_ISOMVA_CUT_BIN4 ) pass = true;
1153     if( MVABin == 5 && mvaval > ELECTRON_LOOSE_ISOMVA_CUT_BIN5 ) pass = true;
1154     if( pass ) status.orStatus(SelectionStatus::LOOSEISO);
1155    
1156     pass = false;
1157     if( MVABin == 0 && mvaval > ELECTRON_TIGHT_ISOMVA_CUT_BIN0 ) pass = true;
1158     if( MVABin == 1 && mvaval > ELECTRON_TIGHT_ISOMVA_CUT_BIN1 ) pass = true;
1159     if( MVABin == 2 && mvaval > ELECTRON_TIGHT_ISOMVA_CUT_BIN2 ) pass = true;
1160     if( MVABin == 3 && mvaval > ELECTRON_TIGHT_ISOMVA_CUT_BIN3 ) pass = true;
1161     if( MVABin == 4 && mvaval > ELECTRON_TIGHT_ISOMVA_CUT_BIN4 ) pass = true;
1162     if( MVABin == 5 && mvaval > ELECTRON_TIGHT_ISOMVA_CUT_BIN5 ) pass = true;
1163     if( pass ) status.orStatus(SelectionStatus::TIGHTISO);
1164 khahn 1.7
1165 khahn 1.4 if(ctrl.debug) cout << "returning status : " << hex << status.getStatus() << dec << endl;
1166     return status;
1167 khahn 1.10
1168 khahn 1.4 }
1169    
1170    
1171     //--------------------------------------------------------------------------------------------------
1172     void initElectronIsoMVA() {
1173     //--------------------------------------------------------------------------------------------------
1174     eleIsoMVA = new mithep::ElectronIDMVA();
1175     vector<string> weightFiles;
1176     weightFiles.push_back("../MitPhysics/data/ElectronMVAWeights/ElectronIso_BDTG_V0_BarrelPt5To10.weights.xml");
1177     weightFiles.push_back("../MitPhysics/data/ElectronMVAWeights/ElectronIso_BDTG_V0_EndcapPt5To10.weights.xml");
1178     weightFiles.push_back("../MitPhysics/data/ElectronMVAWeights/ElectronIso_BDTG_V0_BarrelPt10ToInf.weights.xml");
1179     weightFiles.push_back("../MitPhysics/data/ElectronMVAWeights/ElectronIso_BDTG_V0_EndcapPt10ToInf.weights.xml");
1180     eleIsoMVA->Initialize( "ElectronIsoMVA",
1181     mithep::ElectronIDMVA::kIsoRingsV0,
1182     kTRUE, weightFiles);
1183     }
1184 khahn 1.10
1185    
1186    
1187     //--------------------------------------------------------------------------------------------------
1188     float electronPFIso04(ControlFlags &ctrl,
1189     const mithep::Electron * ele,
1190     const mithep::Vertex & vtx,
1191     const mithep::Array<mithep::PFCandidate> * fPFCandidates,
1192     const mithep::Array<mithep::PileupEnergyDensity> * fPUEnergyDensity,
1193     mithep::ElectronTools::EElectronEffectiveAreaTarget EffectiveAreaVersion,
1194     vector<const mithep::Muon*> muonsToVeto,
1195     vector<const mithep::Electron*> electronsToVeto)
1196     //--------------------------------------------------------------------------------------------------
1197     {
1198    
1199     if( ctrl.debug ) {
1200     cout << "electronIsoMVASelection :: muons to veto " << endl;
1201     for( int i=0; i<muonsToVeto.size(); i++ ) {
1202     const mithep::Muon * vmu = muonsToVeto[i];
1203     cout << "\tpt: " << vmu->Pt()
1204     << "\teta: " << vmu->Eta()
1205     << "\tphi: " << vmu->Phi()
1206     << endl;
1207     }
1208     cout << "electronIsoMVASelection :: electrson to veto " << endl;
1209     for( int i=0; i<electronsToVeto.size(); i++ ) {
1210     const mithep::Electron * vel = electronsToVeto[i];
1211     cout << "\tpt: " << vel->Pt()
1212     << "\teta: " << vel->Eta()
1213     << "\tphi: " << vel->Phi()
1214     << "\ttrk: " << vel->TrackerTrk()
1215     << endl;
1216     }
1217     }
1218    
1219     bool failiso=false;
1220    
1221     //
1222     // tmp iso
1223     //
1224     Double_t tmpChargedIso = 0;
1225     Double_t tmpGammaIso = 0;
1226     Double_t tmpNeutralHadronIso = 0;
1227    
1228     //
1229     // final iso
1230     //
1231     Double_t fChargedIso;
1232     Double_t fGammaIso;
1233     Double_t fNeutralHadronIso;
1234    
1235    
1236     //
1237     //Loop over PF Candidates
1238     //
1239     for(int k=0; k<fPFCandidates->GetEntries(); ++k) {
1240     const mithep::PFCandidate *pf = (mithep::PFCandidate*)((*fPFCandidates)[k]);
1241     Double_t deta = (ele->Eta() - pf->Eta());
1242     Double_t dphi = mithep::MathUtils::DeltaPhi(Double_t(ele->Phi()),Double_t(pf->Phi()));
1243     Double_t dr = mithep::MathUtils::DeltaR(ele->Phi(),ele->Eta(), pf->Phi(), pf->Eta());
1244     if (dr >= 0.4) continue;
1245     if(ctrl.debug) {
1246     cout << "pf :: type: " << pf->PFType() << "\tpt: " << pf->Pt();
1247     if( pf->HasTrackerTrk() ) cout << "\tdZ: " << pf->TrackerTrk()->DzCorrected(vtx);
1248     cout << endl;
1249     }
1250    
1251    
1252     if ( (pf->HasTrackerTrk() && (pf->TrackerTrk() == ele->TrackerTrk())) ||
1253     (pf->HasGsfTrk() && (pf->GsfTrk() == ele->GsfTrk()))) continue;
1254    
1255    
1256     //
1257     // Lepton Footprint Removal
1258     //
1259     Bool_t IsLeptonFootprint = kFALSE;
1260     if (dr < 1.0) {
1261    
1262     //
1263     // Check for electrons
1264     //
1265     for (Int_t q=0; q < electronsToVeto.size(); ++q) {
1266     const mithep::Electron *tmpele = electronsToVeto[q];
1267     // 4l electron
1268     if( pf->HasTrackerTrk() ) {
1269     if( pf->TrackerTrk() == tmpele->TrackerTrk() ) {
1270     if( ctrl.debug) cout << "\tcharged tktrk, matches 4L ele ..." << endl;
1271     IsLeptonFootprint = kTRUE;
1272     }
1273     }
1274     if( pf->HasGsfTrk() ) {
1275     if( pf->GsfTrk() == tmpele->GsfTrk() ) {
1276     if( ctrl.debug) cout << "\tcharged gsftrk, matches 4L ele ..." << endl;
1277     IsLeptonFootprint = kTRUE;
1278     }
1279     }
1280     // PF charged
1281     if (pf->Charge() != 0 && fabs(tmpele->SCluster()->Eta()) > 1.479
1282     && mithep::MathUtils::DeltaR(tmpele->Phi(),tmpele->Eta(), pf->Phi(), pf->Eta()) < 0.015) {
1283     if( ctrl.debug) cout << "\tcharged trk, dR matches 4L ele ..." << endl;
1284     IsLeptonFootprint = kTRUE;
1285     }
1286     // PF gamma
1287     if (abs(pf->PFType()) == PFCandidate::eGamma && fabs(tmpele->SCluster()->Eta()) > 1.479
1288     && mithep::MathUtils::DeltaR(tmpele->Phi(),tmpele->Eta(), pf->Phi(), pf->Eta()) < 0.08) {
1289     if( ctrl.debug) cout << "\tPF gamma, matches 4L ele ..." << endl;
1290     IsLeptonFootprint = kTRUE;
1291     }
1292     } // loop over electrons
1293    
1294     //
1295     // Check for muons
1296     //
1297     for (Int_t q=0; q < muonsToVeto.size(); ++q) {
1298     const mithep::Muon *tmpmu = muonsToVeto[q];
1299     // 4l muon
1300     if( pf->HasTrackerTrk() ) {
1301     if (pf->TrackerTrk() == tmpmu->TrackerTrk() ){
1302     if( ctrl.debug) cout << "\tmatches 4L mu ..." << endl;
1303     IsLeptonFootprint = kTRUE;
1304     }
1305     }
1306     // PF charged
1307     if (pf->Charge() != 0 && mithep::MathUtils::DeltaR(tmpmu->Phi(),tmpmu->Eta(), pf->Phi(), pf->Eta()) < 0.01) {
1308     if( ctrl.debug) cout << "\tcharged trk, dR matches 4L mu ..." << endl;
1309     IsLeptonFootprint = kTRUE;
1310     }
1311     } // loop over muons
1312    
1313    
1314     if (IsLeptonFootprint)
1315     continue;
1316    
1317     //
1318     // Charged Iso Rings
1319     //
1320     if (pf->Charge() != 0 && (pf->HasTrackerTrk()||pf->HasGsfTrk()) ) {
1321    
1322     if( pf->HasTrackerTrk() )
1323     if (abs(pf->TrackerTrk()->DzCorrected(vtx)) > 0.2) continue;
1324     if( pf->HasGsfTrk() )
1325     if (abs(pf->GsfTrk()->DzCorrected(vtx)) > 0.2) continue;
1326    
1327     // Veto any PFmuon, or PFEle
1328     if (abs(pf->PFType()) == PFCandidate::eElectron || abs(pf->PFType()) == PFCandidate::eMuon) continue;
1329    
1330     // Footprint Veto
1331     if (fabs(ele->SCluster()->Eta()) > 1.479 && dr < 0.015) continue;
1332    
1333     if( ctrl.debug) cout << "charged:: pt: " << pf->Pt()
1334     << "\ttype: " << pf->PFType()
1335     << "\ttrk: " << pf->TrackerTrk() << endl;
1336    
1337     tmpChargedIso += pf->Pt();
1338     }
1339    
1340     //
1341     // Gamma Iso
1342     //
1343     else if (abs(pf->PFType()) == PFCandidate::eGamma) {
1344    
1345     if (fabs(ele->SCluster()->Eta()) > 1.479) {
1346     if (mithep::MathUtils::DeltaR(ele->Phi(),ele->Eta(), pf->Phi(), pf->Eta()) < 0.08) continue;
1347     }
1348     if( ctrl.debug) cout << "gamma:: " << pf->Pt() << " "
1349     << dr << endl;
1350     tmpGammaIso += pf->Pt();
1351     }
1352    
1353     //
1354     // Neutral Iso
1355     //
1356     else {
1357     if( ctrl.debug) cout << "neutral:: " << pf->Pt() << " "
1358     << dr << endl;
1359     tmpNeutralHadronIso += pf->Pt();
1360     }
1361    
1362     }
1363    
1364     }
1365    
1366     fChargedIso = ele->Pt()*min((tmpChargedIso)/ele->Pt(), 2.5);
1367    
1368     double rho = 0;
1369     if (!(isnan(fPUEnergyDensity->At(0)->Rho()) || isinf(fPUEnergyDensity->At(0)->Rho())))
1370     rho = fPUEnergyDensity->At(0)->Rho();
1371    
1372     if( ctrl.debug) {
1373     cout << "RHO: " << rho << endl;
1374     cout << "eta: " << ele->SCluster()->Eta() << endl;
1375     cout << "target: " << EffectiveAreaVersion << endl;
1376     cout << "effA 0-1: " << eleT.ElectronEffectiveArea(eleT.kEleNeutralHadronIsoDR0p0To0p1,
1377     ele->SCluster()->Eta(),
1378     EffectiveAreaVersion)
1379     << endl;
1380     cout << "effA 1-2: " << eleT.ElectronEffectiveArea(eleT.kEleNeutralHadronIsoDR0p1To0p2,
1381     ele->SCluster()->Eta(),
1382     EffectiveAreaVersion)
1383     << endl;
1384     cout << "effA 2-3: " << eleT.ElectronEffectiveArea(eleT.kEleNeutralHadronIsoDR0p2To0p3,
1385     ele->SCluster()->Eta(),
1386     EffectiveAreaVersion)
1387     << endl;
1388     cout << "effA 3-4: " << eleT.ElectronEffectiveArea(eleT.kEleNeutralHadronIsoDR0p3To0p4,
1389     ele->SCluster()->Eta(),
1390     EffectiveAreaVersion)
1391     << endl;
1392     }
1393    
1394     fGammaIso = ele->Pt()*max(min((tmpGammaIso
1395     -rho*eleT.ElectronEffectiveArea(eleT.kEleGammaIso04,
1396     ele->SCluster()->Eta(),
1397     EffectiveAreaVersion))/ele->Pt()
1398     ,2.5)
1399     ,0.0);
1400     fNeutralHadronIso = ele->Pt()*max(min((tmpNeutralHadronIso
1401     -rho*eleT.ElectronEffectiveArea(eleT.kEleNeutralHadronIso04,
1402     ele->SCluster()->Eta(),EffectiveAreaVersion))/ele->Pt()
1403     , 2.5)
1404     , 0.0);
1405    
1406     double pfiso = fChargedIso + fGammaIso + fNeutralHadronIso;
1407    
1408     return pfiso;
1409     }
1410    
1411     //--------------------------------------------------------------------------------------------------
1412     SelectionStatus electronIsoReferenceSelection(ControlFlags &ctrl,
1413     const mithep::Electron * ele,
1414     const mithep::Vertex & vtx,
1415     const mithep::Array<mithep::PFCandidate> * fPFCandidates,
1416     const mithep::Array<mithep::PileupEnergyDensity> * fPUEnergyDensity,
1417     mithep::ElectronTools::EElectronEffectiveAreaTarget EffectiveAreaVersion,
1418     vector<const mithep::Muon*> muonsToVeto,
1419     vector<const mithep::Electron*> electronsToVeto)
1420     //--------------------------------------------------------------------------------------------------
1421     {
1422    
1423     SelectionStatus status;
1424    
1425     double pfIso = electronPFIso04( ctrl, ele, vtx, fPFCandidates, fPUEnergyDensity,
1426     EffectiveAreaVersion, muonsToVeto ,electronsToVeto );
1427     bool pass = false;
1428     if( (pfIso/ele->Pt()) < ELECTRON_REFERENCE_PFISO_CUT ) pass = true;
1429    
1430     if( pass ) {
1431     status.orStatus(SelectionStatus::LOOSEISO);
1432     status.orStatus(SelectionStatus::TIGHTISO);
1433     }
1434     if(ctrl.debug) cout << "returning status : " << hex << status.getStatus() << dec << endl;
1435     return status;
1436    
1437     }