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