File indexing completed on 2024-04-06 11:59:02
0001
0002 #include <memory>
0003 #include <string>
0004
0005
0006 #include "DataFormats/Common/interface/Ref.h"
0007 #include "DataFormats/CaloTowers/interface/CaloTowerCollection.h"
0008 #include "DataFormats/CaloTowers/interface/CaloTowerDetId.h"
0009 #include "DataFormats/EcalRecHit/interface/EcalRecHitCollections.h"
0010 #include "DataFormats/HcalRecHit/interface/HcalRecHitCollections.h"
0011 #include "DataFormats/DetId/interface/DetId.h"
0012 #include "DataFormats/GeometryVector/interface/GlobalPoint.h"
0013
0014 #include "CondFormats/HcalObjects/interface/HcalRespCorrs.h"
0015 #include "CondFormats/DataRecord/interface/HcalRespCorrsRcd.h"
0016
0017 #include "FWCore/Framework/interface/Frameworkfwd.h"
0018 #include "FWCore/Framework/interface/one/EDProducer.h"
0019 #include "FWCore/Framework/interface/Event.h"
0020 #include "FWCore/Framework/interface/MakerMacros.h"
0021 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0022 #include "FWCore/Framework/interface/EventSetup.h"
0023 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0024 #include "RecoTracker/TrackProducer/interface/TrackProducerBase.h"
0025 #include "TrackingTools/TransientTrack/interface/TransientTrack.h"
0026 #include "DataFormats/JetReco/interface/CaloJetCollection.h"
0027 #include "FWCore/Utilities/interface/Exception.h"
0028
0029
0030
0031
0032
0033 namespace cms {
0034
0035 class HitReCalibrator : public edm::one::EDProducer<> {
0036 public:
0037 explicit HitReCalibrator(const edm::ParameterSet&);
0038 ~HitReCalibrator() override;
0039
0040 void beginJob() override;
0041
0042 void produce(edm::Event&, const edm::EventSetup&) override;
0043
0044 private:
0045
0046
0047 bool allowMissingInputs_;
0048
0049 const edm::EDGetTokenT<HBHERecHitCollection> tok_hbhe_;
0050 const edm::EDGetTokenT<HORecHitCollection> tok_ho_;
0051 const edm::EDGetTokenT<HFRecHitCollection> tok_hf_;
0052
0053 const edm::ESGetToken<HcalRespCorrs, HcalRespCorrsRcd> tok_resp_;
0054 };
0055 }
0056
0057 namespace cms {
0058
0059 HitReCalibrator::HitReCalibrator(const edm::ParameterSet& iConfig)
0060 : tok_hbhe_(consumes<HBHERecHitCollection>(iConfig.getParameter<edm::InputTag>("hbheInput"))),
0061 tok_ho_(consumes<HORecHitCollection>(iConfig.getParameter<edm::InputTag>("hoInput"))),
0062 tok_hf_(consumes<HFRecHitCollection>(iConfig.getParameter<edm::InputTag>("hfInput"))),
0063 tok_resp_(esConsumes<HcalRespCorrs, HcalRespCorrsRcd>()) {
0064 allowMissingInputs_ = true;
0065
0066
0067
0068 produces<HBHERecHitCollection>("DiJetsHBHEReRecHitCollection");
0069 produces<HORecHitCollection>("DiJetsHOReRecHitCollection");
0070 produces<HFRecHitCollection>("DiJetsHFReRecHitCollection");
0071 }
0072 void HitReCalibrator::beginJob() {}
0073
0074 HitReCalibrator::~HitReCalibrator() {}
0075
0076
0077 void HitReCalibrator::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0078 auto miniDiJetsHBHERecHitCollection = std::make_unique<HBHERecHitCollection>();
0079 auto miniDiJetsHORecHitCollection = std::make_unique<HORecHitCollection>();
0080 auto miniDiJetsHFRecHitCollection = std::make_unique<HFRecHitCollection>();
0081
0082 const HcalRespCorrs* jetRecalib = &iSetup.getData(tok_resp_);
0083
0084 try {
0085 const edm::Handle<HBHERecHitCollection>& hbhe = iEvent.getHandle(tok_hbhe_);
0086 const HBHERecHitCollection Hithbhe = *(hbhe.product());
0087 for (HBHERecHitCollection::const_iterator hbheItr = Hithbhe.begin(); hbheItr != Hithbhe.end(); hbheItr++) {
0088 DetId id = hbheItr->detid();
0089 float recal;
0090 if (jetRecalib->exists(id))
0091 recal = jetRecalib->getValues(id)->getValue();
0092 else
0093 recal = 1.;
0094 float energy = hbheItr->energy();
0095 float time = hbheItr->time();
0096 HBHERecHit* hbhehit = new HBHERecHit(id, recal * energy, time);
0097 miniDiJetsHBHERecHitCollection->push_back(*hbhehit);
0098 }
0099 } catch (cms::Exception& e) {
0100 if (!allowMissingInputs_) {
0101 edm::LogError("HitCalib") << "No HBHE collection ";
0102 throw e;
0103 }
0104 }
0105
0106 try {
0107 const edm::Handle<HORecHitCollection>& ho = iEvent.getHandle(tok_ho_);
0108 const HORecHitCollection Hitho = *(ho.product());
0109 for (HORecHitCollection::const_iterator hoItr = Hitho.begin(); hoItr != Hitho.end(); hoItr++) {
0110 DetId id = hoItr->detid();
0111 float recal;
0112 if (jetRecalib->exists(id))
0113 recal = jetRecalib->getValues(id)->getValue();
0114 else
0115 recal = 1.;
0116 float energy = hoItr->energy();
0117 float time = hoItr->time();
0118 HORecHit* hohit = new HORecHit(id, recal * energy, time);
0119 miniDiJetsHORecHitCollection->push_back(*hohit);
0120 }
0121 } catch (cms::Exception& e) {
0122 if (!allowMissingInputs_) {
0123 edm::LogError("HitCalib") << " No HO collection ";
0124 throw e;
0125 }
0126 }
0127
0128 try {
0129 const edm::Handle<HFRecHitCollection>& hf = iEvent.getHandle(tok_hf_);
0130 const HFRecHitCollection Hithf = *(hf.product());
0131 for (HFRecHitCollection::const_iterator hfItr = Hithf.begin(); hfItr != Hithf.end(); hfItr++) {
0132 DetId id = hfItr->detid();
0133 float recal;
0134 if (jetRecalib->exists(id))
0135 recal = jetRecalib->getValues(id)->getValue();
0136 else
0137 recal = 1.;
0138 float energy = hfItr->energy();
0139 float time = hfItr->time();
0140 HFRecHit* hfhit = new HFRecHit(id, recal * energy, time);
0141 miniDiJetsHFRecHitCollection->push_back(*hfhit);
0142 }
0143 } catch (cms::Exception& e) {
0144 if (!allowMissingInputs_)
0145 throw e;
0146 }
0147
0148
0149
0150 iEvent.put(std::move(miniDiJetsHBHERecHitCollection), "DiJetsHBHEReRecHitCollection");
0151
0152 iEvent.put(std::move(miniDiJetsHORecHitCollection), "DiJetsHOReRecHitCollection");
0153
0154 iEvent.put(std::move(miniDiJetsHFRecHitCollection), "DiJetsHFReRecHitCollection");
0155 }
0156 }
0157
0158 #include "FWCore/Framework/interface/MakerMacros.h"
0159
0160 using cms::HitReCalibrator;
0161 DEFINE_FWK_MODULE(HitReCalibrator);