CMS 3D CMS Logo

List of all members | Public Member Functions | Private Member Functions | Private Attributes | Static Private Attributes
reco::helper::CastorJetIDHelper Class Reference

#include <CastorJetIDHelper.h>

Public Member Functions

void calculate (const edm::Event &event, const reco::BasicJet &jet)
 
 CastorJetIDHelper ()
 
double depth () const
 
double emEnergy () const
 
double fem () const
 
double fhot () const
 
double hadEnergy () const
 
void initValues ()
 
int nTowers () const
 
double sigmaz () const
 
double width () const
 
 ~CastorJetIDHelper ()
 

Private Member Functions

double phiangle (double testphi)
 

Private Attributes

double depth_
 
double emEnergy_
 
double fem_
 
double fhot_
 
double hadEnergy_
 
int nTowers_
 
double sigmaz_
 
double width_
 

Static Private Attributes

static int sanity_checks_left_
 

Detailed Description

Definition at line 13 of file CastorJetIDHelper.h.

Constructor & Destructor Documentation

reco::helper::CastorJetIDHelper::CastorJetIDHelper ( )

Definition at line 14 of file CastorJetIDHelper.cc.

References initValues().

15 {
16 
17  initValues();
18 }
reco::helper::CastorJetIDHelper::~CastorJetIDHelper ( )
inline

Definition at line 18 of file CastorJetIDHelper.h.

References calculate(), initValues(), and metsig::jet.

18 {}

Member Function Documentation

void reco::helper::CastorJetIDHelper::calculate ( const edm::Event event,
const reco::BasicJet jet 
)

Definition at line 33 of file CastorJetIDHelper.cc.

References reco::CastorTower::depth(), depth_, reco::CastorTower::emEnergy(), emEnergy_, reco::LeafCandidate::energy(), fem_, reco::CastorTower::fhot(), fhot_, reco::Jet::getJetConstituents(), reco::CastorTower::hadEnergy(), hadEnergy_, edm::Ref< C, T, F >::id(), initValues(), module::module(), nTowers_, reco::LeafCandidate::phi(), phiangle(), funct::pow(), reco::CastorTower::rechitsBegin(), reco::CastorTower::rechitsEnd(), sigmaz_, mathSSE::sqrt(), and width_.

Referenced by CastorJetIDProducer::produce(), and ~CastorJetIDHelper().

34 {
35  initValues();
36 
37  // calculate Castor JetID properties
38 
39  double zmean = 0.;
40  double z2mean = 0.;
41 
42  std::vector<CandidatePtr> ccp = jet.getJetConstituents();
43  std::vector<CandidatePtr>::const_iterator itParticle;
44  for (itParticle=ccp.begin();itParticle!=ccp.end();++itParticle){
45  const CastorTower* castorcand = dynamic_cast<const CastorTower*>(itParticle->get());
46  emEnergy_ += castorcand->emEnergy();
47  hadEnergy_ += castorcand->hadEnergy();
48  depth_ += castorcand->depth()*castorcand->energy();
49  width_ += pow(phiangle(castorcand->phi() - jet.phi()),2)*castorcand->energy();
50  fhot_ += castorcand->fhot()*castorcand->energy();
51 
52  // loop over rechits
53  for (edm::RefVector<edm::SortedCollection<CastorRecHit> >::iterator it = castorcand->rechitsBegin(); it != castorcand->rechitsEnd(); it++) {
55  double Erechit = rechit_p->energy();
56  HcalCastorDetId id = rechit_p->id();
57  int module = id.module();
58  double zrechit = 0;
59  if (module < 3) zrechit = -14390 - 24.75 - 49.5*(module-1);
60  if (module > 2) zrechit = -14390 - 99 - 49.5 - 99*(module-3);
61  zmean += Erechit*zrechit;
62  z2mean += Erechit*zrechit*zrechit;
63  } // end loop over rechits
64 
65  nTowers_++;
66  }
67  //cout << "" << endl;
68 
69  depth_ = depth_/jet.energy();
70  width_ = sqrt(width_/jet.energy());
71  fhot_ = fhot_/jet.energy();
72  fem_ = emEnergy_/jet.energy();
73 
74  zmean = zmean/jet.energy();
75  z2mean = z2mean/jet.energy();
76  double sigmaz2 = z2mean - zmean*zmean;
77  if(sigmaz2 > 0) sigmaz_ = sqrt(sigmaz2);
78 
79 
80 }
module()
Definition: vlib.cc:994
virtual Constituents getJetConstituents() const
list of constituents
ProductID id() const
Accessor for product ID.
Definition: Ref.h:257
T sqrt(T t)
Definition: SSEVec.h:18
double energy() const final
energy
double phiangle(double testphi)
double phi() const final
momentum azimuthal angle
Definition: vlib.h:208
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:40
double reco::helper::CastorJetIDHelper::depth ( void  ) const
inline

