ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitPhysics/Mods/src/ElectronIDMod.cc
Revision: 1.32
Committed: Wed Aug 12 07:50:31 2009 UTC (15 years, 8 months ago) by loizides
Content type: text/plain
Branch: MAIN
CVS Tags: Mit_011, Mit_010a
Changes since 1.31: +4 -3 lines
Log Message:
Improved comment

File Contents

# Content
1 // $Id: ElectronIDMod.cc,v 1.31 2009/08/11 15:16:17 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 fApplyD0Cut(kTRUE),
31 fD0Cut(0.025),
32 fReverseIsoCut(kFALSE),
33 fReverseD0Cut(kFALSE),
34 fElIdType(kIdUndef),
35 fElIsoType(kIsoUndef),
36 fElectrons(0),
37 fConversions(0),
38 fVertices(0)
39 {
40 // Constructor.
41 }
42
43 //--------------------------------------------------------------------------------------------------
44 Bool_t ElectronIDMod::PassCustomID(const Electron *ele) const
45 {
46 // Based on RecoEgamma/ElectronIdentification/src/CutBasedElectronID.cc.
47
48 Double_t eOverP = ele->ESuperClusterOverP();
49 Double_t fBrem = ele->FBrem();
50
51 if ((eOverP < fCuts[5][0]) && (fBrem < fCuts[5][1]))
52 return kFALSE;
53
54 if (eOverP < fCuts[5][2]*(1-fBrem))
55 return kFALSE;
56
57 Int_t cat = 2;
58 if ((ele->IsEB() && fBrem<0.06) || (ele->IsEE() && fBrem<0.1))
59 cat=1;
60 else if (eOverP < 1.2 && eOverP > 0.8)
61 cat=0;
62
63 Double_t eSeedOverPin = ele->ESeedClusterOverPIn();
64 Double_t hOverE = ele->HadronicOverEm();
65 Double_t sigmaee = ele->CoviEtaiEta();
66 Double_t deltaPhiIn = TMath::Abs(ele->DeltaPhiSuperClusterTrackAtVtx());
67 Double_t deltaEtaIn = TMath::Abs(ele->DeltaEtaSuperClusterTrackAtVtx());
68
69 Int_t eb = 1;
70 if (ele->IsEB())
71 eb = 0;
72
73 if (hOverE>fCuts[0][cat+4*eb])
74 return kFALSE;
75
76 if (sigmaee>fCuts[1][cat+4*eb])
77 return kFALSE;
78
79 if (eOverP<1.5) {
80 if (deltaPhiIn>fCuts[2][cat+4*eb])
81 return kFALSE;
82 } else {
83 if(deltaPhiIn>fCuts[2][3+4*eb])
84 return kFALSE;
85 }
86
87 if(deltaEtaIn>fCuts[3][cat+4*eb])
88 return kFALSE;
89
90 if(eSeedOverPin<fCuts[4][cat+4*eb])
91 return kFALSE;
92
93 return kTRUE;
94 }
95
96 //--------------------------------------------------------------------------------------------------
97 void ElectronIDMod::Process()
98 {
99 // Process entries of the tree.
100
101 LoadEventObject(fElectronBranchName, fElectrons);
102
103 ElectronOArr *GoodElectrons = new ElectronOArr;
104 GoodElectrons->SetName(fGoodElectronsName);
105
106 for (UInt_t i=0; i<fElectrons->GetEntries(); ++i) {
107 const Electron *e = fElectrons->At(i);
108
109 if (e->Pt() <= fElectronPtMin)
110 continue;
111
112 Bool_t idcut = kFALSE;
113 switch (fElIdType) {
114 case kTight:
115 idcut = e->PassTightID();
116 break;
117 case kLoose:
118 idcut = e->PassLooseID();
119 break;
120 case kLikelihood:
121 idcut = (e->IDLikelihood() > fIDLikelihoodCut);
122 break;
123 case kNoId:
124 idcut = kTRUE;
125 break;
126 case kCustomIdLoose:
127 idcut = ElectronIDMod::PassCustomID(e);
128 break;
129 case kCustomIdTight:
130 idcut = ElectronIDMod::PassCustomID(e);
131 break;
132 default:
133 break;
134 }
135
136 if (!idcut)
137 continue;
138
139 Bool_t isocut = kFALSE;
140 switch (fElIsoType) {
141 case kTrackCalo:
142 isocut = (e->TrackIsolationDr03() < fTrackIsolationCut) &&
143 (e->CaloIsolation() < fCaloIsolationCut);
144 break;
145 case kTrackJura:
146 isocut = (e->TrackIsolationDr03() < fTrackIsolationCut) &&
147 (e->EcalRecHitIsoDr04() < fEcalJuraIsoCut) &&
148 (e->HcalIsolation() < fHcalIsolationCut);
149 break;
150 case kTrackJuraSliding:
151 {
152 Double_t totalIso = e->TrackIsolationDr03() + e->EcalRecHitIsoDr04() - 1.5;
153 if ((totalIso < (e->Pt()-10.0)*5.0/15.0 && e->Pt() <= 25) ||
154 (totalIso < 5.0 && e->Pt() > 25) ||
155 totalIso <= 0)
156 isocut = kTRUE;
157
158 if (fReverseIsoCut == kTRUE &&
159 isocut == kFALSE && totalIso < 10)
160 isocut = kTRUE;
161 else if(fReverseIsoCut == kTRUE)
162 isocut = kFALSE;
163 }
164 break;
165 case kNoIso:
166 isocut = kTRUE;
167 break;
168 case kCustomIso:
169 default:
170 break;
171 }
172
173 if (isocut == kFALSE)
174 continue;
175
176 // apply conversion filter
177 Bool_t isGoodConversion = kFALSE;
178 if (fApplyConvFilter) {
179 LoadBranch(fConversionBranchName);
180 for (UInt_t ifc=0; ifc<fConversions->GetEntries(); ifc++) {
181
182 Bool_t ConversionMatchFound = kFALSE;
183 for (UInt_t d=0; d<fConversions->At(ifc)->NDaughters(); d++) {
184 const Track *trk = dynamic_cast<const ChargedParticle*>
185 (fConversions->At(ifc)->Daughter(d))->Trk();
186 if (e->GsfTrk() == trk) {
187 ConversionMatchFound = kTRUE;
188 break;
189 }
190 }
191
192 // if match between the e-track and one of the conversion legs
193 if (ConversionMatchFound == kTRUE){
194 isGoodConversion = (fConversions->At(ifc)->Prob() > 0.0005) &&
195 (fConversions->At(ifc)->Lxy() > 0) &&
196 (fConversions->At(ifc)->Lz() > 0);
197
198 if (isGoodConversion == kTRUE) {
199 for (UInt_t d=0; d<fConversions->At(ifc)->NDaughters(); d++) {
200 const Track *trk = dynamic_cast<const ChargedParticle*>
201 (fConversions->At(ifc)->Daughter(d))->Trk();
202
203 if (trk) {
204 // These requirements are not used for the GSF track
205 if (!(trk->NHits() > 8 && trk->Prob() > 0.005) && trk!=e->GsfTrk())
206 isGoodConversion = kFALSE;
207
208 const StableData *sd = dynamic_cast<const StableData*>
209 (fConversions->At(ifc)->DaughterDat(d));
210 if (sd->NWrongHits() != 0)
211 isGoodConversion = kFALSE;
212
213 } else {
214 isGoodConversion = kFALSE;
215 }
216 }
217 }
218 }
219
220 if (isGoodConversion == kTRUE) break;
221
222 } // loop over all conversions
223
224 }
225 if (isGoodConversion == kTRUE) continue;
226
227 if (fApplyD0Cut) {
228 Bool_t d0cut = kFALSE;
229 LoadBranch(fVertexName);
230 // d0 cut
231 Double_t d0_real = 99999;
232 for(UInt_t i0 = 0; i0 < fVertices->GetEntries(); i0++) {
233 Double_t pD0 = e->GsfTrk()->D0Corrected(*fVertices->At(i0));
234 if(TMath::Abs(pD0) < TMath::Abs(d0_real)) d0_real = TMath::Abs(pD0);
235 }
236 if(d0_real < fD0Cut) d0cut = kTRUE;
237
238 if (fReverseD0Cut == kTRUE &&
239 d0cut == kFALSE && d0_real < 0.05)
240 d0cut = kTRUE;
241 else if(fReverseD0Cut == kTRUE)
242 d0cut = kFALSE;
243
244 if (d0cut == kFALSE)
245 continue;
246 }
247
248 // add good electron
249 GoodElectrons->Add(e);
250 }
251
252 // sort according to pt
253 GoodElectrons->Sort();
254
255 // add to event for other modules to use
256 AddObjThisEvt(GoodElectrons);
257 }
258
259 //--------------------------------------------------------------------------------------------------
260 void ElectronIDMod::SlaveBegin()
261 {
262 // Run startup code on the computer (slave) doing the actual analysis. Here,
263 // we just request the electron collection branch.
264
265 ReqEventObject(fElectronBranchName, fElectrons, kTRUE);
266
267 if (fApplyConvFilter)
268 ReqEventObject(fConversionBranchName, fConversions, kTRUE);
269
270 if (fApplyD0Cut)
271 ReqEventObject(fVertexName, fVertices, kTRUE);
272
273 if (fElectronIDType.CompareTo("Tight") == 0)
274 fElIdType = kTight;
275 else if (fElectronIDType.CompareTo("Loose") == 0)
276 fElIdType = kLoose;
277 else if (fElectronIDType.CompareTo("Likelihood") == 0)
278 fElIdType = kLikelihood;
279 else if (fElectronIDType.CompareTo("NoId") == 0)
280 fElIdType = kNoId;
281 else if (fElectronIDType.CompareTo("CustomLoose") == 0) {
282 fElIdType = kCustomIdLoose;
283 } else if (fElectronIDType.CompareTo("CustomTight") == 0) {
284 fElIdType = kCustomIdTight;
285 }
286 else {
287 SendError(kAbortAnalysis, "SlaveBegin",
288 "The specified electron identification %s is not defined.",
289 fElectronIDType.Data());
290 return;
291 }
292
293 SetCustomIDCuts(fElIdType);
294
295 if (fElectronIsoType.CompareTo("TrackCalo") == 0 )
296 fElIsoType = kTrackCalo;
297 else if (fElectronIsoType.CompareTo("TrackJura") == 0)
298 fElIsoType = kTrackJura;
299 else if(fElectronIsoType.CompareTo("TrackJuraSliding") == 0)
300 fElIsoType = kTrackJuraSliding;
301 else if (fElectronIsoType.CompareTo("NoIso") == 0 )
302 fElIsoType = kNoIso;
303 else if (fElectronIsoType.CompareTo("Custom") == 0 ) {
304 fElIsoType = kCustomIso;
305 SendError(kWarning, "SlaveBegin",
306 "Custom electron isolation is not yet implemented.");
307 } else {
308 SendError(kAbortAnalysis, "SlaveBegin",
309 "The specified electron isolation %s is not defined.",
310 fElectronIsoType.Data());
311 return;
312 }
313 }
314
315 //--------------------------------------------------------------------------------------------------
316 void ElectronIDMod::SetCustomIDCuts(EElIdType idt)
317 {
318 // Set cut values based on RecoEgamma/ElectronIdentification/python/electronIdCutBasedExt_cfi.py.
319 // The following changes are in sigmaetaeta for endcups and deltaetain.
320
321 Double_t tightcuts[6][8]={
322 {0.05, 0.042, 0.045, 0.0, 0.055, 0.037, 0.05, 0.0}, //hovere
323 {0.0125, 0.011, 0.01, 0.0, 0.0295, 0.0292, 0.0283, 0.0}, //sigmaetaeta
324 {0.032, 0.016, 0.0525, 0.09, 0.025, 0.035, 0.065, 0.092}, //deltaphiin
325 {0.0055, 0.003, 0.0075, 0.0, 0.0070, 0.0055, 0.0085, 0.0}, //deltaetain
326 {0.24, 0.94, 0.11, 0.0, 0.32, 0.83, 0.0, 0.0}, //eoverp
327 {0.8,0.2,0.9,0,0,0,0,0}}; //extra cuts fbrem and E_Over_P
328
329 Double_t loosecuts[6][8]={
330 {0.076, 0.033, 0.07, 0.0, 0.083,0.148, 0.033, 0.0}, //hovere
331 {0.0101, 0.0095, 0.0097, 0.0, 0.03, 0.03, 0.03, 0.0}, //sigmaetaeta
332 {0.053, 0.0189, 0.059, 0.099, 0.0278,0.0157, 0.042, 0.080}, //deltaphiin
333 {0.0078, 0.00259, 0.0062, 0.0, 0.0078,0.0061, 0.0061, 0.0}, //deltaetain
334 {0.3, 0.92, 0.211, 0.0, 0.42, 0.88, 0.68, 0.0}, //eoverp
335 {0.8,0.2,0,0,0,0,0,0}}; //extra cuts fbrem and E_Over_P
336
337 switch (idt) {
338 case kCustomIdTight:
339 memcpy(fCuts,tightcuts,sizeof(fCuts));
340 break;
341 case kCustomIdLoose:
342 memcpy(fCuts,loosecuts,sizeof(fCuts));
343 break;
344 default:
345 memset(fCuts,0,sizeof(fCuts));
346 break;
347 }
348 }
349