Back to home page

Project CMSSW displayed by LXR

 
 

    


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

0001 #include "RecoTauTag/HLTProducers/interface/L1HLTJetsMatching.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 //
0008 // class decleration
0009 //
0010 using namespace reco;
0011 using namespace std;
0012 using namespace edm;
0013 using namespace l1extra;
0014 
0015 L1HLTJetsMatching::L1HLTJetsMatching(const edm::ParameterSet& iConfig) {
0016   jetSrc = consumes<edm::View<reco::Candidate> >(iConfig.getParameter<InputTag>("JetSrc"));
0017   tauTrigger = consumes<trigger::TriggerFilterObjectWithRefs>(iConfig.getParameter<InputTag>("L1TauTrigger"));
0018   mEt_Min = iConfig.getParameter<double>("EtMin");
0019 
0020   produces<CaloJetCollection>();
0021 }
0022 
0023 L1HLTJetsMatching::~L1HLTJetsMatching() {}
0024 
0025 void L1HLTJetsMatching::produce(edm::StreamID, 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   typedef std::vector<LeafCandidate> LeafCandidateCollection;
0033 
0034   unique_ptr<CaloJetCollection> tauL2jets(new CaloJetCollection);
0035 
0036   double deltaR = 1.0;
0037   double matchingR = 0.5;
0038   //Getting HLT jets to be matched
0039   edm::Handle<edm::View<Candidate> > tauJets;
0040   iEvent.getByToken(jetSrc, tauJets);
0041 
0042   //        std::cout <<"Size of input jet collection "<<tauJets->size()<<std::endl;
0043 
0044   edm::Handle<trigger::TriggerFilterObjectWithRefs> l1TriggeredTaus;
0045   iEvent.getByToken(tauTrigger, l1TriggeredTaus);
0046 
0047   std::vector<l1extra::L1JetParticleRef> tauCandRefVec;
0048   std::vector<l1extra::L1JetParticleRef> jetCandRefVec;
0049 
0050   l1TriggeredTaus->getObjects(trigger::TriggerL1TauJet, tauCandRefVec);
0051   l1TriggeredTaus->getObjects(trigger::TriggerL1CenJet, jetCandRefVec);
0052   math::XYZPoint a(0., 0., 0.);
0053   CaloJet::Specific f;
0054 
0055   for (unsigned int iL1Tau = 0; iL1Tau < tauCandRefVec.size(); iL1Tau++) {
0056     for (unsigned int iJet = 0; iJet < tauJets->size(); iJet++) {
0057       //Find the relative L2TauJets, to see if it has been reconstructed
0058       const Candidate& myJet = (*tauJets)[iJet];
0059       deltaR = ROOT::Math::VectorUtil::DeltaR(myJet.p4().Vect(), (tauCandRefVec[iL1Tau]->p4()).Vect());
0060       if (deltaR < matchingR) {
0061         //       LeafCandidate myLC(myJet);
0062         CaloJet myCaloJet(myJet.p4(), a, f);
0063         if (myJet.pt() > mEt_Min) {
0064           //          tauL2LC->push_back(myLC);
0065           tauL2jets->push_back(myCaloJet);
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 Candidate& 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         CaloJet myCaloJet(myJet.p4(), a, f);
0080         if (myJet.pt() > mEt_Min) {
0081           //tauL2LC->push_back(myLC);
0082           tauL2jets->push_back(myCaloJet);
0083         }
0084         break;
0085       }
0086     }
0087   }
0088 
0089   //std::cout <<"Size of L1HLT matched jets "<<tauL2jets->size()<<std::endl;
0090 
0091   iEvent.put(std::move(tauL2jets));
0092   // iEvent.put(std::move(tauL2LC));
0093 }
0094 
0095 #include "FWCore/Framework/interface/MakerMacros.h"
0096 DEFINE_FWK_MODULE(L1HLTJetsMatching);