#include <PuppiAlgo.h>

Inheritance diagram for l1tpf_impl::PuppiAlgo:

Public Member Functions
void	doPUGlobals (const std::vector< Region > &rs, float npu, std::vector< float > &globals) const override

const std::vector< std::string > &	puGlobalNames () const override

	PuppiAlgo (const edm::ParameterSet &)

void	runNeutralsPU (Region &r, float npu, const std::vector< float > &globals) const override

	~PuppiAlgo () override

Public Member Functions inherited from l1tpf_impl::PUAlgoBase
virtual void	doVertexing (std::vector< Region > &rs, VertexAlgo algo, float &vz) const

	PUAlgoBase (const edm::ParameterSet &)

virtual void	runChargedPV (Region &r, float z0) const

virtual	~PUAlgoBase ()

Protected Member Functions
virtual void	computePuppiAlphas (const Region &r, std::vector< float > &alphaC, std::vector< float > &alphaF) const

virtual void	computePuppiMedRMS (const std::vector< Region > &rs, float &alphaCMed, float &alphaCRms, float &alphaFMed, float &alphaFRms) const

void	computePuppiWeights (Region &r, const std::vector< float > &alphaC, const std::vector< float > &alphaF, float alphaCMed, float alphaCRms, float alphaFMed, float alphaFRms) const

virtual void	fillPuppi (Region &r) const

Protected Attributes
std::vector< int16_t >	intPuppiEtaCuts_

std::vector< int16_t >	intPuppiPtCuts_

std::vector< int16_t >	intPuppiPtCutsPhotons_

float	puppiDr_

float	puppiDrMin_

std::vector< float >	puppiEtaCuts_

std::vector< float >	puppiPtCuts_

std::vector< float >	puppiPtCutsPhotons_

float	puppiPtMax_

bool	puppiUsingBareTracks_

Protected Attributes inherited from l1tpf_impl::PUAlgoBase
int	debug_

float	etaCharged_

bool	vtxAdaptiveCut_

float	vtxRes_

Additional Inherited Members
Public Types inherited from l1tpf_impl::PUAlgoBase
enum	VertexAlgo { VertexAlgo::Old, VertexAlgo::TP, VertexAlgo::External }
	global operations More...

Detailed Description

Definition at line 8 of file PuppiAlgo.h.

Constructor & Destructor Documentation

PuppiAlgo::PuppiAlgo ( const edm::ParameterSet & iConfig )

Definition at line 19 of file PuppiAlgo.cc.

References l1tpf_impl::PUAlgoBase::debug_, Exception, edm::ParameterSet::getUntrackedParameter(), mps_fire::i, intPuppiEtaCuts_, intPuppiPtCuts_, intPuppiPtCutsPhotons_, dqmiodumpmetadata::n, puppiEtaCuts_, puppiPtCuts_, and puppiPtCutsPhotons_.

     : PUAlgoBase(iConfig),
       puppiDr_(iConfig.getParameter<double>("puppiDr")),
       puppiDrMin_(iConfig.getParameter<double>("puppiDrMin")),
       puppiPtMax_(iConfig.getParameter<double>("puppiPtMax")),
       puppiEtaCuts_(vd2vf(iConfig.getParameter<std::vector<double>>("puppiEtaCuts"))),
       puppiPtCuts_(vd2vf(iConfig.getParameter<std::vector<double>>("puppiPtCuts"))),
       puppiPtCutsPhotons_(vd2vf(iConfig.getParameter<std::vector<double>>("puppiPtCutsPhotons"))),
       puppiUsingBareTracks_(iConfig.getParameter<bool>("puppiUsingBareTracks")) {
   debug_ = iConfig.getUntrackedParameter<int>("puppiDebug", debug_);
   if (puppiEtaCuts_.size() != puppiPtCuts_.size() || puppiPtCuts_.size() != puppiPtCutsPhotons_.size()) {
     throw cms::Exception("Configuration", "Bad PUPPI config");
   }
   for (unsigned int i = 0, n = puppiEtaCuts_.size(); i < n; ++i) {
     intPuppiEtaCuts_.push_back(std::round(puppiEtaCuts_[i] * CaloCluster::ETAPHI_SCALE));
     intPuppiPtCuts_.push_back(std::round(puppiPtCuts_[i] * CaloCluster::PT_SCALE));
     intPuppiPtCutsPhotons_.push_back(std::round(puppiPtCutsPhotons_[i] * CaloCluster::PT_SCALE));
   }
 }

PuppiAlgo::~PuppiAlgo ( )

override

Definition at line 39 of file PuppiAlgo.cc.

39 {}

Member Function Documentation

void PuppiAlgo::computePuppiAlphas	(	const Region &	r,
		std::vector< float > &	alphaC,
		std::vector< float > &	alphaF
	)		const

