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
ElectronMomentumCorrector Class Reference

#include <ElectronMomentumCorrector.h>

Public Member Functions

void correct (reco::GsfElectron &, TrajectoryStateOnSurface &)
 
 ElectronMomentumCorrector ()
 
 ~ElectronMomentumCorrector ()
 

Private Member Functions

float energyError (float E, float *par) const
 

Detailed Description

Definition at line 19 of file ElectronMomentumCorrector.h.

Constructor & Destructor Documentation

ElectronMomentumCorrector::ElectronMomentumCorrector ( )
inline

Definition at line 25 of file ElectronMomentumCorrector.h.

25 {}
ElectronMomentumCorrector::~ElectronMomentumCorrector ( )
inline

Definition at line 27 of file ElectronMomentumCorrector.h.

27 {}

Member Function Documentation

void ElectronMomentumCorrector::correct ( reco::GsfElectron &  electron,
TrajectoryStateOnSurface &  vtxTsos 
)

Definition at line 29 of file ElectronMomentumCorrector.cc.

References reco::GsfElectron::BADTRACK, reco::GsfElectron::BIGBREM, reco::GsfElectron::classification(), reco::GsfElectron::correctedEcalEnergy(), reco::GsfElectron::correctedEcalEnergyError(), reco::GsfElectron::correctMomentum(), funct::false, reco::GsfElectron::GAP, reco::GsfElectron::GOLDEN, reco::GsfElectron::isEB(), reco::GsfElectron::isEE(), GaussianSumUtilities1D::mode(), MultiGaussianStateTransform::multiState1D(), reco::GsfElectron::p4(), reco::GsfElectron::P4_COMBINATION, reco::GsfElectron::p4Error(), hitfit::return, pileupReCalc_HLTpaths::scale, reco::GsfElectron::SHOWERING, mathSSE::sqrt(), reco::btau::trackMomentum, reco::GsfElectron::trackMomentumAtVtx(), funct::true, reco::GsfElectron::UNKNOWN, and SingleGaussianState1D::variance().

Referenced by GsfElectronAlgo::createElectron().

