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 const std::vector<std::string> ElectronDNNEstimator::dnnAvaibleInputs = {
0051     {"pt",
0052      "eta",
0053      "fbrem",
0054      "abs(deltaEtaSuperClusterTrackAtVtx)",
0055      "abs(deltaPhiSuperClusterTrackAtVtx)",
0056      "full5x5_sigmaIetaIeta",
0057      "full5x5_hcalOverEcal",
0058      "eSuperClusterOverP",
0059      "full5x5_e1x5",
0060      "eEleClusterOverPout",
0061      "closestCtfTrackNormChi2",
0062      "closestCtfTrackNLayers",
0063      "gsfTrack.missing_inner_hits",
0064      "dr03TkSumPt",
0065      "dr03EcalRecHitSumEt",
0066      "dr03HcalTowerSumEt",
0067      "gsfTrack.normalizedChi2",
0068      "superCluster.eta",
0069      "ecalPFClusterIso",
0070      "hcalPFClusterIso",
0071      "numberOfBrems",
0072      "abs(deltaEtaSeedClusterTrackAtCalo)",
0073      "hadronicOverEm",
0074      "full5x5_e2x5Max",
0075      "full5x5_e5x5",
0076      "full5x5_sigmaIphiIphi",
0077      "1_minus_full5x5_e1x5_ratio_full5x5_e5x5",
0078      "full5x5_e1x5_ratio_full5x5_e5x5",
0079      "full5x5_r9",
0080      "gsfTrack.trackerLayersWithMeasurement",
0081      "gsfTrack.numberOfValidPixelBarrelHits",
0082      "gsfTrack.numberOfValidPixelEndcapHits",
0083      "superCluster.energy",
0084      "superCluster.rawEnergy",
0085      "superClusterFbrem",
0086      "1_ratio_ecalEnergy_minus_1_ratio_trackMomentumAtVtx.R",
0087      "superCluster.preshowerEnergy_ratio_superCluster.rawEnergy",
0088      "convVtxFitProb",
0089      "superCluster.clustersSize",
0090      "ecalEnergyError_ratio_ecalEnergy",
0091      "superClusterFbrem_plus_superCluster.energy",
0092      "superClusterFbrem_plus_superCluster.rawEnergy",
0093      "trackMomentumError",
0094      "trackMomentumError_ratio_pt",
0095      "full5x5_e5x5_ratio_superCluster.rawEnergy",
0096      "full5x5_e5x5_ratio_superCluster.energy"}};
0097
0098 std::map<std::string, float> ElectronDNNEstimator::getInputsVars(const reco::GsfElectron& ele) const {
0099   // Prepare a map with all the defined variables
0100   std::map<std::string, float> variables;
0101   variables["pt"] = ele.pt();
0102   variables["eta"] = ele.eta();
0103   variables["fbrem"] = ele.fbrem();
0104   variables["abs(deltaEtaSuperClusterTrackAtVtx)"] = std::abs(ele.deltaEtaSuperClusterTrackAtVtx());
0105   variables["abs(deltaPhiSuperClusterTrackAtVtx)"] = std::abs(ele.deltaPhiSuperClusterTrackAtVtx());
0106   variables["full5x5_sigmaIetaIeta"] = ele.full5x5_sigmaIetaIeta();
0107   variables["full5x5_hcalOverEcal"] = ele.full5x5_hcalOverEcalValid() ? ele.full5x5_hcalOverEcal() : 0;
0108   variables["eSuperClusterOverP"] = ele.eSuperClusterOverP();
0109   variables["full5x5_e1x5"] = ele.full5x5_e1x5();
0110   variables["eEleClusterOverPout"] = ele.eEleClusterOverPout();
0111   variables["closestCtfTrackNormChi2"] = ele.closestCtfTrackNormChi2();
0112   variables["closestCtfTrackNLayers"] = ele.closestCtfTrackNLayers();
0113   variables["gsfTrack.missing_inner_hits"] =
0114       ele.gsfTrack()->hitPattern().numberOfLostHits(reco::HitPattern::MISSING_INNER_HITS);
0115   variables["dr03TkSumPt"] = ele.dr03TkSumPt();
0116   variables["dr03EcalRecHitSumEt"] = ele.dr03EcalRecHitSumEt();
0117   variables["dr03HcalTowerSumEt"] = ele.dr03HcalTowerSumEt();
0118   variables["gsfTrack.normalizedChi2"] = ele.gsfTrack()->normalizedChi2();
0119   variables["superCluster.eta"] = ele.superCluster()->eta();
0120   variables["ecalPFClusterIso"] = ele.ecalPFClusterIso();
0121   variables["hcalPFClusterIso"] = ele.hcalPFClusterIso();
0122   variables["numberOfBrems"] = ele.numberOfBrems();
0123   variables["abs(deltaEtaSeedClusterTrackAtCalo)"] = std::abs(ele.deltaEtaSeedClusterTrackAtCalo());
0124   variables["hadronicOverEm"] = ele.hcalOverEcalValid() ? ele.hadronicOverEm() : 0;
0125   variables["full5x5_e2x5Max"] = ele.full5x5_e2x5Max();
0126   variables["full5x5_e5x5"] = ele.full5x5_e5x5();
0127
0128   variables["full5x5_sigmaIphiIphi"] = ele.full5x5_sigmaIphiIphi();
0129   variables["1_minus_full5x5_e1x5_ratio_full5x5_e5x5"] = 1 - ele.full5x5_e1x5() / ele.full5x5_e5x5();
0130   variables["full5x5_e1x5_ratio_full5x5_e5x5"] = ele.full5x5_e1x5() / ele.full5x5_e5x5();
0131   variables["full5x5_r9"] = ele.full5x5_r9();
0132   variables["gsfTrack.trackerLayersWithMeasurement"] = ele.gsfTrack()->hitPattern().trackerLayersWithMeasurement();
0133   variables["gsfTrack.numberOfValidPixelBarrelHits"] = ele.gsfTrack()->hitPattern().numberOfValidPixelBarrelHits();
0134   variables["gsfTrack.numberOfValidPixelEndcapHits"] = ele.gsfTrack()->hitPattern().numberOfValidPixelEndcapHits();
0135   variables["superCluster.energy"] = ele.superCluster()->energy();
0136   variables["superCluster.rawEnergy"] = ele.superCluster()->rawEnergy();
0137   variables["superClusterFbrem"] = ele.superClusterFbrem();
0138   variables["1_ratio_ecalEnergy_minus_1_ratio_trackMomentumAtVtx.R"] =
0139       1 / ele.ecalEnergy() - 1 / ele.trackMomentumAtVtx().R();
0140   variables["superCluster.preshowerEnergy_ratio_superCluster.rawEnergy"] =
0141       ele.superCluster()->preshowerEnergy() / ele.superCluster()->rawEnergy();
0142   variables["convVtxFitProb"] = ele.convVtxFitProb();
0143   variables["superCluster.clustersSize"] = ele.superCluster()->clustersSize();
0144   variables["ecalEnergyError_ratio_ecalEnergy"] = ele.ecalEnergyError() / ele.ecalEnergy();
0145   variables["superClusterFbrem_plus_superCluster.energy"] = ele.superClusterFbrem() + ele.superCluster()->energy();
0146   variables["superClusterFbrem_plus_superCluster.rawEnergy"] =
0147       ele.superClusterFbrem() + ele.superCluster()->rawEnergy();
0148   variables["trackMomentumError"] = ele.trackMomentumError();
0149   variables["trackMomentumError_ratio_pt"] = ele.trackMomentumError() / ele.pt();
0150   variables["full5x5_e5x5_ratio_superCluster.rawEnergy"] = ele.full5x5_e5x5() / ele.superCluster()->rawEnergy();
0151   variables["full5x5_e5x5_ratio_superCluster.energy"] = ele.full5x5_e5x5() / ele.superCluster()->energy();
0152
0153   // Define more variables here and use them directly in the model config!
0154   return variables;
0155 }
0156
0157 std::vector<std::pair<uint, std::vector<float>>> ElectronDNNEstimator::evaluate(
0158     const reco::GsfElectronCollection& electrons) const {
0159   // Collect the map of variables for each candidate and call the dnnHelper
0160   // Scaling, model selection and running is performed in the helper
0161   std::vector<std::map<std::string, float>> inputs;
0162   for (const auto& ele : electrons) {
0163     inputs.push_back(getInputsVars(ele));
0164   }
0165   return dnnHelper_.evaluate(inputs);
0166 }