protectedvirtual

Definition at line 57 of file PuppiAlgo.cc.

References l1tpf_impl::PFParticle::chargedPV, reco::deltaR2(), l1tpf_impl::PropagatedTrack::fromPV, l1tpf_impl::PFParticle::hwId, min(), np, AlCaHLTBitMon_ParallelJobs::p, fireworks::p2, l1tpf_impl::Region::pf, funct::pow(), puppiDr_, puppiDrMin_, puppiPtMax_, puppiUsingBareTracks_, l1tpf_impl::InputRegion::track, and w.

Referenced by runNeutralsPU(), and l1tpf_impl::LinearizedPuppiAlgo::runNeutralsPU().

                                                                                                             {
   alphaC.resize(r.pf.size());
   alphaF.resize(r.pf.size());
   float puppiDr2 = std::pow(puppiDr_, 2), puppiDr2min = std::pow(puppiDrMin_, 2);
   for (unsigned int ip = 0, np = r.pf.size(); ip < np; ++ip) {
     const PFParticle &p = r.pf[ip];
     if (p.hwId <= 1)
       continue;
     // neutral
     alphaC[ip] = 0;
     alphaF[ip] = 0;
     for (const PFParticle &p2 : r.pf) {
       float dr2 = ::deltaR2(p.floatEta(), p.floatPhi(), p2.floatEta(), p2.floatPhi());
       if (dr2 > 0 && dr2 < puppiDr2) {
         float w = std::pow(std::min(p2.floatPt(), puppiPtMax_), 2) / std::max<float>(puppiDr2min, dr2);
         alphaF[ip] += w;
         if (p2.chargedPV)
           alphaC[ip] += w;
       }
     }
     if (puppiUsingBareTracks_) {
       alphaC[ip] = 0;
       for (const PropagatedTrack &p2 : r.track) {
         if (!p2.fromPV)
           continue;
         float dr2 = ::deltaR2(p.floatEta(), p.floatPhi(), p2.floatEta(), p2.floatPhi());
         if (dr2 > 0 && dr2 < puppiDr2) {
           alphaC[ip] += std::pow(std::min(p2.floatPt(), puppiPtMax_), 2) / std::max<float>(puppiDr2min, dr2);
         }
       }
     }
   }
 }

void PuppiAlgo::computePuppiMedRMS	(	const std::vector< Region > &	rs,
		float &	alphaCMed,
		float &	alphaCRms,
		float &	alphaFMed,
		float &	alphaFRms
	)		const

protectedvirtual

Definition at line 157 of file PuppiAlgo.cc.

References funct::abs(), alpha, central, l1tpf_impl::PFParticle::chargedPV, dbgPrintf(), l1tpf_impl::PUAlgoBase::debug_, reco::deltaR2(), l1tpf_impl::PUAlgoBase::etaCharged_, l1tpf_impl::PropagatedTrack::fromPV, l1tpf_impl::PFParticle::hwEta, l1tpf_impl::PFParticle::hwId, log, min(), AlCaHLTBitMon_ParallelJobs::p, fireworks::p2, funct::pow(), puppiDr_, puppiDrMin_, puppiPtMax_, puppiUsingBareTracks_, alignCSCRings::r, mathSSE::sqrt(), and w.

