File indexing completed on 2023-03-17 11:09:39
0001
0002
0003
0004
0005
0006
0007
0008
0009 #include "HLTMuonL1Filter.h"
0010
0011 #include "DataFormats/Common/interface/Handle.h"
0012 #include "DataFormats/HLTReco/interface/TriggerFilterObjectWithRefs.h"
0013 #include "DataFormats/HLTReco/interface/TriggerRefsCollections.h"
0014
0015 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0016
0017 #include "DataFormats/L1GlobalMuonTrigger/interface/L1MuGMTCand.h"
0018 #include "FWCore/Utilities/interface/EDMException.h"
0019
0020 #include "TMath.h"
0021
0022 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0023 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0024 #include "FWCore/Utilities/interface/InputTag.h"
0025
0026 #include <vector>
0027
0028
0029
0030
0031 HLTMuonL1Filter::HLTMuonL1Filter(const edm::ParameterSet& iConfig)
0032 : HLTFilter(iConfig),
0033 l1MuTriggerScalesRcdToken_(esConsumes()),
0034 candTag_(iConfig.getParameter<edm::InputTag>("CandTag")),
0035 candToken_(consumes<l1extra::L1MuonParticleCollection>(candTag_)),
0036 previousCandTag_(iConfig.getParameter<edm::InputTag>("PreviousCandTag")),
0037 previousCandToken_(consumes<trigger::TriggerFilterObjectWithRefs>(previousCandTag_)),
0038 maxEta_(iConfig.getParameter<double>("MaxEta")),
0039 minPt_(iConfig.getParameter<double>("MinPt")),
0040 minN_(iConfig.getParameter<int>("MinN")),
0041 excludeSingleSegmentCSC_(iConfig.getParameter<bool>("ExcludeSingleSegmentCSC")),
0042 csctfTag_(iConfig.getParameter<edm::InputTag>("CSCTFtag")),
0043 csctfToken_(excludeSingleSegmentCSC_ ? consumes<L1CSCTrackCollection>(csctfTag_)
0044 : edm::EDGetTokenT<L1CSCTrackCollection>{}) {
0045 using namespace std;
0046
0047
0048 qualityBitMask_ = 0;
0049 vector<int> selectQualities = iConfig.getParameter<vector<int> >("SelectQualities");
0050 for (int selectQualitie : selectQualities) {
0051 if (selectQualitie > 7) {
0052 throw edm::Exception(edm::errors::Configuration) << "QualityBits must be smaller than 8!";
0053 }
0054 qualityBitMask_ |= 1 << selectQualitie;
0055 }
0056
0057
0058 if (edm::isDebugEnabled()) {
0059 ostringstream ss;
0060 ss << "Constructed with parameters:" << endl;
0061 ss << " CandTag = " << candTag_.encode() << endl;
0062 ss << " PreviousCandTag = " << previousCandTag_.encode() << endl;
0063 ss << " MaxEta = " << maxEta_ << endl;
0064 ss << " MinPt = " << minPt_ << endl;
0065 ss << " SelectQualities =";
0066 for (size_t i = 0; i < 8; i++) {
0067 if ((qualityBitMask_ >> i) % 2)
0068 ss << " " << i;
0069 }
0070 ss << endl;
0071 ss << " MinN = " << minN_ << endl;
0072 ss << " ExcludeSingleSegmentCSC = " << excludeSingleSegmentCSC_ << endl;
0073 ss << " CSCTFtag = " << csctfTag_.encode() << endl;
0074 ss << " saveTags= " << saveTags();
0075 LogDebug("HLTMuonL1Filter") << ss.str();
0076 }
0077 }
0078
0079 HLTMuonL1Filter::~HLTMuonL1Filter() = default;
0080
0081 void HLTMuonL1Filter::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0082 edm::ParameterSetDescription desc;
0083 makeHLTFilterDescription(desc);
0084 desc.add<edm::InputTag>("CandTag", edm::InputTag("hltL1extraParticles"));
0085
0086 desc.add<edm::InputTag>("PreviousCandTag", edm::InputTag(""));
0087 desc.add<double>("MaxEta", 2.5);
0088 desc.add<double>("MinPt", 0.0);
0089 desc.add<int>("MinN", 1);
0090 desc.add<bool>("ExcludeSingleSegmentCSC", false);
0091
0092 desc.add<edm::InputTag>("CSCTFtag", edm::InputTag("csctfDigis"));
0093 {
0094 std::vector<int> temp1;
0095 temp1.reserve(0);
0096 desc.add<std::vector<int> >("SelectQualities", temp1);
0097 }
0098 descriptions.add("hltMuonL1Filter", desc);
0099 }
0100
0101
0102
0103
0104
0105
0106 bool HLTMuonL1Filter::hltFilter(edm::Event& iEvent,
0107 const edm::EventSetup& iSetup,
0108 trigger::TriggerFilterObjectWithRefs& filterproduct) const {
0109 using namespace std;
0110 using namespace edm;
0111 using namespace trigger;
0112 using namespace l1extra;
0113
0114
0115
0116
0117
0118
0119 Handle<L1MuonParticleCollection> allMuons;
0120 iEvent.getByToken(candToken_, allMuons);
0121
0122
0123 const L1CSCTrackCollection* csctfTracks = nullptr;
0124 const L1MuTriggerScales* l1MuTriggerScales = nullptr;
0125
0126
0127 if (excludeSingleSegmentCSC_) {
0128 edm::Handle<L1CSCTrackCollection> csctfTracksHandle;
0129 iEvent.getByToken(csctfToken_, csctfTracksHandle);
0130 csctfTracks = csctfTracksHandle.product();
0131
0132
0133 l1MuTriggerScales = &iSetup.getData(l1MuTriggerScalesRcdToken_);
0134 }
0135
0136
0137 Handle<TriggerFilterObjectWithRefs> previousLevelCands;
0138 iEvent.getByToken(previousCandToken_, previousLevelCands);
0139 vector<L1MuonParticleRef> prevMuons;
0140 previousLevelCands->getObjects(TriggerL1Mu, prevMuons);
0141
0142
0143 int n = 0;
0144 for (size_t i = 0; i < allMuons->size(); i++) {
0145 L1MuonParticleRef muon(allMuons, i);
0146
0147
0148 if (find(prevMuons.begin(), prevMuons.end(), muon) == prevMuons.end())
0149 continue;
0150
0151
0152 if (fabs(muon->eta()) > maxEta_)
0153 continue;
0154
0155
0156 if (muon->pt() < minPt_)
0157 continue;
0158
0159
0160 if (qualityBitMask_) {
0161 int quality = muon->gmtMuonCand().empty() ? 0 : (1 << muon->gmtMuonCand().quality());
0162 if ((quality & qualityBitMask_) == 0)
0163 continue;
0164 }
0165
0166
0167 if (excludeSingleSegmentCSC_ and isSingleSegmentCSC(muon, *csctfTracks, *l1MuTriggerScales))
0168 continue;
0169
0170
0171 n++;
0172 filterproduct.addObject(TriggerL1Mu, muon);
0173 }
0174
0175 if (saveTags())
0176 filterproduct.addCollectionTag(candTag_);
0177
0178
0179 const bool accept(n >= minN_);
0180
0181
0182 if (edm::isDebugEnabled()) {
0183 ostringstream ss;
0184 ss.precision(2);
0185 ss << "L1mu#" << '\t' << "q*pt" << '\t' << '\t' << "eta" << '\t' << "phi" << '\t' << "quality" << '\t' << "isPrev"
0186 << '\t' << "isFired" << '\t' << "isSingleCSC" << endl;
0187 ss << "--------------------------------------------------------------------------" << endl;
0188
0189 vector<L1MuonParticleRef> firedMuons;
0190 filterproduct.getObjects(TriggerL1Mu, firedMuons);
0191 for (size_t i = 0; i < allMuons->size(); i++) {
0192 L1MuonParticleRef mu(allMuons, i);
0193 int quality = mu->gmtMuonCand().empty() ? 0 : mu->gmtMuonCand().quality();
0194 bool isPrev = find(prevMuons.begin(), prevMuons.end(), mu) != prevMuons.end();
0195 bool isFired = find(firedMuons.begin(), firedMuons.end(), mu) != firedMuons.end();
0196 bool isSingleCSC = excludeSingleSegmentCSC_ and isSingleSegmentCSC(mu, *csctfTracks, *l1MuTriggerScales);
0197 ss << i << '\t' << scientific << mu->charge() * mu->pt() << '\t' << fixed << mu->eta() << '\t' << mu->phi()
0198 << '\t' << quality << '\t' << isPrev << '\t' << isFired << '\t' << isSingleCSC << endl;
0199 }
0200 ss << "--------------------------------------------------------------------------" << endl;
0201 LogDebug("HLTMuonL1Filter") << ss.str() << "Decision of filter is " << accept
0202 << ", number of muons passing = " << filterproduct.l1muonSize();
0203 }
0204
0205 return accept;
0206 }
0207
0208 bool HLTMuonL1Filter::isSingleSegmentCSC(l1extra::L1MuonParticleRef const& muon,
0209 L1CSCTrackCollection const& csctfTracks,
0210 L1MuTriggerScales const& l1MuTriggerScales) const {
0211
0212
0213
0214
0215
0216
0217
0218
0219 int csctfMode = -999;
0220
0221
0222 for (auto trk = csctfTracks.begin(); trk < csctfTracks.end(); trk++) {
0223 int trEndcap = (trk->first.endcap() == 2 ? trk->first.endcap() - 3 : trk->first.endcap());
0224 int trSector = 6 * (trk->first.endcap() - 1) + trk->first.sector();
0225
0226
0227
0228 float trEtaScale = l1MuTriggerScales.getRegionalEtaScale(2)->getCenter(trk->first.eta_packed());
0229 float trPhiScale = l1MuTriggerScales.getPhiScale()->getLowEdge(trk->first.localPhi());
0230
0231 double trEta = trEtaScale * trEndcap;
0232
0233 if (trEta < -2.4)
0234 trEta = -2.375;
0235
0236
0237
0238
0239 float trPhi02PI = fmod(trPhiScale + ((trSector - 1) * TMath::Pi() / 3) + (TMath::Pi() / 12), 2 * TMath::Pi());
0240
0241
0242 double trPhi = (trPhi02PI < TMath::Pi() ? trPhi02PI : trPhi02PI - 2 * TMath::Pi());
0243
0244
0245
0246
0247
0248
0249
0250
0251
0252
0253 if (fabs(trEta - muon->eta()) < 0.03 && fabs(trPhi - muon->phi()) < 0.001) {
0254
0255 ptadd thePtAddress(trk->first.ptLUTAddress());
0256 csctfMode = thePtAddress.track_mode;
0257
0258 }
0259 }
0260
0261
0262
0263
0264
0265
0266
0267
0268
0269 return csctfMode == 11;
0270 }
0271
0272
0273 #include "FWCore/Framework/interface/MakerMacros.h"
0274 DEFINE_FWK_MODULE(HLTMuonL1Filter);