PatAlgos/plugins/PATMETProducer.cc

0001 /**
0002   \class    pat::PATMETProducer PATMETProducer.h "PhysicsTools/PatAlgos/interface/PATMETProducer.h"
0003   \brief    Produces the pat::MET
0004
0005    The PATMETProducer produces the analysis-level pat::MET starting from
0006    a collection of objects of METType.
0007
0008   \author   Steven Lowette
0009   \version  $Id: PATMETProducer.h,v 1.10 2009/06/25 23:49:35 gpetrucc Exp $
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     // configurables
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     // tools
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     //MET Significance
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 }  // namespace pat
0082
0083 using namespace pat;
0084
0085 PATMETProducer::PATMETProducer(const edm::ParameterSet& iConfig) : useUserData_(iConfig.exists("userData")) {
0086   // initialize the configurables
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   // Efficiency configurables
0094   addEfficiencies_ = iConfig.getParameter<bool>("addEfficiencies");
0095   if (addEfficiencies_) {
0096     efficiencyLoader_ =
0097         pat::helper::EfficiencyLoader(iConfig.getParameter<edm::ParameterSet>("efficiencies"), consumesCollector());
0098   }
0099
0100   // Resolution configurables
0101   addResolutions_ = iConfig.getParameter<bool>("addResolutions");
0102   if (addResolutions_) {
0103     resolutionLoader_ =
0104         pat::helper::KinResolutionsLoader(iConfig.getParameter<edm::ParameterSet>("resolutions"), consumesCollector());
0105   }
0106
0107   // Check to see if the user wants to add user data
0108   if (useUserData_) {
0109     userDataHelper_ = PATUserDataHelper<MET>(iConfig.getParameter<edm::ParameterSet>("userData"), consumesCollector());
0110   }
0111
0112   // MET Significance
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   // produces vector of mets
0135   produces<std::vector<MET>>();
0136 }
0137
0138 PATMETProducer::~PATMETProducer() {}
0139
0140 void PATMETProducer::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0141   // Get the vector of MET's from the event
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   // Get the vector of generated met from the event if needed
0154   edm::Handle<edm::View<reco::GenMET>> genMETs;
0155   if (addGenMET_) {
0156     iEvent.getByToken(genMETToken_, genMETs);
0157   }
0158
0159   // loop over mets
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     // construct the MET from the ref -> save ref to original object
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     // add the generated MET
0168     if (addGenMET_)
0169       amet.setGenMET((*genMETs)[idx]);
0170
0171     //add the MET significance
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     // correct for muons if demanded... never more: it's now done by JetMETCorrections
0194     // add the MET to the vector of METs
0195     patMETs->push_back(amet);
0196   }
0197
0198   // sort MET in ET .. don't mess with this
0199   //  std::sort(patMETs->begin(), patMETs->end(), eTComparator_);
0200
0201   // put genEvt object in Event
0202   std::unique_ptr<std::vector<MET>> myMETs(patMETs);
0203   iEvent.put(std::move(myMETs));
0204 }
0205
0206 // ParameterSet description for module
0207 void PATMETProducer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0208   edm::ParameterSetDescription iDesc;
0209   iDesc.setComment("PAT MET producer module");
0210
0211   // input source
0212   iDesc.add<edm::InputTag>("metSource", edm::InputTag("no default"))->setComment("input collection");
0213
0214   // MC configurations
0215   iDesc.add<bool>("addGenMET", false);
0216   iDesc.add<edm::InputTag>("genMETSource", edm::InputTag("genMetCalo"));
0217
0218   pat::helper::KinResolutionsLoader::fillDescription(iDesc);
0219
0220   // Efficiency configurables
0221   edm::ParameterSetDescription efficienciesPSet;
0222   efficienciesPSet.setAllowAnything();  // TODO: the pat helper needs to implement a description.
0223   iDesc.add("efficiencies", efficienciesPSet);
0224   iDesc.add<bool>("addEfficiencies", false);
0225
0226   // Check to see if the user wants to add user data
0227   edm::ParameterSetDescription userDataPSet;
0228   PATUserDataHelper<MET>::fillDescription(userDataPSet);
0229   iDesc.addOptional("userData", userDataPSet);
0230
0231   // muon correction
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   // jets
0249   edm::Handle<edm::View<reco::Jet>> inputJets;
0250   event.getByToken(jetToken_, inputJets);
0251
0252   //candidates
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   //Compute the covariance matrix and fill it
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);