CMS 3D CMS Logo

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

#include <TopHLTDQMHelper.h>

Public Member Functions

 CalculateHLT (int maxNJets, double wMass)
 default constructor More...
 
double masslb (reco::RecoCandidate *mu, const reco::MET &met, const reco::Jet &b)
 calculate mlb estimate More...
 
double massTopQuark (const std::vector< reco::Jet > &jets)
 calculate top quark mass estimate More...
 
double massWBoson (const std::vector< reco::Jet > &jets)
 calculate W boson mass estimate More...
 
double tmassTopQuark (reco::RecoCandidate *mu, const reco::MET &met, const reco::Jet &b)
 calculate top quark transverse mass estimate More...
 
double tmassWBoson (reco::RecoCandidate *mu, const reco::MET &met, const reco::Jet &b)
 calculate W boson transverse mass estimate More...
 
 ~CalculateHLT ()
 default destructor More...
 

Private Member Functions

void operator() (const std::vector< reco::Jet > &jets)
 
void operator() (const reco::Jet &bJet, reco::RecoCandidate *lepton, const reco::MET &met)
 

Private Attributes

bool failed_
 indicate failed associations More...
 
double massTopQuark_
 cache of top quark mass estimate More...
 
double massWBoson_
 cache of w boson mass estimate More...
 
int maxNJets_
 max. number of jets to be considered More...
 
double mlb_
 cache of mlb estimate More...
 
double tmassTopQuark_
 cache of top quark transverse mass estimate More...
 
double tmassWBoson_
 cache of W boson transverse mass estimate More...
 
double wMass_
 paramater of the w boson mass More...
 

Detailed Description

Definition at line 67 of file TopHLTDQMHelper.h.

Constructor & Destructor Documentation

CalculateHLT::CalculateHLT ( int  maxNJets,
double  wMass 
)

default constructor

Definition at line 7 of file TopHLTDQMHelper.cc.

7  :
8  failed_(false), maxNJets_(maxNJets), wMass_(wMass), massWBoson_(-1.), massTopQuark_(-1.),tmassWBoson_(-1),tmassTopQuark_(-1),mlb_(-1)
9 {
10 }
bool failed_
indicate failed associations
double wMass_
paramater of the w boson mass
double massWBoson_
cache of w boson mass estimate
double mlb_
cache of mlb estimate
double massTopQuark_
cache of top quark mass estimate
double tmassWBoson_
cache of W boson transverse mass estimate
int maxNJets_
max. number of jets to be considered
double tmassTopQuark_
cache of top quark transverse mass estimate
CalculateHLT::~CalculateHLT ( )
inline

default destructor

Definition at line 72 of file TopHLTDQMHelper.h.

72 {};

Member Function Documentation

double CalculateHLT::masslb ( reco::RecoCandidate mu,
const reco::MET met,
const reco::Jet b 
)

calculate mlb estimate

Definition at line 66 of file TopHLTDQMHelper.cc.

References failed_, mlb_, and operator()().

Referenced by TopHLTSingleLepton::MonitorEnsemble::fill().

67 {
68  if(!failed_&& mlb_<0) operator()(b,mu,met); return mlb_;
69 }
bool failed_
indicate failed associations
void operator()(const std::vector< reco::Jet > &jets)
double mlb_
cache of mlb estimate
double CalculateHLT::massTopQuark ( const std::vector< reco::Jet > &  jets)

calculate top quark mass estimate

Definition at line 21 of file TopHLTDQMHelper.cc.

References failed_, massTopQuark_, and operator()().

Referenced by TopHLTSingleLepton::MonitorEnsemble::fill().

22 {
24 }
bool failed_
indicate failed associations
void operator()(const std::vector< reco::Jet > &jets)
vector< PseudoJet > jets
double massTopQuark_
cache of top quark mass estimate
double CalculateHLT::massWBoson ( const std::vector< reco::Jet > &  jets)

calculate W boson mass estimate

Definition at line 14 of file TopHLTDQMHelper.cc.

References failed_, massWBoson_, and operator()().

Referenced by TopHLTSingleLepton::MonitorEnsemble::fill().

