EgammaTools/plugins/CalibratedElectronProducers.cc

0001 //author: Alan Smithee
0002 //description:
0003 //  this class allows the residual scale and smearing to be applied to electrons
0004 //  the scale and smearing is on the ecal part of the energy
0005 //  hence the E/p combination needs to be re-don, hence the E/p Combination Tools
0006 //  it re-applies the regression with the new corrected ecal energy
0007 //  returns a vector of calibrated energies and correction data, indexed by EGEnergySysIndex
0008 //  a port of EgammaAnalysis/ElectronTools/CalibratedElectronProducerRun2
0009
0010 #include "DataFormats/Common/interface/Handle.h"
0011 #include "FWCore/Framework/interface/stream/EDProducer.h"
0012 #include "FWCore/Framework/interface/Event.h"
0013 #include "FWCore/Framework/interface/EventSetup.h"
0014 #include "FWCore/Framework/interface/ESHandle.h"
0015 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0016 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0017 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0018 #include "FWCore/Utilities/interface/EDGetToken.h"
0019 #include "FWCore/Framework/interface/ConsumesCollector.h"
0020
0021 #include "DataFormats/EgammaCandidates/interface/GsfElectron.h"
0022 #include "DataFormats/EgammaCandidates/interface/GsfElectronFwd.h"
0023 #include "DataFormats/PatCandidates/interface/Electron.h"
0024 #include "DataFormats/Common/interface/ValueMap.h"
0025
0026 #include "CondFormats/DataRecord/interface/GBRDWrapperRcd.h"
0027 #include "CondFormats/GBRForest/interface/GBRForestD.h"
0028 #include "RecoEgamma/EgammaTools/interface/EpCombinationTool.h"
0029 #include "RecoEgamma/EgammaTools/interface/ElectronEnergyCalibrator.h"
0030 #include "RecoEgamma/EgammaTools/interface/EGEnergySysIndex.h"
0031 #include "RecoEgamma/EgammaTools/interface/EgammaRandomSeeds.h"
0032
0033 #include "TRandom2.h"
0034
0035 #include <memory>
0036
0037 #include <vector>
0038
0039 template <typename T>
0040 class CalibratedElectronProducerT : public edm::stream::EDProducer<> {
0041 public:
0042   explicit CalibratedElectronProducerT(const edm::ParameterSet&);
0043   ~CalibratedElectronProducerT() override {}
0044   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0045   void produce(edm::Event&, const edm::EventSetup&) override;
0046
0047 private:
0048   void setSemiDetRandomSeed(const edm::Event& iEvent, const T& obj, size_t nrObjs, size_t objNr);
0049
0050   edm::EDGetTokenT<edm::View<T>> electronToken_;
0051
0052   EpCombinationTool epCombinationTool_;
0053   ElectronEnergyCalibrator energyCorrector_;
0054   std::unique_ptr<TRandom> semiDeterministicRng_;
0055   edm::EDGetTokenT<EcalRecHitCollection> recHitCollectionEBToken_;
0056   edm::EDGetTokenT<EcalRecHitCollection> recHitCollectionEEToken_;
0057   bool produceCalibratedObjs_;
0058   static const std::vector<int> valMapsToStore_;
0059 };
0060
0061 template <typename T>
0062 const std::vector<int> CalibratedElectronProducerT<T>::valMapsToStore_ = {
0063     EGEnergySysIndex::kScaleStatUp,       EGEnergySysIndex::kScaleStatDown,     EGEnergySysIndex::kScaleSystUp,
0064     EGEnergySysIndex::kScaleSystDown,     EGEnergySysIndex::kScaleGainUp,       EGEnergySysIndex::kScaleGainDown,
0065     EGEnergySysIndex::kSmearRhoUp,        EGEnergySysIndex::kSmearRhoDown,      EGEnergySysIndex::kSmearPhiUp,
0066     EGEnergySysIndex::kSmearPhiDown,      EGEnergySysIndex::kScaleUp,           EGEnergySysIndex::kScaleDown,
0067     EGEnergySysIndex::kSmearUp,           EGEnergySysIndex::kSmearDown,         EGEnergySysIndex::kScaleValue,
0068     EGEnergySysIndex::kSmearValue,        EGEnergySysIndex::kSmearNrSigma,      EGEnergySysIndex::kEcalPreCorr,
0069     EGEnergySysIndex::kEcalErrPreCorr,    EGEnergySysIndex::kEcalPostCorr,      EGEnergySysIndex::kEcalErrPostCorr,
0070     EGEnergySysIndex::kEcalTrkPreCorr,    EGEnergySysIndex::kEcalTrkErrPreCorr, EGEnergySysIndex::kEcalTrkPostCorr,
0071     EGEnergySysIndex::kEcalTrkErrPostCorr};
0072
0073 namespace {
0074   template <typename HandleType, typename ValType>
0075   void fillAndStoreValueMap(edm::Event& iEvent,
0076                             HandleType objHandle,
0077                             const std::vector<ValType>& vals,
0078                             const std::string& name) {
0079     auto valMap = std::make_unique<edm::ValueMap<ValType>>();
0080     typename edm::ValueMap<ValType>::Filler filler(*valMap);
0081     filler.insert(objHandle, vals.begin(), vals.end());
0082     filler.fill();
0083     iEvent.put(std::move(valMap), name);
0084   }
0085 }  // namespace
0086
0087 template <typename T>
0088 CalibratedElectronProducerT<T>::CalibratedElectronProducerT(const edm::ParameterSet& conf)
0089     : electronToken_(consumes(conf.getParameter<edm::InputTag>("src"))),
0090       epCombinationTool_{conf.getParameterSet("epCombConfig"), consumesCollector()},
0091       energyCorrector_(epCombinationTool_, conf.getParameter<std::string>("correctionFile")),
0092       recHitCollectionEBToken_(consumes(conf.getParameter<edm::InputTag>("recHitCollectionEB"))),
0093       recHitCollectionEEToken_(consumes(conf.getParameter<edm::InputTag>("recHitCollectionEE"))),
0094       produceCalibratedObjs_(conf.getParameter<bool>("produceCalibratedObjs")) {
0095   energyCorrector_.setMinEt(conf.getParameter<double>("minEtToCalibrate"));
0096
0097   if (conf.getParameter<bool>("semiDeterministic")) {
0098     semiDeterministicRng_ = std::make_unique<TRandom2>();
0099     energyCorrector_.initPrivateRng(semiDeterministicRng_.get());
0100   }
0101
0102   if (produceCalibratedObjs_)
0103     produces<std::vector<T>>();
0104
0105   for (const auto& toStore : valMapsToStore_) {
0106     produces<edm::ValueMap<float>>(EGEnergySysIndex::name(toStore));
0107   }
0108 }
0109
0110 template <typename T>
0111 void CalibratedElectronProducerT<T>::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0112   edm::ParameterSetDescription desc;
0113   desc.add<edm::InputTag>("src", edm::InputTag("gedPhotons"));
0114   desc.add<edm::ParameterSetDescription>("epCombConfig", EpCombinationTool::makePSetDescription());
0115   desc.add<edm::InputTag>("recHitCollectionEB", edm::InputTag("reducedEcalRecHitsEB"));
0116   desc.add<edm::InputTag>("recHitCollectionEE", edm::InputTag("reducedEcalRecHitsEE"));
0117   desc.add<std::string>("correctionFile", std::string());
0118   desc.add<double>("minEtToCalibrate", 5.0);
0119   desc.add<bool>("produceCalibratedObjs", true);
0120   desc.add<bool>("semiDeterministic", true);
0121   std::vector<std::string> valMapsProduced;
0122   valMapsProduced.reserve(valMapsToStore_.size());
0123   for (auto varToStore : valMapsToStore_)
0124     valMapsProduced.push_back(EGEnergySysIndex::name(varToStore));
0125   desc.add<std::vector<std::string>>("valueMapsStored", valMapsProduced)
0126       ->setComment(
0127           "provides to python configs the list of valuemaps stored, can not be overriden in the python config");
0128   descriptions.addWithDefaultLabel(desc);
0129 }
0130
0131 template <typename T>
0132 void CalibratedElectronProducerT<T>::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0133   epCombinationTool_.setEventContent(iSetup);
0134
0135   auto inHandle = iEvent.getHandle(electronToken_);
0136
0137   auto recHitCollectionEBHandle = iEvent.getHandle(recHitCollectionEBToken_);
0138   auto recHitCollectionEEHandle = iEvent.getHandle(recHitCollectionEEToken_);
0139
0140   std::unique_ptr<std::vector<T>> out = std::make_unique<std::vector<T>>();
0141
0142   size_t nrObj = inHandle->size();
0143   std::array<std::vector<float>, EGEnergySysIndex::kNrSysErrs> results;
0144   for (auto& res : results)
0145     res.reserve(nrObj);
0146
0147   const ElectronEnergyCalibrator::EventType evtType =
0148       iEvent.isRealData() ? ElectronEnergyCalibrator::EventType::DATA : ElectronEnergyCalibrator::EventType::MC;
0149
0150   for (const auto& ele : *inHandle) {
0151     out->push_back(ele);
0152
0153     if (semiDeterministicRng_)
0154       setSemiDetRandomSeed(iEvent, ele, nrObj, out->size());
0155
0156     const EcalRecHitCollection* recHits =
0157         (ele.isEB()) ? recHitCollectionEBHandle.product() : recHitCollectionEEHandle.product();
0158     std::array<float, EGEnergySysIndex::kNrSysErrs> uncertainties =
0159         energyCorrector_.calibrate(out->back(), iEvent.id().run(), recHits, iEvent.streamID(), evtType);
0160
0161     for (size_t index = 0; index < EGEnergySysIndex::kNrSysErrs; index++) {
0162       results[index].push_back(uncertainties[index]);
0163     }
0164   }
0165
0166   auto fillAndStore = [&](auto handle) {
0167     for (const auto& mapToStore : valMapsToStore_) {
0168       fillAndStoreValueMap(iEvent, handle, results[mapToStore], EGEnergySysIndex::name(mapToStore));
0169     }
0170   };
0171
0172   if (produceCalibratedObjs_) {
0173     fillAndStore(iEvent.put(std::move(out)));
0174   } else {
0175     fillAndStore(inHandle);
0176   }
0177 }
0178
0179 template <typename T>
0180 void CalibratedElectronProducerT<T>::setSemiDetRandomSeed(const edm::Event& iEvent,
0181                                                           const T& obj,
0182                                                           size_t nrObjs,
0183                                                           size_t objNr) {
0184   if (obj.superCluster().isNonnull()) {
0185     semiDeterministicRng_->SetSeed(egamma::getRandomSeedFromSC(iEvent, obj.superCluster()));
0186   } else {
0187     semiDeterministicRng_->SetSeed(egamma::getRandomSeedFromObj(iEvent, obj, nrObjs, objNr));
0188   }
0189 }
0190
0191 using CalibratedElectronProducer = CalibratedElectronProducerT<reco::GsfElectron>;
0192 using CalibratedPatElectronProducer = CalibratedElectronProducerT<pat::Electron>;
0193
0194 #include "FWCore/Framework/interface/MakerMacros.h"
0195
0196 DEFINE_FWK_MODULE(CalibratedElectronProducer);
0197 DEFINE_FWK_MODULE(CalibratedPatElectronProducer);