Utils/src/MVATools.cc

// $Id: MVATools.cc,v 1.6 2011/12/07 00:45:17 mingyang Exp $

#include "MitPhysics/Utils/interface/PhotonTools.h"
#include "MitPhysics/Utils/interface/MVATools.h"
#include "MitPhysics/Utils/interface/ElectronTools.h"
#include "MitPhysics/Utils/interface/IsolationTools.h"
#include "MitAna/DataTree/interface/StableData.h"
#include <TFile.h>
#include <TRandom3.h>
#include "TMVA/Tools.h"//MVA
#include "TMVA/Reader.h"//MVA

ClassImp(mithep::MVATools)

using namespace mithep;


//--------------------------------------------------------------------------------------------------
MVATools::MVATools():
  fReaderEndcap(0),
  fReaderBarrel(0),
  
  //MVA Variables v4
  HoE(0),
  covIEtaIEta(0),
  tIso1abs(0),
  tIso3abs(0),
  tIso2abs(0),
  R9(0),
  
  absIsoEcal(0),
  absIsoHcal(0),
  RelEMax(0),
  RelETop(0),
  RelEBottom(0),
  RelELeft(0),
  RelERight(0),
  RelE2x5Max(0),
  RelE2x5Top(0),
  RelE2x5Bottom(0),
  RelE2x5Left(0),
  RelE2x5Right(0),
  RelE5x5(0),
  
  EtaWidth(0),
  PhiWidth(0),
  CoviEtaiPhi(0),
  CoviPhiiPhi(0),
  
  NVertexes(0),
  RelPreshowerEnergy(0),

  //MVA v2 and v1
  RelIsoEcal(0),
  RelIsoHcal(0),
  tIso1(0),
  tIso3(0),
  tIso2(0)
{
  // Constructor.
}

