ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitPhysics/Utils/src/MVATools.cc
Revision: 1.18
Committed: Thu Aug 2 13:57:32 2012 UTC (12 years, 9 months ago) by fabstoec
Content type: text/plain
Branch: MAIN
Changes since 1.17: +71 -26 lines
Log Message: 
bug fix

File Contents

# User Rev Content
1 fabstoec 1.18 // $Id: MVATools.cc,v 1.17 2012/08/02 12:30:55 fabstoec Exp $
2 mingyang 1.1
3     #include "MitPhysics/Utils/interface/PhotonTools.h"
4     #include "MitPhysics/Utils/interface/MVATools.h"
5     #include "MitPhysics/Utils/interface/ElectronTools.h"
6     #include "MitPhysics/Utils/interface/IsolationTools.h"
7 mingyang 1.11 #include "MitAna/DataTree/interface/PFCandidateCol.h"
8 mingyang 1.1 #include "MitAna/DataTree/interface/StableData.h"
9 mingyang 1.11 #include <TMath.h>
10 mingyang 1.1 #include <TFile.h>
11     #include <TRandom3.h>
12 fabstoec 1.17 #include <TSystem.h>
13 mingyang 1.1 #include "TMVA/Tools.h"//MVA
14     #include "TMVA/Reader.h"//MVA
15    
16 mingyang 1.11
17 mingyang 1.1 ClassImp(mithep::MVATools)
18    
19     using namespace mithep;
20    
21    
22     //--------------------------------------------------------------------------------------------------
23     MVATools::MVATools():
24     fReaderEndcap(0),
25     fReaderBarrel(0),
26    
27 fabstoec 1.17 fMVAType (MVATools::kNone)
28    
29     // ------------------------------------------------------------------------------
30     // fab: These guys should all go away....
31 mingyang 1.4 //MVA Variables v4
32 fabstoec 1.17 // HoE(0),
33     // covIEtaIEta(0),
34     // tIso1abs(0),
35     // tIso3abs(0),
36     // tIso2abs(0),
37     // R9(0),
38    
39     // absIsoEcal(0),
40     // absIsoHcal(0),
41     // RelEMax(0),
42     // RelETop(0),
43     // RelEBottom(0),
44     // RelELeft(0),
45     // RelERight(0),
46     // RelE2x5Max(0),
47     // RelE2x5Top(0),
48     // RelE2x5Bottom(0),
49     // RelE2x5Left(0),
50     // RelE2x5Right(0),
51     // RelE5x5(0),
52    
53     // EtaWidth(0),
54     // PhiWidth(0),
55     // CoviEtaiPhi(0),
56     // CoviPhiiPhi(0),
57    
58     // NVertexes(0),
59     // RelPreshowerEnergy(0),
60    
61     // //MVA v2 and v1
62     // RelIsoEcal(0),
63     // RelIsoHcal(0),
64     // tIso1(0),
65     // tIso3(0),
66     // tIso2(0),
67     // ScEta(0.)
68    
69 mingyang 1.1 {
70     // Constructor.
71     }
72    
73     //--------------------------------------------------------------------------------------------------
74 fabstoec 1.17 void MVATools::InitializeMVA(MVATools::IdMVAType type) {
75 mingyang 1.1
76 fabstoec 1.17 fMVAType = type;
77    
78 mingyang 1.1 if (fReaderEndcap) delete fReaderEndcap;
79     if (fReaderBarrel) delete fReaderBarrel;
80 fabstoec 1.17
81     // no MVA needed if none requested
82     if( type == kNone ) {
83     return;
84     }
85 mingyang 1.1
86     fReaderEndcap = new TMVA::Reader( "!Color:!Silent:Error" );
87     fReaderBarrel = new TMVA::Reader( "!Color:!Silent:Error" );
88 fabstoec 1.17
89     TString BarrelWeights;
90     TString EndcapWeights;
91    
92     varNames.resize(0);
93 mingyang 1.1
94 fabstoec 1.17 switch (type) {
95 mingyang 1.1
96 fabstoec 1.18 case k2011IdMVA_HZg:
97    
98     EndcapWeights = (gSystem->Getenv("CMSSW_BASE")+
99     TString("/src/MitPhysics/data/")+
100     TString("PhotonId_lowPt_EE_BDTG.")+
101     TString("weights.xml"));
102     BarrelWeights = (gSystem->Getenv("CMSSW_BASE")+
103     TString("/src/MitPhysics/data/")+
104     TString("PhotonId_lowPt_EB_BDTG.")+
105     TString("weights.xml"));
106    
107    
108     varNames.push_back("sigieie");
109     varNames.push_back("covieip");
110     varNames.push_back("s4r" );
111     varNames.push_back("r9" );
112     varNames.push_back("sigeta" );
113     varNames.push_back("sigphi" );
114     varNames.push_back("pfgios" );
115     varNames.push_back("pfciso" );
116     varNames.push_back("rho" );
117     varNames.push_back("sceta" );
118     varNames.push_back("rawe" );
119    
120     mvaVars.resize(varNames.size());
121    
122     std::cout<<" Adding stuff here.... "<<std::endl;
123    
124     for( unsigned int iV = 0; iV < mvaVars.size(); ++iV) {
125     mvaVarMapEB.insert( std::pair<std::string,unsigned int>(varNames[iV], iV) );
126     mvaVarMapEE.insert( std::pair<std::string,unsigned int>(varNames[iV], iV) );
127     }
128    
129     std::cout<<" ... done "<<std::endl;
130    
131     break;
132    
133 fabstoec 1.17 case k2011IdMVA:
134    
135     EndcapWeights = (gSystem->Getenv("CMSSW_BASE")+
136     TString("/src/MitPhysics/data/TMVAClassificationPhotonID_")+
137     TString("Endcap_PassPreSel_Variable_10_BDTnCuts2000_BDT.")+
138     TString("weights.xml"));
139     BarrelWeights = (gSystem->Getenv("CMSSW_BASE")+
140     TString("/src/MitPhysics/data/TMVAClassificationPhotonID_")+
141     TString("Barrel_PassPreSel_Variable_10_BDTnCuts2000_BDT.")+
142     TString("weights.xml"));
143    
144     // set up the variable names
145     mvaVars.resize(12);
146     varNames.push_back("HoE" );
147     varNames.push_back("covIEtaIEta");
148     varNames.push_back("tIso1abs" );
149     varNames.push_back("tIso3abs" );
150     varNames.push_back("tIso2abs" );
151     varNames.push_back("R9" );
152     varNames.push_back("absIsoEcal" );
153     varNames.push_back("absIsoHcal" );
154     varNames.push_back("NVertexes" );
155     varNames.push_back("ScEta" );
156     varNames.push_back("EtaWidth" );
157     varNames.push_back("PhiWidth" );
158    
159     for( unsigned int iV = 0; iV < mvaVars.size(); ++iV) {
160 fabstoec 1.18 mvaVarMapEB.insert( std::pair<std::string,unsigned int>(varNames[iV], iV) );
161     mvaVarMapEE.insert( std::pair<std::string,unsigned int>(varNames[iV], iV) );
162 mingyang 1.4 }
163    
164 fabstoec 1.17 break;
165    
166     case k2012IdMVA:
167     case k2012IdMVA_globe:
168 mingyang 1.4
169 fabstoec 1.17 EndcapWeights = (gSystem->Getenv("CMSSW_BASE")+
170     TString("/src/MitPhysics/data/")+
171     TString("2012ICHEP_PhotonID_Endcap_BDT.")+
172     TString("weights_PSCorr.xml"));
173     BarrelWeights = (gSystem->Getenv("CMSSW_BASE")+
174     TString("/src/MitPhysics/data/")+
175     TString("2012ICHEP_PhotonID_Barrel_BDT.")+
176     TString("weights.xml"));
177    
178     mvaVars.resize(12);
179     varNames.push_back("ph.r9" );
180     varNames.push_back("ph.sigietaieta" );
181     varNames.push_back("ph.scetawidth" );
182     varNames.push_back("ph.scphiwidth" );
183     varNames.push_back("ph.idmva_CoviEtaiPhi" );
184     varNames.push_back("ph.idmva_s4ratio" );
185     varNames.push_back("ph.idmva_GammaIso" );
186     varNames.push_back("ph.idmva_ChargedIso_selvtx" );
187     varNames.push_back("ph.idmva_ChargedIso_worstvtx" );
188     varNames.push_back("ph.sceta" );
189     varNames.push_back("rho" );
190     varNames.push_back("ph.idmva_PsEffWidthSigmaRR" );
191    
192     for( unsigned int iV = 0; iV < mvaVars.size() - 1; ++iV) {
193 fabstoec 1.18 mvaVarMapEB.insert( std::pair<std::string,unsigned int>(varNames[iV], iV) );
194     mvaVarMapEE.insert( std::pair<std::string,unsigned int>(varNames[iV], iV) );
195 mingyang 1.1 }
196    
197 fabstoec 1.17 // pre-shower only used for Endcaps
198 fabstoec 1.18 mvaVarMapEE.insert( std::pair<std::string,unsigned int> ( varNames[mvaVars.size() - 1], mvaVars.size() - 1) );
199 mingyang 1.6
200 fabstoec 1.17 break;
201 mingyang 1.11
202 fabstoec 1.17 default:
203     // no variables... better never called..
204     std::cerr<<" MVATools: Trying to initialize with unknown type."<<std::endl;
205     break;
206 mingyang 1.1 }
207 fabstoec 1.17
208    
209     // looping over both maps and adding Vars to BDT readers
210 fabstoec 1.18 for( unsigned int iV = 0; iV < varNames.size(); ++iV ){
211     std::map<std::string,unsigned int>::const_iterator it = mvaVarMapEB.find(varNames[iV]);
212     if ( it != mvaVarMapEB.end() )
213     fReaderBarrel -> AddVariable( (it->first).c_str(), &(mvaVars[it->second]));
214     it = mvaVarMapEE.find(varNames[iV]);
215     if ( it != mvaVarMapEE.end() )
216     fReaderEndcap -> AddVariable( (it->first).c_str(), &(mvaVars[it->second]));
217     }
218    
219 mingyang 1.1 fReaderEndcap->BookMVA("BDT method",EndcapWeights);
220     fReaderBarrel->BookMVA("BDT method",BarrelWeights);
221 mingyang 1.4
222 mingyang 1.1 assert(fReaderEndcap);
223     assert(fReaderBarrel);
224    
225     }
226    
227 fabstoec 1.17 // //--------------------------------------------------------------------------------------------------
228     // void MVATools::InitializeMVA(int VariableType, TString EndcapWeights,TString BarrelWeights) {
229    
230     // if (fReaderEndcap) delete fReaderEndcap;
231     // if (fReaderBarrel) delete fReaderBarrel;
232    
233     // fReaderEndcap = new TMVA::Reader( "!Color:!Silent:Error" );
234     // fReaderBarrel = new TMVA::Reader( "!Color:!Silent:Error" );
235    
236     // TMVA::Reader *readers[2];
237     // readers[0] = fReaderEndcap;
238     // readers[1] = fReaderBarrel;
239    
240     // for (UInt_t i=0; i<2; ++i) {
241    
242     // if(VariableType==0||VariableType==1||VariableType==2){
243     // readers[i]->AddVariable( "HoE", &HoE );
244     // readers[i]->AddVariable( "covIEtaIEta", &covIEtaIEta );
245     // readers[i]->AddVariable( "tIso1", &tIso1 );
246     // readers[i]->AddVariable( "tIso3", &tIso3 );
247     // readers[i]->AddVariable( "tIso2", &tIso2 );
248     // readers[i]->AddVariable( "R9", &R9 );
249     // }
250    
251     // if(VariableType==3||VariableType==4){
252     // readers[i]->AddVariable( "HoE", &HoE );
253     // readers[i]->AddVariable( "covIEtaIEta", &covIEtaIEta );
254     // readers[i]->AddVariable( "tIso1abs", &tIso1abs );
255     // readers[i]->AddVariable( "tIso3abs", &tIso3abs );
256     // readers[i]->AddVariable( "tIso2abs", &tIso2abs );
257     // readers[i]->AddVariable( "R9", &R9 );
258     // }
259    
260     // if(VariableType==1||VariableType==2){
261     // readers[i]->AddVariable( "RelIsoEcal", &RelIsoEcal );
262     // readers[i]->AddVariable( "RelIsoHcal", &RelIsoHcal );
263     // readers[i]->AddVariable( "RelEMax", &RelEMax );
264     // readers[i]->AddVariable( "RelETop", &RelETop );
265     // readers[i]->AddVariable( "RelEBottom", &RelEBottom );
266     // readers[i]->AddVariable( "RelELeft", &RelELeft );
267     // readers[i]->AddVariable( "RelERight", &RelERight );
268     // readers[i]->AddVariable( "RelE2x5Max", &RelE2x5Max );
269     // readers[i]->AddVariable( "RelE2x5Top", &RelE2x5Top );
270     // readers[i]->AddVariable( "RelE2x5Bottom", &RelE2x5Bottom );
271     // readers[i]->AddVariable( "RelE2x5Left", &RelE2x5Left );
272     // readers[i]->AddVariable( "RelE2x5Right", &RelE2x5Right );
273     // readers[i]->AddVariable( "RelE5x5", &RelE5x5 );
274     // }
275    
276     // if(VariableType==3||VariableType==4){
277     // readers[i]->AddVariable( "absIsoEcal", &absIsoEcal );
278     // readers[i]->AddVariable( "absIsoHcal", &absIsoHcal );
279     // readers[i]->AddVariable( "RelEMax", &RelEMax );
280     // readers[i]->AddVariable( "RelETop", &RelEBottom);
281     // readers[i]->AddVariable( "RelEBottom", &RelEBottom );
282     // readers[i]->AddVariable( "RelELeft", &RelELeft );
283     // readers[i]->AddVariable( "RelERight", &RelERight );
284     // readers[i]->AddVariable( "RelE2x5Max", &RelE2x5Max );
285     // readers[i]->AddVariable( "RelE2x5Top", &RelE2x5Top );
286     // readers[i]->AddVariable( "RelE2x5Bottom", &RelE2x5Bottom );
287     // readers[i]->AddVariable( "RelE2x5Left", &RelE2x5Left );
288     // readers[i]->AddVariable( "RelE2x5Right", &RelE2x5Right );
289     // readers[i]->AddVariable( "RelE5x5", &RelE5x5 );
290     // }
291    
292     // if(VariableType==2||VariableType==3||VariableType==4){
293     // readers[i]->AddVariable( "EtaWidth", &EtaWidth );
294     // readers[i]->AddVariable( "PhiWidth", &PhiWidth );
295     // readers[i]->AddVariable( "CoviEtaiPhi", &CoviEtaiPhi );
296     // readers[i]->AddVariable( "CoviPhiiPhi", &CoviPhiiPhi );
297     // if(VariableType==4){
298     // readers[i]->AddVariable( "NVertexes", &NVertexes );
299     // }
300     // if(i==0){
301     // readers[i]->AddVariable( "RelPreshowerEnergy", &RelPreshowerEnergy );
302     // }
303     // }
304    
305     // if(VariableType==6){
306     // readers[i]->AddVariable( "HoE", &HoE );
307     // readers[i]->AddVariable( "covIEtaIEta", &covIEtaIEta );
308     // readers[i]->AddVariable( "tIso1abs", &tIso1abs );
309     // readers[i]->AddVariable( "tIso3abs", &tIso3abs );
310     // readers[i]->AddVariable( "tIso2abs", &tIso2abs );
311     // readers[i]->AddVariable( "R9", &R9 );
312     // readers[i]->AddVariable( "absIsoEcal", &absIsoEcal );
313     // readers[i]->AddVariable( "absIsoHcal", &absIsoHcal );
314     // readers[i]->AddVariable( "RelE5x5", &RelE5x5 );
315     // readers[i]->AddVariable( "EtaWidth", &EtaWidth );
316     // readers[i]->AddVariable( "PhiWidth", &PhiWidth );
317     // readers[i]->AddVariable( "CoviEtaiPhi", &CoviEtaiPhi );
318     // readers[i]->AddVariable( "CoviPhiiPhi", &CoviPhiiPhi );
319     // readers[i]->AddVariable( "NVertexes", &NVertexes );
320     // if(i==0){
321     // readers[i]->AddVariable( "RelPreshowerEnergy", &RelPreshowerEnergy );
322     // }
323     // }
324    
325     // if(VariableType==7){
326     // readers[i]->AddVariable( "HoE", &HoE );
327     // readers[i]->AddVariable( "covIEtaIEta", &covIEtaIEta );
328     // readers[i]->AddVariable( "tIso1abs", &tIso1abs );
329     // readers[i]->AddVariable( "tIso3abs", &tIso3abs );
330     // readers[i]->AddVariable( "tIso2abs", &tIso2abs );
331     // readers[i]->AddVariable( "R9", &R9 );
332     // readers[i]->AddVariable( "absIsoEcal", &absIsoEcal );
333     // readers[i]->AddVariable( "absIsoHcal", &absIsoHcal );
334     // readers[i]->AddVariable( "NVertexes", &NVertexes );
335     // readers[i]->AddVariable( "ScEta", &ScEta );
336     // }
337    
338     // if(VariableType==10){
339     // readers[i]->AddVariable( "HoE", &HoE );
340     // readers[i]->AddVariable( "covIEtaIEta", &covIEtaIEta );
341     // readers[i]->AddVariable( "tIso1abs", &tIso1abs );
342     // readers[i]->AddVariable( "tIso3abs", &tIso3abs );
343     // readers[i]->AddVariable( "tIso2abs", &tIso2abs );
344     // readers[i]->AddVariable( "R9", &R9 );
345     // readers[i]->AddVariable( "absIsoEcal", &absIsoEcal );
346     // readers[i]->AddVariable( "absIsoHcal", &absIsoHcal );
347     // readers[i]->AddVariable( "NVertexes", &NVertexes );
348     // readers[i]->AddVariable( "ScEta", &ScEta );
349     // readers[i]->AddVariable( "EtaWidth", &EtaWidth );
350     // readers[i]->AddVariable( "PhiWidth", &PhiWidth );
351     // }
352    
353     // if(VariableType==1201){
354     // /*readers[i]->AddVariable( "myphoton_pfchargedisogood03", &myphoton_pfchargedisogood03);
355     // readers[i]->AddVariable( "myphoton_pfchargedisobad03", &myphoton_pfchargedisobad03);
356     // readers[i]->AddVariable( "myphoton_pfphotoniso03", &myphoton_pfphotoniso03 );
357     // readers[i]->AddVariable( "myphoton_sieie", &myphoton_sieie );
358     // readers[i]->AddVariable( "myphoton_sieip", &myphoton_sieip );
359     // readers[i]->AddVariable( "myphoton_etawidth", &myphoton_etawidth );
360     // readers[i]->AddVariable( "myphoton_phiwidth", &myphoton_phiwidth );
361     // readers[i]->AddVariable( "myphoton_r9", &myphoton_r9 );
362     // readers[i]->AddVariable( "myphoton_s4ratio", &myphoton_s4ratio );
363     // readers[i]->AddVariable( "myphoton_SCeta", &myphoton_SCeta );
364     // readers[i]->AddVariable( "event_rho", &event_rho );
365     // if(i==0){
366     // readers[i]->AddVariable( "myphoton_ESEffSigmaRR", &myphoton_ESEffSigmaRR);
367     // }*/
368     // readers[i]->AddVariable( "ph.r9", &myphoton_r9 );
369     // readers[i]->AddVariable( "ph.sigietaieta", &myphoton_sieie );
370     // readers[i]->AddVariable( "ph.scetawidth", &myphoton_etawidth );
371     // readers[i]->AddVariable( "ph.scphiwidth", &myphoton_phiwidth );
372     // readers[i]->AddVariable( "ph.idmva_CoviEtaiPhi", &myphoton_sieip );
373     // readers[i]->AddVariable( "ph.idmva_s4ratio", &myphoton_s4ratio );
374     // readers[i]->AddVariable( "ph.idmva_GammaIso", &myphoton_pfphotoniso03 );
375     // readers[i]->AddVariable( "ph.idmva_ChargedIso_selvtx", &myphoton_pfchargedisogood03);
376     // readers[i]->AddVariable( "ph.idmva_ChargedIso_worstvtx", &myphoton_pfchargedisobad03);
377     // readers[i]->AddVariable( "ph.sceta", &myphoton_SCeta );
378     // readers[i]->AddVariable( "rho", &event_rho );
379     // if(i==0){
380     // //readers[i]->AddVariable( "1.00023*ph.idmva_PsEffWidthSigmaRR + 0.0913", &myphoton_ESEffSigmaRR);
381     // readers[i]->AddVariable( "ph.idmva_PsEffWidthSigmaRR", &myphoton_ESEffSigmaRR);
382     // }
383     // }
384    
385     // }
386    
387     // fReaderEndcap->BookMVA("BDT method",EndcapWeights);
388     // fReaderBarrel->BookMVA("BDT method",BarrelWeights);
389    
390     // assert(fReaderEndcap);
391     // assert(fReaderBarrel);
392    
393     // }
394    
395     // *** REMOVED THIS COMPLETELY. If a module wants to cut on BDT, it should 1,) compute the BDT value (using GetMVAbdtValue(...) ) and then make the cut itself...
396    
397     // 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) {
398    
399     // //initilize the bool value
400     // PassMVA=kFALSE;
401    
402     // Float_t photon_bdt = MVATools::GetMVAbdtValue_2011(p,vtx,trackCol,vtxCol, _tRho, els, applyElectronVeto);
403    
404     // if (isbarrel) {
405     // if(photon_bdt>bdtCutBarrel){
406     // PassMVA=kTRUE;
407     // }
408     // }
409     // else {
410     // if(photon_bdt>bdtCutEndcap){
411     // PassMVA=kTRUE;
412     // }
413     // }
414     // return PassMVA;
415     // }
416    
417     // //---------------------------------------------------------------------------------
418     // Int_t MVATools::PassElectronVetoInt(const Photon* p, const ElectronCol* els) {
419    
420     // // these values are taken from the H2GGlobe code... (actually from Marco/s mail)
421     // float cic4_allcuts_temp_sublead[] = {
422     // 3.8, 2.2, 1.77, 1.29,
423     // 11.7, 3.4, 3.9, 1.84,
424     // 3.5, 2.2, 2.3, 1.45,
425     // 0.0106, 0.0097, 0.028, 0.027,
426     // 0.082, 0.062, 0.065, 0.048,
427     // 0.94, 0.36, 0.94, 0.32,
428     // 1., 0.062, 0.97, 0.97,
429     // 1.5, 1.5, 1.5, 1.5 }; // the last line is PixelmatchVeto and un-used
430    
431     // //initilize the bool value
432     // PassElecVetoInt=0;
433 mingyang 1.1
434 fabstoec 1.17 // dRTrack = PhotonTools::ElectronVetoCiC(p, els);
435 mingyang 1.1
436 fabstoec 1.17 // ScEta_MVA=p->SCluster()->Eta();
437 mingyang 1.1
438 fabstoec 1.17 // isbarrel = (fabs(ScEta_MVA)<1.4442);
439 mingyang 1.1
440 fabstoec 1.17 // R9 = p->R9();
441     // //R9 = p->E33()/p->SCluster()->RawEnergy();
442    
443     // // check which category it is ...
444     // _tCat = 1;
445     // if ( !isbarrel ) _tCat = 3;
446     // if ( R9 < 0.94 ) _tCat++;
447 mingyang 1.1
448 fabstoec 1.17 // //Electron Veto
449     // if(dRTrack > cic4_allcuts_temp_sublead[_tCat-1+6*4]){
450     // PassElecVetoInt=1;
451     // }
452 mingyang 1.1
453 fabstoec 1.17 // return PassElecVetoInt;
454 mingyang 1.1
455 fabstoec 1.17 // }
456 mingyang 1.1
457     //--------------------------------------------------------------------------------------------------
458    
459 fabstoec 1.17
460     Double_t MVATools::GetMVAbdtValue(const Photon* p, const Vertex* vtx, const TrackCol* trackCol, const VertexCol* vtxCol, Double_t _tRho, const PFCandidateCol *fPFCands, const ElectronCol* els, Bool_t applyElectronVeto) {
461    
462     // if there's no reader, or the type is kNone, return the default values of -99.
463     if( ( !fReaderBarrel || !fReaderEndcap ) || fMVAType == kNone ) return -99.;
464    
465     // we compute the variable names... make sure no confilcts when adding new variables...
466 mingyang 1.1
467 fabstoec 1.17 // check if it's a Barrel or EE photon
468     bool isBarrel = ( p->SCluster()->AbsEta() < 1.5 );
469 mingyang 1.11
470 fabstoec 1.18 std::map<std::string,unsigned int>* theVars = ( isBarrel ? &mvaVarMapEB : &mvaVarMapEB );
471    
472 fabstoec 1.17 // loop over all the variables in the map... and keep count (to make sure we have filled all variables)
473     unsigned int varCounter = 0;
474 fabstoec 1.18 for( std::map<std::string,unsigned int>::const_iterator iV = theVars->begin(); iV != theVars->end(); ++iV ) {
475    
476 fabstoec 1.17 TString theVarName = TString(iV->first);
477 fabstoec 1.18 float* theVarValue = &(mvaVars[iV->second]); // pointer to the variable...
478 fabstoec 1.17
479     if(
480     !theVarName.CompareTo("HoE")
481     ) {
482     (*theVarValue) = p->HadOverEm();
483     varCounter++;
484     } else if (
485 fabstoec 1.18 !theVarName.CompareTo("covIEtaIEta") || !theVarName.CompareTo("ph.sigietaieta") || !theVarName.CompareTo("sigieie")
486 fabstoec 1.17 ) {
487     (*theVarValue) = p->CoviEtaiEta();
488     varCounter++;
489     } else if (
490 fabstoec 1.18 !theVarName.CompareTo("R9") || !theVarName.CompareTo("ph.r9") || !theVarName.CompareTo("r9")
491 fabstoec 1.17 ) {
492     (*theVarValue) = p->R9();
493     varCounter++;
494     } else if (
495 fabstoec 1.18 !theVarName.CompareTo("ScEta") || !theVarName.CompareTo("ph.sceta") || !theVarName.CompareTo("sceta")
496 fabstoec 1.17 ) {
497     (*theVarValue) = p->SCluster()->Eta();
498     varCounter++;
499     } else if (
500     !theVarName.CompareTo("rho")
501     ) {
502     (*theVarValue) = _tRho;
503     varCounter++;
504     } else if (
505     !theVarName.CompareTo("tIso1abs")
506     ) {
507     double _ecalIso3 = p->EcalRecHitIsoDr03();
508     double _hcalIso4 = p->HcalTowerSumEtDr04();
509     double _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?
510     (*theVarValue) = _ecalIso3+_hcalIso4+_trackIso1 - 0.17*_tRho;
511     varCounter++;
512     } else if (
513     !theVarName.CompareTo("tIso2abs")
514     ) {
515     double _ecalIso4 = p->EcalRecHitIsoDr04();
516     double _hcalIso4 = p->HcalTowerSumEtDr04();
517     unsigned int wVtxInd = 0;
518     double _trackIso2 = IsolationTools::CiCTrackIsolation(p,vtx, 0.4, 0.02, 0.0, 0.0, 0.1, 1.0, trackCol, &wVtxInd,vtxCol, (!applyElectronVeto ? els : NULL) );
519     (*theVarValue) = _ecalIso4+_hcalIso4+_trackIso2 - 0.52*_tRho;
520     varCounter++;
521     } else if (
522     !theVarName.CompareTo("tIso3abs")
523     ) {
524     (*theVarValue) = IsolationTools::CiCTrackIsolation(p,vtx, 0.3, 0.02, 0.0, 0.0, 0.1, 1.0, trackCol, NULL, NULL, (!applyElectronVeto ? els : NULL) );
525     varCounter++;
526     } else if (
527     !theVarName.CompareTo("absIsoEcal")
528     ) {
529     double _ecalIso3 = p->EcalRecHitIsoDr03();
530     (*theVarValue) = (_ecalIso3-0.17*_tRho);
531     varCounter++;
532     } else if (
533     !theVarName.CompareTo("absIsoHcal")
534     ) {
535     double _hcalIso4 = p->HcalTowerSumEtDr04();
536     (*theVarValue) = (_hcalIso4-0.17*_tRho);
537     varCounter++;
538     } else if (
539     !theVarName.CompareTo("NVertexes")
540     ) {
541    
542     (*theVarValue) = vtxCol->GetEntries();
543     varCounter++;
544     } else if (
545 fabstoec 1.18 !theVarName.CompareTo("EtaWidth") || !theVarName.CompareTo("ph.scetawidth") || !theVarName.CompareTo("sigeta")
546 fabstoec 1.17 ) {
547     (*theVarValue) = p->EtaWidth();
548     varCounter++;
549     } else if (
550 fabstoec 1.18 !theVarName.CompareTo("PhiWidth") || !theVarName.CompareTo("ph.scphiwidth") || !theVarName.CompareTo("sigphi")
551 fabstoec 1.17 ) {
552     (*theVarValue) = p->PhiWidth();
553     varCounter++;
554    
555     } else if (
556 fabstoec 1.18 !theVarName.CompareTo("ph.idmva_CoviEtaiPhi") || !theVarName.CompareTo("covieip")
557 fabstoec 1.17 ) {
558     (*theVarValue) = p->SCluster()->Seed()->CoviEtaiPhi();
559     varCounter++;
560     } else if (
561 fabstoec 1.18 !theVarName.CompareTo("ph.idmva_s4ratio") || !theVarName.CompareTo("s4r")
562 fabstoec 1.17 ) {
563     (*theVarValue) = p->S4Ratio();
564     varCounter++;
565     } else if (
566 fabstoec 1.18 !theVarName.CompareTo("ph.idmva_GammaIso") || !theVarName.CompareTo("pfgiso")
567 fabstoec 1.17 ) {
568     (*theVarValue) = IsolationTools::PFGammaIsolation(p,0.3,0,fPFCands);
569     varCounter++;
570     } else if (
571 fabstoec 1.18 !theVarName.CompareTo("ph.idmva_ChargedIso_selvtx") || !theVarName.CompareTo("pfciso")
572 fabstoec 1.17 ) {
573     (*theVarValue) = IsolationTools::PFChargedIsolation(p,vtx,0.3,0,fPFCands);
574     varCounter++;
575     } else if (
576     !theVarName.CompareTo("ph.idmva_ChargedIso_worstvtx")
577     ) {
578     unsigned int wVtxInd = 0;
579     (*theVarValue) = IsolationTools::PFChargedIsolation(p,vtx,0.3,0,fPFCands,&wVtxInd,vtxCol);
580     varCounter++;
581     } else if (
582     !theVarName.CompareTo("ph.idmva_PsEffWidthSigmaRR")
583     ) {
584     (*theVarValue) = p->EffSigmaRR();
585     varCounter++;
586 fabstoec 1.18 } else if (
587     !theVarName.CompareTo("rawe")
588     ) {
589     (*theVarValue) = p->SCluster()->RawEnergy();
590     varCounter++;
591 fabstoec 1.17 } else {
592     // a variable is not know... copmplain!
593     std::cerr<<" ERROR: MVA Evaluation called with unknown variable name >"<<theVarName<<">."<<std::endl;
594     }
595 mingyang 1.11 }
596 mingyang 1.13
597 fabstoec 1.17 // now all the variables should be filled... check!
598     if( varCounter != theVars->size() )
599     std::cerr<<" ERROR: MVA Evaludation called and not all variables are filled."<<std::endl;
600    
601     // we're ready to compute the MVA value
602     TMVA::Reader* reader = NULL;
603     if (isBarrel)
604 mingyang 1.11 reader = fReaderBarrel;
605 fabstoec 1.17 else
606 mingyang 1.11 reader = fReaderEndcap;
607    
608     assert(reader);
609    
610 fabstoec 1.17 return (reader->EvaluateMVA("BDT method"));
611 mingyang 1.11 }
612    
613 fabstoec 1.17 // Float_t MVATools::GetMVAbdtValue_2012_globe(const Photon* p,const Vertex* vtx,const TrackCol* trackCol,const VertexCol* vtxCol,Double_t _tRho, const PFCandidateCol *fPFCands,const ElectronCol* els,Bool_t applyElectronVeto) {
614    
615     // //get the variables used to compute MVA variables
616     // ecalIso3 = p->EcalRecHitIsoDr03();
617     // ecalIso4 = p->EcalRecHitIsoDr04();
618     // hcalIso4 = p->HcalTowerSumEtDr04();
619    
620     // wVtxInd = 0;
621 mingyang 1.11
622 fabstoec 1.17 // 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?
623    
624     // // track iso worst vtx
625     // trackIso2 = IsolationTools::CiCTrackIsolation(p,vtx, 0.4, 0.02, 0.0, 0.0, 0.1, 1.0, trackCol, &wVtxInd,vtxCol, (!applyElectronVeto ? els : NULL) );
626    
627     // combIso1 = ecalIso3+hcalIso4+trackIso1 - 0.17*_tRho;
628     // combIso2 = ecalIso4+hcalIso4+trackIso2 - 0.52*_tRho;
629    
630     // RawEnergy = p->SCluster()->RawEnergy();
631    
632     // ScEta = p->SCluster()->Eta();
633    
634     // //mva varialbes v1 and v2
635     // tIso1 = (combIso1) *50./p->Et();
636     // tIso3 = (trackIso1)*50./p->Et();
637     // tIso2 = (combIso2) *50./(p->MomVtx(vtxCol->At(wVtxInd)->Position()).Pt());
638     // RelIsoEcal=(ecalIso3-0.17*_tRho)/p->Et();
639     // RelIsoHcal=(hcalIso4-0.17*_tRho)/p->Et();
640    
641     // //compute mva variables for v3
642     // HoE = p->HadOverEm();
643     // covIEtaIEta = p->CoviEtaiEta();
644     // tIso1abs = combIso1;
645     // tIso3abs = trackIso1;
646     // tIso2abs = combIso2;
647     // R9 = p->R9();
648    
649     // absIsoEcal=(ecalIso3-0.17*_tRho);
650     // absIsoHcal=(hcalIso4-0.17*_tRho);
651     // RelEMax=p->SCluster()->Seed()->EMax()/RawEnergy;
652     // RelETop=p->SCluster()->Seed()->ETop()/RawEnergy;
653     // RelEBottom=p->SCluster()->Seed()->EBottom()/RawEnergy;
654     // RelELeft=p->SCluster()->Seed()->ELeft()/RawEnergy;
655     // RelERight=p->SCluster()->Seed()->ERight()/RawEnergy;
656     // RelE2x5Max=p->SCluster()->Seed()->E2x5Max()/RawEnergy;
657     // RelE2x5Top=p->SCluster()->Seed()->E2x5Top()/RawEnergy;
658     // RelE2x5Bottom=p->SCluster()->Seed()->E2x5Bottom()/RawEnergy;
659     // RelE2x5Left=p->SCluster()->Seed()->E2x5Left()/RawEnergy;
660     // RelE2x5Right=p->SCluster()->Seed()->E2x5Right()/RawEnergy;
661     // RelE5x5=p->SCluster()->Seed()->E5x5()/RawEnergy;
662    
663     // EtaWidth=p->EtaWidth();
664     // PhiWidth=p->PhiWidth();
665     // CoviEtaiPhi=p->SCluster()->Seed()->CoviEtaiPhi();
666     // CoviPhiiPhi=p->SCluster()->Seed()->CoviPhiiPhi();
667    
668     // RelPreshowerEnergy=p->SCluster()->PreshowerEnergy()/RawEnergy;
669     // NVertexes=vtxCol->GetEntries();
670    
671     // //spectator variables
672     // Pt_MVA=p->Pt();
673     // ScEta_MVA=p->SCluster()->Eta();
674    
675     // //
676    
677     // isbarrel = (fabs(ScEta_MVA)<1.4442);
678    
679     // //variable 1201
680     // myphoton_pfchargedisogood03=IsolationTools::PFChargedIsolation(p,vtx,0.3,0,fPFCands);
681     // myphoton_pfchargedisobad03=IsolationTools::PFChargedIsolation(p,vtx,0.3,0,fPFCands,&wVtxInd,vtxCol);
682     // myphoton_pfphotoniso03=IsolationTools::PFGammaIsolation(p,0.3,0,fPFCands);
683     // myphoton_sieie=covIEtaIEta;
684     // myphoton_sieip=CoviEtaiPhi;
685     // myphoton_etawidth=EtaWidth;
686     // myphoton_phiwidth=PhiWidth;
687     // myphoton_r9=R9;
688     // myphoton_s4ratio=p->S4Ratio();
689     // myphoton_SCeta=ScEta_MVA;
690     // event_rho= _tRho;
691    
692     // myphoton_ESEffSigmaRR=-99;
693    
694     // if(!isbarrel){
695     // myphoton_ESEffSigmaRR=p->EffSigmaRR();
696     // }
697    
698     // if (isbarrel) {
699     // reader = fReaderBarrel;
700     // }
701     // else {
702     // reader = fReaderEndcap;
703     // }
704    
705     // assert(reader);
706    
707     // double bdt = reader->EvaluateMVA("BDT method");
708    
709     // /* printf("HoE: %f\n",HoE);
710     // printf("covIEtaIEta: %f\n",covIEtaIEta);
711     // printf("tIso1abs: %f\n",tIso1abs);
712     // printf("tIso3abs: %f\n",tIso3abs);
713     // printf("tIso2abs: %f\n",tIso2abs);
714    
715     // printf("absIsoEcal: %f\n",absIsoEcal);
716     // printf("absIsoHcal: %f\n",absIsoHcal);
717     // printf("RelEMax: %f\n",RelEMax);
718     // printf("RelETop: %f\n",RelETop);
719     // printf("RelEBottom: %f\n",RelEBottom);
720     // printf("RelELeft: %f\n",RelELeft);
721     // printf("RelERight: %f\n",RelERight);
722     // printf("RelE2x5Max: %f\n",RelE2x5Max);
723     // printf("RelE2x5Top: %f\n",RelE2x5Top);
724     // printf("RelE2x5Bottom: %f\n",RelE2x5Bottom);
725     // printf("RelE2x5Left: %f\n",RelE2x5Left);
726     // printf("RelE2x5Right;: %f\n",RelE2x5Right);
727     // printf("RelE5x5: %f\n",RelE5x5);
728    
729     // printf("EtaWidth: %f\n",EtaWidth);
730     // printf("PhiWidth: %f\n",PhiWidth);
731     // printf("CoviEtaiPhi: %f\n",CoviEtaiPhi);
732     // printf("CoviPhiiPhi: %f\n",CoviPhiiPhi);
733    
734     // if (!isbarrel) {
735     // printf("RelPreshowerEnergy: %f\n",RelPreshowerEnergy);
736     // }*/
737    
738     // return bdt;
739     // }
740    
741     // Float_t MVATools::GetMVAbdtValue_2011(const Photon* p,const Vertex* vtx,const TrackCol* trackCol,const VertexCol* vtxCol,Double_t _tRho,const ElectronCol* els,Bool_t applyElectronVeto) {
742    
743     // //get the variables used to compute MVA variables
744     // ecalIso3 = p->EcalRecHitIsoDr03();
745     // ecalIso4 = p->EcalRecHitIsoDr04();
746     // hcalIso4 = p->HcalTowerSumEtDr04();
747    
748     // wVtxInd = 0;
749    
750     // 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?
751    
752     // // track iso worst vtx
753     // trackIso2 = IsolationTools::CiCTrackIsolation(p,vtx, 0.4, 0.02, 0.0, 0.0, 0.1, 1.0, trackCol, &wVtxInd,vtxCol, (!applyElectronVeto ? els : NULL) );
754    
755     // combIso1 = ecalIso3+hcalIso4+trackIso1 - 0.17*_tRho;
756     // combIso2 = ecalIso4+hcalIso4+trackIso2 - 0.52*_tRho;
757 mingyang 1.1
758 fabstoec 1.17 // RawEnergy = p->SCluster()->RawEnergy();
759 mingyang 1.1
760 fabstoec 1.17 // ScEta = p->SCluster()->Eta();
761    
762     // //mva varialbes v1 and v2
763     // tIso1 = (combIso1) *50./p->Et();
764     // tIso3 = (trackIso1)*50./p->Et();
765     // tIso2 = (combIso2) *50./(p->MomVtx(vtxCol->At(wVtxInd)->Position()).Pt());
766     // RelIsoEcal=(ecalIso3-0.17*_tRho)/p->Et();
767     // RelIsoHcal=(hcalIso4-0.17*_tRho)/p->Et();
768    
769     // //compute mva variables for v3
770     // HoE = p->HadOverEm();
771     // covIEtaIEta = p->CoviEtaiEta();
772     // tIso1abs = combIso1;
773     // tIso3abs = trackIso1;
774     // tIso2abs = combIso2;
775     // R9 = p->R9();
776    
777     // absIsoEcal=(ecalIso3-0.17*_tRho);
778     // absIsoHcal=(hcalIso4-0.17*_tRho);
779     // RelEMax=p->SCluster()->Seed()->EMax()/RawEnergy;
780     // RelETop=p->SCluster()->Seed()->ETop()/RawEnergy;
781     // RelEBottom=p->SCluster()->Seed()->EBottom()/RawEnergy;
782     // RelELeft=p->SCluster()->Seed()->ELeft()/RawEnergy;
783     // RelERight=p->SCluster()->Seed()->ERight()/RawEnergy;
784     // RelE2x5Max=p->SCluster()->Seed()->E2x5Max()/RawEnergy;
785     // RelE2x5Top=p->SCluster()->Seed()->E2x5Top()/RawEnergy;
786     // RelE2x5Bottom=p->SCluster()->Seed()->E2x5Bottom()/RawEnergy;
787     // RelE2x5Left=p->SCluster()->Seed()->E2x5Left()/RawEnergy;
788     // RelE2x5Right=p->SCluster()->Seed()->E2x5Right()/RawEnergy;
789     // RelE5x5=p->SCluster()->Seed()->E5x5()/RawEnergy;
790    
791     // EtaWidth=p->EtaWidth();
792     // PhiWidth=p->PhiWidth();
793     // CoviEtaiPhi=p->SCluster()->Seed()->CoviEtaiPhi();
794     // CoviPhiiPhi=p->SCluster()->Seed()->CoviPhiiPhi();
795    
796     // RelPreshowerEnergy=p->SCluster()->PreshowerEnergy()/RawEnergy;
797     // NVertexes=vtxCol->GetEntries();
798    
799     // //spectator variables
800     // Pt_MVA=p->Pt();
801     // ScEta_MVA=p->SCluster()->Eta();
802    
803     // //
804    
805     // isbarrel = (fabs(ScEta_MVA)<1.4442);
806    
807     // if (isbarrel) {
808     // reader = fReaderBarrel;
809     // }
810     // else {
811     // reader = fReaderEndcap;
812     // }
813    
814     // assert(reader);
815    
816     // double bdt = reader->EvaluateMVA("BDT method");
817    
818     // /* printf("HoE: %f\n",HoE);
819     // printf("covIEtaIEta: %f\n",covIEtaIEta);
820     // printf("tIso1abs: %f\n",tIso1abs);
821     // printf("tIso3abs: %f\n",tIso3abs);
822     // printf("tIso2abs: %f\n",tIso2abs);
823    
824     // printf("absIsoEcal: %f\n",absIsoEcal);
825     // printf("absIsoHcal: %f\n",absIsoHcal);
826     // printf("RelEMax: %f\n",RelEMax);
827     // printf("RelETop: %f\n",RelETop);
828     // printf("RelEBottom: %f\n",RelEBottom);
829     // printf("RelELeft: %f\n",RelELeft);
830     // printf("RelERight: %f\n",RelERight);
831     // printf("RelE2x5Max: %f\n",RelE2x5Max);
832     // printf("RelE2x5Top: %f\n",RelE2x5Top);
833     // printf("RelE2x5Bottom: %f\n",RelE2x5Bottom);
834     // printf("RelE2x5Left: %f\n",RelE2x5Left);
835     // printf("RelE2x5Right;: %f\n",RelE2x5Right);
836     // printf("RelE5x5: %f\n",RelE5x5);
837    
838     // printf("EtaWidth: %f\n",EtaWidth);
839     // printf("PhiWidth: %f\n",PhiWidth);
840     // printf("CoviEtaiPhi: %f\n",CoviEtaiPhi);
841     // printf("CoviPhiiPhi: %f\n",CoviPhiiPhi);
842    
843     // if (!isbarrel) {
844     // printf("RelPreshowerEnergy: %f\n",RelPreshowerEnergy);
845     // }*/
846 mingyang 1.1
847 fabstoec 1.17 // return bdt;
848     // }