File indexing completed on 2023-03-17 11:16:16
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012 #include "CommonTools/Utils/interface/EtComparator.h"
0013 #include "DataFormats/Candidate/interface/Candidate.h"
0014 #include "DataFormats/Common/interface/View.h"
0015 #include "DataFormats/PatCandidates/interface/MET.h"
0016 #include "DataFormats/PatCandidates/interface/UserData.h"
0017 #include "FWCore/Framework/interface/Event.h"
0018 #include "FWCore/Framework/interface/stream/EDProducer.h"
0019 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0020 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0021 #include "FWCore/ParameterSet/interface/FileInPath.h"
0022 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0023 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0024 #include "PhysicsTools/PatAlgos/interface/EfficiencyLoader.h"
0025 #include "PhysicsTools/PatAlgos/interface/KinResolutionsLoader.h"
0026 #include "PhysicsTools/PatAlgos/interface/PATUserDataHelper.h"
0027 #include "RecoMET/METAlgorithms/interface/METSignificance.h"
0028 #include "CondFormats/DataRecord/interface/JetResolutionRcd.h"
0029 #include "CondFormats/DataRecord/interface/JetResolutionScaleFactorRcd.h"
0030
0031 #include <memory>
0032
0033 namespace pat {
0034
0035 class PATMETProducer : public edm::stream::EDProducer<> {
0036 public:
0037 explicit PATMETProducer(const edm::ParameterSet& iConfig);
0038 ~PATMETProducer() override;
0039
0040 void produce(edm::Event& iEvent, const edm::EventSetup& iSetup) override;
0041
0042 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0043
0044 private:
0045
0046 edm::InputTag metSrc_;
0047 edm::EDGetTokenT<edm::View<reco::MET>> metToken_;
0048 bool addGenMET_;
0049 edm::EDGetTokenT<edm::View<reco::GenMET>> genMETToken_;
0050 bool addResolutions_;
0051 pat::helper::KinResolutionsLoader resolutionLoader_;
0052 bool addMuonCorr_;
0053 edm::InputTag muonSrc_;
0054
0055 GreaterByEt<MET> eTComparator_;
0056
0057 bool addEfficiencies_;
0058 pat::helper::EfficiencyLoader efficiencyLoader_;
0059
0060 bool useUserData_;
0061 pat::PATUserDataHelper<pat::MET> userDataHelper_;
0062
0063
0064 bool calculateMETSignificance_;
0065 metsig::METSignificance* metSigAlgo_;
0066 edm::EDGetTokenT<edm::View<reco::Jet>> jetToken_;
0067 edm::EDGetTokenT<edm::View<reco::Candidate>> pfCandToken_;
0068 std::vector<edm::EDGetTokenT<edm::View<reco::Candidate>>> lepTokens_;
0069 edm::EDGetTokenT<double> rhoToken_;
0070 edm::EDGetTokenT<edm::ValueMap<float>> weightsToken_;
0071 edm::ESGetToken<JME::JetResolutionObject, JetResolutionRcd> jetResPtToken_;
0072 edm::ESGetToken<JME::JetResolutionObject, JetResolutionRcd> jetResPhiToken_;
0073 edm::ESGetToken<JME::JetResolutionObject, JetResolutionScaleFactorRcd> jetSFToken_;
0074
0075 const reco::METCovMatrix getMETCovMatrix(const edm::Event& event,
0076 const edm::EventSetup& iSetup,
0077 const reco::MET& met,
0078 double& sumPtUnclustered) const;
0079 };
0080
0081 }
0082
0083 using namespace pat;
0084
0085 PATMETProducer::PATMETProducer(const edm::ParameterSet& iConfig) : useUserData_(iConfig.exists("userData")) {
0086
0087 metSrc_ = iConfig.getParameter<edm::InputTag>("metSource");
0088 metToken_ = consumes<edm::View<reco::MET>>(metSrc_);
0089 addGenMET_ = iConfig.getParameter<bool>("addGenMET");
0090 genMETToken_ = mayConsume<edm::View<reco::GenMET>>(iConfig.getParameter<edm::InputTag>("genMETSource"));
0091 addResolutions_ = iConfig.getParameter<bool>("addResolutions");
0092
0093
0094 addEfficiencies_ = iConfig.getParameter<bool>("addEfficiencies");
0095 if (addEfficiencies_) {
0096 efficiencyLoader_ =
0097 pat::helper::EfficiencyLoader(iConfig.getParameter<edm::ParameterSet>("efficiencies"), consumesCollector());
0098 }
0099
0100
0101 addResolutions_ = iConfig.getParameter<bool>("addResolutions");
0102 if (addResolutions_) {
0103 resolutionLoader_ =
0104 pat::helper::KinResolutionsLoader(iConfig.getParameter<edm::ParameterSet>("resolutions"), consumesCollector());
0105 }
0106
0107
0108 if (useUserData_) {
0109 userDataHelper_ = PATUserDataHelper<MET>(iConfig.getParameter<edm::ParameterSet>("userData"), consumesCollector());
0110 }
0111
0112
0113 calculateMETSignificance_ = iConfig.getParameter<bool>("computeMETSignificance");
0114 if (calculateMETSignificance_) {
0115 edm::InputTag srcWeights = iConfig.getParameter<edm::InputTag>("srcWeights");
0116 if (!srcWeights.label().empty())
0117 weightsToken_ = consumes<edm::ValueMap<float>>(srcWeights);
0118 metSigAlgo_ = new metsig::METSignificance(iConfig);
0119 rhoToken_ = consumes<double>(iConfig.getParameter<edm::InputTag>("srcRho"));
0120 jetToken_ = consumes<edm::View<reco::Jet>>(iConfig.getParameter<edm::InputTag>("srcJets"));
0121 pfCandToken_ = consumes<edm::View<reco::Candidate>>(iConfig.getParameter<edm::InputTag>("srcPFCands"));
0122 std::vector<edm::InputTag> srcLeptonsTags = iConfig.getParameter<std::vector<edm::InputTag>>("srcLeptons");
0123 for (std::vector<edm::InputTag>::const_iterator it = srcLeptonsTags.begin(); it != srcLeptonsTags.end(); it++) {
0124 lepTokens_.push_back(consumes<edm::View<reco::Candidate>>(*it));
0125 }
0126 auto jetSFType = iConfig.getParameter<std::string>("srcJetSF");
0127 auto jetResPtType = iConfig.getParameter<std::string>("srcJetResPt");
0128 auto jetResPhiType = iConfig.getParameter<std::string>("srcJetResPhi");
0129 jetResPtToken_ = esConsumes(edm::ESInputTag("", jetResPtType));
0130 jetResPhiToken_ = esConsumes(edm::ESInputTag("", jetResPhiType));
0131 jetSFToken_ = esConsumes(edm::ESInputTag("", jetSFType));
0132 }
0133
0134
0135 produces<std::vector<MET>>();
0136 }
0137
0138 PATMETProducer::~PATMETProducer() {}
0139
0140 void PATMETProducer::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0141
0142 edm::Handle<edm::View<reco::MET>> mets;
0143 iEvent.getByToken(metToken_, mets);
0144
0145 if (mets->size() != 1)
0146 throw cms::Exception("Corrupt Data") << "The input MET collection " << metSrc_.encode() << " has size "
0147 << mets->size() << " instead of 1 as it should.\n";
0148 if (efficiencyLoader_.enabled())
0149 efficiencyLoader_.newEvent(iEvent);
0150 if (resolutionLoader_.enabled())
0151 resolutionLoader_.newEvent(iEvent, iSetup);
0152
0153
0154 edm::Handle<edm::View<reco::GenMET>> genMETs;
0155 if (addGenMET_) {
0156 iEvent.getByToken(genMETToken_, genMETs);
0157 }
0158
0159
0160 std::vector<MET>* patMETs = new std::vector<MET>();
0161 for (edm::View<reco::MET>::const_iterator itMET = mets->begin(); itMET != mets->end(); itMET++) {
0162
0163 unsigned int idx = itMET - mets->begin();
0164 edm::RefToBase<reco::MET> metsRef = mets->refAt(idx);
0165 edm::Ptr<reco::MET> metsPtr = mets->ptrAt(idx);
0166 MET amet(metsRef);
0167
0168 if (addGenMET_)
0169 amet.setGenMET((*genMETs)[idx]);
0170
0171
0172 if (calculateMETSignificance_) {
0173 double sumPtUnclustered = 0;
0174 const reco::METCovMatrix& sigcov = getMETCovMatrix(iEvent, iSetup, amet, sumPtUnclustered);
0175 amet.setSignificanceMatrix(sigcov);
0176 double metSig = metSigAlgo_->getSignificance(sigcov, amet);
0177 amet.setMETSignificance(metSig);
0178 amet.setMETSumPtUnclustered(sumPtUnclustered);
0179 }
0180
0181 if (efficiencyLoader_.enabled()) {
0182 efficiencyLoader_.setEfficiencies(amet, metsRef);
0183 }
0184
0185 if (resolutionLoader_.enabled()) {
0186 resolutionLoader_.setResolutions(amet);
0187 }
0188
0189 if (useUserData_) {
0190 userDataHelper_.add(amet, iEvent, iSetup);
0191 }
0192
0193
0194
0195 patMETs->push_back(amet);
0196 }
0197
0198
0199
0200
0201
0202 std::unique_ptr<std::vector<MET>> myMETs(patMETs);
0203 iEvent.put(std::move(myMETs));
0204 }
0205
0206
0207 void PATMETProducer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0208 edm::ParameterSetDescription iDesc;
0209 iDesc.setComment("PAT MET producer module");
0210
0211
0212 iDesc.add<edm::InputTag>("metSource", edm::InputTag("no default"))->setComment("input collection");
0213
0214
0215 iDesc.add<bool>("addGenMET", false);
0216 iDesc.add<edm::InputTag>("genMETSource", edm::InputTag("genMetCalo"));
0217
0218 pat::helper::KinResolutionsLoader::fillDescription(iDesc);
0219
0220
0221 edm::ParameterSetDescription efficienciesPSet;
0222 efficienciesPSet.setAllowAnything();
0223 iDesc.add("efficiencies", efficienciesPSet);
0224 iDesc.add<bool>("addEfficiencies", false);
0225
0226
0227 edm::ParameterSetDescription userDataPSet;
0228 PATUserDataHelper<MET>::fillDescription(userDataPSet);
0229 iDesc.addOptional("userData", userDataPSet);
0230
0231
0232 iDesc.add<bool>("addMuonCorrections", false);
0233 iDesc.add<edm::InputTag>("muonSource", edm::InputTag("muons"));
0234 }
0235
0236 const reco::METCovMatrix PATMETProducer::getMETCovMatrix(const edm::Event& event,
0237 const edm::EventSetup& iSetup,
0238 const reco::MET& met,
0239 double& sumPtUnclustered) const {
0240 std::vector<edm::Handle<reco::CandidateView>> leptons;
0241 for (std::vector<edm::EDGetTokenT<edm::View<reco::Candidate>>>::const_iterator srcLeptons_i = lepTokens_.begin();
0242 srcLeptons_i != lepTokens_.end();
0243 ++srcLeptons_i) {
0244 edm::Handle<reco::CandidateView> leptons_i;
0245 event.getByToken(*srcLeptons_i, leptons_i);
0246 leptons.push_back(leptons_i);
0247 }
0248
0249 edm::Handle<edm::View<reco::Jet>> inputJets;
0250 event.getByToken(jetToken_, inputJets);
0251
0252
0253 edm::Handle<edm::View<reco::Candidate>> inputCands;
0254 event.getByToken(pfCandToken_, inputCands);
0255
0256 edm::Handle<double> rho;
0257 event.getByToken(rhoToken_, rho);
0258
0259 edm::Handle<edm::ValueMap<float>> weights;
0260 if (!weightsToken_.isUninitialized())
0261 event.getByToken(weightsToken_, weights);
0262
0263 JME::JetResolution resPtObj = iSetup.getData(jetResPtToken_);
0264 JME::JetResolution resPhiObj = iSetup.getData(jetResPhiToken_);
0265 JME::JetResolutionScaleFactor resSFObj = iSetup.getData(jetSFToken_);
0266
0267
0268 const edm::ValueMap<float>* weightsPtr = nullptr;
0269 if (met.isWeighted()) {
0270 if (weightsToken_.isUninitialized())
0271 throw cms::Exception("InvalidInput") << "MET is weighted (e.g. PUPPI), but no weights given in PATMETProducer\n";
0272 weightsPtr = &*weights;
0273 }
0274 reco::METCovMatrix cov = metSigAlgo_->getCovariance(*inputJets,
0275 leptons,
0276 inputCands,
0277 *rho,
0278 resPtObj,
0279 resPhiObj,
0280 resSFObj,
0281 event.isRealData(),
0282 sumPtUnclustered,
0283 weightsPtr);
0284
0285 return cov;
0286 }
0287
0288 #include "FWCore/Framework/interface/MakerMacros.h"
0289 DEFINE_FWK_MODULE(PATMETProducer);