File indexing completed on 2024-04-06 12:18:37
0001
0002
0003
0004
0005
0006
0007
0008 #include "HLTMuonL3SimplePreFilter.h"
0009
0010 #include "DataFormats/Common/interface/Handle.h"
0011 #include "DataFormats/Common/interface/RefToBase.h"
0012 #include "DataFormats/TrackReco/interface/Track.h"
0013 #include "DataFormats/HLTReco/interface/TriggerFilterObjectWithRefs.h"
0014 #include "DataFormats/HLTReco/interface/TriggerRefsCollections.h"
0015 #include "DataFormats/RecoCandidate/interface/RecoChargedCandidate.h"
0016 #include "DataFormats/RecoCandidate/interface/RecoChargedCandidateFwd.h"
0017
0018 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0019
0020 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0021 #include <iostream>
0022
0023
0024
0025
0026 using namespace std;
0027 using namespace edm;
0028 using namespace trigger;
0029 using namespace reco;
0030
0031 HLTMuonL3SimplePreFilter::HLTMuonL3SimplePreFilter(const edm::ParameterSet& iConfig)
0032 : HLTFilter(iConfig),
0033 candTag_(iConfig.getParameter<edm::InputTag>("CandTag")),
0034 previousCandTag_(iConfig.getParameter<edm::InputTag>("PreviousCandTag")),
0035 beamspotTag_(iConfig.getParameter<edm::InputTag>("BeamSpotTag")),
0036 min_N_(iConfig.getParameter<int>("MinN")),
0037 max_Eta_(iConfig.getParameter<double>("MaxEta")),
0038 min_Nhits_(iConfig.getParameter<int>("MinNhits")),
0039 max_Dz_(iConfig.getParameter<double>("MaxDz")),
0040 min_DxySig_(iConfig.getParameter<double>("MinDxySig")),
0041 min_Pt_(iConfig.getParameter<double>("MinPt")),
0042 nsigma_Pt_(iConfig.getParameter<double>("NSigmaPt")),
0043 max_NormalizedChi2_(iConfig.getParameter<double>("MaxNormalizedChi2")),
0044 max_DXYBeamSpot_(iConfig.getParameter<double>("MaxDXYBeamSpot")),
0045 min_DXYBeamSpot_(iConfig.getParameter<double>("MinDXYBeamSpot")),
0046 min_NmuonHits_(iConfig.getParameter<int>("MinNmuonHits")),
0047 max_PtDifference_(iConfig.getParameter<double>("MaxPtDifference")),
0048 min_TrackPt_(iConfig.getParameter<double>("MinTrackPt")),
0049 matchPreviousCand_(iConfig.getParameter<bool>("MatchToPreviousCand")) {
0050 candToken_ = consumes<reco::RecoChargedCandidateCollection>(candTag_);
0051 previousCandToken_ = consumes<trigger::TriggerFilterObjectWithRefs>(previousCandTag_);
0052 beamspotToken_ = consumes<reco::BeamSpot>(beamspotTag_);
0053 }
0054
0055 HLTMuonL3SimplePreFilter::~HLTMuonL3SimplePreFilter() = default;
0056
0057
0058
0059
0060 void HLTMuonL3SimplePreFilter::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0061 edm::ParameterSetDescription desc;
0062 makeHLTFilterDescription(desc);
0063 desc.add<edm::InputTag>("BeamSpotTag", edm::InputTag("hltOfflineBeamSpot"));
0064 desc.add<edm::InputTag>("CandTag", edm::InputTag("hltL3MuonCandidates"));
0065 desc.add<edm::InputTag>("PreviousCandTag", edm::InputTag(""));
0066 desc.add<int>("MinN", 1);
0067 desc.add<double>("MaxEta", 2.5);
0068 desc.add<int>("MinNhits", 0);
0069 desc.add<double>("MaxDz", 9999.0);
0070 desc.add<double>("MinDxySig", -1.0);
0071 desc.add<double>("MinPt", 3.0);
0072 desc.add<double>("NSigmaPt", 0.0);
0073 desc.add<double>("MaxNormalizedChi2", 9999.0);
0074 desc.add<double>("MaxDXYBeamSpot", 9999.0);
0075 desc.add<double>("MinDXYBeamSpot", -1.0);
0076 desc.add<int>("MinNmuonHits", 0);
0077 desc.add<double>("MaxPtDifference", 9999.0);
0078 desc.add<double>("MinTrackPt", 0.0);
0079 desc.add<bool>("MatchToPreviousCand", true);
0080 descriptions.add("hltMuonL3SimplePreFilter", desc);
0081 }
0082
0083
0084 bool HLTMuonL3SimplePreFilter::hltFilter(edm::Event& iEvent,
0085 const edm::EventSetup& iSetup,
0086 trigger::TriggerFilterObjectWithRefs& filterproduct) const {
0087
0088
0089
0090
0091 if (saveTags())
0092 filterproduct.addCollectionTag(candTag_);
0093
0094
0095 Handle<RecoChargedCandidateCollection> mucands;
0096 iEvent.getByToken(candToken_, mucands);
0097 Handle<TriggerFilterObjectWithRefs> previousLevelCands;
0098 iEvent.getByToken(previousCandToken_, previousLevelCands);
0099 vector<RecoChargedCandidateRef> vcands;
0100 if (previousLevelCands.isValid()) {
0101 previousLevelCands->getObjects(TriggerMuon, vcands);
0102 }
0103
0104 Handle<BeamSpot> recoBeamSpotHandle;
0105 iEvent.getByToken(beamspotToken_, recoBeamSpotHandle);
0106
0107
0108 int n = 0;
0109 for (unsigned int iMu = 0; iMu < mucands->size(); iMu++) {
0110 RecoChargedCandidateRef cand(mucands, iMu);
0111 LogDebug("HLTMuonL3SimplePreFilter") << "cand isNonnull " << cand.isNonnull();
0112
0113
0114 if (matchPreviousCand_ && !triggerdByPreviousLevel(cand, vcands))
0115 continue;
0116
0117 if (std::abs(cand->eta()) > max_Eta_)
0118 continue;
0119
0120 TrackRef tk = cand->track();
0121 LogDebug("HLTMuonL3SimplePreFilter") << " Muon in loop, q*pt= " << tk->charge() * tk->pt() << " ("
0122 << cand->charge() * cand->pt() << ") "
0123 << ", eta= " << tk->eta() << " (" << cand->eta() << ") "
0124 << ", hits= " << tk->numberOfValidHits() << ", d0= " << tk->d0()
0125 << ", dz= " << tk->dz();
0126
0127
0128 if (tk->numberOfValidHits() < min_Nhits_)
0129 continue;
0130
0131
0132 if (tk->normalizedChi2() > max_NormalizedChi2_)
0133 continue;
0134
0135 if (recoBeamSpotHandle.isValid()) {
0136 const BeamSpot& beamSpot = *recoBeamSpotHandle;
0137
0138
0139 if (std::abs((cand->vz() - beamSpot.z0()) -
0140 ((cand->vx() - beamSpot.x0()) * cand->px() + (cand->vy() - beamSpot.y0()) * cand->py()) /
0141 cand->pt() * cand->pz() / cand->pt()) > max_Dz_)
0142 continue;
0143
0144
0145 if (min_DxySig_ > 0 &&
0146 (tk->dxyError() <= 0 || std::abs(tk->dxy(beamSpot.position()) / tk->dxyError()) < min_DxySig_))
0147 continue;
0148
0149
0150 float absDxy = std::abs(tk->dxy(beamSpot.position()));
0151 if (absDxy > max_DXYBeamSpot_ || absDxy < min_DXYBeamSpot_)
0152 continue;
0153 }
0154
0155
0156 const reco::HitPattern& trackHits = tk->hitPattern();
0157 if (trackHits.numberOfValidMuonHits() < min_NmuonHits_)
0158 continue;
0159
0160
0161 double candPt = cand->pt();
0162 double trackPt = tk->pt();
0163
0164 if (std::abs(candPt - trackPt) > max_PtDifference_)
0165 continue;
0166
0167
0168 if (trackPt < min_TrackPt_)
0169 continue;
0170
0171
0172 double pt = cand->pt();
0173 double err0 = tk->error(0);
0174 double abspar0 = std::abs(tk->parameter(0));
0175 double ptLx = pt;
0176
0177 if (abspar0 > 0)
0178 ptLx += nsigma_Pt_ * err0 / abspar0 * pt;
0179 LogTrace("HLTMuonL3SimplePreFilter") << " ...Muon in loop, trackkRef pt= " << tk->pt() << ", ptLx= " << ptLx
0180 << " cand pT " << cand->pt();
0181 if (ptLx < min_Pt_)
0182 continue;
0183
0184 n++;
0185 filterproduct.addObject(TriggerMuon, cand);
0186 }
0187
0188
0189 const bool accept(n >= min_N_);
0190
0191 LogDebug("HLTMuonL3SimplePreFilter") << " >>>>> Result of HLTMuonL3PreFilter is " << accept
0192 << ", number of muons passing thresholds= " << n;
0193
0194 return accept;
0195 }
0196
0197 bool HLTMuonL3SimplePreFilter::triggerdByPreviousLevel(const reco::RecoChargedCandidateRef& candref,
0198 const std::vector<reco::RecoChargedCandidateRef>& vcands) {
0199 unsigned int i = 0;
0200 unsigned int i_max = vcands.size();
0201 for (; i != i_max; ++i) {
0202 if (candref == vcands[i])
0203 return true;
0204 }
0205
0206 return false;
0207 }
0208
0209
0210 #include "FWCore/Framework/interface/MakerMacros.h"
0211 DEFINE_FWK_MODULE(HLTMuonL3SimplePreFilter);