CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_6_1_2_SLHC2/src/RecoTauTag/RecoTau/src/AntiElectronIDMVA.cc

Go to the documentation of this file.
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   //TMVA::Tools::Instance();
00042 
00043   TMVA::Reader *readerX0BL = new TMVA::Reader( "!Color:Silent:Error" );
00044   readerX0BL->AddVariable("HoP",       &TauLeadPFChargedHadrHoP_);
00045   readerX0BL->AddVariable("EoP",       &TauLeadPFChargedHadrEoP_);
00046   //readerX0BL->AddVariable("emFraction",&TauEmFraction_);
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   //readerX0EC->AddVariable("emFraction",&TauEmFraction_);
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   //GammadEta_ = dEta;
00155   //GammadPhi_ = dPhi;
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   //GammadEta_ = dEta;
00307   //GammadPhi_ = dPhi;
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 }