//--------------------------------------------------------------------------------------------------
void MVATools::InitializeMVA(int VariableType, TString EndcapWeights,TString BarrelWeights) {
  
  if (fReaderEndcap) delete fReaderEndcap;  
  if (fReaderBarrel) delete fReaderBarrel; 
  
  fReaderEndcap = new TMVA::Reader( "!Color:!Silent:Error" );       
  fReaderBarrel = new TMVA::Reader( "!Color:!Silent:Error" );   
  
  TMVA::Reader *readers[2];
  readers[0]  = fReaderEndcap;
  readers[1]  = fReaderBarrel;  
  
  for (UInt_t i=0; i<2; ++i) {

    if(VariableType==0||VariableType==1||VariableType==2){
      readers[i]->AddVariable( "HoE", &HoE );
      readers[i]->AddVariable( "covIEtaIEta", &covIEtaIEta );
      readers[i]->AddVariable( "tIso1", &tIso1 );
      readers[i]->AddVariable( "tIso3", &tIso3 );
      readers[i]->AddVariable( "tIso2", &tIso2 );
      readers[i]->AddVariable( "R9", &R9 );
    }
    
    if(VariableType==3||VariableType==4){
      readers[i]->AddVariable( "HoE", &HoE );
      readers[i]->AddVariable( "covIEtaIEta", &covIEtaIEta );
      readers[i]->AddVariable( "tIso1abs", &tIso1abs );
      readers[i]->AddVariable( "tIso3abs", &tIso3abs );
      readers[i]->AddVariable( "tIso2abs", &tIso2abs );
      readers[i]->AddVariable( "R9", &R9 );
    }
    
    if(VariableType==1||VariableType==2){
      readers[i]->AddVariable( "RelIsoEcal", &RelIsoEcal );
      readers[i]->AddVariable( "RelIsoHcal", &RelIsoHcal );
      readers[i]->AddVariable( "RelEMax", &RelEMax );
      readers[i]->AddVariable( "RelETop", &RelETop );
      readers[i]->AddVariable( "RelEBottom", &RelEBottom );
      readers[i]->AddVariable( "RelELeft", &RelELeft );
      readers[i]->AddVariable( "RelERight", &RelERight );
      readers[i]->AddVariable( "RelE2x5Max", &RelE2x5Max );
      readers[i]->AddVariable( "RelE2x5Top", &RelE2x5Top );
      readers[i]->AddVariable( "RelE2x5Bottom", &RelE2x5Bottom );
      readers[i]->AddVariable( "RelE2x5Left", &RelE2x5Left );
      readers[i]->AddVariable( "RelE2x5Right", &RelE2x5Right );
      readers[i]->AddVariable( "RelE5x5", &RelE5x5 );
    }
    
    if(VariableType==3||VariableType==4){
      readers[i]->AddVariable( "absIsoEcal", &absIsoEcal );
      readers[i]->AddVariable( "absIsoHcal", &absIsoHcal );
      readers[i]->AddVariable( "RelEMax", &RelEMax );
      readers[i]->AddVariable( "RelETop", &RelEBottom);
      readers[i]->AddVariable( "RelEBottom", &RelEBottom );
      readers[i]->AddVariable( "RelELeft", &RelELeft );
      readers[i]->AddVariable( "RelERight", &RelERight );
      readers[i]->AddVariable( "RelE2x5Max", &RelE2x5Max );
      readers[i]->AddVariable( "RelE2x5Top", &RelE2x5Top );
      readers[i]->AddVariable( "RelE2x5Bottom", &RelE2x5Bottom );
      readers[i]->AddVariable( "RelE2x5Left", &RelE2x5Left );
      readers[i]->AddVariable( "RelE2x5Right", &RelE2x5Right );
      readers[i]->AddVariable( "RelE5x5", &RelE5x5 );
    }
    
    if(VariableType==2||VariableType==3||VariableType==4){
      readers[i]->AddVariable( "EtaWidth", &EtaWidth );
      readers[i]->AddVariable( "PhiWidth", &PhiWidth );
      readers[i]->AddVariable( "CoviEtaiPhi", &CoviEtaiPhi );
      readers[i]->AddVariable( "CoviPhiiPhi", &CoviPhiiPhi );
      if(VariableType==4){
        readers[i]->AddVariable( "NVertexes", &NVertexes );
      }
      if(i==0){
        readers[i]->AddVariable( "RelPreshowerEnergy", &RelPreshowerEnergy );
      }
    }

    if(VariableType==6){
      readers[i]->AddVariable( "HoE", &HoE );
      readers[i]->AddVariable( "covIEtaIEta", &covIEtaIEta );
      readers[i]->AddVariable( "tIso1abs", &tIso1abs );
      readers[i]->AddVariable( "tIso3abs", &tIso3abs );
      readers[i]->AddVariable( "tIso2abs", &tIso2abs );
      readers[i]->AddVariable( "R9", &R9 );
      readers[i]->AddVariable( "absIsoEcal", &absIsoEcal );
      readers[i]->AddVariable( "absIsoHcal", &absIsoHcal );
      readers[i]->AddVariable( "RelE5x5", &RelE5x5 );
      readers[i]->AddVariable( "EtaWidth", &EtaWidth );
      readers[i]->AddVariable( "PhiWidth", &PhiWidth );
      readers[i]->AddVariable( "CoviEtaiPhi", &CoviEtaiPhi );
      readers[i]->AddVariable( "CoviPhiiPhi", &CoviPhiiPhi );
      readers[i]->AddVariable( "NVertexes", &NVertexes );
      if(i==0){
      readers[i]->AddVariable( "RelPreshowerEnergy", &RelPreshowerEnergy );
      }
    }

  }
  
  fReaderEndcap->BookMVA("BDT method",EndcapWeights);
  fReaderBarrel->BookMVA("BDT method",BarrelWeights);
  
  assert(fReaderEndcap);
  assert(fReaderBarrel);
  
}

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) {
  
  //initilize the bool value
  PassMVA=kFALSE;
  
  Float_t photon_bdt =  MVATools::GetMVAbdtValue(p,vtx,trackCol,vtxCol, _tRho, els, applyElectronVeto);
  
  if (isbarrel) {
    if(bdt>bdtCutBarrel){
      PassMVA=kTRUE;    
    }
  }
  else {
    if(bdt>bdtCutEndcap){
      PassMVA=kTRUE;    
    }
  } 
  return PassMVA;
}

