CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
List of all members | Classes | Public Types | Public Member Functions | Static Public Member Functions | Private Attributes
reco::Jet Class Reference

Base class for all types of Jets. More...

#include <Jet.h>

Inheritance diagram for reco::Jet:
reco::CompositePtrCandidate reco::LeafCandidate reco::Candidate pat::PATObject< reco::Jet > reco::BasicJet reco::CaloJet reco::GenJet reco::JPTJet reco::PFClusterJet reco::PFJet reco::TrackJet

Classes

class  EtaPhiMoments
 record to store eta-phi first and second moments More...
 

Public Types

typedef edm::Ptr< CandidateConstituent
 
typedef std::vector< ConstituentConstituents
 
- Public Types inherited from reco::CompositePtrCandidate
typedef std::vector< CandidatePtrdaughters
 collection of references to daughters More...
 
typedef std::vector< CandidatePtrmothers
 collection of references to daughters More...
 
- Public Types inherited from reco::LeafCandidate
typedef int Charge
 electric charge type More...
 
typedef CandidateCollection daughters
 collection of daughter candidates More...
 
typedef unsigned int index
 
typedef math::XYZTLorentzVector LorentzVector
 Lorentz vector. More...
 
typedef math::XYZPoint Point
 point in the space More...
 
typedef
math::PtEtaPhiMLorentzVector 
PolarLorentzVector
 Lorentz vector. More...
 
typedef math::XYZVector Vector
 point in the space More...
 
- Public Types inherited from reco::Candidate
enum  { dimension = 3 }
 
enum  { size = dimension * (dimension + 1)/2 }
 matix size More...
 
typedef int Charge
 electric charge type More...
 
typedef candidate::const_iterator const_iterator
 
typedef math::Error< dimension >
::type 
CovarianceMatrix
 covariance error matrix (3x3) More...
 
typedef unsigned int index
 index type More...
 
typedef candidate::iterator iterator
 
typedef math::XYZTLorentzVector LorentzVector
 Lorentz vector. More...
 
typedef math::XYZPoint Point
 point in the space More...
 
typedef
math::PtEtaPhiMLorentzVector 
PolarLorentzVector
 Lorentz vector. More...
 
typedef size_t size_type
 
typedef math::XYZVector Vector
 point in the space More...
 

Public Member Functions

float constituentEtaPhiSpread () const
 
float constituentPtDistribution () const
 
float etaetaMoment () const
 eta-eta second moment, ET weighted More...
 
float etaphiMoment () const
 eta-phi second moment, ET weighted More...
 
EtaPhiMoments etaPhiStatistics () const
 eta-phi statistics, ET weighted More...
 
float etInAnnulus (float fRmin, float fRmax) const
 ET in annulus between rmin and rmax around jet direction. More...
 
virtual Constituents getJetConstituents () const
 list of constituents More...
 
virtual std::vector< const
reco::Candidate * > 
getJetConstituentsQuick () const
 quick list of constituents More...
 
bool isJet () const
 
 Jet ()
 Default constructor. More...
 
 Jet (const LorentzVector &fP4, const Point &fVertex)
 Initiator. More...
 
 Jet (const LorentzVector &fP4, const Point &fVertex, const Constituents &fConstituents)
 
virtual float jetArea () const
 get jet area More...
 
float maxDistance () const
 maximum distance from jet to constituent More...
 
int nCarrying (float fFraction) const
 return # of constituent carrying fraction of energy More...
 
virtual int nConstituents () const
 

of constituents

More...
 
virtual int nPasses () const
 number of passes taken by algorithm More...
 
float phiphiMoment () const
 phi-phi second moment, ET weighted More...
 
virtual float pileup () const
 pileup energy contribution as calculated by algorithm More...
 
virtual std::string print () const
 Print object. More...
 
virtual void scaleEnergy (double fScale)
 scale energy of the jet More...
 
virtual void setJetArea (float fArea)
 set jet area More...
 
virtual void setNPasses (int fPasses)
 Set number of passes taken by algorithm. More...
 
virtual void setPileup (float fEnergy)
 Set pileup energy contribution as calculated by algorithm. More...
 
virtual ~Jet ()
 Destructor. More...
 
- Public Member Functions inherited from reco::CompositePtrCandidate
void addDaughter (const CandidatePtr &)
 add a daughter via a reference More...
 
virtual const_iterator begin () const
 first daughter const_iterator More...
 
virtual iterator begin ()
 first daughter iterator More...
 
void clearDaughters ()
 clear daughter references More...
 
virtual CompositePtrCandidateclone () const
 returns a clone of the candidate More...
 
 CompositePtrCandidate ()
 default constructor More...
 
 CompositePtrCandidate (Charge q, const LorentzVector &p4, const Point &vtx=Point(0, 0, 0), int pdgId=0, int status=0, bool integerCharge=true)
 constructor from values More...
 
 CompositePtrCandidate (Charge q, const PolarLorentzVector &p4, const Point &vtx=Point(0, 0, 0), int pdgId=0, int status=0, bool integerCharge=true)
 constructor from values More...
 
 CompositePtrCandidate (const Candidate &p)
 constructor from a Candidate More...
 
