Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Make Ecal TimeBias correction to use a new GT. #2397

Merged
merged 2 commits into from
Mar 4, 2014
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
4 changes: 0 additions & 4 deletions Configuration/GlobalRuns/python/reco_TLR_311X.py
Original file line number Diff line number Diff line change
Expand Up @@ -49,10 +49,6 @@ def customisePPData(process):
process.hfreco.S8S1stat.flagsToSkip = 18
process.hfreco.S9S1stat.flagsToSkip = 26

##Ecal time bias correction
process.ecalGlobalUncalibRecHit.doEBtimeCorrection = True
process.ecalGlobalUncalibRecHit.doEEtimeCorrection = True

return process


Expand Down
4 changes: 0 additions & 4 deletions Configuration/GlobalRuns/python/reco_TLR_42X.py
Original file line number Diff line number Diff line change
Expand Up @@ -164,10 +164,6 @@ def customisePPData(process):
process.hfreco.S8S1stat.flagsToSkip = 18
process.hfreco.S9S1stat.flagsToSkip = 26

##Ecal time bias correction
process.ecalGlobalUncalibRecHit.doEBtimeCorrection = True
process.ecalGlobalUncalibRecHit.doEEtimeCorrection = True

return process


Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -11,6 +11,3 @@
HcalRemoveAddSevLevel.AddFlag(hcalRecAlgos,"HBHESpikeNoise",12)

CSCHaloData.ExpectedBX = cms.int32(3)

ecalGlobalUncalibRecHit.doEBtimeCorrection = True
ecalGlobalUncalibRecHit.doEEtimeCorrection = True
5 changes: 0 additions & 5 deletions HLTrigger/Configuration/python/customizeHLTforMC.py
Original file line number Diff line number Diff line change
Expand Up @@ -9,11 +9,6 @@ def customizeHLTforMC(process):
#if 'CSCHaloData' in process.__dict__:
# process.CSCHaloData.ExpectedBX = cms.int32(6)

# EcalUncalibRecHitProducer - not used at HLT
#if 'ecalGlobalUncalibRecHit' in process.__dict__:
# process.ecalGlobalUncalibRecHit.doEBtimeCorrection = cms.bool(False)
# process.ecalGlobalUncalibRecHit.doEEtimeCorrection = cms.bool(False)

# HcalRecAlgoESProducer - these flags are not used at HLT (they should stay set to the default value for both data and MC)
#if 'hcalRecAlgos' in process.__dict__:
# import RecoLocalCalo.HcalRecAlgos.RemoveAddSevLevel as HcalRemoveAddSevLevel
Expand Down
165 changes: 64 additions & 101 deletions RecoLocalCalo/EcalRecProducers/plugins/EcalUncalibRecHitWorkerGlobal.cc
Original file line number Diff line number Diff line change
Expand Up @@ -12,6 +12,7 @@
#include "CondFormats/DataRecord/interface/EcalSampleMaskRcd.h"
#include "CondFormats/DataRecord/interface/EcalTimeCalibConstantsRcd.h"
#include "CondFormats/DataRecord/interface/EcalTimeOffsetConstantRcd.h"
#include "CondFormats/DataRecord/interface/EcalTimeBiasCorrectionsRcd.h"

