Skip to content

Commit

Permalink
combination regression
Browse files Browse the repository at this point in the history
  • Loading branch information
@beaudett @lgray
beaudett authored and lgray committed Dec 16, 2013
1 parent 1c96715 commit 91ec9d6
Show file tree
Hide file tree
Showing 4 changed files with 104 additions and 7 deletions.
1 change: 1 addition & 0 deletions RecoEgamma/EgammaElectronAlgos/interface/RegressionHelper.h
View file
Expand Up @@ -47,6 +47,7 @@ class RegressionHelper {
const edm::Handle<EcalRecHitCollection>& rechitsEB,
const edm::Handle<EcalRecHitCollection>& rechitsEE) const;

void applyCombinationRegression(reco::GsfElectron & ele) const;
~RegressionHelper() ;

private:
Expand Down
4 changes: 2 additions & 2 deletions RecoEgamma/EgammaElectronAlgos/src/GsfElectronAlgo.cc
View file
Expand Up @@ -1342,6 +1342,7 @@ void GsfElectronAlgo::createElectron()
tcMatching, tkExtra, ctfInfo,
fiducialFlags,showerShape,
conversionVars ) ;
// Will be overwritten later in the case of the regression
ele->setCorrectedEcalEnergyError(generalData_->superClusterErrorFunction->getValue(*(ele->superCluster()),0)) ;
ele->setP4(GsfElectron::P4_FROM_SUPER_CLUSTER,momentum,0,true) ;

Expand All @@ -1367,7 +1368,6 @@ void GsfElectronAlgo::createElectron()
//====================================================
// classification and corrections
//====================================================

