File indexing completed on 2021-02-14 14:20:42
0001 #include "HLTPixelIsolTrackL1TFilter.h"
0002
0003 #include "DataFormats/Common/interface/Handle.h"
0004 #include "DataFormats/Common/interface/RefToBase.h"
0005
0006 #include "DataFormats/HLTReco/interface/TriggerTypeDefs.h"
0007
0008 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0009 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0010 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0011
0012
0013
0014 HLTPixelIsolTrackL1TFilter::HLTPixelIsolTrackL1TFilter(const edm::ParameterSet& iConfig) : HLTFilter(iConfig) {
0015 candTag_ = iConfig.getParameter<edm::InputTag>("candTag");
0016 hltGTseedlabel_ = iConfig.getParameter<edm::InputTag>("L1GTSeedLabel");
0017 minDeltaPtL1Jet_ = iConfig.getParameter<double>("MinDeltaPtL1Jet");
0018 minpttrack_ = iConfig.getParameter<double>("MinPtTrack");
0019 maxptnearby_ = iConfig.getParameter<double>("MaxPtNearby");
0020 maxetatrack_ = iConfig.getParameter<double>("MaxEtaTrack");
0021 minetatrack_ = iConfig.getParameter<double>("MinEtaTrack");
0022 filterE_ = iConfig.getParameter<bool>("filterTrackEnergy");
0023 minEnergy_ = iConfig.getParameter<double>("MinEnergyTrack");
0024 nMaxTrackCandidates_ = iConfig.getParameter<int>("NMaxTrackCandidates");
0025 dropMultiL2Event_ = iConfig.getParameter<bool>("DropMultiL2Event");
0026 candToken_ = consumes<reco::IsolatedPixelTrackCandidateCollection>(candTag_);
0027 hltGTseedToken_ = consumes<trigger::TriggerFilterObjectWithRefs>(hltGTseedlabel_);
0028 }
0029
0030 HLTPixelIsolTrackL1TFilter::~HLTPixelIsolTrackL1TFilter() = default;
0031
0032 void HLTPixelIsolTrackL1TFilter::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0033 edm::ParameterSetDescription desc;
0034 makeHLTFilterDescription(desc);
0035 desc.add<edm::InputTag>("candTag", edm::InputTag("hltIsolPixelTrackProd"));
0036 desc.add<edm::InputTag>("L1GTSeedLabel", edm::InputTag("hltL1sV0SingleJet60"));
0037 desc.add<double>("MaxPtNearby", 2.0);
0038 desc.add<double>("MinEnergyTrack", 12.0);
0039 desc.add<double>("MinPtTrack", 3.5);
0040 desc.add<double>("MaxEtaTrack", 1.15);
0041 desc.add<double>("MinEtaTrack", 0.0);
0042 desc.add<double>("MinDeltaPtL1Jet", -40000.0);
0043 desc.add<bool>("filterTrackEnergy", true);
0044 desc.add<int>("NMaxTrackCandidates", 10);
0045 desc.add<bool>("DropMultiL2Event", false);
0046 descriptions.add("hltPixelIsolTrackL1TFilter", desc);
0047 }
0048
0049 bool HLTPixelIsolTrackL1TFilter::hltFilter(edm::Event& iEvent,
0050 const edm::EventSetup& iSetup,
0051 trigger::TriggerFilterObjectWithRefs& filterproduct) const {
0052 if (saveTags())
0053 filterproduct.addCollectionTag(candTag_);
0054
0055
0056 edm::Ref<reco::IsolatedPixelTrackCandidateCollection> candref;
0057
0058
0059 edm::Handle<reco::IsolatedPixelTrackCandidateCollection> recotrackcands;
0060 iEvent.getByToken(candToken_, recotrackcands);
0061
0062
0063 int n = 0;
0064 for (unsigned int i = 0; i < recotrackcands->size(); i++) {
0065 candref = edm::Ref<reco::IsolatedPixelTrackCandidateCollection>(recotrackcands, i);
0066
0067
0068 if (!filterE_ && (candref->maxPtPxl() < maxptnearby_) && (candref->pt() > minpttrack_) &&
0069 fabs(candref->track()->eta()) < maxetatrack_ && fabs(candref->track()->eta()) > minetatrack_) {
0070 filterproduct.addObject(trigger::TriggerTrack, candref);
0071 n++;
0072 #ifdef DebugLog
0073 edm::LogInfo("HcalIsoTrack") << "PixelIsolP:Candidate[" << n << "] pt|eta|phi " << candref->pt() << "|"
0074 << candref->track()->pt() << "|" << candref->track()->eta() << "|"
0075 << candref->track()->phi() << "\n";
0076 #endif
0077 }
0078
0079
0080 if (filterE_) {
0081 if ((candref->maxPtPxl() < maxptnearby_) && ((candref->pt()) * cosh(candref->track()->eta()) > minEnergy_) &&
0082 fabs(candref->track()->eta()) < maxetatrack_ && fabs(candref->track()->eta()) > minetatrack_) {
0083 filterproduct.addObject(trigger::TriggerTrack, candref);
0084 n++;
0085 #ifdef DebugLog
0086 edm::LogInfo("HcalIsoTrack") << "PixelIsolE:Candidate[" << n << "] pt|eta|phi " << candref->pt() << "|"
0087 << candref->track()->pt() << "|" << candref->track()->eta() << "|"
0088 << candref->track()->phi() << "\n";
0089 #endif
0090 }
0091 }
0092
0093
0094 if (!dropMultiL2Event_ && n >= nMaxTrackCandidates_)
0095 break;
0096
0097 }
0098
0099 bool accept(n > 0);
0100
0101 if (dropMultiL2Event_ && n > nMaxTrackCandidates_)
0102 accept = false;
0103 #ifdef DebugLog
0104 edm::LogInfo("HcalIsoTrack") << "PixelIsolL1Filter: Tracks " << n << " accept " << accept << "\n";
0105 #endif
0106 return accept;
0107 }
0108
0109
0110 #include "FWCore/Framework/interface/MakerMacros.h"
0111 DEFINE_FWK_MODULE(HLTPixelIsolTrackL1TFilter);