virtual const Candidatedaughter (size_type) const
 return daughter at a given position, i = 0, ... numberOfDaughters() - 1 (read only mode) More...
 
virtual Candidatedaughter (size_type)
 return daughter at a given position, i = 0, ... numberOfDaughters() - 1 More...
 
CandidatePtr daughterPtr (size_type i) const
 reference to daughter at given position More...
 
const daughtersdaughterPtrVector () const
 references to daughtes More...
 
virtual const_iterator end () const
 last daughter const_iterator More...
 
virtual iterator end ()
 last daughter iterator More...
 
virtual const Candidatemother (size_t i=0) const
 return pointer to mother More...
 
virtual size_t numberOfDaughters () const
 number of daughters More...
 
virtual size_t numberOfMothers () const
 number of mothers More...
 
virtual size_type numberOfSourceCandidatePtrs () const
 
virtual CandidatePtr sourceCandidatePtr (size_type i) const
 
virtual ~CompositePtrCandidate ()
 destructor More...
 
- Public Member Functions inherited from reco::LeafCandidate
template<typename S >
daughter_iterator< S >::type beginFilter (const S &s) const
 
virtual Vector boostToCM () const GCC11_FINAL
 
virtual int charge () const GCC11_FINAL
 electric charge More...
 
virtual Candidatedaughter (const std::string &s)
 return daughter with a specified role name More...
 
virtual const Candidatedaughter (const std::string &s) const
 return daughter with a specified role name More...
 
template<typename S >
daughter_iterator< S >::type endFilter (const S &s) const
 
virtual double energy () const GCC11_FINAL
 energy More...
 
virtual double et () const GCC11_FINAL
 transverse energy More...
 
virtual float eta () const GCC11_FINAL
 momentum pseudorapidity More...
 
virtual void fillVertexCovariance (CovarianceMatrix &v) const
 fill SMatrix More...
 
template<typename T >
T get () const
 get a component More...
 
template<typename T , typename Tag >
T get () const
 get a component More...
 
template<typename T >
T get (size_type i) const
 get a component More...
 
template<typename T , typename Tag >
T get (size_type i) const
 get a component More...
 
virtual bool hasMasterClone () const
 
virtual bool hasMasterClonePtr () const
 
virtual bool isCaloMuon () const
 
virtual bool isConvertedPhoton () const
 
virtual bool isElectron () const
 
virtual bool isGlobalMuon () const
 
virtual bool isMuon () const
 
virtual bool isPhoton () const
 
virtual bool isStandAloneMuon () const
 
virtual bool isTrackerMuon () const
 
 LeafCandidate ()
 default constructor More...
 
 LeafCandidate (const Candidate &c)
 
template<typename P4 >
 LeafCandidate (Charge q, const P4 &p4, const Point &vtx=Point(0, 0, 0), int pdgId=0, int status=0, bool integerCharge=true)
 constructor from Any values More...
 
 LeafCandidate (Charge q, const PtEtaPhiMass &p4, const Point &vtx=Point(0, 0, 0), int pdgId=0, int status=0, bool integerCharge=true)
 constructor from values More...
 
 LeafCandidate (Charge q, const LorentzVector &p4, const Point &vtx=Point(0, 0, 0), int pdgId=0, int status=0, bool integerCharge=true)
 constructor from values More...
 
 LeafCandidate (Charge q, const PolarLorentzVector &p4, const Point &vtx=Point(0, 0, 0), int pdgId=0, int status=0, bool integerCharge=true)
 constructor from values More...
 
virtual bool longLived () const GCC11_FINAL
 is long lived? More...
 
virtual float mass () const GCC11_FINAL
 mass More...
 
virtual bool massConstraint () const GCC11_FINAL
 do mass constraint? More...
 
virtual float massSqr () const GCC11_FINAL
 mass squared More...
 
virtual const CandidateBaseRefmasterClone () const
 
virtual const CandidatePtrmasterClonePtr () const
 
template<typename Ref >
Ref masterRef () const
 cast master clone reference to a concrete type More...
 
virtual Vector momentum () const GCC11_FINAL
 spatial momentum vector More...
 
virtual double mt () const GCC11_FINAL
 transverse mass More...
 
virtual double mtSqr () const GCC11_FINAL
 transverse mass squared More...
 
template<typename T >
size_type numberOf () const
 number of components More...
 
template<typename T , typename Tag >
size_type numberOf () const
 number of components More...
 
virtual double p () const GCC11_FINAL
 magnitude of momentum vector More...
 
virtual const LorentzVectorp4 () const GCC11_FINAL
 four-momentum Lorentz vector More...
 
virtual int pdgId () const GCC11_FINAL
 PDG identifier. More...
 
virtual float phi () const GCC11_FINAL
 momentum azimuthal angle More...
 
virtual const PolarLorentzVectorpolarP4 () const GCC11_FINAL
 four-momentum Lorentz vector More...
 
virtual float pt () const GCC11_FINAL
 transverse momentum More...
 
