CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Groups Pages
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)
 
std::vector< PuppiCandidate >
const & 
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 ()
 
std::vector< double > const & puppiWeights ()
 
void setPUProxy (double const iPUProxy)
 
 ~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
 
std::vector< PuppiCandidatefPFParticles
 
std::vector< PuppiCandidatefPFParticlesForVar
 
std::vector< PuppiCandidatefPFParticlesForVarChargedPV
 
double fPtMaxNeutrals
 
double fPtMaxNeutralsStartSlope
 
double fPtMaxPhotons
 
std::vector< PuppiAlgofPuppiAlgo
 
bool fPuppiDiagnostics
 
double fPuppiWeightCut
 
double fPUProxy
 
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 ( const edm::ParameterSet iConfig)

Definition at line 12 of file PuppiContainer.cc.

References edm::ParameterSet::getParameter().

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 }
double fPtMaxNeutralsStartSlope
double fEtaMaxPhotons
double fPtMaxNeutrals
double fPuppiWeightCut
T getParameter(std::string const &) const
Definition: ParameterSet.h:303
double fPtMaxPhotons
std::vector< PuppiAlgo > fPuppiAlgo
PuppiContainer::~PuppiContainer ( )

Definition at line 77 of file PuppiContainer.cc.

77 {}

Member Function Documentation

double PuppiContainer::getChi2FromdZ ( double  iDZ)
protected

Definition at line 232 of file PuppiContainer.cc.

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

232  {
233  //We need to obtain prob of PU + (1-Prob of LV)
234  // Prob(LV) = Gaus(dZ,sigma) where sigma = 1.5mm (its really more like 1mm)
235  //double lProbLV = ROOT::Math::normal_cdf_c(std::abs(iDZ),0.2)*2.; //*2 is to do it double sided
236  //Take iDZ to be corrected by sigma already
237  double lProbLV = ROOT::Math::normal_cdf_c(std::abs(iDZ), 1.) * 2.; //*2 is to do it double sided
238  double lProbPU = 1 - lProbLV;
239  if (lProbPU <= 0)
240  lProbPU = 1e-16; //Quick Trick to through out infs
241  if (lProbPU >= 0)
242  lProbPU = 1 - 1e-16; //Ditto
243  double lChi2PU = TMath::ChisquareQuantile(lProbPU, 1);
244  lChi2PU *= lChi2PU;
245  return lChi2PU;
246 }
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
int PuppiContainer::getPuppiId ( float  iPt,
float  iEta 
)
protected

Definition at line 215 of file PuppiContainer.cc.

References funct::abs(), Puppi_cff::etaMax, Puppi_cff::etaMin, and ptMin.

215  {
216  int lId = -1;
217  for (int i0 = 0; i0 < fNAlgos; i0++) {
218  int nEtaBinsPerAlgo = fPuppiAlgo[i0].etaBins();
219  for (int i1 = 0; i1 < nEtaBinsPerAlgo; i1++) {
220  if ((std::abs(iEta) >= fPuppiAlgo[i0].etaMin(i1)) && (std::abs(iEta) < fPuppiAlgo[i0].etaMax(i1))) {
221  fPuppiAlgo[i0].fixAlgoEtaBin(i1);
222  if (iPt > fPuppiAlgo[i0].ptMin()) {
223  lId = i0;
224  break;
225  }
226  }
227  }
228  }
229  //if(lId == -1) std::cerr << "Error : Full fiducial range is not defined " << std::endl;
230  return lId;
231 }
constexpr float ptMin
tuple etaMin
Definition: Puppi_cff.py:46
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
tuple etaMax
Definition: Puppi_cff.py:47
std::vector< PuppiAlgo > fPuppiAlgo
void PuppiContainer::getRawAlphas ( int  iOpt,
std::vector< PuppiCandidate > const &  iConstits,
std::vector< PuppiCandidate > const &  iParticles,
std::vector< PuppiCandidate > const &  iChargeParticles 
)
protected

Definition at line 193 of file PuppiContainer.cc.

References edm::isFinite(), and LogDebug.

196  {
197  for (int j0 = 0; j0 < fNAlgos; j0++) {
198  for (unsigned int i0 = 0; i0 < iConstits.size(); i0++) {
199  //Get the Puppi Sub Algo (given iteration)
200  int pAlgo = fPuppiAlgo[j0].algoId(iOpt);
201  bool pCharged = fPuppiAlgo[j0].isCharged(iOpt);
202  double pCone = fPuppiAlgo[j0].coneSize(iOpt);
203  //Compute the Puppi Metric
204  double const pVal = goodVar(iConstits[i0], pCharged ? iChargedParticles : iParticles, pAlgo, pCone);
205  fRawAlphas.push_back(pVal);
206  if (!edm::isFinite(pVal)) {
207  LogDebug("NotFound") << "====> Value is Nan " << pVal << " == " << iConstits[i0].pt << " -- "
208  << iConstits[i0].eta << endl;
209  continue;
210  }
211  }
212  }
213 }
constexpr bool isFinite(T x)
std::vector< double > fRawAlphas
double goodVar(PuppiCandidate const &iPart, std::vector< PuppiCandidate > const &iParts, int iOpt, const double iRCone)
std::vector< PuppiAlgo > fPuppiAlgo
#define LogDebug(id)
void PuppiContainer::getRMSAvg ( int  iOpt,
std::vector< PuppiCandidate > const &  iConstits,
std::vector< PuppiCandidate > const &  iParticles,
std::vector< PuppiCandidate > const &  iChargeParticles 
)
protected

Definition at line 137 of file PuppiContainer.cc.

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

140  {
141  for (unsigned int i0 = 0; i0 < iConstits.size(); i0++) {
142  //Calculate the Puppi Algo to use
143  int pPupId = getPuppiId(iConstits[i0].pt, iConstits[i0].eta);
144  if (pPupId == -1 || fPuppiAlgo[pPupId].numAlgos() <= iOpt) {
145  fVals.push_back(-1);
146  continue;
147  }
148  //Get the Puppi Sub Algo (given iteration)
149  int pAlgo = fPuppiAlgo[pPupId].algoId(iOpt);
150  bool pCharged = fPuppiAlgo[pPupId].isCharged(iOpt);
151  double pCone = fPuppiAlgo[pPupId].coneSize(iOpt);
152  // compute the Puppi metric:
153  // - calculate goodVar only for candidates that (1) will not be assigned a predefined weight (e.g 0, 1),
154  // or (2) are required for computations inside puppi-algos (see call to PuppiAlgo::add below)
155  double pVal = -1;
156  bool const getsDefaultWgtIfApplyCHS = iConstits[i0].id == 1 or iConstits[i0].id == 2;
157  if (not((fApplyCHS and getsDefaultWgtIfApplyCHS) or iConstits[i0].id == 3) or
158  (std::abs(iConstits[i0].eta) < fPuppiAlgo[pPupId].etaMaxExtrap() and getsDefaultWgtIfApplyCHS)) {
159  pVal = goodVar(iConstits[i0], pCharged ? iChargedParticles : iParticles, pAlgo, pCone);
160  }
161  fVals.push_back(pVal);
162 
163  if (!edm::isFinite(pVal)) {
164  LogDebug("NotFound") << "====> Value is Nan " << pVal << " == " << iConstits[i0].pt << " -- " << iConstits[i0].eta
165  << endl;
166  continue;
167  }
168 
169  // code added by Nhan: now instead for every algorithm give it all the particles
170  for (int i1 = 0; i1 < fNAlgos; i1++) {
171  // skip cands outside of algo's etaMaxExtrap, as they would anyway be ignored inside PuppiAlgo::add (see end of the block)
172  if (not(std::abs(iConstits[i0].eta) < fPuppiAlgo[i1].etaMaxExtrap() and getsDefaultWgtIfApplyCHS))
173  continue;
174 
175  auto curVal = pVal;
176  // recompute goodVar if algo has changed
177  if (i1 != pPupId) {
178  pAlgo = fPuppiAlgo[i1].algoId(iOpt);
179  pCharged = fPuppiAlgo[i1].isCharged(iOpt);
180  pCone = fPuppiAlgo[i1].coneSize(iOpt);
181  curVal = goodVar(iConstits[i0], pCharged ? iChargedParticles : iParticles, pAlgo, pCone);
182  }
183 
184  fPuppiAlgo[i1].add(iConstits[i0], curVal, iOpt);
185  }
186  }
187 
188  for (int i0 = 0; i0 < fNAlgos; i0++)
189  fPuppiAlgo[i0].computeMedRMS(iOpt);
190 }
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::EventIDconst &, edm::Timestampconst & > We also list in braces which AR_WATCH_USING_METHOD_ is used for those or
Definition: Activities.doc:12
std::vector< double > fVals
constexpr bool isFinite(T x)
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
int getPuppiId(float iPt, float iEta)
double goodVar(PuppiCandidate const &iPart, std::vector< PuppiCandidate > const &iParts, int iOpt, const double iRCone)
std::vector< PuppiAlgo > fPuppiAlgo
#define LogDebug(id)
double PuppiContainer::goodVar ( PuppiCandidate const &  iPart,
std::vector< PuppiCandidate > const &  iParts,
int  iOpt,
const double  iRCone 
)
protected

