File indexing completed on 2023-05-06 02:46:38
0001
0002 #include "HLTPFTauPairLeadTrackDzMatchFilter.h"
0003
0004 #include "DataFormats/Common/interface/Handle.h"
0005 #include "DataFormats/Common/interface/RefToBase.h"
0006 #include "DataFormats/TrackReco/interface/Track.h"
0007 #include "DataFormats/TrackReco/interface/TrackFwd.h"
0008 #include "DataFormats/HLTReco/interface/TriggerFilterObjectWithRefs.h"
0009 #include "DataFormats/HLTReco/interface/TriggerTypeDefs.h"
0010 #include "DataFormats/Math/interface/deltaR.h"
0011
0012 HLTPFTauPairLeadTrackDzMatchFilter::HLTPFTauPairLeadTrackDzMatchFilter(const edm::ParameterSet& conf)
0013 : HLTFilter(conf) {
0014 tauSrc_ = conf.getParameter<edm::InputTag>("tauSrc");
0015 tauSrcToken_ = consumes<reco::PFTauCollection>(tauSrc_);
0016 tauMinPt_ = conf.getParameter<double>("tauMinPt");
0017 tauMaxEta_ = conf.getParameter<double>("tauMaxEta");
0018 tauMinDR_ = conf.getParameter<double>("tauMinDR");
0019 tauLeadTrackMaxDZ_ = conf.getParameter<double>("tauLeadTrackMaxDZ");
0020 triggerType_ = conf.getParameter<int>("triggerType");
0021
0022
0023
0024
0025 if (tauMinDR_ < 0.1)
0026 tauMinDR_ = 0.1;
0027 }
0028
0029 HLTPFTauPairLeadTrackDzMatchFilter::~HLTPFTauPairLeadTrackDzMatchFilter() {}
0030
0031 void HLTPFTauPairLeadTrackDzMatchFilter::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0032 edm::ParameterSetDescription desc;
0033 makeHLTFilterDescription(desc);
0034
0035 desc.add<edm::InputTag>("tauSrc", edm::InputTag("hltPFTaus"));
0036 desc.add<double>("tauMinPt", 5.0);
0037 desc.add<double>("tauMaxEta", 2.5);
0038 desc.add<double>("tauMinDR", 0.1);
0039 desc.add<double>("tauLeadTrackMaxDZ", 0.2);
0040 desc.add<int>("triggerType", trigger::TriggerTau);
0041 descriptions.add("hltPFTauPairLeadTrackDzMatchFilter", desc);
0042 }
0043
0044 bool HLTPFTauPairLeadTrackDzMatchFilter::hltFilter(edm::Event& ev,
0045 const edm::EventSetup& es,
0046 trigger::TriggerFilterObjectWithRefs& filterproduct) const {
0047 using namespace std;
0048 using namespace reco;
0049
0050
0051 if (saveTags())
0052 filterproduct.addCollectionTag(tauSrc_);
0053
0054
0055 PFTauRef ref;
0056
0057
0058 edm::Handle<PFTauCollection> tausHandle;
0059 ev.getByToken(tauSrcToken_, tausHandle);
0060 const PFTauCollection& taus = *tausHandle;
0061 const size_t n_taus = taus.size();
0062
0063
0064 size_t npairs = 0;
0065 if (n_taus > 1)
0066 for (size_t t1 = 0; t1 < n_taus; ++t1) {
0067 if (taus[t1].leadPFChargedHadrCand().isNull() || taus[t1].leadPFChargedHadrCand()->trackRef().isNull() ||
0068 taus[t1].pt() < tauMinPt_ || std::abs(taus[t1].eta()) > tauMaxEta_)
0069 continue;
0070
0071 float mindz = 99.f;
0072 for (size_t t2 = t1 + 1; t2 < n_taus; ++t2) {
0073 if (taus[t2].leadPFChargedHadrCand().isNull() || taus[t2].leadPFChargedHadrCand()->trackRef().isNull() ||
0074 taus[t2].pt() < tauMinPt_ || std::abs(taus[t2].eta()) > tauMaxEta_)
0075 continue;
0076
0077 float dr2 = reco::deltaR2(taus[t1].eta(), taus[t1].phi(), taus[t2].eta(), taus[t2].phi());
0078 float dz =
0079 (taus[t1].leadPFChargedHadrCand()->trackRef()->vz() - taus[t2].leadPFChargedHadrCand()->trackRef()->vz());
0080
0081
0082 if (dr2 < tauMinDR_ * tauMinDR_) {
0083 continue;
0084 }
0085
0086 if (std::abs(dz) < std::abs(mindz))
0087 mindz = dz;
0088
0089
0090 if (std::abs(dz) > tauLeadTrackMaxDZ_) {
0091 continue;
0092 }
0093
0094
0095 ref = PFTauRef(tausHandle, t1);
0096 filterproduct.addObject(triggerType_, ref);
0097
0098 ref = PFTauRef(tausHandle, t2);
0099 filterproduct.addObject(triggerType_, ref);
0100
0101 ++npairs;
0102 }
0103 }
0104
0105
0106 return npairs > 0;
0107 }
0108
0109 #include "FWCore/Framework/interface/MakerMacros.h"
0110 DEFINE_FWK_MODULE(HLTPFTauPairLeadTrackDzMatchFilter);