// classification
ElectronClassification theClassifier ;
theClassifier.classify(*ele) ;
Expand Down Expand Up @@ -1400,7 +1400,7 @@ void GsfElectronAlgo::createElectron()
// momentum
if(generalData_->strategyCfg.useCombinationRegression) // new
{

generalData_->regHelper->applyCombinationRegression(*ele);
}
else if (ele->core()->ecalDrivenSeed()) //original computation
{
Expand Down
92 changes: 91 additions & 1 deletion RecoEgamma/EgammaElectronAlgos/src/RegressionHelper.cc
View file
Expand Up @@ -22,7 +22,8 @@ void RegressionHelper::applyEcalRegression(reco::GsfElectron & ele,
double cor, err;
getEcalRegression(*ele.superCluster(), vertices, rechitsEB, rechitsEE, cor, err);
ele.setCorrectedEcalEnergy(cor * ele.superCluster()->energy());
ele.setCorrectedEcalEnergyError( err * ele.superCluster()->energy());
ele.setCorrectedEcalEnergyError( err * ele.superCluster()->energy());

}

void RegressionHelper::checkSetup(const edm::EventSetup & es) {
Expand Down Expand Up @@ -323,3 +324,92 @@ void RegressionHelper::getEcalRegression(const reco::SuperCluster & sc,
<< "Supercluster seed is either EB nor EE!" << std::endl;
}
}

void RegressionHelper::applyCombinationRegression(reco::GsfElectron & ele) const
{
float energy = ele.correctedEcalEnergy();
float energyError = ele.correctedEcalEnergyError();
float momentum = ele.trackMomentumAtVtx().R();
float momentumError = ele.trackMomentumError();
int elClass = -.1;

switch (ele.classification()) {
case reco::GsfElectron::GOLDEN:
elClass = 0;
break;
case reco::GsfElectron::BIGBREM:
elClass = 1;
break;
case reco::GsfElectron::BADTRACK:
elClass = 2;
break;
case reco::GsfElectron::SHOWERING:
elClass = 3;
break;
case reco::GsfElectron::GAP:
elClass = 4;
break;
default:
elClass = -1;


bool isEcalDriven = ele.ecalDriven();
bool isTrackerDriven = ele.trackerDrivenSeed();
bool isEB = ele.isEB();

// compute relative errors and ratio of errors
float energyRelError = energyError / energy;
float momentumRelError = momentumError / momentum;
float errorRatio = energyRelError / momentumRelError;

// calculate E/p and corresponding error
float eOverP = energy / momentum;
float eOverPerror = sqrt(
(energyError/momentum)*(energyError/momentum) +
(energy*momentumError/momentum/momentum)*
(energy*momentumError/momentum/momentum));

// fill input variables
std::vector<float> regressionInputs ;
regressionInputs.resize(11,0.);

regressionInputs[0] = energy;
regressionInputs[1] = energyRelError;
regressionInputs[2] = momentum;
regressionInputs[3] = momentumRelError;
regressionInputs[4] = errorRatio;
regressionInputs[5] = eOverP;
regressionInputs[6] = eOverPerror;
regressionInputs[7] = static_cast<float>(isEcalDriven);
regressionInputs[8] = static_cast<float>(isTrackerDriven);
regressionInputs[9] = static_cast<float>(elClass);
regressionInputs[10] = static_cast<float>(isEB);

// retrieve combination weight
float weight = 0.;
if(eOverP>0.025
&&fabs(momentum-energy)<15.*sqrt(momentumError*momentumError + energyError*energyError)
) // protection against crazy track measurement
{
weight = combinationReg_->GetResponse(&regressionInputs[0]);
if(weight>1.) weight = 1.;
else if(weight<0.) weight = 0.;
}

float combinedMomentum = weight*momentum + (1.-weight)*energy;
float combinedMomentumError = sqrt(weight*weight*momentumError*momentumError + (1.-weight)*(1.-weight)*energyError*energyError);

// FIXME : pure tracker electrons have track momentum error of 999.
// If the combination try to combine such electrons then the original combined momentum is kept
if(momentumError!=999. || weight==0.)
{
math::XYZTLorentzVector oldMomentum = ele.p4() ;
math::XYZTLorentzVector newMomentum(oldMomentum.x()*combinedMomentum/oldMomentum.t(),
oldMomentum.y()*combinedMomentum/oldMomentum.t(),
oldMomentum.z()*combinedMomentum/oldMomentum.t(),
combinedMomentum);

ele.setP4(reco::GsfElectron::P4_COMBINATION,newMomentum,combinedMomentumError,true);
}
}
}
14 changes: 10 additions & 4 deletions RecoEgamma/EgammaElectronProducers/python/gedGsfElectrons_cfi.py
View file
Expand Up @@ -30,6 +30,11 @@
tag = cms.string('gedelectron_EEUncertainty_offline'),
label = cms.untracked.string('gedelectron_EEUncertainty_offline')
),
cms.PSet(
record = cms.string('GBRWrapperRcd'),
tag = cms.string('gedelectron_p4combination_offline'),
label = cms.untracked.string('gedelectron_p4combination_offline')
),
)
)
gedelectronGBRESSource.connect = cms.string('frontier://FrontierProd/CMS_COND_PAT_000')
Expand All @@ -40,7 +45,8 @@
GBRWrapperRcd = cms.vstring('GBRForest/gedelectron_EBCorrection_offline',
'GBRForest/gedelectron_EECorrection_offline',
'GBRForest/gedelectron_EBUncertainty_offline',
'GBRForest/gedelectron_EEUncertainty_offline')
'GBRForest/gedelectron_EEUncertainty_offline',
'GBRForest/gedelectron_p4combination_offline')
)

gedGsfElectronsTmp = cms.EDProducer("GEDGsfElectronProducer",
Expand Down Expand Up @@ -82,7 +88,7 @@
ambClustersOverlapStrategy = cms.uint32(1),
addPflowElectrons = cms.bool(True), # this one should be transfered to the "core" level
useEcalRegression = cms.bool(True),
useCombinationRegression = cms.bool(False),
useCombinationRegression = cms.bool(True),

# preselection parameters (ecal driven electrons)
minSCEtBarrel = cms.double(4.0),
Expand Down Expand Up @@ -189,12 +195,12 @@
'gedelectron_EECorrection_offline',
'gedelectron_EBUncertainty_offline',
'gedelectron_EEUncertainty_offline'),
combinationRegressionWeightLabels = cms.vstring(),
combinationRegressionWeightLabels = cms.vstring('gedelectron_p4combination_offline'),

ecalWeightsFromDB = cms.bool(True),
# if not from DB. Otherwise, keep empty
ecalRefinedRegressionWeightFiles = cms.vstring(),
combinationWeightsFromDB = cms.bool(False),
combinationWeightsFromDB = cms.bool(True),
# if not from DB. Otherwise, keep empty
combinationRegressionWeightFile = cms.vstring(),

Expand Down

0 comments on commit 91ec9d6

Please sign in to comment.