Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-09-16 00:05:54

0001 #include "RecoTauTag/HLTProducers/interface/L1THLTTauMatching.h"
0002 #include "Math/GenVector/VectorUtil.h"
0003 #include "DataFormats/HLTReco/interface/TriggerTypeDefs.h"
0004 #include "FWCore/Utilities/interface/EDMException.h"
0005 #include "DataFormats/TauReco/interface/PFTau.h"
0006 #include "DataFormats/Common/interface/RefVector.h"
0007 
0008 //
0009 // class decleration
0010 //
0011 using namespace reco;
0012 using namespace std;
0013 using namespace edm;
0014 using namespace trigger;
0015 
0016 L1THLTTauMatching::L1THLTTauMatching(const edm::ParameterSet& iConfig)
0017     : jetSrc(consumes<PFTauCollection>(iConfig.getParameter<InputTag>("JetSrc"))),
0018       tauTrigger(consumes<trigger::TriggerFilterObjectWithRefs>(iConfig.getParameter<InputTag>("L1TauTrigger"))),
0019       mEt_Min(iConfig.getParameter<double>("EtMin")),
0020       reduceTauContent(iConfig.getParameter<bool>("ReduceTauContent")),
0021       keepOriginalVertex(iConfig.getParameter<bool>("KeepOriginalVertex")) {
0022   produces<PFTauCollection>();
0023 }
0024 L1THLTTauMatching::~L1THLTTauMatching() {}
0025 
0026 void L1THLTTauMatching::produce(edm::StreamID iSId, edm::Event& iEvent, const edm::EventSetup& iES) const {
0027   unique_ptr<PFTauCollection> tauL2jets(new PFTauCollection);
0028 
0029   double deltaR = 1.0;
0030   double matchingR = 0.5;
0031 
0032   // Getting HLT jets to be matched
0033   edm::Handle<PFTauCollection> tauJets;
0034   iEvent.getByToken(jetSrc, tauJets);
0035 
0036   edm::Handle<trigger::TriggerFilterObjectWithRefs> l1TriggeredTaus;
0037   iEvent.getByToken(tauTrigger, l1TriggeredTaus);
0038 
0039   l1t::TauVectorRef tauCandRefVec;
0040   l1TriggeredTaus->getObjects(trigger::TriggerL1Tau, tauCandRefVec);
0041 
0042   math::XYZPoint a(0., 0., 0.);
0043 
0044   for (unsigned int iL1Tau = 0; iL1Tau < tauCandRefVec.size(); iL1Tau++) {
0045     for (unsigned int iJet = 0; iJet < tauJets->size(); iJet++) {
0046       // Find the relative L2TauJets, to see if it has been reconstructed
0047       const PFTau& myJet = (*tauJets)[iJet];
0048       deltaR = ROOT::Math::VectorUtil::DeltaR(myJet.p4().Vect(), (tauCandRefVec[iL1Tau]->p4()).Vect());
0049       if (deltaR < matchingR) {
0050         if (myJet.leadChargedHadrCand().isNonnull()) {
0051           a = myJet.leadChargedHadrCand()->vertex();
0052         }
0053 
0054         auto myPFTau =
0055             reduceTauContent ? PFTau(std::numeric_limits<int>::quiet_NaN(), myJet.p4(), myJet.vertex()) : PFTau(myJet);
0056 
0057         if (!keepOriginalVertex) {
0058           myPFTau.setVertex(a);
0059         }
0060 
0061         if (myPFTau.pt() > mEt_Min) {
0062           tauL2jets->push_back(myJet);
0063         }
0064         break;
0065       }
0066     }
0067   }
0068 
0069   iEvent.put(std::move(tauL2jets));
0070 }
0071 
0072 void L1THLTTauMatching::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0073   edm::ParameterSetDescription desc;
0074   desc.add<edm::InputTag>("L1TauTrigger", edm::InputTag("hltL1sDoubleIsoTau40er"))
0075       ->setComment("Name of trigger filter");
0076   desc.add<edm::InputTag>("JetSrc", edm::InputTag("hltSelectedPFTausTrackPt1MediumIsolationReg"))
0077       ->setComment("Input collection of PFTaus");
0078   desc.add<double>("EtMin", 0.0)->setComment("Minimal pT of PFTau to match");
0079   desc.add<bool>("ReduceTauContent", true)->setComment("Should produce taus with reduced content (Only p4 and vertex)");
0080   desc.add<bool>("KeepOriginalVertex", false)
0081       ->setComment("Should use original vertex instead of setting the vertex to that of the leading charged hadron");
0082   descriptions.setComment("This module produces collection of PFTaus matched to L1 Taus / Jets passing a HLT filter.");
0083   descriptions.add("L1THLTTauMatching", desc);
0084 }
0085 
0086 #include "FWCore/Framework/interface/MakerMacros.h"
0087 DEFINE_FWK_MODULE(L1THLTTauMatching);