Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-10-25 09:45:03

0001 #include <TFile.h>
0002 #include "EgammaAnalysis/ElectronTools/interface/EGammaMvaEleEstimator.h"
0003 #include <cmath>
0004 #include <vector>
0005 
0006 #ifndef STANDALONE
0007 #include "DataFormats/TrackReco/interface/Track.h"
0008 #include "DataFormats/GsfTrackReco/interface/GsfTrack.h"
0009 #include "DataFormats/TrackReco/interface/TrackFwd.h"
0010 #include "DataFormats/EgammaCandidates/interface/GsfElectron.h"
0011 #include "DataFormats/EgammaReco/interface/SuperCluster.h"
0012 #include "DataFormats/MuonReco/interface/Muon.h"
0013 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidate.h"
0014 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidateFwd.h"
0015 #include "DataFormats/Common/interface/RefToPtr.h"
0016 #include "DataFormats/VertexReco/interface/Vertex.h"
0017 #include "RecoEcal/EgammaCoreTools/interface/EcalClusterLazyTools.h"
0018 #include "TrackingTools/TransientTrack/interface/TransientTrackBuilder.h"
0019 #include "TrackingTools/IPTools/interface/IPTools.h"
0020 #include "EgammaAnalysis/ElectronTools/interface/ElectronEffectiveArea.h"
0021 #include "DataFormats/Common/interface/RefToPtr.h"
0022 #include "FWCore/Utilities/interface/isFinite.h"
0023 #include <cstdio>
0024 #include <zlib.h>
0025 #endif
0026 
0027 //--------------------------------------------------------------------------------------------------
0028 EGammaMvaEleEstimator::EGammaMvaEleEstimator()
0029     : fMethodname("BDTG method"), fisInitialized(kFALSE), fMVAType(kTrig), fUseBinnedVersion(kTRUE), fNMVABins(0) {
0030   // Constructor.
0031 }
0032 
0033 //--------------------------------------------------------------------------------------------------
0034 EGammaMvaEleEstimator::~EGammaMvaEleEstimator() {
0035   for (unsigned int i = 0; i < fTMVAReader.size(); ++i) {
0036     if (fTMVAReader[i])
0037       delete fTMVAReader[i];
0038   }
0039 }
0040 
0041 //--------------------------------------------------------------------------------------------------
0042 void EGammaMvaEleEstimator::initialize(std::string methodName,
0043                                        std::string weightsfile,
0044                                        EGammaMvaEleEstimator::MVAType type) {
0045   std::vector<std::string> tempWeightFileVector;
0046   tempWeightFileVector.push_back(weightsfile);
0047   initialize(methodName, type, kFALSE, tempWeightFileVector);
0048 }
0049 
0050 //--------------------------------------------------------------------------------------------------
0051 void EGammaMvaEleEstimator::initialize(std::string methodName,
0052                                        EGammaMvaEleEstimator::MVAType type,
0053                                        Bool_t useBinnedVersion,
0054                                        std::vector<std::string> weightsfiles) {
0055   //clean up first
0056   for (unsigned int i = 0; i < fTMVAReader.size(); ++i) {
0057     if (fTMVAReader[i])
0058       delete fTMVAReader[i];
0059   }
0060   fTMVAReader.clear();
0061   fTMVAMethod.clear();
0062 
0063   //initialize
0064   fisInitialized = kTRUE;
0065   fMVAType = type;
0066   fMethodname = methodName;
0067   fUseBinnedVersion = useBinnedVersion;
0068 
0069   //Define expected number of bins
0070   UInt_t ExpectedNBins = 0;
0071   if (!fUseBinnedVersion) {
0072     ExpectedNBins = 1;
0073   } else if (type == kTrig) {
0074     ExpectedNBins = 6;
0075   } else if (type == kTrigNoIP) {
0076     ExpectedNBins = 6;
0077   } else if (type == kNonTrig) {
0078     ExpectedNBins = 6;
0079   } else if (type == kIsoRings) {
0080     ExpectedNBins = 4;
0081   } else if (type == kTrigIDIsoCombined) {
0082     ExpectedNBins = 6;
0083   } else if (type == kTrigIDIsoCombinedPUCorrected) {
0084     ExpectedNBins = 6;
0085   }
0086 
0087   fNMVABins = ExpectedNBins;
0088 
0089   //Check number of weight files given
0090   if (fNMVABins != weightsfiles.size()) {
0091     std::cout << "Error: Expected Number of bins = " << fNMVABins
0092               << " does not equal to weightsfiles.size() = " << weightsfiles.size() << std::endl;
0093 
0094 #ifndef STANDALONE
0095     assert(fNMVABins == weightsfiles.size());
0096 #endif
0097   }
0098 
0099   //Loop over all bins
0100   for (unsigned int i = 0; i < fNMVABins; ++i) {
0101     TMVA::Reader* tmpTMVAReader = new TMVA::Reader("!Color:!Silent:Error");
0102     tmpTMVAReader->SetVerbose(kTRUE);
0103 
0104     if (type == kTrig) {
0105       // Pure tracking variables
0106       tmpTMVAReader->AddVariable("fbrem", &fMVAVar_fbrem);
0107       tmpTMVAReader->AddVariable("kfchi2", &fMVAVar_kfchi2);
0108       tmpTMVAReader->AddVariable("kfhits", &fMVAVar_kfhits);
0109       tmpTMVAReader->AddVariable("gsfchi2", &fMVAVar_gsfchi2);
0110 
0111       // Geometrical matchings
0112       tmpTMVAReader->AddVariable("deta", &fMVAVar_deta);
0113       tmpTMVAReader->AddVariable("dphi", &fMVAVar_dphi);
0114       tmpTMVAReader->AddVariable("detacalo", &fMVAVar_detacalo);
0115 
0116       // Pure ECAL -> shower shapes
0117       tmpTMVAReader->AddVariable("see", &fMVAVar_see);
0118       tmpTMVAReader->AddVariable("spp", &fMVAVar_spp);
0119       tmpTMVAReader->AddVariable("etawidth", &fMVAVar_etawidth);
0120       tmpTMVAReader->AddVariable("phiwidth", &fMVAVar_phiwidth);
0121       tmpTMVAReader->AddVariable("e1x5e5x5", &fMVAVar_OneMinusE1x5E5x5);
0122       tmpTMVAReader->AddVariable("R9", &fMVAVar_R9);
0123 
0124       // Energy matching
0125       tmpTMVAReader->AddVariable("HoE", &fMVAVar_HoE);
0126       tmpTMVAReader->AddVariable("EoP", &fMVAVar_EoP);
0127       tmpTMVAReader->AddVariable("IoEmIoP", &fMVAVar_IoEmIoP);
0128       tmpTMVAReader->AddVariable("eleEoPout", &fMVAVar_eleEoPout);
0129       if (i == 2 || i == 5)
0130         tmpTMVAReader->AddVariable("PreShowerOverRaw", &fMVAVar_PreShowerOverRaw);
0131 
0132       if (!fUseBinnedVersion)
0133         tmpTMVAReader->AddVariable("PreShowerOverRaw", &fMVAVar_PreShowerOverRaw);
0134 
0135       // IP
0136       tmpTMVAReader->AddVariable("d0", &fMVAVar_d0);
0137       tmpTMVAReader->AddVariable("ip3d", &fMVAVar_ip3d);
0138 
0139       tmpTMVAReader->AddSpectator("eta", &fMVAVar_eta);
0140       tmpTMVAReader->AddSpectator("pt", &fMVAVar_pt);
0141     }
0142 
0143     if (type == kTrigNoIP) {
0144       // Pure tracking variables
0145       tmpTMVAReader->AddVariable("fbrem", &fMVAVar_fbrem);
0146       tmpTMVAReader->AddVariable("kfchi2", &fMVAVar_kfchi2);
0147       tmpTMVAReader->AddVariable("kfhits", &fMVAVar_kfhits);
0148       tmpTMVAReader->AddVariable("gsfchi2", &fMVAVar_gsfchi2);
0149 
0150       // Geometrical matchings
0151       tmpTMVAReader->AddVariable("deta", &fMVAVar_deta);
0152       tmpTMVAReader->AddVariable("dphi", &fMVAVar_dphi);
0153       tmpTMVAReader->AddVariable("detacalo", &fMVAVar_detacalo);
0154 
0155       // Pure ECAL -> shower shapes
0156       tmpTMVAReader->AddVariable("see", &fMVAVar_see);
0157       tmpTMVAReader->AddVariable("spp", &fMVAVar_spp);
0158       tmpTMVAReader->AddVariable("etawidth", &fMVAVar_etawidth);
0159       tmpTMVAReader->AddVariable("phiwidth", &fMVAVar_phiwidth);
0160       tmpTMVAReader->AddVariable("e1x5e5x5", &fMVAVar_OneMinusE1x5E5x5);
0161       tmpTMVAReader->AddVariable("R9", &fMVAVar_R9);
0162 
0163       // Energy matching
0164       tmpTMVAReader->AddVariable("HoE", &fMVAVar_HoE);
0165       tmpTMVAReader->AddVariable("EoP", &fMVAVar_EoP);
0166       tmpTMVAReader->AddVariable("IoEmIoP", &fMVAVar_IoEmIoP);
0167       tmpTMVAReader->AddVariable("eleEoPout", &fMVAVar_eleEoPout);
0168       tmpTMVAReader->AddVariable("rho", &fMVAVar_rho);
0169 
0170       if (i == 2 || i == 5)
0171         tmpTMVAReader->AddVariable("PreShowerOverRaw", &fMVAVar_PreShowerOverRaw);
0172 
0173       if (!fUseBinnedVersion)
0174         tmpTMVAReader->AddVariable("PreShowerOverRaw", &fMVAVar_PreShowerOverRaw);
0175 
0176       tmpTMVAReader->AddSpectator("eta", &fMVAVar_eta);
0177       tmpTMVAReader->AddSpectator("pt", &fMVAVar_pt);
0178     }
0179 
0180     if (type == kNonTrig) {
0181       // Pure tracking variables
0182       tmpTMVAReader->AddVariable("fbrem", &fMVAVar_fbrem);
0183       tmpTMVAReader->AddVariable("kfchi2", &fMVAVar_kfchi2);
0184       tmpTMVAReader->AddVariable("kfhits", &fMVAVar_kfhits);
0185       tmpTMVAReader->AddVariable("gsfchi2", &fMVAVar_gsfchi2);
0186 
0187       // Geometrical matchings
0188       tmpTMVAReader->AddVariable("deta", &fMVAVar_deta);
0189       tmpTMVAReader->AddVariable("dphi", &fMVAVar_dphi);
0190       tmpTMVAReader->AddVariable("detacalo", &fMVAVar_detacalo);
0191 
0192       // Pure ECAL -> shower shapes
0193       tmpTMVAReader->AddVariable("see", &fMVAVar_see);
0194       tmpTMVAReader->AddVariable("spp", &fMVAVar_spp);
0195       tmpTMVAReader->AddVariable("etawidth", &fMVAVar_etawidth);
0196       tmpTMVAReader->AddVariable("phiwidth", &fMVAVar_phiwidth);
0197       tmpTMVAReader->AddVariable("e1x5e5x5", &fMVAVar_OneMinusE1x5E5x5);
0198       tmpTMVAReader->AddVariable("R9", &fMVAVar_R9);
0199 
0200       // Energy matching
0201       tmpTMVAReader->AddVariable("HoE", &fMVAVar_HoE);
0202       tmpTMVAReader->AddVariable("EoP", &fMVAVar_EoP);
0203       tmpTMVAReader->AddVariable("IoEmIoP", &fMVAVar_IoEmIoP);
0204       tmpTMVAReader->AddVariable("eleEoPout", &fMVAVar_eleEoPout);
0205       if (i == 2 || i == 5)
0206         tmpTMVAReader->AddVariable("PreShowerOverRaw", &fMVAVar_PreShowerOverRaw);
0207 
0208       if (!fUseBinnedVersion)
0209         tmpTMVAReader->AddVariable("PreShowerOverRaw", &fMVAVar_PreShowerOverRaw);
0210 
0211       tmpTMVAReader->AddSpectator("eta", &fMVAVar_eta);
0212       tmpTMVAReader->AddSpectator("pt", &fMVAVar_pt);
0213     }
0214 
0215     if (type == kIsoRings) {
0216       tmpTMVAReader->AddVariable("ChargedIso_DR0p0To0p1", &fMVAVar_ChargedIso_DR0p0To0p1);
0217       tmpTMVAReader->AddVariable("ChargedIso_DR0p1To0p2", &fMVAVar_ChargedIso_DR0p1To0p2);
0218       tmpTMVAReader->AddVariable("ChargedIso_DR0p2To0p3", &fMVAVar_ChargedIso_DR0p2To0p3);
0219       tmpTMVAReader->AddVariable("ChargedIso_DR0p3To0p4", &fMVAVar_ChargedIso_DR0p3To0p4);
0220       tmpTMVAReader->AddVariable("ChargedIso_DR0p4To0p5", &fMVAVar_ChargedIso_DR0p4To0p5);
0221       tmpTMVAReader->AddVariable("GammaIso_DR0p0To0p1", &fMVAVar_GammaIso_DR0p0To0p1);
0222       tmpTMVAReader->AddVariable("GammaIso_DR0p1To0p2", &fMVAVar_GammaIso_DR0p1To0p2);
0223       tmpTMVAReader->AddVariable("GammaIso_DR0p2To0p3", &fMVAVar_GammaIso_DR0p2To0p3);
0224       tmpTMVAReader->AddVariable("GammaIso_DR0p3To0p4", &fMVAVar_GammaIso_DR0p3To0p4);
0225       tmpTMVAReader->AddVariable("GammaIso_DR0p4To0p5", &fMVAVar_GammaIso_DR0p4To0p5);
0226       tmpTMVAReader->AddVariable("NeutralHadronIso_DR0p0To0p1", &fMVAVar_NeutralHadronIso_DR0p0To0p1);
0227       tmpTMVAReader->AddVariable("NeutralHadronIso_DR0p1To0p2", &fMVAVar_NeutralHadronIso_DR0p1To0p2);
0228       tmpTMVAReader->AddVariable("NeutralHadronIso_DR0p2To0p3", &fMVAVar_NeutralHadronIso_DR0p2To0p3);
0229       tmpTMVAReader->AddVariable("NeutralHadronIso_DR0p3To0p4", &fMVAVar_NeutralHadronIso_DR0p3To0p4);
0230       tmpTMVAReader->AddVariable("NeutralHadronIso_DR0p4To0p5", &fMVAVar_NeutralHadronIso_DR0p4To0p5);
0231       tmpTMVAReader->AddSpectator("eta", &fMVAVar_eta);
0232       tmpTMVAReader->AddSpectator("pt", &fMVAVar_pt);
0233     }
0234 
0235     if (type == kTrigIDIsoCombinedPUCorrected) {
0236       // Pure tracking variables
0237       tmpTMVAReader->AddVariable("fbrem", &fMVAVar_fbrem);
0238       tmpTMVAReader->AddVariable("kfchi2", &fMVAVar_kfchi2);
0239       tmpTMVAReader->AddVariable("kflayers", &fMVAVar_kfhits);
0240       tmpTMVAReader->AddVariable("gsfchi2", &fMVAVar_gsfchi2);
0241 
0242       // Geometrical matchings
0243       tmpTMVAReader->AddVariable("deta", &fMVAVar_deta);
0244       tmpTMVAReader->AddVariable("dphi", &fMVAVar_dphi);
0245       tmpTMVAReader->AddVariable("detacalo", &fMVAVar_detacalo);
0246 
0247       // Pure ECAL -> shower shapes
0248       tmpTMVAReader->AddVariable("see", &fMVAVar_see);
0249       tmpTMVAReader->AddVariable("spp", &fMVAVar_spp);
0250       tmpTMVAReader->AddVariable("etawidth", &fMVAVar_etawidth);
0251       tmpTMVAReader->AddVariable("phiwidth", &fMVAVar_phiwidth);
0252       tmpTMVAReader->AddVariable("OneMinusSeedE1x5OverE5x5", &fMVAVar_OneMinusE1x5E5x5);
0253       tmpTMVAReader->AddVariable("R9", &fMVAVar_R9);
0254 
0255       // Energy matching
0256       tmpTMVAReader->AddVariable("HoE", &fMVAVar_HoE);
0257       tmpTMVAReader->AddVariable("EoP", &fMVAVar_EoP);
0258       tmpTMVAReader->AddVariable("IoEmIoP", &fMVAVar_IoEmIoP);
0259       tmpTMVAReader->AddVariable("EEleoPout", &fMVAVar_eleEoPout);
0260       if (i == 2 || i == 5) {
0261         tmpTMVAReader->AddVariable("PreShowerOverRaw", &fMVAVar_PreShowerOverRaw);
0262       }
0263       if (!fUseBinnedVersion) {
0264         tmpTMVAReader->AddVariable("PreShowerOverRaw", &fMVAVar_PreShowerOverRaw);
0265       }
0266 
0267       // IP
0268       tmpTMVAReader->AddVariable("d0", &fMVAVar_d0);
0269       tmpTMVAReader->AddVariable("ip3d", &fMVAVar_ip3d);
0270 
0271       //isolation variables
0272       tmpTMVAReader->AddVariable("ChargedIso_DR0p0To0p1", &fMVAVar_ChargedIso_DR0p0To0p1);
0273       tmpTMVAReader->AddVariable("ChargedIso_DR0p1To0p2", &fMVAVar_ChargedIso_DR0p1To0p2);
0274       tmpTMVAReader->AddVariable("ChargedIso_DR0p2To0p3", &fMVAVar_ChargedIso_DR0p2To0p3);
0275       tmpTMVAReader->AddVariable("ChargedIso_DR0p3To0p4", &fMVAVar_ChargedIso_DR0p3To0p4);
0276       tmpTMVAReader->AddVariable("ChargedIso_DR0p4To0p5", &fMVAVar_ChargedIso_DR0p4To0p5);
0277       tmpTMVAReader->AddVariable("GammaIso_DR0p0To0p1", &fMVAVar_GammaIso_DR0p0To0p1);
0278       tmpTMVAReader->AddVariable("GammaIso_DR0p1To0p2", &fMVAVar_GammaIso_DR0p1To0p2);
0279       tmpTMVAReader->AddVariable("GammaIso_DR0p2To0p3", &fMVAVar_GammaIso_DR0p2To0p3);
0280       tmpTMVAReader->AddVariable("GammaIso_DR0p3To0p4", &fMVAVar_GammaIso_DR0p3To0p4);
0281       tmpTMVAReader->AddVariable("GammaIso_DR0p4To0p5", &fMVAVar_GammaIso_DR0p4To0p5);
0282       tmpTMVAReader->AddVariable("NeutralHadronIso_DR0p0To0p1", &fMVAVar_NeutralHadronIso_DR0p0To0p1);
0283       tmpTMVAReader->AddVariable("NeutralHadronIso_DR0p1To0p2", &fMVAVar_NeutralHadronIso_DR0p1To0p2);
0284       tmpTMVAReader->AddVariable("NeutralHadronIso_DR0p2To0p3", &fMVAVar_NeutralHadronIso_DR0p2To0p3);
0285       tmpTMVAReader->AddVariable("NeutralHadronIso_DR0p3To0p4", &fMVAVar_NeutralHadronIso_DR0p3To0p4);
0286       tmpTMVAReader->AddVariable("NeutralHadronIso_DR0p4To0p5", &fMVAVar_NeutralHadronIso_DR0p4To0p5);
0287 
0288       //spectators
0289       tmpTMVAReader->AddSpectator("eta", &fMVAVar_eta);
0290       tmpTMVAReader->AddSpectator("pt", &fMVAVar_pt);
0291     }
0292 
0293     if (type == kTrigIDIsoCombined) {
0294       // Pure tracking variables
0295       tmpTMVAReader->AddVariable("fbrem", &fMVAVar_fbrem);
0296       tmpTMVAReader->AddVariable("kfchi2", &fMVAVar_kfchi2);
0297       tmpTMVAReader->AddVariable("kflayers", &fMVAVar_kfhits);
0298       tmpTMVAReader->AddVariable("gsfchi2", &fMVAVar_gsfchi2);
0299 
0300       // Geometrical matchings
0301       tmpTMVAReader->AddVariable("deta", &fMVAVar_deta);
0302       tmpTMVAReader->AddVariable("dphi", &fMVAVar_dphi);
0303       tmpTMVAReader->AddVariable("detacalo", &fMVAVar_detacalo);
0304 
0305       // Pure ECAL -> shower shapes
0306       tmpTMVAReader->AddVariable("see", &fMVAVar_see);
0307       tmpTMVAReader->AddVariable("spp", &fMVAVar_spp);
0308       tmpTMVAReader->AddVariable("etawidth", &fMVAVar_etawidth);
0309       tmpTMVAReader->AddVariable("phiwidth", &fMVAVar_phiwidth);
0310       tmpTMVAReader->AddVariable("OneMinusSeedE1x5OverE5x5", &fMVAVar_OneMinusE1x5E5x5);
0311       tmpTMVAReader->AddVariable("R9", &fMVAVar_R9);
0312 
0313       // Energy matching
0314       tmpTMVAReader->AddVariable("HoE", &fMVAVar_HoE);
0315       tmpTMVAReader->AddVariable("EoP", &fMVAVar_EoP);
0316       tmpTMVAReader->AddVariable("IoEmIoP", &fMVAVar_IoEmIoP);
0317       tmpTMVAReader->AddVariable("EEleoPout", &fMVAVar_eleEoPout);
0318       if (i == 2 || i == 5) {
0319         tmpTMVAReader->AddVariable("PreShowerOverRaw", &fMVAVar_PreShowerOverRaw);
0320       }
0321       if (!fUseBinnedVersion) {
0322         tmpTMVAReader->AddVariable("PreShowerOverRaw", &fMVAVar_PreShowerOverRaw);
0323       }
0324 
0325       // IP
0326       tmpTMVAReader->AddVariable("d0", &fMVAVar_d0);
0327       tmpTMVAReader->AddVariable("ip3d", &fMVAVar_ip3d);
0328 
0329       //isolation variables
0330       tmpTMVAReader->AddVariable("ChargedIso_DR0p0To0p1", &fMVAVar_ChargedIso_DR0p0To0p1);
0331       tmpTMVAReader->AddVariable("ChargedIso_DR0p1To0p2", &fMVAVar_ChargedIso_DR0p1To0p2);
0332       tmpTMVAReader->AddVariable("ChargedIso_DR0p2To0p3", &fMVAVar_ChargedIso_DR0p2To0p3);
0333       tmpTMVAReader->AddVariable("ChargedIso_DR0p3To0p4", &fMVAVar_ChargedIso_DR0p3To0p4);
0334       tmpTMVAReader->AddVariable("ChargedIso_DR0p4To0p5", &fMVAVar_ChargedIso_DR0p4To0p5);
0335       tmpTMVAReader->AddVariable("GammaIso_DR0p0To0p1", &fMVAVar_GammaIso_DR0p0To0p1);
0336       tmpTMVAReader->AddVariable("GammaIso_DR0p1To0p2", &fMVAVar_GammaIso_DR0p1To0p2);
0337       tmpTMVAReader->AddVariable("GammaIso_DR0p2To0p3", &fMVAVar_GammaIso_DR0p2To0p3);
0338       tmpTMVAReader->AddVariable("GammaIso_DR0p3To0p4", &fMVAVar_GammaIso_DR0p3To0p4);
0339       tmpTMVAReader->AddVariable("GammaIso_DR0p4To0p5", &fMVAVar_GammaIso_DR0p4To0p5);
0340       tmpTMVAReader->AddVariable("NeutralHadronIso_DR0p0To0p1", &fMVAVar_NeutralHadronIso_DR0p0To0p1);
0341       tmpTMVAReader->AddVariable("NeutralHadronIso_DR0p1To0p2", &fMVAVar_NeutralHadronIso_DR0p1To0p2);
0342       tmpTMVAReader->AddVariable("NeutralHadronIso_DR0p2To0p3", &fMVAVar_NeutralHadronIso_DR0p2To0p3);
0343       tmpTMVAReader->AddVariable("NeutralHadronIso_DR0p3To0p4", &fMVAVar_NeutralHadronIso_DR0p3To0p4);
0344       tmpTMVAReader->AddVariable("NeutralHadronIso_DR0p4To0p5", &fMVAVar_NeutralHadronIso_DR0p4To0p5);
0345       tmpTMVAReader->AddVariable("rho", &fMVAVar_rho);
0346 
0347       //spectators
0348       tmpTMVAReader->AddSpectator("eta", &fMVAVar_eta);
0349       tmpTMVAReader->AddSpectator("pt", &fMVAVar_pt);
0350     }
0351 
0352 #ifndef STANDALONE
0353     if ((fMethodname == "BDT") && (weightsfiles[i].rfind(".xml.gz") == weightsfiles[i].length() - strlen(".xml.gz"))) {
0354       gzFile file = gzopen(weightsfiles[i].c_str(), "rb");
0355       if (file == nullptr) {
0356         std::cout << "Error opening gzip file associated to " << weightsfiles[i] << std::endl;
0357         throw cms::Exception("Configuration", "Error reading zipped XML file");
0358       }
0359       std::vector<char> data;
0360       data.reserve(1024 * 1024 * 10);
0361       unsigned int bufflen = 32 * 1024;
0362       char* buff = reinterpret_cast<char*>(malloc(bufflen));
0363       if (buff == nullptr) {
0364         std::cout << "Error creating buffer for " << weightsfiles[i] << std::endl;
0365         gzclose(file);
0366         throw cms::Exception("Configuration", "Error reading zipped XML file");
0367       }
0368       int read;
0369       while ((read = gzread(file, buff, bufflen)) != 0) {
0370         if (read == -1) {
0371           std::cout << "Error reading gzip file associated to " << weightsfiles[i] << ": " << gzerror(file, &read)
0372                     << std::endl;
0373           gzclose(file);
0374           free(buff);
0375           throw cms::Exception("Configuration", "Error reading zipped XML file");
0376         }
0377         data.insert(data.end(), buff, buff + read);
0378       }
0379       if (gzclose(file) != Z_OK) {
0380         std::cout << "Error closing gzip file associated to " << weightsfiles[i] << std::endl;
0381       }
0382       free(buff);
0383       data.push_back('\0');  // IMPORTANT
0384       fTMVAMethod.push_back(dynamic_cast<TMVA::MethodBase*>(tmpTMVAReader->BookMVA(TMVA::Types::kBDT, &data[0])));
0385     } else {
0386       fTMVAMethod.push_back(dynamic_cast<TMVA::MethodBase*>(tmpTMVAReader->BookMVA(fMethodname, weightsfiles[i])));
0387     }
0388 #else
0389     fTMVAMethod.push_back(dynamic_cast<TMVA::MethodBase*>(tmpTMVAReader->BookMVA(fMethodname, weightsfiles[i])));
0390 #endif
0391     std::cout << "MVABin " << i << " : MethodName = " << fMethodname << " , type == " << type << " , "
0392               << "Load weights file : " << weightsfiles[i] << std::endl;
0393     fTMVAReader.push_back(tmpTMVAReader);
0394   }
0395   std::cout << "Electron ID MVA Completed\n";
0396 }
0397 
0398 //--------------------------------------------------------------------------------------------------
0399 UInt_t EGammaMvaEleEstimator::GetMVABin(double eta, double pt) const {
0400   //Default is to return the first bin
0401   unsigned int bin = 0;
0402 
0403   if (fMVAType == EGammaMvaEleEstimator::kIsoRings) {
0404     if (pt < 10 && fabs(eta) < 1.479)
0405       bin = 0;
0406     if (pt < 10 && fabs(eta) >= 1.479)
0407       bin = 1;
0408     if (pt >= 10 && fabs(eta) < 1.479)
0409       bin = 2;
0410     if (pt >= 10 && fabs(eta) >= 1.479)
0411       bin = 3;
0412   }
0413 
0414   if (fMVAType == EGammaMvaEleEstimator::kNonTrig) {
0415     bin = 0;
0416     if (pt < 10 && fabs(eta) < 0.8)
0417       bin = 0;
0418     if (pt < 10 && fabs(eta) >= 0.8 && fabs(eta) < 1.479)
0419       bin = 1;
0420     if (pt < 10 && fabs(eta) >= 1.479)
0421       bin = 2;
0422     if (pt >= 10 && fabs(eta) < 0.8)
0423       bin = 3;
0424     if (pt >= 10 && fabs(eta) >= 0.8 && fabs(eta) < 1.479)
0425       bin = 4;
0426     if (pt >= 10 && fabs(eta) >= 1.479)
0427       bin = 5;
0428   }
0429 
0430   if (fMVAType == EGammaMvaEleEstimator::kTrig || fMVAType == EGammaMvaEleEstimator::kTrigIDIsoCombined ||
0431       fMVAType == EGammaMvaEleEstimator::kTrigIDIsoCombinedPUCorrected ||
0432       fMVAType == EGammaMvaEleEstimator::kTrigNoIP) {
0433     bin = 0;
0434     if (pt < 20 && fabs(eta) < 0.8)
0435       bin = 0;
0436     if (pt < 20 && fabs(eta) >= 0.8 && fabs(eta) < 1.479)
0437       bin = 1;
0438     if (pt < 20 && fabs(eta) >= 1.479)
0439       bin = 2;
0440     if (pt >= 20 && fabs(eta) < 0.8)
0441       bin = 3;
0442     if (pt >= 20 && fabs(eta) >= 0.8 && fabs(eta) < 1.479)
0443       bin = 4;
0444     if (pt >= 20 && fabs(eta) >= 1.479)
0445       bin = 5;
0446   }
0447 
0448   return bin;
0449 }
0450 
0451 //--------------------------------------------------------------------------------------------------
0452 // for kTrig algorithm
0453 Double_t EGammaMvaEleEstimator::mvaValue(Double_t fbrem,
0454                                          Double_t kfchi2,
0455                                          Int_t kfhits,
0456                                          Double_t gsfchi2,
0457                                          Double_t deta,
0458                                          Double_t dphi,
0459                                          Double_t detacalo,
0460                                          //Double_t dphicalo,
0461                                          Double_t see,
0462                                          Double_t spp,
0463                                          Double_t etawidth,
0464                                          Double_t phiwidth,
0465                                          Double_t OneMinusE1x5E5x5,
0466                                          Double_t R9,
0467                                          //Int_t    nbrems,
0468                                          Double_t HoE,
0469                                          Double_t EoP,
0470                                          Double_t IoEmIoP,
0471                                          Double_t eleEoPout,
0472                                          Double_t PreShowerOverRaw,
0473                                          //Double_t EoPout,
0474                                          Double_t d0,
0475                                          Double_t ip3d,
0476                                          Double_t eta,
0477                                          Double_t pt,
0478                                          Bool_t printDebug) {
0479   if (!fisInitialized) {
0480     std::cout << "Error: EGammaMvaEleEstimator not properly initialized.\n";
0481     return -9999;
0482   }
0483 
0484   if (fMVAType != EGammaMvaEleEstimator::kTrig) {
0485     std::cout << "Error: This method should be called for kTrig MVA only" << std::endl;
0486     return -9999;
0487   }
0488 
0489   fMVAVar_fbrem = fbrem;
0490   fMVAVar_kfchi2 = kfchi2;
0491   fMVAVar_kfhits = float(kfhits);  // BTD does not support int variables
0492   fMVAVar_gsfchi2 = gsfchi2;
0493 
0494   fMVAVar_deta = deta;
0495   fMVAVar_dphi = dphi;
0496   fMVAVar_detacalo = detacalo;
0497 
0498   fMVAVar_see = see;
0499   fMVAVar_spp = spp;
0500   fMVAVar_etawidth = etawidth;
0501   fMVAVar_phiwidth = phiwidth;
0502   fMVAVar_OneMinusE1x5E5x5 = OneMinusE1x5E5x5;
0503   fMVAVar_R9 = R9;
0504 
0505   fMVAVar_HoE = HoE;
0506   fMVAVar_EoP = EoP;
0507   fMVAVar_IoEmIoP = IoEmIoP;
0508   fMVAVar_eleEoPout = eleEoPout;
0509   fMVAVar_PreShowerOverRaw = PreShowerOverRaw;
0510 
0511   fMVAVar_d0 = d0;
0512   fMVAVar_ip3d = ip3d;
0513   fMVAVar_eta = eta;
0514   fMVAVar_pt = pt;
0515 
0516   bindVariables();
0517   Double_t mva = -9999;
0518   if (fUseBinnedVersion) {
0519     int bin = GetMVABin(fMVAVar_eta, fMVAVar_pt);
0520     mva = fTMVAReader[bin]->EvaluateMVA(fTMVAMethod[bin]);
0521   } else {
0522     mva = fTMVAReader[0]->EvaluateMVA(fTMVAMethod[0]);
0523   }
0524 
0525   if (printDebug) {
0526     std::cout << " *** Inside the class fMethodname " << fMethodname << std::endl;
0527     std::cout << " fbrem " << fMVAVar_fbrem << " kfchi2 " << fMVAVar_kfchi2 << " mykfhits " << fMVAVar_kfhits
0528               << " gsfchi2 " << fMVAVar_gsfchi2 << " deta " << fMVAVar_deta << " dphi " << fMVAVar_dphi << " detacalo "
0529               << fMVAVar_detacalo << " see " << fMVAVar_see << " spp " << fMVAVar_spp << " etawidth "
0530               << fMVAVar_etawidth << " phiwidth " << fMVAVar_phiwidth << " OneMinusE1x5E5x5 "
0531               << fMVAVar_OneMinusE1x5E5x5 << " R9 " << fMVAVar_R9 << " HoE " << fMVAVar_HoE << " EoP " << fMVAVar_EoP
0532               << " IoEmIoP " << fMVAVar_IoEmIoP << " eleEoPout " << fMVAVar_eleEoPout << " PreShowerOverRaw "
0533               << fMVAVar_PreShowerOverRaw << " d0 " << fMVAVar_d0 << " ip3d " << fMVAVar_ip3d << " eta " << fMVAVar_eta
0534               << " pt " << fMVAVar_pt << std::endl;
0535     std::cout << " ### MVA " << mva << std::endl;
0536   }
0537 
0538   return mva;
0539 }
0540 
0541 //--------------------------------------------------------------------------------------------------------
0542 // for kTrigNoIP algorithm
0543 Double_t EGammaMvaEleEstimator::mvaValue(Double_t fbrem,
0544                                          Double_t kfchi2,
0545                                          Int_t kfhits,
0546                                          Double_t gsfchi2,
0547                                          Double_t deta,
0548                                          Double_t dphi,
0549                                          Double_t detacalo,
0550                                          Double_t see,
0551                                          Double_t spp,
0552                                          Double_t etawidth,
0553                                          Double_t phiwidth,
0554                                          Double_t e1x5e5x5,
0555                                          Double_t R9,
0556                                          Double_t HoE,
0557                                          Double_t EoP,
0558                                          Double_t IoEmIoP,
0559                                          Double_t eleEoPout,
0560                                          Double_t rho,
0561                                          Double_t PreShowerOverRaw,
0562                                          Double_t eta,
0563                                          Double_t pt,
0564                                          Bool_t printDebug) {
0565   if (!fisInitialized) {
0566     std::cout << "Error: EGammaMvaEleEstimator not properly initialized.\n";
0567     return -9999;
0568   }
0569 
0570   if (fMVAType != EGammaMvaEleEstimator::kTrigNoIP) {
0571     std::cout << "Error: This method should be called for kTrigNoIP MVA only" << std::endl;
0572     return -9999;
0573   }
0574 
0575   fMVAVar_fbrem = fbrem;
0576   fMVAVar_kfchi2 = kfchi2;
0577   fMVAVar_kfhits = float(kfhits);  // BTD does not support int variables
0578   fMVAVar_gsfchi2 = gsfchi2;
0579 
0580   fMVAVar_deta = deta;
0581   fMVAVar_dphi = dphi;
0582   fMVAVar_detacalo = detacalo;
0583 
0584   fMVAVar_see = see;
0585   fMVAVar_spp = spp;
0586   fMVAVar_etawidth = etawidth;
0587   fMVAVar_phiwidth = phiwidth;
0588   fMVAVar_OneMinusE1x5E5x5 = e1x5e5x5;
0589   fMVAVar_R9 = R9;
0590 
0591   fMVAVar_HoE = HoE;
0592   fMVAVar_EoP = EoP;
0593   fMVAVar_IoEmIoP = IoEmIoP;
0594   fMVAVar_eleEoPout = eleEoPout;
0595   fMVAVar_rho = rho;
0596   fMVAVar_PreShowerOverRaw = PreShowerOverRaw;
0597 
0598   fMVAVar_eta = eta;
0599   fMVAVar_pt = pt;
0600 
0601   bindVariables();
0602   Double_t mva = -9999;
0603   if (fUseBinnedVersion) {
0604     int bin = GetMVABin(fMVAVar_eta, fMVAVar_pt);
0605     mva = fTMVAReader[bin]->EvaluateMVA(fTMVAMethod[bin]);
0606   } else {
0607     mva = fTMVAReader[0]->EvaluateMVA(fTMVAMethod[0]);
0608   }
0609 
0610   if (printDebug) {
0611     std::cout << " *** Inside the class fMethodname " << fMethodname << std::endl;
0612     std::cout << " fbrem " << fMVAVar_fbrem << " kfchi2 " << fMVAVar_kfchi2 << " mykfhits " << fMVAVar_kfhits
0613               << " gsfchi2 " << fMVAVar_gsfchi2 << " deta " << fMVAVar_deta << " dphi " << fMVAVar_dphi << " detacalo "
0614               << fMVAVar_detacalo << " see " << fMVAVar_see << " spp " << fMVAVar_spp << " etawidth "
0615               << fMVAVar_etawidth << " phiwidth " << fMVAVar_phiwidth << " e1x5e5x5 " << fMVAVar_OneMinusE1x5E5x5
0616               << " R9 " << fMVAVar_R9 << " HoE " << fMVAVar_HoE << " EoP " << fMVAVar_EoP << " IoEmIoP "
0617               << fMVAVar_IoEmIoP << " eleEoPout " << fMVAVar_eleEoPout << " rho " << fMVAVar_rho << " PreShowerOverRaw "
0618               << fMVAVar_PreShowerOverRaw << " eta " << fMVAVar_eta << " pt " << fMVAVar_pt << std::endl;
0619     std::cout << " ### MVA " << mva << std::endl;
0620   }
0621 
0622   return mva;
0623 }
0624 
0625 //--------------------------------------------------------------------------------------------------
0626 // for kNonTrig algorithm
0627 Double_t EGammaMvaEleEstimator::mvaValue(Double_t fbrem,
0628                                          Double_t kfchi2,
0629                                          Int_t kfhits,
0630                                          Double_t gsfchi2,
0631                                          Double_t deta,
0632                                          Double_t dphi,
0633                                          Double_t detacalo,
0634                                          //Double_t dphicalo,
0635                                          Double_t see,
0636                                          Double_t spp,
0637                                          Double_t etawidth,
0638                                          Double_t phiwidth,
0639                                          Double_t OneMinusE1x5E5x5,
0640                                          Double_t R9,
0641                                          //Int_t    nbrems,
0642                                          Double_t HoE,
0643                                          Double_t EoP,
0644                                          Double_t IoEmIoP,
0645                                          Double_t eleEoPout,
0646                                          Double_t PreShowerOverRaw,
0647                                          //Double_t EoPout,
0648                                          Double_t eta,
0649                                          Double_t pt,
0650                                          Bool_t printDebug) {
0651   if (!fisInitialized) {
0652     std::cout << "Error: EGammaMvaEleEstimator not properly initialized.\n";
0653     return -9999;
0654   }
0655 
0656   if (fMVAType != EGammaMvaEleEstimator::kNonTrig) {
0657     std::cout << "Error: This method should be called for kNonTrig MVA only" << std::endl;
0658     return -9999;
0659   }
0660 
0661   fMVAVar_fbrem = fbrem;
0662   fMVAVar_kfchi2 = kfchi2;
0663   fMVAVar_kfhits = float(kfhits);  // BTD does not support int variables
0664   fMVAVar_gsfchi2 = gsfchi2;
0665 
0666   fMVAVar_deta = deta;
0667   fMVAVar_dphi = dphi;
0668   fMVAVar_detacalo = detacalo;
0669 
0670   fMVAVar_see = see;
0671   fMVAVar_spp = spp;
0672   fMVAVar_etawidth = etawidth;
0673   fMVAVar_phiwidth = phiwidth;
0674   fMVAVar_OneMinusE1x5E5x5 = OneMinusE1x5E5x5;
0675   fMVAVar_R9 = R9;
0676 
0677   fMVAVar_HoE = HoE;
0678   fMVAVar_EoP = EoP;
0679   fMVAVar_IoEmIoP = IoEmIoP;
0680   fMVAVar_eleEoPout = eleEoPout;
0681   fMVAVar_PreShowerOverRaw = PreShowerOverRaw;
0682 
0683   fMVAVar_eta = eta;
0684   fMVAVar_pt = pt;
0685 
0686   bindVariables();
0687   Double_t mva = -9999;
0688   if (fUseBinnedVersion) {
0689     int bin = GetMVABin(fMVAVar_eta, fMVAVar_pt);
0690     mva = fTMVAReader[bin]->EvaluateMVA(fTMVAMethod[bin]);
0691   } else {
0692     mva = fTMVAReader[0]->EvaluateMVA(fTMVAMethod[0]);
0693   }
0694 
0695   if (printDebug) {
0696     std::cout << " *** Inside the class fMethodname " << fMethodname << std::endl;
0697     std::cout << " fbrem " << fMVAVar_fbrem << " kfchi2 " << fMVAVar_kfchi2 << " mykfhits " << fMVAVar_kfhits
0698               << " gsfchi2 " << fMVAVar_gsfchi2 << " deta " << fMVAVar_deta << " dphi " << fMVAVar_dphi << " detacalo "
0699               << fMVAVar_detacalo << " see " << fMVAVar_see << " spp " << fMVAVar_spp << " etawidth "
0700               << fMVAVar_etawidth << " phiwidth " << fMVAVar_phiwidth << " OneMinusE1x5E5x5 "
0701               << fMVAVar_OneMinusE1x5E5x5 << " R9 " << fMVAVar_R9 << " HoE " << fMVAVar_HoE << " EoP " << fMVAVar_EoP
0702               << " IoEmIoP " << fMVAVar_IoEmIoP << " eleEoPout " << fMVAVar_eleEoPout << " PreShowerOverRaw "
0703               << fMVAVar_PreShowerOverRaw << " eta " << fMVAVar_eta << " pt " << fMVAVar_pt << std::endl;
0704     std::cout << " ### MVA " << mva << std::endl;
0705   }
0706 
0707   return mva;
0708 }
0709 
0710 //--------------------------------------------------------------------------------------------------
0711 Double_t EGammaMvaEleEstimator::IDIsoCombinedMvaValue(Double_t fbrem,
0712                                                       Double_t kfchi2,
0713                                                       Int_t kfhits,
0714                                                       Double_t gsfchi2,
0715                                                       Double_t deta,
0716                                                       Double_t dphi,
0717                                                       Double_t detacalo,
0718                                                       Double_t see,
0719                                                       Double_t spp,
0720                                                       Double_t etawidth,
0721                                                       Double_t phiwidth,
0722                                                       Double_t OneMinusE1x5E5x5,
0723                                                       Double_t R9,
0724                                                       Double_t HoE,
0725                                                       Double_t EoP,
0726                                                       Double_t IoEmIoP,
0727                                                       Double_t eleEoPout,
0728                                                       Double_t PreShowerOverRaw,
0729                                                       Double_t d0,
0730                                                       Double_t ip3d,
0731                                                       Double_t ChargedIso_DR0p0To0p1,
0732                                                       Double_t ChargedIso_DR0p1To0p2,
0733                                                       Double_t ChargedIso_DR0p2To0p3,
0734                                                       Double_t ChargedIso_DR0p3To0p4,
0735                                                       Double_t ChargedIso_DR0p4To0p5,
0736                                                       Double_t GammaIso_DR0p0To0p1,
0737                                                       Double_t GammaIso_DR0p1To0p2,
0738                                                       Double_t GammaIso_DR0p2To0p3,
0739                                                       Double_t GammaIso_DR0p3To0p4,
0740                                                       Double_t GammaIso_DR0p4To0p5,
0741                                                       Double_t NeutralHadronIso_DR0p0To0p1,
0742                                                       Double_t NeutralHadronIso_DR0p1To0p2,
0743                                                       Double_t NeutralHadronIso_DR0p2To0p3,
0744                                                       Double_t NeutralHadronIso_DR0p3To0p4,
0745                                                       Double_t NeutralHadronIso_DR0p4To0p5,
0746                                                       Double_t Rho,
0747                                                       Double_t eta,
0748                                                       Double_t pt,
0749                                                       Bool_t printDebug) {
0750   if (!fisInitialized) {
0751     std::cout << "Error: EGammaMvaEleEstimator not properly initialized.\n";
0752     return -9999;
0753   }
0754 
0755   fMVAVar_fbrem = (fbrem < -1.0) ? -1.0 : fbrem;
0756   fMVAVar_kfchi2 = (kfchi2 > 10) ? 10 : kfchi2;
0757   fMVAVar_kfhits = float(kfhits);  // BTD does not support int variables
0758   fMVAVar_gsfchi2 = (gsfchi2 > 200) ? 200 : gsfchi2;
0759   fMVAVar_deta = (fabs(deta) > 0.06) ? 0.06 : fabs(deta);
0760   fMVAVar_dphi = dphi;
0761   fMVAVar_detacalo = detacalo;
0762 
0763   fMVAVar_see = see;
0764   fMVAVar_spp = spp;
0765   fMVAVar_etawidth = etawidth;
0766   fMVAVar_phiwidth = phiwidth;
0767   fMVAVar_OneMinusE1x5E5x5 = std::max(std::min(double(OneMinusE1x5E5x5), 2.0), -1.0);
0768   fMVAVar_R9 = (R9 > 5) ? 5 : R9;
0769 
0770   fMVAVar_HoE = HoE;
0771   fMVAVar_EoP = (EoP > 20) ? 20 : EoP;
0772   fMVAVar_IoEmIoP = IoEmIoP;
0773   fMVAVar_eleEoPout = (eleEoPout > 20) ? 20 : eleEoPout;
0774   fMVAVar_PreShowerOverRaw = PreShowerOverRaw;
0775 
0776   fMVAVar_d0 = d0;
0777   fMVAVar_ip3d = ip3d;
0778 
0779   fMVAVar_ChargedIso_DR0p0To0p1 = ChargedIso_DR0p0To0p1;
0780   fMVAVar_ChargedIso_DR0p1To0p2 = ChargedIso_DR0p1To0p2;
0781   fMVAVar_ChargedIso_DR0p2To0p3 = ChargedIso_DR0p2To0p3;
0782   fMVAVar_ChargedIso_DR0p3To0p4 = ChargedIso_DR0p3To0p4;
0783   fMVAVar_ChargedIso_DR0p4To0p5 = ChargedIso_DR0p4To0p5;
0784   fMVAVar_GammaIso_DR0p0To0p1 = GammaIso_DR0p0To0p1;
0785   fMVAVar_GammaIso_DR0p1To0p2 = GammaIso_DR0p1To0p2;
0786   fMVAVar_GammaIso_DR0p2To0p3 = GammaIso_DR0p2To0p3;
0787   fMVAVar_GammaIso_DR0p3To0p4 = GammaIso_DR0p3To0p4;
0788   fMVAVar_GammaIso_DR0p4To0p5 = GammaIso_DR0p4To0p5;
0789   fMVAVar_NeutralHadronIso_DR0p0To0p1 = NeutralHadronIso_DR0p0To0p1;
0790   fMVAVar_NeutralHadronIso_DR0p1To0p2 = NeutralHadronIso_DR0p1To0p2;
0791   fMVAVar_NeutralHadronIso_DR0p2To0p3 = NeutralHadronIso_DR0p2To0p3;
0792   fMVAVar_NeutralHadronIso_DR0p3To0p4 = NeutralHadronIso_DR0p3To0p4;
0793   fMVAVar_NeutralHadronIso_DR0p4To0p5 = NeutralHadronIso_DR0p4To0p5;
0794 
0795   fMVAVar_rho = Rho;
0796   fMVAVar_eta = eta;
0797   fMVAVar_pt = pt;
0798 
0799   Double_t mva = -9999;
0800   if (fUseBinnedVersion) {
0801     int bin = GetMVABin(fMVAVar_eta, fMVAVar_pt);
0802     mva = fTMVAReader[bin]->EvaluateMVA(fTMVAMethod[bin]);
0803   } else {
0804     mva = fTMVAReader[0]->EvaluateMVA(fTMVAMethod[0]);
0805   }
0806 
0807   if (printDebug) {
0808     std::cout << " *** Inside the class fMethodname " << fMethodname << std::endl;
0809     std::cout << " fbrem " << fMVAVar_fbrem << " kfchi2 " << fMVAVar_kfchi2 << " mykfhits " << fMVAVar_kfhits
0810               << " gsfchi2 " << fMVAVar_gsfchi2 << " deta " << fMVAVar_deta << " dphi " << fMVAVar_dphi << " detacalo "
0811               << fMVAVar_detacalo << " see " << fMVAVar_see << " spp " << fMVAVar_spp << " etawidth "
0812               << fMVAVar_etawidth << " phiwidth " << fMVAVar_phiwidth << " OneMinusE1x5E5x5 "
0813               << fMVAVar_OneMinusE1x5E5x5 << " R9 " << fMVAVar_R9 << " HoE " << fMVAVar_HoE << " EoP " << fMVAVar_EoP
0814               << " IoEmIoP " << fMVAVar_IoEmIoP << " eleEoPout " << fMVAVar_eleEoPout << " PreShowerOverRaw "
0815               << fMVAVar_PreShowerOverRaw << " d0 " << fMVAVar_d0 << " ip3d " << fMVAVar_ip3d
0816               << " ChargedIso_DR0p0To0p1 " << ChargedIso_DR0p0To0p1 << " ChargedIso_DR0p1To0p2 "
0817               << ChargedIso_DR0p1To0p2 << " ChargedIso_DR0p2To0p3 " << ChargedIso_DR0p2To0p3
0818               << " ChargedIso_DR0p3To0p4 " << ChargedIso_DR0p3To0p4 << " ChargedIso_DR0p4To0p5 "
0819               << ChargedIso_DR0p4To0p5 << " GammaIso_DR0p0To0p1 " << GammaIso_DR0p0To0p1 << " GammaIso_DR0p1To0p2 "
0820               << GammaIso_DR0p1To0p2 << " GammaIso_DR0p2To0p3 " << GammaIso_DR0p2To0p3 << " GammaIso_DR0p3To0p4 "
0821               << GammaIso_DR0p3To0p4 << " GammaIso_DR0p4To0p5 " << GammaIso_DR0p4To0p5
0822               << " NeutralHadronIso_DR0p0To0p1 " << NeutralHadronIso_DR0p0To0p1 << " NeutralHadronIso_DR0p1To0p2 "
0823               << NeutralHadronIso_DR0p1To0p2 << " NeutralHadronIso_DR0p2To0p3 " << NeutralHadronIso_DR0p2To0p3
0824               << " NeutralHadronIso_DR0p3To0p4 " << NeutralHadronIso_DR0p3To0p4 << " NeutralHadronIso_DR0p4To0p5 "
0825               << NeutralHadronIso_DR0p4To0p5 << " Rho " << Rho << " eta " << fMVAVar_eta << " pt " << fMVAVar_pt
0826               << std::endl;
0827     std::cout << " ### MVA " << mva << std::endl;
0828   }
0829 
0830   return mva;
0831 }
0832 
0833 #ifndef STANDALONE
0834 Double_t EGammaMvaEleEstimator::isoMvaValue(Double_t Pt,
0835                                             Double_t Eta,
0836                                             Double_t Rho,
0837                                             ElectronEffectiveArea::ElectronEffectiveAreaTarget EATarget,
0838                                             Double_t ChargedIso_DR0p0To0p1,
0839                                             Double_t ChargedIso_DR0p1To0p2,
0840                                             Double_t ChargedIso_DR0p2To0p3,
0841                                             Double_t ChargedIso_DR0p3To0p4,
0842                                             Double_t ChargedIso_DR0p4To0p5,
0843                                             Double_t GammaIso_DR0p0To0p1,
0844                                             Double_t GammaIso_DR0p1To0p2,
0845                                             Double_t GammaIso_DR0p2To0p3,
0846                                             Double_t GammaIso_DR0p3To0p4,
0847                                             Double_t GammaIso_DR0p4To0p5,
0848                                             Double_t NeutralHadronIso_DR0p0To0p1,
0849                                             Double_t NeutralHadronIso_DR0p1To0p2,
0850                                             Double_t NeutralHadronIso_DR0p2To0p3,
0851                                             Double_t NeutralHadronIso_DR0p3To0p4,
0852                                             Double_t NeutralHadronIso_DR0p4To0p5,
0853                                             Bool_t printDebug) {
0854   if (!fisInitialized) {
0855     std::cout << "Error: EGammaMvaEleEstimator not properly initialized.\n";
0856     return -9999;
0857   }
0858 
0859   fMVAVar_ChargedIso_DR0p0To0p1 = TMath::Min((ChargedIso_DR0p0To0p1) / Pt, 2.5);
0860   fMVAVar_ChargedIso_DR0p1To0p2 = TMath::Min((ChargedIso_DR0p1To0p2) / Pt, 2.5);
0861   fMVAVar_ChargedIso_DR0p2To0p3 = TMath::Min((ChargedIso_DR0p2To0p3) / Pt, 2.5);
0862   fMVAVar_ChargedIso_DR0p3To0p4 = TMath::Min((ChargedIso_DR0p3To0p4) / Pt, 2.5);
0863   fMVAVar_ChargedIso_DR0p4To0p5 = TMath::Min((ChargedIso_DR0p4To0p5) / Pt, 2.5);
0864   fMVAVar_GammaIso_DR0p0To0p1 = TMath::Max(
0865       TMath::Min((GammaIso_DR0p0To0p1 - Rho * ElectronEffectiveArea::GetElectronEffectiveArea(
0866                                                   ElectronEffectiveArea::kEleGammaIsoDR0p0To0p1, Eta, EATarget)) /
0867                      Pt,
0868                  2.5),
0869       0.0);
0870   fMVAVar_GammaIso_DR0p1To0p2 = TMath::Max(
0871       TMath::Min((GammaIso_DR0p1To0p2 - Rho * ElectronEffectiveArea::GetElectronEffectiveArea(
0872                                                   ElectronEffectiveArea::kEleGammaIsoDR0p1To0p2, Eta, EATarget)) /
0873                      Pt,
0874                  2.5),
0875       0.0);
0876   fMVAVar_GammaIso_DR0p2To0p3 = TMath::Max(
0877       TMath::Min((GammaIso_DR0p2To0p3 - Rho * ElectronEffectiveArea::GetElectronEffectiveArea(
0878                                                   ElectronEffectiveArea::kEleGammaIsoDR0p2To0p3, Eta, EATarget)) /
0879                      Pt,
0880                  2.5),
0881       0.0);
0882   fMVAVar_GammaIso_DR0p3To0p4 = TMath::Max(
0883       TMath::Min((GammaIso_DR0p3To0p4 - Rho * ElectronEffectiveArea::GetElectronEffectiveArea(
0884                                                   ElectronEffectiveArea::kEleGammaIsoDR0p3To0p4, Eta, EATarget)) /
0885                      Pt,
0886                  2.5),
0887       0.0);
0888   fMVAVar_GammaIso_DR0p4To0p5 = TMath::Max(
0889       TMath::Min((GammaIso_DR0p4To0p5 - Rho * ElectronEffectiveArea::GetElectronEffectiveArea(
0890                                                   ElectronEffectiveArea::kEleGammaIsoDR0p4To0p5, Eta, EATarget)) /
0891                      Pt,
0892                  2.5),
0893       0.0);
0894   fMVAVar_NeutralHadronIso_DR0p0To0p1 =
0895       TMath::Max(TMath::Min((NeutralHadronIso_DR0p0To0p1 -
0896                              Rho * ElectronEffectiveArea::GetElectronEffectiveArea(
0897                                        ElectronEffectiveArea::kEleNeutralHadronIsoDR0p0To0p1, Eta, EATarget)) /
0898                                 Pt,
0899                             2.5),
0900                  0.0);
0901   fMVAVar_NeutralHadronIso_DR0p1To0p2 =
0902       TMath::Max(TMath::Min((NeutralHadronIso_DR0p1To0p2 -
0903                              Rho * ElectronEffectiveArea::GetElectronEffectiveArea(
0904                                        ElectronEffectiveArea::kEleNeutralHadronIsoDR0p1To0p2, Eta, EATarget)) /
0905                                 Pt,
0906                             2.5),
0907                  0.0);
0908   fMVAVar_NeutralHadronIso_DR0p2To0p3 =
0909       TMath::Max(TMath::Min((NeutralHadronIso_DR0p2To0p3 -
0910                              Rho * ElectronEffectiveArea::GetElectronEffectiveArea(
0911                                        ElectronEffectiveArea::kEleNeutralHadronIsoDR0p2To0p3, Eta, EATarget)) /
0912                                 Pt,
0913                             2.5),
0914                  0.0);
0915   fMVAVar_NeutralHadronIso_DR0p3To0p4 =
0916       TMath::Max(TMath::Min((NeutralHadronIso_DR0p3To0p4 -
0917                              Rho * ElectronEffectiveArea::GetElectronEffectiveArea(
0918                                        ElectronEffectiveArea::kEleNeutralHadronIsoDR0p3To0p4, Eta, EATarget)) /
0919                                 Pt,
0920                             2.5),
0921                  0.0);
0922   fMVAVar_NeutralHadronIso_DR0p4To0p5 =
0923       TMath::Max(TMath::Min((NeutralHadronIso_DR0p4To0p5 -
0924                              Rho * ElectronEffectiveArea::GetElectronEffectiveArea(
0925                                        ElectronEffectiveArea::kEleNeutralHadronIsoDR0p4To0p5, Eta, EATarget)) /
0926                                 Pt,
0927                             2.5),
0928                  0.0);
0929 
0930   // evaluate
0931   int bin = GetMVABin(Eta, Pt);
0932   Double_t mva = fTMVAReader[bin]->EvaluateMVA(fTMVAMethod[bin]);
0933 
0934   if (printDebug) {
0935     std::cout << " *** Inside the class fMethodname " << fMethodname << " fMVAType " << fMVAType << std::endl;
0936     std::cout << "ChargedIso ( 0.0 | 0.1 | 0.2 | 0.3 | 0.4 | 0.5 ): " << fMVAVar_ChargedIso_DR0p0To0p1 << " "
0937               << fMVAVar_ChargedIso_DR0p1To0p2 << " " << fMVAVar_ChargedIso_DR0p2To0p3 << " "
0938               << fMVAVar_ChargedIso_DR0p3To0p4 << " " << fMVAVar_ChargedIso_DR0p4To0p5 << std::endl;
0939     std::cout << "PF Gamma Iso ( 0.0 | 0.1 | 0.2 | 0.3 | 0.4 | 0.5 ): " << fMVAVar_GammaIso_DR0p0To0p1 << " "
0940               << fMVAVar_GammaIso_DR0p1To0p2 << " " << fMVAVar_GammaIso_DR0p2To0p3 << " " << fMVAVar_GammaIso_DR0p3To0p4
0941               << " " << fMVAVar_GammaIso_DR0p4To0p5 << std::endl;
0942     std::cout << "PF Neutral Hadron Iso ( 0.0 | 0.1 | 0.2 | 0.3 | 0.4 | 0.5 ): " << fMVAVar_NeutralHadronIso_DR0p0To0p1
0943               << " " << fMVAVar_NeutralHadronIso_DR0p1To0p2 << " " << fMVAVar_NeutralHadronIso_DR0p2To0p3 << " "
0944               << fMVAVar_NeutralHadronIso_DR0p3To0p4 << " " << fMVAVar_NeutralHadronIso_DR0p4To0p5 << " " << std::endl;
0945     std::cout << " ### MVA " << mva << std::endl;
0946   }
0947 
0948   return mva;
0949 }
0950 
0951 //--------------------------------------------------------------------------------------------------
0952 
0953 // for kTrig and kNonTrig algorithm
0954 Double_t EGammaMvaEleEstimator::mvaValue(const reco::GsfElectron& ele,
0955                                          const reco::Vertex& vertex,
0956                                          const TransientTrackBuilder& transientTrackBuilder,
0957                                          EcalClusterLazyTools const& myEcalCluster,
0958                                          bool printDebug) {
0959   if (!fisInitialized) {
0960     std::cout << "Error: EGammaMvaEleEstimator not properly initialized.\n";
0961     return -9999;
0962   }
0963 
0964   if ((fMVAType != EGammaMvaEleEstimator::kTrig) && (fMVAType != EGammaMvaEleEstimator::kNonTrig)) {
0965     std::cout << "Error: This method should be called for kTrig or kNonTrig MVA only" << std::endl;
0966     return -9999;
0967   }
0968 
0969   bool validKF = false;
0970   reco::TrackRef myTrackRef = ele.closestCtfTrackRef();
0971   validKF = (myTrackRef.isAvailable());
0972   validKF &= (myTrackRef.isNonnull());
0973 
0974   // Pure tracking variables
0975   fMVAVar_fbrem = ele.fbrem();
0976   fMVAVar_kfchi2 = (validKF) ? myTrackRef->normalizedChi2() : 0;
0977   fMVAVar_kfhits = (validKF) ? myTrackRef->hitPattern().trackerLayersWithMeasurement() : -1.;
0978   fMVAVar_kfhitsall =
0979       (validKF) ? myTrackRef->numberOfValidHits() : -1.;  //  save also this in your ntuple as possible alternative
0980   fMVAVar_gsfchi2 = ele.gsfTrack()->normalizedChi2();
0981 
0982   // Geometrical matchings
0983   fMVAVar_deta = ele.deltaEtaSuperClusterTrackAtVtx();
0984   fMVAVar_dphi = ele.deltaPhiSuperClusterTrackAtVtx();
0985   fMVAVar_detacalo = ele.deltaEtaSeedClusterTrackAtCalo();
0986 
0987   // Pure ECAL -> shower shapes
0988   fMVAVar_see = ele.sigmaIetaIeta();  //EleSigmaIEtaIEta
0989   const auto& vCov = myEcalCluster.localCovariances(*(ele.superCluster()->seed()));
0990   if (edm::isFinite(vCov[2]))
0991     fMVAVar_spp = sqrt(vCov[2]);  //EleSigmaIPhiIPhi
0992   else
0993     fMVAVar_spp = 0.;
0994 
0995   fMVAVar_etawidth = ele.superCluster()->etaWidth();
0996   fMVAVar_phiwidth = ele.superCluster()->phiWidth();
0997   fMVAVar_OneMinusE1x5E5x5 = (ele.e5x5()) != 0. ? 1. - (ele.e1x5() / ele.e5x5()) : -1.;
0998   fMVAVar_R9 = myEcalCluster.e3x3(*(ele.superCluster()->seed())) / ele.superCluster()->rawEnergy();
0999 
1000   // Energy matching
1001   fMVAVar_HoE = ele.hadronicOverEm();
1002   fMVAVar_EoP = ele.eSuperClusterOverP();
1003   fMVAVar_IoEmIoP = (1.0 / ele.ecalEnergy()) - (1.0 / ele.p());  // in the future to be changed with ele.gsfTrack()->p()
1004   fMVAVar_eleEoPout = ele.eEleClusterOverPout();
1005   fMVAVar_PreShowerOverRaw = ele.superCluster()->preshowerEnergy() / ele.superCluster()->rawEnergy();
1006 
1007   // Spectators
1008   fMVAVar_eta = ele.superCluster()->eta();
1009   fMVAVar_pt = ele.pt();
1010 
1011   // for triggering electrons get the impact parameteres
1012   if (fMVAType == kTrig) {
1013     //d0
1014     if (ele.gsfTrack().isNonnull()) {
1015       fMVAVar_d0 = (-1.0) * ele.gsfTrack()->dxy(vertex.position());
1016     } else if (ele.closestCtfTrackRef().isNonnull()) {
1017       fMVAVar_d0 = (-1.0) * ele.closestCtfTrackRef()->dxy(vertex.position());
1018     } else {
1019       fMVAVar_d0 = -9999.0;
1020     }
1021 
1022     //default values for IP3D
1023     fMVAVar_ip3d = -999.0;
1024     fMVAVar_ip3dSig = 0.0;
1025     if (ele.gsfTrack().isNonnull()) {
1026       const double gsfsign = ((-ele.gsfTrack()->dxy(vertex.position())) >= 0) ? 1. : -1.;
1027 
1028       const reco::TransientTrack& tt = transientTrackBuilder.build(ele.gsfTrack());
1029       const std::pair<bool, Measurement1D>& ip3dpv = IPTools::absoluteImpactParameter3D(tt, vertex);
1030       if (ip3dpv.first) {
1031         double ip3d = gsfsign * ip3dpv.second.value();
1032         double ip3derr = ip3dpv.second.error();
1033         fMVAVar_ip3d = ip3d;
1034         fMVAVar_ip3dSig = ip3d / ip3derr;
1035       }
1036     }
1037   }
1038 
1039   // evaluate
1040   bindVariables();
1041   Double_t mva = -9999;
1042   if (fUseBinnedVersion) {
1043     int bin = GetMVABin(fMVAVar_eta, fMVAVar_pt);
1044     mva = fTMVAReader[bin]->EvaluateMVA(fTMVAMethod[bin]);
1045   } else {
1046     mva = fTMVAReader[0]->EvaluateMVA(fTMVAMethod[0]);
1047   }
1048 
1049   if (printDebug) {
1050     std::cout << " *** Inside the class fMethodname " << fMethodname << " fMVAType " << fMVAType << std::endl;
1051     std::cout << " fbrem " << fMVAVar_fbrem << " kfchi2 " << fMVAVar_kfchi2 << " mykfhits " << fMVAVar_kfhits
1052               << " gsfchi2 " << fMVAVar_gsfchi2 << " deta " << fMVAVar_deta << " dphi " << fMVAVar_dphi << " detacalo "
1053               << fMVAVar_detacalo << " see " << fMVAVar_see << " spp " << fMVAVar_spp << " etawidth "
1054               << fMVAVar_etawidth << " phiwidth " << fMVAVar_phiwidth << " OneMinusE1x5E5x5 "
1055               << fMVAVar_OneMinusE1x5E5x5 << " R9 " << fMVAVar_R9 << " HoE " << fMVAVar_HoE << " EoP " << fMVAVar_EoP
1056               << " IoEmIoP " << fMVAVar_IoEmIoP << " eleEoPout " << fMVAVar_eleEoPout << " d0 " << fMVAVar_d0
1057               << " ip3d " << fMVAVar_ip3d << " eta " << fMVAVar_eta << " pt " << fMVAVar_pt << std::endl;
1058     std::cout << " ### MVA " << mva << std::endl;
1059   }
1060 
1061   return mva;
1062 }
1063 
1064 // for kTrigNoIP algorithm
1065 Double_t EGammaMvaEleEstimator::mvaValue(const reco::GsfElectron& ele,
1066                                          const reco::Vertex& vertex,
1067                                          double rho,
1068                                          //const TransientTrackBuilder& transientTrackBuilder,
1069                                          EcalClusterLazyTools const& myEcalCluster,
1070                                          bool printDebug) {
1071   if (!fisInitialized) {
1072     std::cout << "Error: EGammaMvaEleEstimator not properly initialized.\n";
1073     return -9999;
1074   }
1075 
1076   if (fMVAType != EGammaMvaEleEstimator::kTrigNoIP) {
1077     std::cout << "Error: This method should be called for kTrigNoIP MVA only" << std::endl;
1078     return -9999;
1079   }
1080 
1081   bool validKF = false;
1082   reco::TrackRef myTrackRef = ele.closestCtfTrackRef();
1083   validKF = (myTrackRef.isAvailable());
1084   validKF &= (myTrackRef.isNonnull());
1085 
1086   // Pure tracking variables
1087   fMVAVar_fbrem = ele.fbrem();
1088   fMVAVar_kfchi2 = (validKF) ? myTrackRef->normalizedChi2() : 0;
1089   fMVAVar_kfhits = (validKF) ? myTrackRef->hitPattern().trackerLayersWithMeasurement() : -1.;
1090   fMVAVar_gsfchi2 = ele.gsfTrack()->normalizedChi2();
1091 
1092   // Geometrical matchings
1093   fMVAVar_deta = ele.deltaEtaSuperClusterTrackAtVtx();
1094   fMVAVar_dphi = ele.deltaPhiSuperClusterTrackAtVtx();
1095   fMVAVar_detacalo = ele.deltaEtaSeedClusterTrackAtCalo();
1096 
1097   // Pure ECAL -> shower shapes
1098   fMVAVar_see = ele.sigmaIetaIeta();  //EleSigmaIEtaIEta
1099   const auto& vCov = myEcalCluster.localCovariances(*(ele.superCluster()->seed()));
1100   if (edm::isFinite(vCov[2]))
1101     fMVAVar_spp = sqrt(vCov[2]);  //EleSigmaIPhiIPhi
1102   else
1103     fMVAVar_spp = 0.;
1104 
1105   fMVAVar_etawidth = ele.superCluster()->etaWidth();
1106   fMVAVar_phiwidth = ele.superCluster()->phiWidth();
1107   fMVAVar_OneMinusE1x5E5x5 = (ele.e5x5()) != 0. ? 1. - (ele.e1x5() / ele.e5x5()) : -1.;
1108   fMVAVar_R9 = myEcalCluster.e3x3(*(ele.superCluster()->seed())) / ele.superCluster()->rawEnergy();
1109 
1110   // Energy matching
1111   fMVAVar_HoE = ele.hadronicOverEm();
1112   fMVAVar_EoP = ele.eSuperClusterOverP();
1113   fMVAVar_IoEmIoP = (1.0 / ele.superCluster()->energy()) -
1114                     (1.0 / ele.gsfTrack()->p());  // in the future to be changed with ele.gsfTrack()->p()
1115   fMVAVar_eleEoPout = ele.eEleClusterOverPout();
1116   fMVAVar_rho = rho;
1117   fMVAVar_PreShowerOverRaw = ele.superCluster()->preshowerEnergy() / ele.superCluster()->rawEnergy();
1118 
1119   // Spectators
1120   fMVAVar_eta = ele.superCluster()->eta();
1121   fMVAVar_pt = ele.pt();
1122 
1123   // evaluate
1124   bindVariables();
1125   Double_t mva = -9999;
1126   if (fUseBinnedVersion) {
1127     int bin = GetMVABin(fMVAVar_eta, fMVAVar_pt);
1128     mva = fTMVAReader[bin]->EvaluateMVA(fTMVAMethod[bin]);
1129   } else {
1130     mva = fTMVAReader[0]->EvaluateMVA(fTMVAMethod[0]);
1131   }
1132 
1133   if (printDebug) {
1134     std::cout << " *** Inside the class fMethodname " << fMethodname << " fMVAType " << fMVAType << std::endl;
1135     std::cout << " fbrem " << fMVAVar_fbrem << " kfchi2 " << fMVAVar_kfchi2 << " mykfhits " << fMVAVar_kfhits
1136               << " gsfchi2 " << fMVAVar_gsfchi2 << " deta " << fMVAVar_deta << " dphi " << fMVAVar_dphi << " detacalo "
1137               << fMVAVar_detacalo
1138               // << " dphicalo " << fMVAVar_dphicalo
1139               << " see " << fMVAVar_see << " spp " << fMVAVar_spp << " etawidth " << fMVAVar_etawidth << " phiwidth "
1140               << fMVAVar_phiwidth << " e1x5e5x5 " << fMVAVar_OneMinusE1x5E5x5 << " R9 "
1141               << fMVAVar_R9
1142               // << " mynbrems " << fMVAVar_nbrems
1143               << " HoE " << fMVAVar_HoE << " EoP " << fMVAVar_EoP << " IoEmIoP " << fMVAVar_IoEmIoP << " eleEoPout "
1144               << fMVAVar_eleEoPout << " rho "
1145               << fMVAVar_rho
1146               // << " EoPout " << fMVAVar_EoPout
1147               << " eta " << fMVAVar_eta << " pt " << fMVAVar_pt << std::endl;
1148     std::cout << " ### MVA " << mva << std::endl;
1149   }
1150 
1151   return mva;
1152 }
1153 
1154 // for kTrig, kNonTrig and kTrigNoIP algorithm
1155 Double_t EGammaMvaEleEstimator::mvaValue(
1156     const pat::Electron& ele, const reco::Vertex& vertex, double rho, bool useFull5x5, bool printDebug) {
1157   if (!fisInitialized) {
1158     std::cout << "Error: EGammaMvaEleEstimator not properly initialized.\n";
1159     return -9999;
1160   }
1161 
1162   bool validKF = false;
1163   reco::TrackRef myTrackRef = ele.closestCtfTrackRef();
1164   validKF = (myTrackRef.isAvailable());
1165   validKF &= (myTrackRef.isNonnull());
1166 
1167   // Pure tracking variables
1168   fMVAVar_fbrem = ele.fbrem();
1169   fMVAVar_kfchi2 = (validKF) ? myTrackRef->normalizedChi2() : 0;
1170   fMVAVar_kfhits = (validKF) ? myTrackRef->hitPattern().trackerLayersWithMeasurement() : -1.;
1171   fMVAVar_gsfchi2 = ele.gsfTrack()->normalizedChi2();
1172 
1173   // Geometrical matchings
1174   fMVAVar_deta = ele.deltaEtaSuperClusterTrackAtVtx();
1175   fMVAVar_dphi = ele.deltaPhiSuperClusterTrackAtVtx();
1176   fMVAVar_detacalo = ele.deltaEtaSeedClusterTrackAtCalo();
1177 
1178   // Pure ECAL -> shower shapes
1179   fMVAVar_see = useFull5x5 ? ele.full5x5_sigmaIetaIeta() : ele.sigmaIetaIeta();  //EleSigmaIEtaIEta
1180   fMVAVar_spp = useFull5x5 ? ele.full5x5_sigmaIphiIphi() : ele.sigmaIphiIphi();  //EleSigmaIEtaIEta
1181 
1182   fMVAVar_etawidth = ele.superCluster()->etaWidth();
1183   fMVAVar_phiwidth = ele.superCluster()->phiWidth();
1184   fMVAVar_OneMinusE1x5E5x5 = useFull5x5
1185                                  ? ((ele.full5x5_e5x5()) != 0. ? 1. - (ele.full5x5_e1x5() / ele.full5x5_e5x5()) : -1.)
1186                                  : ((ele.e5x5()) != 0. ? 1. - (ele.e1x5() / ele.e5x5()) : -1.);
1187   fMVAVar_R9 = useFull5x5 ? ele.full5x5_r9() : ele.r9();
1188 
1189   // Energy matching
1190   fMVAVar_HoE = ele.hadronicOverEm();  // this is identical for both
1191   fMVAVar_EoP = ele.eSuperClusterOverP();
1192 
1193   // unify this in the future?
1194   if (fMVAType == kTrig || fMVAType == kNonTrig) {
1195     fMVAVar_IoEmIoP =
1196         (1.0 / ele.ecalEnergy()) - (1.0 / ele.p());  // in the future to be changed with ele.gsfTrack()->p()
1197   } else {
1198     fMVAVar_IoEmIoP = (1.0 / ele.superCluster()->energy()) -
1199                       (1.0 / ele.gsfTrack()->p());  // in the future to be changed with ele.gsfTrack()->p()
1200   }
1201   fMVAVar_eleEoPout = ele.eEleClusterOverPout();
1202   fMVAVar_rho = rho;
1203   fMVAVar_PreShowerOverRaw = ele.superCluster()->preshowerEnergy() / ele.superCluster()->rawEnergy();
1204 
1205   // for triggering electrons get the impact parameteres
1206   if (fMVAType == kTrig) {
1207     //d0
1208     if (ele.gsfTrack().isNonnull()) {
1209       fMVAVar_d0 = (-1.0) * ele.gsfTrack()->dxy(vertex.position());
1210     } else if (ele.closestCtfTrackRef().isNonnull()) {
1211       fMVAVar_d0 = (-1.0) * ele.closestCtfTrackRef()->dxy(vertex.position());
1212     } else {
1213       fMVAVar_d0 = -9999.0;
1214     }
1215 
1216     //default values for IP3D
1217     fMVAVar_ip3d = -999.0;
1218     fMVAVar_ip3dSig = 0.0;
1219     if (ele.gsfTrack().isNonnull()) {
1220       const double gsfsign = ((-ele.gsfTrack()->dxy(vertex.position())) >= 0) ? 1. : -1.;
1221 
1222       //std::cout << "Warning : if usePV = false when pat-electrons were produced, dB() was calculated with respect to beamspot while original mva uses primary vertex" << std::endl;
1223       double ip3d = gsfsign * ele.dB();
1224       double ip3derr = ele.edB();
1225       fMVAVar_ip3d = ip3d;
1226       fMVAVar_ip3dSig = ip3d / ip3derr;
1227     }
1228   }
1229 
1230   // Spectators
1231   fMVAVar_eta = ele.superCluster()->eta();
1232   fMVAVar_pt = ele.pt();
1233 
1234   // evaluate
1235   bindVariables();
1236   Double_t mva = -9999;
1237   if (fUseBinnedVersion) {
1238     int bin = GetMVABin(fMVAVar_eta, fMVAVar_pt);
1239     mva = fTMVAReader[bin]->EvaluateMVA(fTMVAMethod[bin]);
1240   } else {
1241     mva = fTMVAReader[0]->EvaluateMVA(fTMVAMethod[0]);
1242   }
1243 
1244   if (printDebug) {
1245     std::cout << " *** Inside the class fMethodname " << fMethodname << " fMVAType " << fMVAType << std::endl;
1246     std::cout << " fbrem " << fMVAVar_fbrem << " kfchi2 " << fMVAVar_kfchi2 << " mykfhits " << fMVAVar_kfhits
1247               << " gsfchi2 " << fMVAVar_gsfchi2 << " deta " << fMVAVar_deta << " dphi " << fMVAVar_dphi << " detacalo "
1248               << fMVAVar_detacalo
1249               // << " dphicalo " << fMVAVar_dphicalo
1250               << " see " << fMVAVar_see << " spp " << fMVAVar_spp << " etawidth " << fMVAVar_etawidth << " phiwidth "
1251               << fMVAVar_phiwidth << " e1x5e5x5 " << fMVAVar_OneMinusE1x5E5x5 << " R9 "
1252               << fMVAVar_R9
1253               // << " mynbrems " << fMVAVar_nbrems
1254               << " HoE " << fMVAVar_HoE << " EoP " << fMVAVar_EoP << " IoEmIoP " << fMVAVar_IoEmIoP << " eleEoPout "
1255               << fMVAVar_eleEoPout << " rho "
1256               << fMVAVar_rho
1257               // << " EoPout " << fMVAVar_EoPout
1258               << " d0 " << fMVAVar_d0 << " ip3d " << fMVAVar_ip3d << " eta " << fMVAVar_eta << " pt " << fMVAVar_pt
1259               << std::endl;
1260     std::cout << " ### MVA " << mva << std::endl;
1261   }
1262 
1263   return mva;
1264 }
1265 
1266 // for kTrigNoIP algorithm only.
1267 Double_t EGammaMvaEleEstimator::mvaValue(const pat::Electron& ele, double rho, bool printDebug) {
1268   if (!fisInitialized) {
1269     std::cout << "Error: EGammaMvaEleEstimator not properly initialized.\n";
1270     return -9999;
1271   }
1272 
1273   if ((fMVAType != EGammaMvaEleEstimator::kTrigNoIP)) {
1274     std::cout << "Error: This method should be called for kTrigNoIP mva only" << std::endl;
1275     return -9999;
1276   }
1277 
1278   bool validKF = false;
1279   reco::TrackRef myTrackRef = ele.closestCtfTrackRef();
1280   validKF = (myTrackRef.isAvailable());
1281   validKF &= (myTrackRef.isNonnull());
1282 
1283   // Pure tracking variables
1284   fMVAVar_fbrem = ele.fbrem();
1285   fMVAVar_kfchi2 = (validKF) ? myTrackRef->normalizedChi2() : 0;
1286   fMVAVar_kfhits = (validKF) ? myTrackRef->hitPattern().trackerLayersWithMeasurement() : -1.;
1287   fMVAVar_gsfchi2 = ele.gsfTrack()->normalizedChi2();
1288 
1289   // Geometrical matchings
1290   fMVAVar_deta = ele.deltaEtaSuperClusterTrackAtVtx();
1291   fMVAVar_dphi = ele.deltaPhiSuperClusterTrackAtVtx();
1292   fMVAVar_detacalo = ele.deltaEtaSeedClusterTrackAtCalo();
1293 
1294   // Pure ECAL -> shower shapes
1295   fMVAVar_see = ele.sigmaIetaIeta();  //EleSigmaIEtaIEta
1296 
1297   fMVAVar_spp = ele.sigmaIphiIphi();
1298 
1299   fMVAVar_etawidth = ele.superCluster()->etaWidth();
1300   fMVAVar_phiwidth = ele.superCluster()->phiWidth();
1301   fMVAVar_OneMinusE1x5E5x5 = (ele.e5x5()) != 0. ? 1. - (ele.e1x5() / ele.e5x5()) : -1.;
1302   fMVAVar_R9 = ele.r9();
1303 
1304   // Energy matching
1305   fMVAVar_HoE = ele.hadronicOverEm();
1306   fMVAVar_EoP = ele.eSuperClusterOverP();
1307   fMVAVar_IoEmIoP = (1.0 / ele.superCluster()->energy()) -
1308                     (1.0 / ele.gsfTrack()->p());  // in the future to be changed with ele.gsfTrack()->p()
1309   fMVAVar_eleEoPout = ele.eEleClusterOverPout();
1310   fMVAVar_rho = rho;
1311   fMVAVar_PreShowerOverRaw = ele.superCluster()->preshowerEnergy() / ele.superCluster()->rawEnergy();
1312 
1313   // Spectators
1314   fMVAVar_eta = ele.superCluster()->eta();
1315   fMVAVar_pt = ele.pt();
1316 
1317   // evaluate
1318   bindVariables();
1319   Double_t mva = -9999;
1320   if (fUseBinnedVersion) {
1321     int bin = GetMVABin(fMVAVar_eta, fMVAVar_pt);
1322     mva = fTMVAReader[bin]->EvaluateMVA(fTMVAMethod[bin]);
1323   } else {
1324     mva = fTMVAReader[0]->EvaluateMVA(fTMVAMethod[0]);
1325   }
1326 
1327   if (printDebug) {
1328     std::cout << " *** Inside the class fMethodname " << fMethodname << " fMVAType " << fMVAType << std::endl;
1329     std::cout << " fbrem " << fMVAVar_fbrem << " kfchi2 " << fMVAVar_kfchi2 << " mykfhits " << fMVAVar_kfhits
1330               << " gsfchi2 " << fMVAVar_gsfchi2 << " deta " << fMVAVar_deta << " dphi " << fMVAVar_dphi << " detacalo "
1331               << fMVAVar_detacalo
1332               // << " dphicalo " << fMVAVar_dphicalo
1333               << " see " << fMVAVar_see << " spp " << fMVAVar_spp << " etawidth " << fMVAVar_etawidth << " phiwidth "
1334               << fMVAVar_phiwidth << " e1x5e5x5 " << fMVAVar_OneMinusE1x5E5x5 << " R9 "
1335               << fMVAVar_R9
1336               // << " mynbrems " << fMVAVar_nbrems
1337               << " HoE " << fMVAVar_HoE << " EoP " << fMVAVar_EoP << " IoEmIoP " << fMVAVar_IoEmIoP << " eleEoPout "
1338               << fMVAVar_eleEoPout << " rho "
1339               << fMVAVar_rho
1340               // << " EoPout " << fMVAVar_EoPout
1341               << " eta " << fMVAVar_eta << " pt " << fMVAVar_pt << std::endl;
1342     std::cout << " ### MVA " << mva << std::endl;
1343   }
1344 
1345   return mva;
1346 }
1347 
1348 Double_t EGammaMvaEleEstimator::isoMvaValue(const reco::GsfElectron& ele,
1349                                             const reco::Vertex& vertex,
1350                                             const reco::PFCandidateCollection& PFCandidates,
1351                                             double Rho,
1352                                             ElectronEffectiveArea::ElectronEffectiveAreaTarget EATarget,
1353                                             const reco::GsfElectronCollection& IdentifiedElectrons,
1354                                             const reco::MuonCollection& IdentifiedMuons,
1355                                             bool printDebug) {
1356   if (!fisInitialized) {
1357     std::cout << "Error: EGammaMvaEleEstimator not properly initialized.\n";
1358     return -9999;
1359   }
1360 
1361   // Spectators
1362   fMVAVar_eta = ele.superCluster()->eta();
1363   fMVAVar_pt = ele.pt();
1364 
1365   //**********************************************************
1366   //Isolation variables
1367   //**********************************************************
1368   Double_t tmpChargedIso_DR0p0To0p1 = 0;
1369   Double_t tmpChargedIso_DR0p1To0p2 = 0;
1370   Double_t tmpChargedIso_DR0p2To0p3 = 0;
1371   Double_t tmpChargedIso_DR0p3To0p4 = 0;
1372   Double_t tmpChargedIso_DR0p4To0p5 = 0;
1373   Double_t tmpGammaIso_DR0p0To0p1 = 0;
1374   Double_t tmpGammaIso_DR0p1To0p2 = 0;
1375   Double_t tmpGammaIso_DR0p2To0p3 = 0;
1376   Double_t tmpGammaIso_DR0p3To0p4 = 0;
1377   Double_t tmpGammaIso_DR0p4To0p5 = 0;
1378   Double_t tmpNeutralHadronIso_DR0p0To0p1 = 0;
1379   Double_t tmpNeutralHadronIso_DR0p1To0p2 = 0;
1380   Double_t tmpNeutralHadronIso_DR0p2To0p3 = 0;
1381   Double_t tmpNeutralHadronIso_DR0p3To0p4 = 0;
1382   Double_t tmpNeutralHadronIso_DR0p4To0p5 = 0;
1383 
1384   double electronTrackZ = 0;
1385   if (ele.gsfTrack().isNonnull()) {
1386     electronTrackZ = ele.gsfTrack()->dz(vertex.position());
1387   } else if (ele.closestCtfTrackRef().isNonnull()) {
1388     electronTrackZ = ele.closestCtfTrackRef()->dz(vertex.position());
1389   }
1390 
1391   for (reco::PFCandidateCollection::const_iterator iP = PFCandidates.begin(); iP != PFCandidates.end(); ++iP) {
1392     //exclude the electron itself
1393     if (iP->gsfTrackRef().isNonnull() && ele.gsfTrack().isNonnull() &&
1394         refToPtr(iP->gsfTrackRef()) == refToPtr(ele.gsfTrack()))
1395       continue;
1396     if (iP->trackRef().isNonnull() && ele.closestCtfTrackRef().isNonnull() &&
1397         refToPtr(iP->trackRef()) == refToPtr(ele.closestCtfTrackRef()))
1398       continue;
1399 
1400     //************************************************************
1401     // New Isolation Calculations
1402     //************************************************************
1403     double dr = sqrt(pow(iP->eta() - ele.eta(), 2) + pow(acos(cos(iP->phi() - ele.phi())), 2));
1404     //Double_t deta = (iP->eta() - ele.eta());
1405 
1406     if (dr < 1.0) {
1407       Bool_t IsLeptonFootprint = kFALSE;
1408       //************************************************************
1409       // Lepton Footprint Removal
1410       //************************************************************
1411       for (reco::GsfElectronCollection::const_iterator iE = IdentifiedElectrons.begin();
1412            iE != IdentifiedElectrons.end();
1413            ++iE) {
1414         //if pf candidate matches an electron passing ID cuts, then veto it
1415         if (iP->gsfTrackRef().isNonnull() && iE->gsfTrack().isNonnull() &&
1416             refToPtr(iP->gsfTrackRef()) == refToPtr(iE->gsfTrack()))
1417           IsLeptonFootprint = kTRUE;
1418         if (iP->trackRef().isNonnull() && iE->closestCtfTrackRef().isNonnull() &&
1419             refToPtr(iP->trackRef()) == refToPtr(iE->closestCtfTrackRef()))
1420           IsLeptonFootprint = kTRUE;
1421 
1422         //if pf candidate lies in veto regions of electron passing ID cuts, then veto it
1423         double tmpDR = sqrt(pow(iP->eta() - iE->eta(), 2) + pow(acos(cos(iP->phi() - iE->phi())), 2));
1424         if (iP->trackRef().isNonnull() && fabs(iE->superCluster()->eta()) >= 1.479 && tmpDR < 0.015)
1425           IsLeptonFootprint = kTRUE;
1426         if (iP->particleId() == reco::PFCandidate::gamma && fabs(iE->superCluster()->eta()) >= 1.479 && tmpDR < 0.08)
1427           IsLeptonFootprint = kTRUE;
1428       }
1429       for (reco::MuonCollection::const_iterator iM = IdentifiedMuons.begin(); iM != IdentifiedMuons.end(); ++iM) {
1430         //if pf candidate matches an muon passing ID cuts, then veto it
1431         if (iP->trackRef().isNonnull() && iM->innerTrack().isNonnull() &&
1432             refToPtr(iP->trackRef()) == refToPtr(iM->innerTrack()))
1433           IsLeptonFootprint = kTRUE;
1434 
1435         //if pf candidate lies in veto regions of muon passing ID cuts, then veto it
1436         double tmpDR = sqrt(pow(iP->eta() - iM->eta(), 2) + pow(acos(cos(iP->phi() - iM->phi())), 2));
1437         if (iP->trackRef().isNonnull() && tmpDR < 0.01)
1438           IsLeptonFootprint = kTRUE;
1439       }
1440 
1441       if (!IsLeptonFootprint) {
1442         Bool_t passVeto = kTRUE;
1443         //Charged
1444         if (iP->trackRef().isNonnull()) {
1445           if (!(fabs(iP->trackRef()->dz(vertex.position()) - electronTrackZ) < 0.2))
1446             passVeto = kFALSE;
1447           //************************************************************
1448           // Veto any PFmuon, or PFEle
1449           if (iP->particleId() == reco::PFCandidate::e || iP->particleId() == reco::PFCandidate::mu)
1450             passVeto = kFALSE;
1451           //************************************************************
1452           //************************************************************
1453           // Footprint Veto
1454           if (fabs(fMVAVar_eta) > 1.479 && dr < 0.015)
1455             passVeto = kFALSE;
1456           //************************************************************
1457           if (passVeto) {
1458             if (dr < 0.1)
1459               tmpChargedIso_DR0p0To0p1 += iP->pt();
1460             if (dr >= 0.1 && dr < 0.2)
1461               tmpChargedIso_DR0p1To0p2 += iP->pt();
1462             if (dr >= 0.2 && dr < 0.3)
1463               tmpChargedIso_DR0p2To0p3 += iP->pt();
1464             if (dr >= 0.3 && dr < 0.4)
1465               tmpChargedIso_DR0p3To0p4 += iP->pt();
1466             if (dr >= 0.4 && dr < 0.5)
1467               tmpChargedIso_DR0p4To0p5 += iP->pt();
1468           }  //pass veto
1469         }
1470         //Gamma
1471         else if (iP->particleId() == reco::PFCandidate::gamma) {
1472           //************************************************************
1473           // Footprint Veto
1474           if (fabs(fMVAVar_eta) > 1.479 && dr < 0.08)
1475             passVeto = kFALSE;
1476           //************************************************************
1477           if (passVeto) {
1478             if (dr < 0.1)
1479               tmpGammaIso_DR0p0To0p1 += iP->pt();
1480             if (dr >= 0.1 && dr < 0.2)
1481               tmpGammaIso_DR0p1To0p2 += iP->pt();
1482             if (dr >= 0.2 && dr < 0.3)
1483               tmpGammaIso_DR0p2To0p3 += iP->pt();
1484             if (dr >= 0.3 && dr < 0.4)
1485               tmpGammaIso_DR0p3To0p4 += iP->pt();
1486             if (dr >= 0.4 && dr < 0.5)
1487               tmpGammaIso_DR0p4To0p5 += iP->pt();
1488           }
1489         }
1490         //NeutralHadron
1491         else {
1492           if (dr < 0.1)
1493             tmpNeutralHadronIso_DR0p0To0p1 += iP->pt();
1494           if (dr >= 0.1 && dr < 0.2)
1495             tmpNeutralHadronIso_DR0p1To0p2 += iP->pt();
1496           if (dr >= 0.2 && dr < 0.3)
1497             tmpNeutralHadronIso_DR0p2To0p3 += iP->pt();
1498           if (dr >= 0.3 && dr < 0.4)
1499             tmpNeutralHadronIso_DR0p3To0p4 += iP->pt();
1500           if (dr >= 0.4 && dr < 0.5)
1501             tmpNeutralHadronIso_DR0p4To0p5 += iP->pt();
1502         }
1503       }  //not lepton footprint
1504     }    //in 1.0 dr cone
1505   }      //loop over PF candidates
1506 
1507   fMVAVar_ChargedIso_DR0p0To0p1 = TMath::Min((tmpChargedIso_DR0p0To0p1) / ele.pt(), 2.5);
1508   fMVAVar_ChargedIso_DR0p1To0p2 = TMath::Min((tmpChargedIso_DR0p1To0p2) / ele.pt(), 2.5);
1509   fMVAVar_ChargedIso_DR0p2To0p3 = TMath::Min((tmpChargedIso_DR0p2To0p3) / ele.pt(), 2.5);
1510   fMVAVar_ChargedIso_DR0p3To0p4 = TMath::Min((tmpChargedIso_DR0p3To0p4) / ele.pt(), 2.5);
1511   fMVAVar_ChargedIso_DR0p4To0p5 = TMath::Min((tmpChargedIso_DR0p4To0p5) / ele.pt(), 2.5);
1512   fMVAVar_GammaIso_DR0p0To0p1 = TMath::Max(
1513       TMath::Min(
1514           (tmpGammaIso_DR0p0To0p1 - Rho * ElectronEffectiveArea::GetElectronEffectiveArea(
1515                                               ElectronEffectiveArea::kEleGammaIsoDR0p0To0p1, fMVAVar_eta, EATarget)) /
1516               ele.pt(),
1517           2.5),
1518       0.0);
1519   fMVAVar_GammaIso_DR0p1To0p2 = TMath::Max(
1520       TMath::Min(
1521           (tmpGammaIso_DR0p1To0p2 - Rho * ElectronEffectiveArea::GetElectronEffectiveArea(
1522                                               ElectronEffectiveArea::kEleGammaIsoDR0p1To0p2, fMVAVar_eta, EATarget)) /
1523               ele.pt(),
1524           2.5),
1525       0.0);
1526   fMVAVar_GammaIso_DR0p2To0p3 = TMath::Max(
1527       TMath::Min(
1528           (tmpGammaIso_DR0p2To0p3 - Rho * ElectronEffectiveArea::GetElectronEffectiveArea(
1529                                               ElectronEffectiveArea::kEleGammaIsoDR0p2To0p3, fMVAVar_eta, EATarget)) /
1530               ele.pt(),
1531           2.5),
1532       0.0);
1533   fMVAVar_GammaIso_DR0p3To0p4 = TMath::Max(
1534       TMath::Min(
1535           (tmpGammaIso_DR0p3To0p4 - Rho * ElectronEffectiveArea::GetElectronEffectiveArea(
1536                                               ElectronEffectiveArea::kEleGammaIsoDR0p3To0p4, fMVAVar_eta, EATarget)) /
1537               ele.pt(),
1538           2.5),
1539       0.0);
1540   fMVAVar_GammaIso_DR0p4To0p5 = TMath::Max(
1541       TMath::Min(
1542           (tmpGammaIso_DR0p4To0p5 - Rho * ElectronEffectiveArea::GetElectronEffectiveArea(
1543                                               ElectronEffectiveArea::kEleGammaIsoDR0p4To0p5, fMVAVar_eta, EATarget)) /
1544               ele.pt(),
1545           2.5),
1546       0.0);
1547   fMVAVar_NeutralHadronIso_DR0p0To0p1 =
1548       TMath::Max(TMath::Min((tmpNeutralHadronIso_DR0p0To0p1 -
1549                              Rho * ElectronEffectiveArea::GetElectronEffectiveArea(
1550                                        ElectronEffectiveArea::kEleNeutralHadronIsoDR0p0To0p1, fMVAVar_eta, EATarget)) /
1551                                 ele.pt(),
1552                             2.5),
1553                  0.0);
1554   fMVAVar_NeutralHadronIso_DR0p1To0p2 =
1555       TMath::Max(TMath::Min((tmpNeutralHadronIso_DR0p1To0p2 -
1556                              Rho * ElectronEffectiveArea::GetElectronEffectiveArea(
1557                                        ElectronEffectiveArea::kEleNeutralHadronIsoDR0p1To0p2, fMVAVar_eta, EATarget)) /
1558                                 ele.pt(),
1559                             2.5),
1560                  0.0);
1561   fMVAVar_NeutralHadronIso_DR0p2To0p3 =
1562       TMath::Max(TMath::Min((tmpNeutralHadronIso_DR0p2To0p3 -
1563                              Rho * ElectronEffectiveArea::GetElectronEffectiveArea(
1564                                        ElectronEffectiveArea::kEleNeutralHadronIsoDR0p2To0p3, fMVAVar_eta, EATarget)) /
1565                                 ele.pt(),
1566                             2.5),
1567                  0.0);
1568   fMVAVar_NeutralHadronIso_DR0p3To0p4 =
1569       TMath::Max(TMath::Min((tmpNeutralHadronIso_DR0p3To0p4 -
1570                              Rho * ElectronEffectiveArea::GetElectronEffectiveArea(
1571                                        ElectronEffectiveArea::kEleNeutralHadronIsoDR0p3To0p4, fMVAVar_eta, EATarget)) /
1572                                 ele.pt(),
1573                             2.5),
1574                  0.0);
1575   fMVAVar_NeutralHadronIso_DR0p4To0p5 =
1576       TMath::Max(TMath::Min((tmpNeutralHadronIso_DR0p4To0p5 -
1577                              Rho * ElectronEffectiveArea::GetElectronEffectiveArea(
1578                                        ElectronEffectiveArea::kEleNeutralHadronIsoDR0p4To0p5, fMVAVar_eta, EATarget)) /
1579                                 ele.pt(),
1580                             2.5),
1581                  0.0);
1582 
1583   if (printDebug) {
1584     std::cout << "UseBinnedVersion=" << fUseBinnedVersion << " -> BIN: " << fMVAVar_eta << " " << fMVAVar_pt << " : "
1585               << GetMVABin(fMVAVar_eta, fMVAVar_pt) << std::endl;
1586   }
1587 
1588   // evaluate
1589   bindVariables();
1590   Double_t mva = -9999;
1591 
1592   //   mva = fTMVAReader[0]->EvaluateMVA(fMethodname);
1593   if (fUseBinnedVersion) {
1594     int bin = GetMVABin(fMVAVar_eta, fMVAVar_pt);
1595     mva = fTMVAReader[bin]->EvaluateMVA(fTMVAMethod[bin]);
1596   } else {
1597     mva = fTMVAReader[0]->EvaluateMVA(fTMVAMethod[0]);
1598   }
1599 
1600   if (printDebug) {
1601     std::cout << " *** Inside the class fMethodname " << fMethodname << " fMVAType " << fMVAType << std::endl;
1602     std::cout << "ChargedIso ( 0.0 | 0.1 | 0.2 | 0.3 | 0.4 | 0.5 ): " << fMVAVar_ChargedIso_DR0p0To0p1 << " "
1603               << fMVAVar_ChargedIso_DR0p1To0p2 << " " << fMVAVar_ChargedIso_DR0p2To0p3 << " "
1604               << fMVAVar_ChargedIso_DR0p3To0p4 << " " << fMVAVar_ChargedIso_DR0p4To0p5 << std::endl;
1605     std::cout << "PF Gamma Iso ( 0.0 | 0.1 | 0.2 | 0.3 | 0.4 | 0.5 ): " << fMVAVar_GammaIso_DR0p0To0p1 << " "
1606               << fMVAVar_GammaIso_DR0p1To0p2 << " " << fMVAVar_GammaIso_DR0p2To0p3 << " " << fMVAVar_GammaIso_DR0p3To0p4
1607               << " " << fMVAVar_GammaIso_DR0p4To0p5 << std::endl;
1608     std::cout << "PF Neutral Hadron Iso ( 0.0 | 0.1 | 0.2 | 0.3 | 0.4 | 0.5 ): " << fMVAVar_NeutralHadronIso_DR0p0To0p1
1609               << " " << fMVAVar_NeutralHadronIso_DR0p1To0p2 << " " << fMVAVar_NeutralHadronIso_DR0p2To0p3 << " "
1610               << fMVAVar_NeutralHadronIso_DR0p3To0p4 << " " << fMVAVar_NeutralHadronIso_DR0p4To0p5 << " " << std::endl;
1611     std::cout << " ### MVA " << mva << std::endl;
1612   }
1613 
1614   return mva;
1615 }
1616 
1617 //--------------------------------------------------------------------------------------------------
1618 
1619 Double_t EGammaMvaEleEstimator::IDIsoCombinedMvaValue(const reco::GsfElectron& ele,
1620                                                       const reco::Vertex& vertex,
1621                                                       const TransientTrackBuilder& transientTrackBuilder,
1622                                                       EcalClusterLazyTools const& myEcalCluster,
1623                                                       const reco::PFCandidateCollection& PFCandidates,
1624                                                       double Rho,
1625                                                       ElectronEffectiveArea::ElectronEffectiveAreaTarget EATarget,
1626                                                       bool printDebug) {
1627   if (!fisInitialized) {
1628     std::cout << "Error: EGammaMvaEleEstimator not properly initialized.\n";
1629     return -9999;
1630   }
1631 
1632   bool validKF = false;
1633   reco::TrackRef myTrackRef = ele.closestCtfTrackRef();
1634   validKF = (myTrackRef.isAvailable());
1635   validKF &= (myTrackRef.isNonnull());
1636 
1637   // Pure tracking variables
1638   fMVAVar_fbrem = (ele.fbrem() < -1.) ? -1. : ele.fbrem();
1639   fMVAVar_kfchi2 = (validKF) ? myTrackRef->normalizedChi2() : 0;
1640   if (fMVAVar_kfchi2 > 10)
1641     fMVAVar_kfchi2 = 10;
1642   fMVAVar_kfhits = (validKF) ? myTrackRef->hitPattern().trackerLayersWithMeasurement() : -1.;
1643   fMVAVar_kfhitsall =
1644       (validKF) ? myTrackRef->numberOfValidHits() : -1.;  //  save also this in your ntuple as possible alternative
1645   fMVAVar_gsfchi2 = ele.gsfTrack()->normalizedChi2();
1646   if (fMVAVar_gsfchi2 > 200)
1647     fMVAVar_gsfchi2 = 200;
1648 
1649   // Geometrical matchings
1650   fMVAVar_deta =
1651       (fabs(ele.deltaEtaSuperClusterTrackAtVtx()) > 0.06) ? 0.06 : fabs(ele.deltaEtaSuperClusterTrackAtVtx());
1652   fMVAVar_dphi = ele.deltaPhiSuperClusterTrackAtVtx();
1653   fMVAVar_detacalo = ele.deltaEtaSeedClusterTrackAtCalo();
1654 
1655   // Pure ECAL -> shower shapes
1656   fMVAVar_see = ele.sigmaIetaIeta();  //EleSigmaIEtaIEta
1657   const auto& vCov = myEcalCluster.localCovariances(*(ele.superCluster()->seed()));
1658   if (edm::isFinite(vCov[2]))
1659     fMVAVar_spp = sqrt(vCov[2]);  //EleSigmaIPhiIPhi
1660   else
1661     fMVAVar_spp = 0.;
1662 
1663   fMVAVar_etawidth = ele.superCluster()->etaWidth();
1664   fMVAVar_phiwidth = ele.superCluster()->phiWidth();
1665   fMVAVar_OneMinusE1x5E5x5 = (ele.e5x5()) != 0. ? 1. - (ele.e1x5() / ele.e5x5()) : -1.;
1666   fMVAVar_OneMinusE1x5E5x5 = std::max(std::min(double(fMVAVar_OneMinusE1x5E5x5), 2.0), -1.0);
1667   fMVAVar_R9 = myEcalCluster.e3x3(*(ele.superCluster()->seed())) / ele.superCluster()->rawEnergy();
1668   if (fMVAVar_R9 > 5)
1669     fMVAVar_R9 = 5;
1670 
1671   // Energy matching
1672   fMVAVar_HoE = ele.hadronicOverEm();
1673   fMVAVar_EoP = (ele.eSuperClusterOverP() > 20) ? 20 : ele.eSuperClusterOverP();
1674   fMVAVar_IoEmIoP =
1675       (1.0 / ele.superCluster()->energy()) - (1.0 / ele.trackMomentumAtVtx().R());  //this is the proper variable
1676   fMVAVar_eleEoPout = (ele.eEleClusterOverPout() > 20) ? 20 : ele.eEleClusterOverPout();
1677   fMVAVar_PreShowerOverRaw = ele.superCluster()->preshowerEnergy() / ele.superCluster()->rawEnergy();
1678 
1679   // Spectators
1680   fMVAVar_eta = ele.superCluster()->eta();
1681   fMVAVar_pt = ele.pt();
1682 
1683   // for triggering electrons get the impact parameteres
1684   if (fMVAType == kTrig) {
1685     //d0
1686     if (ele.gsfTrack().isNonnull()) {
1687       fMVAVar_d0 = (-1.0) * ele.gsfTrack()->dxy(vertex.position());
1688     } else if (ele.closestCtfTrackRef().isNonnull()) {
1689       fMVAVar_d0 = (-1.0) * ele.closestCtfTrackRef()->dxy(vertex.position());
1690     } else {
1691       fMVAVar_d0 = -9999.0;
1692     }
1693 
1694     //default values for IP3D
1695     fMVAVar_ip3d = -999.0;
1696     fMVAVar_ip3dSig = 0.0;
1697     if (ele.gsfTrack().isNonnull()) {
1698       const double gsfsign = ((-ele.gsfTrack()->dxy(vertex.position())) >= 0) ? 1. : -1.;
1699 
1700       const reco::TransientTrack& tt = transientTrackBuilder.build(ele.gsfTrack());
1701       const std::pair<bool, Measurement1D>& ip3dpv = IPTools::absoluteImpactParameter3D(tt, vertex);
1702       if (ip3dpv.first) {
1703         double ip3d = gsfsign * ip3dpv.second.value();
1704         double ip3derr = ip3dpv.second.error();
1705         fMVAVar_ip3d = ip3d;
1706         fMVAVar_ip3dSig = ip3d / ip3derr;
1707       }
1708     }
1709   }
1710 
1711   //**********************************************************
1712   //Isolation variables
1713   //**********************************************************
1714   Double_t tmpChargedIso_DR0p0To0p1 = 0;
1715   Double_t tmpChargedIso_DR0p1To0p2 = 0;
1716   Double_t tmpChargedIso_DR0p2To0p3 = 0;
1717   Double_t tmpChargedIso_DR0p3To0p4 = 0;
1718   Double_t tmpChargedIso_DR0p4To0p5 = 0;
1719   Double_t tmpGammaIso_DR0p0To0p1 = 0;
1720   Double_t tmpGammaIso_DR0p1To0p2 = 0;
1721   Double_t tmpGammaIso_DR0p2To0p3 = 0;
1722   Double_t tmpGammaIso_DR0p3To0p4 = 0;
1723   Double_t tmpGammaIso_DR0p4To0p5 = 0;
1724   Double_t tmpNeutralHadronIso_DR0p0To0p1 = 0;
1725   Double_t tmpNeutralHadronIso_DR0p1To0p2 = 0;
1726   Double_t tmpNeutralHadronIso_DR0p2To0p3 = 0;
1727   Double_t tmpNeutralHadronIso_DR0p3To0p4 = 0;
1728   Double_t tmpNeutralHadronIso_DR0p4To0p5 = 0;
1729 
1730   for (reco::PFCandidateCollection::const_iterator iP = PFCandidates.begin(); iP != PFCandidates.end(); ++iP) {
1731     double dr = sqrt(pow(iP->eta() - ele.eta(), 2) + pow(acos(cos(iP->phi() - ele.phi())), 2));
1732 
1733     Bool_t passVeto = kTRUE;
1734     //Charged
1735     if (iP->trackRef().isNonnull()) {
1736       //make sure charged pf candidates pass the PFNoPU condition (assumed)
1737 
1738       //************************************************************
1739       // Veto any PFmuon, or PFEle
1740       if (iP->particleId() == reco::PFCandidate::e || iP->particleId() == reco::PFCandidate::mu)
1741         passVeto = kFALSE;
1742       //************************************************************
1743       //************************************************************
1744       // Footprint Veto
1745       if (fabs(fMVAVar_eta) > 1.479 && dr < 0.015)
1746         passVeto = kFALSE;
1747       //************************************************************
1748       if (passVeto) {
1749         if (dr < 0.1)
1750           tmpChargedIso_DR0p0To0p1 += iP->pt();
1751         if (dr >= 0.1 && dr < 0.2)
1752           tmpChargedIso_DR0p1To0p2 += iP->pt();
1753         if (dr >= 0.2 && dr < 0.3)
1754           tmpChargedIso_DR0p2To0p3 += iP->pt();
1755         if (dr >= 0.3 && dr < 0.4)
1756           tmpChargedIso_DR0p3To0p4 += iP->pt();
1757         if (dr >= 0.4 && dr < 0.5)
1758           tmpChargedIso_DR0p4To0p5 += iP->pt();
1759       }  //pass veto
1760     }
1761     //Gamma
1762     else if (iP->particleId() == reco::PFCandidate::gamma) {
1763       //************************************************************
1764       // Footprint Veto
1765       if (fabs(fMVAVar_eta) > 1.479 && dr < 0.08)
1766         passVeto = kFALSE;
1767       //************************************************************
1768       if (passVeto) {
1769         if (dr < 0.1)
1770           tmpGammaIso_DR0p0To0p1 += iP->pt();
1771         if (dr >= 0.1 && dr < 0.2)
1772           tmpGammaIso_DR0p1To0p2 += iP->pt();
1773         if (dr >= 0.2 && dr < 0.3)
1774           tmpGammaIso_DR0p2To0p3 += iP->pt();
1775         if (dr >= 0.3 && dr < 0.4)
1776           tmpGammaIso_DR0p3To0p4 += iP->pt();
1777         if (dr >= 0.4 && dr < 0.5)
1778           tmpGammaIso_DR0p4To0p5 += iP->pt();
1779       }
1780     }
1781     //NeutralHadron
1782     else {
1783       if (dr < 0.1)
1784         tmpNeutralHadronIso_DR0p0To0p1 += iP->pt();
1785       if (dr >= 0.1 && dr < 0.2)
1786         tmpNeutralHadronIso_DR0p1To0p2 += iP->pt();
1787       if (dr >= 0.2 && dr < 0.3)
1788         tmpNeutralHadronIso_DR0p2To0p3 += iP->pt();
1789       if (dr >= 0.3 && dr < 0.4)
1790         tmpNeutralHadronIso_DR0p3To0p4 += iP->pt();
1791       if (dr >= 0.4 && dr < 0.5)
1792         tmpNeutralHadronIso_DR0p4To0p5 += iP->pt();
1793     }
1794   }  //loop over PF candidates
1795 
1796   if (fMVAType == kTrigIDIsoCombinedPUCorrected) {
1797     fMVAVar_ChargedIso_DR0p0To0p1 = TMath::Min((tmpChargedIso_DR0p0To0p1) / ele.pt(), 2.5);
1798     fMVAVar_ChargedIso_DR0p1To0p2 = TMath::Min((tmpChargedIso_DR0p1To0p2) / ele.pt(), 2.5);
1799     fMVAVar_ChargedIso_DR0p2To0p3 = TMath::Min((tmpChargedIso_DR0p2To0p3) / ele.pt(), 2.5);
1800     fMVAVar_ChargedIso_DR0p3To0p4 = TMath::Min((tmpChargedIso_DR0p3To0p4) / ele.pt(), 2.5);
1801     fMVAVar_ChargedIso_DR0p4To0p5 = TMath::Min((tmpChargedIso_DR0p4To0p5) / ele.pt(), 2.5);
1802     fMVAVar_GammaIso_DR0p0To0p1 = TMath::Max(
1803         TMath::Min(
1804             (tmpGammaIso_DR0p0To0p1 - Rho * ElectronEffectiveArea::GetElectronEffectiveArea(
1805                                                 ElectronEffectiveArea::kEleGammaIsoDR0p0To0p1, fMVAVar_eta, EATarget)) /
1806                 ele.pt(),
1807             2.5),
1808         0.0);
1809     fMVAVar_GammaIso_DR0p1To0p2 = TMath::Max(
1810         TMath::Min(
1811             (tmpGammaIso_DR0p1To0p2 - Rho * ElectronEffectiveArea::GetElectronEffectiveArea(
1812                                                 ElectronEffectiveArea::kEleGammaIsoDR0p1To0p2, fMVAVar_eta, EATarget)) /
1813                 ele.pt(),
1814             2.5),
1815         0.0);
1816     fMVAVar_GammaIso_DR0p2To0p3 = TMath::Max(
1817         TMath::Min(
1818             (tmpGammaIso_DR0p2To0p3 - Rho * ElectronEffectiveArea::GetElectronEffectiveArea(
1819                                                 ElectronEffectiveArea::kEleGammaIsoDR0p2To0p3, fMVAVar_eta, EATarget)) /
1820                 ele.pt(),
1821             2.5),
1822         0.0);
1823     fMVAVar_GammaIso_DR0p3To0p4 = TMath::Max(
1824         TMath::Min(
1825             (tmpGammaIso_DR0p3To0p4 - Rho * ElectronEffectiveArea::GetElectronEffectiveArea(
1826                                                 ElectronEffectiveArea::kEleGammaIsoDR0p3To0p4, fMVAVar_eta, EATarget)) /
1827                 ele.pt(),
1828             2.5),
1829         0.0);
1830     fMVAVar_GammaIso_DR0p4To0p5 = TMath::Max(
1831         TMath::Min(
1832             (tmpGammaIso_DR0p4To0p5 - Rho * ElectronEffectiveArea::GetElectronEffectiveArea(
1833                                                 ElectronEffectiveArea::kEleGammaIsoDR0p4To0p5, fMVAVar_eta, EATarget)) /
1834                 ele.pt(),
1835             2.5),
1836         0.0);
1837     fMVAVar_NeutralHadronIso_DR0p0To0p1 = TMath::Max(
1838         TMath::Min((tmpNeutralHadronIso_DR0p0To0p1 -
1839                     Rho * ElectronEffectiveArea::GetElectronEffectiveArea(
1840                               ElectronEffectiveArea::kEleNeutralHadronIsoDR0p0To0p1, fMVAVar_eta, EATarget)) /
1841                        ele.pt(),
1842                    2.5),
1843         0.0);
1844     fMVAVar_NeutralHadronIso_DR0p1To0p2 = TMath::Max(
1845         TMath::Min((tmpNeutralHadronIso_DR0p1To0p2 -
1846                     Rho * ElectronEffectiveArea::GetElectronEffectiveArea(
1847                               ElectronEffectiveArea::kEleNeutralHadronIsoDR0p1To0p2, fMVAVar_eta, EATarget)) /
1848                        ele.pt(),
1849                    2.5),
1850         0.0);
1851     fMVAVar_NeutralHadronIso_DR0p2To0p3 = TMath::Max(
1852         TMath::Min((tmpNeutralHadronIso_DR0p2To0p3 -
1853                     Rho * ElectronEffectiveArea::GetElectronEffectiveArea(
1854                               ElectronEffectiveArea::kEleNeutralHadronIsoDR0p2To0p3, fMVAVar_eta, EATarget)) /
1855                        ele.pt(),
1856                    2.5),
1857         0.0);
1858     fMVAVar_NeutralHadronIso_DR0p3To0p4 = TMath::Max(
1859         TMath::Min((tmpNeutralHadronIso_DR0p3To0p4 -
1860                     Rho * ElectronEffectiveArea::GetElectronEffectiveArea(
1861                               ElectronEffectiveArea::kEleNeutralHadronIsoDR0p3To0p4, fMVAVar_eta, EATarget)) /
1862                        ele.pt(),
1863                    2.5),
1864         0.0);
1865     fMVAVar_NeutralHadronIso_DR0p4To0p5 = TMath::Max(
1866         TMath::Min((tmpNeutralHadronIso_DR0p4To0p5 -
1867                     Rho * ElectronEffectiveArea::GetElectronEffectiveArea(
1868                               ElectronEffectiveArea::kEleNeutralHadronIsoDR0p4To0p5, fMVAVar_eta, EATarget)) /
1869                        ele.pt(),
1870                    2.5),
1871         0.0);
1872   } else if (fMVAType == kTrigIDIsoCombined) {
1873     fMVAVar_ChargedIso_DR0p0To0p1 = TMath::Min((tmpChargedIso_DR0p0To0p1) / ele.pt(), 2.5);
1874     fMVAVar_ChargedIso_DR0p1To0p2 = TMath::Min((tmpChargedIso_DR0p1To0p2) / ele.pt(), 2.5);
1875     fMVAVar_ChargedIso_DR0p2To0p3 = TMath::Min((tmpChargedIso_DR0p2To0p3) / ele.pt(), 2.5);
1876     fMVAVar_ChargedIso_DR0p3To0p4 = TMath::Min((tmpChargedIso_DR0p3To0p4) / ele.pt(), 2.5);
1877     fMVAVar_ChargedIso_DR0p4To0p5 = TMath::Min((tmpChargedIso_DR0p4To0p5) / ele.pt(), 2.5);
1878     fMVAVar_GammaIso_DR0p0To0p1 = TMath::Max(TMath::Min((tmpGammaIso_DR0p0To0p1) / ele.pt(), 2.5), 0.0);
1879     fMVAVar_GammaIso_DR0p1To0p2 = TMath::Max(TMath::Min((tmpGammaIso_DR0p1To0p2) / ele.pt(), 2.5), 0.0);
1880     fMVAVar_GammaIso_DR0p2To0p3 = TMath::Max(TMath::Min((tmpGammaIso_DR0p2To0p3) / ele.pt(), 2.5), 0.0);
1881     fMVAVar_GammaIso_DR0p3To0p4 = TMath::Max(TMath::Min((tmpGammaIso_DR0p3To0p4) / ele.pt(), 2.5), 0.0);
1882     fMVAVar_GammaIso_DR0p4To0p5 = TMath::Max(TMath::Min((tmpGammaIso_DR0p4To0p5) / ele.pt(), 2.5), 0.0);
1883     fMVAVar_NeutralHadronIso_DR0p0To0p1 = TMath::Max(TMath::Min((tmpNeutralHadronIso_DR0p0To0p1) / ele.pt(), 2.5), 0.0);
1884     fMVAVar_NeutralHadronIso_DR0p1To0p2 = TMath::Max(TMath::Min((tmpNeutralHadronIso_DR0p1To0p2) / ele.pt(), 2.5), 0.0);
1885     fMVAVar_NeutralHadronIso_DR0p2To0p3 = TMath::Max(TMath::Min((tmpNeutralHadronIso_DR0p2To0p3) / ele.pt(), 2.5), 0.0);
1886     fMVAVar_NeutralHadronIso_DR0p3To0p4 = TMath::Max(TMath::Min((tmpNeutralHadronIso_DR0p3To0p4) / ele.pt(), 2.5), 0.0);
1887     fMVAVar_NeutralHadronIso_DR0p4To0p5 = TMath::Max(TMath::Min((tmpNeutralHadronIso_DR0p4To0p5) / ele.pt(), 2.5), 0.0);
1888     fMVAVar_rho = Rho;
1889   } else {
1890     std::cout << "Warning: Type " << fMVAType << " is not supported.\n";
1891   }
1892 
1893   // evaluate
1894   Double_t mva = -9999;
1895   if (fUseBinnedVersion) {
1896     int bin = GetMVABin(fMVAVar_eta, fMVAVar_pt);
1897     mva = fTMVAReader[bin]->EvaluateMVA(fTMVAMethod[bin]);
1898   } else {
1899     mva = fTMVAReader[0]->EvaluateMVA(fTMVAMethod[0]);
1900   }
1901 
1902   if (printDebug) {
1903     std::cout << " *** Inside the class fMethodname " << fMethodname << " fMVAType " << fMVAType << std::endl;
1904     std::cout << " fbrem " << fMVAVar_fbrem << " kfchi2 " << fMVAVar_kfchi2 << " mykfhits " << fMVAVar_kfhits
1905               << " gsfchi2 " << fMVAVar_gsfchi2 << " deta " << fMVAVar_deta << " dphi " << fMVAVar_dphi << " detacalo "
1906               << fMVAVar_detacalo << " see " << fMVAVar_see << " spp " << fMVAVar_spp << " etawidth "
1907               << fMVAVar_etawidth << " phiwidth " << fMVAVar_phiwidth << " OneMinusE1x5E5x5 "
1908               << fMVAVar_OneMinusE1x5E5x5 << " R9 " << fMVAVar_R9 << " HoE " << fMVAVar_HoE << " EoP " << fMVAVar_EoP
1909               << " IoEmIoP " << fMVAVar_IoEmIoP << " eleEoPout " << fMVAVar_eleEoPout << " d0 " << fMVAVar_d0
1910               << " ip3d " << fMVAVar_ip3d << " eta " << fMVAVar_eta << " pt " << fMVAVar_pt << std::endl;
1911     std::cout << "ChargedIso ( 0.0 | 0.1 | 0.2 | 0.3 | 0.4 | 0.5 ): " << fMVAVar_ChargedIso_DR0p0To0p1 << " "
1912               << fMVAVar_ChargedIso_DR0p1To0p2 << " " << fMVAVar_ChargedIso_DR0p2To0p3 << " "
1913               << fMVAVar_ChargedIso_DR0p3To0p4 << " " << fMVAVar_ChargedIso_DR0p4To0p5 << std::endl;
1914     std::cout << "PF Gamma Iso ( 0.0 | 0.1 | 0.2 | 0.3 | 0.4 | 0.5 ): " << fMVAVar_GammaIso_DR0p0To0p1 << " "
1915               << fMVAVar_GammaIso_DR0p1To0p2 << " " << fMVAVar_GammaIso_DR0p2To0p3 << " " << fMVAVar_GammaIso_DR0p3To0p4
1916               << " " << fMVAVar_GammaIso_DR0p4To0p5 << std::endl;
1917     std::cout << "PF Neutral Hadron Iso ( 0.0 | 0.1 | 0.2 | 0.3 | 0.4 | 0.5 ): " << fMVAVar_NeutralHadronIso_DR0p0To0p1
1918               << " " << fMVAVar_NeutralHadronIso_DR0p1To0p2 << " " << fMVAVar_NeutralHadronIso_DR0p2To0p3 << " "
1919               << fMVAVar_NeutralHadronIso_DR0p3To0p4 << " " << fMVAVar_NeutralHadronIso_DR0p4To0p5 << " " << std::endl;
1920     std::cout << "Rho : " << Rho << std::endl;
1921     std::cout << " ### MVA " << mva << std::endl;
1922   }
1923 
1924   return mva;
1925 }
1926 
1927 #endif
1928 
1929 void EGammaMvaEleEstimator::bindVariables() {
1930   // this binding is needed for variables that sometime diverge.
1931 
1932   if (fMVAVar_fbrem < -1.)
1933     fMVAVar_fbrem = -1.;
1934 
1935   fMVAVar_deta = fabs(fMVAVar_deta);
1936   if (fMVAVar_deta > 0.06)
1937     fMVAVar_deta = 0.06;
1938 
1939   fMVAVar_dphi = fabs(fMVAVar_dphi);
1940   if (fMVAVar_dphi > 0.6)
1941     fMVAVar_dphi = 0.6;
1942 
1943   if (fMVAVar_EoP > 20.)
1944     fMVAVar_EoP = 20.;
1945 
1946   if (fMVAVar_eleEoPout > 20.)
1947     fMVAVar_eleEoPout = 20.;
1948 
1949   fMVAVar_detacalo = fabs(fMVAVar_detacalo);
1950   if (fMVAVar_detacalo > 0.2)
1951     fMVAVar_detacalo = 0.2;
1952 
1953   if (fMVAVar_OneMinusE1x5E5x5 < -1.)
1954     fMVAVar_OneMinusE1x5E5x5 = -1;
1955 
1956   if (fMVAVar_OneMinusE1x5E5x5 > 2.)
1957     fMVAVar_OneMinusE1x5E5x5 = 2.;
1958 
1959   if (fMVAVar_R9 > 5)
1960     fMVAVar_R9 = 5;
1961 
1962   if (fMVAVar_gsfchi2 > 200.)
1963     fMVAVar_gsfchi2 = 200;
1964 
1965   if (fMVAVar_kfchi2 > 10.)
1966     fMVAVar_kfchi2 = 10.;
1967 
1968     // Needed for a bug in CMSSW_420, fixed in more recent CMSSW versions
1969 #ifndef STANDALONE
1970   if (edm::isNotFinite(fMVAVar_spp))
1971 #else
1972   if (std::isnan(fMVAVar_spp))
1973 #endif
1974     fMVAVar_spp = 0.;
1975 
1976   return;
1977 }