ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitHzz4l/LeptonSelection/src/IsolationSelection.cc
(Generate patch)

Comparing UserCode/MitHzz4l/LeptonSelection/src/IsolationSelection.cc (file contents):
Revision 1.2 by khahn, Fri Feb 17 14:49:00 2012 UTC vs.
Revision 1.10 by khahn, Sat May 5 13:09:05 2012 UTC

# Line 3 | Line 3
3   #include "IsolationSelection.h"
4   #include "IsolationSelectionDefs.h"
5  
6 < bool pairwiseIsoSelection( ControlFlags &ctrl, vector<SimpleLepton> &lepvec, float rho ) {
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 > //--------------------------------------------------------------------------------------------------
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  
105    bool passiso=true;
106  
# Line 75 | Line 172 | bool pairwiseIsoSelection( ControlFlags
172    return passiso;
173   }
174  
175 <
176 < SelectionStatus passMuonIsoSelection( ControlFlags &ctrl, const mithep::TMuon * mu ) {
177 <
178 <  float reliso = mu->pfIso03/mu->pt;
179 <  bool isEB = (fabs(mu->eta) < 1.479 ? 1 : 0 );  
175 > //--------------------------------------------------------------------------------------------------
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    bool failiso = false;
185 <  if( isEB && mu->pt > 20 && reliso > PFISO_MU_LOOSE_EB_HIGHPT ) {  
185 >  if( isEB && mu->Pt() > 20 && reliso > PFISO_MU_LOOSE_EB_HIGHPT ) {  
186      failiso = true;
187    }
188 <  if( isEB && mu->pt < 20 && reliso > PFISO_MU_LOOSE_EB_LOWPT ) {
188 >  if( isEB && mu->Pt() < 20 && reliso > PFISO_MU_LOOSE_EB_LOWPT ) {
189      failiso = true;
190    }
191 <  if( !(isEB) && mu->pt > 20 && reliso > PFISO_MU_LOOSE_EE_HIGHPT ) {
191 >  if( !(isEB) && mu->Pt() > 20 && reliso > PFISO_MU_LOOSE_EE_HIGHPT ) {
192      failiso = true;
193    }
194 <  if( !(isEB) && mu->pt < 20 && reliso > PFISO_MU_LOOSE_EE_LOWPT ) {
194 >  if( !(isEB) && mu->Pt() < 20 && reliso > PFISO_MU_LOOSE_EE_LOWPT ) {
195      failiso = true;
196    }
197  
# Line 101 | Line 202 | SelectionStatus passMuonIsoSelection( Co
202  
203   };
204  
205 <
206 < SelectionStatus failEleIso(ControlFlags &ctrl, const mithep::TElectron * ele) {
205 > //--------------------------------------------------------------------------------------------------
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  
213    bool failiso=false;
214  
215 <  float reliso = ele->pfIso04/ele->pt;
216 <  bool isEB = (fabs(ele->eta) < 1.479 ? 1 : 0 );  
217 <  if( isEB && ele->pt > 20 && reliso > PFISO_ELE_LOOSE_EB_HIGHPT ) {
215 >  float reliso = computePFEleIso(ele,fVertex,fPFCandidates,0.4)/ele->Pt();
216 >
217 >  if( ele->IsEB() && ele->Pt() > 20 && reliso > PFISO_ELE_LOOSE_EB_HIGHPT ) {
218      failiso = true;
219    }
220 <  if( isEB && ele->pt < 20 && reliso > PFISO_ELE_LOOSE_EB_LOWPT ) {
220 >  if( ele->IsEB() && ele->Pt() < 20 && reliso > PFISO_ELE_LOOSE_EB_LOWPT ) {
221      failiso = true;
222    }
223    if(ctrl.debug) cout << "before iso check ..." << endl;
224 <  if( !(isEB) && ele->pt > 20 && reliso > PFISO_ELE_LOOSE_EE_HIGHPT ) {
224 >  if( !(ele->IsEB()) && ele->Pt() > 20 && reliso > PFISO_ELE_LOOSE_EE_HIGHPT ) {
225      if(ctrl.debug) cout << "\tit fails ..." << endl;
226      failiso = true;
227    }
228 <  if( !(isEB) && ele->pt < 20 && reliso > PFISO_ELE_LOOSE_EE_LOWPT ) {
228 >  if( !(ele->IsEB()) && ele->Pt() < 20 && reliso > PFISO_ELE_LOOSE_EE_LOWPT ) {
229      failiso = true;
230    }
231  
232    SelectionStatus status;
233    if( !failiso ) {
234 <    status.setStatus(SelectionStatus::LOOSEISO);
235 <    status.setStatus(SelectionStatus::TIGHTISO);
234 >    status.orStatus(SelectionStatus::LOOSEISO);
235 >    status.orStatus(SelectionStatus::TIGHTISO);
236 >  }
237 >  if(ctrl.debug) cout << "returning status : " << hex << status.getStatus() << dec << endl;
238 >  return status;
239 >
240 > }
241 >
242 >
243 > bool noIso(ControlFlags &, vector<SimpleLepton> &, float rho) {
244 >
245 >        return true;
246 > }
247 >
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 >  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 >  bool failiso=false;
279 >
280 >  //
281 >  // tmp iso rings
282 >  //
283 >  Double_t tmpChargedIso_DR0p0To0p1  = 0;
284 >  Double_t tmpChargedIso_DR0p1To0p2  = 0;
285 >  Double_t tmpChargedIso_DR0p2To0p3  = 0;
286 >  Double_t tmpChargedIso_DR0p3To0p4  = 0;
287 >  Double_t tmpChargedIso_DR0p4To0p5  = 0;
288 >  Double_t tmpChargedIso_DR0p5To0p7  = 0;
289 >
290 >  Double_t tmpGammaIso_DR0p0To0p1  = 0;
291 >  Double_t tmpGammaIso_DR0p1To0p2  = 0;
292 >  Double_t tmpGammaIso_DR0p2To0p3  = 0;
293 >  Double_t tmpGammaIso_DR0p3To0p4  = 0;
294 >  Double_t tmpGammaIso_DR0p4To0p5  = 0;
295 >  Double_t tmpGammaIso_DR0p5To0p7  = 0;
296 >
297 >  Double_t tmpNeutralHadronIso_DR0p0To0p1  = 0;
298 >  Double_t tmpNeutralHadronIso_DR0p1To0p2  = 0;
299 >  Double_t tmpNeutralHadronIso_DR0p2To0p3  = 0;
300 >  Double_t tmpNeutralHadronIso_DR0p3To0p4  = 0;
301 >  Double_t tmpNeutralHadronIso_DR0p4To0p5  = 0;
302 >  Double_t tmpNeutralHadronIso_DR0p5To0p7  = 0;
303 >
304 >        
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 >  Double_t fChargedIso_DR0p5To0p7;
315 >
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 >  Double_t fGammaIso_DR0p5To0p7;
322 >
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 >  Double_t fNeutralHadronIso_DR0p5To0p7;
329 >
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 >    if (pf->HasTrackerTrk() && (pf->TrackerTrk() == mu->TrackerTrk()) ) continue;
343 >
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 >        // 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 >        // 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 >        if (abs(pf->PFType()) == PFCandidate::eGamma && fabs(tmpele->SCluster()->Eta()) > 1.479
370 >            && 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 >        // 4l muon
380 >        if( pf->HasTrackerTrk() ) {
381 >          if( pf->TrackerTrk() == tmpmu->TrackerTrk() )
382 >            IsLeptonFootprint = kTRUE;
383 >        }
384 >        // 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 >    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 >
401 >      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 >
414 >      // Footprint Veto
415 >      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 >    }
422 >
423 >    //
424 >    // Gamma Iso Rings
425 >    //
426 >    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 >    }
434 >
435 >    //
436 >    // Other Neutral Iso Rings
437 >    //
438 >    else {
439 >      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 >    }
446 >
447 >    }
448 >
449 >  }
450 >
451 >  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 >  fChargedIso_DR0p3To0p4   = min((tmpChargedIso_DR0p3To0p4)/mu->Pt(), 2.5);
455 >  fChargedIso_DR0p4To0p5   = min((tmpChargedIso_DR0p4To0p5)/mu->Pt(), 2.5);
456 >
457 >
458 >  double rho = 0;
459 >  if (!(isnan(fPUEnergyDensity->At(0)->Rho()) || isinf(fPUEnergyDensity->At(0)->Rho())))
460 >    rho = fPUEnergyDensity->At(0)->Rho();
461 >  
462 >
463 >  fGammaIso_DR0p0To0p1 = max(min((tmpGammaIso_DR0p0To0p1
464 >                                  -rho*muT.MuonEffectiveArea(muT.kMuGammaIsoDR0p0To0p1,mu->Eta(),EffectiveAreaVersion))/mu->Pt()
465 >                                 ,2.5)
466 >                             ,0.0);
467 >  fGammaIso_DR0p1To0p2 = max(min((tmpGammaIso_DR0p1To0p2
468 >                                  -rho*muT.MuonEffectiveArea(muT.kMuGammaIsoDR0p1To0p2,mu->Eta(),EffectiveAreaVersion))/mu->Pt()
469 >                                 ,2.5)
470 >                             ,0.0);
471 >  fGammaIso_DR0p2To0p3 = max(min((tmpGammaIso_DR0p2To0p3
472 >                                  -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 >
486 >  fNeutralHadronIso_DR0p0To0p1 = max(min((tmpNeutralHadronIso_DR0p0To0p1
487 >                                          -rho*muT.MuonEffectiveArea(muT.kMuNeutralHadronIsoDR0p0To0p1,
488 >                                                                 mu->Eta(),EffectiveAreaVersion))/mu->Pt()
489 >                                         , 2.5)
490 >                                     , 0.0);
491 >  fNeutralHadronIso_DR0p1To0p2 = max(min((tmpNeutralHadronIso_DR0p1To0p2
492 >                                            -rho*muT.MuonEffectiveArea(muT.kMuNeutralHadronIsoDR0p1To0p2,
493 >                                                                   mu->Eta(),EffectiveAreaVersion))/mu->Pt()
494 >                                           , 2.5)
495 >                                       , 0.0);
496 >  fNeutralHadronIso_DR0p2To0p3 = max(min((tmpNeutralHadronIso_DR0p2To0p3
497 >                                          -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 >
513 >  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 >  bool pass;
534 >
535 >  pass = false;
536 >  if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
537 >      && fabs(mu->Eta()) <= 1.5 && mu->Pt() <= 10 && mvaval >= MUON_ISOMVA_LOOSE_FORPFID_CUT_BIN0)   pass = true;
538 >  else if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
539 >           && fabs(mu->Eta()) <= 1.5 && mu->Pt() > 10 && mvaval >= MUON_ISOMVA_LOOSE_FORPFID_CUT_BIN1)  pass = true;
540 >  else if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
541 >           && fabs(mu->Eta()) > 1.5 && mu->Pt() <= 10 && mvaval >= MUON_ISOMVA_LOOSE_FORPFID_CUT_BIN2)  pass = true;
542 >  else if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
543 >           && 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 >
548 >  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 >
561 >  //  pass &= (fChargedIso_DR0p0To0p1 + fChargedIso_DR0p1To0p2 + fChargedIso_DR0p2To0p3 < 0.7);
562 >
563 >  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 >  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 >  muIsoMVA->Initialize( "MuonIsoMVA",
580 >                        mithep::MuonIDMVA::kIsoRingsV0,
581 >                        kTRUE, weightFiles);
582 > }
583 >
584 >
585 > //--------------------------------------------------------------------------------------------------
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 >
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 >  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 >  bool failiso=false;
822 >
823 >  //
824 >  // tmp iso rings
825 >  //
826 >  Double_t tmpChargedIso_DR0p0To0p1  = 0;
827 >  Double_t tmpChargedIso_DR0p1To0p2  = 0;
828 >  Double_t tmpChargedIso_DR0p2To0p3  = 0;
829 >  Double_t tmpChargedIso_DR0p3To0p4  = 0;
830 >  Double_t tmpChargedIso_DR0p4To0p5  = 0;
831 >  Double_t tmpChargedIso_DR0p5To0p7  = 0;
832 >
833 >  Double_t tmpGammaIso_DR0p0To0p1  = 0;
834 >  Double_t tmpGammaIso_DR0p1To0p2  = 0;
835 >  Double_t tmpGammaIso_DR0p2To0p3  = 0;
836 >  Double_t tmpGammaIso_DR0p3To0p4  = 0;
837 >  Double_t tmpGammaIso_DR0p4To0p5  = 0;
838 >  Double_t tmpGammaIso_DR0p5To0p7  = 0;
839 >
840 >  Double_t tmpNeutralHadronIso_DR0p0To0p1  = 0;
841 >  Double_t tmpNeutralHadronIso_DR0p1To0p2  = 0;
842 >  Double_t tmpNeutralHadronIso_DR0p2To0p3  = 0;
843 >  Double_t tmpNeutralHadronIso_DR0p3To0p4  = 0;
844 >  Double_t tmpNeutralHadronIso_DR0p4To0p5  = 0;
845 >  Double_t tmpNeutralHadronIso_DR0p5To0p7  = 0;
846 >
847 >        
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 >    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 >
887 >    if ( (pf->HasTrackerTrk() && (pf->TrackerTrk() == ele->TrackerTrk())) ||
888 >         (pf->HasGsfTrk() && (pf->GsfTrk() == ele->GsfTrk()))) continue;
889 >    
890 >
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 >        // 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 >        // PF charged
916 >        if (pf->Charge() != 0 && fabs(tmpele->SCluster()->Eta()) > 1.479
917 >            && 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 >          IsLeptonFootprint = kTRUE;
920 >        }
921 >        // PF gamma
922 >        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 >          IsLeptonFootprint = kTRUE;
926 >        }
927 >      } // 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 >        // 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 >        // PF charged
942 >        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 >          IsLeptonFootprint = kTRUE;
945 >        }
946 >      } // loop over muons
947 >
948 >
949 >    if (IsLeptonFootprint)
950 >      continue;
951 >
952 >    //
953 >    // Charged Iso Rings
954 >    //
955 >    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 >
962 >      // Veto any PFmuon, or PFEle
963 >      if (abs(pf->PFType()) == PFCandidate::eElectron || abs(pf->PFType()) == PFCandidate::eMuon) continue;
964 >
965 >      // Footprint Veto
966 >      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 >
979 >    }
980 >
981 >    //
982 >    // Gamma Iso Rings
983 >    //
984 >    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 >
1000 >    }
1001 >
1002 >    //
1003 >    // Other Neutral Iso Rings
1004 >    //
1005 >    else {
1006 >      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 >    }
1015 >
1016 >    }
1017 >
1018 >  }
1019 >
1020 >  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 >  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 >
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 >                                  -rho*eleT.ElectronEffectiveArea(eleT.kEleGammaIsoDR0p0To0p1,
1054 >                                                              ele->SCluster()->Eta(),
1055 >                                                              EffectiveAreaVersion))/ele->Pt()
1056 >                                 ,2.5)
1057 >                             ,0.0);
1058 >  fGammaIso_DR0p1To0p2 = max(min((tmpGammaIso_DR0p1To0p2
1059 >                                  -rho*eleT.ElectronEffectiveArea(eleT.kEleGammaIsoDR0p1To0p2,
1060 >                                                              ele->SCluster()->Eta(),
1061 >                                                              EffectiveAreaVersion))/ele->Pt()
1062 >                                 ,2.5)
1063 >                             ,0.0);
1064 >  fGammaIso_DR0p2To0p3 = max(min((tmpGammaIso_DR0p2To0p3
1065 >                                  -rho*eleT.ElectronEffectiveArea(eleT.kEleGammaIsoDR0p2To0p3,
1066 >                                                              ele->SCluster()->Eta()
1067 >                                                              ,EffectiveAreaVersion))/ele->Pt()
1068 >                                 ,2.5)
1069 >                             ,0.0);
1070 >  fGammaIso_DR0p3To0p4 = max(min((tmpGammaIso_DR0p3To0p4
1071 >                                  -rho*eleT.ElectronEffectiveArea(eleT.kEleGammaIsoDR0p3To0p4,
1072 >                                                              ele->SCluster()->Eta(),
1073 >                                                              EffectiveAreaVersion))/ele->Pt()
1074 >                                 ,2.5)
1075 >                             ,0.0);
1076 >  fGammaIso_DR0p4To0p5 = max(min((tmpGammaIso_DR0p4To0p5
1077 >                                  -rho*eleT.ElectronEffectiveArea(eleT.kEleGammaIsoDR0p4To0p5,
1078 >                                                              ele->SCluster()->Eta(),
1079 >                                                              EffectiveAreaVersion))/ele->Pt()
1080 >                                 ,2.5)
1081 >                             ,0.0);
1082 >
1083 >
1084 >  fNeutralHadronIso_DR0p0To0p1 = max(min((tmpNeutralHadronIso_DR0p0To0p1
1085 >                                          -rho*eleT.ElectronEffectiveArea(eleT.kEleNeutralHadronIsoDR0p0To0p1,
1086 >                                                                 ele->SCluster()->Eta(),EffectiveAreaVersion))/ele->Pt()
1087 >                                         , 2.5)
1088 >                                     , 0.0);
1089 >  fNeutralHadronIso_DR0p1To0p2 = max(min((tmpNeutralHadronIso_DR0p1To0p2
1090 >                                            -rho*eleT.ElectronEffectiveArea(eleT.kEleNeutralHadronIsoDR0p1To0p2,
1091 >                                                                   ele->SCluster()->Eta(),EffectiveAreaVersion))/ele->Pt()
1092 >                                           , 2.5)
1093 >                                       , 0.0);
1094 >  fNeutralHadronIso_DR0p2To0p3 = max(min((tmpNeutralHadronIso_DR0p2To0p3
1095 >                                          -rho*eleT.ElectronEffectiveArea(eleT.kEleNeutralHadronIsoDR0p2To0p3,
1096 >                                                                 ele->SCluster()->Eta(),EffectiveAreaVersion))/ele->Pt()
1097 >                                         , 2.5)
1098 >                                     , 0.0);
1099 >  fNeutralHadronIso_DR0p3To0p4 = max(min((tmpNeutralHadronIso_DR0p3To0p4
1100 >                                          -rho*eleT.ElectronEffectiveArea(eleT.kEleNeutralHadronIsoDR0p3To0p4,
1101 >                                                                 ele->SCluster()->Eta(), EffectiveAreaVersion))/ele->Pt()
1102 >                                         , 2.5)
1103 >                                     , 0.0);
1104 >  fNeutralHadronIso_DR0p4To0p5 = max(min((tmpNeutralHadronIso_DR0p4To0p5
1105 >                                          -rho*eleT.ElectronEffectiveArea(eleT.kEleNeutralHadronIsoDR0p4To0p5,
1106 >                                                                 ele->SCluster()->Eta(), EffectiveAreaVersion))/ele->Pt()
1107 >                                         , 2.5)
1108 >                                     , 0.0);
1109 >
1110 >  double mvaval = eleIsoMVA->MVAValue_IsoRings( ele->Pt(),
1111 >                                                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 >
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 >  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 >
1165 >  if(ctrl.debug) cout << "returning status : " << hex << status.getStatus() << dec << endl;
1166 >  return status;
1167 >  
1168 > }
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 >
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;

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines