ElectronTools/src/ElectronEnergyCalibratorRun2.cc

0001 #include "EgammaAnalysis/ElectronTools/interface/ElectronEnergyCalibratorRun2.h"
0002 #include <CLHEP/Random/RandGaussQ.h>
0003 #include "FWCore/ServiceRegistry/interface/Service.h"
0004 #include "FWCore/Utilities/interface/RandomNumberGenerator.h"
0005 #include "FWCore/Utilities/interface/Exception.h"
0006
0007 ElectronEnergyCalibratorRun2::ElectronEnergyCalibratorRun2(EpCombinationTool &combinator,
0008                                                            bool isMC,
0009                                                            bool synchronization,
0010                                                            std::string correctionFile)
0011     : epCombinationTool_(&combinator),
0012       isMC_(isMC),
0013       synchronization_(synchronization),
0014       rng_(nullptr),
0015       _correctionRetriever(correctionFile)  // here is opening the files and reading the corrections
0016 {
0017   if (isMC_) {
0018     _correctionRetriever.doScale = false;
0019     _correctionRetriever.doSmearings = true;
0020   } else {
0021     _correctionRetriever.doScale = true;
0022     _correctionRetriever.doSmearings = false;
0023   }
0024 }
0025
0026 ElectronEnergyCalibratorRun2::~ElectronEnergyCalibratorRun2() {}
0027
0028 void ElectronEnergyCalibratorRun2::initPrivateRng(TRandom *rnd) { rng_ = rnd; }
0029
0030 void ElectronEnergyCalibratorRun2::calibrate(reco::GsfElectron &electron,
0031                                              unsigned int runNumber,
0032                                              edm::StreamID const &id) const {
0033   SimpleElectron simple(electron, runNumber, isMC_);
0034   calibrate(simple, id);
0035   simple.writeTo(electron);
0036 }
0037 void ElectronEnergyCalibratorRun2::calibrate(SimpleElectron &electron, edm::StreamID const &id) const {
0038   assert(isMC_ == electron.isMC());
0039   float smear = 0.0, scale = 1.0;
0040   float aeta = std::abs(electron.getEta());  //, r9 = electron.getR9();
0041   float et = electron.getNewEnergy() / cosh(aeta);
0042
0043   scale = _correctionRetriever.ScaleCorrection(electron.getRunNumber(), electron.isEB(), electron.getR9(), aeta, et);
0044   smear = _correctionRetriever.getSmearingSigma(
0045       electron.getRunNumber(), electron.isEB(), electron.getR9(), aeta, et, 0., 0.);
0046
0047   double newEcalEnergy, newEcalEnergyError;
0048   if (isMC_) {
0049     double corr = 1.0 + smear * gauss(id);
0050     newEcalEnergy = electron.getNewEnergy() * corr;
0051     newEcalEnergyError = std::hypot(electron.getNewEnergyError() * corr, smear * newEcalEnergy);
0052   } else {
0053     newEcalEnergy = electron.getNewEnergy() * scale;
0054     newEcalEnergyError = std::hypot(electron.getNewEnergyError() * scale, smear * newEcalEnergy);
0055   }
0056   electron.setNewEnergy(newEcalEnergy);
0057   electron.setNewEnergyError(newEcalEnergyError);
0058   epCombinationTool_->combine(electron);
0059 }
0060
0061 double ElectronEnergyCalibratorRun2::gauss(edm::StreamID const &id) const {
0062   if (synchronization_)
0063     return 1.0;
0064   if (rng_) {
0065     return rng_->Gaus();
0066   } else {
0067     edm::Service<edm::RandomNumberGenerator> rng;
0068     if (!rng.isAvailable()) {
0069       throw cms::Exception("Configuration")
0070           << "XXXXXXX requires the RandomNumberGeneratorService\n"
0071              "which is not present in the configuration file.  You must add the service\n"
0072              "in the configuration file or remove the modules that require it.";
0073     }
0074     CLHEP::RandGaussQ gaussDistribution(rng->getEngine(id), 0.0, 1.0);
0075     return gaussDistribution.fire();
0076   }
0077 }