CMS 3D CMS Logo

List of all members | Public Member Functions | Protected Member Functions | Protected Attributes
PuppiContainer Class Reference

#include <PuppiContainer.h>

Public Member Functions

void initialize (const std::vector< RecoObj > &iRecoObjects)
 
const std::vector< PuppiCandidate > & pfParticles () const
 
const std::vector< double > & puppiAlphas ()
 
const std::vector< double > & puppiAlphasMed ()
 
const std::vector< double > & puppiAlphasRMS ()
 
 PuppiContainer (const edm::ParameterSet &iConfig)
 
int puppiNAlgos ()
 
const std::vector< double > & puppiRawAlphas ()
 
const std::vector< double > & puppiWeights ()
 
void setNPV (int iNPV)
 
 ~PuppiContainer ()
 

Protected Member Functions

double getChi2FromdZ (double iDZ)
 
int getPuppiId (float iPt, float iEta)
 
void getRawAlphas (int iOpt, std::vector< PuppiCandidate > const &iConstits, std::vector< PuppiCandidate > const &iParticles, std::vector< PuppiCandidate > const &iChargeParticles)
 
void getRMSAvg (int iOpt, std::vector< PuppiCandidate > const &iConstits, std::vector< PuppiCandidate > const &iParticles, std::vector< PuppiCandidate > const &iChargeParticles)
 
double goodVar (PuppiCandidate const &iPart, std::vector< PuppiCandidate > const &iParts, int iOpt, const double iRCone)
 
double var_within_R (int iId, const std::vector< PuppiCandidate > &particles, const PuppiCandidate &centre, const double R)
 

Protected Attributes

std::vector< double > fAlphaMed
 
std::vector< double > fAlphaRMS
 
bool fApplyCHS
 
double fEtaMaxPhotons
 
bool fInvert
 
int fNAlgos
 
double fNeutralMinPt
 
double fNeutralSlope
 
int fNPV
 
std::vector< PuppiCandidatefPFParticles
 
std::vector< PuppiCandidatefPFParticlesForVar
 
std::vector< PuppiCandidatefPFParticlesForVarChargedPV
 
double fPtMaxNeutrals
 
double fPtMaxNeutralsStartSlope
 
double fPtMaxPhotons
 
std::vector< PuppiAlgofPuppiAlgo
 
bool fPuppiDiagnostics
 
double fPuppiWeightCut
 
std::vector< double > fRawAlphas
 
const std::vector< RecoObj > * fRecoParticles
 
bool fUseExp
 
std::vector< double > fVals
 
std::vector< double > fWeights
 

Detailed Description

Definition at line 8 of file PuppiContainer.h.

Constructor & Destructor Documentation

◆ PuppiContainer()

PuppiContainer::PuppiContainer ( const edm::ParameterSet iConfig)

Definition at line 12 of file PuppiContainer.cc.

12  {
13  fPuppiDiagnostics = iConfig.getParameter<bool>("puppiDiagnostics");
14  fApplyCHS = iConfig.getParameter<bool>("applyCHS");
15  fInvert = iConfig.getParameter<bool>("invertPuppi");
16  fUseExp = iConfig.getParameter<bool>("useExp");
17  fPuppiWeightCut = iConfig.getParameter<double>("MinPuppiWeight");
18  fPtMaxPhotons = iConfig.getParameter<double>("PtMaxPhotons");
19  fEtaMaxPhotons = iConfig.getParameter<double>("EtaMaxPhotons");
20  fPtMaxNeutrals = iConfig.getParameter<double>("PtMaxNeutrals");
21  fPtMaxNeutralsStartSlope = iConfig.getParameter<double>("PtMaxNeutralsStartSlope");
22  std::vector<edm::ParameterSet> lAlgos = iConfig.getParameter<std::vector<edm::ParameterSet> >("algos");
23  fNAlgos = lAlgos.size();
24  for (unsigned int i0 = 0; i0 < lAlgos.size(); i0++) {
25  PuppiAlgo pPuppiConfig(lAlgos[i0]);
26  fPuppiAlgo.push_back(pPuppiConfig);
27  }
28 }

References edm::ParameterSet::getParameter().

◆ ~PuppiContainer()

PuppiContainer::~PuppiContainer ( )

Definition at line 78 of file PuppiContainer.cc.

78 {}

Member Function Documentation

◆ getChi2FromdZ()

double PuppiContainer::getChi2FromdZ ( double  iDZ)
protected

Definition at line 225 of file PuppiContainer.cc.

225  {
226  //We need to obtain prob of PU + (1-Prob of LV)
227  // Prob(LV) = Gaus(dZ,sigma) where sigma = 1.5mm (its really more like 1mm)
228  //double lProbLV = ROOT::Math::normal_cdf_c(std::abs(iDZ),0.2)*2.; //*2 is to do it double sided
229  //Take iDZ to be corrected by sigma already
230  double lProbLV = ROOT::Math::normal_cdf_c(std::abs(iDZ), 1.) * 2.; //*2 is to do it double sided
231  double lProbPU = 1 - lProbLV;
232  if (lProbPU <= 0)
233  lProbPU = 1e-16; //Quick Trick to through out infs
234  if (lProbPU >= 0)
235  lProbPU = 1 - 1e-16; //Ditto
236  double lChi2PU = TMath::ChisquareQuantile(lProbPU, 1);
237  lChi2PU *= lChi2PU;
238  return lChi2PU;
239 }

References funct::abs(), and MillePedeFileConverter_cfg::e.

◆ getPuppiId()

int PuppiContainer::getPuppiId ( float  iPt,
float  iEta 
)
protected

Definition at line 208 of file PuppiContainer.cc.

208  {
209  int lId = -1;
210  for (int i0 = 0; i0 < fNAlgos; i0++) {
211  int nEtaBinsPerAlgo = fPuppiAlgo[i0].etaBins();
212  for (int i1 = 0; i1 < nEtaBinsPerAlgo; i1++) {
213  if ((std::abs(iEta) >= fPuppiAlgo[i0].etaMin(i1)) && (std::abs(iEta) < fPuppiAlgo[i0].etaMax(i1))) {
214  fPuppiAlgo[i0].fixAlgoEtaBin(i1);
215  if (iPt > fPuppiAlgo[i0].ptMin()) {
216  lId = i0;
217  break;
218  }
219  }
220  }
221  }
222  //if(lId == -1) std::cerr << "Error : Full fiducial range is not defined " << std::endl;
223  return lId;
224 }

References funct::abs(), ALCARECOTkAlBeamHalo_cff::etaMax, ALCARECOTkAlBeamHalo_cff::etaMin, testProducerWithPsetDescEmpty_cfi::i1, L1TowerCalibrationProducer_cfi::iEta, and ptMin.

◆ getRawAlphas()

void PuppiContainer::getRawAlphas ( int  iOpt,
std::vector< PuppiCandidate > const &  iConstits,
std::vector< PuppiCandidate > const &  iParticles,
std::vector< PuppiCandidate > const &  iChargeParticles 
)
protected

Definition at line 186 of file PuppiContainer.cc.

189  {
190  for (int j0 = 0; j0 < fNAlgos; j0++) {
191  for (unsigned int i0 = 0; i0 < iConstits.size(); i0++) {
192  //Get the Puppi Sub Algo (given iteration)
193  int pAlgo = fPuppiAlgo[j0].algoId(iOpt);
194  bool pCharged = fPuppiAlgo[j0].isCharged(iOpt);
195  double pCone = fPuppiAlgo[j0].coneSize(iOpt);
196  //Compute the Puppi Metric
197  double const pVal = goodVar(iConstits[i0], pCharged ? iChargedParticles : iParticles, pAlgo, pCone);
198  fRawAlphas.push_back(pVal);
199  if (!edm::isFinite(pVal)) {
200  LogDebug("NotFound") << "====> Value is Nan " << pVal << " == " << iConstits[i0].pt << " -- "
201  << iConstits[i0].eta << endl;
202  continue;
203  }
204  }
205  }
206 }

References edm::isFinite(), and LogDebug.

◆ getRMSAvg()

void PuppiContainer::getRMSAvg ( int  iOpt,
std::vector< PuppiCandidate > const &  iConstits,
std::vector< PuppiCandidate > const &  iParticles,
std::vector< PuppiCandidate > const &  iChargeParticles 
)
protected

Definition at line 130 of file PuppiContainer.cc.

133  {
134  for (unsigned int i0 = 0; i0 < iConstits.size(); i0++) {
135  //Calculate the Puppi Algo to use
136  int pPupId = getPuppiId(iConstits[i0].pt, iConstits[i0].eta);
137  if (pPupId == -1 || fPuppiAlgo[pPupId].numAlgos() <= iOpt) {
138  fVals.push_back(-1);
139  continue;
140  }
141  //Get the Puppi Sub Algo (given iteration)
142  int pAlgo = fPuppiAlgo[pPupId].algoId(iOpt);
143  bool pCharged = fPuppiAlgo[pPupId].isCharged(iOpt);
144  double pCone = fPuppiAlgo[pPupId].coneSize(iOpt);
145  // compute the Puppi metric:
146  // - calculate goodVar only for candidates that (1) will not be assigned a predefined weight (e.g 0, 1),
147  // or (2) are required for computations inside puppi-algos (see call to PuppiAlgo::add below)
148  double pVal = -1;
149  bool const getsDefaultWgtIfApplyCHS = iConstits[i0].id == 1 or iConstits[i0].id == 2;
150  if (not((fApplyCHS and getsDefaultWgtIfApplyCHS) or iConstits[i0].id == 3) or
151  (std::abs(iConstits[i0].eta) < fPuppiAlgo[pPupId].etaMaxExtrap() and getsDefaultWgtIfApplyCHS)) {
152  pVal = goodVar(iConstits[i0], pCharged ? iChargedParticles : iParticles, pAlgo, pCone);
153  }
154  fVals.push_back(pVal);
155 
156  if (!edm::isFinite(pVal)) {
157  LogDebug("NotFound") << "====> Value is Nan " << pVal << " == " << iConstits[i0].pt << " -- " << iConstits[i0].eta
158  << endl;
159  continue;
160  }
161 
162  // code added by Nhan: now instead for every algorithm give it all the particles
163  for (int i1 = 0; i1 < fNAlgos; i1++) {
164  // skip cands outside of algo's etaMaxExtrap, as they would anyway be ignored inside PuppiAlgo::add (see end of the block)
165  if (not(std::abs(iConstits[i0].eta) < fPuppiAlgo[i1].etaMaxExtrap() and getsDefaultWgtIfApplyCHS))
166  continue;
167 
168  auto curVal = pVal;
169  // recompute goodVar if algo has changed
170  if (i1 != pPupId) {
171  pAlgo = fPuppiAlgo[i1].algoId(iOpt);
172  pCharged = fPuppiAlgo[i1].isCharged(iOpt);
173  pCone = fPuppiAlgo[i1].coneSize(iOpt);
174  curVal = goodVar(iConstits[i0], pCharged ? iChargedParticles : iParticles, pAlgo, pCone);
175  }
176 
177  fPuppiAlgo[i1].add(iConstits[i0], curVal, iOpt);
178  }
179  }
180 
181  for (int i0 = 0; i0 < fNAlgos; i0++)
182  fPuppiAlgo[i0].computeMedRMS(iOpt);
183 }

References funct::abs(), PVValHelper::eta, testProducerWithPsetDescEmpty_cfi::i1, edm::isFinite(), LogDebug, or, and DiDispStaMuonMonitor_cfi::pt.

◆ goodVar()

double PuppiContainer::goodVar ( PuppiCandidate const &  iPart,
std::vector< PuppiCandidate > const &  iParts,
int  iOpt,
const double  iRCone 
)
protected

Definition at line 80 of file PuppiContainer.cc.

83  {
84  return var_within_R(iOpt, iParts, iPart, iRCone);
85 }

◆ initialize()

void PuppiContainer::initialize ( const std::vector< RecoObj > &  iRecoObjects)

Definition at line 30 of file PuppiContainer.cc.

30  {
31  //Clear everything
32  fPFParticles.resize(0);
33  fPFParticlesForVar.resize(0);
35  fWeights.resize(0);
36  fVals.resize(0);
37  fRawAlphas.resize(0);
38  fAlphaMed.resize(0);
39  fAlphaRMS.resize(0);
40  fNPV = 1.;
41  //Link to the RecoObjects
42  fRecoParticles = &iRecoObjects;
43  fPFParticles.reserve(iRecoObjects.size());
44  fPFParticlesForVar.reserve(iRecoObjects.size());
45  fPFParticlesForVarChargedPV.reserve(iRecoObjects.size());
46  for (auto const &rParticle : *fRecoParticles) {
47  PuppiCandidate pCand;
48  pCand.id = rParticle.id;
49  if (edm::isFinite(rParticle.rapidity)) {
50  pCand.pt = rParticle.pt;
51  pCand.eta = rParticle.eta;
52  pCand.rapidity = rParticle.rapidity;
53  pCand.phi = rParticle.phi;
54  pCand.m = rParticle.m;
55  } else {
56  pCand.pt = 0.;
57  pCand.eta = 99.;
58  pCand.rapidity = 99.;
59  pCand.phi = 0.;
60  pCand.m = 0.;
61  }
62 
63  fPFParticles.push_back(pCand);
64 
65  // skip candidates to be ignored in the computation
66  // of PUPPI's alphas (e.g. electrons and muons if puppiNoLep=True)
67  if (std::abs(rParticle.id) == 3)
68  continue;
69 
70  fPFParticlesForVar.push_back(pCand);
71  // charged candidates assigned to LV
72  if (std::abs(rParticle.id) == 1)
73  fPFParticlesForVarChargedPV.push_back(pCand);
74  // if(fNPV < rParticle.vtxId) fNPV = rParticle.vtxId;
75  }
76 }

References funct::abs(), PuppiCandidate::eta, PuppiCandidate::id, edm::isFinite(), PuppiCandidate::m, PuppiCandidate::phi, PuppiCandidate::pt, and PuppiCandidate::rapidity.

◆ pfParticles()

const std::vector<PuppiCandidate>& PuppiContainer::pfParticles ( ) const
inline

Definition at line 15 of file PuppiContainer.h.

15 { return fPFParticles; }

References fPFParticles.

◆ puppiAlphas()

const std::vector<double>& PuppiContainer::puppiAlphas ( )
inline

Definition at line 18 of file PuppiContainer.h.

18 { return fVals; }

References fVals.

◆ puppiAlphasMed()

const std::vector<double>& PuppiContainer::puppiAlphasMed ( )
inline

Definition at line 20 of file PuppiContainer.h.

20 { return fAlphaMed; }

References fAlphaMed.

◆ puppiAlphasRMS()

const std::vector<double>& PuppiContainer::puppiAlphasRMS ( )
inline

Definition at line 21 of file PuppiContainer.h.

21 { return fAlphaRMS; }

References fAlphaRMS.

◆ puppiNAlgos()

int PuppiContainer::puppiNAlgos ( )
inline

Definition at line 23 of file PuppiContainer.h.

23 { return fNAlgos; }

References fNAlgos.

◆ puppiRawAlphas()

const std::vector<double>& PuppiContainer::puppiRawAlphas ( )
inline

Definition at line 17 of file PuppiContainer.h.

17 { return fRawAlphas; }

References fRawAlphas.

◆ puppiWeights()

const std::vector< double > & PuppiContainer::puppiWeights ( )

Definition at line 240 of file PuppiContainer.cc.