virtual double px () const GCC11_FINAL
 x coordinate of momentum vector More...
 
virtual double py () const GCC11_FINAL
 y coordinate of momentum vector More...
 
virtual double pz () const GCC11_FINAL
 z coordinate of momentum vector More...
 
virtual double rapidity () const GCC11_FINAL
 rapidity More...
 
virtual void setCharge (Charge q) GCC11_FINAL
 set electric charge More...
 
virtual void setLongLived () GCC11_FINAL
 set long lived flag More...
 
virtual void setMass (double m) GCC11_FINAL
 set particle mass More...
 
virtual void setMassConstraint () GCC11_FINAL
 set mass constraint flag More...
 
virtual void setP4 (const LorentzVector &p4) GCC11_FINAL
 set 4-momentum More...
 
virtual void setP4 (const PolarLorentzVector &p4) GCC11_FINAL
 set 4-momentum More...
 
virtual void setPdgId (int pdgId) GCC11_FINAL
 
virtual void setPz (double pz) GCC11_FINAL
 
virtual void setStatus (int status) GCC11_FINAL
 set status word More...
 
virtual void setThreeCharge (Charge qx3) GCC11_FINAL
 set electric charge More...
 
virtual void setVertex (const Point &vertex)
 set vertex More...
 
virtual int status () const GCC11_FINAL
 status word More...
 
virtual double theta () const GCC11_FINAL
 momentum polar angle More...
 
virtual int threeCharge () const GCC11_FINAL
 electric charge More...
 
virtual const Pointvertex () const
 vertex position (overwritten by PF...) More...
 
virtual double vertexChi2 () const
 chi-squares More...
 
virtual double vertexCovariance (int i, int j) const
 (i, j)-th element of error matrix, i, j = 0, ... 2 More...
 
CovarianceMatrix vertexCovariance () const GCC11_FINAL
 return SMatrix More...
 
virtual double vertexNdof () const
 
virtual double vertexNormalizedChi2 () const
 chi-squared divided by n.d.o.f. More...
 
virtual double vx () const
 x coordinate of vertex position More...
 
virtual double vy () const
 y coordinate of vertex position More...
 
virtual double vz () const
 z coordinate of vertex position More...
 
virtual double y () const GCC11_FINAL
 rapidity More...
 
virtual ~LeafCandidate ()
 destructor More...
 
- Public Member Functions inherited from reco::Candidate
template<typename S >
daughter_iterator< S >::type beginFilter (const S &s) const
 
 Candidate ()
 default constructor More...
 
template<typename S >
daughter_iterator< S >::type endFilter (const S &s) const
 
template<typename T >
T get () const
 get a component More...
 
template<typename T , typename Tag >
T get () const
 get a component More...
 
template<typename T >
T get (size_type i) const
 get a component More...
 
template<typename T , typename Tag >
T get (size_type i) const
 get a component More...
 
template<typename Ref >
Ref masterRef () const
 cast master clone reference to a concrete type More...
 
template<typename T >
size_type numberOf () const
 number of components More...
 
template<typename T , typename Tag >
size_type numberOf () const
 number of components More...
 
virtual void setSourceCandidatePtr (const CandidatePtr &ptr)
 Set the ptr to the source Candidate. More...
 
virtual ~Candidate ()
 destructor More...
 

Static Public Member Functions

static float detectorEta (float fZVertex, float fPhysicsEta)
 static function to convert physics eta to detector eta More...
 
static Candidate::LorentzVector detectorP4 (const Candidate::Point &vertex, const Candidate &inParticle)
 
static float physicsEta (float fZVertex, float fDetectorEta)
 static function to convert detector eta to physics eta More...
 
static Candidate::LorentzVector physicsP4 (const Candidate::Point &newVertex, const Candidate &inParticle, const Candidate::Point &oldVertex=Candidate::Point(0, 0, 0))
 

Private Attributes

float mJetArea
 
int mPassNumber
 
float mPileupEnergy
 

Additional Inherited Members

- Static Public Attributes inherited from reco::LeafCandidate
static const unsigned int longLivedTag = 65536
 long lived flag More...
 
static const unsigned int massConstraintTag = 131072
 do mass constraint flag More...
 
- Protected Member Functions inherited from reco::LeafCandidate
void cacheCartesian () const
 set internal cache More...
 
void cachePolar () const
 set internal cache More...
 
void clearCache () const
 clear internal cache More...
 
- Protected Attributes inherited from reco::LeafCandidate
bool cacheCartesianFixed_
 
bool cachePolarFixed_
 has cache been set? More...
 
float eta_
 
float mass_
 
LorentzVector p4Cartesian_
 internal cache for p4 More...
 
PolarLorentzVector p4Polar_
 internal cache for p4 More...
 
int pdgId_
 PDG identifier. More...
 
float phi_
 
float pt_
 four-momentum Lorentz vector More...
 
Charge qx3_
 electric charge More...
 
int status_
 status word More...
 
Point vertex_
 vertex position More...
 

Detailed Description

Base class for all types of Jets.

