File indexing completed on 2025-01-31 02:19:45
0001 #ifndef CorrectedJetProducer_h
0002 #define CorrectedJetProducer_h
0003
0004 #include <sstream>
0005 #include <string>
0006 #include <vector>
0007
0008 #include "CommonTools/Utils/interface/PtComparator.h"
0009 #include "DataFormats/Common/interface/Handle.h"
0010 #include "DataFormats/Common/interface/Ref.h"
0011 #include "DataFormats/Common/interface/RefToBase.h"
0012 #include "FWCore/Framework/interface/Event.h"
0013 #include "FWCore/Framework/interface/EventSetup.h"
0014 #include "FWCore/Framework/interface/global/EDProducer.h"
0015 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0016 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0017 #include "FWCore/Utilities/interface/InputTag.h"
0018 #include "FWCore/Utilities/interface/transform.h"
0019 #include "JetMETCorrections/JetCorrector/interface/JetCorrector.h"
0020
0021 namespace edm {
0022 class ParameterSet;
0023 }
0024
0025 namespace reco {
0026 template <class T>
0027 class CorrectedJetProducer : public edm::global::EDProducer<> {
0028 public:
0029 typedef std::vector<T> JetCollection;
0030 explicit CorrectedJetProducer(const edm::ParameterSet& fParameters);
0031 ~CorrectedJetProducer() override {}
0032 void produce(edm::StreamID, edm::Event&, const edm::EventSetup&) const override;
0033 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0034
0035 private:
0036 const edm::EDGetTokenT<JetCollection> mInput;
0037 const std::vector<edm::EDGetTokenT<reco::JetCorrector> > mCorrectorTokens;
0038 const bool mVerbose;
0039 };
0040 }
0041
0042
0043
0044 namespace reco {
0045
0046 template <class T>
0047 CorrectedJetProducer<T>::CorrectedJetProducer(const edm::ParameterSet& fConfig)
0048 : mInput(consumes<JetCollection>(fConfig.getParameter<edm::InputTag>("src"))),
0049 mCorrectorTokens(
0050 edm::vector_transform(fConfig.getParameter<std::vector<edm::InputTag> >("correctors"),
0051 [this](edm::InputTag const& tag) { return consumes<reco::JetCorrector>(tag); })),
0052 mVerbose(fConfig.getUntrackedParameter<bool>("verbose", false)) {
0053 std::string alias = fConfig.getUntrackedParameter<std::string>("alias", "");
0054 if (alias.empty())
0055 produces<JetCollection>();
0056 else
0057 produces<JetCollection>().setBranchAlias(alias);
0058 }
0059
0060 template <class T>
0061 void CorrectedJetProducer<T>::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0062 edm::ParameterSetDescription desc;
0063 desc.add<edm::InputTag>("src", edm::InputTag(""));
0064 desc.add<std::vector<edm::InputTag> >("correctors", {});
0065 desc.addUntracked<bool>("verbose", false);
0066 desc.addUntracked<std::string>("alias", "");
0067 descriptions.addWithDefaultLabel(desc);
0068 }
0069
0070 template <class T>
0071 void CorrectedJetProducer<T>::produce(edm::StreamID, edm::Event& fEvent, const edm::EventSetup& fSetup) const {
0072
0073 std::vector<reco::JetCorrector const*> correctors(mCorrectorTokens.size(), nullptr);
0074
0075
0076 for (unsigned i = 0; i < mCorrectorTokens.size(); i++) {
0077 edm::Handle<reco::JetCorrector> handle;
0078 fEvent.getByToken(mCorrectorTokens[i], handle);
0079 correctors[i] = handle.product();
0080 }
0081 edm::Handle<JetCollection> jets;
0082 fEvent.getByToken(mInput, jets);
0083 std::unique_ptr<JetCollection> result(new JetCollection);
0084 typename JetCollection::const_iterator jet;
0085 for (jet = jets->begin(); jet != jets->end(); jet++) {
0086 const T* referenceJet = &*jet;
0087 int index = jet - jets->begin();
0088 edm::RefToBase<reco::Jet> jetRef(edm::Ref<JetCollection>(jets, index));
0089 T correctedJet = *jet;
0090 if (mVerbose)
0091 std::cout << "CorrectedJetProducer::produce-> original jet: " << jet->print() << std::endl;
0092 for (unsigned i = 0; i < mCorrectorTokens.size(); ++i) {
0093 if (!(correctors[i]->vectorialCorrection())) {
0094
0095 double scale = 1.;
0096 if (!(correctors[i]->refRequired()))
0097 scale = correctors[i]->correction(*referenceJet);
0098 else
0099 scale = correctors[i]->correction(*referenceJet, jetRef);
0100 if (mVerbose)
0101 std::cout << "CorrectedJetProducer::produce-> Corrector # " << i << ", correction factor: " << scale
0102 << std::endl;
0103 correctedJet.scaleEnergy(scale);
0104 referenceJet = &correctedJet;
0105 } else {
0106
0107 reco::JetCorrector::LorentzVector corrected;
0108 double scale = correctors[i]->correction(*referenceJet, jetRef, corrected);
0109 if (mVerbose)
0110 std::cout << "CorrectedJetProducer::produce-> Corrector # " << i << ", correction factor: " << scale
0111 << std::endl;
0112 correctedJet.setP4(corrected);
0113 referenceJet = &correctedJet;
0114 }
0115 }
0116 if (mVerbose)
0117 std::cout << "CorrectedJetProducer::produce-> corrected jet: " << correctedJet.print() << std::endl;
0118 result->push_back(correctedJet);
0119 }
0120 NumericSafeGreaterByPt<T> compJets;
0121
0122 std::sort(result->begin(), result->end(), compJets);
0123
0124 fEvent.put(std::move(result));
0125 }
0126
0127 }
0128
0129 #endif