ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitPhysics/Mods/src/ElectronIDMod.cc
Revision: 1.48
Committed: Sun Nov 8 13:52:29 2009 UTC (15 years, 5 months ago) by sixie
Content type: text/plain
Branch: MAIN
CVS Tags: Mit_012d, Mit_012c, Mit_012b
Changes since 1.47: +3 -1 lines
Log Message:
fix bug with conversion filter. was failing incorrectly if we did not use conversion filter.

File Contents

# Content
1 // $Id: ElectronIDMod.cc,v 1.47 2009/11/02 13:52:55 ceballos Exp $
2
3 #include "MitPhysics/Mods/interface/ElectronIDMod.h"
4 #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 #include "MitPhysics/Init/interface/ModNames.h"
9
10 using namespace mithep;
11
12 ClassImp(mithep::ElectronIDMod)
13
14 //--------------------------------------------------------------------------------------------------
15 ElectronIDMod::ElectronIDMod(const char *name, const char *title) :
16 BaseMod(name,title),
17 fElectronBranchName(Names::gkElectronBrn),
18 fConversionBranchName(Names::gkMvfConversionBrn),
19 fGoodElectronsName(ModNames::gkGoodElectronsName),
20 fVertexName(string("PrimaryVertexesBeamSpot").c_str()),
21 fElectronIDType("CustomTight"),
22 fElectronIsoType("TrackJuraSliding"),
23 fElectronPtMin(10),
24 fIDLikelihoodCut(0.9),
25 fTrackIsolationCut(5.0),
26 fCaloIsolationCut(5.0),
27 fEcalJuraIsoCut(5.0),
28 fHcalIsolationCut(5.0),
29 fApplyConvFilter(kTRUE),
30 fWrongHitsRequirement(kTRUE),
31 fApplyD0Cut(kTRUE),
32 fChargeFilter(kTRUE),
33 fD0Cut(0.025),
34 fReverseIsoCut(kFALSE),
35 fReverseD0Cut(kFALSE),
36 fElIdType(kIdUndef),
37 fElIsoType(kIsoUndef),
38 fElectrons(0),
39 fConversions(0),
40 fVertices(0)
41 {
42 // Constructor.
43 }
44
45 //--------------------------------------------------------------------------------------------------
46 Bool_t ElectronIDMod::PassCustomID(const Electron *ele) const
47 {
48 // Based on RecoEgamma/ElectronIdentification/src/CutBasedElectronID.cc.
49
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 if (deltaPhiIn>fCuts[2][cat+4*eb])
82 return kFALSE;
83
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 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 return passChargeFilter;
262 }
263
264 //--------------------------------------------------------------------------------------------------
265 void ElectronIDMod::Process()
266 {
267 // Process entries of the tree.
268
269 LoadEventObject(fElectronBranchName, fElectrons);
270
271 ElectronOArr *GoodElectrons = new ElectronOArr;
272 GoodElectrons->SetName(fGoodElectronsName);
273
274 for (UInt_t i=0; i<fElectrons->GetEntries(); ++i) {
275 const Electron *e = fElectrons->At(i);
276
277 if (e->Pt() <= fElectronPtMin)
278 continue;
279
280 //apply id cut
281 Bool_t idcut = PassIDCut(e, fElIdType);
282 if (!idcut)
283 continue;
284
285 //apply Isolation Cut
286 Bool_t isocut = PassIsolationCut(e, fElIsoType);
287 if (!isocut)
288 continue;
289
290 // apply conversion filter
291 Bool_t passConvVeto = kFALSE;
292 if (fApplyConvFilter) {
293 LoadEventObject(fConversionBranchName, fConversions);
294 passConvVeto = PassConversionFilter(e, fConversions);
295 } else {
296 passConvVeto = kTRUE;
297 }
298 if (passConvVeto == kFALSE) continue;
299
300 // apply d0 cut
301 if (fApplyD0Cut) {
302 LoadEventObject(fVertexName, fVertices);
303 Bool_t passD0cut = PassD0Cut(e, fVertices);
304 if (!passD0cut)
305 continue;
306 }
307
308 //apply charge filter
309 if(fChargeFilter == kTRUE) {
310 Bool_t passChargeFilter = PassChargeFilter(e);
311 if (!passChargeFilter) continue;
312 }
313
314 // add good electron
315 GoodElectrons->Add(e);
316 }
317
318 // sort according to pt
319 GoodElectrons->Sort();
320
321 // add to event for other modules to use
322 AddObjThisEvt(GoodElectrons);
323 }
324
325 //--------------------------------------------------------------------------------------------------
326 void ElectronIDMod::SlaveBegin()
327 {
328 // Run startup code on the computer (slave) doing the actual analysis. Here,
329 // we just request the electron collection branch.
330
331 ReqEventObject(fElectronBranchName, fElectrons, kTRUE);
332
333 if (fApplyConvFilter)
334 ReqEventObject(fConversionBranchName, fConversions, kTRUE);
335
336 if (fApplyD0Cut)
337 ReqEventObject(fVertexName, fVertices, kTRUE);
338
339 Setup();
340 }
341
342 //--------------------------------------------------------------------------------------------------
343 void ElectronIDMod::Setup()
344 {
345 // Set all options properly before execution.
346
347 if (fElectronIDType.CompareTo("Tight") == 0)
348 fElIdType = kTight;
349 else if (fElectronIDType.CompareTo("Loose") == 0)
350 fElIdType = kLoose;
351 else if (fElectronIDType.CompareTo("Likelihood") == 0)
352 fElIdType = kLikelihood;
353 else if (fElectronIDType.CompareTo("NoId") == 0)
354 fElIdType = kNoId;
355 else if (fElectronIDType.CompareTo("ZeeId") == 0)
356 fElIdType = kZeeId;
357 else if (fElectronIDType.CompareTo("CustomLoose") == 0) {
358 fElIdType = kCustomIdLoose;
359 } else if (fElectronIDType.CompareTo("CustomTight") == 0) {
360 fElIdType = kCustomIdTight;
361 }
362 else {
363 SendError(kAbortAnalysis, "SlaveBegin",
364 "The specified electron identification %s is not defined.",
365 fElectronIDType.Data());
366 return;
367 }
368
369 SetCustomIDCuts(fElIdType);
370
371 if (fElectronIsoType.CompareTo("TrackCalo") == 0 )
372 fElIsoType = kTrackCalo;
373 else if (fElectronIsoType.CompareTo("TrackJura") == 0)
374 fElIsoType = kTrackJura;
375 else if(fElectronIsoType.CompareTo("TrackJuraSliding") == 0)
376 fElIsoType = kTrackJuraSliding;
377 else if (fElectronIsoType.CompareTo("NoIso") == 0 )
378 fElIsoType = kNoIso;
379 else if (fElectronIsoType.CompareTo("ZeeIso") == 0 )
380 fElIsoType = kZeeIso;
381 else if (fElectronIsoType.CompareTo("Custom") == 0 ) {
382 fElIsoType = kCustomIso;
383 SendError(kWarning, "SlaveBegin",
384 "Custom electron isolation is not yet implemented.");
385 } else {
386 SendError(kAbortAnalysis, "SlaveBegin",
387 "The specified electron isolation %s is not defined.",
388 fElectronIsoType.Data());
389 return;
390 }
391
392 }
393
394 //--------------------------------------------------------------------------------------------------
395 void ElectronIDMod::SetCustomIDCuts(EElIdType idt)
396 {
397 // Set cut values based on RecoEgamma/ElectronIdentification/python/electronIdCutBasedExt_cfi.py.
398 // The following changes are in sigmaetaeta for endcups and deltaetain.
399
400 Double_t tightcuts[6][8]={
401 {0.086, 0.1, 0.052, 0.0, 0.050, 0.059, 0.061, 0.0}, //hovere
402 {0.011, 0.011, 0.011, 0.0, 0.033, 0.029, 0.030, 0.0}, //sigmaetaeta
403 {0.038, 0.024, 0.045, 0.0, 0.034, 0.017, 0.026, 0.0}, //deltaphiin
404 {0.0081, 0.0029, 0.0051, 0.0, 0.0070, 0.0062, 0.0088, 0.0}, //deltaetain
405 {0.0, 0.9, 0.0, 0.0, 0.0, 0.78, 0.0, 0.0}, //eoverp
406 {0.8,0.2,0.9,0,0,0,0,0}}; //extra cuts fbrem and E_Over_P
407
408 Double_t loosecuts[6][8]={
409 {0.076, 0.033, 0.07, 0.0, 0.083,0.148, 0.033, 0.0}, //hovere
410 {0.0101, 0.0095, 0.0097, 0.0, 0.03, 0.03, 0.03, 0.0}, //sigmaetaeta
411 {0.053, 0.0189, 0.059, 0.099, 0.0278,0.0157, 0.042, 0.080}, //deltaphiin
412 {0.0078, 0.00259, 0.0062, 0.0, 0.0078,0.0061, 0.0061, 0.0}, //deltaetain
413 {0.3, 0.92, 0.211, 0.0, 0.42, 0.88, 0.68, 0.0}, //eoverp
414 {0.8,0.2,0,0,0,0,0,0}}; //extra cuts fbrem and E_Over_P
415
416
417 switch (idt) {
418 case kCustomIdTight:
419 memcpy(fCuts,tightcuts,sizeof(fCuts));
420 break;
421 case kCustomIdLoose:
422 memcpy(fCuts,loosecuts,sizeof(fCuts));
423 break;
424 default:
425 memset(fCuts,0,sizeof(fCuts));
426 break;
427 }
428 }