Referenced by doPUGlobals().

                                                                                                                {
   std::vector<float> alphaFs;
   std::vector<float> alphaCs;
   int16_t ietacut = std::round(etaCharged_ * CaloCluster::ETAPHI_SCALE);
   float puppiDr2 = std::pow(puppiDr_, 2), puppiDr2min = std::pow(puppiDrMin_, 2);
   for (const Region &r : rs) {
     for (const PFParticle &p : r.pf) {
       bool central = std::abs(p.hwEta) < ietacut;
       if (r.relativeCoordinates)
         central = (r.globalAbsEta(p.floatEta()) < etaCharged_);  // FIXME could make a better integer implementation
       if (central) {
         if (p.hwId > 1 || p.chargedPV)
           continue;
       }
       float alphaC = 0, alphaF = 0;
       for (const PFParticle &p2 : r.pf) {
         float dr2 = ::deltaR2(p.floatEta(), p.floatPhi(), p2.floatEta(), p2.floatPhi());
         if (dr2 > 0 && dr2 < puppiDr2) {
           float w = std::pow(std::min(p2.floatPt(), puppiPtMax_), 2) / std::max<float>(puppiDr2min, dr2);
           alphaF += w;
           if (p2.chargedPV)
             alphaC += w;
         }
       }
       if (puppiUsingBareTracks_) {
         alphaC = 0;
         for (const PropagatedTrack &p2 : r.track) {
           if (!p2.fromPV)
             continue;
           float dr2 = ::deltaR2(p.floatEta(), p.floatPhi(), p2.floatEta(), p2.floatPhi());
           if (dr2 > 0 && dr2 < puppiDr2) {
             alphaC += std::pow(std::min(p2.floatPt(), puppiPtMax_), 2) / std::max<float>(puppiDr2min, dr2);
           }
         }
       }
       if (central) {
         if (alphaC > 0)
           alphaCs.push_back(std::log(alphaC));
       } else {
         if (alphaF > 0)
           alphaFs.push_back(std::log(alphaF));
       }
     }
   }
   std::sort(alphaCs.begin(), alphaCs.end());
   std::sort(alphaFs.begin(), alphaFs.end());
 
   if (alphaCs.size() > 1) {
     alphaCMed = alphaCs[alphaCs.size() / 2 + 1];
     double sum = 0.0;
     for (float alpha : alphaCs)
       sum += std::pow(alpha - alphaCMed, 2);
     alphaCRms = std::sqrt(float(sum) / alphaCs.size());
   } else {
     alphaCMed = 8.;
     alphaCRms = 8.;
   }
 
   if (alphaFs.size() > 1) {
     alphaFMed = alphaFs[alphaFs.size() / 2 + 1];
     double sum = 0.0;
     for (float alpha : alphaFs)
       sum += std::pow(alpha - alphaFMed, 2);
     alphaFRms = std::sqrt(float(sum) / alphaFs.size());
   } else {
     alphaFMed = 6.;
     alphaFRms = 6.;
   }
   if (debug_)
     dbgPrintf("PUPPI \t alphaC = %+6.2f +- %6.2f (%4lu), alphaF = %+6.2f +- %6.2f (%4lu)\n",
               alphaCMed,
               alphaCRms,
               alphaCs.size(),
               alphaFMed,
               alphaFRms,
               alphaFs.size());
 }

void PuppiAlgo::computePuppiWeights	(	Region &	r,
		const std::vector< float > &	alphaC,
		const std::vector< float > &	alphaF,
		float	alphaCMed,
		float	alphaCRms,
		float	alphaFMed,
		float	alphaFRms
	)		const

protected

Definition at line 91 of file PuppiAlgo.cc.

References funct::abs(), alpha, central, l1t::PFCandidate::ChargedHadron, l1tpf_impl::PFParticle::chargedPV, dbgPrintf(), l1tpf_impl::PUAlgoBase::debug_, l1t::PFCandidate::Electron, l1tpf_impl::PUAlgoBase::etaCharged_, l1tpf_impl::Region::globalAbsEta(), l1tpf_impl::PFParticle::hwEta, l1tpf_impl::PFParticle::hwId, log, l1t::PFCandidate::Muon, np, AlCaHLTBitMon_ParallelJobs::p, l1tpf_impl::Region::pf, funct::pow(), and l1tpf_impl::Region::relativeCoordinates.

Referenced by runNeutralsPU().

                                                            {
   int16_t ietacut = std::round(etaCharged_ * CaloCluster::ETAPHI_SCALE);
   for (unsigned int ip = 0, np = r.pf.size(); ip < np; ++ip) {
     PFParticle &p = r.pf[ip];
     // charged
     if (p.hwId == l1t::PFCandidate::ChargedHadron || p.hwId == l1t::PFCandidate::Electron ||
         p.hwId == l1t::PFCandidate::Muon) {
       p.setPuppiW(p.chargedPV || p.hwId == l1t::PFCandidate::Muon ? 1.0 : 0);
       if (debug_)
         dbgPrintf(
             "PUPPI \t charged id %1d pt %7.2f eta %+5.2f phi %+5.2f  alpha %+7.2f x2 %+7.2f --> puppi weight %.3f   "
             "puppi pt %7.2f \n",
             p.hwId,
             p.floatPt(),
             p.floatEta(),
             p.floatPhi(),
             0.,
             0.,
             p.floatPuppiW(),
             p.floatPt() * p.floatPuppiW());
       continue;
     }
     // neutral
     float alpha = -99, x2 = -99;
     bool central = std::abs(p.hwEta) < ietacut;
     if (r.relativeCoordinates)
       central =
           (std::abs(r.globalAbsEta(p.floatEta())) < etaCharged_);  // FIXME could make a better integer implementation
     if (central) {
       if (alphaC[ip] > 0) {
         alpha = std::log(alphaC[ip]);
         x2 = (alpha - alphaCMed) * std::abs(alpha - alphaCMed) / std::pow(alphaCRms, 2);
         p.setPuppiW(ROOT::Math::chisquared_cdf(x2, 1));
       } else {
         p.setPuppiW(0);
       }
     } else {
       if (alphaF[ip] > 0) {
         alpha = std::log(alphaF[ip]);
         x2 = (alpha - alphaFMed) * std::abs(alpha - alphaFMed) / std::pow(alphaFRms, 2);
         p.setPuppiW(ROOT::Math::chisquared_cdf(x2, 1));
       } else {
         p.setPuppiW(0);
       }
     }
     if (debug_)
       dbgPrintf(
           "PUPPI \t neutral id %1d pt %7.2f eta %+5.2f phi %+5.2f  alpha %+7.2f x2 %+7.2f --> puppi weight %.3f   "
           "puppi pt %7.2f \n",
           p.hwId,
           p.floatPt(),
           p.floatEta(),
           p.floatPhi(),
           alpha,
           x2,
           p.floatPuppiW(),
           p.floatPt() * p.floatPuppiW());
   }
 }

