File indexing completed on 2024-04-06 12:19:26
0001
0002
0003
0004 #include "FWCore/Framework/interface/Frameworkfwd.h"
0005 #include "FWCore/Framework/interface/stream/EDProducer.h"
0006 #include "FWCore/Framework/interface/Event.h"
0007 #include "FWCore/Framework/interface/EventSetup.h"
0008 #include "FWCore/Framework/interface/MakerMacros.h"
0009 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0010 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0011 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0012 #include "FWCore/Utilities/interface/InputTag.h"
0013
0014 #include "DataFormats/Candidate/interface/Candidate.h"
0015 #include "DataFormats/METReco/interface/PFMET.h"
0016 #include "DataFormats/METReco/interface/CorrMETData.h"
0017
0018 #include "JetMETCorrections/Type1MET/interface/AddCorrectionsToGenericMET.h"
0019
0020 #include "HLTrigger/HLTcore/interface/defaultModuleLabel.h"
0021
0022 #include <vector>
0023
0024
0025 class CorrectedPFMETProducer : public edm::stream::EDProducer<> {
0026 public:
0027 explicit CorrectedPFMETProducer(const edm::ParameterSet& cfg)
0028 : corrector(), token_(consumes<METCollection>(cfg.getParameter<edm::InputTag>("src"))) {
0029 std::vector<edm::InputTag> corrInputTags = cfg.getParameter<std::vector<edm::InputTag> >("srcCorrections");
0030 std::vector<edm::EDGetTokenT<CorrMETData> > corrTokens;
0031 for (std::vector<edm::InputTag>::const_iterator inputTag = corrInputTags.begin(); inputTag != corrInputTags.end();
0032 ++inputTag) {
0033 corrTokens.push_back(consumes<CorrMETData>(*inputTag));
0034 }
0035
0036 corrector.setCorTokens(corrTokens);
0037
0038 produces<METCollection>("");
0039 }
0040
0041 ~CorrectedPFMETProducer() override {}
0042
0043 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0044 edm::ParameterSetDescription desc;
0045 desc.add<edm::InputTag>("src", edm::InputTag("pfMet", ""));
0046 std::vector<edm::InputTag> tmpv;
0047 tmpv.push_back(edm::InputTag("corrPfMetType1", "type1"));
0048 desc.add<std::vector<edm::InputTag> >("srcCorrections", tmpv);
0049 descriptions.add(defaultModuleLabel<CorrectedPFMETProducer>(), desc);
0050 }
0051
0052 private:
0053 AddCorrectionsToGenericMET corrector;
0054
0055 typedef std::vector<reco::PFMET> METCollection;
0056
0057 edm::EDGetTokenT<METCollection> token_;
0058
0059 void produce(edm::Event& evt, const edm::EventSetup& es) override {
0060 edm::Handle<METCollection> srcMETCollection;
0061 evt.getByToken(token_, srcMETCollection);
0062
0063 const reco::PFMET& srcMET = (*srcMETCollection)[0];
0064
0065 reco::PFMET outMET = corrector.getCorrectedPFMET(srcMET, evt);
0066
0067 std::unique_ptr<METCollection> product(new METCollection);
0068 product->push_back(outMET);
0069 evt.put(std::move(product));
0070 }
0071 };
0072
0073
0074
0075 DEFINE_FWK_MODULE(CorrectedPFMETProducer);