ElectronIdentification/src/ElectronDNNEstimator.cc

0001 #include "RecoEgamma/ElectronIdentification/interface/ElectronDNNEstimator.h"
0002 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0003 #include "FWCore/Utilities/interface/FileInPath.h"
0004 #include "DataFormats/GsfTrackReco/interface/GsfTrack.h"
0005
0006 #include <iostream>
0007 #include <fstream>
0008 #include <memory>
0009
0010 using namespace std::placeholders;
0011
0012 inline uint electronModelSelector(
0013     const std::map<std::string, float>& vars, float ptThr, float etaThr, float endcapBoundary, float extEtaBoundary) {
0014   /*
0015   Selection of the model to be applied on the electron based on pt/eta cuts or whatever selection
0016   */
0017   const auto pt = vars.at("pt");
0018   const auto absEta = std::abs(vars.at("superCluster.eta"));
0019   if (absEta <= endcapBoundary) {
0020     if (pt < ptThr)
0021       return 0;
0022     else {
0023       if (absEta <= etaThr) {
0024         return 1;
0025       } else {
0026         return 2;
0027       }
0028     }
0029   } else {
0030     if (absEta < extEtaBoundary)
0031       return 3;
0032     else {
0033       return 4;
0034     }
0035   }
0036 }
0037
0038 ElectronDNNEstimator::ElectronDNNEstimator(const egammaTools::DNNConfiguration& cfg, const bool useEBModelInGap)
0039     : dnnHelper_(cfg,
0040                  std::bind(electronModelSelector,
0041                            _1,
0042                            ElectronDNNEstimator::ptThreshold,
0043                            (useEBModelInGap) ? ElectronDNNEstimator::ecalBarrelMaxEtaWithGap
0044                                              : ElectronDNNEstimator::ecalBarrelMaxEtaNoGap,
0045                            ElectronDNNEstimator::endcapBoundary,
0046                            ElectronDNNEstimator::extEtaBoundary),
0047                  ElectronDNNEstimator::dnnAvaibleInputs),
0048       useEBModelInGap_(useEBModelInGap) {}
0049
0050 std::vector<tensorflow::Session*> ElectronDNNEstimator::getSessions() const { return dnnHelper_.getSessions(); };
0051
0052 const std::vector<std::string> ElectronDNNEstimator::dnnAvaibleInputs = {
0053     {"pt",
0054      "eta",
0055      "fbrem",
0056      "abs(deltaEtaSuperClusterTrackAtVtx)",
0057      "abs(deltaPhiSuperClusterTrackAtVtx)",
0058      "full5x5_sigmaIetaIeta",
0059      "full5x5_hcalOverEcal",
0060      "eSuperClusterOverP",
0061      "full5x5_e1x5",
0062      "eEleClusterOverPout",
0063      "closestCtfTrackNormChi2",
0064      "closestCtfTrackNLayers",
0065      "gsfTrack.missing_inner_hits",
0066      "dr03TkSumPt",
0067      "dr03EcalRecHitSumEt",
0068      "dr03HcalTowerSumEt",
0069      "gsfTrack.normalizedChi2",
0070      "superCluster.eta",
0071      "ecalPFClusterIso",
0072      "hcalPFClusterIso",
0073      "numberOfBrems",
0074      "abs(deltaEtaSeedClusterTrackAtCalo)",
0075      "hadronicOverEm",
0076      "full5x5_e2x5Max",
0077      "full5x5_e5x5",
0078      "full5x5_sigmaIphiIphi",
0079      "1_minus_full5x5_e1x5_ratio_full5x5_e5x5",
0080      "full5x5_e1x5_ratio_full5x5_e5x5",
0081      "full5x5_r9",
0082      "gsfTrack.trackerLayersWithMeasurement",
0083      "gsfTrack.numberOfValidPixelBarrelHits",
0084      "gsfTrack.numberOfValidPixelEndcapHits",
0085      "superCluster.energy",
0086      "superCluster.rawEnergy",
0087      "superClusterFbrem",
0088      "1_ratio_ecalEnergy_minus_1_ratio_trackMomentumAtVtx.R",
0089      "superCluster.preshowerEnergy_ratio_superCluster.rawEnergy",
0090      "convVtxFitProb",
0091      "superCluster.clustersSize",
0092      "ecalEnergyError_ratio_ecalEnergy",
0093      "superClusterFbrem_plus_superCluster.energy",
0094      "superClusterFbrem_plus_superCluster.rawEnergy",
0095      "trackMomentumError",
0096      "trackMomentumError_ratio_pt",
0097      "full5x5_e5x5_ratio_superCluster.rawEnergy",
0098      "full5x5_e5x5_ratio_superCluster.energy"}};
0099
0100 std::map<std::string, float> ElectronDNNEstimator::getInputsVars(const reco::GsfElectron& ele) const {
0101   // Prepare a map with all the defined variables
0102   std::map<std::string, float> variables;
0103   variables["pt"] = ele.pt();
0104   variables["eta"] = ele.eta();
0105   variables["fbrem"] = ele.fbrem();
0106   variables["abs(deltaEtaSuperClusterTrackAtVtx)"] = std::abs(ele.deltaEtaSuperClusterTrackAtVtx());
0107   variables["abs(deltaPhiSuperClusterTrackAtVtx)"] = std::abs(ele.deltaPhiSuperClusterTrackAtVtx());
0108   variables["full5x5_sigmaIetaIeta"] = ele.full5x5_sigmaIetaIeta();
0109   variables["full5x5_hcalOverEcal"] = ele.full5x5_hcalOverEcalValid() ? ele.full5x5_hcalOverEcal() : 0;
0110   variables["eSuperClusterOverP"] = ele.eSuperClusterOverP();
0111   variables["full5x5_e1x5"] = ele.full5x5_e1x5();
0112   variables["eEleClusterOverPout"] = ele.eEleClusterOverPout();
0113   variables["closestCtfTrackNormChi2"] = ele.closestCtfTrackNormChi2();
0114   variables["closestCtfTrackNLayers"] = ele.closestCtfTrackNLayers();
0115   variables["gsfTrack.missing_inner_hits"] =
0116       ele.gsfTrack()->hitPattern().numberOfLostHits(reco::HitPattern::MISSING_INNER_HITS);
0117   variables["dr03TkSumPt"] = ele.dr03TkSumPt();
0118   variables["dr03EcalRecHitSumEt"] = ele.dr03EcalRecHitSumEt();
0119   variables["dr03HcalTowerSumEt"] = ele.dr03HcalTowerSumEt();
0120   variables["gsfTrack.normalizedChi2"] = ele.gsfTrack()->normalizedChi2();
0121   variables["superCluster.eta"] = ele.superCluster()->eta();
0122   variables["ecalPFClusterIso"] = ele.ecalPFClusterIso();
0123   variables["hcalPFClusterIso"] = ele.hcalPFClusterIso();
0124   variables["numberOfBrems"] = ele.numberOfBrems();
0125   variables["abs(deltaEtaSeedClusterTrackAtCalo)"] = std::abs(ele.deltaEtaSeedClusterTrackAtCalo());
0126   variables["hadronicOverEm"] = ele.hcalOverEcalValid() ? ele.hadronicOverEm() : 0;
0127   variables["full5x5_e2x5Max"] = ele.full5x5_e2x5Max();
0128   variables["full5x5_e5x5"] = ele.full5x5_e5x5();
0129
0130   variables["full5x5_sigmaIphiIphi"] = ele.full5x5_sigmaIphiIphi();
0131   variables["1_minus_full5x5_e1x5_ratio_full5x5_e5x5"] = 1 - ele.full5x5_e1x5() / ele.full5x5_e5x5();
0132   variables["full5x5_e1x5_ratio_full5x5_e5x5"] = ele.full5x5_e1x5() / ele.full5x5_e5x5();
0133   variables["full5x5_r9"] = ele.full5x5_r9();
0134   variables["gsfTrack.trackerLayersWithMeasurement"] = ele.gsfTrack()->hitPattern().trackerLayersWithMeasurement();
0135   variables["gsfTrack.numberOfValidPixelBarrelHits"] = ele.gsfTrack()->hitPattern().numberOfValidPixelBarrelHits();
0136   variables["gsfTrack.numberOfValidPixelEndcapHits"] = ele.gsfTrack()->hitPattern().numberOfValidPixelEndcapHits();
0137   variables["superCluster.energy"] = ele.superCluster()->energy();
0138   variables["superCluster.rawEnergy"] = ele.superCluster()->rawEnergy();
0139   variables["superClusterFbrem"] = ele.superClusterFbrem();
0140   variables["1_ratio_ecalEnergy_minus_1_ratio_trackMomentumAtVtx.R"] =
0141       1 / ele.ecalEnergy() - 1 / ele.trackMomentumAtVtx().R();
0142   variables["superCluster.preshowerEnergy_ratio_superCluster.rawEnergy"] =
0143       ele.superCluster()->preshowerEnergy() / ele.superCluster()->rawEnergy();
0144   variables["convVtxFitProb"] = ele.convVtxFitProb();
0145   variables["superCluster.clustersSize"] = ele.superCluster()->clustersSize();
0146   variables["ecalEnergyError_ratio_ecalEnergy"] = ele.ecalEnergyError() / ele.ecalEnergy();
0147   variables["superClusterFbrem_plus_superCluster.energy"] = ele.superClusterFbrem() + ele.superCluster()->energy();
0148   variables["superClusterFbrem_plus_superCluster.rawEnergy"] =
0149       ele.superClusterFbrem() + ele.superCluster()->rawEnergy();
0150   variables["trackMomentumError"] = ele.trackMomentumError();
0151   variables["trackMomentumError_ratio_pt"] = ele.trackMomentumError() / ele.pt();
0152   variables["full5x5_e5x5_ratio_superCluster.rawEnergy"] = ele.full5x5_e5x5() / ele.superCluster()->rawEnergy();
0153   variables["full5x5_e5x5_ratio_superCluster.energy"] = ele.full5x5_e5x5() / ele.superCluster()->energy();
0154
0155   // Define more variables here and use them directly in the model config!
0156   return variables;
0157 }
0158
0159 std::vector<std::pair<uint, std::vector<float>>> ElectronDNNEstimator::evaluate(
0160     const reco::GsfElectronCollection& electrons, const std::vector<tensorflow::Session*>& sessions) const {
0161   // Collect the map of variables for each candidate and call the dnnHelper
0162   // Scaling, model selection and running is performed in the helper
0163   std::vector<std::map<std::string, float>> inputs;
0164   for (const auto& ele : electrons) {
0165     inputs.push_back(getInputsVars(ele));
0166   }
0167   return dnnHelper_.evaluate(inputs, sessions);
0168 }