Definition at line 79 of file PuppiContainer.cc.

82  {
83  return var_within_R(iOpt, iParts, iPart, iRCone);
84 }
double var_within_R(int iId, const std::vector< PuppiCandidate > &particles, const PuppiCandidate &centre, const double R)
void PuppiContainer::initialize ( const std::vector< RecoObj > &  iRecoObjects)

Definition at line 30 of file PuppiContainer.cc.

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

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  fPUProxy = 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  }
75 }
std::vector< PuppiCandidate > fPFParticlesForVar
std::vector< PuppiCandidate > fPFParticles
std::vector< double > fVals
constexpr bool isFinite(T x)
std::vector< double > fRawAlphas
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
std::vector< double > fAlphaMed
std::vector< PuppiCandidate > fPFParticlesForVarChargedPV
const std::vector< RecoObj > * fRecoParticles
std::vector< double > fWeights
std::vector< double > fAlphaRMS
std::vector<PuppiCandidate> const& PuppiContainer::pfParticles ( ) const
inline

Definition at line 15 of file PuppiContainer.h.

References fPFParticles.

15 { return fPFParticles; }
std::vector< PuppiCandidate > fPFParticles
const std::vector<double>& PuppiContainer::puppiAlphas ( )
inline

Definition at line 18 of file PuppiContainer.h.

References fVals.

18 { return fVals; }
std::vector< double > fVals
const std::vector<double>& PuppiContainer::puppiAlphasMed ( )
inline

Definition at line 20 of file PuppiContainer.h.

References fAlphaMed.

20 { return fAlphaMed; }
std::vector< double > fAlphaMed
const std::vector<double>& PuppiContainer::puppiAlphasRMS ( )
inline

Definition at line 21 of file PuppiContainer.h.

References fAlphaRMS.

21 { return fAlphaRMS; }
std::vector< double > fAlphaRMS
int PuppiContainer::puppiNAlgos ( )
inline

Definition at line 23 of file PuppiContainer.h.

References fNAlgos.

23 { return fNAlgos; }
const std::vector<double>& PuppiContainer::puppiRawAlphas ( )
inline

Definition at line 17 of file PuppiContainer.h.

References fRawAlphas.

17 { return fRawAlphas; }
std::vector< double > fRawAlphas
std::vector< double > const & PuppiContainer::puppiWeights ( )

Definition at line 247 of file PuppiContainer.cc.

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

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

Definition at line 13 of file PuppiContainer.h.

References fPUProxy.

13 { fPUProxy = iPUProxy; }
double PuppiContainer::var_within_R ( int  iId,
const std::vector< PuppiCandidate > &  particles,
const PuppiCandidate centre,
const double  R 
)
protected

Definition at line 86 of file PuppiContainer.cc.

