diff --git a/CommonTools/PileupAlgos/interface/PuppiAlgo.h b/CommonTools/PileupAlgos/interface/PuppiAlgo.h index b3c9ede5eb2c2..7aa3868e2b794 100644 --- a/CommonTools/PileupAlgos/interface/PuppiAlgo.h +++ b/CommonTools/PileupAlgos/interface/PuppiAlgo.h @@ -9,7 +9,7 @@ class PuppiAlgo { public: - PuppiAlgo(edm::ParameterSet &iConfig); + PuppiAlgo(const edm::ParameterSet &iConfig); ~PuppiAlgo(); static void fillDescriptionsPuppiAlgo(edm::ParameterSetDescription &desc); //Computing Mean and RMS @@ -19,7 +19,7 @@ class PuppiAlgo { void computeMedRMS(const unsigned int &iAlgo); //Get the Weight double compute(std::vector const &iVals, double iChi2) const; - const std::vector &alphas() { return fPups; } + const std::vector &alphas() const { return fPups; } //Helpers inline int etaBins() const { return fEtaMin.size(); } inline double etaMin(int i) const { return fEtaMin[i]; } diff --git a/CommonTools/PileupAlgos/interface/PuppiContainer.h b/CommonTools/PileupAlgos/interface/PuppiContainer.h index 479cf542a4585..4ad1dd951a6cd 100644 --- a/CommonTools/PileupAlgos/interface/PuppiContainer.h +++ b/CommonTools/PileupAlgos/interface/PuppiContainer.h @@ -8,51 +8,45 @@ class PuppiContainer { public: PuppiContainer(const edm::ParameterSet &iConfig); - ~PuppiContainer(); - void initialize(const std::vector &iRecoObjects); - void setPUProxy(double const iPUProxy) { fPUProxy = iPUProxy; } - std::vector const &pfParticles() const { return fPFParticles; } - std::vector const &puppiWeights(); - const std::vector &puppiRawAlphas() { return fRawAlphas; } - const std::vector &puppiAlphas() { return fVals; } - // const std::vector puppiAlpha () {return fAlpha;} - const std::vector &puppiAlphasMed() { return fAlphaMed; } - const std::vector &puppiAlphasRMS() { return fAlphaRMS; } + struct Weights { + std::vector weights; + std::vector puppiRawAlphas; + std::vector puppiAlphas; + std::vector puppiAlphasMed; + std::vector puppiAlphasRMS; + }; - int puppiNAlgos() { return fNAlgos; } + Weights calculatePuppiWeights(const std::vector &iRecoObjects, double iPUProxy); -protected: - double goodVar(PuppiCandidate const &iPart, std::vector const &iParts, int iOpt, const double iRCone); + int puppiNAlgos() const { return fPuppiAlgo.size(); } + +private: + void initialize(const std::vector &iRecoObjects, + std::vector &fPFParticles, + std::vector &fPFParticlesForVar, + std::vector &fPFParticlesForVarChargedPV) const; + + double goodVar(PuppiCandidate const &iPart, + std::vector const &iParts, + int iOpt, + const double iRCone) const; void getRMSAvg(int iOpt, std::vector const &iConstits, std::vector const &iParticles, - std::vector const &iChargeParticles); - void getRawAlphas(int iOpt, - std::vector const &iConstits, - std::vector const &iParticles, - std::vector const &iChargeParticles); - double getChi2FromdZ(double iDZ); + std::vector const &iChargeParticles, + std::vector &oVals); + std::vector getRawAlphas(int iOpt, + std::vector const &iConstits, + std::vector const &iParticles, + std::vector const &iChargeParticles) const; + double getChi2FromdZ(double iDZ) const; int getPuppiId(float iPt, float iEta); double var_within_R(int iId, const std::vector &particles, const PuppiCandidate ¢re, - const double R); - - bool fPuppiDiagnostics; - const std::vector *fRecoParticles; - std::vector fPFParticles; - std::vector fPFParticlesForVar; - std::vector fPFParticlesForVarChargedPV; - std::vector fWeights; - std::vector fVals; - std::vector fRawAlphas; - std::vector fAlphaMed; - std::vector fAlphaRMS; + const double R) const; - bool fApplyCHS; - bool fInvert; - bool fUseExp; double fNeutralMinPt; double fNeutralSlope; double fPuppiWeightCut; @@ -60,8 +54,11 @@ class PuppiContainer { double fEtaMaxPhotons; double fPtMaxNeutrals; double fPtMaxNeutralsStartSlope; - int fNAlgos; - double fPUProxy; std::vector fPuppiAlgo; + + bool fPuppiDiagnostics; + bool fApplyCHS; + bool fInvert; + bool fUseExp; }; #endif diff --git a/CommonTools/PileupAlgos/plugins/PuppiProducer.cc b/CommonTools/PileupAlgos/plugins/PuppiProducer.cc index 2ec716126bf75..165ed313dfb8a 100644 --- a/CommonTools/PileupAlgos/plugins/PuppiProducer.cc +++ b/CommonTools/PileupAlgos/plugins/PuppiProducer.cc @@ -25,7 +25,6 @@ class PuppiProducer : public edm::stream::EDProducer<> { public: explicit PuppiProducer(const edm::ParameterSet&); - ~PuppiProducer() override; static void fillDescriptions(edm::ConfigurationDescriptions& descriptions); typedef math::XYZTLorentzVector LorentzVector; @@ -39,9 +38,7 @@ class PuppiProducer : public edm::stream::EDProducer<> { typedef edm::Association CandToVertex; private: - virtual void beginJob(); void produce(edm::Event&, const edm::EventSetup&) override; - virtual void endJob(); edm::EDGetTokenT tokenPFCandidates_; edm::EDGetTokenT tokenVertices_; @@ -87,12 +84,11 @@ class PuppiProducer : public edm::stream::EDProducer<> { int fVtxNdofCut; double fVtxZCut; bool fUsePUProxyValue; - std::unique_ptr fPuppiContainer; - std::vector fRecoObjCollection; + PuppiContainer fPuppiContainer; }; // ------------------------------------------------------------------------------------------ -PuppiProducer::PuppiProducer(const edm::ParameterSet& iConfig) { +PuppiProducer::PuppiProducer(const edm::ParameterSet& iConfig) : fPuppiContainer(iConfig) { fPuppiDiagnostics = iConfig.getParameter("puppiDiagnostics"); fPuppiNoLep = iConfig.getParameter("puppiNoLep"); fUseFromPVLooseTight = iConfig.getParameter("UseFromPVLooseTight"); @@ -113,7 +109,6 @@ PuppiProducer::PuppiProducer(const edm::ParameterSet& iConfig) { fClonePackedCands = iConfig.getParameter("clonePackedCands"); fVtxNdofCut = iConfig.getParameter("vtxNdofCut"); fVtxZCut = iConfig.getParameter("vtxZCut"); - fPuppiContainer = std::make_unique(iConfig); tokenPFCandidates_ = consumes(iConfig.getParameter("candName")); tokenVertices_ = consumes(iConfig.getParameter("vertexName")); @@ -150,8 +145,6 @@ PuppiProducer::PuppiProducer(const edm::ParameterSet& iConfig) { } } // ------------------------------------------------------------------------------------------ -PuppiProducer::~PuppiProducer() {} -// ------------------------------------------------------------------------------------------ void PuppiProducer::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) { // Get PFCandidate Collection edm::Handle hPFProduct; @@ -180,11 +173,13 @@ void PuppiProducer::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) { } } + std::vector recoObjCollection; + + PuppiContainer::Weights weightsInfo; std::vector lWeights; if (!fUseExistingWeights) { //Fill the reco objects - fRecoObjCollection.clear(); - fRecoObjCollection.reserve(pfCol->size()); + recoObjCollection.reserve(pfCol->size()); int iCand = 0; for (auto const& aPF : *pfCol) { RecoObj pReco; @@ -347,15 +342,13 @@ void PuppiProducer::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) { } } - fRecoObjCollection.push_back(pReco); + recoObjCollection.push_back(pReco); iCand++; } - fPuppiContainer->initialize(fRecoObjCollection); - fPuppiContainer->setPUProxy(puProxyValue); - //Compute the weights and get the particles - lWeights = fPuppiContainer->puppiWeights(); + weightsInfo = fPuppiContainer.calculatePuppiWeights(recoObjCollection, puProxyValue); + lWeights = std::move(weightsInfo.weights); } else { //Use the existing weights lWeights.reserve(pfCol->size()); @@ -457,8 +450,8 @@ void PuppiProducer::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) { p4PupFiller.insert(hPFProduct, puppiP4s.begin(), puppiP4s.end()); p4PupFiller.fill(); - iEvent.emplace(ptokenPupOut_, lPupOut); - iEvent.emplace(ptokenP4PupOut_, p4PupOut); + iEvent.emplace(ptokenPupOut_, std::move(lPupOut)); + iEvent.emplace(ptokenP4PupOut_, std::move(p4PupOut)); if (fUseExistingWeights || fClonePackedCands) { edm::OrphanHandle oh = iEvent.emplace(ptokenPackedPuppiCandidates_, fPackedPuppiCandidates); @@ -477,31 +470,22 @@ void PuppiProducer::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) { edm::ValueMap::Filler filler(pfMap_p); filler.insert(hPFProduct, values.begin(), values.end()); filler.fill(); - iEvent.emplace(ptokenValues_, pfMap_p); + iEvent.emplace(ptokenValues_, std::move(pfMap_p)); ////////////////////////////////////////////// if (fPuppiDiagnostics && !fUseExistingWeights) { // all the different alphas per particle // THE alpha per particle - std::vector theAlphas(fPuppiContainer->puppiAlphas()); - std::vector theAlphasMed(fPuppiContainer->puppiAlphasMed()); - std::vector theAlphasRms(fPuppiContainer->puppiAlphasRMS()); - std::vector alphas(fPuppiContainer->puppiRawAlphas()); - double nalgos(fPuppiContainer->puppiNAlgos()); + double nalgos(fPuppiContainer.puppiNAlgos()); - iEvent.emplace(ptokenRawAlphas_, alphas); + iEvent.emplace(ptokenRawAlphas_, std::move(weightsInfo.puppiRawAlphas)); iEvent.emplace(ptokenNalgos_, nalgos); - iEvent.emplace(ptokenAlphas_, theAlphas); - iEvent.emplace(ptokenAlphasMed_, theAlphasMed); - iEvent.emplace(ptokenAlphasRms_, theAlphasRms); + iEvent.emplace(ptokenAlphas_, std::move(weightsInfo.puppiAlphas)); + iEvent.emplace(ptokenAlphasMed_, std::move(weightsInfo.puppiAlphasMed)); + iEvent.emplace(ptokenAlphasRms_, std::move(weightsInfo.puppiAlphasRMS)); } } -// ------------------------------------------------------------------------------------------ -void PuppiProducer::beginJob() {} -// ------------------------------------------------------------------------------------------ -void PuppiProducer::endJob() {} -// ------------------------------------------------------------------------------------------ void PuppiProducer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) { edm::ParameterSetDescription desc; desc.add("puppiDiagnostics", false); diff --git a/CommonTools/PileupAlgos/src/PuppiAlgo.cc b/CommonTools/PileupAlgos/src/PuppiAlgo.cc index d5a9f47e8cd20..b0f79f1df64d8 100644 --- a/CommonTools/PileupAlgos/src/PuppiAlgo.cc +++ b/CommonTools/PileupAlgos/src/PuppiAlgo.cc @@ -5,7 +5,7 @@ #include "Math/ProbFunc.h" #include "TMath.h" -PuppiAlgo::PuppiAlgo(edm::ParameterSet &iConfig) { +PuppiAlgo::PuppiAlgo(edm::ParameterSet const &iConfig) { fEtaMin = iConfig.getParameter>("etaMin"); fEtaMax = iConfig.getParameter>("etaMax"); fPtMin = iConfig.getParameter>("ptMin"); diff --git a/CommonTools/PileupAlgos/src/PuppiContainer.cc b/CommonTools/PileupAlgos/src/PuppiContainer.cc index ceb08185e561a..da1b7a6f44e06 100644 --- a/CommonTools/PileupAlgos/src/PuppiContainer.cc +++ b/CommonTools/PileupAlgos/src/PuppiContainer.cc @@ -20,30 +20,20 @@ PuppiContainer::PuppiContainer(const edm::ParameterSet &iConfig) { fPtMaxNeutrals = iConfig.getParameter("PtMaxNeutrals"); fPtMaxNeutralsStartSlope = iConfig.getParameter("PtMaxNeutralsStartSlope"); std::vector lAlgos = iConfig.getParameter >("algos"); - fNAlgos = lAlgos.size(); - for (unsigned int i0 = 0; i0 < lAlgos.size(); i0++) { - PuppiAlgo pPuppiConfig(lAlgos[i0]); - fPuppiAlgo.push_back(pPuppiConfig); + fPuppiAlgo.reserve(lAlgos.size()); + for (auto const &algos : lAlgos) { + fPuppiAlgo.emplace_back(algos); } } -void PuppiContainer::initialize(const std::vector &iRecoObjects) { - //Clear everything - fPFParticles.resize(0); - fPFParticlesForVar.resize(0); - fPFParticlesForVarChargedPV.resize(0); - fWeights.resize(0); - fVals.resize(0); - fRawAlphas.resize(0); - fAlphaMed.resize(0); - fAlphaRMS.resize(0); - fPUProxy = 1.; - //Link to the RecoObjects - fRecoParticles = &iRecoObjects; - fPFParticles.reserve(iRecoObjects.size()); - fPFParticlesForVar.reserve(iRecoObjects.size()); - fPFParticlesForVarChargedPV.reserve(iRecoObjects.size()); - for (auto const &rParticle : *fRecoParticles) { +void PuppiContainer::initialize(const std::vector &iRecoObjects, + std::vector &pfParticles, + std::vector &pfParticlesForVar, + std::vector &pfParticlesForVarChargedPV) const { + pfParticles.reserve(iRecoObjects.size()); + pfParticlesForVar.reserve(iRecoObjects.size()); + pfParticlesForVarChargedPV.reserve(iRecoObjects.size()); + for (auto const &rParticle : iRecoObjects) { PuppiCandidate pCand; pCand.id = rParticle.id; if (edm::isFinite(rParticle.rapidity)) { @@ -60,33 +50,31 @@ void PuppiContainer::initialize(const std::vector &iRecoObjects) { pCand.m = 0.; } - fPFParticles.push_back(pCand); + pfParticles.push_back(pCand); // skip candidates to be ignored in the computation // of PUPPI's alphas (e.g. electrons and muons if puppiNoLep=True) if (std::abs(rParticle.id) == 3) continue; - fPFParticlesForVar.push_back(pCand); + pfParticlesForVar.push_back(pCand); // charged candidates assigned to LV if (std::abs(rParticle.id) == 1) - fPFParticlesForVarChargedPV.push_back(pCand); + pfParticlesForVarChargedPV.push_back(pCand); } } -PuppiContainer::~PuppiContainer() {} - double PuppiContainer::goodVar(PuppiCandidate const &iPart, std::vector const &iParts, int iOpt, - const double iRCone) { + const double iRCone) const { return var_within_R(iOpt, iParts, iPart, iRCone); } double PuppiContainer::var_within_R(int iId, const vector &particles, const PuppiCandidate ¢re, - const double R) { + const double R) const { if (iId == -1) return 1.; @@ -137,12 +125,13 @@ double PuppiContainer::var_within_R(int iId, void PuppiContainer::getRMSAvg(int iOpt, std::vector const &iConstits, std::vector const &iParticles, - std::vector const &iChargedParticles) { + std::vector const &iChargedParticles, + std::vector &oVals) { for (unsigned int i0 = 0; i0 < iConstits.size(); i0++) { //Calculate the Puppi Algo to use int pPupId = getPuppiId(iConstits[i0].pt, iConstits[i0].eta); if (pPupId == -1 || fPuppiAlgo[pPupId].numAlgos() <= iOpt) { - fVals.push_back(-1); + oVals.push_back(-1); continue; } //Get the Puppi Sub Algo (given iteration) @@ -158,7 +147,7 @@ void PuppiContainer::getRMSAvg(int iOpt, (std::abs(iConstits[i0].eta) < fPuppiAlgo[pPupId].etaMaxExtrap() and getsDefaultWgtIfApplyCHS)) { pVal = goodVar(iConstits[i0], pCharged ? iChargedParticles : iParticles, pAlgo, pCone); } - fVals.push_back(pVal); + oVals.push_back(pVal); if (!edm::isFinite(pVal)) { LogDebug("NotFound") << "====> Value is Nan " << pVal << " == " << iConstits[i0].pt << " -- " << iConstits[i0].eta @@ -167,60 +156,66 @@ void PuppiContainer::getRMSAvg(int iOpt, } // code added by Nhan: now instead for every algorithm give it all the particles - for (int i1 = 0; i1 < fNAlgos; i1++) { + int count = 0; + for (auto &algo : fPuppiAlgo) { + int index = count++; // skip cands outside of algo's etaMaxExtrap, as they would anyway be ignored inside PuppiAlgo::add (see end of the block) - if (not(std::abs(iConstits[i0].eta) < fPuppiAlgo[i1].etaMaxExtrap() and getsDefaultWgtIfApplyCHS)) + if (not(std::abs(iConstits[i0].eta) < algo.etaMaxExtrap() and getsDefaultWgtIfApplyCHS)) continue; auto curVal = pVal; // recompute goodVar if algo has changed - if (i1 != pPupId) { - pAlgo = fPuppiAlgo[i1].algoId(iOpt); - pCharged = fPuppiAlgo[i1].isCharged(iOpt); - pCone = fPuppiAlgo[i1].coneSize(iOpt); + if (index != pPupId) { + pAlgo = algo.algoId(iOpt); + pCharged = algo.isCharged(iOpt); + pCone = algo.coneSize(iOpt); curVal = goodVar(iConstits[i0], pCharged ? iChargedParticles : iParticles, pAlgo, pCone); } - fPuppiAlgo[i1].add(iConstits[i0], curVal, iOpt); + algo.add(iConstits[i0], curVal, iOpt); } } - for (int i0 = 0; i0 < fNAlgos; i0++) - fPuppiAlgo[i0].computeMedRMS(iOpt); + for (auto &algo : fPuppiAlgo) + algo.computeMedRMS(iOpt); } //In fact takes the median not the average -void PuppiContainer::getRawAlphas(int iOpt, - std::vector const &iConstits, - std::vector const &iParticles, - std::vector const &iChargedParticles) { - for (int j0 = 0; j0 < fNAlgos; j0++) { - for (unsigned int i0 = 0; i0 < iConstits.size(); i0++) { +std::vector PuppiContainer::getRawAlphas(int iOpt, + std::vector const &iConstits, + std::vector const &iParticles, + std::vector const &iChargedParticles) const { + std::vector oRawAlphas; + oRawAlphas.reserve(fPuppiAlgo.size() * iConstits.size()); + for (auto &algo : fPuppiAlgo) { + for (auto const &constit : iConstits) { //Get the Puppi Sub Algo (given iteration) - int pAlgo = fPuppiAlgo[j0].algoId(iOpt); - bool pCharged = fPuppiAlgo[j0].isCharged(iOpt); - double pCone = fPuppiAlgo[j0].coneSize(iOpt); + int pAlgo = algo.algoId(iOpt); + bool pCharged = algo.isCharged(iOpt); + double pCone = algo.coneSize(iOpt); //Compute the Puppi Metric - double const pVal = goodVar(iConstits[i0], pCharged ? iChargedParticles : iParticles, pAlgo, pCone); - fRawAlphas.push_back(pVal); + double const pVal = goodVar(constit, pCharged ? iChargedParticles : iParticles, pAlgo, pCone); + oRawAlphas.push_back(pVal); if (!edm::isFinite(pVal)) { - LogDebug("NotFound") << "====> Value is Nan " << pVal << " == " << iConstits[i0].pt << " -- " - << iConstits[i0].eta << endl; + LogDebug("NotFound") << "====> Value is Nan " << pVal << " == " << constit.pt << " -- " << constit.eta << endl; continue; } } } + return oRawAlphas; } int PuppiContainer::getPuppiId(float iPt, float iEta) { int lId = -1; - for (int i0 = 0; i0 < fNAlgos; i0++) { - int nEtaBinsPerAlgo = fPuppiAlgo[i0].etaBins(); + int count = 0; + for (auto &algo : fPuppiAlgo) { + int index = count++; + int nEtaBinsPerAlgo = algo.etaBins(); for (int i1 = 0; i1 < nEtaBinsPerAlgo; i1++) { - if ((std::abs(iEta) >= fPuppiAlgo[i0].etaMin(i1)) && (std::abs(iEta) < fPuppiAlgo[i0].etaMax(i1))) { - fPuppiAlgo[i0].fixAlgoEtaBin(i1); - if (iPt > fPuppiAlgo[i0].ptMin()) { - lId = i0; + if ((std::abs(iEta) >= algo.etaMin(i1)) && (std::abs(iEta) < algo.etaMax(i1))) { + algo.fixAlgoEtaBin(i1); + if (iPt > algo.ptMin()) { + lId = index; break; } } @@ -229,7 +224,7 @@ int PuppiContainer::getPuppiId(float iPt, float iEta) { //if(lId == -1) std::cerr << "Error : Full fiducial range is not defined " << std::endl; return lId; } -double PuppiContainer::getChi2FromdZ(double iDZ) { +double PuppiContainer::getChi2FromdZ(double iDZ) const { //We need to obtain prob of PU + (1-Prob of LV) // Prob(LV) = Gaus(dZ,sigma) where sigma = 1.5mm (its really more like 1mm) //double lProbLV = ROOT::Math::normal_cdf_c(std::abs(iDZ),0.2)*2.; //*2 is to do it double sided @@ -244,39 +239,49 @@ double PuppiContainer::getChi2FromdZ(double iDZ) { lChi2PU *= lChi2PU; return lChi2PU; } -std::vector const &PuppiContainer::puppiWeights() { - int lNParticles = fRecoParticles->size(); +PuppiContainer::Weights PuppiContainer::calculatePuppiWeights(const std::vector &iRecoObjects, + double iPUProxy) { + std::vector pfParticles; + std::vector pfParticlesForVar; + std::vector pfParticlesForVarChargedPV; + + initialize(iRecoObjects, pfParticles, pfParticlesForVar, pfParticlesForVarChargedPV); + + int lNParticles = iRecoObjects.size(); + + Weights returnValue; + returnValue.weights.reserve(lNParticles); + returnValue.puppiAlphas.reserve(lNParticles); - fWeights.clear(); - fWeights.reserve(lNParticles); - fVals.clear(); - fVals.reserve(lNParticles); - for (int i0 = 0; i0 < fNAlgos; i0++) - fPuppiAlgo[i0].reset(); + //guarantee all algos are rest before leaving this function + auto doReset = [this](void *) { + for (auto &algo : fPuppiAlgo) + algo.reset(); + }; + std::unique_ptr guard(&fPuppiAlgo, doReset); int lNMaxAlgo = 1; - for (int i0 = 0; i0 < fNAlgos; i0++) - lNMaxAlgo = std::max(fPuppiAlgo[i0].numAlgos(), lNMaxAlgo); + for (auto &algo : fPuppiAlgo) + lNMaxAlgo = std::max(algo.numAlgos(), lNMaxAlgo); //Run through all compute mean and RMS for (int i0 = 0; i0 < lNMaxAlgo; i0++) { - getRMSAvg(i0, fPFParticles, fPFParticlesForVar, fPFParticlesForVarChargedPV); + getRMSAvg(i0, pfParticles, pfParticlesForVar, pfParticlesForVarChargedPV, returnValue.puppiAlphas); } if (fPuppiDiagnostics) - getRawAlphas(0, fPFParticles, fPFParticlesForVar, fPFParticlesForVarChargedPV); + returnValue.puppiRawAlphas = getRawAlphas(0, pfParticles, pfParticlesForVar, pfParticlesForVarChargedPV); std::vector pVals; pVals.reserve(lNParticles); for (int i0 = 0; i0 < lNParticles; i0++) { //Refresh pVals.clear(); - double pWeight = 1; //Get the Puppi Id and if ill defined move on - const auto &rParticle = (*fRecoParticles)[i0]; + const auto &rParticle = iRecoObjects[i0]; int pPupId = getPuppiId(rParticle.pt, rParticle.eta); if (pPupId == -1) { - fWeights.push_back(0); - fAlphaMed.push_back(-10); - fAlphaRMS.push_back(-10); + returnValue.weights.push_back(0); + returnValue.puppiAlphasMed.push_back(-10); + returnValue.puppiAlphasRMS.push_back(-10); continue; } @@ -292,9 +297,9 @@ std::vector const &PuppiContainer::puppiWeights() { //Fill and compute the PuppiWeight int lNAlgos = fPuppiAlgo[pPupId].numAlgos(); for (int i1 = 0; i1 < lNAlgos; i1++) - pVals.push_back(fVals[lNParticles * i1 + i0]); + pVals.push_back(returnValue.puppiAlphas[lNParticles * i1 + i0]); - pWeight = fPuppiAlgo[pPupId].compute(pVals, pChi2); + double pWeight = fPuppiAlgo[pPupId].compute(pVals, pChi2); //Apply the CHS weights if (rParticle.id == 1 && fApplyCHS) pWeight = 1; @@ -307,34 +312,34 @@ std::vector const &PuppiContainer::puppiWeights() { if (!edm::isFinite(pWeight)) { pWeight = 0.0; LogDebug("PuppiWeightError") << "====> Weight is nan : " << pWeight << " : pt " << rParticle.pt - << " -- eta : " << rParticle.eta << " -- Value" << fVals[i0] + << " -- eta : " << rParticle.eta << " -- Value" << returnValue.puppiAlphas[i0] << " -- id : " << rParticle.id << " -- NAlgos: " << lNAlgos << std::endl; } //Basic Cuts - if (pWeight * fPFParticles[i0].pt < fPuppiAlgo[pPupId].neutralPt(fPUProxy) && rParticle.id == 0) + if (pWeight * pfParticles[i0].pt < fPuppiAlgo[pPupId].neutralPt(iPUProxy) && rParticle.id == 0) pWeight = 0; //threshold cut on the neutral Pt // Protect high pT photons (important for gamma to hadronic recoil balance) - if (fPtMaxPhotons > 0 && rParticle.pdgId == 22 && std::abs(fPFParticles[i0].eta) < fEtaMaxPhotons && - fPFParticles[i0].pt > fPtMaxPhotons) + if (fPtMaxPhotons > 0 && rParticle.pdgId == 22 && std::abs(pfParticles[i0].eta) < fEtaMaxPhotons && + pfParticles[i0].pt > fPtMaxPhotons) pWeight = 1.; // Protect high pT neutrals else if ((fPtMaxNeutrals > 0) && (rParticle.id == 0)) pWeight = std::clamp( - (fPFParticles[i0].pt - fPtMaxNeutralsStartSlope) / (fPtMaxNeutrals - fPtMaxNeutralsStartSlope), pWeight, 1.); + (pfParticles[i0].pt - fPtMaxNeutralsStartSlope) / (fPtMaxNeutrals - fPtMaxNeutralsStartSlope), pWeight, 1.); if (pWeight < fPuppiWeightCut) pWeight = 0; //==> Elminate the low Weight stuff if (fInvert) pWeight = 1. - pWeight; //std::cout << "rParticle.pt = " << rParticle.pt << ", rParticle.charge = " << rParticle.charge << ", rParticle.id = " << rParticle.id << ", weight = " << pWeight << std::endl; - fWeights.push_back(pWeight); - fAlphaMed.push_back(fPuppiAlgo[pPupId].median()); - fAlphaRMS.push_back(fPuppiAlgo[pPupId].rms()); - //Now get rid of the thrown out weights for the particle collection + returnValue.weights.push_back(pWeight); + returnValue.puppiAlphasMed.push_back(fPuppiAlgo[pPupId].median()); + returnValue.puppiAlphasRMS.push_back(fPuppiAlgo[pPupId].rms()); + //Now get rid of the thrown out returnValue.weights for the particle collection // leave these lines in, in case want to move eventually to having no 1-to-1 correspondence between puppi and pf cands // if( std::abs(pWeight) < std::numeric_limits::denorm_min() ) continue; // this line seems not to work like it's supposed to... // if(std::abs(pWeight) <= 0. ) continue; } - return fWeights; + return returnValue; }