GenericJet describes jets made from arbitrary constituents, No direct constituents reference is stored for now

Author
Fedor Ratnikov, UMd
Version
Mar 23, 2007 by F.R.
Id:
GenericJet.h,v 1.8 2007/08/24 17:35:23 fedor Exp

Jet describes properties common for all kinds of jets, essentially kinematics. Base class for all types of Jets.

Author
Fedor Ratnikov, UMd
Version
Original: April 22, 2005 by Fernando Varela Rodriguez.
May 23, 2006 by F.R.
Id:
Jet.h,v 1.32 2012/02/01 14:51:10 pandolf Exp

Definition at line 21 of file Jet.h.

Member Typedef Documentation

Definition at line 23 of file Jet.h.

typedef std::vector<Constituent> reco::Jet::Constituents

Definition at line 24 of file Jet.h.

Constructor & Destructor Documentation

reco::Jet::Jet ( )
inline

Default constructor.

Definition at line 37 of file Jet.h.

37 : mJetArea (0), mPileupEnergy (0), mPassNumber (0) {}
int mPassNumber
Definition: Jet.h:123
float mJetArea
Definition: Jet.h:121
float mPileupEnergy
Definition: Jet.h:122
Jet::Jet ( const LorentzVector fP4,
const Point fVertex 
)

Initiator.

Definition at line 175 of file Jet.cc.

177  : CompositePtrCandidate (0, fP4, fVertex),
178  mJetArea (0),
179  mPileupEnergy (0),
180  mPassNumber (0)
181 {}
int mPassNumber
Definition: Jet.h:123
CompositePtrCandidate()
default constructor
float mJetArea
Definition: Jet.h:121
float mPileupEnergy
Definition: Jet.h:122
Jet::Jet ( const LorentzVector fP4,
const Point fVertex,
const Constituents fConstituents 
)

Definition at line 162 of file Jet.cc.

References reco::CompositePtrCandidate::addDaughter(), and i.

165  : CompositePtrCandidate (0, fP4, fVertex),
166  mJetArea (0),
167  mPileupEnergy (0),
168  mPassNumber (0)
169 {
170  for (unsigned i = 0; i < fConstituents.size (); i++) {
171  addDaughter (fConstituents [i]);
172  }
173 }
int i
Definition: DBlmapReader.cc:9
int mPassNumber
Definition: Jet.h:123
CompositePtrCandidate()
default constructor
float mJetArea
Definition: Jet.h:121
void addDaughter(const CandidatePtr &)
add a daughter via a reference
float mPileupEnergy
Definition: Jet.h:122
virtual reco::Jet::~Jet ( )
inlinevirtual

Destructor.

Reimplemented in pat::Jet.

Definition at line 42 of file Jet.h.

42 {}

Member Function Documentation

float Jet::constituentEtaPhiSpread ( ) const

Definition at line 395 of file Jet.cc.

References reco::deltaR(), PFRecoTauDiscriminationAgainstElectronDeadECAL_cfi::dR, getJetConstituents(), and reco::LeafCandidate::pt().

Referenced by FFTGenericScaleCalculator::mapFFTJet().

395  {
396 
397  Jet::Constituents constituents = this->getJetConstituents();
398 
399 
400  float sum_pt2 = 0.;
401  float sum_pt2deltaR2 = 0.;
402 
403  for( unsigned iConst=0; iConst<constituents.size(); ++iConst ) {
404 
405  LorentzVector thisConstituent = constituents[iConst]->p4();
406 
407  float pt = thisConstituent.Pt();
408  float pt2 = pt*pt;
409  double dR = deltaR (*this, *(constituents[iConst]));
410  float pt2deltaR2 = pt*pt*dR*dR;
411 
412  sum_pt2 += pt2;
413  sum_pt2deltaR2 += pt2deltaR2;
414 
415  } //for constituents
416 
417  float rmsCand_value = (sum_pt2>0.) ? sum_pt2deltaR2/sum_pt2 : 0.;
418 
419  return rmsCand_value;
420 
421 } //constituentEtaPhiSpread
math::XYZTLorentzVector LorentzVector
std::vector< Constituent > Constituents
Definition: Jet.h:24
auto deltaR(const T1 &t1, const T2 &t2) -> decltype(t1.eta())
Definition: deltaR.h:30
virtual float pt() const GCC11_FINAL
transverse momentum
virtual Constituents getJetConstituents() const
list of constituents
Definition: Jet.cc:351
float Jet::constituentPtDistribution ( ) const

Definition at line 370 of file Jet.cc.

References getJetConstituents(), reco::LeafCandidate::pt(), and mathSSE::sqrt().

Referenced by FFTGenericScaleCalculator::mapFFTJet().

