PatAlgos/plugins/PATJetUpdater.cc

0001 /**
0002   \class    pat::PATJetUpdater PATJetUpdater.h "PhysicsTools/PatAlgos/interface/PATJetUpdater.h"
0003   \brief    Produces pat::Jet's
0004
0005    The PATJetUpdater produces analysis-level pat::Jet's starting from
0006    a collection of pat::Jet's and updates information.
0007
0008   \author   Andreas Hinzmann
0009   \version  $Id: PATJetUpdater.h,v 1.00 2014/03/11 18:13:54 srappocc Exp $
0010 */
0011
0012 #include "CommonTools/Utils/interface/PtComparator.h"
0013 #include "DataFormats/Candidate/interface/CandAssociation.h"
0014 #include "DataFormats/Common/interface/Association.h"
0015 #include "DataFormats/Common/interface/ValueMap.h"
0016 #include "DataFormats/Common/interface/View.h"
0017 #include "DataFormats/PatCandidates/interface/Jet.h"
0018 #include "DataFormats/PatCandidates/interface/JetCorrFactors.h"
0019 #include "DataFormats/PatCandidates/interface/UserData.h"
0020 #include "FWCore/Framework/interface/Event.h"
0021 #include "FWCore/Framework/interface/stream/EDProducer.h"
0022 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0023 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0024 #include "FWCore/ParameterSet/interface/FileInPath.h"
0025 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0026 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0027 #include "FWCore/Utilities/interface/InputTag.h"
0028 #include "FWCore/Utilities/interface/transform.h"
0029 #include "PhysicsTools/PatAlgos/interface/PATUserDataHelper.h"
0030
0031 #include <algorithm>
0032 #include <memory>
0033 #include <vector>
0034
0035 namespace pat {
0036
0037   class PATJetUpdater : public edm::stream::EDProducer<> {
0038   public:
0039     explicit PATJetUpdater(const edm::ParameterSet& iConfig);
0040     ~PATJetUpdater() override;
0041
0042     void produce(edm::Event& iEvent, const edm::EventSetup& iSetup) override;
0043
0044     static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0045
0046   private:
0047     // configurables
0048     edm::EDGetTokenT<edm::View<reco::Jet>> jetsToken_;
0049     bool sort_;
0050     bool addJetCorrFactors_;
0051     std::vector<edm::EDGetTokenT<edm::ValueMap<JetCorrFactors>>> jetCorrFactorsTokens_;
0052
0053     bool addBTagInfo_;
0054     bool addDiscriminators_;
0055     std::vector<edm::InputTag> discriminatorTags_;
0056     std::vector<edm::EDGetTokenT<reco::JetFloatAssociation::Container>> discriminatorTokens_;
0057     std::vector<std::string> discriminatorLabels_;
0058     bool addTagInfos_;
0059     std::vector<edm::InputTag> tagInfoTags_;
0060     std::vector<edm::EDGetTokenT<edm::View<reco::BaseTagInfo>>> tagInfoTokens_;
0061     std::vector<std::string> tagInfoLabels_;
0062
0063     GreaterByPt<Jet> pTComparator_;
0064
0065     bool useUserData_;
0066     pat::PATUserDataHelper<pat::Jet> userDataHelper_;
0067     //
0068     bool printWarning_;  // this is introduced to issue warnings only once per job
0069   };
0070
0071 }  // namespace pat
0072
0073 using namespace pat;
0074
0075 PATJetUpdater::PATJetUpdater(const edm::ParameterSet& iConfig)
0076     : useUserData_(iConfig.exists("userData")), printWarning_(iConfig.getParameter<bool>("printWarning")) {
0077   // initialize configurables
0078   jetsToken_ = consumes<edm::View<reco::Jet>>(iConfig.getParameter<edm::InputTag>("jetSource"));
0079   sort_ = iConfig.getParameter<bool>("sort");
0080   addJetCorrFactors_ = iConfig.getParameter<bool>("addJetCorrFactors");
0081   if (addJetCorrFactors_) {
0082     jetCorrFactorsTokens_ = edm::vector_transform(
0083         iConfig.getParameter<std::vector<edm::InputTag>>("jetCorrFactorsSource"),
0084         [this](edm::InputTag const& tag) { return mayConsume<edm::ValueMap<JetCorrFactors>>(tag); });
0085   }
0086   addBTagInfo_ = iConfig.getParameter<bool>("addBTagInfo");
0087   addDiscriminators_ = iConfig.getParameter<bool>("addDiscriminators");
0088   discriminatorTags_ = iConfig.getParameter<std::vector<edm::InputTag>>("discriminatorSources");
0089   discriminatorTokens_ = edm::vector_transform(discriminatorTags_, [this](edm::InputTag const& tag) {
0090     return mayConsume<reco::JetFloatAssociation::Container>(tag);
0091   });
0092   addTagInfos_ = iConfig.getParameter<bool>("addTagInfos");
0093   tagInfoTags_ = iConfig.getParameter<std::vector<edm::InputTag>>("tagInfoSources");
0094   tagInfoTokens_ = edm::vector_transform(
0095       tagInfoTags_, [this](edm::InputTag const& tag) { return mayConsume<edm::View<reco::BaseTagInfo>>(tag); });
0096   if (discriminatorTags_.empty()) {
0097     addDiscriminators_ = false;
0098   } else {
0099     for (std::vector<edm::InputTag>::const_iterator it = discriminatorTags_.begin(), ed = discriminatorTags_.end();
0100          it != ed;
0101          ++it) {
0102       std::string label = it->label();
0103       std::string::size_type pos = label.find("JetTags");
0104       if ((pos != std::string::npos) && (pos != label.length() - 7)) {
0105         label.erase(pos + 7);  // trim a tail after "JetTags"
0106       }
0107       if (!it->instance().empty()) {
0108         label = (label + std::string(":") + it->instance());
0109       }
0110       discriminatorLabels_.push_back(label);
0111     }
0112   }
0113   if (tagInfoTags_.empty()) {
0114     addTagInfos_ = false;
0115   } else {
0116     for (std::vector<edm::InputTag>::const_iterator it = tagInfoTags_.begin(), ed = tagInfoTags_.end(); it != ed;
0117          ++it) {
0118       std::string label = it->label();
0119       std::string::size_type pos = label.find("TagInfos");
0120       if ((pos != std::string::npos) && (pos != label.length() - 8)) {
0121         label.erase(pos + 8);  // trim a tail after "TagInfos"
0122       }
0123       tagInfoLabels_.push_back(label);
0124     }
0125   }
0126   if (!addBTagInfo_) {
0127     addDiscriminators_ = false;
0128     addTagInfos_ = false;
0129   }
0130   // Check to see if the user wants to add user data
0131   if (useUserData_) {
0132     userDataHelper_ = PATUserDataHelper<Jet>(iConfig.getParameter<edm::ParameterSet>("userData"), consumesCollector());
0133   }
0134   // produces vector of jets
0135   produces<std::vector<Jet>>();
0136   produces<edm::OwnVector<reco::BaseTagInfo>>("tagInfos");
0137 }
0138
0139 PATJetUpdater::~PATJetUpdater() {}
0140
0141 void PATJetUpdater::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0142   // Get the vector of jets
0143   edm::Handle<edm::View<reco::Jet>> jets;
0144   iEvent.getByToken(jetsToken_, jets);
0145
0146   // read in the jet correction factors ValueMap
0147   std::vector<edm::ValueMap<JetCorrFactors>> jetCorrs;
0148   if (addJetCorrFactors_) {
0149     for (size_t i = 0; i < jetCorrFactorsTokens_.size(); ++i) {
0150       edm::Handle<edm::ValueMap<JetCorrFactors>> jetCorr;
0151       iEvent.getByToken(jetCorrFactorsTokens_[i], jetCorr);
0152       jetCorrs.push_back(*jetCorr);
0153     }
0154   }
0155
0156   // Get the vector of jet tags with b-tagging info
0157   std::vector<edm::Handle<reco::JetFloatAssociation::Container>> jetDiscriminators;
0158   if (addBTagInfo_ && addDiscriminators_) {
0159     jetDiscriminators.resize(discriminatorTokens_.size());
0160     for (size_t i = 0; i < discriminatorTokens_.size(); ++i) {
0161       iEvent.getByToken(discriminatorTokens_[i], jetDiscriminators[i]);
0162     }
0163   }
0164   std::vector<edm::Handle<edm::View<reco::BaseTagInfo>>> jetTagInfos;
0165   if (addBTagInfo_ && addTagInfos_) {
0166     jetTagInfos.resize(tagInfoTokens_.size());
0167     for (size_t i = 0; i < tagInfoTokens_.size(); ++i) {
0168       iEvent.getByToken(tagInfoTokens_[i], jetTagInfos[i]);
0169     }
0170   }
0171
0172   // loop over jets
0173   auto patJets = std::make_unique<std::vector<Jet>>();
0174
0175   auto tagInfosOut = std::make_unique<edm::OwnVector<reco::BaseTagInfo>>();
0176
0177   edm::RefProd<edm::OwnVector<reco::BaseTagInfo>> h_tagInfosOut =
0178       iEvent.getRefBeforePut<edm::OwnVector<reco::BaseTagInfo>>("tagInfos");
0179
0180   for (edm::View<reco::Jet>::const_iterator itJet = jets->begin(); itJet != jets->end(); itJet++) {
0181     // construct the Jet from the ref -> save ref to original object
0182     unsigned int idx = itJet - jets->begin();
0183     const edm::RefToBase<reco::Jet> jetRef = jets->refAt(idx);
0184     const edm::RefToBase<Jet> patJetRef(jetRef.castTo<JetRef>());
0185     Jet ajet(patJetRef);
0186
0187     if (addJetCorrFactors_) {
0188       // undo previous jet energy corrections
0189       ajet.setP4(ajet.correctedP4(0));
0190       // clear previous JetCorrFactors
0191       ajet.jec_.clear();
0192       // add additional JetCorrs to the jet
0193       for (unsigned int i = 0; i < jetCorrFactorsTokens_.size(); ++i) {
0194         const JetCorrFactors& jcf = jetCorrs[i][jetRef];
0195         // uncomment for debugging
0196         // jcf.print();
0197         ajet.addJECFactors(jcf);
0198       }
0199       std::vector<std::string> levels = jetCorrs[0][jetRef].correctionLabels();
0200       if (std::find(levels.begin(), levels.end(), "L2L3Residual") != levels.end()) {
0201         ajet.initializeJEC(jetCorrs[0][jetRef].jecLevel("L2L3Residual"));
0202       } else if (std::find(levels.begin(), levels.end(), "L3Absolute") != levels.end()) {
0203         ajet.initializeJEC(jetCorrs[0][jetRef].jecLevel("L3Absolute"));
0204       } else {
0205         ajet.initializeJEC(jetCorrs[0][jetRef].jecLevel("Uncorrected"));
0206         if (printWarning_) {
0207           edm::LogWarning("L3Absolute not found")
0208               << "L2L3Residual and L3Absolute are not part of the jetCorrFactors\n"
0209               << "of module " << jetCorrs[0][jetRef].jecSet() << ". Jets will remain"
0210               << " uncorrected.";
0211           printWarning_ = false;
0212         }
0213       }
0214     }
0215
0216     // add b-tag info if available & required
0217     if (addBTagInfo_) {
0218       if (addDiscriminators_) {
0219         for (size_t k = 0; k < jetDiscriminators.size(); ++k) {
0220           float value = (*jetDiscriminators[k])[jetRef];
0221           ajet.addBDiscriminatorPair(std::make_pair(discriminatorLabels_[k], value));
0222         }
0223       }
0224       if (addTagInfos_) {
0225         for (size_t k = 0; k < jetTagInfos.size(); ++k) {
0226           const edm::View<reco::BaseTagInfo>& taginfos = *jetTagInfos[k];
0227           // This is not associative, so we have to search the jet
0228           edm::Ptr<reco::BaseTagInfo> match;
0229           // Try first by 'same index'
0230           if ((idx < taginfos.size()) && (taginfos[idx].jet() == jetRef)) {
0231             match = taginfos.ptrAt(idx);
0232           } else {
0233             // otherwise fail back to a simple search
0234             for (edm::View<reco::BaseTagInfo>::const_iterator itTI = taginfos.begin(), edTI = taginfos.end();
0235                  itTI != edTI;
0236                  ++itTI) {
0237               if (itTI->jet() == jetRef) {
0238                 match = taginfos.ptrAt(itTI - taginfos.begin());
0239                 break;
0240               }
0241             }
0242           }
0243           if (match.isNonnull()) {
0244             tagInfosOut->push_back(match->clone());
0245             // set the "forward" ptr to the thinned collection
0246             edm::Ptr<reco::BaseTagInfo> tagInfoForwardPtr(
0247                 h_tagInfosOut.id(), &tagInfosOut->back(), tagInfosOut->size() - 1);
0248             // set the "backward" ptr to the original collection for association
0249             const edm::Ptr<reco::BaseTagInfo>& tagInfoBackPtr(match);
0250             // make FwdPtr
0251             TagInfoFwdPtrCollection::value_type tagInfoFwdPtr(tagInfoForwardPtr, tagInfoBackPtr);
0252             ajet.addTagInfo(tagInfoLabels_[k], tagInfoFwdPtr);
0253           }
0254         }
0255       }
0256     }
0257
0258     if (useUserData_) {
0259       userDataHelper_.add(ajet, iEvent, iSetup);
0260     }
0261
0262     // reassign the original object reference to preserve reference to the original jet the input PAT jet was derived from
0263     // (this needs to be done at the end since cloning the input PAT jet would interfere with adding UserData)
0264     ajet.refToOrig_ = patJetRef->originalObjectRef();
0265
0266     patJets->push_back(ajet);
0267   }
0268
0269   // sort jets in pt
0270   if (sort_) {
0271     std::sort(patJets->begin(), patJets->end(), pTComparator_);
0272   }
0273
0274   // put genEvt  in Event
0275   iEvent.put(std::move(patJets));
0276
0277   iEvent.put(std::move(tagInfosOut), "tagInfos");
0278 }
0279
0280 // ParameterSet description for module
0281 void PATJetUpdater::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0282   edm::ParameterSetDescription iDesc;
0283   iDesc.setComment("PAT jet producer module");
0284
0285   // input source
0286   iDesc.add<edm::InputTag>("jetSource", edm::InputTag("no default"))->setComment("input collection");
0287
0288   // sort inputs (by pt)
0289   iDesc.add<bool>("sort", true);
0290
0291   // tag info
0292   iDesc.add<bool>("addTagInfos", true);
0293   std::vector<edm::InputTag> emptyVInputTags;
0294   iDesc.add<std::vector<edm::InputTag>>("tagInfoSources", emptyVInputTags);
0295
0296   // jet energy corrections
0297   iDesc.add<bool>("addJetCorrFactors", true);
0298   iDesc.add<std::vector<edm::InputTag>>("jetCorrFactorsSource", emptyVInputTags);
0299
0300   // btag discriminator tags
0301   iDesc.add<bool>("addBTagInfo", true);
0302   iDesc.add<bool>("addDiscriminators", true);
0303   iDesc.add<std::vector<edm::InputTag>>("discriminatorSources", emptyVInputTags);
0304
0305   // silent warning if false
0306   iDesc.add<bool>("printWarning", true);
0307
0308   // Check to see if the user wants to add user data
0309   edm::ParameterSetDescription userDataPSet;
0310   PATUserDataHelper<Jet>::fillDescription(userDataPSet);
0311   iDesc.addOptional("userData", userDataPSet);
0312
0313   descriptions.add("PATJetUpdater", iDesc);
0314 }
0315
0316 #include "FWCore/Framework/interface/MakerMacros.h"
0317
0318 DEFINE_FWK_MODULE(PATJetUpdater);