File indexing completed on 2025-03-14 23:36:30
0001 #include "L1Trigger/Phase2L1GMT/interface/L1TPhase2GMTEndcapStubProcessor.h"
0002 #include "L1Trigger/L1TMuon/interface/MuonTriggerPrimitive.h"
0003 #include <cmath>
0004 #include <iostream>
0005 #include <string>
0006 #include <sstream>
0007
0008 L1TPhase2GMTEndcapStubProcessor::L1TPhase2GMTEndcapStubProcessor() : minBX_(-3), maxBX_(3) {}
0009
0010 L1TPhase2GMTEndcapStubProcessor::L1TPhase2GMTEndcapStubProcessor(const edm::ParameterSet& iConfig)
0011 : minBX_(iConfig.getParameter<int>("minBX")),
0012 maxBX_(iConfig.getParameter<int>("maxBX")),
0013 coord1LSB_(iConfig.getParameter<double>("coord1LSB")),
0014 coord2LSB_(iConfig.getParameter<double>("coord2LSB")),
0015 eta1LSB_(iConfig.getParameter<double>("eta1LSB")),
0016 eta2LSB_(iConfig.getParameter<double>("eta2LSB")),
0017 etaMatch_(iConfig.getParameter<double>("etaMatch")),
0018 phiMatch_(iConfig.getParameter<double>("phiMatch")),
0019 verbose_(iConfig.getParameter<unsigned int>("verbose")) {}
0020
0021 L1TPhase2GMTEndcapStubProcessor::~L1TPhase2GMTEndcapStubProcessor() {}
0022
0023 l1t::MuonStub L1TPhase2GMTEndcapStubProcessor::buildCSCOnlyStub(const CSCDetId& detid,
0024 const CSCCorrelatedLCTDigi& digi,
0025 const L1TMuon::GeometryTranslator* translator,
0026 unsigned int tag) {
0027 int endcap = detid.endcap();
0028 int station = detid.station();
0029 int chamber = detid.chamber();
0030 int ring = detid.ring();
0031
0032 L1TMuon::TriggerPrimitive primitive(detid, digi);
0033
0034 const GlobalPoint& gp = translator->getGlobalPoint(primitive);
0035
0036 int phi = int(gp.phi().value() / coord1LSB_);
0037 int eta1 = int(gp.eta() / eta1LSB_);
0038
0039 int wheel = 0;
0040
0041 int sign = endcap == 2 ? -1 : 1;
0042
0043 if (ring == 3)
0044 wheel = sign * 3;
0045 else if (ring == 2)
0046 wheel = sign * 4;
0047 else if (ring == 1)
0048 wheel = sign * 5;
0049
0050 int sector = fabs(chamber);
0051
0052 int bx = digi.getBX() - 8;
0053 int quality = 1;
0054
0055 uint tfLayer = 0;
0056 if ((ring == 3 || ring == 2) && station == 1)
0057 tfLayer = 4;
0058 else if (ring == 1 && station == 1)
0059 tfLayer = 0;
0060 else if (station == 2)
0061 tfLayer = 2;
0062 else if (station == 3)
0063 tfLayer = 1;
0064 else if (station == 4)
0065 tfLayer = 3;
0066
0067 l1t::MuonStub stub(wheel, sector, station, tfLayer, phi, 0, tag, bx, quality, eta1, 0, 1, 0);
0068
0069 stub.setOfflineQuantities(gp.phi().value(), 0.0, gp.eta(), 0.0);
0070 return stub;
0071 }
0072
0073 l1t::MuonStub L1TPhase2GMTEndcapStubProcessor::buildRPCOnlyStub(const RPCDetId& detid,
0074 const RPCDigi& digi,
0075 const L1TMuon::GeometryTranslator* translator) {
0076 L1TMuon::TriggerPrimitive primitive(detid, digi);
0077 const GlobalPoint& gp = translator->getGlobalPoint(primitive);
0078
0079 int phi2 = int(gp.phi().value() / coord2LSB_);
0080 int eta2 = int(gp.eta() / eta2LSB_);
0081
0082 int wheel = -(6 - detid.ring()) * detid.region();
0083 int sector = (detid.sector() - 1) * 6 + detid.subsector();
0084 int station = detid.station();
0085 bool tag = detid.trIndex();
0086 int bx = digi.bx();
0087 int quality = 2;
0088
0089 int ring = detid.ring();
0090
0091 uint tfLayer = 0;
0092
0093 if ((ring == 3 || ring == 2) && station == 1)
0094 tfLayer = 4;
0095
0096
0097 else if (station == 2)
0098 tfLayer = 2;
0099 else if (station == 3)
0100 tfLayer = 1;
0101 else if (station == 4)
0102 tfLayer = 3;
0103
0104 l1t::MuonStub stub(wheel, sector, station, tfLayer, 0, phi2, tag, bx, quality, 0, eta2, 2, 0);
0105 stub.setOfflineQuantities(gp.phi().value(), gp.phi().value(), gp.eta(), gp.eta());
0106 return stub;
0107 }
0108
0109 l1t::MuonStubCollection L1TPhase2GMTEndcapStubProcessor::combineStubs(const l1t::MuonStubCollection& cscStubs,
0110 const l1t::MuonStubCollection& rpcStubs) {
0111 l1t::MuonStubCollection out;
0112 l1t::MuonStubCollection usedRPC;
0113 l1t::MuonStubCollection cleanedRPC;
0114
0115 l1t::MuonStubCollection cleanedCSC;
0116
0117
0118 l1t::MuonStubCollection allCSC = cscStubs;
0119
0120 while (!allCSC.empty()) {
0121 l1t::MuonStub stub = allCSC[0];
0122 l1t::MuonStubCollection freeCSC;
0123 for (uint i = 1; i < allCSC.size(); ++i) {
0124 if ((stub.etaRegion() == allCSC[i].etaRegion()) && (stub.depthRegion() == allCSC[i].depthRegion()) &&
0125 (fabs(deltaPhi(stub.offline_coord1(), allCSC[i].offline_coord1())) < 0.001)) {
0126 if (fabs(stub.offline_eta1() - allCSC[i].offline_eta1()) > 0.001) {
0127 stub.setEta(stub.eta1(), allCSC[i].eta1(), 3);
0128 stub.setOfflineQuantities(stub.offline_coord1(), 0.0, stub.offline_eta1(), allCSC[i].offline_eta1());
0129 }
0130 } else {
0131 freeCSC.push_back(allCSC[i]);
0132 }
0133 }
0134 cleanedCSC.push_back(stub);
0135 allCSC = freeCSC;
0136 }
0137
0138 for (const auto& csc : cleanedCSC) {
0139 int nRPC = 0;
0140 float phiF = 0.0;
0141 float etaF = 0.0;
0142 int phi = 0;
0143 int eta = 0;
0144 for (const auto& rpc : rpcStubs) {
0145 if (csc.depthRegion() != rpc.depthRegion())
0146 continue;
0147 if (fabs(deltaPhi(csc.offline_coord1(), rpc.offline_coord2())) < phiMatch_ &&
0148 fabs(csc.offline_eta1() - rpc.offline_eta2()) < etaMatch_ && csc.bxNum() == rpc.bxNum()) {
0149 phiF += rpc.offline_coord2();
0150 etaF += rpc.offline_eta2();
0151 phi += rpc.coord2();
0152 eta += rpc.eta2();
0153 nRPC++;
0154 usedRPC.push_back(rpc);
0155 }
0156 }
0157
0158 int finalRPCPhi = 0;
0159 int finalRPCEta = 0;
0160 double offline_finalRPCPhi = 0;
0161 double offline_finalRPCEta = 0;
0162 if (nRPC != 0) {
0163 finalRPCPhi = phi / nRPC;
0164 finalRPCEta = eta / nRPC;
0165 offline_finalRPCPhi = phiF / nRPC;
0166 offline_finalRPCEta = etaF / nRPC;
0167 l1t::MuonStub stub(csc.etaRegion(),
0168 csc.phiRegion(),
0169 csc.depthRegion(),
0170 csc.tfLayer(),
0171 csc.coord1(),
0172 finalRPCPhi,
0173 0,
0174 csc.bxNum(),
0175 3,
0176 csc.eta1(),
0177 finalRPCEta,
0178 3,
0179 0);
0180 stub.setOfflineQuantities(csc.offline_coord1(), offline_finalRPCPhi, csc.offline_eta1(), offline_finalRPCEta);
0181 out.push_back(stub);
0182 } else {
0183 out.push_back(csc);
0184 }
0185 }
0186
0187
0188
0189 for (const auto& rpc : rpcStubs) {
0190 bool keep = true;
0191 for (const auto& rpc2 : usedRPC) {
0192 if (rpc == rpc2) {
0193 keep = false;
0194 break;
0195 }
0196 }
0197 if (keep)
0198 cleanedRPC.push_back(rpc);
0199 }
0200
0201 while (!cleanedRPC.empty()) {
0202 l1t::MuonStubCollection freeRPC;
0203
0204 int nRPC = 1;
0205 float phiF = cleanedRPC[0].offline_coord2();
0206 float etaF = cleanedRPC[0].offline_eta2();
0207 int phi = cleanedRPC[0].coord2();
0208 int eta = cleanedRPC[0].eta2();
0209
0210 for (unsigned i = 1; i < cleanedRPC.size(); ++i) {
0211 if (fabs(deltaPhi(cleanedRPC[0].offline_coord2(), cleanedRPC[i].offline_coord2())) < phiMatch_ &&
0212 cleanedRPC[0].depthRegion() == cleanedRPC[i].depthRegion() &&
0213 fabs(cleanedRPC[0].offline_eta2() - cleanedRPC[i].offline_eta2()) < etaMatch_ &&
0214 cleanedRPC[0].bxNum() == cleanedRPC[i].bxNum()) {
0215 phiF += cleanedRPC[i].offline_coord2();
0216 etaF += cleanedRPC[i].offline_eta2();
0217 phi += cleanedRPC[i].coord2();
0218 eta += cleanedRPC[i].eta2();
0219 nRPC++;
0220 } else {
0221 freeRPC.push_back(cleanedRPC[i]);
0222 }
0223 }
0224 l1t::MuonStub stub(cleanedRPC[0].etaRegion(),
0225 cleanedRPC[0].phiRegion(),
0226 cleanedRPC[0].depthRegion(),
0227 cleanedRPC[0].tfLayer(),
0228 0,
0229 phi / nRPC,
0230 0,
0231 cleanedRPC[0].bxNum(),
0232 2,
0233 0,
0234 eta / nRPC,
0235 2,
0236 0);
0237 stub.setOfflineQuantities(phiF / nRPC, phiF / nRPC, etaF / nRPC, etaF / nRPC);
0238 out.push_back(stub);
0239 cleanedRPC = freeRPC;
0240 };
0241 return out;
0242 }
0243
0244 l1t::MuonStubCollection L1TPhase2GMTEndcapStubProcessor::makeStubs(
0245 const MuonDigiCollection<CSCDetId, CSCCorrelatedLCTDigi>& csc,
0246 const MuonDigiCollection<RPCDetId, RPCDigi>& cleaned,
0247 const L1TMuon::GeometryTranslator* t,
0248 const edm::EventSetup& iSetup) {
0249 l1t::MuonStubCollection cscStubs;
0250 auto chamber = csc.begin();
0251 auto chend = csc.end();
0252 for (; chamber != chend; ++chamber) {
0253 auto digi = (*chamber).second.first;
0254 auto dend = (*chamber).second.second;
0255 unsigned int tag = 0;
0256 for (; digi != dend; ++digi) {
0257 l1t::MuonStub stub = buildCSCOnlyStub((*chamber).first, *digi, t, tag);
0258 tag = tag + 1;
0259 if (stub.bxNum() >= minBX_ && stub.bxNum() <= maxBX_)
0260 cscStubs.push_back(stub);
0261 }
0262 }
0263
0264 l1t::MuonStubCollection rpcStubs;
0265
0266 auto rpcchamber = cleaned.begin();
0267 auto rpcchend = cleaned.end();
0268 for (; rpcchamber != rpcchend; ++rpcchamber) {
0269 if ((*rpcchamber).first.region() == 0)
0270 continue;
0271 auto digi = (*rpcchamber).second.first;
0272 auto dend = (*rpcchamber).second.second;
0273 for (; digi != dend; ++digi) {
0274 l1t::MuonStub stub = buildRPCOnlyStub((*rpcchamber).first, *digi, t);
0275 if (stub.bxNum() >= minBX_ && stub.bxNum() <= maxBX_)
0276 rpcStubs.push_back(stub);
0277 }
0278 }
0279
0280 l1t::MuonStubCollection combinedStubs = combineStubs(cscStubs, rpcStubs);
0281
0282 if (verbose_) {
0283 edm::LogInfo("EndcapStub") << "CSC Stubs";
0284 for (const auto& stub : cscStubs)
0285 edm::LogInfo("EndcapStub") << "CSC Stub bx=" << stub.bxNum() << " TF=" << stub.tfLayer()
0286 << " etaRegion=" << stub.etaRegion() << " phiRegion=" << stub.phiRegion()
0287 << " depthRegion=" << stub.depthRegion() << " coord1=" << stub.offline_coord1() << ","
0288 << stub.coord1() << " coord2=" << stub.offline_coord2() << "," << stub.coord2()
0289 << " eta1=" << stub.offline_eta1() << "," << stub.eta1()
0290 << " eta2=" << stub.offline_eta2() << "," << stub.eta2()
0291 << " quality=" << stub.quality() << " etaQuality=" << stub.etaQuality();
0292
0293 edm::LogInfo("EndcapStub") << "RPC Stubs";
0294 for (const auto& stub : rpcStubs)
0295 edm::LogInfo("EndcapStub") << "RPC Stub bx=" << stub.bxNum() << " TF=" << stub.tfLayer()
0296 << " etaRegion=" << stub.etaRegion() << " phiRegion=" << stub.phiRegion()
0297 << " depthRegion=" << stub.depthRegion() << " coord1=" << stub.offline_coord1() << ","
0298 << stub.coord1() << " coord2=" << stub.offline_coord2() << "," << stub.coord2()
0299 << " eta1=" << stub.offline_eta1() << "," << stub.eta1()
0300 << " eta2=" << stub.offline_eta2() << "," << stub.eta2()
0301 << " quality=" << stub.quality() << " etaQuality=" << stub.etaQuality();
0302
0303 for (const auto& stub : combinedStubs)
0304 edm::LogInfo("EndcapStub") << "Combined Stub bx=" << stub.bxNum() << " TF=" << stub.tfLayer()
0305 << " etaRegion=" << stub.etaRegion() << " phiRegion=" << stub.phiRegion()
0306 << " depthRegion=" << stub.depthRegion() << " coord1=" << stub.offline_coord1() << ","
0307 << stub.coord1() << " coord2=" << stub.offline_coord2() << "," << stub.coord2()
0308 << " eta1=" << stub.offline_eta1() << "," << stub.eta1()
0309 << " eta2=" << stub.offline_eta2() << "," << stub.eta2()
0310 << " quality=" << stub.quality() << " etaQuality=" << stub.etaQuality();
0311 }
0312
0313 return combinedStubs;
0314 }