File indexing completed on 2024-04-06 12:24:00
0001 #include "FWCore/Framework/interface/Event.h"
0002 #include "FWCore/Framework/interface/global/EDProducer.h"
0003 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0004
0005 #include "DataFormats/PatCandidates/interface/Jet.h"
0006 #include "DataFormats/PatCandidates/interface/MET.h"
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020
0021
0022
0023
0024
0025
0026
0027
0028
0029
0030
0031
0032
0033
0034 class JetEnergyShift : public edm::global::EDProducer<> {
0035 public:
0036
0037 explicit JetEnergyShift(const edm::ParameterSet&);
0038
0039 private:
0040
0041 void produce(edm::StreamID, edm::Event&, const edm::EventSetup&) const override;
0042
0043 private:
0044
0045 edm::EDGetTokenT<std::vector<pat::Jet>> inputJetsToken_;
0046
0047 edm::EDGetTokenT<std::vector<pat::MET>> inputMETsToken_;
0048
0049 std::string outputJets_;
0050
0051 std::string outputMETs_;
0052
0053 double scaleFactor_;
0054
0055 double jetPTThresholdForMET_;
0056
0057 double jetEMLimitForMET_;
0058 };
0059
0060 JetEnergyShift::JetEnergyShift(const edm::ParameterSet& cfg)
0061 : inputJetsToken_(consumes<std::vector<pat::Jet>>(cfg.getParameter<edm::InputTag>("inputJets"))),
0062 inputMETsToken_(consumes<std::vector<pat::MET>>(cfg.getParameter<edm::InputTag>("inputMETs"))),
0063 scaleFactor_(cfg.getParameter<double>("scaleFactor")),
0064 jetPTThresholdForMET_(cfg.getParameter<double>("jetPTThresholdForMET")),
0065 jetEMLimitForMET_(cfg.getParameter<double>("jetEMLimitForMET")) {
0066
0067 outputJets_ = cfg.getParameter<edm::InputTag>("inputJets").label();
0068 outputMETs_ = cfg.getParameter<edm::InputTag>("inputMETs").label();
0069
0070 produces<std::vector<pat::Jet>>(outputJets_);
0071 produces<std::vector<pat::MET>>(outputMETs_);
0072 }
0073
0074 void JetEnergyShift::produce(edm::StreamID, edm::Event& event, const edm::EventSetup& setup) const {
0075 edm::Handle<std::vector<pat::Jet>> jets;
0076 event.getByToken(inputJetsToken_, jets);
0077
0078 edm::Handle<std::vector<pat::MET>> mets;
0079 event.getByToken(inputMETsToken_, mets);
0080
0081 auto pJets = std::make_unique<std::vector<pat::Jet>>();
0082 auto pMETs = std::make_unique<std::vector<pat::MET>>();
0083
0084 double dPx = 0.;
0085 double dPy = 0.;
0086 double dSumEt = 0.;
0087
0088 for (std::vector<pat::Jet>::const_iterator jet = jets->begin(); jet != jets->end(); ++jet) {
0089 pat::Jet scaledJet = *jet;
0090 scaledJet.scaleEnergy(scaleFactor_);
0091 pJets->push_back(scaledJet);
0092
0093
0094 if (jet->correctedJet("raw").pt() > jetPTThresholdForMET_ && jet->emEnergyFraction() < jetEMLimitForMET_) {
0095 dPx += scaledJet.px() - jet->px();
0096 dPy += scaledJet.py() - jet->py();
0097 dSumEt += scaledJet.et() - jet->et();
0098 }
0099 }
0100
0101
0102 pat::MET met = *(mets->begin());
0103 double scaledMETPx = met.px() - dPx;
0104 double scaledMETPy = met.py() - dPy;
0105 pat::MET scaledMET(
0106 reco::MET(met.sumEt() + dSumEt,
0107 reco::MET::LorentzVector(
0108 scaledMETPx, scaledMETPy, 0, sqrt(scaledMETPx * scaledMETPx + scaledMETPy * scaledMETPy)),
0109 reco::MET::Point(0, 0, 0)));
0110 pMETs->push_back(scaledMET);
0111 event.put(std::move(pJets), outputJets_);
0112 event.put(std::move(pMETs), outputMETs_);
0113 }
0114
0115 #include "FWCore/Framework/interface/MakerMacros.h"
0116 DEFINE_FWK_MODULE(JetEnergyShift);