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
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
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
0064 fisInitialized = kTRUE;
0065 fMVAType = type;
0066 fMethodname = methodName;
0067 fUseBinnedVersion = useBinnedVersion;
0068
0069
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
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
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
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
0112 tmpTMVAReader->AddVariable("deta", &fMVAVar_deta);
0113 tmpTMVAReader->AddVariable("dphi", &fMVAVar_dphi);
0114 tmpTMVAReader->AddVariable("detacalo", &fMVAVar_detacalo);
0115
0116
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
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
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
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
0151 tmpTMVAReader->AddVariable("deta", &fMVAVar_deta);
0152 tmpTMVAReader->AddVariable("dphi", &fMVAVar_dphi);
0153 tmpTMVAReader->AddVariable("detacalo", &fMVAVar_detacalo);
0154
0155
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
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
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
0188 tmpTMVAReader->AddVariable("deta", &fMVAVar_deta);
0189 tmpTMVAReader->AddVariable("dphi", &fMVAVar_dphi);
0190 tmpTMVAReader->AddVariable("detacalo", &fMVAVar_detacalo);
0191
0192
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
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
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
0243 tmpTMVAReader->AddVariable("deta", &fMVAVar_deta);
0244 tmpTMVAReader->AddVariable("dphi", &fMVAVar_dphi);
0245 tmpTMVAReader->AddVariable("detacalo", &fMVAVar_detacalo);
0246
0247
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
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
0268 tmpTMVAReader->AddVariable("d0", &fMVAVar_d0);
0269 tmpTMVAReader->AddVariable("ip3d", &fMVAVar_ip3d);
0270
0271
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
0289 tmpTMVAReader->AddSpectator("eta", &fMVAVar_eta);
0290 tmpTMVAReader->AddSpectator("pt", &fMVAVar_pt);
0291 }
0292
0293 if (type == kTrigIDIsoCombined) {
0294
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
0301 tmpTMVAReader->AddVariable("deta", &fMVAVar_deta);
0302 tmpTMVAReader->AddVariable("dphi", &fMVAVar_dphi);
0303 tmpTMVAReader->AddVariable("detacalo", &fMVAVar_detacalo);
0304
0305
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
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
0326 tmpTMVAReader->AddVariable("d0", &fMVAVar_d0);
0327 tmpTMVAReader->AddVariable("ip3d", &fMVAVar_ip3d);
0328
0329
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
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');
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
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
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
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
0468 Double_t HoE,
0469 Double_t EoP,
0470 Double_t IoEmIoP,
0471 Double_t eleEoPout,
0472 Double_t PreShowerOverRaw,
0473
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);
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
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);
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
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
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
0642 Double_t HoE,
0643 Double_t EoP,
0644 Double_t IoEmIoP,
0645 Double_t eleEoPout,
0646 Double_t PreShowerOverRaw,
0647
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);
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);
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
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
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
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.;
0980 fMVAVar_gsfchi2 = ele.gsfTrack()->normalizedChi2();
0981
0982
0983 fMVAVar_deta = ele.deltaEtaSuperClusterTrackAtVtx();
0984 fMVAVar_dphi = ele.deltaPhiSuperClusterTrackAtVtx();
0985 fMVAVar_detacalo = ele.deltaEtaSeedClusterTrackAtCalo();
0986
0987
0988 fMVAVar_see = ele.sigmaIetaIeta();
0989 const auto& vCov = myEcalCluster.localCovariances(*(ele.superCluster()->seed()));
0990 if (edm::isFinite(vCov[2]))
0991 fMVAVar_spp = sqrt(vCov[2]);
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
1001 fMVAVar_HoE = ele.hadronicOverEm();
1002 fMVAVar_EoP = ele.eSuperClusterOverP();
1003 fMVAVar_IoEmIoP = (1.0 / ele.ecalEnergy()) - (1.0 / ele.p());
1004 fMVAVar_eleEoPout = ele.eEleClusterOverPout();
1005 fMVAVar_PreShowerOverRaw = ele.superCluster()->preshowerEnergy() / ele.superCluster()->rawEnergy();
1006
1007
1008 fMVAVar_eta = ele.superCluster()->eta();
1009 fMVAVar_pt = ele.pt();
1010
1011
1012 if (fMVAType == kTrig) {
1013
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
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
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
1065 Double_t EGammaMvaEleEstimator::mvaValue(const reco::GsfElectron& ele,
1066 const reco::Vertex& vertex,
1067 double rho,
1068
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
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
1093 fMVAVar_deta = ele.deltaEtaSuperClusterTrackAtVtx();
1094 fMVAVar_dphi = ele.deltaPhiSuperClusterTrackAtVtx();
1095 fMVAVar_detacalo = ele.deltaEtaSeedClusterTrackAtCalo();
1096
1097
1098 fMVAVar_see = ele.sigmaIetaIeta();
1099 const auto& vCov = myEcalCluster.localCovariances(*(ele.superCluster()->seed()));
1100 if (edm::isFinite(vCov[2]))
1101 fMVAVar_spp = sqrt(vCov[2]);
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
1111 fMVAVar_HoE = ele.hadronicOverEm();
1112 fMVAVar_EoP = ele.eSuperClusterOverP();
1113 fMVAVar_IoEmIoP = (1.0 / ele.superCluster()->energy()) -
1114 (1.0 / ele.gsfTrack()->p());
1115 fMVAVar_eleEoPout = ele.eEleClusterOverPout();
1116 fMVAVar_rho = rho;
1117 fMVAVar_PreShowerOverRaw = ele.superCluster()->preshowerEnergy() / ele.superCluster()->rawEnergy();
1118
1119
1120 fMVAVar_eta = ele.superCluster()->eta();
1121 fMVAVar_pt = ele.pt();
1122
1123
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
1139 << " see " << fMVAVar_see << " spp " << fMVAVar_spp << " etawidth " << fMVAVar_etawidth << " phiwidth "
1140 << fMVAVar_phiwidth << " e1x5e5x5 " << fMVAVar_OneMinusE1x5E5x5 << " R9 "
1141 << fMVAVar_R9
1142
1143 << " HoE " << fMVAVar_HoE << " EoP " << fMVAVar_EoP << " IoEmIoP " << fMVAVar_IoEmIoP << " eleEoPout "
1144 << fMVAVar_eleEoPout << " rho "
1145 << fMVAVar_rho
1146
1147 << " eta " << fMVAVar_eta << " pt " << fMVAVar_pt << std::endl;
1148 std::cout << " ### MVA " << mva << std::endl;
1149 }
1150
1151 return mva;
1152 }
1153
1154
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
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
1174 fMVAVar_deta = ele.deltaEtaSuperClusterTrackAtVtx();
1175 fMVAVar_dphi = ele.deltaPhiSuperClusterTrackAtVtx();
1176 fMVAVar_detacalo = ele.deltaEtaSeedClusterTrackAtCalo();
1177
1178
1179 fMVAVar_see = useFull5x5 ? ele.full5x5_sigmaIetaIeta() : ele.sigmaIetaIeta();
1180 fMVAVar_spp = useFull5x5 ? ele.full5x5_sigmaIphiIphi() : ele.sigmaIphiIphi();
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
1190 fMVAVar_HoE = ele.hadronicOverEm();
1191 fMVAVar_EoP = ele.eSuperClusterOverP();
1192
1193
1194 if (fMVAType == kTrig || fMVAType == kNonTrig) {
1195 fMVAVar_IoEmIoP =
1196 (1.0 / ele.ecalEnergy()) - (1.0 / ele.p());
1197 } else {
1198 fMVAVar_IoEmIoP = (1.0 / ele.superCluster()->energy()) -
1199 (1.0 / ele.gsfTrack()->p());
1200 }
1201 fMVAVar_eleEoPout = ele.eEleClusterOverPout();
1202 fMVAVar_rho = rho;
1203 fMVAVar_PreShowerOverRaw = ele.superCluster()->preshowerEnergy() / ele.superCluster()->rawEnergy();
1204
1205
1206 if (fMVAType == kTrig) {
1207
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
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
1223 double ip3d = gsfsign * ele.dB();
1224 double ip3derr = ele.edB();
1225 fMVAVar_ip3d = ip3d;
1226 fMVAVar_ip3dSig = ip3d / ip3derr;
1227 }
1228 }
1229
1230
1231 fMVAVar_eta = ele.superCluster()->eta();
1232 fMVAVar_pt = ele.pt();
1233
1234
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
1250 << " see " << fMVAVar_see << " spp " << fMVAVar_spp << " etawidth " << fMVAVar_etawidth << " phiwidth "
1251 << fMVAVar_phiwidth << " e1x5e5x5 " << fMVAVar_OneMinusE1x5E5x5 << " R9 "
1252 << fMVAVar_R9
1253
1254 << " HoE " << fMVAVar_HoE << " EoP " << fMVAVar_EoP << " IoEmIoP " << fMVAVar_IoEmIoP << " eleEoPout "
1255 << fMVAVar_eleEoPout << " rho "
1256 << fMVAVar_rho
1257
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
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
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
1290 fMVAVar_deta = ele.deltaEtaSuperClusterTrackAtVtx();
1291 fMVAVar_dphi = ele.deltaPhiSuperClusterTrackAtVtx();
1292 fMVAVar_detacalo = ele.deltaEtaSeedClusterTrackAtCalo();
1293
1294
1295 fMVAVar_see = ele.sigmaIetaIeta();
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
1305 fMVAVar_HoE = ele.hadronicOverEm();
1306 fMVAVar_EoP = ele.eSuperClusterOverP();
1307 fMVAVar_IoEmIoP = (1.0 / ele.superCluster()->energy()) -
1308 (1.0 / ele.gsfTrack()->p());
1309 fMVAVar_eleEoPout = ele.eEleClusterOverPout();
1310 fMVAVar_rho = rho;
1311 fMVAVar_PreShowerOverRaw = ele.superCluster()->preshowerEnergy() / ele.superCluster()->rawEnergy();
1312
1313
1314 fMVAVar_eta = ele.superCluster()->eta();
1315 fMVAVar_pt = ele.pt();
1316
1317
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
1333 << " see " << fMVAVar_see << " spp " << fMVAVar_spp << " etawidth " << fMVAVar_etawidth << " phiwidth "
1334 << fMVAVar_phiwidth << " e1x5e5x5 " << fMVAVar_OneMinusE1x5E5x5 << " R9 "
1335 << fMVAVar_R9
1336
1337 << " HoE " << fMVAVar_HoE << " EoP " << fMVAVar_EoP << " IoEmIoP " << fMVAVar_IoEmIoP << " eleEoPout "
1338 << fMVAVar_eleEoPout << " rho "
1339 << fMVAVar_rho
1340
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
1362 fMVAVar_eta = ele.superCluster()->eta();
1363 fMVAVar_pt = ele.pt();
1364
1365
1366
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
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
1402
1403 double dr = sqrt(pow(iP->eta() - ele.eta(), 2) + pow(acos(cos(iP->phi() - ele.phi())), 2));
1404
1405
1406 if (dr < 1.0) {
1407 Bool_t IsLeptonFootprint = kFALSE;
1408
1409
1410
1411 for (reco::GsfElectronCollection::const_iterator iE = IdentifiedElectrons.begin();
1412 iE != IdentifiedElectrons.end();
1413 ++iE) {
1414
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
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
1431 if (iP->trackRef().isNonnull() && iM->innerTrack().isNonnull() &&
1432 refToPtr(iP->trackRef()) == refToPtr(iM->innerTrack()))
1433 IsLeptonFootprint = kTRUE;
1434
1435
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
1444 if (iP->trackRef().isNonnull()) {
1445 if (!(fabs(iP->trackRef()->dz(vertex.position()) - electronTrackZ) < 0.2))
1446 passVeto = kFALSE;
1447
1448
1449 if (iP->particleId() == reco::PFCandidate::e || iP->particleId() == reco::PFCandidate::mu)
1450 passVeto = kFALSE;
1451
1452
1453
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 }
1469 }
1470
1471 else if (iP->particleId() == reco::PFCandidate::gamma) {
1472
1473
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
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 }
1504 }
1505 }
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
1589 bindVariables();
1590 Double_t mva = -9999;
1591
1592
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
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.;
1645 fMVAVar_gsfchi2 = ele.gsfTrack()->normalizedChi2();
1646 if (fMVAVar_gsfchi2 > 200)
1647 fMVAVar_gsfchi2 = 200;
1648
1649
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
1656 fMVAVar_see = ele.sigmaIetaIeta();
1657 const auto& vCov = myEcalCluster.localCovariances(*(ele.superCluster()->seed()));
1658 if (edm::isFinite(vCov[2]))
1659 fMVAVar_spp = sqrt(vCov[2]);
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
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());
1676 fMVAVar_eleEoPout = (ele.eEleClusterOverPout() > 20) ? 20 : ele.eEleClusterOverPout();
1677 fMVAVar_PreShowerOverRaw = ele.superCluster()->preshowerEnergy() / ele.superCluster()->rawEnergy();
1678
1679
1680 fMVAVar_eta = ele.superCluster()->eta();
1681 fMVAVar_pt = ele.pt();
1682
1683
1684 if (fMVAType == kTrig) {
1685
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
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
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
1735 if (iP->trackRef().isNonnull()) {
1736
1737
1738
1739
1740 if (iP->particleId() == reco::PFCandidate::e || iP->particleId() == reco::PFCandidate::mu)
1741 passVeto = kFALSE;
1742
1743
1744
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 }
1760 }
1761
1762 else if (iP->particleId() == reco::PFCandidate::gamma) {
1763
1764
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
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 }
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
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
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
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 }