370  {
371 
372  Jet::Constituents constituents = this->getJetConstituents();
373 
374  float sum_pt2 = 0.;
375  float sum_pt = 0.;
376 
377  for( unsigned iConst=0; iConst<constituents.size(); ++iConst ) {
378 
379  float pt = constituents[iConst]->p4().Pt();
380  float pt2 = pt*pt;
381 
382  sum_pt += pt;
383  sum_pt2 += pt2;
384 
385  } //for constituents
386 
387  float ptD_value = (sum_pt>0.) ? sqrt( sum_pt2 / (sum_pt*sum_pt) ) : 0.;
388 
389  return ptD_value;
390 
391 } //constituentPtDistribution
std::vector< Constituent > Constituents
Definition: Jet.h:24
T sqrt(T t)
Definition: SSEVec.h:48
virtual float pt() const GCC11_FINAL
transverse momentum
virtual Constituents getJetConstituents() const
list of constituents
Definition: Jet.cc:351
float Jet::detectorEta ( float  fZVertex,
float  fPhysicsEta 
)
static

static function to convert physics eta to detector eta

Definition at line 326 of file Jet.cc.

Referenced by reco::GenJet::detectorEta().

326  {
327  CaloPointZ refPoint (fZVertex, fPhysicsEta);
328  return refPoint.etaReference (0.);
329 }
Candidate::LorentzVector Jet::detectorP4 ( const Candidate::Point vertex,
const Candidate inParticle 
)
static

Definition at line 340 of file Jet.cc.

References reco::Candidate::energy(), reco::Candidate::momentum(), np, reco::LeafCandidate::p(), and p3.

Referenced by reco::CaloJet::detectorP4().

340  {
341  CaloPoint3D<Point> caloPoint(vertex,inParticle.momentum()); // Jet position in Calo.
342  static const Point np(0,0,0);
343  Vector detectorDir = caloPoint.caloPoint() - np;
344  double p = inParticle.momentum().r();
345  Vector p3 = p * detectorDir.unit();
346  LorentzVector returnVector(p3.x(), p3.y(), p3.z(), inParticle.energy());
347  return returnVector;
348 }
virtual double energy() const =0
energy
virtual double p() const GCC11_FINAL
magnitude of momentum vector
virtual const Point & vertex() const
vertex position (overwritten by PF...)
math::XYZTLorentzVector LorentzVector
virtual Vector momentum() const =0
spatial momentum vector
int np
Definition: AMPTWrapper.h:33
fixed size vector
Definition: Vector.h:23
math::XYZPoint Point
double p3[4]
Definition: TauolaWrapper.h:91
float Jet::etaetaMoment ( ) const

eta-eta second moment, ET weighted

eta-eta second moment

Definition at line 224 of file Jet.cc.

References getJetConstituentsQuick(), i, relativeConstraints::value, and histoStyle::weight.

Referenced by fireworks::makeEveJetCone(), and JetIDSelectionFunctor::operator()().

224  {
225  std::vector<const Candidate*> towers = getJetConstituentsQuick ();
226  double sumw = 0;
227  double sum = 0;
228  double sum2 = 0;
229  int i = towers.size();
230  while (--i >= 0) {
231  double value = towers[i]->eta();
232  double weight = towers[i]->et();
233  sumw += weight;
234  sum += value * weight;
235  sum2 += value * value * weight;
236  }
237  return sumw > 0 ? (sum2 - sum*sum/sumw ) / sumw : 0;
238 }
int i
Definition: DBlmapReader.cc:9
virtual std::vector< const reco::Candidate * > getJetConstituentsQuick() const
quick list of constituents
Definition: Jet.cc:359
int weight
Definition: histoStyle.py:50
float Jet::etaphiMoment ( ) const

eta-phi second moment, ET weighted

eta-phi second moment

Definition at line 259 of file Jet.cc.

References reco::deltaPhi(), getJetConstituentsQuick(), i, reco::LeafCandidate::phi(), and histoStyle::weight.

259  {
260  double phiRef = phi();
261  std::vector<const Candidate*> towers = getJetConstituentsQuick ();
262  double sumw = 0;
263  double sumA = 0;
264  double sumB = 0;
265  double sumAB = 0;
266  int i = towers.size();
267  while (--i >= 0) {
268  double valueA = towers[i]->eta();
269  double valueB = deltaPhi (towers[i]->phi(), phiRef);
270  double weight = towers[i]->et();
271  sumw += weight;
272  sumA += valueA * weight;
273  sumB += valueB * weight;
274  sumAB += valueA * valueB * weight;
275  }
276  return sumw > 0 ? (sumAB - sumA*sumB/sumw ) / sumw : 0;
277 }
int i
Definition: DBlmapReader.cc:9
virtual float phi() const GCC11_FINAL
momentum azimuthal angle
double deltaPhi(double phi1, double phi2)
Definition: deltaPhi.h:12
virtual std::vector< const reco::Candidate * > getJetConstituentsQuick() const
quick list of constituents
Definition: Jet.cc:359
int weight
Definition: histoStyle.py:50
Jet::EtaPhiMoments Jet::etaPhiStatistics ( ) const

eta-phi statistics, ET weighted

eta-phi statistics

Definition at line 184 of file Jet.cc.

