ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitPhysics/Utils/src/MVATools.cc
(Generate patch)

Comparing UserCode/MitPhysics/Utils/src/MVATools.cc (file contents):
Revision 1.2 by mingyang, Thu Oct 13 18:34:50 2011 UTC vs.
Revision 1.10 by bendavid, Mon Dec 19 23:45:00 2011 UTC

# Line 20 | Line 20 | MVATools::MVATools():
20    fReaderEndcap(0),
21    fReaderBarrel(0),
22    
23 <  //MVA variables 0
23 >  //MVA Variables v4
24    HoE(0),
25    covIEtaIEta(0),
26 <  tIso1(0),
27 <  tIso3(0),
28 <  tIso2(0),
26 >  tIso1abs(0),
27 >  tIso3abs(0),
28 >  tIso2abs(0),
29    R9(0),
30    
31 <  //MVA variables 1
32 <  RelIsoEcal(0),
33 <  RelIsoHcal(0),
34 <  
31 >  absIsoEcal(0),
32 >  absIsoHcal(0),
33    RelEMax(0),
34    RelETop(0),
35    RelEBottom(0),
# Line 44 | Line 42 | MVATools::MVATools():
42    RelE2x5Right(0),
43    RelE5x5(0),
44    
47 –  //MVA variables 2
45    EtaWidth(0),
46    PhiWidth(0),
47 <  CoviEtaiPhi(0),
47 >  CoviEtaiPhi(0),
48    CoviPhiiPhi(0),
52 –  RelPreshowerEnergy(0)
49    
50 +  NVertexes(0),
51 +  RelPreshowerEnergy(0),
52 +
53 +  //MVA v2 and v1
54 +  RelIsoEcal(0),
55 +  RelIsoHcal(0),
56 +  tIso1(0),
57 +  tIso3(0),
58 +  tIso2(0),
59 +  ScEta(0.)
60   {
61    // Constructor.
62   }
# Line 69 | Line 75 | void MVATools::InitializeMVA(int Variabl
75    readers[1]  = fReaderBarrel;  
76    
77    for (UInt_t i=0; i<2; ++i) {
78 <    
78 >
79      if(VariableType==0||VariableType==1||VariableType==2){
80        readers[i]->AddVariable( "HoE", &HoE );
81        readers[i]->AddVariable( "covIEtaIEta", &covIEtaIEta );
# Line 79 | Line 85 | void MVATools::InitializeMVA(int Variabl
85        readers[i]->AddVariable( "R9", &R9 );
86      }
87      
88 +    if(VariableType==3||VariableType==4){
89 +      readers[i]->AddVariable( "HoE", &HoE );
90 +      readers[i]->AddVariable( "covIEtaIEta", &covIEtaIEta );
91 +      readers[i]->AddVariable( "tIso1abs", &tIso1abs );
92 +      readers[i]->AddVariable( "tIso3abs", &tIso3abs );
93 +      readers[i]->AddVariable( "tIso2abs", &tIso2abs );
94 +      readers[i]->AddVariable( "R9", &R9 );
95 +    }
96 +    
97      if(VariableType==1||VariableType==2){
98 <      readers[i]->AddVariable( "RelIsoEcal", &RelIsoEcal);
99 <      readers[i]->AddVariable( "RelIsoHcal", &RelIsoHcal);
100 <      readers[i]->AddVariable( "RelEMax", &RelEMax);
101 <      readers[i]->AddVariable( "RelETop",&RelETop);
102 <      readers[i]->AddVariable( "RelEBottom", &RelEBottom);
98 >      readers[i]->AddVariable( "RelIsoEcal", &RelIsoEcal );
99 >      readers[i]->AddVariable( "RelIsoHcal", &RelIsoHcal );
100 >      readers[i]->AddVariable( "RelEMax", &RelEMax );
101 >      readers[i]->AddVariable( "RelETop", &RelETop );
102 >      readers[i]->AddVariable( "RelEBottom", &RelEBottom );
103        readers[i]->AddVariable( "RelELeft", &RelELeft );
104 <      readers[i]->AddVariable( "RelERight", &RelERight);
104 >      readers[i]->AddVariable( "RelERight", &RelERight );
105        readers[i]->AddVariable( "RelE2x5Max", &RelE2x5Max );
106        readers[i]->AddVariable( "RelE2x5Top", &RelE2x5Top );
107        readers[i]->AddVariable( "RelE2x5Bottom", &RelE2x5Bottom );
# Line 95 | Line 110 | void MVATools::InitializeMVA(int Variabl
110        readers[i]->AddVariable( "RelE5x5", &RelE5x5 );
111      }
112      
113 <    if(VariableType==2){
113 >    if(VariableType==3||VariableType==4){
114 >      readers[i]->AddVariable( "absIsoEcal", &absIsoEcal );
115 >      readers[i]->AddVariable( "absIsoHcal", &absIsoHcal );
116 >      readers[i]->AddVariable( "RelEMax", &RelEMax );
117 >      readers[i]->AddVariable( "RelETop", &RelEBottom);
118 >      readers[i]->AddVariable( "RelEBottom", &RelEBottom );
119 >      readers[i]->AddVariable( "RelELeft", &RelELeft );
120 >      readers[i]->AddVariable( "RelERight", &RelERight );
121 >      readers[i]->AddVariable( "RelE2x5Max", &RelE2x5Max );
122 >      readers[i]->AddVariable( "RelE2x5Top", &RelE2x5Top );
123 >      readers[i]->AddVariable( "RelE2x5Bottom", &RelE2x5Bottom );
124 >      readers[i]->AddVariable( "RelE2x5Left", &RelE2x5Left );
125 >      readers[i]->AddVariable( "RelE2x5Right", &RelE2x5Right );
126 >      readers[i]->AddVariable( "RelE5x5", &RelE5x5 );
127 >    }
128 >    
129 >    if(VariableType==2||VariableType==3||VariableType==4){
130        readers[i]->AddVariable( "EtaWidth", &EtaWidth );
131        readers[i]->AddVariable( "PhiWidth", &PhiWidth );
132        readers[i]->AddVariable( "CoviEtaiPhi", &CoviEtaiPhi );
133        readers[i]->AddVariable( "CoviPhiiPhi", &CoviPhiiPhi );
134 +      if(VariableType==4){
135 +        readers[i]->AddVariable( "NVertexes", &NVertexes );
136 +      }
137        if(i==0){
138          readers[i]->AddVariable( "RelPreshowerEnergy", &RelPreshowerEnergy );
139        }
140      }
141 +
142 +    if(VariableType==6){
143 +      readers[i]->AddVariable( "HoE", &HoE );
144 +      readers[i]->AddVariable( "covIEtaIEta", &covIEtaIEta );
145 +      readers[i]->AddVariable( "tIso1abs", &tIso1abs );
146 +      readers[i]->AddVariable( "tIso3abs", &tIso3abs );
147 +      readers[i]->AddVariable( "tIso2abs", &tIso2abs );
148 +      readers[i]->AddVariable( "R9", &R9 );
149 +      readers[i]->AddVariable( "absIsoEcal", &absIsoEcal );
150 +      readers[i]->AddVariable( "absIsoHcal", &absIsoHcal );
151 +      readers[i]->AddVariable( "RelE5x5", &RelE5x5 );
152 +      readers[i]->AddVariable( "EtaWidth", &EtaWidth );
153 +      readers[i]->AddVariable( "PhiWidth", &PhiWidth );
154 +      readers[i]->AddVariable( "CoviEtaiPhi", &CoviEtaiPhi );
155 +      readers[i]->AddVariable( "CoviPhiiPhi", &CoviPhiiPhi );
156 +      readers[i]->AddVariable( "NVertexes", &NVertexes );
157 +      if(i==0){
158 +      readers[i]->AddVariable( "RelPreshowerEnergy", &RelPreshowerEnergy );
159 +      }
160 +    }
161 +    
162 +    if(VariableType==7){
163 +      readers[i]->AddVariable( "HoE", &HoE );
164 +      readers[i]->AddVariable( "covIEtaIEta", &covIEtaIEta );
165 +      readers[i]->AddVariable( "tIso1abs", &tIso1abs );
166 +      readers[i]->AddVariable( "tIso3abs", &tIso3abs );
167 +      readers[i]->AddVariable( "tIso2abs", &tIso2abs );
168 +      readers[i]->AddVariable( "R9", &R9 );
169 +      readers[i]->AddVariable( "absIsoEcal", &absIsoEcal );
170 +      readers[i]->AddVariable( "absIsoHcal", &absIsoHcal );
171 +      readers[i]->AddVariable( "NVertexes", &NVertexes );
172 +      readers[i]->AddVariable( "ScEta", &ScEta );
173 +    }    
174 +
175 +    if(VariableType==10){
176 +      readers[i]->AddVariable( "HoE", &HoE );
177 +      readers[i]->AddVariable( "covIEtaIEta", &covIEtaIEta );
178 +      readers[i]->AddVariable( "tIso1abs", &tIso1abs );
179 +      readers[i]->AddVariable( "tIso3abs", &tIso3abs );
180 +      readers[i]->AddVariable( "tIso2abs", &tIso2abs );
181 +      readers[i]->AddVariable( "R9", &R9 );
182 +      readers[i]->AddVariable( "absIsoEcal", &absIsoEcal );
183 +      readers[i]->AddVariable( "absIsoHcal", &absIsoHcal );
184 +      readers[i]->AddVariable( "NVertexes", &NVertexes );
185 +      readers[i]->AddVariable( "ScEta", &ScEta );
186 +      readers[i]->AddVariable( "EtaWidth", &EtaWidth );
187 +      readers[i]->AddVariable( "PhiWidth", &PhiWidth );      
188 +    }    
189 +
190    }
191    
192    fReaderEndcap->BookMVA("BDT method",EndcapWeights);
193    fReaderBarrel->BookMVA("BDT method",BarrelWeights);
194 <
194 >  
195    assert(fReaderEndcap);
196    assert(fReaderBarrel);
197    
198   }
199  
200 < //--------------------------------------------------------------------------------------------------
118 <
119 < Bool_t MVATools::PassMVASelection(const Photon* p,const Vertex* vtx,const TrackCol* trackCol,const VertexCol* vtxCol,Double_t _tRho,const ElectronCol* els,double MVAPtMin, Float_t bdtCutBarrel, Float_t bdtCutEndcap) {
120 <  
121 <  // these values are taken from the H2GGlobe code... (actually from Marco/s mail)
122 <  float cic4_allcuts_temp_sublead[] = {
123 <    3.8,         2.2,         1.77,        1.29,
124 <    11.7,        3.4,         3.9,         1.84,
125 <    3.5,         2.2,         2.3,         1.45,
126 <    0.0106,      0.0097,      0.028,       0.027,
127 <    0.082,       0.062,       0.065,       0.048,
128 <    0.94,        0.36,        0.94,        0.32,
129 <    1.,          0.062,       0.97,        0.97,
130 <    1.5,         1.5,         1.5,         1.5 };  // the last line is PixelmatchVeto and un-used
131 <  
200 > Bool_t MVATools::PassMVASelection(const Photon* p,const Vertex* vtx,const TrackCol* trackCol,const VertexCol* vtxCol,Double_t _tRho,Float_t bdtCutBarrel, Float_t bdtCutEndcap, const ElectronCol* els, Bool_t applyElectronVeto) {
201    
202    //initilize the bool value
203    PassMVA=kFALSE;
135 –  PassElecVeto=kFALSE;
136 –  //Bool_t PassElecVeto=PhotonTools::PassElectronVeto(p,els);
204    
205 +  Float_t photon_bdt =  MVATools::GetMVAbdtValue(p,vtx,trackCol,vtxCol, _tRho, els, applyElectronVeto);
206    
207 <  //get the variables used to compute MVA variables
208 <  ecalIso3 = p->EcalRecHitIsoDr03();
209 <  ecalIso4 = p->EcalRecHitIsoDr04();
142 <  hcalIso4 = p->HcalTowerSumEtDr04();
143 <  
144 <  wVtxInd = 0;
145 <  
146 <  trackIso1 = IsolationTools::CiCTrackIsolation(p,vtx, 0.3, 0.02, 0.0, 0.0, 0.1, 1.0, trackCol);//Question Ming:whyfPV->At(0) instead of selected vertex using ranking method?
147 <  
148 <  // track iso only
149 <  trackIso3 = IsolationTools::CiCTrackIsolation(p,vtx, 0.3, 0.02, 0.0, 0.0, 0.1, 1.0, trackCol);
150 <  
151 <  // track iso worst vtx
152 <  trackIso2 = IsolationTools::CiCTrackIsolation(p,vtx, 0.4, 0.02, 0.0, 0.0, 0.1, 1.0, trackCol, &wVtxInd,vtxCol);
153 <  
154 <  combIso1 = ecalIso3+hcalIso4+trackIso1 - 0.17*_tRho;
155 <  combIso2 = ecalIso4+hcalIso4+trackIso2 - 0.52*_tRho;
156 <  
157 <  RawEnergy = p->SCluster()->RawEnergy();
158 <  
159 <  dRTrack = PhotonTools::ElectronVetoCiC(p, els);
160 <  
161 <  //compute MVA variables 0
162 <  HoE = p->HadOverEm();
163 <  covIEtaIEta = p->CoviEtaiEta();
164 <  tIso1 = (combIso1) *50./p->Et();
165 <  tIso3 = (trackIso3)*50./p->Et();
166 <  tIso2 = (combIso2) *50./(p->MomVtx(vtxCol->At(wVtxInd)->Position()).Pt());
167 <  R9 = p->R9();
168 <  
169 <  //newly added MVA variables 1
170 <  RelIsoEcal=(ecalIso3-0.17*_tRho)/p->Et();
171 <  RelIsoHcal=(hcalIso4-0.17*_tRho)/p->Et();
172 <  
173 <  RelEMax=p->SCluster()->Seed()->EMax()/RawEnergy;
174 <  RelETop=p->SCluster()->Seed()->ETop()/RawEnergy;
175 <  RelEBottom=p->SCluster()->Seed()->EBottom()/RawEnergy;
176 <  RelELeft=p->SCluster()->Seed()->ELeft()/RawEnergy;
177 <  RelERight=p->SCluster()->Seed()->ERight()/RawEnergy;
178 <  RelE2x5Max=p->SCluster()->Seed()->E2x5Max()/RawEnergy;
179 <  RelE2x5Top=p->SCluster()->Seed()->E2x5Top()/RawEnergy;
180 <  RelE2x5Bottom=p->SCluster()->Seed()->E2x5Bottom()/RawEnergy;
181 <  RelE2x5Left=p->SCluster()->Seed()->E2x5Left()/RawEnergy;
182 <  RelE2x5Right=p->SCluster()->Seed()->E2x5Right()/RawEnergy;
183 <  RelE5x5=p->SCluster()->Seed()->E5x5()/RawEnergy;
184 <  
185 <  //newly added MVA variables 2
186 <  EtaWidth=p->SCluster()->EtaWidth();
187 <  PhiWidth=p->SCluster()->PhiWidth();
188 <  CoviEtaiPhi=p->SCluster()->Seed()->CoviEtaiPhi();
189 <  CoviPhiiPhi=p->SCluster()->Seed()->CoviPhiiPhi();
190 <  RelPreshowerEnergy=p->SCluster()->PreshowerEnergy()/RawEnergy;
191 <  
192 <  //spectator variables
193 <  Pt_MVA=p->Pt();
194 <  ScEta_MVA=p->SCluster()->Eta();
195 <  
196 <  isbarrel = (fabs(ScEta_MVA)<1.4442);
197 <  
198 <  // check which category it is ...
199 <  _tCat = 1;
200 <  if ( !isbarrel ) _tCat = 3;
201 <  if ( R9 < 0.94 ) _tCat++;
202 <  
203 <  //Electron Veto
204 <  if(dRTrack > cic4_allcuts_temp_sublead[_tCat-1+6*4]){
205 <    PassElecVeto=kTRUE;
206 <  }
207 <  
208 <  if(PassElecVeto && Pt_MVA>MVAPtMin && ((fabs(ScEta_MVA)<1.4442)||(fabs(ScEta_MVA)>1.566 && fabs(ScEta_MVA)<2.5)) && tIso1<250 && tIso2<250 && tIso3<250){
209 <    
210 <    if (isbarrel) {
211 <      reader = fReaderBarrel;
212 <    }
213 <    else {
214 <      reader = fReaderEndcap;
215 <    }
216 <    
217 <    bdt = reader->EvaluateMVA( "BDT method" );
218 <    
219 <    if (isbarrel) {
220 <      if(bdt>bdtCutBarrel){
221 <        PassMVA=kTRUE;  
222 <      }
207 >  if (isbarrel) {
208 >    if(bdt>bdtCutBarrel){
209 >      PassMVA=kTRUE;    
210      }
224 –    else {
225 –      if(bdt>bdtCutEndcap){
226 –        PassMVA=kTRUE;  
227 –      }
228 –    }
229 –
211    }
212 <
212 >  else {
213 >    if(bdt>bdtCutEndcap){
214 >      PassMVA=kTRUE;    
215 >    }
216 >  }
217    return PassMVA;
218   }
219  
# Line 256 | Line 241 | Int_t MVATools::PassElectronVetoInt(cons
241    isbarrel = (fabs(ScEta_MVA)<1.4442);
242  
243    R9 = p->R9();
244 +  //R9 = p->E33()/p->SCluster()->RawEnergy();
245    
246    // check which category it is ...
247    _tCat = 1;
# Line 273 | Line 259 | Int_t MVATools::PassElectronVetoInt(cons
259  
260   //--------------------------------------------------------------------------------------------------
261  
262 < Float_t MVATools::GetMVAbdtValue(const Photon* p,const Vertex* vtx,const TrackCol* trackCol,const VertexCol* vtxCol,Double_t _tRho,const ElectronCol* els) {
262 > Float_t MVATools::GetMVAbdtValue(const Photon* p,const Vertex* vtx,const TrackCol* trackCol,const VertexCol* vtxCol,Double_t _tRho, const ElectronCol* els, Bool_t applyElectronVeto) {
263    
264    //get the variables used to compute MVA variables
265    ecalIso3 = p->EcalRecHitIsoDr03();
# Line 282 | Line 268 | Float_t MVATools::GetMVAbdtValue(const P
268    
269    wVtxInd = 0;
270    
271 <  trackIso1 = IsolationTools::CiCTrackIsolation(p,vtx, 0.3, 0.02, 0.0, 0.0, 0.1, 1.0, trackCol);//Question Ming:whyfPV->At(0) instead of selected vertex using ranking method?
272 <  
287 <  // track iso only
288 <  trackIso3 = IsolationTools::CiCTrackIsolation(p,vtx, 0.3, 0.02, 0.0, 0.0, 0.1, 1.0, trackCol);
289 <  
271 >  trackIso1 = IsolationTools::CiCTrackIsolation(p,vtx, 0.3, 0.02, 0.0, 0.0, 0.1, 1.0, trackCol, NULL, NULL, (!applyElectronVeto ? els : NULL) );//Question Ming:whyfPV->At(0) instead of selected vertex using ranking method?
272 >    
273    // track iso worst vtx
274 <  trackIso2 = IsolationTools::CiCTrackIsolation(p,vtx, 0.4, 0.02, 0.0, 0.0, 0.1, 1.0, trackCol, &wVtxInd,vtxCol);
274 >  trackIso2 = IsolationTools::CiCTrackIsolation(p,vtx, 0.4, 0.02, 0.0, 0.0, 0.1, 1.0, trackCol, &wVtxInd,vtxCol, (!applyElectronVeto ? els : NULL) );
275    
276    combIso1 = ecalIso3+hcalIso4+trackIso1 - 0.17*_tRho;
277    combIso2 = ecalIso4+hcalIso4+trackIso2 - 0.52*_tRho;
278    
279    RawEnergy = p->SCluster()->RawEnergy();
280    
281 <  dRTrack = PhotonTools::ElectronVetoCiC(p, els);
281 >  ScEta = p->SCluster()->Eta();
282    
283 <  //compute MVA variables 0
301 <  HoE = p->HadOverEm();
302 <  covIEtaIEta = p->CoviEtaiEta();
283 >  //mva varialbes v1 and v2
284    tIso1 = (combIso1) *50./p->Et();
285 <  tIso3 = (trackIso3)*50./p->Et();
285 >  tIso3 = (trackIso1)*50./p->Et();
286    tIso2 = (combIso2) *50./(p->MomVtx(vtxCol->At(wVtxInd)->Position()).Pt());
306 –  R9 = p->R9();
307 –  
308 –  //newly added MVA variables 1
287    RelIsoEcal=(ecalIso3-0.17*_tRho)/p->Et();
288    RelIsoHcal=(hcalIso4-0.17*_tRho)/p->Et();
289 <  
289 >
290 >  //compute mva variables for v3
291 >  HoE = p->HadOverEm();
292 >  covIEtaIEta = p->CoviEtaiEta();
293 >  tIso1abs = combIso1;
294 >  tIso3abs = trackIso1;
295 >  tIso2abs = combIso2;
296 >  R9 = p->R9();
297 >
298 >  absIsoEcal=(ecalIso3-0.17*_tRho);
299 >  absIsoHcal=(hcalIso4-0.17*_tRho);
300    RelEMax=p->SCluster()->Seed()->EMax()/RawEnergy;
301    RelETop=p->SCluster()->Seed()->ETop()/RawEnergy;
302    RelEBottom=p->SCluster()->Seed()->EBottom()/RawEnergy;
# Line 321 | Line 309 | Float_t MVATools::GetMVAbdtValue(const P
309    RelE2x5Right=p->SCluster()->Seed()->E2x5Right()/RawEnergy;
310    RelE5x5=p->SCluster()->Seed()->E5x5()/RawEnergy;
311    
312 <  //newly added MVA variables 2
313 <  EtaWidth=p->SCluster()->EtaWidth();
326 <  PhiWidth=p->SCluster()->PhiWidth();
312 >  EtaWidth=p->EtaWidth();
313 >  PhiWidth=p->PhiWidth();
314    CoviEtaiPhi=p->SCluster()->Seed()->CoviEtaiPhi();
315    CoviPhiiPhi=p->SCluster()->Seed()->CoviPhiiPhi();
316 +
317    RelPreshowerEnergy=p->SCluster()->PreshowerEnergy()/RawEnergy;
318 +  NVertexes=vtxCol->GetEntries();
319    
320    //spectator variables
321    Pt_MVA=p->Pt();
# Line 345 | Line 334 | Float_t MVATools::GetMVAbdtValue(const P
334  
335    bdt = reader->EvaluateMVA("BDT method");
336  
337 <  /*  printf("HoE: %f\n",HoE);
337 >  /* printf("HoE: %f\n",HoE);
338    printf("covIEtaIEta: %f\n",covIEtaIEta);
339 <  printf("tIso1: %f\n",tIso1);
340 <  printf("tIso3: %f\n",tIso3);
341 <  printf("tIso2: %f\n",tIso2);
353 <  printf("tIso1: %f\n",tIso1);
354 <  printf("R9: %f\n",R9);
355 <  printf("RelIsoEcal: %f\n",RelIsoEcal);
356 <  printf("RelIsoHcal: %f\n",RelIsoHcal);
339 >  printf("tIso1abs: %f\n",tIso1abs);
340 >  printf("tIso3abs: %f\n",tIso3abs);
341 >  printf("tIso2abs: %f\n",tIso2abs);
342    
343 +  printf("absIsoEcal: %f\n",absIsoEcal);
344 +  printf("absIsoHcal: %f\n",absIsoHcal);
345    printf("RelEMax: %f\n",RelEMax);
346    printf("RelETop: %f\n",RelETop);
347    printf("RelEBottom: %f\n",RelEBottom);
# Line 372 | Line 359 | Float_t MVATools::GetMVAbdtValue(const P
359    printf("CoviEtaiPhi: %f\n",CoviEtaiPhi);
360    printf("CoviPhiiPhi: %f\n",CoviPhiiPhi);
361    
362 <  if (isbarrel) {
362 >  if (!isbarrel) {
363      printf("RelPreshowerEnergy: %f\n",RelPreshowerEnergy);
364      }*/
365    

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines