File indexing completed on 2024-04-06 12:18:36
0001
0002
0003
0004
0005
0006
0007
0008
0009 #include "HLTMuonL1toL3TkPreFilter.h"
0010
0011 #include "DataFormats/Common/interface/Handle.h"
0012 #include "DataFormats/Common/interface/RefToBase.h"
0013
0014 #include "DataFormats/HLTReco/interface/TriggerFilterObjectWithRefs.h"
0015 #include "DataFormats/HLTReco/interface/TriggerRefsCollections.h"
0016
0017 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0018 #include "DataFormats/TrackReco/interface/Track.h"
0019 #include "DataFormats/RecoCandidate/interface/RecoChargedCandidate.h"
0020 #include "DataFormats/RecoCandidate/interface/RecoChargedCandidateFwd.h"
0021 #include "DataFormats/MuonReco/interface/MuonTrackLinks.h"
0022 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
0023 #include "DataFormats/MuonSeed/interface/L3MuonTrajectorySeed.h"
0024 #include "DataFormats/MuonSeed/interface/L3MuonTrajectorySeedCollection.h"
0025
0026 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0027
0028
0029
0030
0031 using namespace std;
0032 using namespace edm;
0033 using namespace reco;
0034 using namespace trigger;
0035
0036 HLTMuonL1toL3TkPreFilter::HLTMuonL1toL3TkPreFilter(const ParameterSet& iConfig)
0037 : HLTFilter(iConfig),
0038 beamspotTag_(iConfig.getParameter<edm::InputTag>("BeamSpotTag")),
0039 beamspotToken_(consumes<reco::BeamSpot>(beamspotTag_)),
0040 candTag_(iConfig.getParameter<InputTag>("CandTag")),
0041 candToken_(consumes<reco::RecoChargedCandidateCollection>(candTag_)),
0042 previousCandTag_(iConfig.getParameter<InputTag>("PreviousCandTag")),
0043 previousCandToken_(consumes<trigger::TriggerFilterObjectWithRefs>(previousCandTag_)),
0044 min_N_(iConfig.getParameter<int>("MinN")),
0045 max_Eta_(iConfig.getParameter<double>("MaxEta")),
0046 min_Nhits_(iConfig.getParameter<int>("MinNhits")),
0047 max_Dr_(iConfig.getParameter<double>("MaxDr")),
0048 max_Dz_(iConfig.getParameter<double>("MaxDz")),
0049 min_Pt_(iConfig.getParameter<double>("MinPt")),
0050 nsigma_Pt_(iConfig.getParameter<double>("NSigmaPt")) {
0051 LogDebug("HLTMuonL1toL3TkPreFilter") << " CandTag/MinN/MaxEta/MinNhits/MaxDr/MaxDz/MinPt/NSigmaPt : "
0052 << candTag_.encode() << " " << min_N_ << " " << max_Eta_ << " " << min_Nhits_
0053 << " " << max_Dr_ << " " << max_Dz_ << " " << min_Pt_ << " " << nsigma_Pt_;
0054
0055
0056 produces<TriggerFilterObjectWithRefs>();
0057 }
0058
0059 HLTMuonL1toL3TkPreFilter::~HLTMuonL1toL3TkPreFilter() = default;
0060
0061
0062
0063
0064 void HLTMuonL1toL3TkPreFilter::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0065 edm::ParameterSetDescription desc;
0066 makeHLTFilterDescription(desc);
0067 desc.add<edm::InputTag>("BeamSpotTag", edm::InputTag("hltBeamSpotTag"));
0068 desc.add<edm::InputTag>("CandTag", edm::InputTag("hltCandTag"));
0069 desc.add<edm::InputTag>("PreviousCandTag", edm::InputTag("hltPreviousCandTag"));
0070 desc.add<int>("MinN", 0);
0071 desc.add<double>("MaxEta", 9999.0);
0072 desc.add<int>("MinNhits", 0);
0073 desc.add<double>("MaxDr", 9999.0);
0074 desc.add<double>("MaxDz", 9999.0);
0075 desc.add<double>("MinPt", 0.0);
0076 desc.add<double>("NSigmaPt", 9999.0);
0077 descriptions.add("hltMuonL1toL3TkPreFilter", desc);
0078 }
0079
0080
0081 bool HLTMuonL1toL3TkPreFilter::hltFilter(Event& iEvent,
0082 const EventSetup& iSetup,
0083 trigger::TriggerFilterObjectWithRefs& filterproduct) const {
0084
0085
0086
0087
0088
0089 Handle<RecoChargedCandidateCollection> mucands;
0090 iEvent.getByToken(candToken_, mucands);
0091 if (saveTags())
0092 filterproduct.addCollectionTag(candTag_);
0093
0094 std::map<l1extra::L1MuonParticleRef, std::vector<RecoChargedCandidateRef> > L1toL3s;
0095 unsigned int n = 0;
0096 unsigned int maxN = mucands->size();
0097 for (; n != maxN; n++) {
0098 TrackRef tk = (*mucands)[n].track();
0099 edm::Ref<L3MuonTrajectorySeedCollection> l3seedRef =
0100 tk->seedRef().castTo<edm::Ref<L3MuonTrajectorySeedCollection> >();
0101 l1extra::L1MuonParticleRef l1mu = l3seedRef->l1Particle();
0102 L1toL3s[l1mu].push_back(RecoChargedCandidateRef(mucands, n));
0103 }
0104
0105
0106 Handle<TriggerFilterObjectWithRefs> previousLevelCands;
0107 iEvent.getByToken(previousCandToken_, previousLevelCands);
0108 BeamSpot beamSpot;
0109 Handle<BeamSpot> recoBeamSpotHandle;
0110 iEvent.getByToken(beamspotToken_, recoBeamSpotHandle);
0111 beamSpot = *recoBeamSpotHandle;
0112
0113
0114 vector<l1extra::L1MuonParticleRef> vl1cands;
0115 previousLevelCands->getObjects(TriggerL1Mu, vl1cands);
0116
0117 auto L1toL3s_it = L1toL3s.begin();
0118 auto L1toL3s_end = L1toL3s.end();
0119 for (; L1toL3s_it != L1toL3s_end; ++L1toL3s_it) {
0120 if (!triggeredAtL1(L1toL3s_it->first, vl1cands))
0121 continue;
0122
0123
0124 unsigned int iTk = 0;
0125 unsigned int maxItk = L1toL3s_it->second.size();
0126 for (; iTk != maxItk; iTk++) {
0127 RecoChargedCandidateRef& cand = L1toL3s_it->second[iTk];
0128 TrackRef tk = cand->track();
0129
0130 if (fabs(tk->eta()) > max_Eta_)
0131 continue;
0132
0133
0134 if (tk->numberOfValidHits() < min_Nhits_)
0135 continue;
0136
0137
0138
0139 if (fabs(tk->dxy(beamSpot.position())) > max_Dr_)
0140 continue;
0141
0142
0143 if (fabs(tk->dz()) > max_Dz_)
0144 continue;
0145
0146
0147 double pt = tk->pt();
0148 double err0 = tk->error(0);
0149 double abspar0 = fabs(tk->parameter(0));
0150 double ptLx = pt;
0151
0152 if (abspar0 > 0)
0153 ptLx += nsigma_Pt_ * err0 / abspar0 * pt;
0154 LogTrace("HLTMuonL1toL3TkPreFilter") << " ...Muon in loop, pt= " << pt << ", ptLx= " << ptLx;
0155 if (ptLx < min_Pt_)
0156 continue;
0157
0158
0159 filterproduct.addObject(TriggerMuon, cand);
0160 break;
0161 }
0162
0163 }
0164
0165 vector<RecoChargedCandidateRef> vref;
0166 filterproduct.getObjects(TriggerMuon, vref);
0167 for (auto& i : vref) {
0168 TrackRef tk = i->track();
0169 LogDebug("HLTMuonL1toL3TkPreFilter") << " Track passing filter: pt= " << tk->pt() << ", eta: " << tk->eta();
0170 }
0171
0172
0173 const bool accept((int)n >= min_N_);
0174
0175 LogDebug("HLTMuonL1toL3TkPreFilter") << " >>>>> Result of HLTMuonL1toL3TkPreFilter is " << accept
0176 << ", number of muons passing thresholds= " << n;
0177
0178 return accept;
0179 }
0180
0181 bool HLTMuonL1toL3TkPreFilter::triggeredAtL1(const l1extra::L1MuonParticleRef& l1mu,
0182 std::vector<l1extra::L1MuonParticleRef>& vcands) const {
0183 bool ok = false;
0184
0185
0186 for (auto& vcand : vcands) {
0187
0188 if (vcand == l1mu) {
0189 ok = true;
0190 LogDebug("HLTMuonL1toL3TkPreFilter") << "The L1 mu triggered";
0191 break;
0192 }
0193 }
0194 return ok;
0195 }
0196
0197
0198 #include "FWCore/Framework/interface/MakerMacros.h"
0199 DEFINE_FWK_MODULE(HLTMuonL1toL3TkPreFilter);