References reco::deltaPhi(), reco::LeafCandidate::eta(), reco::Jet::EtaPhiMoments::etaEtaMoment, reco::Jet::EtaPhiMoments::etaMean, reco::Jet::EtaPhiMoments::etaPhiMoment, getJetConstituentsQuick(), i, reco::LeafCandidate::phi(), reco::Jet::EtaPhiMoments::phiMean, reco::Jet::EtaPhiMoments::phiPhiMoment, query::result, and histoStyle::weight.

184  {
185  std::vector<const Candidate*> towers = getJetConstituentsQuick ();
186  double phiRef = phi();
187  double sumw = 0;
188  double sumEta = 0;
189  double sumEta2 = 0;
190  double sumPhi = 0;
191  double sumPhi2 = 0;
192  double sumEtaPhi = 0;
193  int i = towers.size();
194  while (--i >= 0) {
195  double eta = towers[i]->eta();
196  double phi = deltaPhi (towers[i]->phi(), phiRef);
197  double weight = towers[i]->et();
198  sumw += weight;
199  sumEta += eta * weight;
200  sumEta2 += eta * eta * weight;
201  sumPhi += phi * weight;
202  sumPhi2 += phi * phi * weight;
203  sumEtaPhi += eta * phi * weight;
204  }
206  if (sumw > 0) {
207  result.etaMean = sumEta / sumw;
208  result.phiMean = deltaPhi (phiRef + sumPhi, 0.);
209  result.etaEtaMoment = (sumEta2 - sumEta * sumEta / sumw) / sumw;
210  result.phiPhiMoment = (sumPhi2 - sumPhi * sumPhi / sumw) / sumw;
211  result.etaPhiMoment = (sumEtaPhi - sumEta * sumPhi / sumw) / sumw;
212  }
213  else {
214  result.etaMean = 0;
215  result.phiMean = 0;
216  result.etaEtaMoment = 0;
217  result.phiPhiMoment = 0;
218  result.etaPhiMoment = 0;
219  }
220  return result;
221 }
int i
Definition: DBlmapReader.cc:9
virtual float phi() const GCC11_FINAL
momentum azimuthal angle
record to store eta-phi first and second moments
Definition: Jet.h:27
tuple result
Definition: query.py:137
virtual float eta() const GCC11_FINAL
momentum pseudorapidity
double deltaPhi(double phi1, double phi2)
Definition: deltaPhi.h:12
virtual std::vector< const reco::Candidate * > getJetConstituentsQuick() const
quick list of constituents
Definition: Jet.cc:359
int weight
Definition: histoStyle.py:50
float Jet::etInAnnulus ( float  fRmin,
float  fRmax 
) const

ET in annulus between rmin and rmax around jet direction.

et in annulus between rmin and rmax around jet direction

Definition at line 280 of file Jet.cc.

References reco::deltaR(), getJetConstituentsQuick(), i, and query::result.

280  {
281  float result = 0;
282  std::vector<const Candidate*> towers = getJetConstituentsQuick ();
283  int i = towers.size ();
284  while (--i >= 0) {
285  double r = deltaR (*this, *(towers[i]));
286  if (r >= fRmin && r < fRmax) result += towers[i]->et ();
287  }
288  return result;
289 }
int i
Definition: DBlmapReader.cc:9
tuple result
Definition: query.py:137
auto deltaR(const T1 &t1, const T2 &t2) -> decltype(t1.eta())
Definition: deltaR.h:30
virtual std::vector< const reco::Candidate * > getJetConstituentsQuick() const
quick list of constituents
Definition: Jet.cc:359
Jet::Constituents Jet::getJetConstituents ( ) const
virtual
std::vector< const Candidate * > Jet::getJetConstituentsQuick ( ) const
virtual

quick list of constituents

Definition at line 359 of file Jet.cc.

References reco::CompositePtrCandidate::daughter(), i, reco::CompositePtrCandidate::numberOfDaughters(), and query::result.

Referenced by FWPFTauProxyBuilder::buildViewType(), etaetaMoment(), etaphiMoment(), etaPhiStatistics(), etInAnnulus(), maxDistance(), nCarrying(), and phiphiMoment().

359  {
360  std::vector<const Candidate*> result;
361  int nDaughters = numberOfDaughters();
362  for (int i = 0; i < nDaughters; ++i) {
363  result.push_back (daughter (i));
364  }
365  return result;
366 }
int i
Definition: DBlmapReader.cc:9
virtual size_t numberOfDaughters() const
number of daughters
tuple result
Definition: query.py:137
virtual const Candidate * daughter(size_type) const
return daughter at a given position, i = 0, ... numberOfDaughters() - 1 (read only mode) ...
bool Jet::isJet ( ) const
virtual

Reimplemented from reco::LeafCandidate.

Definition at line 449 of file Jet.cc.

449  {
450  return true;
451 }
virtual float reco::Jet::jetArea ( ) const
inlinevirtual

get jet area

Definition at line 106 of file Jet.h.

References mJetArea.

Referenced by L1FastjetCorrector::correction(), PileupJetIdProducer::produce(), and HiL1Subtractor::produce().

106 {return mJetArea;}
float mJetArea
Definition: Jet.h:121
float Jet::maxDistance ( ) const

maximum distance from jet to constituent

Definition at line 307 of file Jet.cc.

