File indexing completed on 2023-03-17 11:09:39
0001
0002 #include "HLTMuonL1TRegionalFilter.h"
0003 #include "DataFormats/Common/interface/Handle.h"
0004 #include "DataFormats/HLTReco/interface/TriggerFilterObjectWithRefs.h"
0005 #include "DataFormats/HLTReco/interface/TriggerRefsCollections.h"
0006 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0007 #include "FWCore/Utilities/interface/EDMException.h"
0008 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0009
0010 HLTMuonL1TRegionalFilter::HLTMuonL1TRegionalFilter(const edm::ParameterSet& iConfig)
0011 : HLTFilter(iConfig),
0012 candTag_(iConfig.getParameter<edm::InputTag>("CandTag")),
0013 candToken_(consumes<l1t::MuonBxCollection>(candTag_)),
0014 previousCandTag_(iConfig.getParameter<edm::InputTag>("PreviousCandTag")),
0015 previousCandToken_(consumes<trigger::TriggerFilterObjectWithRefs>(previousCandTag_)),
0016 minN_(iConfig.getParameter<int>("MinN")),
0017 centralBxOnly_(iConfig.getParameter<bool>("CentralBxOnly")) {
0018 using namespace std;
0019 using namespace edm;
0020
0021
0022 const vector<ParameterSet> cuts = iConfig.getParameter<vector<ParameterSet> >("Cuts");
0023 size_t ranges = cuts.size();
0024 if (ranges == 0) {
0025 throw edm::Exception(errors::Configuration) << "Please provide at least one PSet in the Cuts VPSet!";
0026 }
0027 etaBoundaries_.reserve(ranges + 1);
0028 minPts_.reserve(ranges);
0029 qualityBitMasks_.reserve(ranges);
0030 for (size_t i = 0; i < ranges; i++) {
0031
0032 vector<double> etaRange = cuts[i].getParameter<vector<double> >("EtaRange");
0033 if (etaRange.size() != 2 || etaRange[0] >= etaRange[1]) {
0034 throw edm::Exception(errors::Configuration) << "EtaRange must have two non-equal values in increasing order!";
0035 }
0036 if (i == 0) {
0037 etaBoundaries_.push_back(etaRange[0]);
0038 } else if (etaBoundaries_[i] != etaRange[0]) {
0039 throw edm::Exception(errors::Configuration)
0040 << "EtaRanges must be disjoint without gaps and listed in increasing eta order!";
0041 }
0042 etaBoundaries_.push_back(etaRange[1]);
0043
0044
0045 minPts_.push_back(cuts[i].getParameter<double>("MinPt"));
0046
0047
0048 qualityBitMasks_.push_back(0);
0049 vector<unsigned int> qualities = cuts[i].getParameter<vector<unsigned int> >("QualityBits");
0050 for (unsigned int qualitie : qualities) {
0051
0052
0053
0054 qualityBitMasks_[i] |= 1 << qualitie;
0055 }
0056 }
0057
0058
0059 if (edm::isDebugEnabled()) {
0060 ostringstream ss;
0061 ss << "Constructed with parameters:" << endl;
0062 ss << " CandTag = " << candTag_.encode() << endl;
0063 ss << " PreviousCandTag = " << previousCandTag_.encode() << endl;
0064 ss << " EtaBoundaries = \t" << etaBoundaries_[0];
0065 for (size_t i = 1; i < etaBoundaries_.size(); i++) {
0066 ss << '\t' << etaBoundaries_[i];
0067 }
0068 ss << endl;
0069 ss << " MinPts = \t " << minPts_[0];
0070 for (size_t i = 1; i < minPts_.size(); i++) {
0071 ss << "\t " << minPts_[i];
0072 }
0073 ss << endl;
0074 ss << " QualityBitMasks = \t " << qualityBitMasks_[0];
0075 for (size_t i = 1; i < qualityBitMasks_.size(); i++) {
0076 ss << "\t " << qualityBitMasks_[i];
0077 }
0078 ss << endl;
0079 ss << " MinN = " << minN_ << endl;
0080 ss << " saveTags= " << saveTags();
0081 LogDebug("HLTMuonL1TRegionalFilter") << ss.str();
0082 }
0083 }
0084
0085 HLTMuonL1TRegionalFilter::~HLTMuonL1TRegionalFilter() = default;
0086
0087 void HLTMuonL1TRegionalFilter::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0088 edm::ParameterSetDescription desc;
0089 makeHLTFilterDescription(desc);
0090 desc.add<edm::InputTag>("CandTag", edm::InputTag("hltGmtStage2Digis"));
0091 desc.add<edm::InputTag>("PreviousCandTag", edm::InputTag("hltL1sL1SingleMu20"));
0092 desc.add<int>("MinN", 1);
0093 desc.add<bool>("CentralBxOnly", true);
0094
0095 edm::ParameterSetDescription validator;
0096 std::vector<edm::ParameterSet> defaults(3);
0097
0098 std::vector<double> etaRange;
0099 double minPt;
0100 std::vector<unsigned int> qualityBits;
0101
0102 etaRange.clear();
0103 etaRange.push_back(-2.5);
0104 etaRange.push_back(+2.5);
0105 minPt = 20.0;
0106 qualityBits.clear();
0107 qualityBits.push_back(6);
0108 qualityBits.push_back(7);
0109 validator.add<std::vector<double> >("EtaRange", etaRange);
0110 validator.add<double>("MinPt", minPt);
0111 validator.add<std::vector<unsigned int> >("QualityBits", qualityBits);
0112
0113 etaRange.clear();
0114 etaRange.push_back(-2.5);
0115 etaRange.push_back(-1.6);
0116 minPt = 20.0;
0117 qualityBits.clear();
0118 qualityBits.push_back(6);
0119 qualityBits.push_back(7);
0120 defaults[0].addParameter<std::vector<double> >("EtaRange", etaRange);
0121 defaults[0].addParameter<double>("MinPt", minPt);
0122 defaults[0].addParameter<std::vector<unsigned int> >("QualityBits", qualityBits);
0123
0124 etaRange.clear();
0125 etaRange.push_back(-1.6);
0126 etaRange.push_back(+1.6);
0127 minPt = 20.0;
0128 qualityBits.clear();
0129 qualityBits.push_back(7);
0130 defaults[1].addParameter<std::vector<double> >("EtaRange", etaRange);
0131 defaults[1].addParameter<double>("MinPt", minPt);
0132 defaults[1].addParameter<std::vector<unsigned int> >("QualityBits", qualityBits);
0133
0134 etaRange.clear();
0135 etaRange.push_back(+1.6);
0136 etaRange.push_back(+2.5);
0137 minPt = 20.0;
0138 qualityBits.clear();
0139 qualityBits.push_back(6);
0140 qualityBits.push_back(7);
0141 edm::ParameterSetDescription element2;
0142 defaults[2].addParameter<std::vector<double> >("EtaRange", etaRange);
0143 defaults[2].addParameter<double>("MinPt", minPt);
0144 defaults[2].addParameter<std::vector<unsigned int> >("QualityBits", qualityBits);
0145
0146 desc.addVPSet("Cuts", validator, defaults);
0147
0148 descriptions.add("hltMuonL1TRegionalFilter", desc);
0149 }
0150
0151 bool HLTMuonL1TRegionalFilter::hltFilter(edm::Event& iEvent,
0152 const edm::EventSetup& iSetup,
0153 trigger::TriggerFilterObjectWithRefs& filterproduct) const {
0154 using namespace std;
0155 using namespace edm;
0156 using namespace trigger;
0157 using namespace l1t;
0158
0159
0160
0161
0162
0163
0164 Handle<l1t::MuonBxCollection> allMuons;
0165 iEvent.getByToken(candToken_, allMuons);
0166
0167
0168 Handle<TriggerFilterObjectWithRefs> previousLevelCands;
0169 iEvent.getByToken(previousCandToken_, previousLevelCands);
0170
0171 vector<MuonRef> prevMuons;
0172 previousLevelCands->getObjects(TriggerL1Mu, prevMuons);
0173
0174
0175 int n = 0;
0176
0177 for (int ibx = allMuons->getFirstBX(); ibx <= allMuons->getLastBX(); ++ibx) {
0178 if (centralBxOnly_ && (ibx != 0))
0179 continue;
0180 for (auto it = allMuons->begin(ibx); it != allMuons->end(ibx); it++) {
0181 MuonRef muon(allMuons, distance(allMuons->begin(allMuons->getFirstBX()), it));
0182
0183
0184 if (find(prevMuons.begin(), prevMuons.end(), muon) == prevMuons.end())
0185 continue;
0186
0187
0188 float eta = muon->eta();
0189 int region = -1;
0190 for (size_t r = 0; r < etaBoundaries_.size() - 1; r++) {
0191 if (etaBoundaries_[r] <= eta && eta <= etaBoundaries_[r + 1]) {
0192 region = r;
0193 break;
0194 }
0195 }
0196 if (region == -1)
0197 continue;
0198
0199
0200 if (muon->pt() < minPts_[region])
0201 continue;
0202
0203
0204 if (qualityBitMasks_[region]) {
0205 int quality = (it->hwQual() == 0 ? 0 : (1 << it->hwQual()));
0206 if ((quality & qualityBitMasks_[region]) == 0)
0207 continue;
0208 }
0209
0210
0211 n++;
0212 filterproduct.addObject(TriggerL1Mu, muon);
0213 }
0214 }
0215
0216 if (saveTags())
0217 filterproduct.addCollectionTag(candTag_);
0218
0219
0220 const bool accept(n >= minN_);
0221
0222
0223 if (edm::isDebugEnabled()) {
0224 LogTrace("HLTMuonL1TRegionalFilter") << "\nHLTMuonL1TRegionalFilter -----------------------------------------------"
0225 << endl;
0226 LogTrace("HLTMuonL1TRegionalFilter") << "L1mu#" << '\t' << "q*pt" << '\t' << '\t' << "eta" << '\t' << "phi" << '\t'
0227 << "quality" << '\t' << "isPrev\t " << endl;
0228 LogTrace("HLTMuonL1TRegionalFilter") << "--------------------------------------------------------------------------"
0229 << endl;
0230
0231 vector<MuonRef> firedMuons;
0232 filterproduct.getObjects(TriggerL1Mu, firedMuons);
0233 for (size_t i = 0; i < firedMuons.size(); i++) {
0234 l1t::MuonRef mu = firedMuons[i];
0235 bool isPrev = find(prevMuons.begin(), prevMuons.end(), mu) != prevMuons.end();
0236 LogTrace("HLTMuonL1TRegionalFilter")
0237 << i << '\t' << setprecision(2) << scientific << mu->charge() * mu->pt() << '\t' << fixed << mu->eta() << '\t'
0238 << mu->phi() << '\t' << mu->hwQual() << '\t' << isPrev << endl;
0239 }
0240 LogTrace("HLTMuonL1TRegionalFilter") << "--------------------------------------------------------------------------"
0241 << endl;
0242 LogTrace("HLTMuonL1TRegionalFilter") << "Decision of this filter is " << accept
0243 << ", number of muons passing = " << filterproduct.l1tmuonSize();
0244 }
0245
0246 return accept;
0247 }
0248
0249
0250 #include "FWCore/Framework/interface/MakerMacros.h"
0251 DEFINE_FWK_MODULE(HLTMuonL1TRegionalFilter);