From b2bdbef6284897f11f90256b7526319ce3234b96 Mon Sep 17 00:00:00 2001 From: Srecko Date: Tue, 5 Dec 2023 11:54:03 +0100 Subject: [PATCH] * e2x2 produced by the EGammaClusterShapeProducer * Photon XGBoost MVA estimator and producer. Implemented using native XGBoost C API for photon MVA score Currently diphoton score filter is not yet submitted. Model files will be added to the cms-data XGBoost Double photon combination filter fix module name code-format fix comments from maintainers Update HLTrigger/Egamma/plugins/HLTEgammaDoubleXGBoostCombFilter.cc Co-authored-by: Marco Musich Update HLTrigger/Egamma/plugins/HLTEgammaDoubleXGBoostCombFilter.cc Co-authored-by: Marco Musich Update HLTrigger/Egamma/plugins/HLTEgammaDoubleXGBoostCombFilter.cc Co-authored-by: Marco Musich Update HLTrigger/Egamma/plugins/HLTEgammaDoubleXGBoostCombFilter.cc Co-authored-by: Marco Musich Update HLTrigger/Egamma/plugins/HLTEgammaDoubleXGBoostCombFilter.cc Co-authored-by: Marco Musich Update HLTrigger/Egamma/plugins/HLTEgammaDoubleXGBoostCombFilter.cc Co-authored-by: Marco Musich Update RecoEgamma/PhotonIdentification/plugins/PhotonXGBoostProducer.cc Co-authored-by: Marco Musich Update RecoEgamma/PhotonIdentification/plugins/PhotonXGBoostProducer.cc Co-authored-by: Marco Musich Update RecoEgamma/PhotonIdentification/plugins/PhotonXGBoostProducer.cc Co-authored-by: Marco Musich Update RecoEgamma/PhotonIdentification/plugins/PhotonXGBoostProducer.cc Co-authored-by: Marco Musich Update RecoEgamma/PhotonIdentification/plugins/PhotonXGBoostProducer.cc Co-authored-by: Marco Musich tidy up and sort headers code-format Update RecoEgamma/PhotonIdentification/plugins/PhotonXGBoostProducer.cc Co-authored-by: Marco Musich --- .../HLTEgammaDoubleXGBoostCombFilter.cc | 157 ++++++++++++++++++ .../plugins/EgammaHLTClusterShapeProducer.cc | 10 ++ RecoEgamma/PhotonIdentification/BuildFile.xml | 1 + .../interface/PhotonXGBoostEstimator.h | 28 ++++ .../plugins/PhotonXGBoostProducer.cc | 137 +++++++++++++++ .../src/PhotonXGBoostEstimator.cc | 60 +++++++ 6 files changed, 393 insertions(+) create mode 100644 HLTrigger/Egamma/plugins/HLTEgammaDoubleXGBoostCombFilter.cc create mode 100644 RecoEgamma/PhotonIdentification/interface/PhotonXGBoostEstimator.h create mode 100644 RecoEgamma/PhotonIdentification/plugins/PhotonXGBoostProducer.cc create mode 100644 RecoEgamma/PhotonIdentification/src/PhotonXGBoostEstimator.cc diff --git a/HLTrigger/Egamma/plugins/HLTEgammaDoubleXGBoostCombFilter.cc b/HLTrigger/Egamma/plugins/HLTEgammaDoubleXGBoostCombFilter.cc new file mode 100644 index 0000000000000..8d3c22ea7fe95 --- /dev/null +++ b/HLTrigger/Egamma/plugins/HLTEgammaDoubleXGBoostCombFilter.cc @@ -0,0 +1,157 @@ +#include "DataFormats/HLTReco/interface/TriggerFilterObjectWithRefs.h" +#include "DataFormats/RecoCandidate/interface/RecoEcalCandidate.h" +#include "DataFormats/RecoCandidate/interface/RecoEcalCandidateIsolation.h" +#include "FWCore/Framework/interface/Event.h" +#include "FWCore/MessageLogger/interface/MessageLogger.h" +#include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h" +#include "FWCore/ParameterSet/interface/ParameterSet.h" +#include "FWCore/ParameterSet/interface/ParameterSetDescription.h" +#include "FWCore/Utilities/interface/InputTag.h" +#include "HLTrigger/HLTcore/interface/HLTFilter.h" + +#include +#include + +namespace edm { + class ConfigurationDescriptions; +} + +class HLTEgammaDoubleXGBoostCombFilter : public HLTFilter { +public: + explicit HLTEgammaDoubleXGBoostCombFilter(edm::ParameterSet const&); + ~HLTEgammaDoubleXGBoostCombFilter() override {} + + static void fillDescriptions(edm::ConfigurationDescriptions& descriptions); + +private: + bool hltFilter(edm::Event& event, + const edm::EventSetup& setup, + trigger::TriggerFilterObjectWithRefs& filterproduct) const override; + + const double highMassCut_; + const std::vector leadCutHighMass1_; + const std::vector subCutHighMass1_; + const std::vector leadCutHighMass2_; + const std::vector subCutHighMass2_; + const std::vector leadCutHighMass3_; + const std::vector subCutHighMass3_; + + const double lowMassCut_; + const std::vector leadCutLowMass1_; + const std::vector subCutLowMass1_; + const std::vector leadCutLowMass2_; + const std::vector subCutLowMass2_; + const std::vector leadCutLowMass3_; + const std::vector subCutLowMass3_; + + const edm::EDGetTokenT candToken_; + const edm::EDGetTokenT mvaToken_; +}; + +HLTEgammaDoubleXGBoostCombFilter::HLTEgammaDoubleXGBoostCombFilter(edm::ParameterSet const& config) + : HLTFilter(config), + highMassCut_(config.getParameter("highMassCut")), + leadCutHighMass1_(config.getParameter>("leadCutHighMass1")), + subCutHighMass1_(config.getParameter>("subCutHighMass1")), + leadCutHighMass2_(config.getParameter>("leadCutHighMass2")), + subCutHighMass2_(config.getParameter>("subCutHighMass2")), + leadCutHighMass3_(config.getParameter>("leadCutHighMass3")), + subCutHighMass3_(config.getParameter>("subCutHighMass3")), + + lowMassCut_(config.getParameter("lowMassCut")), + leadCutLowMass1_(config.getParameter>("leadCutLowMass1")), + subCutLowMass1_(config.getParameter>("subCutLowMass1")), + leadCutLowMass2_(config.getParameter>("leadCutLowMass2")), + subCutLowMass2_(config.getParameter>("subCutLowMass2")), + leadCutLowMass3_(config.getParameter>("leadCutLowMass3")), + subCutLowMass3_(config.getParameter>("subCutLowMass3")), + + candToken_(consumes(config.getParameter("candTag"))), + mvaToken_(consumes(config.getParameter("mvaPhotonTag"))) {} + +void HLTEgammaDoubleXGBoostCombFilter::fillDescriptions(edm::ConfigurationDescriptions& descriptions) { + edm::ParameterSetDescription desc; + makeHLTFilterDescription(desc); + + desc.add("highMassCut", 95.0); + desc.add>("leadCutHighMass1", {0.98, 0.95}); + desc.add>("subCutHighMass1", {0.0, 0.04}); + desc.add>("leadCutHighMass2", {0.85, 0.85}); + desc.add>("subCutHighMass2", {0.04, 0.08}); + desc.add>("leadCutHighMass3", {0.30, 0.50}); + desc.add>("subCutHighMass3", {0.15, 0.20}); + + desc.add("lowMassCut", 60.0); + desc.add>("leadCutLowMass1", {0.98, 0.90}); + desc.add>("subCutLowMass1", {0.04, 0.05}); + desc.add>("leadCutLowMass2", {0.90, 0.80}); + desc.add>("subCutLowMass2", {0.10, 0.10}); + desc.add>("leadCutLowMass3", {0.60, 0.60}); + desc.add>("subCutLowMass3", {0.30, 0.30}); + + desc.add("candTag", edm::InputTag("hltEgammaCandidatesUnseeded")); + desc.add("mvaPhotonTag", edm::InputTag("PhotonXGBoostProducer")); + + descriptions.addWithDefaultLabel(desc); +} + +bool HLTEgammaDoubleXGBoostCombFilter::hltFilter(edm::Event& event, + const edm::EventSetup& setup, + trigger::TriggerFilterObjectWithRefs& filterproduct) const { + //producer collection (hltEgammaCandidates(Unseeded)) + const auto& recCollection = event.getHandle(candToken_); + + //get hold of photon MVA association map + const auto& mvaMap = event.getHandle(mvaToken_); + + std::vector p4s(recCollection->size()); + std::vector isTight(recCollection->size()); + + bool accept = false; + + for (size_t i = 0; i < recCollection->size(); i++) { + edm::Ref refi(recCollection, i); + float EtaSCi = refi->eta(); + int eta1 = (std::abs(EtaSCi) < 1.5) ? 0 : 1; + float mvaScorei = (*mvaMap).find(refi)->val; + math::XYZTLorentzVector p4i = refi->p4(); + for (size_t j = i + 1; j < recCollection->size(); j++) { + edm::Ref refj(recCollection, j); + float EtaSCj = refj->eta(); + int eta2 = (std::abs(EtaSCj) < 1.5) ? 0 : 1; + float mvaScorej = (*mvaMap).find(refj)->val; + math::XYZTLorentzVector p4j = refj->p4(); + math::XYZTLorentzVector pairP4 = p4i + p4j; + double mass = pairP4.M(); + if (mass >= highMassCut_) { + if (mvaScorei >= mvaScorej && ((mvaScorei > leadCutHighMass1_[eta1] && mvaScorej > subCutHighMass1_[eta2]) || + (mvaScorei > leadCutHighMass2_[eta1] && mvaScorej > subCutHighMass2_[eta2]) || + (mvaScorei > leadCutHighMass3_[eta1] && mvaScorej > subCutHighMass3_[eta2]))) { + accept = true; + } //if scoreI > scoreJ + else if (mvaScorej > mvaScorei && + ((mvaScorej > leadCutHighMass1_[eta1] && mvaScorei > subCutHighMass1_[eta2]) || + (mvaScorej > leadCutHighMass2_[eta1] && mvaScorei > subCutHighMass2_[eta2]) || + (mvaScorej > leadCutHighMass3_[eta1] && mvaScorei > subCutHighMass3_[eta2]))) { + accept = true; + } // if scoreJ > scoreI + } //If high mass + else if (mass > lowMassCut_ && mass < highMassCut_) { + if (mvaScorei >= mvaScorej && ((mvaScorei > leadCutLowMass1_[eta1] && mvaScorej > subCutLowMass1_[eta2]) || + (mvaScorei > leadCutLowMass2_[eta1] && mvaScorej > subCutLowMass2_[eta2]) || + (mvaScorei > leadCutLowMass3_[eta1] && mvaScorej > subCutLowMass3_[eta2]))) { + accept = true; + } //if scoreI > scoreJ + else if (mvaScorej > mvaScorei && ((mvaScorej > leadCutLowMass1_[eta1] && mvaScorei > subCutLowMass1_[eta2]) || + (mvaScorej > leadCutLowMass2_[eta1] && mvaScorei > subCutLowMass2_[eta2]) || + (mvaScorej > leadCutLowMass3_[eta1] && mvaScorei > subCutLowMass3_[eta2]))) { + accept = true; + } //if scoreJ > scoreI + } //If low mass + } //j loop + } //i loop + return accept; +} //Definition + +#include "FWCore/Framework/interface/MakerMacros.h" +DEFINE_FWK_MODULE(HLTEgammaDoubleXGBoostCombFilter); diff --git a/RecoEgamma/EgammaHLTProducers/plugins/EgammaHLTClusterShapeProducer.cc b/RecoEgamma/EgammaHLTProducers/plugins/EgammaHLTClusterShapeProducer.cc index 09dda363550c5..58f0e4e4c42d0 100644 --- a/RecoEgamma/EgammaHLTProducers/plugins/EgammaHLTClusterShapeProducer.cc +++ b/RecoEgamma/EgammaHLTProducers/plugins/EgammaHLTClusterShapeProducer.cc @@ -63,6 +63,7 @@ EgammaHLTClusterShapeProducer::EgammaHLTClusterShapeProducer(const edm::Paramete produces("sigmaIPhiIPhi5x5NoiseCleaned"); produces("sMajor"); produces("sMinor"); + produces("e2x2"); } EgammaHLTClusterShapeProducer::~EgammaHLTClusterShapeProducer() {} @@ -108,6 +109,8 @@ void EgammaHLTClusterShapeProducer::produce(edm::StreamID sid, reco::RecoEcalCandidateIsolationMap clshSMajorMap(recoecalcandHandle); reco::RecoEcalCandidateIsolationMap clshSMinorMap(recoecalcandHandle); + reco::RecoEcalCandidateIsolationMap e2x2Map(recoecalcandHandle); + for (unsigned int iRecoEcalCand = 0; iRecoEcalCand < recoecalcandHandle->size(); iRecoEcalCand++) { reco::RecoEcalCandidateRef recoecalcandref(recoecalcandHandle, iRecoEcalCand); if (recoecalcandref->superCluster()->seed()->seed().det() != DetId::Ecal) { //HGCAL, skip for now @@ -122,6 +125,8 @@ void EgammaHLTClusterShapeProducer::produce(edm::StreamID sid, clshSMajorMap.insert(recoecalcandref, 0); clshSMinorMap.insert(recoecalcandref, 0); + e2x2Map.insert(recoecalcandref, 0); + continue; } @@ -161,6 +166,9 @@ void EgammaHLTClusterShapeProducer::produce(edm::StreamID sid, float sMin = moments.sMin; clshSMajorMap.insert(recoecalcandref, sMaj); clshSMinorMap.insert(recoecalcandref, sMin); + + auto const e2x2 = lazyTools.e2x2(*(recoecalcandref->superCluster()->seed())); + e2x2Map.insert(recoecalcandref, e2x2); } iEvent.put(std::make_unique(clshMap)); @@ -175,6 +183,8 @@ void EgammaHLTClusterShapeProducer::produce(edm::StreamID sid, iEvent.put(std::make_unique(clshSMajorMap), "sMajor"); iEvent.put(std::make_unique(clshSMinorMap), "sMinor"); + + iEvent.put(std::make_unique(e2x2Map), "e2x2"); } #include "FWCore/Framework/interface/MakerMacros.h" diff --git a/RecoEgamma/PhotonIdentification/BuildFile.xml b/RecoEgamma/PhotonIdentification/BuildFile.xml index 37b3d83765254..5e23ba0376cda 100644 --- a/RecoEgamma/PhotonIdentification/BuildFile.xml +++ b/RecoEgamma/PhotonIdentification/BuildFile.xml @@ -22,6 +22,7 @@ + diff --git a/RecoEgamma/PhotonIdentification/interface/PhotonXGBoostEstimator.h b/RecoEgamma/PhotonIdentification/interface/PhotonXGBoostEstimator.h new file mode 100644 index 0000000000000..fbacb9e5729aa --- /dev/null +++ b/RecoEgamma/PhotonIdentification/interface/PhotonXGBoostEstimator.h @@ -0,0 +1,28 @@ +#ifndef ReciEgamma_PhotonIdentification_PhotonXGBoostEstimator_h +#define ReciEgamma_PhotonIdentification_PhotonXGBoostEstimator_h + +#include "FWCore/ParameterSet/interface/FileInPath.h" +#include "xgboost/c_api.h" + +class PhotonXGBoostEstimator { +public: + PhotonXGBoostEstimator(const edm::FileInPath& weightsFile, int best_ntree_limit); + ~PhotonXGBoostEstimator(); + + float computeMva(float rawEnergyIn, + float r9In, + float sigmaIEtaIEtaIn, + float etaWidthIn, + float phiWidthIn, + float e2x2In, + float etaIn, + float hOvrEIn, + float ecalPFIsoIn) const; + +private: + BoosterHandle booster_; + int best_ntree_limit_ = -1; + std::string config_; +}; + +#endif diff --git a/RecoEgamma/PhotonIdentification/plugins/PhotonXGBoostProducer.cc b/RecoEgamma/PhotonIdentification/plugins/PhotonXGBoostProducer.cc new file mode 100644 index 0000000000000..720b50fd1ee65 --- /dev/null +++ b/RecoEgamma/PhotonIdentification/plugins/PhotonXGBoostProducer.cc @@ -0,0 +1,137 @@ +#include "DataFormats/Common/interface/AssociationMap.h" +#include "DataFormats/EgammaReco/interface/SuperCluster.h" +#include "DataFormats/EgammaReco/interface/SuperClusterFwd.h" +#include "DataFormats/RecoCandidate/interface/RecoEcalCandidate.h" +#include "DataFormats/RecoCandidate/interface/RecoEcalCandidateIsolation.h" +#include "FWCore/Framework/interface/global/EDProducer.h" +#include "FWCore/Framework/interface/one/EDProducer.h" +#include "FWCore/Framework/interface/Event.h" +#include "FWCore/MessageLogger/interface/MessageLogger.h" +#include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h" +#include "FWCore/ParameterSet/interface/FileInPath.h" +#include "FWCore/ParameterSet/interface/ParameterSet.h" +#include "FWCore/ParameterSet/interface/ParameterSetDescription.h" +#include "FWCore/Utilities/interface/InputTag.h" +#include "RecoEgamma/PhotonIdentification/interface/PhotonXGBoostEstimator.h" + +#include +#include + +class PhotonXGBoostProducer : public edm::global::EDProducer<> { +public: + explicit PhotonXGBoostProducer(edm::ParameterSet const&); + ~PhotonXGBoostProducer() = default; + + static void fillDescriptions(edm::ConfigurationDescriptions& descriptions); + +private: + void produce(edm::StreamID, edm::Event&, edm::EventSetup const&) const override; + + const edm::EDGetTokenT candToken_; + const edm::EDGetTokenT tokenR9_; + const edm::EDGetTokenT tokenHoE_; + const edm::EDGetTokenT tokenSigmaiEtaiEta_; + const edm::EDGetTokenT tokenE2x2_; + const edm::EDGetTokenT tokenIso_; + const edm::FileInPath mvaFileXgbB_; + const edm::FileInPath mvaFileXgbE_; + const unsigned mvaNTreeLimitB_; + const unsigned mvaNTreeLimitE_; + const double mvaThresholdEt_; + std::unique_ptr mvaEstimatorB_; + std::unique_ptr mvaEstimatorE_; +}; + +PhotonXGBoostProducer::PhotonXGBoostProducer(edm::ParameterSet const& config) + : candToken_(consumes(config.getParameter("candTag"))), + tokenR9_(consumes(config.getParameter("inputTagR9"))), + tokenHoE_(consumes(config.getParameter("inputTagHoE"))), + tokenSigmaiEtaiEta_( + consumes(config.getParameter("inputTagSigmaiEtaiEta"))), + tokenE2x2_(consumes(config.getParameter("inputTagE2x2"))), + tokenIso_(consumes(config.getParameter("inputTagIso"))), + mvaFileXgbB_(config.getParameter("mvaFileXgbB")), + mvaFileXgbE_(config.getParameter("mvaFileXgbE")), + mvaNTreeLimitB_(config.getParameter("mvaNTreeLimitB")), + mvaNTreeLimitE_(config.getParameter("mvaNTreeLimitE")), + mvaThresholdEt_(config.getParameter("mvaThresholdEt")) { + mvaEstimatorB_ = std::make_unique(mvaFileXgbB_, mvaNTreeLimitB_); + mvaEstimatorE_ = std::make_unique(mvaFileXgbE_, mvaNTreeLimitE_); + produces(); +} + +void PhotonXGBoostProducer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) { + edm::ParameterSetDescription desc; + desc.add("candTag", edm::InputTag("hltEgammaCandidatesUnseeded")); + desc.add("inputTagR9", edm::InputTag("hltEgammaR9IDUnseeded", "r95x5")); + desc.add("inputTagHoE", edm::InputTag("hltEgammaHoverEUnseeded")); + desc.add("inputTagSigmaiEtaiEta", + edm::InputTag("hltEgammaClusterShapeUnseeded", "sigmaIEtaIEta5x5NoiseCleaned")); + desc.add("inputTagE2x2", edm::InputTag("hltEgammaClusterShapeUnseeded", "e2x2")); + desc.add("inputTagIso", edm::InputTag("hltEgammaEcalPFClusterIsoUnseeded")); + desc.add("mvaFileXgbB", + edm::FileInPath("RecoEgamma/PhotonIdentification/data/xgb_photonmva_barrel_v1.bin")); + desc.add("mvaFileXgbE", + edm::FileInPath("RecoEgamma/PhotonIdentification/data/xgb_photonmva_endcap_v1.bin")); + desc.add("mvaNTreeLimitB", 55); + desc.add("mvaNTreeLimitE", 48); + desc.add("mvaThresholdEt", 0); + descriptions.addWithDefaultLabel(desc); +} + +void PhotonXGBoostProducer::produce(edm::StreamID, edm::Event& event, edm::EventSetup const& setup) const { + const auto& recCollection = event.getHandle(candToken_); + + //get hold of r9 association map + const auto& r9Map = event.getHandle(tokenR9_); + + //get hold of HoE association map + const auto& hoEMap = event.getHandle(tokenHoE_); + + //get hold of isolated association map + const auto& sigmaiEtaiEtaMap = event.getHandle(tokenSigmaiEtaiEta_); + + //get hold of e2x2 (s4) association map + const auto& e2x2Map = event.getHandle(tokenE2x2_); + + //get hold of Ecal isolation association map + const auto& isoMap = event.getHandle(tokenIso_); + + //output + reco::RecoEcalCandidateIsolationMap mvaScoreMap(recCollection); + + for (size_t i = 0; i < recCollection->size(); i++) { + edm::Ref ref(recCollection, i); + + float etaSC = ref->eta(); + + float scEnergy = ref->superCluster()->energy(); + float r9 = (*r9Map).find(ref)->val; + float hoe = (*hoEMap).find(ref)->val / scEnergy; + float siEtaiEta = (*sigmaiEtaiEtaMap).find(ref)->val; + float e2x2 = (*e2x2Map).find(ref)->val; + float iso = (*isoMap).find(ref)->val; + + float rawEnergy = ref->superCluster()->rawEnergy(); + float etaW = ref->superCluster()->etaWidth(); + float phiW = ref->superCluster()->phiWidth(); + + float scEt = scEnergy / cosh(etaSC); + if (scEt < 0.) + scEt = 0.; /* first and second order terms assume non-negative energies */ + + float xgbScore = -100.; + //compute only above threshold used for training and cand filter, else store negative value into the association map. + if (scEt >= mvaThresholdEt_) { + if (std::abs(etaSC) < 1.5) + xgbScore = mvaEstimatorB_->computeMva(rawEnergy, r9, siEtaiEta, etaW, phiW, e2x2, etaSC, hoe, iso); + else + xgbScore = mvaEstimatorE_->computeMva(rawEnergy, r9, siEtaiEta, etaW, phiW, e2x2, etaSC, hoe, iso); + } + mvaScoreMap.insert(ref, xgbScore); + } + event.put(std::make_unique(mvaScoreMap)); +} + +#include "FWCore/Framework/interface/MakerMacros.h" +DEFINE_FWK_MODULE(PhotonXGBoostProducer); diff --git a/RecoEgamma/PhotonIdentification/src/PhotonXGBoostEstimator.cc b/RecoEgamma/PhotonIdentification/src/PhotonXGBoostEstimator.cc new file mode 100644 index 0000000000000..4cc34b3f3077a --- /dev/null +++ b/RecoEgamma/PhotonIdentification/src/PhotonXGBoostEstimator.cc @@ -0,0 +1,60 @@ +#include "RecoEgamma/PhotonIdentification/interface/PhotonXGBoostEstimator.h" +#include + +PhotonXGBoostEstimator::PhotonXGBoostEstimator(const edm::FileInPath& weightsFile, int best_ntree_limit) { + XGBoosterCreate(NULL, 0, &booster_); + XGBoosterLoadModel(booster_, weightsFile.fullPath().c_str()); + best_ntree_limit_ = best_ntree_limit; + + std::stringstream config; + config << "{\"training\": false, \"type\": 0, \"iteration_begin\": 0, \"iteration_end\": " << best_ntree_limit_ + << ", \"strict_shape\": false}"; + config_ = config.str(); +} + +PhotonXGBoostEstimator::~PhotonXGBoostEstimator() { XGBoosterFree(booster_); } + +namespace { + enum inputIndexes { + rawEnergy = 0, // 0 + r9 = 1, // 1 + sigmaIEtaIEta = 2, // 2 + etaWidth = 3, // 3 + phiWidth = 4, // 4 + e2x2 = 5, // 5 + eta = 6, // 6 + hOvrE = 7, // 7 + ecalPFIso = 8, // 8 + }; +} // namespace + +float PhotonXGBoostEstimator::computeMva(float rawEnergyIn, + float r9In, + float sigmaIEtaIEtaIn, + float etaWidthIn, + float phiWidthIn, + float e2x2In, + float etaIn, + float hOvrEIn, + float ecalPFIsoIn) const { + float var[9]; + var[rawEnergy] = rawEnergyIn; + var[r9] = r9In; + var[sigmaIEtaIEta] = sigmaIEtaIEtaIn; + var[etaWidth] = etaWidthIn; + var[phiWidth] = phiWidthIn; + var[e2x2] = e2x2In; + var[eta] = etaIn; + var[hOvrE] = hOvrEIn; + var[ecalPFIso] = ecalPFIsoIn; + + DMatrixHandle dmat; + XGDMatrixCreateFromMat(var, 1, 9, -999.9f, &dmat); + uint64_t const* out_shape; + uint64_t out_dim; + const float* out_result = NULL; + XGBoosterPredictFromDMatrix(booster_, dmat, config_.c_str(), &out_shape, &out_dim, &out_result); + float ret = out_result[0]; + XGDMatrixFree(dmat); + return ret; +}