References reco::deltaR(), PFRecoTauDiscriminationAgainstElectronDeadECAL_cfi::dR, getJetConstituentsQuick(), i, and query::result.

307  {
308  float result = 0;
309  std::vector<const Candidate*> towers = getJetConstituentsQuick ();
310  for (unsigned i = 0; i < towers.size(); ++i) {
311  float dR = deltaR (*(towers[i]), *this);
312  if (dR > result) result = dR;
313  }
314  return result;
315 }
int i
Definition: DBlmapReader.cc:9
tuple result
Definition: query.py:137
auto deltaR(const T1 &t1, const T2 &t2) -> decltype(t1.eta())
Definition: deltaR.h:30
virtual std::vector< const reco::Candidate * > getJetConstituentsQuick() const
quick list of constituents
Definition: Jet.cc:359
int Jet::nCarrying ( float  fFraction) const

return # of constituent carrying fraction of energy

return # of constituent carring fraction of energy. Assume ordered towers

Definition at line 292 of file Jet.cc.

References reco::LeafCandidate::et(), getJetConstituentsQuick(), i, and query::result.

Referenced by reco::CaloJet::n60(), pat::Jet::n60(), reco::CaloJet::n90(), and pat::Jet::n90().

292  {
293  std::vector<const Candidate*> towers = getJetConstituentsQuick ();
294  if (fFraction >= 1) return towers.size();
295  double totalEt = 0;
296  for (unsigned i = 0; i < towers.size(); ++i) totalEt += towers[i]->et();
297  double fractionEnergy = totalEt * fFraction;
298  unsigned result = 0;
299  for (; result < towers.size(); ++result) {
300  fractionEnergy -= towers[result]->et();
301  if (fractionEnergy <= 0) return result+1;
302  }
303  return 0;
304 }
virtual double et() const GCC11_FINAL
transverse energy
int i
Definition: DBlmapReader.cc:9
tuple result
Definition: query.py:137
virtual std::vector< const reco::Candidate * > getJetConstituentsQuick() const
quick list of constituents
Definition: Jet.cc:359
virtual int reco::Jet::nConstituents ( ) const
inlinevirtual

of constituents

Definition at line 66 of file Jet.h.

References reco::CompositePtrCandidate::numberOfDaughters().

Referenced by FFTGenericScaleCalculator::mapFFTJet(), print(), and JetIdSelector< T >::produce().

66 {return numberOfDaughters();}
virtual size_t numberOfDaughters() const
number of daughters
virtual int reco::Jet::nPasses ( ) const
inlinevirtual

number of passes taken by algorithm

Definition at line 116 of file Jet.h.

References mPassNumber.

116 {return mPassNumber;}
int mPassNumber
Definition: Jet.h:123
float Jet::phiphiMoment ( ) const

phi-phi second moment, ET weighted

phi-phi second moment

Definition at line 241 of file Jet.cc.

References reco::deltaPhi(), getJetConstituentsQuick(), i, reco::LeafCandidate::phi(), relativeConstraints::value, and histoStyle::weight.

Referenced by fireworks::makeEveJetCone(), and JetIDSelectionFunctor::operator()().

241  {
242  double phiRef = phi();
243  std::vector<const Candidate*> towers = getJetConstituentsQuick ();
244  double sumw = 0;
245  double sum = 0;
246  double sum2 = 0;
247  int i = towers.size();
248  while (--i >= 0) {
249  double value = deltaPhi (towers[i]->phi(), phiRef);
250  double weight = towers[i]->et();
251  sumw += weight;
252  sum += value * weight;
253  sum2 += value * value * weight;
254  }
255  return sumw > 0 ? (sum2 - sum*sum/sumw ) / sumw : 0;
256 }
int i
Definition: DBlmapReader.cc:9
virtual float phi() const GCC11_FINAL
momentum azimuthal angle
double deltaPhi(double phi1, double phi2)
Definition: deltaPhi.h:12
virtual std::vector< const reco::Candidate * > getJetConstituentsQuick() const
quick list of constituents
Definition: Jet.cc:359
int weight
Definition: histoStyle.py:50
float Jet::physicsEta ( float  fZVertex,
float  fDetectorEta 
)
static

static function to convert detector eta to physics eta

Definition at line 319 of file Jet.cc.

Referenced by JetTracksAssociationDRVertexAssigned::produce().

319  {
320  CaloPointZ refPoint (0., fDetectorEta);
321  return refPoint.etaReference (fZVertex);
322 }
Candidate::LorentzVector Jet::physicsP4 ( const Candidate::Point newVertex,
const Candidate inParticle,
const Candidate::Point oldVertex = Candidate::Point(0,0,0) 
)
static

Definition at line 331 of file Jet.cc.

References reco::Candidate::energy(), reco::Candidate::momentum(), reco::LeafCandidate::p(), and p3.

Referenced by reco::CaloJet::physicsP4().

