00001 #include <TFile.h>
00002 #include <TMath.h>
00003 #include "RecoTauTag/RecoTau//interface/AntiElectronIDMVA.h"
00004 #include "RecoTauTag/RecoTau/interface/TMVAZipReader.h"
00005
00006 AntiElectronIDMVA::AntiElectronIDMVA():
00007 isInitialized_(kFALSE),
00008 methodName_("BDT")
00009 {
00010 for(UInt_t i=0; i<6; ++i) {
00011 fTMVAReader_[i] = 0;
00012 }
00013
00014 verbosity_ = 1;
00015 }
00016
00017
00018 AntiElectronIDMVA::~AntiElectronIDMVA()
00019 {
00020 for(UInt_t i=0; i<6; ++i) {
00021 if (fTMVAReader_[i]) delete fTMVAReader_[i];
00022 }
00023 }
00024
00025 void AntiElectronIDMVA::Initialize(std::string methodName,
00026 std::string oneProng0Pi0_BL,
00027 std::string oneProng1pi0wGSF_BL,
00028 std::string oneProng1pi0woGSF_BL,
00029 std::string oneProng0Pi0_EC,
00030 std::string oneProng1pi0wGSF_EC,
00031 std::string oneProng1pi0woGSF_EC
00032 ){
00033
00034 for(UInt_t i=0; i<6; ++i) {
00035 if (fTMVAReader_[i]) delete fTMVAReader_[i];
00036 }
00037
00038 isInitialized_ = kTRUE;
00039 methodName_ = methodName;
00040
00041
00042
00043 TMVA::Reader *readerX0BL = new TMVA::Reader( "!Color:Silent:Error" );
00044 readerX0BL->AddVariable("HoP", &TauLeadPFChargedHadrHoP_);
00045 readerX0BL->AddVariable("EoP", &TauLeadPFChargedHadrEoP_);
00046
00047 readerX0BL->SetVerbose(verbosity_);
00048 reco::details::loadTMVAWeights(readerX0BL, methodName_, oneProng0Pi0_BL );
00049
00050 TMVA::Reader *reader11BL = new TMVA::Reader( "!Color:Silent:Error" );
00051 reader11BL->AddVariable("mva", &TauLeadPFChargedHadrMva_);
00052 reader11BL->AddVariable("visMass", &TauVisMass_);
00053 reader11BL->AddVariable("etaMom2*TMath::Sqrt(gammaFrac)*pt", &GammadEta_);
00054 reader11BL->AddVariable("phiMom2*TMath::Sqrt(gammaFrac)*pt", &GammadPhi_);
00055 reader11BL->AddVariable("gammaFrac", &GammadPt_);
00056 reader11BL->SetVerbose(verbosity_);
00057 reco::details::loadTMVAWeights(reader11BL, methodName_, oneProng1pi0wGSF_BL );
00058
00059 TMVA::Reader *reader01BL = new TMVA::Reader( "!Color:Silent:Error" );
00060 reader01BL->AddVariable("visMass", &TauVisMass_);
00061 reader01BL->AddVariable("etaMom2*TMath::Sqrt(gammaFrac)*pt", &GammadEta_);
00062 reader01BL->AddVariable("phiMom2*TMath::Sqrt(gammaFrac)*pt", &GammadPhi_);
00063 reader01BL->AddVariable("gammaFrac", &GammadPt_);
00064 reader01BL->SetVerbose(verbosity_);
00065 reco::details::loadTMVAWeights(reader01BL, methodName_, oneProng1pi0woGSF_BL );
00066
00068
00069 TMVA::Reader *readerX0EC = new TMVA::Reader( "!Color:Silent:Error" );
00070 readerX0EC->AddVariable("HoP", &TauLeadPFChargedHadrHoP_);
00071 readerX0EC->AddVariable("EoP", &TauLeadPFChargedHadrEoP_);
00072
00073 readerX0EC->SetVerbose(verbosity_);
00074 reco::details::loadTMVAWeights(readerX0EC, methodName_, oneProng0Pi0_EC );
00075
00076 TMVA::Reader *reader11EC = new TMVA::Reader( "!Color:Silent:Error" );
00077 reader11EC->AddVariable("mva", &TauLeadPFChargedHadrMva_);
00078 reader11EC->AddVariable("visMass", &TauVisMass_);
00079 reader11EC->AddVariable("etaMom2*TMath::Sqrt(gammaFrac)*pt", &GammadEta_);
00080 reader11EC->AddVariable("phiMom2*TMath::Sqrt(gammaFrac)*pt", &GammadPhi_);
00081 reader11EC->AddVariable("gammaFrac", &GammadPt_);
00082 reader11EC->SetVerbose(verbosity_);
00083 reco::details::loadTMVAWeights(reader11EC, methodName_, oneProng1pi0wGSF_EC );
00084
00085 TMVA::Reader *reader01EC = new TMVA::Reader( "!Color:Silent:Error" );
00086 reader01EC->AddVariable("visMass", &TauVisMass_);
00087 reader01EC->AddVariable("etaMom2*TMath::Sqrt(gammaFrac)*pt", &GammadEta_);
00088 reader01EC->AddVariable("phiMom2*TMath::Sqrt(gammaFrac)*pt", &GammadPhi_);
00089 reader01EC->AddVariable("gammaFrac", &GammadPt_);
00090 reader01EC->SetVerbose(verbosity_);
00091 reco::details::loadTMVAWeights(reader01EC, methodName_, oneProng1pi0woGSF_EC );
00092
00093 fTMVAReader_[0] = readerX0BL;
00094 fTMVAReader_[1] = reader11BL;
00095 fTMVAReader_[2] = reader01BL;
00096 fTMVAReader_[3] = readerX0EC;
00097 fTMVAReader_[4] = reader11EC;
00098 fTMVAReader_[5] = reader01EC;
00099 }
00100
00101
00102 double AntiElectronIDMVA::MVAValue(Float_t TauEta, Float_t TauPt,
00103 Float_t TauSignalPFChargedCands, Float_t TauSignalPFGammaCands,
00104 Float_t TauLeadPFChargedHadrMva,
00105 Float_t TauLeadPFChargedHadrHoP, Float_t TauLeadPFChargedHadrEoP,
00106 Float_t TauHasGsf, Float_t TauVisMass, Float_t TauEmFraction,
00107 std::vector<Float_t>* GammasdEta, std::vector<Float_t>* GammasdPhi, std::vector<Float_t>* GammasPt
00108 ){
00109
00110 if (!isInitialized_) {
00111 std::cout << "Error: AntiElectronMVA with method 1 not properly initialized.\n";
00112 return -999;
00113 }
00114
00115 double mva;
00116
00117 TauVisMass_ = TauVisMass;
00118 TauLeadPFChargedHadrMva_ = TMath::Max(TauLeadPFChargedHadrMva,float(-1.0));
00119 TauLeadPFChargedHadrHoP_ = TauLeadPFChargedHadrHoP;
00120 TauLeadPFChargedHadrEoP_ = TauLeadPFChargedHadrEoP;
00121 TauEmFraction_ = TMath::Max(TauEmFraction,float(0.0));
00122
00123 float sumPt = 0;
00124 float dEta = 0;
00125 float dEta2 = 0;
00126 float dPhi = 0;
00127 float dPhi2 = 0;
00128 float sumPt2 = 0;
00129
00130 for(unsigned int k = 0 ; k < GammasPt->size() ; k++){
00131 float pt_k = (*GammasPt)[k];
00132 float phi_k = (*GammasdPhi)[k];
00133 if ((*GammasdPhi)[k] > TMath::Pi()) phi_k = (*GammasdPhi)[k] - 2*TMath::Pi();
00134 else if((*GammasdPhi)[k] < -TMath::Pi()) phi_k = (*GammasdPhi)[k] + 2*TMath::Pi();
00135 float eta_k = (*GammasdEta)[k];
00136 sumPt += pt_k;
00137 sumPt2 += (pt_k*pt_k);
00138 dEta += (pt_k*eta_k);
00139 dEta2 += (pt_k*eta_k*eta_k);
00140 dPhi += (pt_k*phi_k);
00141 dPhi2 += (pt_k*phi_k*phi_k);
00142 }
00143
00144 GammadPt_ = sumPt/TauPt;
00145
00146 if(sumPt>0){
00147 dEta /= sumPt;
00148 dPhi /= sumPt;
00149 dEta2 /= sumPt;
00150 dPhi2 /= sumPt;
00151
00152 }
00153
00154
00155
00156
00157 GammadEta_ = TMath::Sqrt(dEta2)*TMath::Sqrt(GammadPt_)*TauPt;
00158 GammadPhi_ = TMath::Sqrt(dPhi2)*TMath::Sqrt(GammadPt_)*TauPt;
00159
00160
00161 if( TauSignalPFChargedCands==3 )
00162 mva = 1.0;
00163 else if( TauSignalPFChargedCands==1 && TauSignalPFGammaCands==0){
00164 if(TMath::Abs(TauEta)<1.5)
00165 mva = fTMVAReader_[0]->EvaluateMVA( methodName_ );
00166 else
00167 mva = fTMVAReader_[3]->EvaluateMVA( methodName_ );
00168 }
00169 else if( TauSignalPFChargedCands==1 && TauSignalPFGammaCands>0 && TauHasGsf>0.5){
00170 if(TMath::Abs(TauEta)<1.5)
00171 mva = fTMVAReader_[1]->EvaluateMVA( methodName_ );
00172 else
00173 mva = fTMVAReader_[4]->EvaluateMVA( methodName_ );
00174 }
00175 else if( TauSignalPFChargedCands==1 && TauSignalPFGammaCands>0 && TauHasGsf<0.5){
00176 if(TMath::Abs(TauEta)<1.5)
00177 mva = fTMVAReader_[2]->EvaluateMVA( methodName_ );
00178 else
00179 mva = fTMVAReader_[5]->EvaluateMVA( methodName_ );
00180 }
00181 else{
00182 mva = -99;
00183 }
00184
00185 return mva;
00186
00187 }
00188
00189 double AntiElectronIDMVA::MVAValue(Float_t TauEta, Float_t TauPt,
00190 Float_t TauSignalPFChargedCands, Float_t TauSignalPFGammaCands,
00191 Float_t TauLeadPFChargedHadrMva,
00192 Float_t TauLeadPFChargedHadrHoP, Float_t TauLeadPFChargedHadrEoP,
00193 Float_t TauHasGsf, Float_t TauVisMass, Float_t TauEmFraction,
00194 Float_t GammaEtaMom, Float_t GammaPhiMom, Float_t GammaEnFrac
00195 ){
00196
00197 if (!isInitialized_) {
00198 std::cout << "Error: AntiElectronMVA with method 2 not properly initialized.\n";
00199 return -999;
00200 }
00201
00202
00203 double mva;
00204
00205 TauVisMass_ = TauVisMass;
00206 TauLeadPFChargedHadrMva_ = TMath::Max(TauLeadPFChargedHadrMva,float(-1.0));
00207 TauLeadPFChargedHadrHoP_ = TauLeadPFChargedHadrHoP;
00208 TauLeadPFChargedHadrEoP_ = TauLeadPFChargedHadrEoP;
00209 TauEmFraction_ = TMath::Max(TauEmFraction,float(0.0));
00210 GammadPt_ = GammaEnFrac;
00211 GammadEta_ = GammaEtaMom;
00212 GammadPhi_ = GammaPhiMom;
00213
00214 if( TauSignalPFChargedCands==3 )
00215 mva = 1.0;
00216 else if( TauSignalPFChargedCands==1 && TauSignalPFGammaCands==0){
00217 if(TMath::Abs(TauEta)<1.5)
00218 mva = fTMVAReader_[0]->EvaluateMVA( methodName_ );
00219 else
00220 mva = fTMVAReader_[3]->EvaluateMVA( methodName_ );
00221 }
00222 else if( TauSignalPFChargedCands==1 && TauSignalPFGammaCands>0 && TauHasGsf>0.5){
00223 if(TMath::Abs(TauEta)<1.5)
00224 mva = fTMVAReader_[1]->EvaluateMVA( methodName_ );
00225 else
00226 mva = fTMVAReader_[4]->EvaluateMVA( methodName_ );
00227 }
00228 else if( TauSignalPFChargedCands==1 && TauSignalPFGammaCands>0 && TauHasGsf<0.5){
00229 if(TMath::Abs(TauEta)<1.5)
00230 mva = fTMVAReader_[2]->EvaluateMVA( methodName_ );
00231 else
00232 mva = fTMVAReader_[5]->EvaluateMVA( methodName_ );
00233 }
00234 else{
00235 mva = -99.;
00236 }
00237
00238 return mva;
00239
00240 }
00241
00242
00243
00244 double AntiElectronIDMVA::MVAValue(const reco::PFTauRef& thePFTauRef){
00245
00246 if (!isInitialized_) {
00247 std::cout << "Error: AntiElectronMVA with method 3 not properly initialized.\n";
00248 return -999;
00249 }
00250
00251 double mva;
00252
00253 TauVisMass_ = (*thePFTauRef).mass();
00254 TauLeadPFChargedHadrMva_ = TMath::Max((*thePFTauRef).electronPreIDOutput(),float(-1.0));
00255 TauLeadPFChargedHadrHoP_ = ((*thePFTauRef).leadPFChargedHadrCand())->hcalEnergy()/(*thePFTauRef).leadPFChargedHadrCand()->p();
00256 TauLeadPFChargedHadrEoP_ = ((*thePFTauRef).leadPFChargedHadrCand())->ecalEnergy()/(*thePFTauRef).leadPFChargedHadrCand()->p();
00257 TauEmFraction_ = TMath::Max((*thePFTauRef).emFraction(),float(0.0));
00258
00259 std::vector<float> GammasdEta;
00260 std::vector<float> GammasdPhi;
00261 std::vector<float> GammasPt;
00262
00263 for(unsigned int k = 0 ; k < ((*thePFTauRef).signalPFGammaCands()).size() ; k++){
00264 reco::PFCandidateRef gamma = ((*thePFTauRef).signalPFGammaCands()).at(k);
00265 if( ((*thePFTauRef).leadPFChargedHadrCand()).isNonnull() ){
00266 GammasdEta.push_back( gamma->eta() - (*thePFTauRef).leadPFChargedHadrCand()->eta() );
00267 GammasdPhi.push_back( gamma->phi() - (*thePFTauRef).leadPFChargedHadrCand()->phi() );
00268 }
00269 else{
00270 GammasdEta.push_back( gamma->eta() - (*thePFTauRef).eta() );
00271 GammasdPhi.push_back( gamma->phi() - (*thePFTauRef).phi() );
00272 }
00273 GammasPt.push_back( gamma->pt() );
00274 }
00275
00276 float sumPt = 0;
00277 float dEta = 0;
00278 float dEta2 = 0;
00279 float dPhi = 0;
00280 float dPhi2 = 0;
00281 float sumPt2 = 0;
00282
00283 for(unsigned int k = 0 ; k < GammasPt.size() ; k++){
00284 float pt_k = GammasPt[k];
00285 float phi_k = GammasdPhi[k];
00286 if (GammasdPhi[k] > TMath::Pi()) phi_k = GammasdPhi[k] - 2*TMath::Pi();
00287 else if(GammasdPhi[k] < -TMath::Pi()) phi_k = GammasdPhi[k] + 2*TMath::Pi();
00288 float eta_k = GammasdEta[k];
00289 sumPt += pt_k;
00290 sumPt2 += (pt_k*pt_k);
00291 dEta += (pt_k*eta_k);
00292 dEta2 += (pt_k*eta_k*eta_k);
00293 dPhi += (pt_k*phi_k);
00294 dPhi2 += (pt_k*phi_k*phi_k);
00295 }
00296
00297 GammadPt_ = sumPt/(*thePFTauRef).pt();
00298
00299 if(sumPt>0){
00300 dEta /= sumPt;
00301 dPhi /= sumPt;
00302 dEta2 /= sumPt;
00303 dPhi2 /= sumPt;
00304 }
00305
00306
00307
00308
00309 GammadEta_ = TMath::Sqrt(dEta2)*TMath::Sqrt(GammadPt_)*(*thePFTauRef).pt();
00310 GammadPhi_ = TMath::Sqrt(dPhi2)*TMath::Sqrt(GammadPt_)*(*thePFTauRef).pt();
00311
00312 if( ((*thePFTauRef).signalPFChargedHadrCands()).size() == 3)
00313 mva = 1.0;
00314 else if( ((*thePFTauRef).signalPFChargedHadrCands()).size()==1 && ((*thePFTauRef).signalPFGammaCands()).size()==0){
00315 if(TMath::Abs((*thePFTauRef).eta())<1.5)
00316 mva = fTMVAReader_[0]->EvaluateMVA( methodName_ );
00317 else
00318 mva = fTMVAReader_[3]->EvaluateMVA( methodName_ );
00319 }
00320 else if( ((*thePFTauRef).signalPFChargedHadrCands()).size()==1 && ((*thePFTauRef).signalPFGammaCands()).size()>0 && (((*thePFTauRef).leadPFChargedHadrCand())->gsfTrackRef()).isNonnull()){
00321 if(TMath::Abs((*thePFTauRef).eta())<1.5)
00322 mva = fTMVAReader_[1]->EvaluateMVA( methodName_ );
00323 else
00324 mva = fTMVAReader_[4]->EvaluateMVA( methodName_ );
00325 }
00326 else if( ((*thePFTauRef).signalPFChargedHadrCands()).size()==1 && ((*thePFTauRef).signalPFGammaCands()).size()>0 && !(((*thePFTauRef).leadPFChargedHadrCand())->gsfTrackRef()).isNonnull()){
00327 if(TMath::Abs((*thePFTauRef).eta())<1.5)
00328 mva = fTMVAReader_[2]->EvaluateMVA( methodName_ );
00329 else
00330 mva = fTMVAReader_[5]->EvaluateMVA( methodName_ );
00331 }
00332 else{
00333 mva = -99;
00334 }
00335
00336 return mva;
00337
00338 }