File indexing completed on 2024-04-06 12:21:21
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013 #include <iostream>
0014 #include <iomanip>
0015 #include <iterator>
0016 #include <cmath>
0017 #include <map>
0018
0019 #include "L1Trigger/L1TTwinMux/interface/RPCHitCleaner.h"
0020 #include "L1Trigger/L1TTwinMux/interface/RPCtoDTTranslator.h"
0021 #include "Geometry/RPCGeometry/interface/RPCRoll.h"
0022 #include "Geometry/RPCGeometry/interface/RPCGeometry.h"
0023 #include "CondFormats/DataRecord/interface/RPCEMapRcd.h"
0024 #include "Geometry/Records/interface/MuonGeometryRecord.h"
0025 #include "L1Trigger/L1TCommon/interface/BitShift.h"
0026
0027 using namespace std;
0028
0029 RPCtoDTTranslator::RPCtoDTTranslator(RPCDigiCollection const& inrpcDigis) : m_rpcDigis{inrpcDigis} {}
0030
0031 namespace {
0032 constexpr int max_rpc_bx = 2;
0033 constexpr int min_rpc_bx = -2;
0034
0035 struct rpc_hit {
0036 int bx;
0037 int station;
0038 int sector;
0039 int wheel;
0040 RPCDetId detid;
0041 int strip;
0042 int roll;
0043 int layer;
0044 rpc_hit(int pbx, int pstation, int psector, int pwheel, RPCDetId pdet, int pstrip, int proll, int player)
0045 : bx(pbx),
0046 station(pstation),
0047 sector(psector),
0048 wheel(pwheel),
0049 detid(pdet),
0050 strip(pstrip),
0051 roll(proll),
0052 layer(player) {}
0053 };
0054
0055
0056
0057 class BxToHit {
0058 public:
0059 BxToHit() : m_hits{} {}
0060
0061 static bool outOfRange(int iBX) { return (iBX > max_rpc_bx or iBX < min_rpc_bx); }
0062
0063 int& operator[](int iBX) { return m_hits[iBX - min_rpc_bx]; }
0064
0065 size_t size() const { return m_hits.size(); }
0066
0067 private:
0068 std::array<int, max_rpc_bx - min_rpc_bx + 1> m_hits;
0069 };
0070 }
0071
0072 void RPCtoDTTranslator::run(const RPCGeometry& rpcGeometry) {
0073 std::vector<L1MuDTChambPhDigi> l1ttma_out;
0074 std::vector<L1MuDTChambPhDigi> l1ttma_hits_out;
0075
0076 std::vector<rpc_hit> vrpc_hit_layer1, vrpc_hit_layer2, vrpc_hit_st3, vrpc_hit_st4;
0077
0078
0079 for (auto chamber = m_rpcDigis.begin(); chamber != m_rpcDigis.end(); ++chamber) {
0080 RPCDetId detid = (*chamber).first;
0081 for (auto digi = (*chamber).second.first; digi != (*chamber).second.second; ++digi) {
0082 if (detid.region() != 0)
0083 continue;
0084 if (BxToHit::outOfRange(digi->bx()))
0085 continue;
0086 if (detid.layer() == 1)
0087 vrpc_hit_layer1.emplace_back(digi->bx(),
0088 detid.station(),
0089 detid.sector(),
0090 detid.ring(),
0091 detid,
0092 digi->strip(),
0093 detid.roll(),
0094 detid.layer());
0095 if (detid.station() == 3)
0096 vrpc_hit_st3.emplace_back(digi->bx(),
0097 detid.station(),
0098 detid.sector(),
0099 detid.ring(),
0100 detid,
0101 digi->strip(),
0102 detid.roll(),
0103 detid.layer());
0104 if (detid.layer() == 2)
0105 vrpc_hit_layer2.emplace_back(digi->bx(),
0106 detid.station(),
0107 detid.sector(),
0108 detid.ring(),
0109 detid,
0110 digi->strip(),
0111 detid.roll(),
0112 detid.layer());
0113 if (detid.station() == 4)
0114 vrpc_hit_st4.emplace_back(digi->bx(),
0115 detid.station(),
0116 detid.sector(),
0117 detid.ring(),
0118 detid,
0119 digi->strip(),
0120 detid.roll(),
0121 detid.layer());
0122 }
0123 }
0124
0125 vector<int> vcluster_size;
0126 int cluster_id = -1;
0127 int itr = 0;
0128
0129 std::map<RPCHitCleaner::detId_Ext, int> hits;
0130 int cluster_size = 0;
0131 for (auto chamber = m_rpcDigis.begin(); chamber != m_rpcDigis.end(); ++chamber) {
0132 RPCDetId detid = (*chamber).first;
0133 int strip_n1 = -10000;
0134 int bx_n1 = -10000;
0135 if (detid.region() != 0)
0136 continue;
0137 for (auto digi = (*chamber).second.first; digi != (*chamber).second.second; ++digi) {
0138 if (fabs(digi->bx()) > 3)
0139 continue;
0140
0141
0142 if (abs(digi->strip() - strip_n1) != 1 || digi->bx() != bx_n1) {
0143 if (itr != 0)
0144 vcluster_size.push_back(cluster_size);
0145 cluster_size = 0;
0146 cluster_id++;
0147 }
0148 itr++;
0149 cluster_size++;
0150
0151
0152 RPCHitCleaner::detId_Ext tmp{detid, digi->bx(), digi->strip()};
0153 hits[tmp] = cluster_id;
0154
0155 strip_n1 = digi->strip();
0156 bx_n1 = digi->bx();
0157 }
0158 }
0159 vcluster_size.push_back(cluster_size);
0160
0161 for (int wh = -2; wh <= 2; wh++) {
0162 for (int sec = 1; sec <= 12; sec++) {
0163 for (int st = 1; st <= 4; st++) {
0164 int rpcbx = 0;
0165 std::vector<int> delta_phib;
0166 bool found_hits = false;
0167 std::vector<int> rpc2dt_phi, rpc2dt_phib;
0168
0169 int itr1 = 0;
0170 for (unsigned int l1 = 0; l1 < vrpc_hit_layer1.size(); l1++) {
0171 RPCHitCleaner::detId_Ext tmp{vrpc_hit_layer1[l1].detid, vrpc_hit_layer1[l1].bx, vrpc_hit_layer1[l1].strip};
0172 int id = hits[tmp];
0173 int phi1 = radialAngle(vrpc_hit_layer1[l1].detid, rpcGeometry, vrpc_hit_layer1[l1].strip);
0174 if (vcluster_size[id] == 2 && itr1 == 0) {
0175 itr1++;
0176 continue;
0177 }
0178 if (vcluster_size[id] == 2 && itr1 == 1) {
0179 itr1 = 0;
0180 phi1 = phi1 + (radialAngle(vrpc_hit_layer1[l1 - 1].detid, rpcGeometry, vrpc_hit_layer1[l1 - 1].strip));
0181 phi1 /= 2;
0182 }
0183 int itr2 = 0;
0184 for (unsigned int l2 = 0; l2 < vrpc_hit_layer2.size(); l2++) {
0185 if (vrpc_hit_layer1[l1].station != st || vrpc_hit_layer2[l2].station != st)
0186 continue;
0187 if (vrpc_hit_layer1[l1].sector != sec || vrpc_hit_layer2[l2].sector != sec)
0188 continue;
0189 if (vrpc_hit_layer1[l1].wheel != wh || vrpc_hit_layer2[l2].wheel != wh)
0190 continue;
0191 if (vrpc_hit_layer1[l1].bx != vrpc_hit_layer2[l2].bx)
0192 continue;
0193 RPCHitCleaner::detId_Ext tmp{vrpc_hit_layer2[l2].detid, vrpc_hit_layer2[l2].bx, vrpc_hit_layer2[l2].strip};
0194 int id = hits[tmp];
0195
0196 if (vcluster_size[id] == 2 && itr2 == 0) {
0197 itr2++;
0198 continue;
0199 }
0200
0201
0202 int phi2 = radialAngle(vrpc_hit_layer2[l2].detid, rpcGeometry, vrpc_hit_layer2[l2].strip);
0203 if (vcluster_size[id] == 2 && itr2 == 1) {
0204 itr2 = 0;
0205 phi2 = phi2 + (radialAngle(vrpc_hit_layer2[l2 - 1].detid, rpcGeometry, vrpc_hit_layer2[l2 - 1].strip));
0206 phi2 /= 2;
0207 }
0208 int average = l1t::bitShift(((phi1 + phi2) / 2), 2);
0209 rpc2dt_phi.push_back(average);
0210
0211
0212
0213
0214 int xin = localXX(l1t::bitShift(phi1, 2), 1, vrpc_hit_layer1[l1].station);
0215 int xout = localXX(l1t::bitShift(phi2, 2), 2, vrpc_hit_layer2[l2].station);
0216 if (vcluster_size[id] == 2 && itr2 == 1) {
0217 int phi1_n1 = radialAngle(vrpc_hit_layer1[l1 - 1].detid, rpcGeometry, vrpc_hit_layer1[l1 - 1].strip);
0218 int phi2_n1 = radialAngle(vrpc_hit_layer2[l2 - 1].detid, rpcGeometry, vrpc_hit_layer2[l2 - 1].strip);
0219 xin += localXX(l1t::bitShift(phi1_n1, 2), 1, vrpc_hit_layer1[l1].station);
0220 xout += localXX(l1t::bitShift(phi2_n1, 2), 2, vrpc_hit_layer2[l2].station);
0221 xin /= 2;
0222 xout /= 2;
0223 }
0224
0225 int phi_b = bendingAngle(xin, xout, average);
0226
0227 rpc2dt_phib.push_back(phi_b);
0228
0229 delta_phib.push_back(abs(phi_b));
0230 found_hits = true;
0231 rpcbx = vrpc_hit_layer1[l1].bx;
0232 }
0233 }
0234 if (found_hits) {
0235
0236 int min_index = std::distance(delta_phib.begin(), std::min_element(delta_phib.begin(), delta_phib.end())) + 0;
0237
0238 l1ttma_out.emplace_back(rpcbx, wh, sec - 1, st, rpc2dt_phi[min_index], rpc2dt_phib[min_index], 3, 0, 0, 2);
0239 }
0240
0241 BxToHit hit;
0242 itr1 = 0;
0243 for (unsigned int l1 = 0; l1 < vrpc_hit_layer1.size(); l1++) {
0244 if (vrpc_hit_layer1[l1].station != st || st > 2 || vrpc_hit_layer1[l1].sector != sec ||
0245 vrpc_hit_layer1[l1].wheel != wh)
0246 continue;
0247
0248 RPCHitCleaner::detId_Ext tmp{vrpc_hit_layer1[l1].detid, vrpc_hit_layer1[l1].bx, vrpc_hit_layer1[l1].strip};
0249 int id = hits[tmp];
0250 if (vcluster_size[id] == 2 && itr1 == 0) {
0251 itr1++;
0252 continue;
0253 }
0254 int phi2 = radialAngle(vrpc_hit_layer1[l1].detid, rpcGeometry, vrpc_hit_layer1[l1].strip);
0255 phi2 = l1t::bitShift(phi2, 2);
0256 if (vcluster_size[id] == 2 && itr1 == 1) {
0257 itr1 = 0;
0258 phi2 =
0259 phi2 + l1t::bitShift(
0260 radialAngle(vrpc_hit_layer1[l1 - 1].detid, rpcGeometry, vrpc_hit_layer1[l1 - 1].strip), 2);
0261 phi2 /= 2;
0262 }
0263
0264 l1ttma_hits_out.emplace_back(
0265 vrpc_hit_layer1[l1].bx, wh, sec - 1, st, phi2, 0, 3, hit[vrpc_hit_layer1[l1].bx], 0, 2);
0266 hit[vrpc_hit_layer1[l1].bx]++;
0267 }
0268 itr1 = 0;
0269 for (unsigned int l2 = 0; l2 < vrpc_hit_layer2.size(); l2++) {
0270 if (vrpc_hit_layer2[l2].station != st || st > 2 || vrpc_hit_layer2[l2].sector != sec ||
0271 vrpc_hit_layer2[l2].wheel != wh)
0272 continue;
0273 RPCHitCleaner::detId_Ext tmp{vrpc_hit_layer2[l2].detid, vrpc_hit_layer2[l2].bx, vrpc_hit_layer2[l2].strip};
0274 int id = hits[tmp];
0275
0276 if (vcluster_size[id] == 2 && itr1 == 0) {
0277 itr1++;
0278 continue;
0279 }
0280 int phi2 = radialAngle(vrpc_hit_layer2[l2].detid, rpcGeometry, vrpc_hit_layer2[l2].strip);
0281 phi2 = l1t::bitShift(phi2, 2);
0282 if (vcluster_size[id] == 2 && itr1 == 1) {
0283 itr1 = 0;
0284 phi2 =
0285 phi2 + l1t::bitShift(
0286 radialAngle(vrpc_hit_layer2[l2 - 1].detid, rpcGeometry, vrpc_hit_layer2[l2 - 1].strip), 2);
0287 phi2 /= 2;
0288 }
0289 l1ttma_hits_out.emplace_back(
0290 vrpc_hit_layer2[l2].bx, wh, sec - 1, st, phi2, 0, 3, hit[vrpc_hit_layer2[l2].bx], 0, 2);
0291 hit[vrpc_hit_layer2[l2].bx]++;
0292 }
0293 itr1 = 0;
0294
0295 for (unsigned int l1 = 0; l1 < vrpc_hit_st3.size(); l1++) {
0296 if (st != 3 || vrpc_hit_st3[l1].station != 3 || vrpc_hit_st3[l1].wheel != wh ||
0297 vrpc_hit_st3[l1].sector != sec)
0298 continue;
0299 RPCHitCleaner::detId_Ext tmp{vrpc_hit_st3[l1].detid, vrpc_hit_st3[l1].bx, vrpc_hit_st3[l1].strip};
0300 int id = hits[tmp];
0301
0302 if (vcluster_size[id] == 2 && itr1 == 0) {
0303 itr1++;
0304 continue;
0305 }
0306 int phi2 = radialAngle(vrpc_hit_st3[l1].detid, rpcGeometry, vrpc_hit_st3[l1].strip);
0307 phi2 = l1t::bitShift(phi2, 2);
0308 if (vcluster_size[id] == 2 && itr1 == 1) {
0309 itr1 = 0;
0310 phi2 = phi2 +
0311 l1t::bitShift(radialAngle(vrpc_hit_st3[l1 - 1].detid, rpcGeometry, vrpc_hit_st3[l1 - 1].strip), 2);
0312 phi2 /= 2;
0313 }
0314 l1ttma_hits_out.emplace_back(
0315 vrpc_hit_st3[l1].bx, wh, sec - 1, st, phi2, 0, 3, hit[vrpc_hit_st3[l1].bx], 0, 2);
0316 hit[vrpc_hit_st3[l1].bx]++;
0317 }
0318 itr1 = 0;
0319
0320 for (unsigned int l1 = 0; l1 < vrpc_hit_st4.size(); l1++) {
0321 if (st != 4 || vrpc_hit_st4[l1].station != 4 || vrpc_hit_st4[l1].wheel != wh ||
0322 vrpc_hit_st4[l1].sector != sec)
0323 continue;
0324
0325 RPCHitCleaner::detId_Ext tmp{vrpc_hit_st4[l1].detid, vrpc_hit_st4[l1].bx, vrpc_hit_st4[l1].strip};
0326 int id = hits[tmp];
0327 if (vcluster_size[id] == 2 && itr1 == 0) {
0328 itr1++;
0329 continue;
0330 }
0331 int phi2 = radialAngle(vrpc_hit_st4[l1].detid, rpcGeometry, vrpc_hit_st4[l1].strip);
0332 phi2 = l1t::bitShift(phi2, 2);
0333 if (vcluster_size[id] == 2 && itr1 == 1) {
0334 itr1 = 0;
0335 phi2 = phi2 +
0336 l1t::bitShift(radialAngle(vrpc_hit_st4[l1 - 1].detid, rpcGeometry, vrpc_hit_st4[l1 - 1].strip), 2);
0337 phi2 /= 2;
0338 }
0339 l1ttma_hits_out.emplace_back(
0340 vrpc_hit_st4[l1].bx, wh, sec - 1, st, phi2, 0, 3, hit[vrpc_hit_st4[l1].bx], 0, 2);
0341 hit[vrpc_hit_st4[l1].bx]++;
0342
0343
0344
0345 }
0346 }
0347 }
0348 }
0349
0350 m_rpcdt_translated.setContainer(l1ttma_out);
0351
0352 m_rpchitsdt_translated.setContainer(l1ttma_hits_out);
0353 }
0354
0355
0356 int RPCtoDTTranslator::radialAngle(RPCDetId detid, const RPCGeometry& rpcGeometry, int strip) {
0357 int radialAngle;
0358
0359
0360 const RPCRoll* roll = rpcGeometry.roll(detid);
0361 GlobalPoint stripPosition = roll->toGlobal(roll->centreOfStrip(strip));
0362
0363 double globalphi = stripPosition.phi();
0364 int sector = (roll->id()).sector();
0365 if (sector == 1)
0366 radialAngle = int(globalphi * 1024);
0367 else {
0368 if (globalphi >= 0)
0369 radialAngle = int((globalphi - (sector - 1) * Geom::pi() / 6.) * 1024);
0370 else
0371 radialAngle = int((globalphi + (13 - sector) * Geom::pi() / 6.) * 1024);
0372 }
0373 return radialAngle;
0374 }
0375
0376
0377 int RPCtoDTTranslator::localX(RPCDetId detid, const RPCGeometry& rpcGeometry, int strip) {
0378 const RPCRoll* roll = rpcGeometry.roll(detid);
0379
0380
0381 GlobalPoint p1cmPRG = roll->toGlobal(LocalPoint(1, 0, 0));
0382 GlobalPoint m1cmPRG = roll->toGlobal(LocalPoint(-1, 0, 0));
0383 float phip1cm = p1cmPRG.phi();
0384 float phim1cm = m1cmPRG.phi();
0385 int direction = (phip1cm - phim1cm) / abs(phip1cm - phim1cm);
0386
0387
0388 return direction * roll->centreOfStrip(strip).x();
0389 }
0390
0391 int RPCtoDTTranslator::bendingAngle(int xin, int xout, int phi) {
0392
0393 int atanv = (int)(atan((xout - xin) / 34.6) * 512);
0394 if (atanv > 512)
0395 return 512;
0396 int rvalue = atanv - phi / 8;
0397 return rvalue;
0398 }
0399
0400 int RPCtoDTTranslator::localXX(int phi, int layer, int station) {
0401
0402 double R[2][2] = {{410.0, 444.8}, {492.7, 527.3}};
0403 double rvalue = R[station - 1][layer - 1] * tan(phi / 4096.);
0404 return rvalue;
0405 }