Back to home page

Project CMSSW displayed by LXR

 
 

    


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

0001 #include "RecoTauTag/HLTProducers/interface/L1HLTTauMatching.h"
0002 #include "Math/GenVector/VectorUtil.h"
0003 #include "DataFormats/L1Trigger/interface/L1JetParticle.h"
0004 #include "DataFormats/L1Trigger/interface/L1JetParticleFwd.h"
0005 #include "DataFormats/HLTReco/interface/TriggerTypeDefs.h"
0006 #include "FWCore/Utilities/interface/EDMException.h"
0007 #include "DataFormats/TauReco/interface/PFTau.h"
0008 
0009 //
0010 // class decleration
0011 //
0012 using namespace reco;
0013 using namespace std;
0014 using namespace edm;
0015 using namespace l1extra;
0016 
0017 L1HLTTauMatching::L1HLTTauMatching(const edm::ParameterSet& iConfig)
0018     : jetSrc(consumes<PFTauCollection>(iConfig.getParameter<InputTag>("JetSrc"))),
0019       tauTrigger(consumes<trigger::TriggerFilterObjectWithRefs>(iConfig.getParameter<InputTag>("L1TauTrigger"))),
0020       mEt_Min(iConfig.getParameter<double>("EtMin")) {
0021   produces<PFTauCollection>();
0022 }
0023 L1HLTTauMatching::~L1HLTTauMatching() {}
0024 
0025 void L1HLTTauMatching::produce(edm::StreamID iSId, edm::Event& iEvent, const edm::EventSetup& iES) const {
0026   using namespace edm;
0027   using namespace std;
0028   using namespace reco;
0029   using namespace trigger;
0030   using namespace l1extra;
0031 
0032   unique_ptr<PFTauCollection> tauL2jets(new PFTauCollection);
0033 
0034   double deltaR = 1.0;
0035   double matchingR = 0.5;
0036   //Getting HLT jets to be matched
0037   edm::Handle<PFTauCollection> tauJets;
0038   iEvent.getByToken(jetSrc, tauJets);
0039 
0040   //        std::cout <<"Size of input jet collection "<<tauJets->size()<<std::endl;
0041 
0042   edm::Handle<trigger::TriggerFilterObjectWithRefs> l1TriggeredTaus;
0043   iEvent.getByToken(tauTrigger, l1TriggeredTaus);
0044 
0045   vector<l1extra::L1JetParticleRef> tauCandRefVec;
0046   vector<l1extra::L1JetParticleRef> jetCandRefVec;
0047   l1TriggeredTaus->getObjects(trigger::TriggerL1TauJet, tauCandRefVec);
0048   l1TriggeredTaus->getObjects(trigger::TriggerL1CenJet, jetCandRefVec);
0049 
0050   math::XYZPoint a(0., 0., 0.);
0051 
0052   for (unsigned int iL1Tau = 0; iL1Tau < tauCandRefVec.size(); iL1Tau++) {
0053     for (unsigned int iJet = 0; iJet < tauJets->size(); iJet++) {
0054       //Find the relative L2TauJets, to see if it has been reconstructed
0055       const PFTau& myJet = (*tauJets)[iJet];
0056       deltaR = ROOT::Math::VectorUtil::DeltaR(myJet.p4().Vect(), (tauCandRefVec[iL1Tau]->p4()).Vect());
0057       if (deltaR < matchingR) {
0058         //       LeafCandidate myLC(myJet);
0059         if (myJet.leadChargedHadrCand().isNonnull()) {
0060           a = myJet.leadChargedHadrCand()->vertex();
0061         }
0062         PFTau myPFTau(std::numeric_limits<int>::quiet_NaN(), myJet.p4(), a);
0063         if (myJet.pt() > mEt_Min) {
0064           //          tauL2LC->push_back(myLC);
0065           tauL2jets->push_back(myPFTau);
0066         }
0067         break;
0068       }
0069     }
0070   }
0071 
0072   for (unsigned int iL1Tau = 0; iL1Tau < jetCandRefVec.size(); iL1Tau++) {
0073     for (unsigned int iJet = 0; iJet < tauJets->size(); iJet++) {
0074       const PFTau& myJet = (*tauJets)[iJet];
0075       //Find the relative L2TauJets, to see if it has been reconstructed
0076       deltaR = ROOT::Math::VectorUtil::DeltaR(myJet.p4().Vect(), (jetCandRefVec[iL1Tau]->p4()).Vect());
0077       if (deltaR < matchingR) {
0078         //       LeafCandidate myLC(myJet);
0079         if (myJet.leadChargedHadrCand().isNonnull()) {
0080           a = myJet.leadChargedHadrCand()->vertex();
0081         }
0082 
0083         PFTau myPFTau(std::numeric_limits<int>::quiet_NaN(), myJet.p4(), a);
0084         if (myJet.pt() > mEt_Min) {
0085           //          tauL2LC->push_back(myLC);
0086           tauL2jets->push_back(myPFTau);
0087         }
0088         break;
0089       }
0090     }
0091   }
0092 
0093   //std::cout <<"Size of L1HLT matched jets "<<tauL2jets->size()<<std::endl;
0094 
0095   iEvent.put(std::move(tauL2jets));
0096   // iEvent.put(std::move(tauL2LC));
0097 }
0098 
0099 void L1HLTTauMatching::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0100   edm::ParameterSetDescription desc;
0101   desc.add<edm::InputTag>("L1TauTrigger", edm::InputTag("hltL1sDoubleIsoTau40er"))
0102       ->setComment("Name of trigger filter");
0103   desc.add<edm::InputTag>("JetSrc", edm::InputTag("hltSelectedPFTausTrackPt1MediumIsolationReg"))
0104       ->setComment("Input collection of PFTaus");
0105   desc.add<double>("EtMin", 0.0)->setComment("Minimal pT of PFTau to match");
0106   descriptions.setComment(
0107       "This module produces collection of PFTaus matched to L1ExtraTaus/Jets passing a HLT filter (Only p4 and vertex "
0108       "of returned PFTaus are set).");
0109   descriptions.add("L1HLTJetsMatching", desc);
0110 }
0111 
0112 #include "FWCore/Framework/interface/MakerMacros.h"
0113 DEFINE_FWK_MODULE(L1HLTTauMatching);