CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
MET.cc
Go to the documentation of this file.
1 // File: MET.cc
2 // Description: see MET.h
3 // Author: Michael Schmitt, R. Cavanaugh University of Florida
4 // Creation Date: MHS MAY 30, 2005 initial version
5 
7 #include "TVector.h"
8 #include "TMatrix.h"
9 
10 using namespace std;
11 using namespace reco;
12 
13 //---------------------------------------------------------------------------
14 // Default Constructor;
15 //-----------------------------------
16 MET::MET()
17 {
18  sumet = 0.0;
19  elongit = 0.0;
20  signif_dxx=signif_dyy=signif_dyx=signif_dxy=0.;
21 }
22 //---------------------------------------------------------------------------
23 
24 //---------------------------------------------------------------------------
25 // Constructer for the case when only p4_ = (mEx, mEy, 0, mEt) is known.
26 // The vertex information is currently not used (but may be in the future)
27 // and is required by the RecoCandidate constructer.
28 //-----------------------------------
29 MET::MET( const LorentzVector& p4_, const Point& vtx_ ) :
30  RecoCandidate( 0, p4_, vtx_ )
31 {
32  sumet = 0.0;
33  elongit = 0.0;
35 }
36 //---------------------------------------------------------------------------
37 
38 //---------------------------------------------------------------------------
39 // Constructer for the case when the SumET is known in addition to
40 // p4_ = (mEx, mEy, 0, mEt). The vertex information is currently not used
41 // (but see above).
42 //-----------------------------------
43 MET::MET( double sumet_, const LorentzVector& p4_, const Point& vtx_ ) :
44  RecoCandidate( 0, p4_, vtx_ )
45 {
46  sumet = sumet_;
47  elongit = 0.0;
49 }
50 //---------------------------------------------------------------------------
51 
52 //---------------------------------------------------------------------------
53 // Constructor for the case when the SumET, the corrections which
54 // were applied to the MET, as well the MET itself p4_ = (mEx, mEy, 0, mEt)
55 // are all known. See above concerning the vertex information.
56 //-----------------------------------
57 MET::MET( double sumet_, std::vector<CorrMETData> corr_,
58  const LorentzVector& p4_, const Point& vtx_ ) :
59  RecoCandidate( 0, p4_, vtx_ )
60 {
61  sumet = sumet_;
62  elongit = 0.0;
64  //-----------------------------------
65  // Fill the vector containing the corrections (corr) with vector of
66  // known corrections (corr_) passed in via the constructor.
67  std::vector<CorrMETData>::const_iterator i;
68  for( i = corr_.begin(); i != corr_.end(); i++ )
69  {
70  corr.push_back( *i );
71  }
72 }
73 //---------------------------------------------------------------------------
74 
75 
76 //-----------------------------------------------------------------
77 //explicit clone function
78 MET * MET::clone() const {
79  return new MET( * this );
80 }
81 //-----------------------------------------------------------------
82 
83 // function that calculates the MET significance from the vector information.
84 double MET::significance() const {
85  if(signif_dxx==0 && signif_dyy==0 && signif_dxy==0 && signif_dyx==0)
86  return -1;
87  TMatrixD metmat = getSignificanceMatrix();
88  TVectorD metvec(2);
89  metvec(0)=this->px();
90  metvec(1)=this->py();
91  double signif = -1;
92  if(std::fabs(metmat.Determinant())>0.000001){
93  metmat.Invert();
94  signif = metvec * (metmat * metvec);
95  }
96  return signif;
97 }
98 
99 //---------------------------------------------------------------------------
100 // Returns the vector of all corrections applied to the x component of the
101 // missing transverse momentum, mEx
102 //-----------------------------------
103 std::vector<double> MET::dmEx() const
104 {
105  std::vector<double> deltas;
106  std::vector<CorrMETData>::const_iterator i;
107  for( i = corr.begin(); i != corr.end(); i++ )
108  {
109  deltas.push_back( i->mex );
110  }
111  return deltas;
112 }
113 //---------------------------------------------------------------------------
114 
115 //---------------------------------------------------------------------------
116 // Returns the vector of all corrections applied to the y component of the
117 // missing transverse momentum, mEy
118 //-----------------------------------
119 std::vector<double> MET::dmEy() const
120 {
121  std::vector<double> deltas;
122  std::vector<CorrMETData>::const_iterator i;
123  for( i = corr.begin(); i != corr.end(); i++ )
124  {
125  deltas.push_back( i->mey );
126  }
127  return deltas;
128 }
129 //---------------------------------------------------------------------------
130 
131 //---------------------------------------------------------------------------
132 // Returns the vector of all corrections applied to the scalar sum of the
133 // transverse energy (over all objects)
134 //-----------------------------------
135 std::vector<double> MET::dsumEt() const
136 {
137  std::vector<double> deltas;
138  std::vector<CorrMETData>::const_iterator i;
139  for( i = corr.begin(); i != corr.end(); i++ )
140  {
141  deltas.push_back( i->sumet );
142  }
143  return deltas;
144 }
145 
146 //---------------------------------------------------------------------------
147 // returns the significance matrix
148 //---------------------------------------------------------------------------
149 
150 TMatrixD MET::getSignificanceMatrix(void) const
151 {
152  TMatrixD result(2,2);
153  result(0,0)=signif_dxx;
154  result(0,1)=signif_dxy;
155  result(1,0)=signif_dyx;
156  result(1,1)=signif_dyy;
157  return result;
158 }
159 
160 //---------------------------------------------------------------------------
161 // Returns the vector of all corrections applied to the significance of the
162 // transverse energy (over all objects)
163 //-----------------------------------
164 std::vector<double> MET::dSignificance() const
165 {
166  std::vector<double> deltas;
167  std::vector<CorrMETData>::const_iterator i;
168  for( i = corr.begin(); i != corr.end(); i++ )
169  {
170  deltas.push_back( i->significance );
171  }
172  return deltas;
173 }
174 //---------------------------------------------------------------------------
175 
176 //---------------------------------------------------------------------------
177 // Required RecoCandidate polymorphism
178 //-----------------------------------
179 bool MET::overlap( const Candidate & ) const
180 {
181  return false;
182 }
183 //---------------------------------------------------------------------------
184 
185 void
186 MET::setSignificanceMatrix(const TMatrixD &inmatrix){
187  signif_dxx=inmatrix(0,0);
188  signif_dxy=inmatrix(0,1);
189  signif_dyx=inmatrix(1,0);
190  signif_dyy=inmatrix(1,1);
191  return;
192 }
std::vector< double > dsumEt() const
Definition: MET.cc:135
int i
Definition: DBlmapReader.cc:9
double signif_dyx
Definition: MET.h:74
MET * clone() const
returns a clone of the Candidate object
Definition: MET.cc:78
double signif_dxy
Definition: MET.h:75
double signif_dxx
Definition: MET.h:72
std::vector< double > dmEx() const
Definition: MET.cc:103
MET()
Definition: MET.cc:16
TMatrixD getSignificanceMatrix(void) const
Definition: MET.cc:150
Definition: MET.h:32
tuple result
Definition: query.py:137
std::vector< CorrMETData > corr
Definition: MET.h:76
void setSignificanceMatrix(const TMatrixD &matrix)
Definition: MET.cc:186
std::vector< double > dmEy() const
Definition: MET.cc:119
double signif_dyy
Definition: MET.h:73
virtual double px() const
x coordinate of momentum vector
math::XYZTLorentzVector LorentzVector
Lorentz vector.
Definition: Candidate.h:39
virtual bool overlap(const Candidate &) const
check overlap with another candidate
Definition: MET.cc:179
std::vector< double > dSignificance() const
Definition: MET.cc:164
math::XYZPoint Point
point in the space
Definition: Candidate.h:43
double elongit
Definition: MET.h:70
double significance() const
Definition: MET.cc:84
double sumet
Definition: MET.h:69
virtual double py() const
y coordinate of momentum vector