File indexing completed on 2023-03-17 10:59:16
0001 #include "EgammaAnalysis/ElectronTools/interface/PhotonEnergyCalibratorRun2.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 PhotonEnergyCalibratorRun2::PhotonEnergyCalibratorRun2(bool isMC, bool synchronization, std::string correctionFile)
0008 : isMC_(isMC),
0009 synchronization_(synchronization),
0010 rng_(nullptr),
0011 _correctionRetriever(correctionFile)
0012 {
0013 if (isMC_) {
0014 _correctionRetriever.doScale = false;
0015 _correctionRetriever.doSmearings = true;
0016 } else {
0017 _correctionRetriever.doScale = true;
0018 _correctionRetriever.doSmearings = false;
0019 }
0020 }
0021
0022 PhotonEnergyCalibratorRun2::~PhotonEnergyCalibratorRun2() {}
0023
0024 void PhotonEnergyCalibratorRun2::initPrivateRng(TRandom *rnd) { rng_ = rnd; }
0025
0026 void PhotonEnergyCalibratorRun2::calibrate(reco::Photon &photon,
0027 unsigned int runNumber,
0028 edm::StreamID const &id) const {
0029 SimplePhoton simple(photon, runNumber, isMC_);
0030 calibrate(simple, id);
0031 simple.writeTo(photon);
0032 }
0033
0034 void PhotonEnergyCalibratorRun2::calibrate(SimplePhoton &photon, edm::StreamID const &id) const {
0035 assert(isMC_ == photon.isMC());
0036 float smear = 0.0, scale = 1.0;
0037 float aeta = std::abs(photon.getEta());
0038 float et = photon.getNewEnergy() / cosh(aeta);
0039
0040 scale = _correctionRetriever.ScaleCorrection(photon.getRunNumber(), photon.isEB(), photon.getR9(), aeta, et);
0041 smear = _correctionRetriever.getSmearingSigma(photon.getRunNumber(), photon.isEB(), photon.getR9(), aeta, et, 0., 0.);
0042
0043 double newEcalEnergy, newEcalEnergyError;
0044 if (isMC_) {
0045 double corr = 1.0 + smear * gauss(id);
0046 newEcalEnergy = photon.getNewEnergy() * corr;
0047 newEcalEnergyError = std::hypot(photon.getNewEnergyError() * corr, smear * newEcalEnergy);
0048 } else {
0049 newEcalEnergy = photon.getNewEnergy() * scale;
0050 newEcalEnergyError = std::hypot(photon.getNewEnergyError() * scale, smear * newEcalEnergy);
0051 }
0052 photon.setNewEnergy(newEcalEnergy);
0053 photon.setNewEnergyError(newEcalEnergyError);
0054 }
0055
0056 double PhotonEnergyCalibratorRun2::gauss(edm::StreamID const &id) const {
0057 if (synchronization_)
0058 return 1.0;
0059 if (rng_) {
0060 return rng_->Gaus();
0061 } else {
0062 edm::Service<edm::RandomNumberGenerator> rng;
0063 if (!rng.isAvailable()) {
0064 throw cms::Exception("Configuration")
0065 << "XXXXXXX requires the RandomNumberGeneratorService\n"
0066 "which is not present in the configuration file. You must add the service\n"
0067 "in the configuration file or remove the modules that require it.";
0068 }
0069 CLHEP::RandGaussQ gaussDistribution(rng->getEngine(id), 0.0, 1.0);
0070 return gaussDistribution.fire();
0071 }
0072 }