331  {
332  CaloPoint3D<Point> caloPoint(oldVertex,inParticle.momentum()); // Jet position in Calo.
333  Vector physicsDir = caloPoint.caloPoint() - newVertex;
334  double p = inParticle.momentum().r();
335  Vector p3 = p * physicsDir.unit();
336  LorentzVector returnVector(p3.x(), p3.y(), p3.z(), inParticle.energy());
337  return returnVector;
338 }
virtual double energy() const =0
energy
virtual double p() const GCC11_FINAL
magnitude of momentum vector
math::XYZTLorentzVector LorentzVector
virtual Vector momentum() const =0
spatial momentum vector
fixed size vector
Definition: Vector.h:23
double p3[4]
Definition: TauolaWrapper.h:91
virtual float reco::Jet::pileup ( ) const
inlinevirtual

pileup energy contribution as calculated by algorithm

Definition at line 111 of file Jet.h.

References mPileupEnergy.

111 {return mPileupEnergy;}
float mPileupEnergy
Definition: Jet.h:122
std::string Jet::print ( void  ) const
virtual

Print object.

Reimplemented in reco::PFJet, reco::CaloJet, reco::JPTJet, reco::GenJet, reco::TrackJet, reco::PFClusterJet, and reco::BasicJet.

Definition at line 425 of file Jet.cc.

References reco::CompositePtrCandidate::daughterPtr(), reco::LeafCandidate::eta(), edm::Ptr< T >::id(), edm::Ptr< T >::isNonnull(), edm::Ptr< T >::key(), nConstituents(), reco::CompositePtrCandidate::numberOfDaughters(), dbtoconf::out, reco::LeafCandidate::p(), reco::LeafCandidate::phi(), reco::LeafCandidate::pt(), reco::LeafCandidate::px(), reco::LeafCandidate::py(), and reco::LeafCandidate::pz().

Referenced by reco::BasicJet::print(), reco::PFClusterJet::print(), reco::TrackJet::print(), reco::GenJet::print(), reco::JPTJet::print(), reco::CaloJet::print(), and reco::PFJet::print().

425  {
426  std::ostringstream out;
427  out << "Jet p/px/py/pz/pt: " << p() << '/' << px () << '/' << py() << '/' << pz() << '/' << pt() << std::endl
428  << " eta/phi: " << eta () << '/' << phi () << std::endl
429  << " # of constituents: " << nConstituents () << std::endl;
430  out << " Constituents:" << std::endl;
431  for (unsigned index = 0; index < numberOfDaughters(); index++) {
432  Constituent constituent = daughterPtr (index); // deref
433  if (constituent.isNonnull()) {
434  out << " #" << index << " p/pt/eta/phi: "
435  << constituent->p() << '/' << constituent->pt() << '/' << constituent->eta() << '/' << constituent->phi()
436  << " productId/index: " << constituent.id() << '/' << constituent.key() << std::endl;
437  }
438  else {
439  out << " #" << index << " constituent is not available in the event" << std::endl;
440  }
441  }
442  return out.str ();
443 }
virtual double p() const GCC11_FINAL
magnitude of momentum vector
CandidatePtr daughterPtr(size_type i) const
reference to daughter at given position
virtual double pz() const GCC11_FINAL
z coordinate of momentum vector
virtual double py() const GCC11_FINAL
y coordinate of momentum vector
virtual float phi() const GCC11_FINAL
momentum azimuthal angle
virtual size_t numberOfDaughters() const
number of daughters
virtual double px() const GCC11_FINAL
x coordinate of momentum vector
unsigned int index
Definition: LeafCandidate.h:34
virtual float eta() const GCC11_FINAL
momentum pseudorapidity
tuple out
Definition: dbtoconf.py:99
edm::Ptr< Candidate > Constituent
Definition: Jet.h:23
virtual int nConstituents() const
of constituents
Definition: Jet.h:66
virtual float pt() const GCC11_FINAL
transverse momentum
void Jet::scaleEnergy ( double  fScale)
virtual
virtual void reco::Jet::setJetArea ( float  fArea)
inlinevirtual
virtual void reco::Jet::setNPasses ( int  fPasses)
inlinevirtual

Set number of passes taken by algorithm.

Definition at line 114 of file Jet.h.

References mPassNumber.

114 {mPassNumber = fPasses;}
int mPassNumber
Definition: Jet.h:123
virtual void reco::Jet::setPileup ( float  fEnergy)
inlinevirtual

Set pileup energy contribution as calculated by algorithm.

Definition at line 109 of file Jet.h.

References mPileupEnergy.

Referenced by HiL1Subtractor::produce(), FastjetJetProducer::produceTrackJets(), and cms::SubEventGenJetProducer::runAlgorithm().

109 {mPileupEnergy = fEnergy;}
float mPileupEnergy
Definition: Jet.h:122

Member Data Documentation

float reco::Jet::mJetArea
private

Definition at line 121 of file Jet.h.

Referenced by jetArea(), and setJetArea().

int reco::Jet::mPassNumber
private

Definition at line 123 of file Jet.h.

Referenced by nPasses(), and setNPasses().

float reco::Jet::mPileupEnergy
private

Definition at line 122 of file Jet.h.

Referenced by pileup(), and setPileup().