//---------------------------------------------------------------------------------
Int_t MVATools::PassElectronVetoInt(const Photon* p, const ElectronCol* els) {
  
  // these values are taken from the H2GGlobe code... (actually from Marco/s mail)
  float cic4_allcuts_temp_sublead[] = { 
    3.8,         2.2,         1.77,        1.29,
    11.7,        3.4,         3.9,         1.84,
    3.5,         2.2,         2.3,         1.45,
    0.0106,      0.0097,      0.028,       0.027,
    0.082,       0.062,       0.065,       0.048,
    0.94,        0.36,        0.94,        0.32,
    1.,          0.062,       0.97,        0.97,
    1.5,         1.5,         1.5,         1.5 };  // the last line is PixelmatchVeto and un-used
  
  //initilize the bool value
  PassElecVetoInt=0;
  
  dRTrack = PhotonTools::ElectronVetoCiC(p, els);
  
  ScEta_MVA=p->SCluster()->Eta();
  
  isbarrel = (fabs(ScEta_MVA)<1.4442);

  //R9 = p->R9();
  R9 = p->E33()/p->SCluster()->RawEnergy();
  
  // check which category it is ...
  _tCat = 1;
  if ( !isbarrel ) _tCat = 3;
  if ( R9 < 0.94 ) _tCat++;
  
  //Electron Veto
  if(dRTrack > cic4_allcuts_temp_sublead[_tCat-1+6*4]){
    PassElecVetoInt=1;
  }
  
  return  PassElecVetoInt;
  
}

