Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 13:32:46

0001 // -*- C++ -*-
0002 //
0003 // Package:    PhysicsTools/NanoAOD
0004 // Class:      LeptonJetVarProducer
0005 //
0006 /**\class LeptonJetVarProducer LeptonJetVarProducer.cc PhysicsTools/NanoAOD/plugins/LeptonJetVarProducer.cc
0007 
0008  Description: [one line class summary]
0009 
0010  Implementation:
0011      [Notes on implementation]
0012 */
0013 //
0014 // Original Author:  Marco Peruzzi
0015 //         Created:  Tue, 05 Sep 2017 12:24:38 GMT
0016 //
0017 //
0018 
0019 // system include files
0020 #include <memory>
0021 
0022 // user include files
0023 #include "FWCore/Framework/interface/Frameworkfwd.h"
0024 #include "FWCore/Framework/interface/global/EDProducer.h"
0025 
0026 #include "FWCore/Framework/interface/Event.h"
0027 #include "FWCore/Framework/interface/MakerMacros.h"
0028 
0029 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0030 #include "FWCore/Utilities/interface/StreamID.h"
0031 
0032 #include "DataFormats/PatCandidates/interface/Muon.h"
0033 #include "DataFormats/PatCandidates/interface/Jet.h"
0034 #include "DataFormats/PatCandidates/interface/Electron.h"
0035 #include "DataFormats/VertexReco/interface/Vertex.h"
0036 
0037 #include "DataFormats/Common/interface/View.h"
0038 #include "DataFormats/Candidate/interface/VertexCompositePtrCandidate.h"
0039 
0040 #include "PhysicsTools/NanoAOD/interface/MatchingUtils.h"
0041 
0042 //
0043 // class declaration
0044 //
0045 
0046 template <typename T>
0047 class LeptonJetVarProducer : public edm::global::EDProducer<> {
0048 public:
0049   explicit LeptonJetVarProducer(const edm::ParameterSet& iConfig)
0050       : srcJet_(consumes<edm::View<pat::Jet>>(iConfig.getParameter<edm::InputTag>("srcJet"))),
0051         srcLep_(consumes<edm::View<T>>(iConfig.getParameter<edm::InputTag>("srcLep"))),
0052         srcVtx_(consumes<std::vector<reco::Vertex>>(iConfig.getParameter<edm::InputTag>("srcVtx"))) {
0053     produces<edm::ValueMap<float>>("ptRatio");
0054     produces<edm::ValueMap<float>>("ptRel");
0055     produces<edm::ValueMap<float>>("jetNDauChargedMVASel");
0056     produces<edm::ValueMap<reco::CandidatePtr>>("jetForLepJetVar");
0057   }
0058   ~LeptonJetVarProducer() override{};
0059 
0060   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0061 
0062 private:
0063   void produce(edm::StreamID, edm::Event&, const edm::EventSetup&) const override;
0064 
0065   std::tuple<float, float, float> calculatePtRatioRel(edm::Ptr<reco::Candidate> lep,
0066                                                       edm::Ptr<pat::Jet> jet,
0067                                                       const reco::Vertex& vtx) const;
0068 
0069   // ----------member data ---------------------------
0070 
0071   edm::EDGetTokenT<edm::View<pat::Jet>> srcJet_;
0072   edm::EDGetTokenT<edm::View<T>> srcLep_;
0073   edm::EDGetTokenT<std::vector<reco::Vertex>> srcVtx_;
0074 };
0075 
0076 //
0077 // constants, enums and typedefs
0078 //
0079 
0080 //
0081 // static data member definitions
0082 //
0083 
0084 //
0085 // member functions
0086 //
0087 
0088 // ------------ method called to produce the data  ------------
0089 template <typename T>
0090 void LeptonJetVarProducer<T>::produce(edm::StreamID streamID, edm::Event& iEvent, const edm::EventSetup& iSetup) const {
0091   edm::Handle<edm::View<pat::Jet>> srcJet;
0092   iEvent.getByToken(srcJet_, srcJet);
0093   edm::Handle<edm::View<T>> srcLep;
0094   iEvent.getByToken(srcLep_, srcLep);
0095   edm::Handle<std::vector<reco::Vertex>> srcVtx;
0096   iEvent.getByToken(srcVtx_, srcVtx);
0097 
0098   unsigned int nJet = srcJet->size();
0099   unsigned int nLep = srcLep->size();
0100 
0101   std::vector<float> ptRatio(nLep, -1);
0102   std::vector<float> ptRel(nLep, -1);
0103   std::vector<float> jetNDauChargedMVASel(nLep, 0);
0104   std::vector<reco::CandidatePtr> jetForLepJetVar(nLep, reco::CandidatePtr());
0105 
0106   const auto& pv = (*srcVtx)[0];
0107 
0108   for (unsigned int il = 0; il < nLep; il++) {
0109     for (unsigned int ij = 0; ij < nJet; ij++) {
0110       auto lep = srcLep->ptrAt(il);
0111       auto jet = srcJet->ptrAt(ij);
0112       if (matchByCommonSourceCandidatePtr(*lep, *jet)) {
0113         auto res = calculatePtRatioRel(lep, jet, pv);
0114         ptRatio[il] = std::get<0>(res);
0115         ptRel[il] = std::get<1>(res);
0116         jetNDauChargedMVASel[il] = std::get<2>(res);
0117         jetForLepJetVar[il] = jet;
0118         break;  // take leading jet with shared source candidates
0119       }
0120     }
0121   }
0122 
0123   std::unique_ptr<edm::ValueMap<float>> ptRatioV(new edm::ValueMap<float>());
0124   edm::ValueMap<float>::Filler fillerRatio(*ptRatioV);
0125   fillerRatio.insert(srcLep, ptRatio.begin(), ptRatio.end());
0126   fillerRatio.fill();
0127   iEvent.put(std::move(ptRatioV), "ptRatio");
0128 
0129   std::unique_ptr<edm::ValueMap<float>> ptRelV(new edm::ValueMap<float>());
0130   edm::ValueMap<float>::Filler fillerRel(*ptRelV);
0131   fillerRel.insert(srcLep, ptRel.begin(), ptRel.end());
0132   fillerRel.fill();
0133   iEvent.put(std::move(ptRelV), "ptRel");
0134 
0135   std::unique_ptr<edm::ValueMap<float>> jetNDauChargedMVASelV(new edm::ValueMap<float>());
0136   edm::ValueMap<float>::Filler fillerNDau(*jetNDauChargedMVASelV);
0137   fillerNDau.insert(srcLep, jetNDauChargedMVASel.begin(), jetNDauChargedMVASel.end());
0138   fillerNDau.fill();
0139   iEvent.put(std::move(jetNDauChargedMVASelV), "jetNDauChargedMVASel");
0140 
0141   std::unique_ptr<edm::ValueMap<reco::CandidatePtr>> jetForLepJetVarV(new edm::ValueMap<reco::CandidatePtr>());
0142   edm::ValueMap<reco::CandidatePtr>::Filler fillerjetForLepJetVar(*jetForLepJetVarV);
0143   fillerjetForLepJetVar.insert(srcLep, jetForLepJetVar.begin(), jetForLepJetVar.end());
0144   fillerjetForLepJetVar.fill();
0145   iEvent.put(std::move(jetForLepJetVarV), "jetForLepJetVar");
0146 }
0147 
0148 template <typename T>
0149 std::tuple<float, float, float> LeptonJetVarProducer<T>::calculatePtRatioRel(edm::Ptr<reco::Candidate> lep,
0150                                                                              edm::Ptr<pat::Jet> jet,
0151                                                                              const reco::Vertex& vtx) const {
0152   auto rawp4 = jet->correctedP4("Uncorrected");
0153   auto lepp4 = lep->p4();
0154 
0155   if ((rawp4 - lepp4).R() < 1e-4)
0156     return std::tuple<float, float, float>(1.0, 0.0, 0.0);
0157 
0158   auto l1corrFactor = jet->jecFactor("L1FastJet") / jet->jecFactor("Uncorrected");
0159 
0160   auto jetp4 = (rawp4 - lepp4 * (1.0 / l1corrFactor)) * (jet->pt() / rawp4.pt()) + lepp4;
0161   auto ptratio = lepp4.pt() / jetp4.pt();
0162   auto ptrel = lepp4.Vect().Cross((jetp4 - lepp4).Vect().Unit()).R();
0163 
0164   unsigned int jndau = 0;
0165   for (const auto& _d : jet->daughterPtrVector()) {
0166     const auto d = dynamic_cast<const pat::PackedCandidate*>(_d.get());
0167     if (d->charge() == 0)
0168       continue;
0169     if (d->fromPV() <= 1)
0170       continue;
0171     if (deltaR(*d, *lep) > 0.4)
0172       continue;
0173     if (!(d->hasTrackDetails()))
0174       continue;
0175     auto tk = d->pseudoTrack();
0176     if (tk.pt() > 1 && tk.hitPattern().numberOfValidHits() >= 8 && tk.hitPattern().numberOfValidPixelHits() >= 2 &&
0177         tk.normalizedChi2() < 5 && fabs(tk.dxy(vtx.position())) < 0.2 && fabs(tk.dz(vtx.position())) < 17)
0178       jndau++;
0179   }
0180 
0181   return std::tuple<float, float, float>(ptratio, ptrel, float(jndau));
0182 }
0183 
0184 // ------------ method fills 'descriptions' with the allowed parameters for the module  ------------
0185 template <typename T>
0186 void LeptonJetVarProducer<T>::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0187   //The following says we do not know what parameters are allowed so do no validation
0188   // Please change this to state exactly what you do use, even if it is no parameters
0189   edm::ParameterSetDescription desc;
0190   desc.add<edm::InputTag>("srcJet")->setComment("jet input collection");
0191   desc.add<edm::InputTag>("srcLep")->setComment("lepton input collection");
0192   desc.add<edm::InputTag>("srcVtx")->setComment("primary vertex input collection");
0193   std::string modname;
0194   if (typeid(T) == typeid(pat::Muon))
0195     modname += "Muon";
0196   else if (typeid(T) == typeid(pat::Electron))
0197     modname += "Electron";
0198   modname += "JetVarProducer";
0199   descriptions.add(modname, desc);
0200 }
0201 
0202 typedef LeptonJetVarProducer<pat::Muon> MuonJetVarProducer;
0203 typedef LeptonJetVarProducer<pat::Electron> ElectronJetVarProducer;
0204 
0205 //define this as a plug-in
0206 DEFINE_FWK_MODULE(MuonJetVarProducer);
0207 DEFINE_FWK_MODULE(ElectronJetVarProducer);