EcalUncalibRecHitWorkerGlobal::EcalUncalibRecHitWorkerGlobal(const edm::ParameterSet&ps,edm::ConsumesCollector& c) :
EcalUncalibRecHitWorkerBaseClass(ps,c)
Expand Down Expand Up @@ -39,21 +40,6 @@ EcalUncalibRecHitWorkerGlobal::EcalUncalibRecHitWorkerGlobal(const edm::Paramete
outOfTimeThreshG61mEE_ = ps.getParameter<double>("outOfTimeThresholdGain61mEE");
amplitudeThreshEB_ = ps.getParameter<double>("amplitudeThresholdEB");
amplitudeThreshEE_ = ps.getParameter<double>("amplitudeThresholdEE");
// amplitude-dependent correction of time
doEBtimeCorrection_ = ps.getParameter<bool>("doEBtimeCorrection");
doEEtimeCorrection_ = ps.getParameter<bool>("doEEtimeCorrection");
EBtimeCorrAmplitudeBins_ = ps.getParameter<std::vector<double> >("EBtimeCorrAmplitudeBins");
EBtimeCorrShiftBins_ = ps.getParameter<std::vector<double> >("EBtimeCorrShiftBins");
EEtimeCorrAmplitudeBins_ = ps.getParameter<std::vector<double> >("EEtimeCorrAmplitudeBins");
EEtimeCorrShiftBins_ = ps.getParameter<std::vector<double> >("EEtimeCorrShiftBins");
if(EBtimeCorrAmplitudeBins_.size() != EBtimeCorrShiftBins_.size()) {
doEBtimeCorrection_ = false;
edm::LogError("EcalRecHitError") << "Size of EBtimeCorrAmplitudeBins different from EBtimeCorrShiftBins. Forcing no time corrections for EB. ";
}
if(EEtimeCorrAmplitudeBins_.size() != EEtimeCorrShiftBins_.size()) {
doEEtimeCorrection_ = false;
edm::LogError("EcalRecHitError") << "Size of EEtimeCorrAmplitudeBins different from EEtimeCorrShiftBins. Forcing no time corrections for EE. ";
}

// spike threshold
ebSpikeThresh_ = ps.getParameter<double>("ebSpikeThreshold");
Expand Down Expand Up @@ -97,21 +83,6 @@ EcalUncalibRecHitWorkerGlobal::EcalUncalibRecHitWorkerGlobal(const edm::Paramete
outOfTimeThreshG61mEE_ = ps.getParameter<double>("outOfTimeThresholdGain61mEE");
amplitudeThreshEB_ = ps.getParameter<double>("amplitudeThresholdEB");
amplitudeThreshEE_ = ps.getParameter<double>("amplitudeThresholdEE");
// amplitude-dependent correction of time
doEBtimeCorrection_ = ps.getParameter<bool>("doEBtimeCorrection");
doEEtimeCorrection_ = ps.getParameter<bool>("doEEtimeCorrection");
EBtimeCorrAmplitudeBins_ = ps.getParameter<std::vector<double> >("EBtimeCorrAmplitudeBins");
EBtimeCorrShiftBins_ = ps.getParameter<std::vector<double> >("EBtimeCorrShiftBins");
EEtimeCorrAmplitudeBins_ = ps.getParameter<std::vector<double> >("EEtimeCorrAmplitudeBins");
EEtimeCorrShiftBins_ = ps.getParameter<std::vector<double> >("EEtimeCorrShiftBins");
if(EBtimeCorrAmplitudeBins_.size() != EBtimeCorrShiftBins_.size()) {
doEBtimeCorrection_ = false;
edm::LogError("EcalRecHitError") << "Size of EBtimeCorrAmplitudeBins different from EBtimeCorrShiftBins. Forcing no time corrections for EB. ";
}
if(EEtimeCorrAmplitudeBins_.size() != EEtimeCorrShiftBins_.size()) {
doEEtimeCorrection_ = false;
edm::LogError("EcalRecHitError") << "Size of EEtimeCorrAmplitudeBins different from EEtimeCorrShiftBins. Forcing no time corrections for EE. ";
}

// spike threshold
ebSpikeThresh_ = ps.getParameter<double>("ebSpikeThreshold");
Expand Down Expand Up @@ -158,6 +129,9 @@ EcalUncalibRecHitWorkerGlobal::set(const edm::EventSetup& es)
// for the leading edge method
es.get<EcalTimeCalibConstantsRcd>().get(itime);
es.get<EcalTimeOffsetConstantRcd>().get(offtime);

// for the time correction methods
es.get<EcalTimeBiasCorrectionsRcd>().get(timeCorrBias_);
}


Expand All @@ -177,84 +151,71 @@ int EcalUncalibRecHitWorkerGlobal::isSaturated(const C & dataFrame)
return -1; // no saturation found
}


double EcalUncalibRecHitWorkerGlobal::timeCorrectionEB(float ampliEB){
// computed initially in ns. Than turned in the BX's, as EcalUncalibratedRecHit need be.
double theCorrection=0;


int myBin = -1;
for (int bin=0; bin<(int)EBtimeCorrAmplitudeBins_.size(); bin++ ){
if(ampliEB > EBtimeCorrAmplitudeBins_.at(bin)) {
myBin = bin; }
else break;
/**
* Amplitude-dependent time corrections; EE and EB have separate corrections:
* EXtimeCorrAmplitudes (ADC) and EXtimeCorrShifts (ns) need to have the same number of elements
* Bins must be ordered in amplitude. First-last bins take care of under-overflows.
*
* The algorithm is the same for EE and EB, only the correction vectors are different.
*
* @return Jitter (in clock cycles) which will be added to UncalibRechit.setJitter(), 0 if no correction is applied.
*/
double EcalUncalibRecHitWorkerGlobal::timeCorrection(
float ampli,
const std::vector<float>& amplitudeBins,
const std::vector<float>& shiftBins) {

// computed initially in ns. Than turned in the BX's, as
// EcalUncalibratedRecHit need be.
double theCorrection = 0;

// sanity check for arrays
if (amplitudeBins.size() == 0) {
edm::LogError("EcalRecHitError")
<< "timeCorrAmplitudeBins is empty, forcing no time bias corrections.";

return 0;
}

if (myBin == -1)
{
theCorrection = EBtimeCorrShiftBins_.at(0);
}
else if ( myBin == ((int)(EBtimeCorrAmplitudeBins_.size()-1)) )
{
theCorrection = EBtimeCorrShiftBins_.at( myBin );
}
else if ( -1 < myBin && myBin < ((int)EBtimeCorrAmplitudeBins_.size()-1) )
{
// interpolate linearly between two assingned points
theCorrection = ( EBtimeCorrShiftBins_.at(myBin+1) - EBtimeCorrShiftBins_.at(myBin) );
theCorrection *= ( ((double)ampliEB) - EBtimeCorrAmplitudeBins_.at(myBin) ) / ( EBtimeCorrAmplitudeBins_.at(myBin+1) - EBtimeCorrAmplitudeBins_.at(myBin) );
theCorrection += EBtimeCorrShiftBins_.at(myBin);
}
else
{
edm::LogError("EcalRecHitError") << "Assigning time correction impossible. Setting it to 0 ";
theCorrection = 0.;
}

// convert ns into clocks
return theCorrection/25.;
}
if (amplitudeBins.size() != shiftBins.size()) {
edm::LogError("EcalRecHitError")
<< "Size of timeCorrAmplitudeBins different from "
"timeCorrShiftBins. Forcing no time bias corrections. ";

return 0;
}

double EcalUncalibRecHitWorkerGlobal::timeCorrectionEE(float ampliEE){
// computed initially in ns. Than turned in the BX's, as EcalUncalibratedRecHit need be.
double theCorrection=0;

int myBin = -1;
for (int bin=0; bin<(int)EEtimeCorrAmplitudeBins_.size(); bin++ ){
if(ampliEE > EEtimeCorrAmplitudeBins_.at(bin)) {
myBin = bin; }
else break;
for (int bin = 0; bin < (int) amplitudeBins.size(); bin++) {
if (ampli > amplitudeBins.at(bin)) {
myBin = bin;
} else {
break;
}
}

if (myBin == -1) {
theCorrection = shiftBins.at(0);
} else if (myBin == ((int)(amplitudeBins.size() - 1))) {
theCorrection = shiftBins.at(myBin);
} else if (-1 < myBin && myBin < ((int) amplitudeBins.size() - 1)) {
// interpolate linearly between two assingned points
theCorrection = (shiftBins.at(myBin + 1) - shiftBins.at(myBin));
theCorrection *= (((double) ampli) - amplitudeBins.at(myBin)) /
(amplitudeBins.at(myBin + 1) - amplitudeBins.at(myBin));
theCorrection += shiftBins.at(myBin);
} else {
edm::LogError("EcalRecHitError")
<< "Assigning time correction impossible. Setting it to 0 ";
theCorrection = 0.;
}

if (myBin == -1)
{
theCorrection = EEtimeCorrShiftBins_.at(0);
}
else if ( myBin == ((int)(EEtimeCorrAmplitudeBins_.size()-1)) )
{
theCorrection = EEtimeCorrShiftBins_.at( myBin );
}
else if ( -1 < myBin && myBin < ((int)EEtimeCorrAmplitudeBins_.size()-1) )
{
// interpolate linearly between two assingned points
theCorrection = ( EEtimeCorrShiftBins_.at(myBin+1) - EEtimeCorrShiftBins_.at(myBin) );
theCorrection *= ( ((double)ampliEE) - EEtimeCorrAmplitudeBins_.at(myBin) ) / ( EEtimeCorrAmplitudeBins_.at(myBin+1) - EEtimeCorrAmplitudeBins_.at(myBin) );
theCorrection += EEtimeCorrShiftBins_.at(myBin);
}
else
{
edm::LogError("EcalRecHitError") << "Assigning time correction impossible. Setting it to 0 ";
theCorrection = 0.;
}


// convert ns into clocks
return theCorrection/25.;
return theCorrection / 25.;
}




bool
EcalUncalibRecHitWorkerGlobal::run( const edm::Event & evt,
const EcalDigiCollection::const_iterator & itdg,
Expand Down Expand Up @@ -393,8 +354,9 @@ EcalUncalibRecHitWorkerGlobal::run( const edm::Event & evt,
ratioMethod_endcap_.computeTime( EEtimeFitParameters_, EEtimeFitLimits_, EEamplitudeFitParameters_ );
ratioMethod_endcap_.computeAmplitude( EEamplitudeFitParameters_);
EcalUncalibRecHitRatioMethodAlgo<EEDataFrame>::CalculatedRecHit crh = ratioMethod_endcap_.getCalculatedRecHit();
double theTimeCorrectionEE=0;
if(doEEtimeCorrection_) theTimeCorrectionEE = timeCorrectionEE( uncalibRecHit.amplitude() );
double theTimeCorrectionEE = timeCorrection(uncalibRecHit.amplitude(),
timeCorrBias_->EETimeCorrAmplitudeBins, timeCorrBias_->EETimeCorrShiftBins);

uncalibRecHit.setJitter( crh.timeMax - 5 + theTimeCorrectionEE);
uncalibRecHit.setJitterError( std::sqrt(pow(crh.timeError,2) + std::pow(EEtimeConstantTerm_,2)/std::pow(clockToNsConstant,2)) );
uncalibRecHit.setOutOfTimeEnergy( crh.amplitudeMax );
Expand Down Expand Up @@ -430,8 +392,9 @@ EcalUncalibRecHitWorkerGlobal::run( const edm::Event & evt,
ratioMethod_barrel_.computeTime( EBtimeFitParameters_, EBtimeFitLimits_, EBamplitudeFitParameters_ );
ratioMethod_barrel_.computeAmplitude( EBamplitudeFitParameters_);
EcalUncalibRecHitRatioMethodAlgo<EBDataFrame>::CalculatedRecHit crh = ratioMethod_barrel_.getCalculatedRecHit();
double theTimeCorrectionEB=0;
if(doEBtimeCorrection_) theTimeCorrectionEB = timeCorrectionEB( uncalibRecHit.amplitude() );

double theTimeCorrectionEB = timeCorrection(uncalibRecHit.amplitude(),
timeCorrBias_->EBTimeCorrAmplitudeBins, timeCorrBias_->EBTimeCorrShiftBins);

uncalibRecHit.setJitter( crh.timeMax - 5 + theTimeCorrectionEB);

Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -22,6 +22,7 @@
#include "CondFormats/EcalObjects/interface/EcalWeightXtalGroups.h"
#include "CondFormats/EcalObjects/interface/EcalTBWeights.h"
#include "CondFormats/EcalObjects/interface/EcalSampleMask.h"
#include "CondFormats/EcalObjects/interface/EcalTimeBiasCorrections.h"
#include "SimCalorimetry/EcalSimAlgos/interface/EBShape.h"
#include "SimCalorimetry/EcalSimAlgos/interface/EEShape.h"

Expand Down Expand Up @@ -53,8 +54,8 @@ class EcalUncalibRecHitWorkerGlobal : public EcalUncalibRecHitWorkerBaseClass {

template < class C > int isSaturated(const C & digi);

double timeCorrectionEB(float ampliEB);
double timeCorrectionEE(float ampliEE);
double timeCorrection(float ampli,
const std::vector<float>& amplitudeBins, const std::vector<float>& shiftBins);

// weights method
edm::ESHandle<EcalWeightXtalGroups> grps;
Expand All @@ -76,12 +77,7 @@ class EcalUncalibRecHitWorkerGlobal : public EcalUncalibRecHitWorkerBaseClass {
std::vector<double> EEamplitudeFitParameters_;
std::pair<double,double> EBtimeFitLimits_;
std::pair<double,double> EEtimeFitLimits_;
bool doEBtimeCorrection_;
bool doEEtimeCorrection_;
std::vector<double> EBtimeCorrAmplitudeBins_;
std::vector<double> EBtimeCorrShiftBins_;
std::vector<double> EEtimeCorrAmplitudeBins_;
std::vector<double> EEtimeCorrShiftBins_;

EcalUncalibRecHitRatioMethodAlgo<EBDataFrame> ratioMethod_barrel_;
EcalUncalibRecHitRatioMethodAlgo<EEDataFrame> ratioMethod_endcap_;

Expand All @@ -101,6 +97,8 @@ class EcalUncalibRecHitWorkerGlobal : public EcalUncalibRecHitWorkerBaseClass {
double amplitudeThreshEE_;
double ebSpikeThresh_;

edm::ESHandle<EcalTimeBiasCorrections> timeCorrBias_;

// leading edge method
edm::ESHandle<EcalTimeCalibConstants> itime;
edm::ESHandle<EcalTimeOffsetConstant> offtime;
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -30,37 +30,6 @@
outOfTimeThresholdGain61mEE = cms.double(10), # times estimated precision
amplitudeThresholdEB = cms.double(10),
amplitudeThresholdEE = cms.double(10),
#amplitude-dependent time corrections; EE and EB have separate corrections
#EXtimeCorrAmplitudes (ADC) and EXtimeCorrShifts (ns) need to have the same number of elements
#Bins must be ordered in amplitude. First-last bins take care of under-overflows.
doEBtimeCorrection = cms.bool(False),
doEEtimeCorrection = cms.bool(False),

EBtimeCorrAmplitudeBins = cms.vdouble(
7.9, 8.9, 10, 11.2, 12.5, 14.1, 15.8, 17.7, 19.9, 22.3, 25, 28.1, 31.5, 35.3, 39.7,
44.5, 49.9, 56, 62.8, 70.5, 79.1, 88.8, 99.6, 111.7, 125.4, 140.7, 157.9, 177.1, 198.7, 223,
250.2, 280.7, 315, 353.4, 396.5, 444.9, 499.2, 560.1, 628.4, 705.1, 791.1, 887.7, 996, 1117.5, 1253.9,
1406.8, 1578.5, 1771.1, 1987.2, 2229.7, 2501.8, 2807, 3149.5, 3533.8, 3895.9, 3896, 4311.8, 4837.9, 5428.2, 6090.6,
6833.7, 7667.5, 8603.1, 9652.9, 10830, 12152, 13635, 15298, 17165, 19260, 21610),

EBtimeCorrShiftBins = cms.vdouble(
-1.770, -1.770, -1.770, -1.770, -1.666, -1.430, -1.233, -1.012, -0.866, -0.736, -0.640, -0.561, -0.505, -0.452, -0.405,
-0.363, -0.335, -0.305, -0.279, -0.260, -0.239, -0.220, -0.204, -0.191, -0.186, -0.177, -0.158, -0.137, -0.126, -0.115,
-0.104, -0.096, -0.085, -0.064, -0.056, -0.036, -0.020, -0.006, -0.020, -0.009, -0.020, 0.005, 0.053, 0.076, 0.093,
0.137, 0.143, 0.171, 0.222, 0.229, 0.271, 0.298, 0.312, 0.307, 0.254 , -0.997 ,-0.859 , -0.819, -0.775, -0.589,
-0.428, -0.288, -0.434, -0.277, -0.210, -0.179, -0.134, 0.362, 0.152, -0.282, -0.382),

EEtimeCorrAmplitudeBins = cms.vdouble(
15.7, 17.6, 19.7, 22.1, 24.8, 27.9, 31.3, 35.1, 39.4, 44.2, 49.6, 55.6, 62.4, 70, 78.6,
88.1, 98.9, 111, 124.5, 139.7, 156.7, 175.9, 197.3, 221.4, 248.4, 278.7, 312.7, 350.9, 393.7, 441.7,
495.6, 556.1, 624, 700.1, 785.5, 881.4, 988.9, 1109.6, 1245, 1396.9, 1567.3, 1758.6, 1973.1, 2213.9, 2484,
2787.1, 3127.2, 3508.8, 3936.9, 4417.3, 4956.3, 5561.1, 6239.6, 7001, 7522.8, 8440.7, 9470.6, 10626),

EEtimeCorrShiftBins = cms.vdouble(
-0.896, -0.896, -0.896, -0.896, -0.563, -0.392, -0.287, -0.203, -0.135, -0.100, -0.068, -0.050, -0.060, -0.052, -0.055,
-0.050, -0.052, -0.056, -0.055, -0.056, -0.048, -0.037, -0.038, -0.037, -0.025, -0.026, -0.024, -0.013, -0.003, 0.005,
0.020, 0.026, 0.008, 0.007, -0.006, 0.024, 0.045, 0.062, 0.085, 0.088 , 0.111 , 0.139, 0.156, 0.176, 0.210,
0.242, 0.267, 0.301, 0.318, 0.278, 0.287, 0.218, 0.305, 0.245, 0.184, -0.159, -0.095, 0.037),

ebSpikeThreshold = cms.double(1.042),

Expand Down