void PuppiAlgo::doPUGlobals	(	const std::vector< Region > &	rs,
		float	npu,
		std::vector< float > &	globals
	)		const

overridevirtual

Implements l1tpf_impl::PUAlgoBase.

Definition at line 45 of file PuppiAlgo.cc.

References computePuppiMedRMS().

                                                                                                    {
   globals.resize(4);
   computePuppiMedRMS(rs, globals[0], globals[1], globals[2], globals[3]);
 }

void PuppiAlgo::fillPuppi ( Region & r ) const

protectedvirtual

Definition at line 236 of file PuppiAlgo.cc.

References funct::abs(), l1t::PFCandidate::ChargedHadron, l1t::PFCandidate::Electron, l1tpf_impl::Region::globalAbsEta(), l1tpf_impl::PFParticle::hwEta, l1tpf_impl::PFParticle::hwId, l1tpf_impl::PFParticle::hwPt, l1tpf_impl::PFParticle::hwPuppiWeight, intPuppiEtaCuts_, intPuppiPtCuts_, intPuppiPtCutsPhotons_, l1t::PFCandidate::Muon, AlCaHLTBitMon_ParallelJobs::p, l1tpf_impl::Region::pf, l1t::PFCandidate::Photon, l1tpf_impl::Region::puppi, and l1tpf_impl::Region::relativeCoordinates.

Referenced by runNeutralsPU(), and l1tpf_impl::LinearizedPuppiAlgo::runNeutralsPU().

                                          {
   uint16_t PUPPIW_0p01 = std::round(0.01 * PFParticle::PUPPI_SCALE);
   r.puppi.clear();
   for (PFParticle &p : r.pf) {
     if (p.hwId == l1t::PFCandidate::ChargedHadron || p.hwId == l1t::PFCandidate::Electron ||
         p.hwId == l1t::PFCandidate::Muon) {  // charged
       if (p.hwPuppiWeight > 0) {
         r.puppi.push_back(p);
       }
     } else {  // neutral
       if (p.hwPuppiWeight > PUPPIW_0p01) {
         // FIXME would work better with PUPPI_SCALE being a power of two, to do the shift
         // FIXME done with floats
         int16_t hwPt = (float(p.hwPt) * float(p.hwPuppiWeight) / float(PFParticle::PUPPI_SCALE));
         int16_t hwPtCut = 0, hwAbsEta = r.relativeCoordinates
                                             ? round(r.globalAbsEta(p.floatEta()) * CaloCluster::ETAPHI_SCALE)
                                             : std::abs(p.hwEta);
         for (unsigned int ietaBin = 0, nBins = intPuppiEtaCuts_.size(); ietaBin < nBins; ++ietaBin) {
           if (hwAbsEta < intPuppiEtaCuts_[ietaBin]) {
             hwPtCut = (p.hwId == l1t::PFCandidate::Photon ? intPuppiPtCutsPhotons_[ietaBin] : intPuppiPtCuts_[ietaBin]);
             break;
           }
         }
         if (hwPt > hwPtCut) {
           r.puppi.push_back(p);
           r.puppi.back().hwPt = hwPt;
         }
       }
     }
   }
 }

const std::vector< std::string > & PuppiAlgo::puGlobalNames ( ) const

overridevirtual

Reimplemented from l1tpf_impl::PUAlgoBase.

Definition at line 41 of file PuppiAlgo.cc.

                                                            {
   static const std::vector<std::string> names_{"alphaCMed", "alphaCRms", "alphaFMed", "alphaFRms"};
   return names_;
 }

void PuppiAlgo::runNeutralsPU	(	Region &	r,
		float	npu,
		const std::vector< float > &	globals
	)		const

overridevirtual

Implements l1tpf_impl::PUAlgoBase.

Definition at line 50 of file PuppiAlgo.cc.

References computePuppiAlphas(), computePuppiWeights(), and fillPuppi().

                                                                                          {
   std::vector<float> alphaC, alphaF;
   computePuppiAlphas(r, alphaC, alphaF);
   computePuppiWeights(r, alphaC, alphaF, globals[0], globals[1], globals[2], globals[3]);
   fillPuppi(r);
 }