Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-03-05 02:39:29

0001 /** Description: Check correlation between PFJet pairs and filtered PFTau pairs and store the PFJet pairs.
0002 For (j1, j2, t1, t2) where j1, j2 from the PFJet collection and t1, t2 from the filtered PFTau collection,
0003 the module checks if there is no overlap (within dRmin) between j1, j2, t1, t2, i.e. they are 4 different objects.
0004 In addition, the module imposes the following cuts:
0005 * mjjMin: the min invariant mass cut on (j1, j2)
0006 * extraTauPtCut: the leading tau pt cut on (t1, t2) (under the assumption t1, t2 are products of a subleading pt filter with minN = 2)
0007 The module stores j1, j2 of any (j1, j2, t1, t2) that satisfies the conditions above. */
0008 
0009 // user include files
0010 #include "FWCore/Framework/interface/Frameworkfwd.h"
0011 #include "FWCore/Framework/interface/global/EDProducer.h"
0012 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0013 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0014 #include "FWCore/Framework/interface/Event.h"
0015 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0016 #include "FWCore/Utilities/interface/InputTag.h"
0017 #include "FWCore/Utilities/interface/Exception.h"
0018 #include "DataFormats/HLTReco/interface/TriggerFilterObjectWithRefs.h"
0019 #include "DataFormats/HLTReco/interface/TriggerTypeDefs.h"
0020 #include "DataFormats/TauReco/interface/PFTau.h"
0021 #include "DataFormats/Math/interface/deltaR.h"
0022 #include "CommonTools/Utils/interface/PtComparator.h"
0023 
0024 //
0025 // class definition
0026 //
0027 class HLTPFDiJetCorrCheckerWithDiTau : public edm::global::EDProducer<> {
0028 public:
0029   explicit HLTPFDiJetCorrCheckerWithDiTau(const edm::ParameterSet&);
0030   void produce(edm::StreamID, edm::Event&, const edm::EventSetup&) const override;
0031   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0032 
0033 private:
0034   const edm::EDGetTokenT<trigger::TriggerFilterObjectWithRefs> tauSrc_;
0035   const edm::EDGetTokenT<reco::PFJetCollection> pfJetSrc_;
0036   const double extraTauPtCut_;
0037   const double mjjMin_;
0038   const double dRmin_, dRmin2_;
0039   // pt comparator
0040   GreaterByPt<reco::PFJet> pTComparator_;
0041 };
0042 
0043 //
0044 // class declaration
0045 //
0046 HLTPFDiJetCorrCheckerWithDiTau::HLTPFDiJetCorrCheckerWithDiTau(const edm::ParameterSet& iConfig)
0047     : tauSrc_(consumes(iConfig.getParameter<edm::InputTag>("tauSrc"))),
0048       pfJetSrc_(consumes(iConfig.getParameter<edm::InputTag>("pfJetSrc"))),
0049       extraTauPtCut_(iConfig.getParameter<double>("extraTauPtCut")),
0050       mjjMin_(iConfig.getParameter<double>("mjjMin")),
0051       dRmin_(iConfig.getParameter<double>("dRmin")),
0052       dRmin2_(dRmin_ * dRmin_) {
0053   if (dRmin_ <= 0.) {
0054     throw cms::Exception("HLTPFDiJetCorrCheckerWithDiTau")
0055         << "invalid value for parameter \"dRmin\" (must be > 0): " << dRmin_;
0056   }
0057   produces<reco::PFJetCollection>();
0058 }
0059 
0060 void HLTPFDiJetCorrCheckerWithDiTau::produce(edm::StreamID iSId, edm::Event& iEvent, const edm::EventSetup& iES) const {
0061   auto const& pfJets = iEvent.get(pfJetSrc_);
0062 
0063   auto cleanedPFJets = std::make_unique<reco::PFJetCollection>();
0064   cleanedPFJets->reserve(pfJets.size());
0065 
0066   trigger::VRpftau taus;
0067   iEvent.get(tauSrc_).getObjects(trigger::TriggerTau, taus);
0068 
0069   std::set<unsigned int> indices;
0070 
0071   if (pfJets.size() > 1 && taus.size() > 1) {
0072     for (unsigned int iJet1 = 0; iJet1 < pfJets.size(); iJet1++) {
0073       for (unsigned int iJet2 = iJet1 + 1; iJet2 < pfJets.size(); iJet2++) {
0074         bool correctComb = false;
0075         const reco::PFJet& myPFJet1 = pfJets[iJet1];
0076         const reco::PFJet& myPFJet2 = pfJets[iJet2];
0077 
0078         if ((myPFJet1.p4() + myPFJet2.p4()).M() < mjjMin_)
0079           continue;
0080 
0081         for (unsigned int iTau1 = 0; iTau1 < taus.size(); iTau1++) {
0082           if (reco::deltaR2(taus[iTau1]->p4(), myPFJet1.p4()) < dRmin2_)
0083             continue;
0084           if (reco::deltaR2(taus[iTau1]->p4(), myPFJet2.p4()) < dRmin2_)
0085             continue;
0086 
0087           for (unsigned int iTau2 = iTau1 + 1; iTau2 < taus.size(); iTau2++) {
0088             if (taus[iTau1]->pt() < extraTauPtCut_ && taus[iTau2]->pt() < extraTauPtCut_)
0089               continue;
0090 
0091             if (reco::deltaR2(taus[iTau2]->p4(), myPFJet1.p4()) < dRmin2_)
0092               continue;
0093             if (reco::deltaR2(taus[iTau2]->p4(), myPFJet2.p4()) < dRmin2_)
0094               continue;
0095 
0096             correctComb = true;
0097             break;
0098           }
0099           if (correctComb)
0100             break;
0101         }
0102 
0103         if (correctComb) {
0104           indices.insert(iJet1);
0105           indices.insert(iJet2);
0106         }
0107       }
0108     }
0109 
0110     for (const auto& i : indices)
0111       cleanedPFJets->emplace_back(pfJets[i]);
0112   }
0113   // sort jets in pt
0114   std::sort(cleanedPFJets->begin(), cleanedPFJets->end(), pTComparator_);
0115   iEvent.put(std::move(cleanedPFJets));
0116 }
0117 
0118 void HLTPFDiJetCorrCheckerWithDiTau::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0119   edm::ParameterSetDescription desc;
0120   desc.add<edm::InputTag>("pfJetSrc", edm::InputTag("hltAK4PFJetsCorrected"))->setComment("Input collection of PFJets");
0121   desc.add<edm::InputTag>("tauSrc", edm::InputTag("hltSinglePFTau20TrackPt1LooseChargedIsolationReg"))
0122       ->setComment("Input collection of PFTaus that have passed ID and isolation requirements");
0123   desc.add<double>("extraTauPtCut", 45)->setComment("In case of asymmetric tau pt cuts");
0124   desc.add<double>("mjjMin", 500)->setComment("VBF dijet mass condition");
0125   desc.add<double>("dRmin", 0.5)->setComment("Minimum dR between PFJets and filtered PFTaus");
0126   descriptions.setComment(
0127       "This module produces a collection of PFJets that are cross-cleaned with respect to PFTaus passing a HLT "
0128       "filter.");
0129   descriptions.addWithDefaultLabel(desc);
0130 }
0131 
0132 //
0133 // module registration
0134 //
0135 #include "FWCore/Framework/interface/MakerMacros.h"
0136 DEFINE_FWK_MODULE(HLTPFDiJetCorrCheckerWithDiTau);