#include <ElectronMomentumCorrector.h>
Public Member Functions | |
void | correct (reco::GsfElectron &, TrajectoryStateOnSurface &) |
ElectronMomentumCorrector () | |
~ElectronMomentumCorrector () | |
Private Member Functions | |
float | energyError (float E, float *par) const |
Definition at line 19 of file ElectronMomentumCorrector.h.
ElectronMomentumCorrector::ElectronMomentumCorrector | ( | ) | [inline] |
Definition at line 25 of file ElectronMomentumCorrector.h.
{}
ElectronMomentumCorrector::~ElectronMomentumCorrector | ( | ) | [inline] |
Definition at line 27 of file ElectronMomentumCorrector.h.
{}
void ElectronMomentumCorrector::correct | ( | reco::GsfElectron & | electron, |
TrajectoryStateOnSurface & | vtxTsos | ||
) |
Definition at line 33 of file ElectronMomentumCorrector.cc.
References reco::GsfElectron::classification(), reco::GsfElectron::correctMomentum(), reco::GsfElectron::ecalEnergy(), reco::GsfElectron::ecalEnergyError(), reco::GsfElectron::GAP, reco::GsfElectron::GOLDEN, reco::GsfElectron::isEB(), reco::GsfElectron::isEE(), reco::GsfElectron::isGap(), GaussianSumUtilities1D::mode(), MultiGaussianStateTransform::multiState1D(), reco::GsfElectron::p4(), reco::GsfElectron::P4_COMBINATION, reco::GsfElectron::p4Error(), reco::GsfElectron::SHOWERING, mathSSE::sqrt(), reco::btau::trackMomentum, reco::GsfElectron::trackMomentumAtVtx(), reco::GsfElectron::UNKNOWN, and SingleGaussianState1D::variance().
Referenced by GsfElectronAlgo::createElectron().
{ if (electron.p4Error(reco::GsfElectron::P4_COMBINATION)!= 999.) { edm::LogWarning("ElectronMomentumCorrector::correct")<<"already done" ; return ; } math::XYZTLorentzVector newMomentum = electron.p4() ; // default int elClass = electron.classification() ; // irrelevant classification if ( (elClass <= reco::GsfElectron::UNKNOWN) || (elClass>reco::GsfElectron::GAP) ) { edm::LogWarning("ElectronMomentumCorrector::correct")<<"unexpected classification" ; return ; } float scEnergy = electron.ecalEnergy() ; float errorEnergy_ = electron.ecalEnergyError() ; float trackMomentum = electron.trackMomentumAtVtx().R() ; float errorTrackMomentum = 999. ; // retreive momentum error MultiGaussianState1D qpState(MultiGaussianStateTransform::multiState1D(vtxTsos,0)); GaussianSumUtilities1D qpUtils(qpState); errorTrackMomentum = trackMomentum*trackMomentum*sqrt(qpUtils.mode().variance()); float finalMomentum = electron.p4().t(); // initial float finalMomentumError = 999.; // first check for large errors if (errorTrackMomentum/trackMomentum > 0.5 && errorEnergy_/scEnergy <= 0.5) { finalMomentum = scEnergy; finalMomentumError = errorEnergy_; } else if (errorTrackMomentum/trackMomentum <= 0.5 && errorEnergy_/scEnergy > 0.5){ finalMomentum = trackMomentum; finalMomentumError = errorTrackMomentum; } else if (errorTrackMomentum/trackMomentum > 0.5 && errorEnergy_/scEnergy > 0.5){ if (errorTrackMomentum/trackMomentum < errorEnergy_/scEnergy) { finalMomentum = trackMomentum; finalMomentumError = errorTrackMomentum; } else{ finalMomentum = scEnergy; finalMomentumError = errorEnergy_; } } // then apply the combination algorithm else { // calculate E/p and corresponding error float eOverP = scEnergy / trackMomentum; float errorEOverP = sqrt( (errorEnergy_/trackMomentum)*(errorEnergy_/trackMomentum) + (scEnergy*errorTrackMomentum/trackMomentum/trackMomentum)* (scEnergy*errorTrackMomentum/trackMomentum/trackMomentum)); if ( eOverP > 1 + 2.5*errorEOverP ) { finalMomentum = scEnergy; finalMomentumError = errorEnergy_; if ((elClass==reco::GsfElectron::GOLDEN) && electron.isEB() && (eOverP<1.15)) { if (scEnergy<15) {finalMomentum = trackMomentum ; finalMomentumError = errorTrackMomentum;} } } else if ( eOverP < 1 - 2.5*errorEOverP ) { finalMomentum = scEnergy; finalMomentumError = errorEnergy_; if (elClass==reco::GsfElectron::SHOWERING) { if (electron.isEB()) { if (scEnergy<18) { finalMomentum = trackMomentum; finalMomentumError = errorTrackMomentum; } } else if (electron.isEE()) { if (scEnergy<13) {finalMomentum = trackMomentum; finalMomentumError = errorTrackMomentum;} } else { edm::LogWarning("ElectronMomentumCorrector::correct")<<"nor barrel neither endcap electron ?!" ; } } else if (electron.isGap()) { if (scEnergy<60) { finalMomentum = trackMomentum; finalMomentumError = errorTrackMomentum ; } } } else { // combination finalMomentum = (scEnergy/errorEnergy_/errorEnergy_ + trackMomentum/errorTrackMomentum/errorTrackMomentum) / (1/errorEnergy_/errorEnergy_ + 1/errorTrackMomentum/errorTrackMomentum); float finalMomentumVariance = 1 / (1/errorEnergy_/errorEnergy_ + 1/errorTrackMomentum/errorTrackMomentum); finalMomentumError = sqrt(finalMomentumVariance); } } math::XYZTLorentzVector oldMomentum = electron.p4() ; newMomentum = math::XYZTLorentzVector ( oldMomentum.x()*finalMomentum/oldMomentum.t(), oldMomentum.y()*finalMomentum/oldMomentum.t(), oldMomentum.z()*finalMomentum/oldMomentum.t(), finalMomentum ) ; // final set electron.correctMomentum(newMomentum,errorTrackMomentum,finalMomentumError); }
float ElectronMomentumCorrector::energyError | ( | float | E, |
float * | par | ||
) | const [private] |