ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitPhysics/Utils/src/ElectronTools.cc
Revision: 1.2
Committed: Sat Apr 10 19:51:58 2010 UTC (15 years ago) by sixie
Content type: text/plain
Branch: MAIN
CVS Tags: Mit_013d
Changes since 1.1: +23 -1 lines
Log Message:
Add possibility for D0 to beamspot

File Contents

# User Rev Content
1 sixie 1.2 // $Id: ElectronTools.cc,v 1.1 2010/04/10 18:06:51 sixie Exp $
2 sixie 1.1
3     #include "MitPhysics/Utils/interface/ElectronTools.h"
4     #include "MitAna/DataTree/interface/StableData.h"
5     #include <TFile.h>
6    
7     ClassImp(mithep::ElectronTools)
8    
9     using namespace mithep;
10    
11     //--------------------------------------------------------------------------------------------------
12     ElectronTools::ElectronTools()
13     {
14     // Constructor.
15     }
16    
17     //--------------------------------------------------------------------------------------------------
18     Bool_t ElectronTools::PassCustomID(const Electron *ele, EElIdType idType) {
19    
20     Double_t fCuts[6][8]; //!custom id cuts
21    
22     Double_t tightcuts[6][8]={
23     {0.086, 0.1, 0.052, 0.0, 0.050, 0.059, 0.061, 0.0}, //hovere
24     {0.011, 0.011, 0.011, 0.0, 0.033, 0.029, 0.030, 0.0}, //sigmaetaeta
25     {0.038, 0.024, 0.045, 0.0, 0.034, 0.017, 0.026, 0.0}, //deltaphiin
26     {0.0081, 0.0029, 0.0051, 0.0, 0.0070, 0.0062, 0.0088, 0.0}, //deltaetain
27     {0.0, 0.9, 0.0, 0.0, 0.0, 0.78, 0.0, 0.0}, //eoverp
28     {0.8,0.2,0.9,0,0,0,0,0}}; //extra cuts fbrem and E_Over_P
29    
30     Double_t loosecuts[6][8]={
31     {0.076, 0.033, 0.07, 0.0, 0.083,0.148, 0.033, 0.0}, //hovere
32     {0.0101, 0.0095, 0.0097, 0.0, 0.03, 0.03, 0.03, 0.0}, //sigmaetaeta
33     {0.053, 0.0189, 0.059, 0.099, 0.0278,0.0157, 0.042, 0.080}, //deltaphiin
34     {0.0078, 0.00259, 0.0062, 0.0, 0.0078,0.0061, 0.0061, 0.0}, //deltaetain
35     {0.3, 0.92, 0.211, 0.0, 0.42, 0.88, 0.68, 0.0}, //eoverp
36     {0.8,0.2,0,0,0,0,0,0}}; //extra cuts fbrem and E_Over_P
37    
38     Double_t VBTFWorkingPoint95[6][8] = {
39     {0.05, 0.05, 0.05, 0.05, 0.04, 0.04, 0.04, 0.04 }, //hovere
40     {0.01, 0.01, 0.01, 0.01, 0.03, 0.03, 0.03, 0.03 }, //sigmaetaeta
41     {0.8, 0.8, 0.8, 0.8, 0.7, 0.7, 0.7, 0.7 }, //deltaphiin
42     {0.006, 0.006, 0.006, 0.006, 0.008, 0.008, 0.008, 0.008 }, //deltaetain
43     {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 }, //eoverp
44     {0.0, 0.0, 0, 0, 0, 0, 0, 0 } //extra cuts fbrem and E_Over_P
45     };
46    
47    
48     Double_t VBTFWorkingPoint90[6][8] = {
49     {0.05, 0.05, 0.05, 0.05, 0.025, 0.025, 0.025, 0.025 }, //hovere
50     {0.01, 0.01, 0.01, 0.01, 0.03, 0.03, 0.03, 0.03 }, //sigmaetaeta
51     {0.04, 0.04, 0.04, 0.04, 0.025, 0.025, 0.025, 0.025 }, //deltaphiin
52     {0.006, 0.006, 0.006, 0.006, 0.008, 0.008, 0.008, 0.008 }, //deltaetain
53     {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 }, //eoverp
54     {0.0, 0.0, 0, 0, 0, 0, 0, 0 } //extra cuts fbrem and E_Over_P
55     };
56    
57     Double_t VBTFWorkingPoint80[6][8] = {
58     {0.05, 0.05, 0.05, 0.05, 0.025, 0.025, 0.025, 0.025}, //hovere
59     {0.01, 0.01, 0.01, 0.01, 0.03, 0.03, 0.03, 0.03 }, //sigmaetaeta
60     {0.02, 0.02, 0.02, 0.02, 0.02, 0.02, 0.02, 0.02 }, //deltaphiin
61     {0.006, 0.006, 0.006, 0.006, 0.006, 0.006, 0.006, 0.006}, //deltaetain
62     {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 }, //eoverp
63     {0.0, 0.0, 0, 0, 0, 0, 0, 0 } //extra cuts fbrem and E_Over_P
64     };
65    
66     Double_t VBTFWorkingPoint70[6][8] = {
67     {0.02, 0.02, 0.02, 0.02, 0.025, 0.025, 0.025, 0.025}, //hovere
68     {0.01, 0.01, 0.01, 0.01, 0.03, 0.03, 0.03, 0.03 }, //sigmaetaeta
69     {0.02, 0.02, 0.02, 0.02, 0.02, 0.02, 0.02, 0.02 }, //deltaphiin
70     {0.006, 0.006, 0.006, 0.006, 0.003, 0.003, 0.003, 0.003}, //deltaetain
71     {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 }, //eoverp
72     {0.0, 0.0, 0, 0, 0, 0, 0, 0 } //extra cuts fbrem and E_Over_P
73     };
74    
75     switch (idType) {
76     case kCustomIdTight:
77     memcpy(fCuts,tightcuts,sizeof(fCuts));
78     break;
79     case kCustomIdLoose:
80     memcpy(fCuts,loosecuts,sizeof(fCuts));
81     break;
82     case kVBTFWorkingPoint95Id:
83     memcpy(fCuts,VBTFWorkingPoint95,sizeof(fCuts));
84     break;
85     case kVBTFWorkingPoint90Id:
86     memcpy(fCuts,VBTFWorkingPoint90,sizeof(fCuts));
87     break;
88     case kVBTFWorkingPoint80Id:
89     memcpy(fCuts,VBTFWorkingPoint80,sizeof(fCuts));
90     break;
91     case kVBTFWorkingPoint70Id:
92     memcpy(fCuts,VBTFWorkingPoint70,sizeof(fCuts));
93     break;
94     default:
95     memset(fCuts,0,sizeof(fCuts));
96     break;
97     }
98    
99    
100     // Based on RecoEgamma/ElectronIdentification/src/CutBasedElectronID.cc.
101     Double_t eOverP = ele->ESuperClusterOverP();
102     Double_t fBrem = ele->FBrem();
103    
104     if ((eOverP < fCuts[5][0]) && (fBrem < fCuts[5][1]))
105     return kFALSE;
106    
107     if (eOverP < fCuts[5][2]*(1-fBrem))
108     return kFALSE;
109    
110     Int_t cat = 2;
111     if ((ele->IsEB() && fBrem<0.06) || (ele->IsEE() && fBrem<0.1))
112     cat=1;
113     else if (eOverP < 1.2 && eOverP > 0.8)
114     cat=0;
115    
116     Double_t eSeedOverPin = ele->ESeedClusterOverPIn();
117     Double_t hOverE = ele->HadronicOverEm();
118     Double_t sigmaee = ele->CoviEtaiEta();
119     Double_t deltaPhiIn = TMath::Abs(ele->DeltaPhiSuperClusterTrackAtVtx());
120     Double_t deltaEtaIn = TMath::Abs(ele->DeltaEtaSuperClusterTrackAtVtx());
121    
122     Int_t eb = 1;
123     if (ele->IsEB())
124     eb = 0;
125    
126     if (hOverE>fCuts[0][cat+4*eb])
127     return kFALSE;
128    
129     if (sigmaee>fCuts[1][cat+4*eb])
130     return kFALSE;
131    
132     if (deltaPhiIn>fCuts[2][cat+4*eb])
133     return kFALSE;
134    
135     if(deltaEtaIn>fCuts[3][cat+4*eb])
136     return kFALSE;
137    
138     if(eSeedOverPin<fCuts[4][cat+4*eb])
139     return kFALSE;
140    
141     return kTRUE;
142     }
143    
144     //--------------------------------------------------------------------------------------------------
145     Bool_t ElectronTools::PassCustomIso(const Electron *ele, EElIsoType isoType)
146     {
147     Bool_t pass = kTRUE;
148     Double_t fIsoCuts[4][2]; //!custom isolation cuts
149    
150     Double_t VBTFWorkingPoint95[4][2] = {
151     {7.0 , 8.0 }, //TrkIso
152     {5.0 , 3.0 }, //ECALIso
153     {5.0 , 2.0 }, //HCALIso
154     {9999, 9999 } //Combined
155     };
156    
157     Double_t VBTFWorkingPoint90[4][2] = {
158     {6.0 , 6.0 }, //TrkIso
159     {5.0 , 2.5 }, //ECALIso
160     {5.0 , 1.5 }, //HCALIso
161     {9999, 9999 } //Combined
162     };
163    
164     Double_t VBTFWorkingPoint80[4][2] = {
165     {3.0 , 1.5 }, //TrkIso
166     {4.0 , 2.5 }, //ECALIso
167     {5.0 , 0.7 }, //HCALIso
168     {9999, 9999 } //Combined
169     };
170    
171     Double_t VBTFWorkingPoint70[4][2] = {
172     {2.5 , 0.8 }, //TrkIso
173     {3.0 , 2.5 }, //ECALIso
174     {5.0 , 0.25 }, //HCALIso
175     {9999, 9999 } //Combined
176     };
177    
178     switch (isoType) {
179     case kVBTFWorkingPoint95Iso:
180     memcpy(fIsoCuts,VBTFWorkingPoint95,sizeof(fIsoCuts));
181     break;
182     case kVBTFWorkingPoint90Iso:
183     memcpy(fIsoCuts,VBTFWorkingPoint90,sizeof(fIsoCuts));
184     break;
185     case kVBTFWorkingPoint80Iso:
186     memcpy(fIsoCuts,VBTFWorkingPoint80,sizeof(fIsoCuts));
187     break;
188     case kVBTFWorkingPoint70Iso:
189     memcpy(fIsoCuts,VBTFWorkingPoint70,sizeof(fIsoCuts));
190     break;
191     default:
192     memset(fIsoCuts,0,sizeof(fIsoCuts));
193     break;
194     }
195    
196     Double_t trkIso = ele->TrackIsolationDr03();
197     Double_t ecalIso = ele->EcalRecHitIsoDr04();
198     Double_t hcalIso = ele->HcalIsolation();
199     Double_t combinedIso = ele->TrackIsolationDr03() + ele->EcalRecHitIsoDr04() - 1.5;
200    
201     Int_t eb = 1;
202     if (ele->IsEB())
203     eb = 0;
204    
205     if (trkIso>fIsoCuts[0][eb])
206     pass = kFALSE;
207     if (ecalIso>fIsoCuts[1][eb])
208     pass = kFALSE;
209     if (hcalIso>fIsoCuts[2][eb])
210     pass = kFALSE;
211     if (combinedIso>fIsoCuts[3][eb])
212     pass = kFALSE;
213    
214    
215     return pass;
216     }
217    
218    
219     //--------------------------------------------------------------------------------------------------
220     Bool_t ElectronTools::PassConversionFilter(const Electron *ele,
221     const DecayParticleCol *conversions,
222     Bool_t WrongHitsRequirement)
223     {
224     Bool_t isGoodConversion = kFALSE;
225    
226     for (UInt_t ifc=0; ifc<conversions->GetEntries(); ifc++) {
227     Bool_t ConversionMatchFound = kFALSE;
228     for (UInt_t d=0; d<conversions->At(ifc)->NDaughters(); d++) {
229     const Track *trk = dynamic_cast<const ChargedParticle*>
230     (conversions->At(ifc)->Daughter(d))->Trk();
231     if (ele->GsfTrk() == trk) {
232     ConversionMatchFound = kTRUE;
233     break;
234     }
235     }
236    
237     // if match between the e-track and one of the conversion legs
238     if (ConversionMatchFound == kTRUE){
239     isGoodConversion = (conversions->At(ifc)->Prob() > 1e-6) &&
240     (conversions->At(ifc)->Lxy() > 0) &&
241     (conversions->At(ifc)->Position().Rho() > 2.0);
242    
243     if (isGoodConversion == kTRUE) {
244     for (UInt_t d=0; d<conversions->At(ifc)->NDaughters(); d++) {
245     const Track *trk = dynamic_cast<const ChargedParticle*>
246     (conversions->At(ifc)->Daughter(d))->Trk();
247     if (trk) {
248     // These requirements are not used for the GSF track
249     if (!(trk->NHits() >= 3 && trk->Prob() > 1e-6) && trk!=ele->GsfTrk())
250     isGoodConversion = kFALSE;
251     const StableData *sd = dynamic_cast<const StableData*>
252     (conversions->At(ifc)->DaughterDat(d));
253     if (WrongHitsRequirement && sd->NWrongHits() != 0)
254     isGoodConversion = kFALSE;
255     } else {
256     isGoodConversion = kFALSE;
257     }
258     }
259     }
260     }
261     if (isGoodConversion == kTRUE) break;
262    
263     } // loop over all conversions
264    
265     return !isGoodConversion;
266     }
267    
268     //--------------------------------------------------------------------------------------------------
269     Bool_t ElectronTools::PassD0Cut(const Electron *ele, const VertexCol *vertices, Double_t fD0Cut,
270     Bool_t fReverseD0Cut)
271     {
272     Bool_t d0cut = kFALSE;
273     // d0 cut
274     Double_t d0_real = 99999;
275     for(UInt_t i0 = 0; i0 < vertices->GetEntries(); i0++) {
276     Double_t pD0 = ele->GsfTrk()->D0Corrected(*vertices->At(i0));
277     if(TMath::Abs(pD0) < TMath::Abs(d0_real)) d0_real = TMath::Abs(pD0);
278     }
279     if(d0_real < fD0Cut) d0cut = kTRUE;
280    
281     if (fReverseD0Cut == kTRUE &&
282     d0cut == kFALSE && d0_real < 0.05)
283     d0cut = kTRUE;
284     else if(fReverseD0Cut == kTRUE)
285     d0cut = kFALSE;
286    
287     return d0cut;
288     }
289    
290     //--------------------------------------------------------------------------------------------------
291 sixie 1.2 Bool_t ElectronTools::PassD0Cut(const Electron *ele, const BeamSpotCol *beamspots, Double_t fD0Cut,
292     Bool_t fReverseD0Cut)
293     {
294     Bool_t d0cut = kFALSE;
295     // d0 cut
296     Double_t d0_real = 99999;
297     for(UInt_t i0 = 0; i0 < beamspots->GetEntries(); i0++) {
298     Double_t pD0 = ele->GsfTrk()->D0Corrected(*beamspots->At(i0));
299     if(TMath::Abs(pD0) < TMath::Abs(d0_real)) d0_real = TMath::Abs(pD0);
300     }
301     if(d0_real < fD0Cut) d0cut = kTRUE;
302    
303     if (fReverseD0Cut == kTRUE &&
304     d0cut == kFALSE && d0_real < 0.05)
305     d0cut = kTRUE;
306     else if(fReverseD0Cut == kTRUE)
307     d0cut = kFALSE;
308    
309     return d0cut;
310     }
311    
312     //--------------------------------------------------------------------------------------------------
313 sixie 1.1 Bool_t ElectronTools::PassChargeFilter(const Electron *ele)
314     {
315     Bool_t passChargeFilter = kTRUE;
316     if(ele->TrackerTrk() &&
317     ele->TrackerTrk()->Charge() != ele->Charge()) passChargeFilter = kFALSE;
318    
319     return passChargeFilter;
320     }
321    
322     //--------------------------------------------------------------------------------------------------
323     Bool_t ElectronTools::PassSpikeRemovalFilter(const Electron *ele)
324     {
325     Bool_t passSpikeRemovalFilter = kTRUE;
326     if(ele->SCluster()->Seed()->Energy() > 5.0 &&
327     ele->SCluster()->Seed()->EMax() / ele->SCluster()->Seed()->E3x3() > 0.95
328     ) {
329     passSpikeRemovalFilter = kFALSE;
330     }
331    
332     // For Now Only use the EMax/E3x3 prescription.
333     // if(ele->SCluster()->Seed()->Energy() > 5.0 &&
334     // (1 - (ele->SCluster()->Seed()->E1x3() + ele->SCluster()->Seed()->E3x1() - 2*ele->SCluster()->Seed()->EMax())) > 0.95
335     // ) {
336     // passSpikeRemovalFilter = kFALSE;
337     // }
338    
339     return passSpikeRemovalFilter;
340     }