File indexing completed on 2024-04-06 12:19:26
0001 #include "JetMETCorrections/Type1MET/plugins/MultShiftMETcorrInputProducer.h"
0002
0003 #include "FWCore/Utilities/interface/Exception.h"
0004
0005 #include "DataFormats/METReco/interface/CorrMETData.h"
0006 #include "DataFormats/VertexReco/interface/Vertex.h"
0007 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidate.h"
0008 #include "DataFormats/Common/interface/View.h"
0009 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0010
0011 #include <TString.h>
0012
0013 #include <memory>
0014
0015 int MultShiftMETcorrInputProducer::translateTypeToAbsPdgId(reco::PFCandidate::ParticleType type) {
0016 switch (type) {
0017 case reco::PFCandidate::ParticleType::h:
0018 return 211;
0019 case reco::PFCandidate::ParticleType::e:
0020 return 11;
0021 case reco::PFCandidate::ParticleType::mu:
0022 return 13;
0023 case reco::PFCandidate::ParticleType::gamma:
0024 return 22;
0025 case reco::PFCandidate::ParticleType::h0:
0026 return 130;
0027 case reco::PFCandidate::ParticleType::h_HF:
0028 return 1;
0029 case reco::PFCandidate::ParticleType::egamma_HF:
0030 return 2;
0031 case reco::PFCandidate::ParticleType::X:
0032 default:
0033 return 0;
0034 }
0035 }
0036
0037 MultShiftMETcorrInputProducer::MultShiftMETcorrInputProducer(const edm::ParameterSet& cfg)
0038 : moduleLabel_(cfg.getParameter<std::string>("@module_label")) {
0039 pflow_ = consumes<edm::View<reco::Candidate>>(cfg.getParameter<edm::InputTag>("srcPFlow"));
0040 vertices_ = consumes<edm::View<reco::Vertex>>(cfg.getParameter<edm::InputTag>("vertexCollection"));
0041
0042 edm::InputTag srcWeights = cfg.getParameter<edm::InputTag>("srcWeights");
0043 if (!srcWeights.label().empty())
0044 weightsToken_ = consumes<edm::ValueMap<float>>(srcWeights);
0045
0046 cfgCorrParameters_ = cfg.getParameter<std::vector<edm::ParameterSet>>("parameters");
0047 etaMin_.clear();
0048 etaMax_.clear();
0049 type_.clear();
0050 varType_.clear();
0051
0052 produces<CorrMETData>();
0053
0054 for (std::vector<edm::ParameterSet>::const_iterator v = cfgCorrParameters_.begin(); v != cfgCorrParameters_.end();
0055 v++) {
0056 TString corrPxFormula = v->getParameter<std::string>("fx");
0057 TString corrPyFormula = v->getParameter<std::string>("fy");
0058 std::vector<double> corrPxParams = v->getParameter<std::vector<double>>("px");
0059 std::vector<double> corrPyParams = v->getParameter<std::vector<double>>("py");
0060
0061 formula_x_.push_back(std::make_unique<TF1>(
0062 std::string(moduleLabel_).append("_").append(v->getParameter<std::string>("name")).append("_corrPx").c_str(),
0063 v->getParameter<std::string>("fx").c_str()));
0064 formula_y_.push_back(std::make_unique<TF1>(
0065 std::string(moduleLabel_).append("_").append(v->getParameter<std::string>("name")).append("_corrPy").c_str(),
0066 v->getParameter<std::string>("fy").c_str()));
0067
0068 for (unsigned i = 0; i < corrPxParams.size(); i++)
0069 formula_x_.back()->SetParameter(i, corrPxParams[i]);
0070 for (unsigned i = 0; i < corrPyParams.size(); i++)
0071 formula_y_.back()->SetParameter(i, corrPyParams[i]);
0072
0073 counts_.push_back(0);
0074 sumPt_.push_back(0.);
0075 etaMin_.push_back(v->getParameter<double>("etaMin"));
0076 etaMax_.push_back(v->getParameter<double>("etaMax"));
0077 type_.push_back(v->getParameter<int>("type"));
0078 varType_.push_back(v->getParameter<int>("varType"));
0079 }
0080 }
0081
0082 MultShiftMETcorrInputProducer::~MultShiftMETcorrInputProducer() {}
0083
0084 void MultShiftMETcorrInputProducer::produce(edm::Event& evt, const edm::EventSetup& es) {
0085
0086 edm::Handle<edm::View<reco::Vertex>> hpv;
0087 evt.getByToken(vertices_, hpv);
0088 if (!hpv.isValid()) {
0089 edm::LogError("MultShiftMETcorrInputProducer::produce") << "could not find vertex collection ";
0090 }
0091 uint ngoodVertices = 0;
0092 for (unsigned i = 0; i < hpv->size(); i++) {
0093 if ((*hpv)[i].ndof() > 4 && (fabs((*hpv)[i].z()) <= 24.) && (fabs((*hpv)[i].position().rho()) <= 2.0))
0094 ngoodVertices += 1;
0095 }
0096 uint nVertices = hpv->size();
0097
0098 for (unsigned i = 0; i < counts_.size(); i++)
0099 counts_[i] = 0;
0100 for (unsigned i = 0; i < sumPt_.size(); i++)
0101 sumPt_[i] = 0.;
0102
0103 edm::Handle<edm::View<reco::Candidate>> particleFlow;
0104 evt.getByToken(pflow_, particleFlow);
0105
0106 edm::Handle<edm::ValueMap<float>> weights;
0107 if (!weightsToken_.isUninitialized())
0108 evt.getByToken(weightsToken_, weights);
0109 if (std::find(varType_.begin(), varType_.end(), 0) != varType_.end() ||
0110 std::find(varType_.begin(), varType_.end(), 2) != varType_.end()) {
0111 for (unsigned i = 0; i < particleFlow->size(); ++i) {
0112 const reco::Candidate& c = particleFlow->at(i);
0113 for (unsigned j = 0; j < type_.size(); j++) {
0114 if (abs(c.pdgId()) == translateTypeToAbsPdgId(reco::PFCandidate::ParticleType(type_[j]))) {
0115 if ((c.eta() > etaMin_[j]) and (c.eta() < etaMax_[j])) {
0116 float weight = (!weightsToken_.isUninitialized()) ? (*weights)[particleFlow->ptrAt(i)] : 1.0;
0117 counts_[j] += (weight > 0);
0118 sumPt_[j] += c.pt() * weight;
0119 continue;
0120 }
0121 }
0122 }
0123 }
0124 }
0125
0126
0127 std::unique_ptr<CorrMETData> metCorr(new CorrMETData());
0128
0129 double corx = 0.;
0130 double cory = 0.;
0131
0132 for (std::vector<edm::ParameterSet>::const_iterator v = cfgCorrParameters_.begin(); v != cfgCorrParameters_.end();
0133 v++) {
0134 unsigned j = v - cfgCorrParameters_.begin();
0135
0136 double val(0.);
0137 if (varType_[j] == 0) {
0138 val = counts_[j];
0139 } else if (varType_[j] == 1) {
0140 val = ngoodVertices;
0141 } else if (varType_[j] == 2) {
0142 val = sumPt_[j];
0143 } else if (varType_[j] == 3) {
0144 val = nVertices;
0145 }
0146
0147 corx -= formula_x_[j]->Eval(val);
0148 cory -= formula_y_[j]->Eval(val);
0149
0150 }
0151
0152 metCorr->mex = corx;
0153 metCorr->mey = cory;
0154 evt.put(std::move(metCorr), "");
0155 }
0156
0157 #include "FWCore/Framework/interface/MakerMacros.h"
0158
0159 DEFINE_FWK_MODULE(MultShiftMETcorrInputProducer);
