ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitPhysics/Mods/src/ElectronIDMod.cc
Revision: 1.42
Committed: Mon Oct 26 14:33:20 2009 UTC (15 years, 6 months ago) by sixie
Content type: text/plain
Branch: MAIN
Changes since 1.41: +206 -139 lines
Log Message:
Restructure electron id module to be more modular. add different accessors for determining whether it passes specific parts of the electron ID, eg. conversion veto, d0cut, ID, isolation.

File Contents

# User Rev Content
1 sixie 1.42 // $Id: ElectronIDMod.cc,v 1.41 2009/10/21 13:45:44 ceballos Exp $
2 loizides 1.1
3     #include "MitPhysics/Mods/interface/ElectronIDMod.h"
4 loizides 1.26 #include "MitAna/DataTree/interface/StableData.h"
5     #include "MitAna/DataTree/interface/ElectronCol.h"
6     #include "MitAna/DataTree/interface/VertexCol.h"
7     #include "MitAna/DataTree/interface/DecayParticleCol.h"
8 loizides 1.5 #include "MitPhysics/Init/interface/ModNames.h"
9 loizides 1.1
10     using namespace mithep;
11    
12     ClassImp(mithep::ElectronIDMod)
13    
14     //--------------------------------------------------------------------------------------------------
15 loizides 1.5 ElectronIDMod::ElectronIDMod(const char *name, const char *title) :
16 loizides 1.1 BaseMod(name,title),
17 loizides 1.5 fElectronBranchName(Names::gkElectronBrn),
18 ceballos 1.12 fConversionBranchName(Names::gkMvfConversionBrn),
19 loizides 1.5 fGoodElectronsName(ModNames::gkGoodElectronsName),
20 ceballos 1.12 fVertexName(string("PrimaryVertexesBeamSpot").c_str()),
21 ceballos 1.31 fElectronIDType("CustomTight"),
22 ceballos 1.8 fElectronIsoType("TrackJuraSliding"),
23 loizides 1.1 fElectronPtMin(10),
24     fIDLikelihoodCut(0.9),
25     fTrackIsolationCut(5.0),
26     fCaloIsolationCut(5.0),
27 loizides 1.5 fEcalJuraIsoCut(5.0),
28     fHcalIsolationCut(5.0),
29 loizides 1.14 fApplyConvFilter(kTRUE),
30 bendavid 1.38 fWrongHitsRequirement(kTRUE),
31 ceballos 1.18 fApplyD0Cut(kTRUE),
32 loizides 1.14 fD0Cut(0.025),
33 ceballos 1.40 fChargeFilter(kTRUE),
34 loizides 1.25 fReverseIsoCut(kFALSE),
35     fReverseD0Cut(kFALSE),
36 loizides 1.5 fElIdType(kIdUndef),
37 ceballos 1.12 fElIsoType(kIsoUndef),
38 loizides 1.14 fElectrons(0),
39     fConversions(0),
40 loizides 1.25 fVertices(0)
41 loizides 1.1 {
42     // Constructor.
43     }
44    
45     //--------------------------------------------------------------------------------------------------
46 loizides 1.30 Bool_t ElectronIDMod::PassCustomID(const Electron *ele) const
47     {
48 loizides 1.32 // Based on RecoEgamma/ElectronIdentification/src/CutBasedElectronID.cc.
49 loizides 1.30
50     Double_t eOverP = ele->ESuperClusterOverP();
51     Double_t fBrem = ele->FBrem();
52    
53     if ((eOverP < fCuts[5][0]) && (fBrem < fCuts[5][1]))
54     return kFALSE;
55    
56     if (eOverP < fCuts[5][2]*(1-fBrem))
57     return kFALSE;
58    
59     Int_t cat = 2;
60     if ((ele->IsEB() && fBrem<0.06) || (ele->IsEE() && fBrem<0.1))
61     cat=1;
62     else if (eOverP < 1.2 && eOverP > 0.8)
63     cat=0;
64    
65     Double_t eSeedOverPin = ele->ESeedClusterOverPIn();
66     Double_t hOverE = ele->HadronicOverEm();
67     Double_t sigmaee = ele->CoviEtaiEta();
68     Double_t deltaPhiIn = TMath::Abs(ele->DeltaPhiSuperClusterTrackAtVtx());
69     Double_t deltaEtaIn = TMath::Abs(ele->DeltaEtaSuperClusterTrackAtVtx());
70    
71     Int_t eb = 1;
72     if (ele->IsEB())
73     eb = 0;
74    
75     if (hOverE>fCuts[0][cat+4*eb])
76     return kFALSE;
77    
78     if (sigmaee>fCuts[1][cat+4*eb])
79     return kFALSE;
80    
81 ceballos 1.40 if (deltaPhiIn>fCuts[2][cat+4*eb])
82     return kFALSE;
83 loizides 1.30
84     if(deltaEtaIn>fCuts[3][cat+4*eb])
85     return kFALSE;
86    
87     if(eSeedOverPin<fCuts[4][cat+4*eb])
88     return kFALSE;
89    
90     return kTRUE;
91     }
92    
93     //--------------------------------------------------------------------------------------------------
94 sixie 1.42 Bool_t ElectronIDMod::PassIDCut(const Electron *ele, EElIdType idType) const
95     {
96    
97     Bool_t idcut = kFALSE;
98     switch (idType) {
99     case kTight:
100     idcut = ele->PassTightID();
101     break;
102     case kLoose:
103     idcut = ele->PassLooseID();
104     break;
105     case kLikelihood:
106     idcut = (ele->IDLikelihood() > fIDLikelihoodCut);
107     break;
108     case kNoId:
109     idcut = kTRUE;
110     break;
111     case kCustomIdLoose:
112     idcut = ElectronIDMod::PassCustomID(ele);
113     break;
114     case kCustomIdTight:
115     idcut = ElectronIDMod::PassCustomID(ele);
116     break;
117     case kZeeId:
118     if (ele->IsEB()) {
119     idcut = (ele->CoviEtaiEta() < 0.01 && ele->DeltaEtaSuperClusterTrackAtVtx() < 0.0071);
120     } else {
121     idcut = (ele->CoviEtaiEta() < 0.028 && ele->DeltaEtaSuperClusterTrackAtVtx() < 0.0066);
122     }
123     break;
124     default:
125     break;
126     }
127    
128     return idcut;
129     }
130    
131    
132     //--------------------------------------------------------------------------------------------------
133     Bool_t ElectronIDMod::PassIsolationCut(const Electron *ele, EElIsoType isoType) const
134     {
135    
136     Bool_t isocut = kFALSE;
137     switch (isoType) {
138     case kTrackCalo:
139     isocut = (ele->TrackIsolationDr03() < fTrackIsolationCut) &&
140     (ele->CaloIsolation() < fCaloIsolationCut);
141     break;
142     case kTrackJura:
143     isocut = (ele->TrackIsolationDr03() < fTrackIsolationCut) &&
144     (ele->EcalRecHitIsoDr04() < fEcalJuraIsoCut) &&
145     (ele->HcalIsolation() < fHcalIsolationCut);
146     break;
147     case kTrackJuraSliding:
148     {
149     Double_t totalIso = ele->TrackIsolationDr03() + ele->EcalRecHitIsoDr04() - 1.5;
150     if (totalIso < (ele->Pt()-10.0)*4.5/20.0 ||
151     totalIso <= 0)
152     isocut = kTRUE;
153    
154     if (fReverseIsoCut == kTRUE &&
155     isocut == kFALSE && totalIso < 10)
156     isocut = kTRUE;
157     else if(fReverseIsoCut == kTRUE)
158     isocut = kFALSE;
159     }
160     break;
161     case kNoIso:
162     isocut = kTRUE;
163     break;
164     case kZeeIso:
165     if (ele->IsEB()) {
166     isocut = (ele->TrackIsolationDr04() < 7.2 && ele->EcalRecHitIsoDr04() < 5.7 && ele->HcalTowerSumEtDr04() < 8.1);
167     } else {
168     isocut = (ele->TrackIsolationDr04() < 5.1 && ele->EcalRecHitIsoDr04() < 5.0 && ele->HcalTowerSumEtDr04() < 3.4);
169     }
170     break;
171     case kCustomIso:
172     default:
173     break;
174     }
175    
176     return isocut;
177     }
178    
179    
180     //--------------------------------------------------------------------------------------------------
181     Bool_t ElectronIDMod::PassConversionFilter(const Electron *ele, const DecayParticleCol *conversions) const
182     {
183     Bool_t isGoodConversion = kFALSE;
184    
185     for (UInt_t ifc=0; ifc<conversions->GetEntries(); ifc++) {
186    
187     Bool_t ConversionMatchFound = kFALSE;
188     for (UInt_t d=0; d<conversions->At(ifc)->NDaughters(); d++) {
189     const Track *trk = dynamic_cast<const ChargedParticle*>
190     (conversions->At(ifc)->Daughter(d))->Trk();
191     if (ele->GsfTrk() == trk) {
192     ConversionMatchFound = kTRUE;
193     break;
194     }
195     }
196    
197     // if match between the e-track and one of the conversion legs
198     if (ConversionMatchFound == kTRUE){
199     isGoodConversion = (conversions->At(ifc)->Prob() > 1e-6) &&
200     (conversions->At(ifc)->Lxy() > 0) &&
201     (conversions->At(ifc)->Lz() > 0) &&
202     (conversions->At(ifc)->Position().Rho() > 2.0);
203    
204     if (isGoodConversion == kTRUE) {
205     for (UInt_t d=0; d<conversions->At(ifc)->NDaughters(); d++) {
206     const Track *trk = dynamic_cast<const ChargedParticle*>
207     (conversions->At(ifc)->Daughter(d))->Trk();
208    
209     if (trk) {
210     // These requirements are not used for the GSF track
211     if (!(trk->NHits() >= 3 && trk->Prob() > 1e-6) && trk!=ele->GsfTrk())
212     isGoodConversion = kFALSE;
213    
214     const StableData *sd = dynamic_cast<const StableData*>
215     (conversions->At(ifc)->DaughterDat(d));
216     if (fWrongHitsRequirement && sd->NWrongHits() != 0)
217     isGoodConversion = kFALSE;
218    
219     } else {
220     isGoodConversion = kFALSE;
221     }
222     }
223     }
224     }
225    
226     if (isGoodConversion == kTRUE) break;
227    
228     } // loop over all conversions
229    
230     return isGoodConversion;
231     }
232    
233     //--------------------------------------------------------------------------------------------------
234     Bool_t ElectronIDMod::PassD0Cut(const Electron *ele, const VertexCol *vertices) const
235     {
236     Bool_t d0cut = kFALSE;
237     // d0 cut
238     Double_t d0_real = 99999;
239     for(UInt_t i0 = 0; i0 < vertices->GetEntries(); i0++) {
240     Double_t pD0 = ele->GsfTrk()->D0Corrected(*vertices->At(i0));
241     if(TMath::Abs(pD0) < TMath::Abs(d0_real)) d0_real = TMath::Abs(pD0);
242     }
243     if(d0_real < fD0Cut) d0cut = kTRUE;
244    
245     if (fReverseD0Cut == kTRUE &&
246     d0cut == kFALSE && d0_real < 0.05)
247     d0cut = kTRUE;
248     else if(fReverseD0Cut == kTRUE)
249     d0cut = kFALSE;
250    
251     return d0cut;
252     }
253    
254     //--------------------------------------------------------------------------------------------------
255     Bool_t ElectronIDMod::PassChargeFilter(const Electron *ele) const
256     {
257     Bool_t passChargeFilter = kTRUE;
258     if(ele->TrackerTrk() &&
259     ele->TrackerTrk()->Charge() != ele->Charge()) passChargeFilter = kFALSE;
260    
261    
262     return passChargeFilter;
263     }
264    
265    
266     //--------------------------------------------------------------------------------------------------
267 loizides 1.1 void ElectronIDMod::Process()
268     {
269     // Process entries of the tree.
270    
271 loizides 1.23 LoadEventObject(fElectronBranchName, fElectrons);
272 loizides 1.1
273 loizides 1.6 ElectronOArr *GoodElectrons = new ElectronOArr;
274     GoodElectrons->SetName(fGoodElectronsName);
275 loizides 1.1
276 ceballos 1.18 for (UInt_t i=0; i<fElectrons->GetEntries(); ++i) {
277 loizides 1.5 const Electron *e = fElectrons->At(i);
278 loizides 1.1
279 loizides 1.5 if (e->Pt() <= fElectronPtMin)
280     continue;
281 loizides 1.1
282 sixie 1.42 //apply id cut
283     Bool_t idcut = PassIDCut(e, fElIdType);
284 loizides 1.5 if (!idcut)
285     continue;
286    
287 sixie 1.42 //apply Isolation Cut
288     Bool_t isocut = PassIsolationCut(e, fElIsoType);
289     if (!isocut)
290 loizides 1.5 continue;
291    
292 loizides 1.14 // apply conversion filter
293 ceballos 1.12 Bool_t isGoodConversion = kFALSE;
294 loizides 1.14 if (fApplyConvFilter) {
295 sixie 1.42 LoadEventObject(fConversionBranchName, fConversions);
296     isGoodConversion = PassConversionFilter(e, fConversions);
297 ceballos 1.12 }
298 sixie 1.42 if (isGoodConversion) continue;
299    
300     // apply d0 cut
301 ceballos 1.15 if (fApplyD0Cut) {
302 sixie 1.42 LoadEventObject(fVertexName, fVertices);
303     Bool_t passD0cut = PassD0Cut(e, *&fVertices);
304     if (!passD0cut)
305 ceballos 1.24 continue;
306 ceballos 1.12 }
307    
308 sixie 1.42 //apply charge filter
309     if(fChargeFilter == kTRUE) {
310     Bool_t passChargeFilter = PassChargeFilter(e);
311     if (!passChargeFilter) continue;
312     }
313 loizides 1.5 // add good electron
314 ceballos 1.12 GoodElectrons->Add(e);
315 loizides 1.5 }
316 loizides 1.1
317 loizides 1.9 // sort according to pt
318     GoodElectrons->Sort();
319    
320 loizides 1.5 // add to event for other modules to use
321 loizides 1.6 AddObjThisEvt(GoodElectrons);
322 loizides 1.1 }
323    
324     //--------------------------------------------------------------------------------------------------
325     void ElectronIDMod::SlaveBegin()
326     {
327     // Run startup code on the computer (slave) doing the actual analysis. Here,
328 loizides 1.5 // we just request the electron collection branch.
329 loizides 1.1
330 loizides 1.23 ReqEventObject(fElectronBranchName, fElectrons, kTRUE);
331 loizides 1.17
332 ceballos 1.16 if (fApplyConvFilter)
333 loizides 1.23 ReqEventObject(fConversionBranchName, fConversions, kTRUE);
334 loizides 1.17
335 ceballos 1.15 if (fApplyD0Cut)
336 loizides 1.23 ReqEventObject(fVertexName, fVertices, kTRUE);
337 loizides 1.1
338 sixie 1.42 Setup();
339     }
340    
341     //--------------------------------------------------------------------------------------------------
342     void ElectronIDMod::Setup()
343     {
344     // Set all options properly before execution.
345    
346 loizides 1.5 if (fElectronIDType.CompareTo("Tight") == 0)
347     fElIdType = kTight;
348     else if (fElectronIDType.CompareTo("Loose") == 0)
349     fElIdType = kLoose;
350     else if (fElectronIDType.CompareTo("Likelihood") == 0)
351     fElIdType = kLikelihood;
352 loizides 1.10 else if (fElectronIDType.CompareTo("NoId") == 0)
353     fElIdType = kNoId;
354 sixie 1.42 else if (fElectronIDType.CompareTo("ZeeId") == 0)
355     fElIdType = kZeeId;
356 loizides 1.30 else if (fElectronIDType.CompareTo("CustomLoose") == 0) {
357 peveraer 1.29 fElIdType = kCustomIdLoose;
358 loizides 1.30 } else if (fElectronIDType.CompareTo("CustomTight") == 0) {
359 peveraer 1.29 fElIdType = kCustomIdTight;
360 loizides 1.30 }
361 peveraer 1.29 else {
362 loizides 1.5 SendError(kAbortAnalysis, "SlaveBegin",
363     "The specified electron identification %s is not defined.",
364     fElectronIDType.Data());
365     return;
366     }
367    
368 loizides 1.30 SetCustomIDCuts(fElIdType);
369    
370 loizides 1.5 if (fElectronIsoType.CompareTo("TrackCalo") == 0 )
371     fElIsoType = kTrackCalo;
372     else if (fElectronIsoType.CompareTo("TrackJura") == 0)
373     fElIsoType = kTrackJura;
374     else if(fElectronIsoType.CompareTo("TrackJuraSliding") == 0)
375     fElIsoType = kTrackJuraSliding;
376     else if (fElectronIsoType.CompareTo("NoIso") == 0 )
377     fElIsoType = kNoIso;
378 sixie 1.42 else if (fElectronIsoType.CompareTo("ZeeIso") == 0 )
379     fElIsoType = kZeeIso;
380 loizides 1.5 else if (fElectronIsoType.CompareTo("Custom") == 0 ) {
381     fElIsoType = kCustomIso;
382     SendError(kWarning, "SlaveBegin",
383     "Custom electron isolation is not yet implemented.");
384     } else {
385     SendError(kAbortAnalysis, "SlaveBegin",
386     "The specified electron isolation %s is not defined.",
387     fElectronIsoType.Data());
388     return;
389     }
390 sixie 1.42
391 loizides 1.1 }
392 loizides 1.30
393     //--------------------------------------------------------------------------------------------------
394     void ElectronIDMod::SetCustomIDCuts(EElIdType idt)
395     {
396 loizides 1.32 // Set cut values based on RecoEgamma/ElectronIdentification/python/electronIdCutBasedExt_cfi.py.
397     // The following changes are in sigmaetaeta for endcups and deltaetain.
398 loizides 1.30
399     Double_t tightcuts[6][8]={
400 ceballos 1.40 {0.086, 0.1, 0.052, 0.0, 0.050, 0.059, 0.061, 0.0}, //hovere
401     {0.01, 0.01, 0.01, 0.0, 0.033, 0.029, 0.028, 0.0}, //sigmaetaeta
402     {0.038, 0.024, 0.038, 0.0, 0.034, 0.011, 0.023, 0.0}, //deltaphiin
403     {0.0081, 0.0029, 0.0051, 0.0, 0.0056, 0.0062, 0.0088, 0.0}, //deltaetain
404     {0.0, 0.9, 0.0, 0.0, 0.0, 0.78, 0.0, 0.0}, //eoverp
405     {0.8,0.2,0.9,0,0,0,0,0}}; //extra cuts fbrem and E_Over_P
406 loizides 1.30
407     Double_t loosecuts[6][8]={
408 loizides 1.39 {0.076, 0.033, 0.07, 0.0, 0.083,0.148, 0.033, 0.0}, //hovere
409     {0.0101, 0.0095, 0.0097, 0.0, 0.03, 0.03, 0.03, 0.0}, //sigmaetaeta
410     {0.053, 0.0189, 0.059, 0.099, 0.0278,0.0157, 0.042, 0.080}, //deltaphiin
411     {0.0078, 0.00259, 0.0062, 0.0, 0.0078,0.0061, 0.0061, 0.0}, //deltaetain
412     {0.3, 0.92, 0.211, 0.0, 0.42, 0.88, 0.68, 0.0}, //eoverp
413     {0.8,0.2,0,0,0,0,0,0}}; //extra cuts fbrem and E_Over_P
414 loizides 1.30
415 sixie 1.42
416 loizides 1.30 switch (idt) {
417     case kCustomIdTight:
418     memcpy(fCuts,tightcuts,sizeof(fCuts));
419     break;
420     case kCustomIdLoose:
421     memcpy(fCuts,loosecuts,sizeof(fCuts));
422     break;
423     default:
424     memset(fCuts,0,sizeof(fCuts));
425     break;
426     }
427     }