240  {
241  int lNParticles = fRecoParticles->size();
242 
243  fWeights.clear();
244  fWeights.reserve(lNParticles);
245  fVals.clear();
246  fVals.reserve(lNParticles);
247  for (int i0 = 0; i0 < fNAlgos; i0++)
248  fPuppiAlgo[i0].reset();
249 
250  int lNMaxAlgo = 1;
251  for (int i0 = 0; i0 < fNAlgos; i0++)
252  lNMaxAlgo = std::max(fPuppiAlgo[i0].numAlgos(), lNMaxAlgo);
253  //Run through all compute mean and RMS
254  for (int i0 = 0; i0 < lNMaxAlgo; i0++) {
256  }
257  if (fPuppiDiagnostics)
259 
260  std::vector<double> pVals;
261  pVals.reserve(lNParticles);
262  for (int i0 = 0; i0 < lNParticles; i0++) {
263  //Refresh
264  pVals.clear();
265  double pWeight = 1;
266  //Get the Puppi Id and if ill defined move on
267  const auto &rParticle = (*fRecoParticles)[i0];
268  int pPupId = getPuppiId(rParticle.pt, rParticle.eta);
269  if (pPupId == -1) {
270  fWeights.push_back(0);
271  fAlphaMed.push_back(-10);
272  fAlphaRMS.push_back(-10);
273  continue;
274  }
275 
276  // fill the p-values
277  double pChi2 = 0;
278  if (fUseExp) {
279  //Compute an Experimental Puppi Weight with delta Z info (very simple example)
280  pChi2 = getChi2FromdZ(rParticle.dZ);
281  //Now make sure Neutrals are not set
282  if ((std::abs(rParticle.pdgId) == 22) || (std::abs(rParticle.pdgId) == 130))
283  pChi2 = 0;
284  }
285  //Fill and compute the PuppiWeight
286  int lNAlgos = fPuppiAlgo[pPupId].numAlgos();
287  for (int i1 = 0; i1 < lNAlgos; i1++)
288  pVals.push_back(fVals[lNParticles * i1 + i0]);
289 
290  pWeight = fPuppiAlgo[pPupId].compute(pVals, pChi2);
291  //Apply the CHS weights
292  if (rParticle.id == 1 && fApplyCHS)
293  pWeight = 1;
294  if (rParticle.id == 2 && fApplyCHS)
295  pWeight = 0;
296  //Apply weight of 1 for leptons if puppiNoLep
297  if (rParticle.id == 3)
298  pWeight = 1;
299  //Basic Weight Checks
300  if (!edm::isFinite(pWeight)) {
301  pWeight = 0.0;
302  LogDebug("PuppiWeightError") << "====> Weight is nan : " << pWeight << " : pt " << rParticle.pt
303  << " -- eta : " << rParticle.eta << " -- Value" << fVals[i0]
304  << " -- id : " << rParticle.id << " -- NAlgos: " << lNAlgos << std::endl;
305  }
306  //Basic Cuts
307  if (pWeight * fPFParticles[i0].pt < fPuppiAlgo[pPupId].neutralPt(fNPV) && rParticle.id == 0)
308  pWeight = 0; //threshold cut on the neutral Pt
309  // Protect high pT photons (important for gamma to hadronic recoil balance)
310  if (fPtMaxPhotons > 0 && rParticle.pdgId == 22 && std::abs(fPFParticles[i0].eta) < fEtaMaxPhotons &&
311  fPFParticles[i0].pt > fPtMaxPhotons)
312  pWeight = 1.;
313  // Protect high pT neutrals
314  else if ((fPtMaxNeutrals > 0) && (rParticle.id == 0))
315  pWeight = std::clamp(
317  if (pWeight < fPuppiWeightCut)
318  pWeight = 0; //==> Elminate the low Weight stuff
319  if (fInvert)
320  pWeight = 1. - pWeight;
321  //std::cout << "rParticle.pt = " << rParticle.pt << ", rParticle.charge = " << rParticle.charge << ", rParticle.id = " << rParticle.id << ", weight = " << pWeight << std::endl;
322 
323  fWeights.push_back(pWeight);
324  fAlphaMed.push_back(fPuppiAlgo[pPupId].median());
325  fAlphaRMS.push_back(fPuppiAlgo[pPupId].rms());
326  //Now get rid of the thrown out weights for the particle collection
327 
328  // leave these lines in, in case want to move eventually to having no 1-to-1 correspondence between puppi and pf cands
329  // if( std::abs(pWeight) < std::numeric_limits<double>::denorm_min() ) continue; // this line seems not to work like it's supposed to...
330  // if(std::abs(pWeight) <= 0. ) continue;
331  }
332  return fWeights;
333 }

References funct::abs(), testProducerWithPsetDescEmpty_cfi::i1, edm::isFinite(), LogDebug, SiStripPI::max, pfDeepBoostedJetPreprocessParams_cfi::median, DiDispStaMuonMonitor_cfi::pt, reset(), and SiStripPI::rms.

◆ setNPV()

void PuppiContainer::setNPV ( int  iNPV)
inline

Definition at line 13 of file PuppiContainer.h.

13 { fNPV = iNPV; }

References fNPV.

◆ var_within_R()

double PuppiContainer::var_within_R ( int  iId,
const std::vector< PuppiCandidate > &  particles,
const PuppiCandidate centre,
const double  R 
)
protected

Definition at line 87 of file PuppiContainer.cc.

90  {
91  if (iId == -1)
92  return 1.;
93 
94  double const r2 = R * R;
95  double var = 0.;
96 
97  for (auto const &cand : particles) {
98  if (std::abs(cand.rapidity - centre.rapidity) < R) {
99  auto const dr2y = reco::deltaR2(cand.rapidity, cand.phi, centre.rapidity, centre.phi);
100  if (dr2y < r2) {
101  auto const dr2 = reco::deltaR2(cand.eta, cand.phi, centre.eta, centre.phi);
102  if (dr2 < 0.0001)
103  continue;
104  auto const pt = cand.pt;
105  if (iId == 5)
106  var += (pt * pt / dr2);
107  else if (iId == 4)
108  var += pt;
109  else if (iId == 3)
110  var += (1. / dr2);
111  else if (iId == 2)
112  var += (1. / dr2);
113  else if (iId == 1)
114  var += pt;
115  else if (iId == 0)
116  var += (pt / dr2);
117  }
118  }
119  }
120 
121  if ((var != 0.) and ((iId == 0) or (iId == 3) or (iId == 5)))
122  var = log(var);
123  else if (iId == 1)
124  var += centre.pt;
125 
126  return var;
127 }

References funct::abs(), reco::deltaR2(), PuppiCandidate::eta, dqm-mbProfile::log, or, ecalTrigSettings_cff::particles, PuppiCandidate::phi, PuppiCandidate::pt, DiDispStaMuonMonitor_cfi::pt, dttmaxenums::R, diffTwoXMLs::r2, PuppiCandidate::rapidity, and trigObjTnPSource_cfi::var.

Member Data Documentation

◆ fAlphaMed

std::vector<double> PuppiContainer::fAlphaMed
protected

Definition at line 50 of file PuppiContainer.h.

Referenced by puppiAlphasMed().

◆ fAlphaRMS

std::vector<double> PuppiContainer::fAlphaRMS
protected

Definition at line 51 of file PuppiContainer.h.

Referenced by puppiAlphasRMS().

◆ fApplyCHS

bool PuppiContainer::fApplyCHS
protected

Definition at line 53 of file PuppiContainer.h.

◆ fEtaMaxPhotons

double PuppiContainer::fEtaMaxPhotons
protected

Definition at line 60 of file PuppiContainer.h.

◆ fInvert

bool PuppiContainer::fInvert
protected

Definition at line 54 of file PuppiContainer.h.

◆ fNAlgos

int PuppiContainer::fNAlgos
protected

Definition at line 63 of file PuppiContainer.h.

Referenced by puppiNAlgos().

◆ fNeutralMinPt

double PuppiContainer::fNeutralMinPt
protected

Definition at line 56 of file PuppiContainer.h.

◆ fNeutralSlope

double PuppiContainer::fNeutralSlope
protected

Definition at line 57 of file PuppiContainer.h.

◆ fNPV

int PuppiContainer::fNPV
protected

Definition at line 64 of file PuppiContainer.h.

Referenced by setNPV().

◆ fPFParticles

std::vector<PuppiCandidate> PuppiContainer::fPFParticles
protected

Definition at line 44 of file PuppiContainer.h.

Referenced by pfParticles().

◆ fPFParticlesForVar

std::vector<PuppiCandidate> PuppiContainer::fPFParticlesForVar
protected

Definition at line 45 of file PuppiContainer.h.

◆ fPFParticlesForVarChargedPV

std::vector<PuppiCandidate> PuppiContainer::fPFParticlesForVarChargedPV
protected

Definition at line 46 of file PuppiContainer.h.

◆ fPtMaxNeutrals

double PuppiContainer::fPtMaxNeutrals
protected

Definition at line 61 of file PuppiContainer.h.

◆ fPtMaxNeutralsStartSlope

double PuppiContainer::fPtMaxNeutralsStartSlope
protected

Definition at line 62 of file PuppiContainer.h.

◆ fPtMaxPhotons

double PuppiContainer::fPtMaxPhotons
protected

Definition at line 59 of file PuppiContainer.h.

◆ fPuppiAlgo

std::vector<PuppiAlgo> PuppiContainer::fPuppiAlgo
protected

Definition at line 65 of file PuppiContainer.h.

◆ fPuppiDiagnostics

bool PuppiContainer::fPuppiDiagnostics
protected

Definition at line 42 of file PuppiContainer.h.

◆ fPuppiWeightCut

double PuppiContainer::fPuppiWeightCut
protected

Definition at line 58 of file PuppiContainer.h.

◆ fRawAlphas

std::vector<double> PuppiContainer::fRawAlphas
protected

Definition at line 49 of file PuppiContainer.h.

Referenced by puppiRawAlphas().

◆ fRecoParticles

const std::vector<RecoObj>* PuppiContainer::fRecoParticles
protected

Definition at line 43 of file PuppiContainer.h.

◆ fUseExp

bool PuppiContainer::fUseExp
protected

Definition at line 55 of file PuppiContainer.h.

◆ fVals

std::vector<double> PuppiContainer::fVals
protected

Definition at line 48 of file PuppiContainer.h.

Referenced by puppiAlphas().

◆ fWeights

std::vector<double> PuppiContainer::fWeights
protected

Definition at line 47 of file PuppiContainer.h.

PuppiCandidate::pt
double pt
Definition: PuppiCandidate.h:5
PuppiContainer::fPuppiWeightCut
double fPuppiWeightCut
Definition: PuppiContainer.h:58
ALCARECOTkAlBeamHalo_cff.etaMin
etaMin
GeV.
Definition: ALCARECOTkAlBeamHalo_cff.py:32
PuppiContainer::getChi2FromdZ
double getChi2FromdZ(double iDZ)
Definition: PuppiContainer.cc:225
PuppiContainer::fPtMaxNeutrals
double fPtMaxNeutrals
Definition: PuppiContainer.h:61
PuppiCandidate::eta
double eta
Definition: PuppiCandidate.h:6
DiDispStaMuonMonitor_cfi.pt
pt
Definition: DiDispStaMuonMonitor_cfi.py:39
pfDeepBoostedJetPreprocessParams_cfi.median
median
Definition: pfDeepBoostedJetPreprocessParams_cfi.py:12
testProducerWithPsetDescEmpty_cfi.i1
i1
Definition: testProducerWithPsetDescEmpty_cfi.py:45
ptMin
constexpr float ptMin
Definition: PhotonIDValueMapProducer.cc:155
SiStripPI::rms
Definition: SiStripPayloadInspectorHelper.h:169
PuppiContainer::fRecoParticles
const std::vector< RecoObj > * fRecoParticles
Definition: PuppiContainer.h:43
ecalTrigSettings_cff.particles
particles
Definition: ecalTrigSettings_cff.py:11
PuppiCandidate
Definition: PuppiCandidate.h:4
PuppiContainer::fRawAlphas
std::vector< double > fRawAlphas
Definition: PuppiContainer.h:49
PuppiContainer::fPFParticlesForVarChargedPV
std::vector< PuppiCandidate > fPFParticlesForVarChargedPV
Definition: PuppiContainer.h:46
trigObjTnPSource_cfi.var
var
Definition: trigObjTnPSource_cfi.py:21
PVValHelper::eta
Definition: PVValidationHelpers.h:69
PuppiContainer::fPtMaxPhotons
double fPtMaxPhotons
Definition: PuppiContainer.h:59
PuppiContainer::getPuppiId
int getPuppiId(float iPt, float iEta)
Definition: PuppiContainer.cc:208
PuppiContainer::goodVar
double goodVar(PuppiCandidate const &iPart, std::vector< PuppiCandidate > const &iParts, int iOpt, const double iRCone)
Definition: PuppiContainer.cc:80
PuppiAlgo
Definition: PuppiAlgo.h:10
PuppiContainer::fVals
std::vector< double > fVals
Definition: PuppiContainer.h:48
PuppiCandidate::id
int id
Definition: PuppiCandidate.h:10
PuppiContainer::fInvert
bool fInvert
Definition: PuppiContainer.h:54
PuppiContainer::fPuppiAlgo
std::vector< PuppiAlgo > fPuppiAlgo
Definition: PuppiContainer.h:65
LogDebug
#define LogDebug(id)
Definition: MessageLogger.h:223
PuppiContainer::fPuppiDiagnostics
bool fPuppiDiagnostics
Definition: PuppiContainer.h:42
PuppiContainer::getRawAlphas
void getRawAlphas(int iOpt, std::vector< PuppiCandidate > const &iConstits, std::vector< PuppiCandidate > const &iParticles, std::vector< PuppiCandidate > const &iChargeParticles)
Definition: PuppiContainer.cc:186
PuppiContainer::var_within_R
double var_within_R(int iId, const std::vector< PuppiCandidate > &particles, const PuppiCandidate &centre, const double R)
Definition: PuppiContainer.cc:87
SiStripPI::max
Definition: SiStripPayloadInspectorHelper.h:169
PuppiContainer::fApplyCHS
bool fApplyCHS
Definition: PuppiContainer.h:53
reco::deltaR2
constexpr auto deltaR2(const T1 &t1, const T2 &t2) -> decltype(t1.eta())
Definition: deltaR.h:16
diffTwoXMLs.r2
r2
Definition: diffTwoXMLs.py:73
cand
Definition: decayParser.h:32
PuppiContainer::fAlphaRMS
std::vector< double > fAlphaRMS
Definition: PuppiContainer.h:51
PuppiCandidate::rapidity
double rapidity
Definition: PuppiCandidate.h:9
PuppiContainer::fAlphaMed
std::vector< double > fAlphaMed
Definition: PuppiContainer.h:50
PuppiContainer::fPFParticlesForVar
std::vector< PuppiCandidate > fPFParticlesForVar
Definition: PuppiContainer.h:45
PuppiContainer::fEtaMaxPhotons
double fEtaMaxPhotons
Definition: PuppiContainer.h:60
PuppiContainer::fPFParticles
std::vector< PuppiCandidate > fPFParticles
Definition: PuppiContainer.h:44
PuppiContainer::fPtMaxNeutralsStartSlope
double fPtMaxNeutralsStartSlope
Definition: PuppiContainer.h:62
ALCARECOTkAlBeamHalo_cff.etaMax
etaMax
Definition: ALCARECOTkAlBeamHalo_cff.py:33
or
The Signals That Services Can Subscribe To This is based on ActivityRegistry and is current per Services can connect to the signals distributed by the ActivityRegistry in order to monitor the activity of the application Each possible callback has some defined which we here list in angle e< void, edm::EventID const &, edm::Timestamp const & > We also list in braces which AR_WATCH_USING_METHOD_ is used for those or
Definition: Activities.doc:12
edm::ParameterSet::getParameter
T getParameter(std::string const &) const
Definition: ParameterSet.h:303
PuppiCandidate::phi
double phi
Definition: PuppiCandidate.h:7
dqm-mbProfile.log
log
Definition: dqm-mbProfile.py:17
PuppiContainer::fWeights
std::vector< double > fWeights
Definition: PuppiContainer.h:47
PuppiContainer::getRMSAvg
void getRMSAvg(int iOpt, std::vector< PuppiCandidate > const &iConstits, std::vector< PuppiCandidate > const &iParticles, std::vector< PuppiCandidate > const &iChargeParticles)
Definition: PuppiContainer.cc:130
reset
void reset(double vett[256])
Definition: TPedValues.cc:11
L1TowerCalibrationProducer_cfi.iEta
iEta
Definition: L1TowerCalibrationProducer_cfi.py:60
PuppiContainer::fUseExp
bool fUseExp
Definition: PuppiContainer.h:55
funct::abs
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
edm::isFinite
constexpr bool isFinite(T x)
PuppiCandidate::m
double m
Definition: PuppiCandidate.h:8
dttmaxenums::R
Definition: DTTMax.h:29
PuppiContainer::fNPV
int fNPV
Definition: PuppiContainer.h:64
PuppiContainer::fNAlgos
int fNAlgos
Definition: PuppiContainer.h:63
MillePedeFileConverter_cfg.e
e
Definition: MillePedeFileConverter_cfg.py:37