Definition at line 32 of file CastorJetIDHelper.h.

References depth_.

Referenced by CastorJetIDProducer::produce().

double reco::helper::CastorJetIDHelper::emEnergy ( ) const
inline

Definition at line 28 of file CastorJetIDHelper.h.

References emEnergy_.

Referenced by CastorJetIDProducer::produce().

double reco::helper::CastorJetIDHelper::fem ( ) const
inline

Definition at line 30 of file CastorJetIDHelper.h.

References fem_.

Referenced by CastorJetIDProducer::produce().

double reco::helper::CastorJetIDHelper::fhot ( ) const
inline

Definition at line 33 of file CastorJetIDHelper.h.

References fhot_.

Referenced by CastorJetIDProducer::produce().

double reco::helper::CastorJetIDHelper::hadEnergy ( ) const
inline

Definition at line 29 of file CastorJetIDHelper.h.

References hadEnergy_.

Referenced by CastorJetIDProducer::produce().

void reco::helper::CastorJetIDHelper::initValues ( )
int reco::helper::CastorJetIDHelper::nTowers ( ) const
inline

Definition at line 35 of file CastorJetIDHelper.h.

References nTowers_, and phiangle().

Referenced by CastorJetIDProducer::produce().

double reco::helper::CastorJetIDHelper::phiangle ( double  testphi)
private

Definition at line 83 of file CastorJetIDHelper.cc.

References M_PI.

Referenced by calculate(), and nTowers().

83  {
84  double phi = testphi;
85  while (phi>M_PI) phi -= (2*M_PI);
86  while (phi<-M_PI) phi += (2*M_PI);
87  return phi;
88 }
#define M_PI
double reco::helper::CastorJetIDHelper::sigmaz ( ) const
inline

Definition at line 34 of file CastorJetIDHelper.h.

References sigmaz_.

Referenced by CastorJetIDProducer::produce().

double reco::helper::CastorJetIDHelper::width ( ) const
inline

Member Data Documentation

double reco::helper::CastorJetIDHelper::depth_
private

Definition at line 48 of file CastorJetIDHelper.h.

Referenced by calculate(), depth(), and initValues().

double reco::helper::CastorJetIDHelper::emEnergy_
private

Definition at line 44 of file CastorJetIDHelper.h.

Referenced by calculate(), emEnergy(), and initValues().

double reco::helper::CastorJetIDHelper::fem_
private

Definition at line 46 of file CastorJetIDHelper.h.

Referenced by calculate(), fem(), and initValues().

double reco::helper::CastorJetIDHelper::fhot_
private

Definition at line 49 of file CastorJetIDHelper.h.

Referenced by calculate(), fhot(), and initValues().

double reco::helper::CastorJetIDHelper::hadEnergy_
private

Definition at line 45 of file CastorJetIDHelper.h.

Referenced by calculate(), hadEnergy(), and initValues().

int reco::helper::CastorJetIDHelper::nTowers_
private

Definition at line 51 of file CastorJetIDHelper.h.

Referenced by calculate(), initValues(), and nTowers().

int reco::helper::CastorJetIDHelper::sanity_checks_left_
staticprivate

Definition at line 53 of file CastorJetIDHelper.h.

double reco::helper::CastorJetIDHelper::sigmaz_
private

Definition at line 50 of file CastorJetIDHelper.h.

Referenced by calculate(), initValues(), and sigmaz().

double reco::helper::CastorJetIDHelper::width_
private

Definition at line 47 of file CastorJetIDHelper.h.

Referenced by calculate(), initValues(), and width().