Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 14:26:58

0001 // -*- C++ -*-
0002 //
0003 // Package:    RecoTauTag/HLTProducers
0004 // Class:      L1TJetsMatching
0005 //
0006 /**\class L1TJetsMatching L1TJetsMatching.h 
0007  RecoTauTag/HLTProducers/interface/L1TJetsMatching.h
0008  Description: 
0009  Matching L1 to PF/Calo Jets. Used for HLT_VBF paths.
0010     *Matches PF/Calo Jets to L1 jets from the dedicated seed
0011     *Adds selection criteria to the leading/subleading jets as well as the maximum dijet mass
0012     *Separates collections of PF/Calo jets into two categories
0013  
0014  
0015 */
0016 //
0017 // Original Author:  Vukasin Milosevic
0018 //         Created:  Thu, 01 Jun 2017 17:23:00 GMT
0019 //
0020 //
0021 
0022 #ifndef RecoTauTag_HLTProducers_L1TJetsMatching_h
0023 #define RecoTauTag_HLTProducers_L1TJetsMatching_h
0024 
0025 // user include files
0026 #include "FWCore/Framework/interface/Frameworkfwd.h"
0027 #include "FWCore/Framework/interface/global/EDProducer.h"
0028 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0029 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0030 #include "FWCore/Framework/interface/Event.h"
0031 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0032 #include "FWCore/Utilities/interface/InputTag.h"
0033 #include "DataFormats/Common/interface/Handle.h"
0034 #include "DataFormats/JetReco/interface/PFJetCollection.h"
0035 #include "DataFormats/TauReco/interface/PFTauFwd.h"
0036 #include "DataFormats/HLTReco/interface/TriggerFilterObjectWithRefs.h"
0037 
0038 #include "Math/GenVector/VectorUtil.h"
0039 #include "DataFormats/HLTReco/interface/TriggerTypeDefs.h"
0040 #include "FWCore/Utilities/interface/EDMException.h"
0041 #include "DataFormats/JetReco/interface/PFJet.h"
0042 
0043 #include "HLTrigger/HLTcore/interface/defaultModuleLabel.h"
0044 #include "DataFormats/Math/interface/deltaR.h"
0045 
0046 #include <map>
0047 #include <vector>
0048 
0049 template <typename T>
0050 class L1TJetsMatching : public edm::global::EDProducer<> {
0051 public:
0052   explicit L1TJetsMatching(const edm::ParameterSet&);
0053   ~L1TJetsMatching() override;
0054   void produce(edm::StreamID, edm::Event&, const edm::EventSetup&) const override;
0055   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0056 
0057 private:
0058   const edm::EDGetTokenT<std::vector<T>> jetSrc_;
0059   const edm::EDGetTokenT<trigger::TriggerFilterObjectWithRefs> jetTrigger_;
0060   const double pt1Min_;
0061   const double pt2Min_;
0062   const double pt3Min_;
0063   const double mjjMin_;
0064   const double matchingR_;
0065   const double matchingR2_;
0066 };
0067 //
0068 // class decleration
0069 //
0070 template <typename T>
0071 std::pair<std::vector<T>, std::vector<T>> categorise(
0072     const std::vector<T>& pfMatchedJets, double pt1, double pt2, double pt3, double Mjj) {
0073   std::pair<std::vector<T>, std::vector<T>> output;
0074   unsigned int i1 = 0;
0075   unsigned int i2 = 0;
0076   double mjj = 0;
0077   if (pfMatchedJets.size() > 1) {
0078     for (unsigned int i = 0; i < pfMatchedJets.size() - 1; i++) {
0079       const T& myJet1 = (pfMatchedJets)[i];
0080 
0081       for (unsigned int j = i + 1; j < pfMatchedJets.size(); j++) {
0082         const T& myJet2 = (pfMatchedJets)[j];
0083 
0084         const double mjj_test = (myJet1.p4() + myJet2.p4()).M();
0085 
0086         if (mjj_test > mjj) {
0087           mjj = mjj_test;
0088           i1 = i;
0089           i2 = j;
0090         }
0091       }
0092     }
0093 
0094     const T& myJet1 = (pfMatchedJets)[i1];
0095     const T& myJet2 = (pfMatchedJets)[i2];
0096 
0097     if ((mjj > Mjj) && (myJet1.pt() >= pt1) && (myJet2.pt() > pt2)) {
0098       output.first.push_back(myJet1);
0099       output.first.push_back(myJet2);
0100     }
0101 
0102     if ((mjj > Mjj) && (myJet1.pt() < pt3) && (myJet1.pt() > pt2) && (myJet2.pt() > pt2)) {
0103       const T& myJetTest = (pfMatchedJets)[0];
0104       if (myJetTest.pt() > pt3) {
0105         output.second.push_back(myJet1);
0106         output.second.push_back(myJet2);
0107         output.second.push_back(myJetTest);
0108       }
0109     }
0110   }
0111 
0112   return output;
0113 }
0114 template <typename T>
0115 L1TJetsMatching<T>::L1TJetsMatching(const edm::ParameterSet& iConfig)
0116     : jetSrc_(consumes<std::vector<T>>(iConfig.getParameter<edm::InputTag>("JetSrc"))),
0117       jetTrigger_(consumes<trigger::TriggerFilterObjectWithRefs>(iConfig.getParameter<edm::InputTag>("L1JetTrigger"))),
0118       pt1Min_(iConfig.getParameter<double>("pt1Min")),
0119       pt2Min_(iConfig.getParameter<double>("pt2Min")),
0120       pt3Min_(iConfig.getParameter<double>("pt3Min")),
0121       mjjMin_(iConfig.getParameter<double>("mjjMin")),
0122       matchingR_(iConfig.getParameter<double>("matchingR")),
0123       matchingR2_(matchingR_ * matchingR_) {
0124   produces<std::vector<T>>("TwoJets");
0125   produces<std::vector<T>>("ThreeJets");
0126 }
0127 template <typename T>
0128 L1TJetsMatching<T>::~L1TJetsMatching() {}
0129 
0130 template <typename T>
0131 void L1TJetsMatching<T>::produce(edm::StreamID iSId, edm::Event& iEvent, const edm::EventSetup& iES) const {
0132   unique_ptr<std::vector<T>> pfMatchedJets(new std::vector<T>);
0133   std::pair<std::vector<T>, std::vector<T>> output;
0134 
0135   // Getting HLT jets to be matched
0136   edm::Handle<std::vector<T>> pfJets;
0137   iEvent.getByToken(jetSrc_, pfJets);
0138 
0139   edm::Handle<trigger::TriggerFilterObjectWithRefs> l1TriggeredJets;
0140   iEvent.getByToken(jetTrigger_, l1TriggeredJets);
0141 
0142   //l1t::TauVectorRef jetCandRefVec;
0143   l1t::JetVectorRef jetCandRefVec;
0144   l1TriggeredJets->getObjects(trigger::TriggerL1Jet, jetCandRefVec);
0145 
0146   math::XYZPoint a(0., 0., 0.);
0147 
0148   //std::cout<<"PFsize= "<<pfJets->size()<<endl<<" L1size= "<<jetCandRefVec.size()<<std::endl;
0149   for (unsigned int iJet = 0; iJet < pfJets->size(); iJet++) {
0150     const T& myJet = (*pfJets)[iJet];
0151     for (unsigned int iL1Jet = 0; iL1Jet < jetCandRefVec.size(); iL1Jet++) {
0152       // Find the relative L2pfJets, to see if it has been reconstructed
0153       //  if ((iJet<3) && (iL1Jet==0))  std::cout<<myJet.p4().Pt()<<" ";
0154       if ((reco::deltaR2(myJet.p4(), jetCandRefVec[iL1Jet]->p4()) < matchingR2_) && (myJet.pt() > pt2Min_)) {
0155         pfMatchedJets->push_back(myJet);
0156         break;
0157       }
0158     }
0159   }
0160 
0161   output = categorise(*pfMatchedJets, pt1Min_, pt2Min_, pt3Min_, mjjMin_);
0162   unique_ptr<std::vector<T>> output1(new std::vector<T>(output.first));
0163   unique_ptr<std::vector<T>> output2(new std::vector<T>(output.second));
0164 
0165   iEvent.put(std::move(output1), "TwoJets");
0166   iEvent.put(std::move(output2), "ThreeJets");
0167 }
0168 template <typename T>
0169 void L1TJetsMatching<T>::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0170   edm::ParameterSetDescription desc;
0171   desc.add<edm::InputTag>("L1JetTrigger", edm::InputTag("hltL1DiJetVBF"))->setComment("Name of trigger filter");
0172   desc.add<edm::InputTag>("JetSrc", edm::InputTag("hltAK4PFJetsTightIDCorrected"))
0173       ->setComment("Input collection of PFJets");
0174   desc.add<double>("pt1Min", 110.0)->setComment("Minimal pT1 of PFJets to match");
0175   desc.add<double>("pt2Min", 35.0)->setComment("Minimal pT2 of PFJets to match");
0176   desc.add<double>("pt3Min", 110.0)->setComment("Minimum pT3 of PFJets to match");
0177   desc.add<double>("mjjMin", 650.0)->setComment("Minimal mjj of matched PFjets");
0178   desc.add<double>("matchingR", 0.5)->setComment("dR value used for matching");
0179   descriptions.setComment(
0180       "This module produces collection of PFJetss matched to L1 Taus / Jets passing a HLT filter (Only p4 and vertex "
0181       "of returned PFJetss are set).");
0182   descriptions.add(defaultModuleLabel<L1TJetsMatching<T>>(), desc);
0183 }
0184 
0185 #endif