References funct::abs(), reco::deltaR2(), PuppiCandidate::eta, log, or, PuppiCandidate::phi, PuppiCandidate::pt, DiDispStaMuonMonitor_cfi::pt, dttmaxenums::R, diffTwoXMLs::r2, PuppiCandidate::rapidity, and ALCARECOEcalPhiSym_cff::var.

89  {
90  if (iId == -1)
91  return 1.;
92 
93  double const r2 = R * R;
94  double var = 0.;
95 
96  for (auto const &cand : particles) {
97  if (std::abs(cand.rapidity - centre.rapidity) < R) {
98  auto const dr2y = reco::deltaR2(cand.rapidity, cand.phi, centre.rapidity, centre.phi);
99  if (dr2y < r2) {
100  auto const dr2 = reco::deltaR2(cand.eta, cand.phi, centre.eta, centre.phi);
101  if (dr2 < 0.0001)
102  continue;
103  auto const pt = cand.pt;
104  switch (iId) {
105  case 5:
106  var += (pt * pt / dr2);
107  break;
108  case 4:
109  var += pt;
110  break;
111  case 3:
112  var += (1. / dr2);
113  break;
114  case 2:
115  var += (1. / dr2);
116  break;
117  case 1:
118  var += pt;
119  break;
120  case 0:
121  var += (pt / dr2);
122  break;
123  }
124  }
125  }
126  }
127 
128  if ((var != 0.) and ((iId == 0) or (iId == 3) or (iId == 5)))
129  var = log(var);
130  else if (iId == 1)
131  var += centre.pt;
132 
133  return var;
134 }
static std::vector< std::string > checklist log
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::EventIDconst &, edm::Timestampconst & > We also list in braces which AR_WATCH_USING_METHOD_ is used for those or
Definition: Activities.doc:12
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
constexpr auto deltaR2(const T1 &t1, const T2 &t2) -> decltype(t1.eta())
Definition: deltaR.h:16

Member Data Documentation

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

Definition at line 50 of file PuppiContainer.h.

Referenced by puppiAlphasMed().

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

Definition at line 51 of file PuppiContainer.h.

Referenced by puppiAlphasRMS().

bool PuppiContainer::fApplyCHS
protected

Definition at line 53 of file PuppiContainer.h.

double PuppiContainer::fEtaMaxPhotons
protected

Definition at line 60 of file PuppiContainer.h.

bool PuppiContainer::fInvert
protected

Definition at line 54 of file PuppiContainer.h.

int PuppiContainer::fNAlgos
protected

Definition at line 63 of file PuppiContainer.h.

Referenced by puppiNAlgos().

double PuppiContainer::fNeutralMinPt
protected

Definition at line 56 of file PuppiContainer.h.

double PuppiContainer::fNeutralSlope
protected

Definition at line 57 of file PuppiContainer.h.

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

Definition at line 44 of file PuppiContainer.h.

Referenced by pfParticles().

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

Definition at line 45 of file PuppiContainer.h.

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

Definition at line 46 of file PuppiContainer.h.

double PuppiContainer::fPtMaxNeutrals
protected

Definition at line 61 of file PuppiContainer.h.

double PuppiContainer::fPtMaxNeutralsStartSlope
protected

Definition at line 62 of file PuppiContainer.h.

double PuppiContainer::fPtMaxPhotons
protected

Definition at line 59 of file PuppiContainer.h.

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

Definition at line 65 of file PuppiContainer.h.

bool PuppiContainer::fPuppiDiagnostics
protected

Definition at line 42 of file PuppiContainer.h.

double PuppiContainer::fPuppiWeightCut
protected

Definition at line 58 of file PuppiContainer.h.

double PuppiContainer::fPUProxy
protected

Definition at line 64 of file PuppiContainer.h.

Referenced by setPUProxy().

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

Definition at line 49 of file PuppiContainer.h.

Referenced by puppiRawAlphas().

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

Definition at line 43 of file PuppiContainer.h.

bool PuppiContainer::fUseExp
protected

Definition at line 55 of file PuppiContainer.h.

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

Definition at line 48 of file PuppiContainer.h.

Referenced by puppiAlphas().

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

Definition at line 47 of file PuppiContainer.h.