//--------------------------------------------------------------------------------------------------

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) {
  
  //get the variables used to compute MVA variables
  ecalIso3 = p->EcalRecHitIsoDr03();
  ecalIso4 = p->EcalRecHitIsoDr04();
  hcalIso4 = p->HcalTowerSumEtDr04();
  
  wVtxInd = 0;
  
  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?
  
  // track iso only
  trackIso3 = IsolationTools::CiCTrackIsolation(p,vtx, 0.3, 0.02, 0.0, 0.0, 0.1, 1.0, trackCol, NULL, NULL, (!applyElectronVeto ? els : NULL) );
  
  // track iso worst vtx
  trackIso2 = IsolationTools::CiCTrackIsolation(p,vtx, 0.4, 0.02, 0.0, 0.0, 0.1, 1.0, trackCol, &wVtxInd,vtxCol, (!applyElectronVeto ? els : NULL) );
  
  combIso1 = ecalIso3+hcalIso4+trackIso1 - 0.17*_tRho;
  combIso2 = ecalIso4+hcalIso4+trackIso2 - 0.52*_tRho;
  
  RawEnergy = p->SCluster()->RawEnergy();
  
  //mva varialbes v1 and v2
  tIso1 = (combIso1) *50./p->Et();
  tIso3 = (trackIso3)*50./p->Et();
  tIso2 = (combIso2) *50./(p->MomVtx(vtxCol->At(wVtxInd)->Position()).Pt());
  RelIsoEcal=(ecalIso3-0.17*_tRho)/p->Et();
  RelIsoHcal=(hcalIso4-0.17*_tRho)/p->Et();

  //compute mva variables for v3
  HoE = p->HadOverEm();
  covIEtaIEta = p->CoviEtaiEta();
  tIso1abs = combIso1;
  tIso3abs = trackIso3;
  tIso2abs = combIso2;
  R9 = p->R9();

  absIsoEcal=(ecalIso3-0.17*_tRho);
  absIsoHcal=(hcalIso4-0.17*_tRho);
  RelEMax=p->SCluster()->Seed()->EMax()/RawEnergy;
  RelETop=p->SCluster()->Seed()->ETop()/RawEnergy;
  RelEBottom=p->SCluster()->Seed()->EBottom()/RawEnergy;
  RelELeft=p->SCluster()->Seed()->ELeft()/RawEnergy;
  RelERight=p->SCluster()->Seed()->ERight()/RawEnergy;
  RelE2x5Max=p->SCluster()->Seed()->E2x5Max()/RawEnergy;
  RelE2x5Top=p->SCluster()->Seed()->E2x5Top()/RawEnergy;
  RelE2x5Bottom=p->SCluster()->Seed()->E2x5Bottom()/RawEnergy;
  RelE2x5Left=p->SCluster()->Seed()->E2x5Left()/RawEnergy;
  RelE2x5Right=p->SCluster()->Seed()->E2x5Right()/RawEnergy;
  RelE5x5=p->SCluster()->Seed()->E5x5()/RawEnergy;
  
  EtaWidth=p->SCluster()->EtaWidth();
  PhiWidth=p->SCluster()->PhiWidth();
  CoviEtaiPhi=p->SCluster()->Seed()->CoviEtaiPhi(); 
  CoviPhiiPhi=p->SCluster()->Seed()->CoviPhiiPhi();

  RelPreshowerEnergy=p->SCluster()->PreshowerEnergy()/RawEnergy;
  NVertexes=vtxCol->GetEntries();
  
  //spectator variables
  Pt_MVA=p->Pt();
  ScEta_MVA=p->SCluster()->Eta();
  
  isbarrel = (fabs(ScEta_MVA)<1.4442);
  
  if (isbarrel) {
    reader = fReaderBarrel;
  }
  else {
    reader = fReaderEndcap;
  }
  
  assert(reader); 

  bdt = reader->EvaluateMVA("BDT method");

  /* printf("HoE: %f\n",HoE);
  printf("covIEtaIEta: %f\n",covIEtaIEta);
  printf("tIso1abs: %f\n",tIso1abs);
  printf("tIso3abs: %f\n",tIso3abs);
  printf("tIso2abs: %f\n",tIso2abs);
  
  printf("absIsoEcal: %f\n",absIsoEcal);
  printf("absIsoHcal: %f\n",absIsoHcal);
  printf("RelEMax: %f\n",RelEMax);
  printf("RelETop: %f\n",RelETop);
  printf("RelEBottom: %f\n",RelEBottom);
  printf("RelELeft: %f\n",RelELeft);
  printf("RelERight: %f\n",RelERight);
  printf("RelE2x5Max: %f\n",RelE2x5Max);
  printf("RelE2x5Top: %f\n",RelE2x5Top);
  printf("RelE2x5Bottom: %f\n",RelE2x5Bottom);
  printf("RelE2x5Left: %f\n",RelE2x5Left);
  printf("RelE2x5Right;: %f\n",RelE2x5Right);
  printf("RelE5x5: %f\n",RelE5x5);
  
  printf("EtaWidth: %f\n",EtaWidth);
  printf("PhiWidth: %f\n",PhiWidth);
  printf("CoviEtaiPhi: %f\n",CoviEtaiPhi);
  printf("CoviPhiiPhi: %f\n",CoviPhiiPhi);
  
  if (!isbarrel) {
    printf("RelPreshowerEnergy: %f\n",RelPreshowerEnergy);
    }*/
  
  return bdt;
}
Revision:	1.7
Committed:	Sun Dec 11 00:03:05 2011 UTC (13 years, 4 months ago) by bendavid
Content type:	text/plain
Branch:	MAIN
Changes since 1.6:	+9 -8 lines
Log Message:	more photon updates for mva plus id analysis
#	User	Rev	Content
1	bendavid	1.7	// $Id: MVATools.cc,v 1.6 2011/12/07 00:45:17 mingyang 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			#include "MitAna/DataTree/interface/StableData.h"
8			#include <TFile.h>
9			#include <TRandom3.h>
10			#include "TMVA/Tools.h"//MVA
11			#include "TMVA/Reader.h"//MVA
12
13			ClassImp(mithep::MVATools)
14
15			using namespace mithep;
16
17
18			//--------------------------------------------------------------------------------------------------
19			MVATools::MVATools():
20			fReaderEndcap(0),
21			fReaderBarrel(0),
22
23	mingyang	1.4	//MVA Variables v4
24	mingyang	1.1	HoE(0),
25			covIEtaIEta(0),
26	mingyang	1.4	tIso1abs(0),
27			tIso3abs(0),
28			tIso2abs(0),
29	mingyang	1.1	R9(0),
30
31	mingyang	1.4	absIsoEcal(0),
32			absIsoHcal(0),
33	mingyang	1.1	RelEMax(0),
34			RelETop(0),
35			RelEBottom(0),
36			RelELeft(0),
37			RelERight(0),
38			RelE2x5Max(0),
39			RelE2x5Top(0),
40			RelE2x5Bottom(0),
41			RelE2x5Left(0),
42			RelE2x5Right(0),
43			RelE5x5(0),
44
45			EtaWidth(0),
46			PhiWidth(0),
47	mingyang	1.4	CoviEtaiPhi(0),
48	mingyang	1.1	CoviPhiiPhi(0),
49
50	mingyang	1.4	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	mingyang	1.1	{
60			// Constructor.
61			}
62
63			//--------------------------------------------------------------------------------------------------
64			void MVATools::InitializeMVA(int VariableType, TString EndcapWeights,TString BarrelWeights) {
65
66			if (fReaderEndcap) delete fReaderEndcap;
67			if (fReaderBarrel) delete fReaderBarrel;
68
69			fReaderEndcap = new TMVA::Reader( "!Color:!Silent:Error" );
70			fReaderBarrel = new TMVA::Reader( "!Color:!Silent:Error" );
71
72			TMVA::Reader *readers[2];
73			readers[0] = fReaderEndcap;
74			readers[1] = fReaderBarrel;
75
76			for (UInt_t i=0; i<2; ++i) {
77	mingyang	1.4
78	mingyang	1.1	if(VariableType==0\|\|VariableType==1\|\|VariableType==2){
79			readers[i]->AddVariable( "HoE", &HoE );
80			readers[i]->AddVariable( "covIEtaIEta", &covIEtaIEta );
81			readers[i]->AddVariable( "tIso1", &tIso1 );
82			readers[i]->AddVariable( "tIso3", &tIso3 );
83			readers[i]->AddVariable( "tIso2", &tIso2 );
84			readers[i]->AddVariable( "R9", &R9 );
85			}
86
87	mingyang	1.4	if(VariableType==3\|\|VariableType==4){
88			readers[i]->AddVariable( "HoE", &HoE );
89			readers[i]->AddVariable( "covIEtaIEta", &covIEtaIEta );
90			readers[i]->AddVariable( "tIso1abs", &tIso1abs );
91			readers[i]->AddVariable( "tIso3abs", &tIso3abs );
92			readers[i]->AddVariable( "tIso2abs", &tIso2abs );
93			readers[i]->AddVariable( "R9", &R9 );
94			}
95
96	mingyang	1.1	if(VariableType==1\|\|VariableType==2){
97	mingyang	1.4	readers[i]->AddVariable( "RelIsoEcal", &RelIsoEcal );
98			readers[i]->AddVariable( "RelIsoHcal", &RelIsoHcal );
99			readers[i]->AddVariable( "RelEMax", &RelEMax );
100			readers[i]->AddVariable( "RelETop", &RelETop );
101			readers[i]->AddVariable( "RelEBottom", &RelEBottom );
102			readers[i]->AddVariable( "RelELeft", &RelELeft );
103			readers[i]->AddVariable( "RelERight", &RelERight );
104			readers[i]->AddVariable( "RelE2x5Max", &RelE2x5Max );
105			readers[i]->AddVariable( "RelE2x5Top", &RelE2x5Top );
106			readers[i]->AddVariable( "RelE2x5Bottom", &RelE2x5Bottom );
107			readers[i]->AddVariable( "RelE2x5Left", &RelE2x5Left );
108			readers[i]->AddVariable( "RelE2x5Right", &RelE2x5Right );
109			readers[i]->AddVariable( "RelE5x5", &RelE5x5 );
110			}
111
112			if(VariableType==3\|\|VariableType==4){
113			readers[i]->AddVariable( "absIsoEcal", &absIsoEcal );
114			readers[i]->AddVariable( "absIsoHcal", &absIsoHcal );
115			readers[i]->AddVariable( "RelEMax", &RelEMax );
116			readers[i]->AddVariable( "RelETop", &RelEBottom);
117			readers[i]->AddVariable( "RelEBottom", &RelEBottom );
118	mingyang	1.1	readers[i]->AddVariable( "RelELeft", &RelELeft );
119	mingyang	1.4	readers[i]->AddVariable( "RelERight", &RelERight );
120	mingyang	1.1	readers[i]->AddVariable( "RelE2x5Max", &RelE2x5Max );
121			readers[i]->AddVariable( "RelE2x5Top", &RelE2x5Top );
122			readers[i]->AddVariable( "RelE2x5Bottom", &RelE2x5Bottom );
123			readers[i]->AddVariable( "RelE2x5Left", &RelE2x5Left );
124			readers[i]->AddVariable( "RelE2x5Right", &RelE2x5Right );
125			readers[i]->AddVariable( "RelE5x5", &RelE5x5 );
126			}
127
128	mingyang	1.4	if(VariableType==2\|\|VariableType==3\|\|VariableType==4){
129	mingyang	1.1	readers[i]->AddVariable( "EtaWidth", &EtaWidth );
130			readers[i]->AddVariable( "PhiWidth", &PhiWidth );
131			readers[i]->AddVariable( "CoviEtaiPhi", &CoviEtaiPhi );
132			readers[i]->AddVariable( "CoviPhiiPhi", &CoviPhiiPhi );
133	mingyang	1.4	if(VariableType==4){
134			readers[i]->AddVariable( "NVertexes", &NVertexes );
135			}
136	mingyang	1.1	if(i==0){
137			readers[i]->AddVariable( "RelPreshowerEnergy", &RelPreshowerEnergy );
138			}
139			}
140	mingyang	1.6
141			if(VariableType==6){
142			readers[i]->AddVariable( "HoE", &HoE );
143			readers[i]->AddVariable( "covIEtaIEta", &covIEtaIEta );
144			readers[i]->AddVariable( "tIso1abs", &tIso1abs );
145			readers[i]->AddVariable( "tIso3abs", &tIso3abs );
146			readers[i]->AddVariable( "tIso2abs", &tIso2abs );
147			readers[i]->AddVariable( "R9", &R9 );
148			readers[i]->AddVariable( "absIsoEcal", &absIsoEcal );
149			readers[i]->AddVariable( "absIsoHcal", &absIsoHcal );
150			readers[i]->AddVariable( "RelE5x5", &RelE5x5 );
151			readers[i]->AddVariable( "EtaWidth", &EtaWidth );
152			readers[i]->AddVariable( "PhiWidth", &PhiWidth );
153			readers[i]->AddVariable( "CoviEtaiPhi", &CoviEtaiPhi );
154			readers[i]->AddVariable( "CoviPhiiPhi", &CoviPhiiPhi );
155			readers[i]->AddVariable( "NVertexes", &NVertexes );
156			if(i==0){
157			readers[i]->AddVariable( "RelPreshowerEnergy", &RelPreshowerEnergy );
158			}
159			}
160
161	mingyang	1.1	}
162
163			fReaderEndcap->BookMVA("BDT method",EndcapWeights);
164			fReaderBarrel->BookMVA("BDT method",BarrelWeights);
165	mingyang	1.4
166	mingyang	1.1	assert(fReaderEndcap);
167			assert(fReaderBarrel);
168
169			}
170
171	bendavid	1.7	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) {
172	mingyang	1.1
173			//initilize the bool value
174			PassMVA=kFALSE;
175
176	bendavid	1.7	Float_t photon_bdt = MVATools::GetMVAbdtValue(p,vtx,trackCol,vtxCol, _tRho, els, applyElectronVeto);
177	mingyang	1.1
178	mingyang	1.4	if (isbarrel) {
179			if(bdt>bdtCutBarrel){
180			PassMVA=kTRUE;
181			}
182	mingyang	1.1	}
183	mingyang	1.4	else {
184			if(bdt>bdtCutEndcap){
185			PassMVA=kTRUE;
186	mingyang	1.1	}
187	mingyang	1.4	}
188	mingyang	1.1	return PassMVA;
189			}
190
191			//---------------------------------------------------------------------------------
192			Int_t MVATools::PassElectronVetoInt(const Photon* p, const ElectronCol* els) {
193
194			// these values are taken from the H2GGlobe code... (actually from Marco/s mail)
195			float cic4_allcuts_temp_sublead[] = {
196			3.8, 2.2, 1.77, 1.29,
197			11.7, 3.4, 3.9, 1.84,
198			3.5, 2.2, 2.3, 1.45,
199			0.0106, 0.0097, 0.028, 0.027,
200			0.082, 0.062, 0.065, 0.048,
201			0.94, 0.36, 0.94, 0.32,
202			1., 0.062, 0.97, 0.97,
203			1.5, 1.5, 1.5, 1.5 }; // the last line is PixelmatchVeto and un-used
204
205			//initilize the bool value
206			PassElecVetoInt=0;
207
208			dRTrack = PhotonTools::ElectronVetoCiC(p, els);
209
210			ScEta_MVA=p->SCluster()->Eta();
211
212			isbarrel = (fabs(ScEta_MVA)<1.4442);
213
214	bendavid	1.7	//R9 = p->R9();
215			R9 = p->E33()/p->SCluster()->RawEnergy();
216	mingyang	1.1
217			// check which category it is ...
218			_tCat = 1;
219			if ( !isbarrel ) _tCat = 3;
220			if ( R9 < 0.94 ) _tCat++;
221
222			//Electron Veto
223			if(dRTrack > cic4_allcuts_temp_sublead[_tCat-1+6*4]){
224			PassElecVetoInt=1;
225			}
226
227			return PassElecVetoInt;
228
229			}
230
231			//--------------------------------------------------------------------------------------------------
232
233	bendavid	1.7	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) {
234	mingyang	1.1
235			//get the variables used to compute MVA variables
236			ecalIso3 = p->EcalRecHitIsoDr03();
237			ecalIso4 = p->EcalRecHitIsoDr04();
238			hcalIso4 = p->HcalTowerSumEtDr04();
239
240			wVtxInd = 0;
241
242	bendavid	1.7	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?
243	mingyang	1.1
244			// track iso only
245	bendavid	1.7	trackIso3 = IsolationTools::CiCTrackIsolation(p,vtx, 0.3, 0.02, 0.0, 0.0, 0.1, 1.0, trackCol, NULL, NULL, (!applyElectronVeto ? els : NULL) );
246	mingyang	1.1
247			// track iso worst vtx
248	bendavid	1.7	trackIso2 = IsolationTools::CiCTrackIsolation(p,vtx, 0.4, 0.02, 0.0, 0.0, 0.1, 1.0, trackCol, &wVtxInd,vtxCol, (!applyElectronVeto ? els : NULL) );
249	mingyang	1.1
250			combIso1 = ecalIso3+hcalIso4+trackIso1 - 0.17*_tRho;
251			combIso2 = ecalIso4+hcalIso4+trackIso2 - 0.52*_tRho;
252
253			RawEnergy = p->SCluster()->RawEnergy();
254
255	mingyang	1.4	//mva varialbes v1 and v2
256	mingyang	1.1	tIso1 = (combIso1) *50./p->Et();
257			tIso3 = (trackIso3)*50./p->Et();
258			tIso2 = (combIso2) *50./(p->MomVtx(vtxCol->At(wVtxInd)->Position()).Pt());
259			RelIsoEcal=(ecalIso3-0.17*_tRho)/p->Et();
260			RelIsoHcal=(hcalIso4-0.17*_tRho)/p->Et();
261	mingyang	1.4
262			//compute mva variables for v3
263			HoE = p->HadOverEm();
264			covIEtaIEta = p->CoviEtaiEta();
265			tIso1abs = combIso1;
266			tIso3abs = trackIso3;
267			tIso2abs = combIso2;
268			R9 = p->R9();
269
270			absIsoEcal=(ecalIso3-0.17*_tRho);
271			absIsoHcal=(hcalIso4-0.17*_tRho);
272	mingyang	1.1	RelEMax=p->SCluster()->Seed()->EMax()/RawEnergy;
273			RelETop=p->SCluster()->Seed()->ETop()/RawEnergy;
274			RelEBottom=p->SCluster()->Seed()->EBottom()/RawEnergy;
275			RelELeft=p->SCluster()->Seed()->ELeft()/RawEnergy;
276			RelERight=p->SCluster()->Seed()->ERight()/RawEnergy;
277			RelE2x5Max=p->SCluster()->Seed()->E2x5Max()/RawEnergy;
278			RelE2x5Top=p->SCluster()->Seed()->E2x5Top()/RawEnergy;
279			RelE2x5Bottom=p->SCluster()->Seed()->E2x5Bottom()/RawEnergy;
280			RelE2x5Left=p->SCluster()->Seed()->E2x5Left()/RawEnergy;
281			RelE2x5Right=p->SCluster()->Seed()->E2x5Right()/RawEnergy;
282			RelE5x5=p->SCluster()->Seed()->E5x5()/RawEnergy;
283
284			EtaWidth=p->SCluster()->EtaWidth();
285			PhiWidth=p->SCluster()->PhiWidth();
286			CoviEtaiPhi=p->SCluster()->Seed()->CoviEtaiPhi();
287			CoviPhiiPhi=p->SCluster()->Seed()->CoviPhiiPhi();
288	mingyang	1.4
289	mingyang	1.1	RelPreshowerEnergy=p->SCluster()->PreshowerEnergy()/RawEnergy;
290	mingyang	1.4	NVertexes=vtxCol->GetEntries();
291	mingyang	1.1
292			//spectator variables
293			Pt_MVA=p->Pt();
294			ScEta_MVA=p->SCluster()->Eta();
295
296			isbarrel = (fabs(ScEta_MVA)<1.4442);
297
298			if (isbarrel) {
299			reader = fReaderBarrel;
300			}
301			else {
302			reader = fReaderEndcap;
303			}
304
305			assert(reader);
306
307			bdt = reader->EvaluateMVA("BDT method");
308
309	mingyang	1.5	/* printf("HoE: %f\n",HoE);
310	mingyang	1.1	printf("covIEtaIEta: %f\n",covIEtaIEta);
311	mingyang	1.4	printf("tIso1abs: %f\n",tIso1abs);
312			printf("tIso3abs: %f\n",tIso3abs);
313			printf("tIso2abs: %f\n",tIso2abs);
314	mingyang	1.1
315	mingyang	1.4	printf("absIsoEcal: %f\n",absIsoEcal);
316			printf("absIsoHcal: %f\n",absIsoHcal);
317	mingyang	1.1	printf("RelEMax: %f\n",RelEMax);
318			printf("RelETop: %f\n",RelETop);
319			printf("RelEBottom: %f\n",RelEBottom);
320			printf("RelELeft: %f\n",RelELeft);
321			printf("RelERight: %f\n",RelERight);
322			printf("RelE2x5Max: %f\n",RelE2x5Max);
323			printf("RelE2x5Top: %f\n",RelE2x5Top);
324			printf("RelE2x5Bottom: %f\n",RelE2x5Bottom);
325			printf("RelE2x5Left: %f\n",RelE2x5Left);
326			printf("RelE2x5Right;: %f\n",RelE2x5Right);
327			printf("RelE5x5: %f\n",RelE5x5);
328
329			printf("EtaWidth: %f\n",EtaWidth);
330			printf("PhiWidth: %f\n",PhiWidth);
331			printf("CoviEtaiPhi: %f\n",CoviEtaiPhi);
332			printf("CoviPhiiPhi: %f\n",CoviPhiiPhi);
333
334	mingyang	1.4	if (!isbarrel) {
335	mingyang	1.1	printf("RelPreshowerEnergy: %f\n",RelPreshowerEnergy);
336	mingyang	1.5	}*/
337	mingyang	1.1
338			return bdt;
339			}