15 {
16  if(!failed_&& massWBoson_<0) operator()(jets); return massWBoson_;
17 }
bool failed_
indicate failed associations
void operator()(const std::vector< reco::Jet > &jets)
double massWBoson_
cache of w boson mass estimate
vector< PseudoJet > jets
void CalculateHLT::operator() ( const std::vector< reco::Jet > &  jets)
private

do the calculation; this is called only once per event by the first function call to return a mass estimate. The once calculated values are cached afterwards

Definition at line 94 of file TopHLTDQMHelper.cc.

References failed_, scaleCards::mass, massTopQuark_, massWBoson_, maxNJets_, and wMass_.

Referenced by masslb(), massTopQuark(), massWBoson(), tmassTopQuark(), and tmassWBoson().

95 {
96 
97  if(maxNJets_<0) maxNJets_=jets.size();
98  failed_= jets.size()<(unsigned int) maxNJets_;
99  if( failed_){ return; }
100 
101  // associate those jets with maximum pt of the vectorial
102  // sum to the hadronic decay chain
103  double maxPt=-1.;
104  std::vector<int> maxPtIndices;
105  maxPtIndices.push_back(-1);
106  maxPtIndices.push_back(-1);
107  maxPtIndices.push_back(-1);
108  for(int idx=0; idx<maxNJets_; ++idx){
109  for(int jdx=0; jdx<maxNJets_; ++jdx){ if(jdx<=idx) continue;
110  for(int kdx=0; kdx<maxNJets_; ++kdx){ if(kdx==idx || kdx==jdx) continue;
111  reco::Particle::LorentzVector sum = jets[idx].p4()+jets[jdx].p4()+jets[kdx].p4();
112  if( maxPt<0. || maxPt<sum.pt() ){
113  maxPt=sum.pt();
114  maxPtIndices.clear();
115  maxPtIndices.push_back(idx);
116  maxPtIndices.push_back(jdx);
117  maxPtIndices.push_back(kdx);
118  }
119  }
120  }
121  }
122  massTopQuark_= (jets[maxPtIndices[0]].p4()+
123  jets[maxPtIndices[1]].p4()+
124  jets[maxPtIndices[2]].p4()).mass();
125 
126  // associate those jets that get closest to the W mass
127  // with their invariant mass to the W boson
128  double wDist =-1.;
129  std::vector<int> wMassIndices;
130  wMassIndices.push_back(-1);
131  wMassIndices.push_back(-1);
132  for(unsigned idx=0; idx<maxPtIndices.size(); ++idx){
133  for(unsigned jdx=0; jdx<maxPtIndices.size(); ++jdx){
134  if( jdx==idx || maxPtIndices[idx]>maxPtIndices[jdx] ) continue;
135  reco::Particle::LorentzVector sum = jets[maxPtIndices[idx]].p4()+jets[maxPtIndices[jdx]].p4();
136  if( wDist<0. || wDist>fabs(sum.mass()-wMass_) ){
137  wDist=fabs(sum.mass()-wMass_);
138  wMassIndices.clear();
139  wMassIndices.push_back(maxPtIndices[idx]);
140  wMassIndices.push_back(maxPtIndices[jdx]);
141  }
142  }
143  }
144  massWBoson_= (jets[wMassIndices[0]].p4()+
145  jets[wMassIndices[1]].p4()).mass();
146 }
bool failed_
indicate failed associations
double wMass_
paramater of the w boson mass
double massWBoson_
cache of w boson mass estimate
vector< PseudoJet > jets
double massTopQuark_
cache of top quark mass estimate
int maxNJets_
max. number of jets to be considered
tuple mass
Definition: scaleCards.py:27
math::XYZTLorentzVector LorentzVector
Lorentz vector.
Definition: Particle.h:25
void CalculateHLT::operator() ( const reco::Jet bJet,
reco::RecoCandidate lepton,
const reco::MET met 
)
private

Definition at line 79 of file TopHLTDQMHelper.cc.

References gather_cfg::cout, mlb_, reco::LeafCandidate::p4(), funct::pow(), reco::LeafCandidate::pt(), reco::LeafCandidate::px(), reco::LeafCandidate::py(), mathSSE::sqrt(), tmassTopQuark_, and tmassWBoson_.