30  {
31  if (electron.p4Error(reco::GsfElectron::P4_COMBINATION)!= 999.)
32  {
33  edm::LogWarning("ElectronMomentumCorrector::correct")<<"already done" ;
34  return ;
35  }
36 
37  //math::XYZTLorentzVector newMomentum = electron.p4() ; // default
38  int elClass = electron.classification() ;
39 
40  // irrelevant classification
41  if ( (elClass <= reco::GsfElectron::UNKNOWN) ||
42  (elClass>reco::GsfElectron::GAP) )
43  {
44  edm::LogWarning("ElectronMomentumCorrector::correct")<<"unexpected classification" ;
45  return ;
46  }
47 
48 
49  //=======================================================================================
50  // cluster energy
51  //=======================================================================================
52 
53  float scEnergy = electron.correctedEcalEnergy() ;
54  float errorEnergy = electron.correctedEcalEnergyError() ;
55 
56 
57  //=======================================================================================
58  // track momentum
59  //=======================================================================================
60 
61  // basic values
62  float trackMomentum = electron.trackMomentumAtVtx().R() ;
63  //float errorTrackMomentum = 999. ;
64 
65  // tracker momentum scale corrections (Mykhailo Dalchenko)
66  double scale=1.;
67  if (electron.isEB()){
68  if (elClass==0) {scale = 1./(0.00104*sqrt(trackMomentum)+1);}
69  if (elClass==1) {scale = 1./(0.0017*sqrt(trackMomentum)+0.9986);}
70  if (elClass==3) {scale = 1./(1.004 - 0.00021*trackMomentum);}
71  if (elClass==4) {scale = 0.995;}
72  } else if (electron.isEE()){
73  if (elClass==3){scale = 1./(1.01432-0.00201872*trackMomentum+0.0000142621*trackMomentum*trackMomentum);}
74  if (elClass==4){scale = 1./(0.996859-0.000345347*trackMomentum);}
75  }
76  if (scale<0.) scale = 1.; // CC added protection
77  trackMomentum = trackMomentum*scale ;
78 
79  // error (must be done after trackMomentum rescaling)
80  MultiGaussianState1D qpState(MultiGaussianStateTransform::multiState1D(vtxTsos,0)) ;
81  GaussianSumUtilities1D qpUtils(qpState) ;
82  float errorTrackMomentum = trackMomentum*trackMomentum*sqrt(qpUtils.mode().variance()) ;
83 
84 
85  //=======================================================================================
86  // combination
87  //=======================================================================================
88 
89  float finalMomentum = electron.p4().t(); // initial
90  float finalMomentumError = 999.;
91 
92  // first check for large errors
93 
94  if (errorTrackMomentum/trackMomentum > 0.5 && errorEnergy/scEnergy <= 0.5) {
95  finalMomentum = scEnergy; finalMomentumError = errorEnergy;
96  }
97  else if (errorTrackMomentum/trackMomentum <= 0.5 && errorEnergy/scEnergy > 0.5){
98  finalMomentum = trackMomentum; finalMomentumError = errorTrackMomentum;
99  }
100  else if (errorTrackMomentum/trackMomentum > 0.5 && errorEnergy/scEnergy > 0.5){
101  if (errorTrackMomentum/trackMomentum < errorEnergy/scEnergy) {
102  finalMomentum = trackMomentum; finalMomentumError = errorTrackMomentum;
103  }
104  else{
105  finalMomentum = scEnergy; finalMomentumError = errorEnergy;
106  }
107  }
108 
109  // then apply the combination algorithm
110  else {
111 
112  // calculate E/p and corresponding error
113  float eOverP = scEnergy / trackMomentum;
114  float errorEOverP = sqrt(
115  (errorEnergy/trackMomentum)*(errorEnergy/trackMomentum) +
116  (scEnergy*errorTrackMomentum/trackMomentum/trackMomentum)*
117  (scEnergy*errorTrackMomentum/trackMomentum/trackMomentum));
118 
119 // if ( eOverP > 1 + 2.5*errorEOverP )
120 // {
121 // finalMomentum = scEnergy; finalMomentumError = errorEnergy;
122 // if ((elClass==reco::GsfElectron::GOLDEN) && electron.isEB() && (eOverP<1.15))
123 // {
124 // if (scEnergy<15) {finalMomentum = trackMomentum ; finalMomentumError = errorTrackMomentum;}
125 // }
126 // }
127 // else if ( eOverP < 1 - 2.5*errorEOverP )
128 // {
129 // finalMomentum = scEnergy; finalMomentumError = errorEnergy;
130 // if (elClass==reco::GsfElectron::SHOWERING)
131 // {
132 // if (electron.isEB())
133 // {
134 // if (scEnergy<18)
135 // { finalMomentum = trackMomentum; finalMomentumError = errorTrackMomentum; }
136 // }
137 // else if (electron.isEE())
138 // {
139 // if (scEnergy<13)
140 // {finalMomentum = trackMomentum; finalMomentumError = errorTrackMomentum;}
141 // }
142 // else
143 // { edm::LogWarning("ElectronMomentumCorrector::correct")<<"nor barrel neither endcap electron ?!" ; }
144 // }
145 // else if (electron.isGap())
146 // {
147 // if (scEnergy<60)
148 // { finalMomentum = trackMomentum; finalMomentumError = errorTrackMomentum ; }
149 // }
150 // }
151 
152  bool eleIsNotInCombination = false ;
153  if ( (eOverP > 1 + 2.5*errorEOverP) || (eOverP < 1 - 2.5*errorEOverP) || (eOverP < 0.8) || (eOverP > 1.3) )
154  { eleIsNotInCombination = true ; }
155  if (eleIsNotInCombination)
156  {
157  if (eOverP > 1)
158  { finalMomentum = scEnergy ; finalMomentumError = errorEnergy ; }
159  else
160  {
161  if (elClass == reco::GsfElectron::GOLDEN)
162  { finalMomentum = scEnergy; finalMomentumError = errorEnergy; }
163  if (elClass == reco::GsfElectron::BIGBREM)
164  {
165  if (scEnergy<36)
166  { finalMomentum = trackMomentum ; finalMomentumError = errorTrackMomentum ; }
167  else
168  { finalMomentum = scEnergy ; finalMomentumError = errorEnergy ; }
169  }
170  if (elClass == reco::GsfElectron::BADTRACK)
171  { finalMomentum = scEnergy; finalMomentumError = errorEnergy ; }
172  if (elClass == reco::GsfElectron::SHOWERING)
173  {
174  if (scEnergy<30)
175  { finalMomentum = trackMomentum ; finalMomentumError = errorTrackMomentum; }
176  else
177  { finalMomentum = scEnergy; finalMomentumError = errorEnergy;}
178  }
179  if (elClass == reco::GsfElectron::GAP)
180  {
181  if (scEnergy<60)
182  { finalMomentum = trackMomentum ; finalMomentumError = errorTrackMomentum ; }
183  else
184  { finalMomentum = scEnergy; finalMomentumError = errorEnergy ; }
185  }
186  }
187  }
188 
189  else
190  {
191  // combination
192  finalMomentum = (scEnergy/errorEnergy/errorEnergy + trackMomentum/errorTrackMomentum/errorTrackMomentum) /
193  (1/errorEnergy/errorEnergy + 1/errorTrackMomentum/errorTrackMomentum);
194  float finalMomentumVariance = 1 / (1/errorEnergy/errorEnergy + 1/errorTrackMomentum/errorTrackMomentum);
195  finalMomentumError = sqrt(finalMomentumVariance);
196  }
197 
198  }
199 
200 
201  //=======================================================================================
202  // final set
203  //=======================================================================================
204 
205  math::XYZTLorentzVector oldMomentum = electron.p4() ;
206  math::XYZTLorentzVector newMomentum = math::XYZTLorentzVector
207  ( oldMomentum.x()*finalMomentum/oldMomentum.t(),
208  oldMomentum.y()*finalMomentum/oldMomentum.t(),
209  oldMomentum.z()*finalMomentum/oldMomentum.t(),
210  finalMomentum ) ;
211 
212  electron.correctMomentum(newMomentum,errorTrackMomentum,finalMomentumError);
213 
214  }
reco::GsfElectron::p4
const LorentzVector & p4(P4Kind kind) const
Definition: GsfElectron.cc:204
reco::GsfElectron::correctMomentum
void correctMomentum(const LorentzVector &p4, float trackMomentumError, float p4Error)
Definition: GsfElectron.h:734
reco::GsfElectron::trackMomentumAtVtx
math::XYZVectorF trackMomentumAtVtx() const
Definition: GsfElectron.h:272
reco::GsfElectron::p4Error
float p4Error(P4Kind kind) const
Definition: GsfElectron.cc:216
reco::GsfElectron::isEE
bool isEE() const
Definition: GsfElectron.h:331
reco::GsfElectron::isEB
bool isEB() const
Definition: GsfElectron.h:330
math::XYZTLorentzVector
XYZTLorentzVectorD XYZTLorentzVector
Lorentz vector with cylindrical internal representation using pseudorapidity.
Definition: LorentzVector.h:30
MultiGaussianStateTransform::multiState1D
MultiGaussianState1D multiState1D(const std::vector< MultiGaussianState< N >::Vector > &, const std::vector< MultiGaussianState< N >::Matrix > &, const std::vector< double > &, unsigned int)
Definition: MultiGaussianStateTransform.cc:71
pileupReCalc_HLTpaths.scale
int scale
Definition: pileupReCalc_HLTpaths.py:111
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:48
reco::GsfElectron::correctedEcalEnergy
float correctedEcalEnergy() const
Definition: GsfElectron.h:713
reco::GsfElectron::classification
Classification classification() const
Definition: GsfElectron.h:635
reco::GsfElectron::correctedEcalEnergyError
float correctedEcalEnergyError() const
Definition: GsfElectron.h:714
float ElectronMomentumCorrector::energyError ( float  E,
float *  par 
) const
private
Generated for CMSSW Reference Manual by doxygen 1.8.5