ElectronTools/src/EGammaMvaEleEstimatorCSA14.cc

0001 #include <TFile.h>
0002 #include "EgammaAnalysis/ElectronTools/interface/EGammaMvaEleEstimatorCSA14.h"
0003 #include <cmath>
0004 #include <vector>
0005 #include <cstdio>
0006 #include <zlib.h>
0007 #include "TMVA/MethodBase.h"
0008 #include "FWCore/Utilities/interface/isFinite.h"
0009
0010 //--------------------------------------------------------------------------------------------------
0011 EGammaMvaEleEstimatorCSA14::EGammaMvaEleEstimatorCSA14()
0012     : fMethodname("BDTG method"), fisInitialized(kFALSE), fMVAType(kTrig), fUseBinnedVersion(kTRUE), fNMVABins(0) {
0013   // Constructor.
0014 }
0015
0016 //--------------------------------------------------------------------------------------------------
0017 EGammaMvaEleEstimatorCSA14::~EGammaMvaEleEstimatorCSA14() {
0018   for (unsigned int i = 0; i < fTMVAReader.size(); ++i) {
0019     if (fTMVAReader[i])
0020       delete fTMVAReader[i];
0021   }
0022 }
0023
0024 //--------------------------------------------------------------------------------------------------
0025 void EGammaMvaEleEstimatorCSA14::initialize(std::string methodName,
0026                                             std::string weightsfile,
0027                                             EGammaMvaEleEstimatorCSA14::MVAType type) {
0028   std::vector<std::string> tempWeightFileVector;
0029   tempWeightFileVector.push_back(weightsfile);
0030   initialize(methodName, type, kFALSE, tempWeightFileVector);
0031 }
0032
0033 //--------------------------------------------------------------------------------------------------
0034 void EGammaMvaEleEstimatorCSA14::initialize(std::string methodName,
0035                                             EGammaMvaEleEstimatorCSA14::MVAType type,
0036                                             Bool_t useBinnedVersion,
0037                                             std::vector<std::string> weightsfiles) {
0038   //clean up first
0039   for (unsigned int i = 0; i < fTMVAReader.size(); ++i) {
0040     if (fTMVAReader[i])
0041       delete fTMVAReader[i];
0042   }
0043   fTMVAReader.clear();
0044   fTMVAMethod.clear();
0045   //initialize
0046   fisInitialized = kTRUE;
0047   fMVAType = type;
0048   fMethodname = methodName;
0049   fUseBinnedVersion = useBinnedVersion;
0050
0051   //Define expected number of bins
0052   UInt_t ExpectedNBins = 0;
0053   if (type == kTrig) {
0054     ExpectedNBins = 2;
0055   } else if (type == kNonTrig) {
0056     ExpectedNBins = 4;
0057   } else if (type == kNonTrigPhys14) {
0058     ExpectedNBins = 6;
0059   }
0060
0061   fNMVABins = ExpectedNBins;
0062
0063   //Check number of weight files given
0064   if (fNMVABins != weightsfiles.size()) {
0065     std::cout << "Error: Expected Number of bins = " << fNMVABins
0066               << " does not equal to weightsfiles.size() = " << weightsfiles.size() << std::endl;
0067
0068 #ifndef STANDALONE
0069     assert(fNMVABins == weightsfiles.size());
0070 #endif
0071   }
0072
0073   //Loop over all bins
0074   for (unsigned int i = 0; i < fNMVABins; ++i) {
0075     TMVA::Reader *tmpTMVAReader = new TMVA::Reader("!Color:!Silent:Error");
0076     tmpTMVAReader->SetVerbose(kTRUE);
0077
0078     if (type == kTrig) {
0079       // Pure tracking variables
0080       // Pure tracking variables
0081       tmpTMVAReader->AddVariable("fBrem", &fMVAVar_fbrem);
0082       tmpTMVAReader->AddVariable("kfchi2", &fMVAVar_kfchi2);
0083       tmpTMVAReader->AddVariable("kfhits", &fMVAVar_kfhits);
0084       tmpTMVAReader->AddVariable("gsfChi2", &fMVAVar_gsfchi2);
0085
0086       // Geometrical matchings
0087       tmpTMVAReader->AddVariable("eledeta", &fMVAVar_deta);
0088       tmpTMVAReader->AddVariable("eledphi", &fMVAVar_dphi);
0089       tmpTMVAReader->AddVariable("detacalo", &fMVAVar_detacalo);
0090
0091       // Pure ECAL -> shower shapes
0092       tmpTMVAReader->AddVariable("noZSsee", &fMVAVar_see);
0093       tmpTMVAReader->AddVariable("noZSspp", &fMVAVar_spp);
0094       tmpTMVAReader->AddVariable("etawidth", &fMVAVar_etawidth);
0095       tmpTMVAReader->AddVariable("phiwidth", &fMVAVar_phiwidth);
0096       tmpTMVAReader->AddVariable("noZSe1x5e5x5", &fMVAVar_OneMinusE1x5E5x5);
0097       tmpTMVAReader->AddVariable("noZSr9", &fMVAVar_R9);
0098
0099       // Energy matching
0100       tmpTMVAReader->AddVariable("HtoE", &fMVAVar_HoE);
0101       tmpTMVAReader->AddVariable("EoP", &fMVAVar_EoP);
0102       tmpTMVAReader->AddVariable("IoEmIoP", &fMVAVar_IoEmIoP);
0103       tmpTMVAReader->AddVariable("EEleoPout", &fMVAVar_eleEoPout);
0104       if (i == 1)
0105         tmpTMVAReader->AddVariable("PreShowerOverRaw", &fMVAVar_PreShowerOverRaw);
0106
0107       tmpTMVAReader->AddSpectator("pt", &fMVAVar_pt);
0108       tmpTMVAReader->AddSpectator("absEta", &fMVAVar_abseta);
0109     }
0110
0111     if ((type == kNonTrig) || (type == kNonTrigPhys14)) {
0112       tmpTMVAReader->AddVariable("ele_kfhits", &fMVAVar_kfhits);
0113       // Pure ECAL -> shower shapes
0114       tmpTMVAReader->AddVariable("ele_oldsigmaietaieta", &fMVAVar_see);
0115       tmpTMVAReader->AddVariable("ele_oldsigmaiphiiphi", &fMVAVar_spp);
0116       tmpTMVAReader->AddVariable("ele_oldcircularity", &fMVAVar_OneMinusE1x5E5x5);
0117       tmpTMVAReader->AddVariable("ele_oldr9", &fMVAVar_R9);
0118       tmpTMVAReader->AddVariable("ele_scletawidth", &fMVAVar_etawidth);
0119       tmpTMVAReader->AddVariable("ele_sclphiwidth", &fMVAVar_phiwidth);
0120       tmpTMVAReader->AddVariable("ele_he", &fMVAVar_HoE);
0121       if ((type == kNonTrig) && (i == 1 || i == 3))
0122         tmpTMVAReader->AddVariable("ele_psEoverEraw", &fMVAVar_PreShowerOverRaw);
0123       if ((type == kNonTrigPhys14) && (i == 2 || i == 5))
0124         tmpTMVAReader->AddVariable("ele_psEoverEraw", &fMVAVar_PreShowerOverRaw);
0125
0126       //Pure tracking variables
0127       tmpTMVAReader->AddVariable("ele_kfchi2", &fMVAVar_kfchi2);
0128       tmpTMVAReader->AddVariable("ele_chi2_hits", &fMVAVar_gsfchi2);
0129       // Energy matching
0130       tmpTMVAReader->AddVariable("ele_fbrem", &fMVAVar_fbrem);
0131       tmpTMVAReader->AddVariable("ele_ep", &fMVAVar_EoP);
0132       tmpTMVAReader->AddVariable("ele_eelepout", &fMVAVar_eleEoPout);
0133       tmpTMVAReader->AddVariable("ele_IoEmIop", &fMVAVar_IoEmIoP);
0134
0135       // Geometrical matchings
0136       tmpTMVAReader->AddVariable("ele_deltaetain", &fMVAVar_deta);
0137       tmpTMVAReader->AddVariable("ele_deltaphiin", &fMVAVar_dphi);
0138       tmpTMVAReader->AddVariable("ele_deltaetaseed", &fMVAVar_detacalo);
0139
0140       tmpTMVAReader->AddSpectator("ele_pT", &fMVAVar_pt);
0141       tmpTMVAReader->AddSpectator("ele_isbarrel", &fMVAVar_isBarrel);
0142       tmpTMVAReader->AddSpectator("ele_isendcap", &fMVAVar_isEndcap);
0143       if (type == kNonTrigPhys14)
0144         tmpTMVAReader->AddSpectator("scl_eta", &fMVAVar_SCeta);
0145     }
0146
0147 #ifndef STANDALONE
0148     if ((fMethodname.find("BDT") == 0) &&
0149         (weightsfiles[i].rfind(".xml.gz") == weightsfiles[i].length() - strlen(".xml.gz"))) {
0150       gzFile file = gzopen(weightsfiles[i].c_str(), "rb");
0151       if (file == nullptr) {
0152         std::cout << "Error opening gzip file associated to " << weightsfiles[i] << std::endl;
0153         throw cms::Exception("Configuration", "Error reading zipped XML file");
0154       }
0155       std::vector<char> data;
0156       data.reserve(1024 * 1024 * 10);
0157       unsigned int bufflen = 32 * 1024;
0158       char *buff = reinterpret_cast<char *>(malloc(bufflen));
0159       if (buff == nullptr) {
0160         std::cout << "Error creating buffer for " << weightsfiles[i] << std::endl;
0161         gzclose(file);
0162         throw cms::Exception("Configuration", "Error reading zipped XML file");
0163       }
0164       int read;
0165       while ((read = gzread(file, buff, bufflen)) != 0) {
0166         if (read == -1) {
0167           std::cout << "Error reading gzip file associated to " << weightsfiles[i] << ": " << gzerror(file, &read)
0168                     << std::endl;
0169           gzclose(file);
0170           free(buff);
0171           throw cms::Exception("Configuration", "Error reading zipped XML file");
0172         }
0173         data.insert(data.end(), buff, buff + read);
0174       }
0175       if (gzclose(file) != Z_OK) {
0176         std::cout << "Error closing gzip file associated to " << weightsfiles[i] << std::endl;
0177       }
0178       free(buff);
0179       data.push_back('\0');  // IMPORTANT
0180       fTMVAMethod.push_back(dynamic_cast<TMVA::MethodBase *>(tmpTMVAReader->BookMVA(TMVA::Types::kBDT, &data[0])));
0181     } else {
0182       if (weightsfiles[i].rfind(".xml.gz") == weightsfiles[i].length() - strlen(".xml.gz")) {
0183         std::cout << "Error: xml.gz unsupported for method " << fMethodname << ", weight file " << weightsfiles[i]
0184                   << std::endl;
0185         throw cms::Exception("Configuration", "Error reading zipped XML file");
0186       }
0187       fTMVAMethod.push_back(dynamic_cast<TMVA::MethodBase *>(tmpTMVAReader->BookMVA(fMethodname, weightsfiles[i])));
0188     }
0189 #else
0190     if (weightsfiles[i].rfind(".xml.gz") == weightsfiles[i].length() - strlen(".xml.gz")) {
0191       std::cout << "Error: xml.gz unsupported for method " << fMethodname << ", weight file " << weightsfiles[i]
0192                 << std::endl;
0193       abort();
0194     }
0195     fTMVAMethod.push_back(dynamic_cast<TMVA::MethodBase *>(tmpTMVAReader->BookMVA(fMethodname, weightsfiles[i])));
0196 #endif
0197     std::cout << "MVABin " << i << " : MethodName = " << fMethodname << " , type == " << type << " , "
0198               << "Load weights file : " << weightsfiles[i] << std::endl;
0199     fTMVAReader.push_back(tmpTMVAReader);
0200   }
0201   std::cout << "Electron ID MVA Completed\n";
0202 }
0203
0204 //--------------------------------------------------------------------------------------------------
0205 UInt_t EGammaMvaEleEstimatorCSA14::GetMVABin(double eta, double pt) const {
0206   //Default is to return the first bin
0207   unsigned int bin = 0;
0208
0209   if (fMVAType == EGammaMvaEleEstimatorCSA14::kNonTrig) {
0210     bin = 0;
0211     if (pt < 10 && fabs(eta) < 1.479)
0212       bin = 0;
0213     if (pt < 10 && fabs(eta) >= 1.479)
0214       bin = 1;
0215     if (pt >= 10 && fabs(eta) < 1.479)
0216       bin = 2;
0217     if (pt >= 10 && fabs(eta) >= 1.479)
0218       bin = 3;
0219   }
0220
0221   if (fMVAType == EGammaMvaEleEstimatorCSA14::kTrig) {
0222     bin = 0;
0223     if (pt >= 10 && fabs(eta) < 1.479)
0224       bin = 0;
0225     if (pt >= 10 && fabs(eta) >= 1.479)
0226       bin = 1;
0227   }
0228
0229   if (fMVAType == EGammaMvaEleEstimatorCSA14::kNonTrigPhys14) {
0230     bin = 0;
0231     if (pt < 10 && fabs(eta) < 0.8)
0232       bin = 0;
0233     if (pt < 10 && fabs(eta) >= 0.8 && fabs(eta) < 1.479)
0234       bin = 1;
0235     if (pt < 10 && fabs(eta) >= 1.479)
0236       bin = 2;
0237     if (pt >= 10 && fabs(eta) < 0.8)
0238       bin = 3;
0239     if (pt >= 10 && fabs(eta) >= 0.8 && fabs(eta) < 1.479)
0240       bin = 4;
0241     if (pt >= 10 && fabs(eta) >= 1.479)
0242       bin = 5;
0243   }
0244
0245   return bin;
0246 }
0247
0248 //--------------------------------------------------------------------------------------------------
0249
0250 // for kTrig and kNonTrig algorithm
0251 Double_t EGammaMvaEleEstimatorCSA14::mvaValue(const reco::GsfElectron &ele,
0252                                               const reco::Vertex &vertex,
0253                                               const TransientTrackBuilder &transientTrackBuilder,
0254                                               noZS::EcalClusterLazyTools myEcalCluster,
0255                                               bool printDebug) {
0256   if (!fisInitialized) {
0257     std::cout << "Error: EGammaMvaEleEstimatorCSA14 not properly initialized.\n";
0258     return -9999;
0259   }
0260
0261   if ((fMVAType != EGammaMvaEleEstimatorCSA14::kTrig) && (fMVAType != EGammaMvaEleEstimatorCSA14::kNonTrig) &&
0262       (fMVAType != EGammaMvaEleEstimatorCSA14::kNonTrigPhys14)) {
0263     std::cout << "Error: This method should be called for kTrig or kNonTrig or kNonTrigPhys14 MVA only" << std::endl;
0264     return -9999;
0265   }
0266
0267   bool validKF = false;
0268   reco::TrackRef myTrackRef = ele.closestCtfTrackRef();
0269   validKF = (myTrackRef.isAvailable());
0270   validKF &= (myTrackRef.isNonnull());
0271
0272   // Pure tracking variables
0273   fMVAVar_fbrem = ele.fbrem();
0274   fMVAVar_kfchi2 = (validKF) ? myTrackRef->normalizedChi2() : 0;
0275   fMVAVar_kfhits = (validKF) ? myTrackRef->hitPattern().trackerLayersWithMeasurement() : -1.;
0276   fMVAVar_kfhitsall =
0277       (validKF) ? myTrackRef->numberOfValidHits() : -1.;  //  save also this in your ntuple as possible alternative
0278   fMVAVar_gsfchi2 = ele.gsfTrack()->normalizedChi2();
0279
0280   // Geometrical matchings
0281   fMVAVar_deta = ele.deltaEtaSuperClusterTrackAtVtx();
0282   fMVAVar_dphi = ele.deltaPhiSuperClusterTrackAtVtx();
0283   fMVAVar_detacalo = ele.deltaEtaSeedClusterTrackAtCalo();
0284
0285   // Pure ECAL -> shower shapes
0286   const auto &vCov = myEcalCluster.localCovariances(*(ele.superCluster()->seed()));
0287   if (edm::isFinite(vCov[0]))
0288     fMVAVar_see = sqrt(vCov[0]);  //EleSigmaIEtaIEta
0289   else
0290     fMVAVar_see = 0.;
0291   if (edm::isFinite(vCov[2]))
0292     fMVAVar_spp = sqrt(vCov[2]);  //EleSigmaIPhiIPhi
0293   else
0294     fMVAVar_spp = 0.;
0295
0296   fMVAVar_etawidth = ele.superCluster()->etaWidth();
0297   fMVAVar_phiwidth = ele.superCluster()->phiWidth();
0298   fMVAVar_OneMinusE1x5E5x5 =
0299       (ele.e5x5()) != 0.
0300           ? 1. - (myEcalCluster.e1x5(*(ele.superCluster()->seed())) / myEcalCluster.e5x5(*(ele.superCluster()->seed())))
0301           : -1.;
0302   fMVAVar_R9 = myEcalCluster.e3x3(*(ele.superCluster()->seed())) / ele.superCluster()->rawEnergy();
0303
0304   // Energy matching
0305   fMVAVar_HoE = ele.hadronicOverEm();
0306   fMVAVar_EoP = ele.eSuperClusterOverP();
0307   fMVAVar_IoEmIoP = (1.0 / ele.ecalEnergy()) - (1.0 / ele.p());  // in the future to be changed with ele.gsfTrack()->p()
0308   fMVAVar_eleEoPout = ele.eEleClusterOverPout();
0309   fMVAVar_PreShowerOverRaw = ele.superCluster()->preshowerEnergy() / ele.superCluster()->rawEnergy();
0310
0311   // Spectators
0312   fMVAVar_eta = ele.superCluster()->eta();
0313   fMVAVar_abseta = fabs(ele.superCluster()->eta());
0314   fMVAVar_pt = ele.pt();
0315   fMVAVar_isBarrel = (ele.superCluster()->eta() < 1.479);
0316   fMVAVar_isEndcap = (ele.superCluster()->eta() > 1.479);
0317   fMVAVar_SCeta = ele.superCluster()->eta();
0318
0319   // for triggering electrons get the impact parameteres
0320   if (fMVAType == kTrig) {
0321     //d0
0322     if (ele.gsfTrack().isNonnull()) {
0323       fMVAVar_d0 = (-1.0) * ele.gsfTrack()->dxy(vertex.position());
0324     } else if (ele.closestCtfTrackRef().isNonnull()) {
0325       fMVAVar_d0 = (-1.0) * ele.closestCtfTrackRef()->dxy(vertex.position());
0326     } else {
0327       fMVAVar_d0 = -9999.0;
0328     }
0329
0330     //default values for IP3D
0331     fMVAVar_ip3d = -999.0;
0332     fMVAVar_ip3dSig = 0.0;
0333     if (ele.gsfTrack().isNonnull()) {
0334       const double gsfsign = ((-ele.gsfTrack()->dxy(vertex.position())) >= 0) ? 1. : -1.;
0335
0336       const reco::TransientTrack &tt = transientTrackBuilder.build(ele.gsfTrack());
0337       const std::pair<bool, Measurement1D> &ip3dpv = IPTools::absoluteImpactParameter3D(tt, vertex);
0338       if (ip3dpv.first) {
0339         double ip3d = gsfsign * ip3dpv.second.value();
0340         double ip3derr = ip3dpv.second.error();
0341         fMVAVar_ip3d = ip3d;
0342         fMVAVar_ip3dSig = ip3d / ip3derr;
0343       }
0344     }
0345   }
0346
0347   // evaluate
0348   bindVariables();
0349   Double_t mva = -9999;
0350   if (fUseBinnedVersion) {
0351     int bin = GetMVABin(fMVAVar_eta, fMVAVar_pt);
0352     mva = fTMVAReader[bin]->EvaluateMVA(fTMVAMethod[bin]);
0353   } else {
0354     mva = fTMVAReader[0]->EvaluateMVA(fTMVAMethod[0]);
0355   }
0356
0357   if (printDebug) {
0358     std::cout << " *** Inside the class fMethodname " << fMethodname << " fMVAType " << fMVAType << std::endl;
0359     std::cout << " fbrem " << fMVAVar_fbrem << " kfchi2 " << fMVAVar_kfchi2 << " mykfhits " << fMVAVar_kfhits
0360               << " gsfchi2 " << fMVAVar_gsfchi2 << " deta " << fMVAVar_deta << " dphi " << fMVAVar_dphi << " detacalo "
0361               << fMVAVar_detacalo << " see " << fMVAVar_see << " spp " << fMVAVar_spp << " etawidth "
0362               << fMVAVar_etawidth << " phiwidth " << fMVAVar_phiwidth << " OneMinusE1x5E5x5 "
0363               << fMVAVar_OneMinusE1x5E5x5 << " R9 " << fMVAVar_R9 << " HoE " << fMVAVar_HoE << " EoP " << fMVAVar_EoP
0364               << " IoEmIoP " << fMVAVar_IoEmIoP << " eleEoPout " << fMVAVar_eleEoPout << " d0 " << fMVAVar_d0
0365               << " ip3d " << fMVAVar_ip3d << " eta " << fMVAVar_eta << " pt " << fMVAVar_pt << std::endl;
0366     std::cout << " ### MVA " << mva << std::endl;
0367   }
0368
0369   return mva;
0370 }
0371
0372 Double_t EGammaMvaEleEstimatorCSA14::mvaValue(const pat::Electron &ele, bool printDebug) {
0373   if (!fisInitialized) {
0374     std::cout << "Error: EGammaMvaEleEstimatorCSA14 not properly initialized.\n";
0375     return -9999;
0376   }
0377
0378   if ((fMVAType != EGammaMvaEleEstimatorCSA14::kTrig) && (fMVAType != EGammaMvaEleEstimatorCSA14::kNonTrig) &&
0379       (fMVAType != EGammaMvaEleEstimatorCSA14::kNonTrigPhys14)) {
0380     std::cout << "Error: This method should be called for kTrig or kNonTrig or kNonTrigPhys14 MVA only" << std::endl;
0381     return -9999;
0382   }
0383
0384   bool validKF = false;
0385   reco::TrackRef myTrackRef = ele.closestCtfTrackRef();
0386   validKF = (myTrackRef.isAvailable());
0387   validKF &= (myTrackRef.isNonnull());
0388
0389   // Pure tracking variables
0390   fMVAVar_fbrem = ele.fbrem();
0391   fMVAVar_kfchi2 = (validKF) ? myTrackRef->normalizedChi2() : 0;
0392   fMVAVar_kfhits = (validKF) ? myTrackRef->hitPattern().trackerLayersWithMeasurement() : -1.;
0393   fMVAVar_kfhitsall =
0394       (validKF) ? myTrackRef->numberOfValidHits() : -1.;  //  save also this in your ntuple as possible alternative
0395   fMVAVar_gsfchi2 = ele.gsfTrack()->normalizedChi2();
0396
0397   // Geometrical matchings
0398   fMVAVar_deta = ele.deltaEtaSuperClusterTrackAtVtx();
0399   fMVAVar_dphi = ele.deltaPhiSuperClusterTrackAtVtx();
0400   fMVAVar_detacalo = ele.deltaEtaSeedClusterTrackAtCalo();
0401
0402   // Pure ECAL -> shower shapes
0403   fMVAVar_see = ele.full5x5_sigmaIetaIeta();  //EleSigmaIEtaIEta
0404   fMVAVar_spp = ele.full5x5_sigmaIphiIphi();  //EleSigmaIPhiIPhi
0405
0406   fMVAVar_etawidth = ele.superCluster()->etaWidth();
0407   fMVAVar_phiwidth = ele.superCluster()->phiWidth();
0408   fMVAVar_OneMinusE1x5E5x5 = (ele.full5x5_e5x5()) != 0. ? 1. - (ele.full5x5_e1x5() / ele.full5x5_e5x5()) : -1.;
0409   fMVAVar_R9 = ele.full5x5_r9();
0410
0411   // Energy matching
0412   fMVAVar_HoE = ele.hadronicOverEm();
0413   fMVAVar_EoP = ele.eSuperClusterOverP();
0414   fMVAVar_IoEmIoP = (1.0 / ele.ecalEnergy()) - (1.0 / ele.p());  // in the future to be changed with ele.gsfTrack()->p()
0415   fMVAVar_eleEoPout = ele.eEleClusterOverPout();
0416   fMVAVar_PreShowerOverRaw = ele.superCluster()->preshowerEnergy() / ele.superCluster()->rawEnergy();
0417
0418   // Spectators
0419   fMVAVar_eta = ele.superCluster()->eta();
0420   fMVAVar_abseta = fabs(ele.superCluster()->eta());
0421   fMVAVar_pt = ele.pt();
0422   fMVAVar_isBarrel = (ele.superCluster()->eta() < 1.479);
0423   fMVAVar_isEndcap = (ele.superCluster()->eta() > 1.479);
0424   fMVAVar_SCeta = ele.superCluster()->eta();
0425
0426   // evaluate
0427   bindVariables();
0428   Double_t mva = -9999;
0429   if (fUseBinnedVersion) {
0430     int bin = GetMVABin(fMVAVar_eta, fMVAVar_pt);
0431     mva = fTMVAReader[bin]->EvaluateMVA(fTMVAMethod[bin]);
0432   } else {
0433     mva = fTMVAReader[0]->EvaluateMVA(fTMVAMethod[0]);
0434   }
0435
0436   if (printDebug) {
0437     std::cout << " *** Inside the class fMethodname " << fMethodname << " fMVAType " << fMVAType << std::endl;
0438     std::cout << " fbrem " << fMVAVar_fbrem << " kfchi2 " << fMVAVar_kfchi2 << " mykfhits " << fMVAVar_kfhits
0439               << " gsfchi2 " << fMVAVar_gsfchi2 << " deta " << fMVAVar_deta << " dphi " << fMVAVar_dphi << " detacalo "
0440               << fMVAVar_detacalo << " see " << fMVAVar_see << " spp " << fMVAVar_spp << " etawidth "
0441               << fMVAVar_etawidth << " phiwidth " << fMVAVar_phiwidth << " OneMinusE1x5E5x5 "
0442               << fMVAVar_OneMinusE1x5E5x5 << " R9 " << fMVAVar_R9 << " HoE " << fMVAVar_HoE << " EoP " << fMVAVar_EoP
0443               << " IoEmIoP " << fMVAVar_IoEmIoP << " eleEoPout " << fMVAVar_eleEoPout << " eta " << fMVAVar_eta
0444               << " pt " << fMVAVar_pt << std::endl;
0445     std::cout << " ### MVA " << mva << std::endl;
0446   }
0447
0448   return mva;
0449 }
0450
0451 void EGammaMvaEleEstimatorCSA14::bindVariables() {
0452   // this binding is needed for variables that sometime diverge.
0453
0454   if (fMVAVar_fbrem < -1.)
0455     fMVAVar_fbrem = -1.;
0456
0457   fMVAVar_deta = fabs(fMVAVar_deta);
0458   if (fMVAVar_deta > 0.06)
0459     fMVAVar_deta = 0.06;
0460
0461   fMVAVar_dphi = fabs(fMVAVar_dphi);
0462   if (fMVAVar_dphi > 0.6)
0463     fMVAVar_dphi = 0.6;
0464
0465   if (fMVAVar_EoP > 20.)
0466     fMVAVar_EoP = 20.;
0467
0468   if (fMVAVar_eleEoPout > 20.)
0469     fMVAVar_eleEoPout = 20.;
0470
0471   fMVAVar_detacalo = fabs(fMVAVar_detacalo);
0472   if (fMVAVar_detacalo > 0.2)
0473     fMVAVar_detacalo = 0.2;
0474
0475   if (fMVAVar_OneMinusE1x5E5x5 < -1.)
0476     fMVAVar_OneMinusE1x5E5x5 = -1;
0477
0478   if (fMVAVar_OneMinusE1x5E5x5 > 2.)
0479     fMVAVar_OneMinusE1x5E5x5 = 2.;
0480
0481   if (fMVAVar_R9 > 5)
0482     fMVAVar_R9 = 5;
0483
0484   if (fMVAVar_gsfchi2 > 200.)
0485     fMVAVar_gsfchi2 = 200;
0486
0487   if (fMVAVar_kfchi2 > 10.)
0488     fMVAVar_kfchi2 = 10.;
0489
0490     // Needed for a bug in CMSSW_420, fixed in more recent CMSSW versions
0491 #ifndef STANDALONE
0492   if (edm::isNotFinite(fMVAVar_spp))
0493 #else
0494   if (std::isnan(fMVAVar_spp))
0495 #endif
0496     fMVAVar_spp = 0.;
0497
0498   return;
0499 }