79  {
80  double metT = sqrt(pow(met.px(),2) + pow(met.py(),2));
81  double lepT = sqrt(pow(lepton->px(),2) + pow(lepton->py(),2));
82  double bT = sqrt(pow(bJet.px(),2) + pow(bJet.py(),2));
83  reco::Particle::LorentzVector WT = lepton->p4() + met.p4();
84  cout<<"in calculate:\n\t"<<bJet.pt()<<"\t"<<lepton->pt()<<"\t"<<met.pt()<<endl;
85  tmassWBoson_ = sqrt(pow(metT+lepT,2) - (WT.px()*WT.px()) - (WT.py()*WT.py()));
86  reco::Particle::LorentzVector topT = WT + bJet.p4();
87  tmassTopQuark_ = sqrt(pow((metT+lepT+bT),2) - (topT.px()*topT.px()) - (topT.py()*topT.py()));
88  reco::Particle::LorentzVector lb = bJet.p4() + lepton->p4();
89  mlb_ = lb.mass();
90 }
double mlb_
cache of mlb estimate
T sqrt(T t)
Definition: SSEVec.h:46
virtual double px() const
x coordinate of momentum vector
virtual double pt() const
transverse momentum
double tmassWBoson_
cache of W boson transverse mass estimate
tuple cout
Definition: gather_cfg.py:121
double tmassTopQuark_
cache of top quark transverse mass estimate
virtual const LorentzVector & p4() const
four-momentum Lorentz vector
math::XYZTLorentzVector LorentzVector
Lorentz vector.
Definition: Particle.h:25
virtual double py() const
y coordinate of momentum vector
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:40
double CalculateHLT::tmassTopQuark ( reco::RecoCandidate mu,
const reco::MET met,
const reco::Jet b 
)

calculate top quark transverse mass estimate

Definition at line 73 of file TopHLTDQMHelper.cc.

References failed_, operator()(), and tmassTopQuark_.

Referenced by TopHLTSingleLepton::MonitorEnsemble::fill().

74 {
75  if(!failed_&& tmassTopQuark_<0) operator()(b,lepton,met); return tmassTopQuark_;
76 }
bool failed_
indicate failed associations
void operator()(const std::vector< reco::Jet > &jets)
double tmassTopQuark_
cache of top quark transverse mass estimate
double CalculateHLT::tmassWBoson ( reco::RecoCandidate mu,
const reco::MET met,
const reco::Jet b 
)

calculate W boson transverse mass estimate

calculate W boson transverse mass estimate

Definition at line 59 of file TopHLTDQMHelper.cc.

References failed_, operator()(), and tmassWBoson_.

Referenced by TopHLTSingleLepton::MonitorEnsemble::fill().

60 {
61  if(!failed_&& tmassWBoson_<0) operator()(b,mu,met); return tmassWBoson_;
62 }
bool failed_
indicate failed associations
void operator()(const std::vector< reco::Jet > &jets)
double tmassWBoson_
cache of W boson transverse mass estimate

Member Data Documentation

bool CalculateHLT::failed_
private

indicate failed associations

Definition at line 100 of file TopHLTDQMHelper.h.

Referenced by masslb(), massTopQuark(), massWBoson(), operator()(), tmassTopQuark(), and tmassWBoson().

double CalculateHLT::massTopQuark_
private

cache of top quark mass estimate

Definition at line 108 of file TopHLTDQMHelper.h.

Referenced by massTopQuark(), and operator()().

double CalculateHLT::massWBoson_
private

cache of w boson mass estimate

Definition at line 106 of file TopHLTDQMHelper.h.

Referenced by massWBoson(), and operator()().

int CalculateHLT::maxNJets_
private

max. number of jets to be considered

Definition at line 102 of file TopHLTDQMHelper.h.

Referenced by operator()().

double CalculateHLT::mlb_
private

cache of mlb estimate

Definition at line 114 of file TopHLTDQMHelper.h.

Referenced by masslb(), and operator()().

double CalculateHLT::tmassTopQuark_
private

cache of top quark transverse mass estimate

Definition at line 112 of file TopHLTDQMHelper.h.

Referenced by operator()(), and tmassTopQuark().

double CalculateHLT::tmassWBoson_
private

cache of W boson transverse mass estimate

Definition at line 110 of file TopHLTDQMHelper.h.

Referenced by operator()(), and tmassWBoson().

double CalculateHLT::wMass_
private

paramater of the w boson mass

Definition at line 104 of file TopHLTDQMHelper.h.

Referenced by operator()().