File indexing completed on 2024-09-07 04:36:57
0001 #include "L1Trigger/L1TMuonEndCap/interface/AngleCalculation.h"
0002 #include "DataFormats/L1TMuon/interface/L1TMuonSubsystems.h"
0003 #include "helper.h" // to_hex, to_binary
0004
0005 namespace {
0006 const int bw_fph = 13;
0007 const int bw_th = 7;
0008 const int invalid_dtheta = (1 << bw_th) - 1;
0009 const int invalid_dphi = (1 << bw_fph) - 1;
0010 }
0011
0012 void AngleCalculation::configure(int verbose,
0013 int endcap,
0014 int sector,
0015 int bx,
0016 int bxWindow,
0017 int thetaWindow,
0018 int thetaWindowZone0,
0019 bool bugME11Dupes,
0020 bool bugAmbigThetaWin,
0021 bool twoStationSameBX) {
0022 verbose_ = verbose;
0023 endcap_ = endcap;
0024 sector_ = sector;
0025 bx_ = bx;
0026
0027 bxWindow_ = bxWindow;
0028 thetaWindow_ = thetaWindow;
0029 thetaWindowZone0_ = thetaWindowZone0;
0030 bugME11Dupes_ = bugME11Dupes;
0031 bugAmbigThetaWin_ = bugAmbigThetaWin;
0032 twoStationSameBX_ = twoStationSameBX;
0033 }
0034
0035 void AngleCalculation::process(emtf::zone_array<EMTFTrackCollection>& zone_tracks) const {
0036 for (int izone = 0; izone < emtf::NUM_ZONES; ++izone) {
0037 EMTFTrackCollection& tracks = zone_tracks.at(izone);
0038
0039 EMTFTrackCollection::iterator tracks_it = tracks.begin();
0040 EMTFTrackCollection::iterator tracks_end = tracks.end();
0041
0042
0043 for (; tracks_it != tracks_end; ++tracks_it) {
0044 calculate_angles(*tracks_it, izone);
0045 }
0046
0047
0048
0049
0050 erase_tracks(tracks);
0051
0052 tracks_it = tracks.begin();
0053 tracks_end = tracks.end();
0054
0055
0056
0057 for (; tracks_it != tracks_end; ++tracks_it) {
0058 calculate_bx(*tracks_it);
0059 }
0060 }
0061
0062 if (verbose_ > 0) {
0063 for (const auto& tracks : zone_tracks) {
0064 for (const auto& track : tracks) {
0065 std::cout << "deltas: z: " << track.Zone() - 1 << " pat: " << track.Winner()
0066 << " rank: " << to_hex(track.Rank()) << " delta_ph: " << array_as_string(track.PtLUT().delta_ph)
0067 << " delta_th: " << array_as_string(track.PtLUT().delta_th)
0068 << " sign_ph: " << array_as_string(track.PtLUT().sign_ph)
0069 << " sign_th: " << array_as_string(track.PtLUT().sign_th) << " phi: " << track.Phi_fp()
0070 << " theta: " << track.Theta_fp() << " cpat: " << array_as_string(track.PtLUT().cpattern)
0071 << " v: " << array_as_string(track.PtLUT().bt_vi) << " h: " << array_as_string(track.PtLUT().bt_hi)
0072 << " c: " << array_as_string(track.PtLUT().bt_ci) << " s: " << array_as_string(track.PtLUT().bt_si)
0073 << std::endl;
0074 }
0075 }
0076 }
0077 }
0078
0079 void AngleCalculation::calculate_angles(EMTFTrack& track, const int izone) const {
0080
0081 std::array<EMTFHitCollection, emtf::NUM_STATIONS> st_conv_hits;
0082
0083 for (int istation = 0; istation < emtf::NUM_STATIONS; ++istation) {
0084 for (const auto& conv_hit : track.Hits()) {
0085 if ((conv_hit.Station() - 1) == istation) {
0086 st_conv_hits.at(istation).push_back(conv_hit);
0087 }
0088 }
0089
0090 if (bugME11Dupes_) {
0091 emtf_assert(st_conv_hits.at(istation).size() <= 4);
0092 } else {
0093 emtf_assert(st_conv_hits.at(istation).size() <= 2);
0094 }
0095 }
0096 emtf_assert(st_conv_hits.size() == emtf::NUM_STATIONS);
0097
0098
0099
0100 std::array<int, emtf::NUM_STATION_PAIRS>
0101 best_dtheta_arr;
0102 std::array<int, emtf::NUM_STATION_PAIRS> best_dtheta_sign_arr;
0103 std::array<int, emtf::NUM_STATION_PAIRS>
0104 best_dphi_arr;
0105 std::array<int, emtf::NUM_STATION_PAIRS> best_dphi_sign_arr;
0106
0107
0108
0109
0110 std::array<int, emtf::NUM_STATION_PAIRS> best_theta_arr;
0111 std::array<int, emtf::NUM_STATION_PAIRS> best_phi_arr;
0112
0113
0114 std::array<bool, emtf::NUM_STATION_PAIRS> best_dtheta_valid_arr;
0115
0116
0117
0118 best_dtheta_arr.fill(invalid_dtheta);
0119 best_dtheta_sign_arr.fill(1);
0120 best_dphi_arr.fill(invalid_dphi);
0121 best_dphi_sign_arr.fill(1);
0122 best_phi_arr.fill(0);
0123 best_theta_arr.fill(0);
0124 best_dtheta_valid_arr.fill(false);
0125
0126
0127 auto abs_diff = [](int a, int b) { return std::abs(a - b); };
0128
0129
0130 int ipair = 0;
0131
0132 for (int ist1 = 0; ist1 + 1 < emtf::NUM_STATIONS; ++ist1) {
0133 for (int ist2 = ist1 + 1; ist2 < emtf::NUM_STATIONS; ++ist2) {
0134 const EMTFHitCollection& conv_hitsA = st_conv_hits.at(ist1);
0135 const EMTFHitCollection& conv_hitsB = st_conv_hits.at(ist2);
0136
0137
0138 for (const auto& conv_hitA : conv_hitsA) {
0139 for (const auto& conv_hitB : conv_hitsB) {
0140
0141
0142
0143 int thA = conv_hitA.Theta_fp();
0144 int thB = conv_hitB.Theta_fp();
0145 int dth = abs_diff(thA, thB);
0146 int dth_sign = (thA <= thB);
0147 emtf_assert(thA != 0 && thB != 0);
0148 emtf_assert(dth < invalid_dtheta);
0149
0150 if (best_dtheta_arr.at(ipair) >= dth) {
0151 best_dtheta_arr.at(ipair) = dth;
0152 best_dtheta_sign_arr.at(ipair) = dth_sign;
0153 best_dtheta_valid_arr.at(ipair) = true;
0154
0155
0156
0157
0158 best_theta_arr.at(ipair) = (ipair < 3) ? thB : thA;
0159 }
0160
0161
0162 int phA = conv_hitA.Phi_fp();
0163 int phB = conv_hitB.Phi_fp();
0164 int dph = abs_diff(phA, phB);
0165 int dph_sign = (phA <= phB);
0166
0167 if (best_dphi_arr.at(ipair) >= dph) {
0168 best_dphi_arr.at(ipair) = dph;
0169 best_dphi_sign_arr.at(ipair) = dph_sign;
0170
0171
0172
0173 best_phi_arr.at(ipair) = (ipair < 3) ? phB : phA;
0174 }
0175 }
0176 }
0177
0178 ++ipair;
0179 }
0180 }
0181 emtf_assert(ipair == emtf::NUM_STATION_PAIRS);
0182
0183
0184
0185
0186
0187
0188 std::array<bool, emtf::NUM_STATION_PAIRS> best_dtheta_valid_arr_1;
0189
0190 for (int ipair = 0; ipair < emtf::NUM_STATION_PAIRS; ++ipair) {
0191 if (izone == 0)
0192 best_dtheta_valid_arr_1.at(ipair) =
0193 best_dtheta_valid_arr.at(ipair) && (best_dtheta_arr.at(ipair) <= thetaWindowZone0_);
0194 else
0195 best_dtheta_valid_arr_1.at(ipair) =
0196 best_dtheta_valid_arr.at(ipair) && (best_dtheta_arr.at(ipair) <= thetaWindow_);
0197 }
0198
0199
0200
0201 int vmask1 = 0, vmask2 = 0, vmask3 = 0;
0202
0203 if (best_dtheta_valid_arr_1.at(0)) {
0204 vmask1 |= 0b0011;
0205 }
0206 if (best_dtheta_valid_arr_1.at(1)) {
0207 vmask1 |= 0b0101;
0208 }
0209 if (best_dtheta_valid_arr_1.at(2)) {
0210 vmask1 |= 0b1001;
0211 }
0212 if (best_dtheta_valid_arr_1.at(3)) {
0213 vmask2 |= 0b0110;
0214 }
0215 if (best_dtheta_valid_arr_1.at(4)) {
0216 vmask2 |= 0b1010;
0217 }
0218 if (best_dtheta_valid_arr_1.at(5)) {
0219 vmask3 |= 0b1100;
0220 }
0221
0222
0223
0224
0225
0226
0227
0228
0229
0230
0231 int vstat = vmask1;
0232 if ((vstat & vmask2) != 0 || vstat == 0)
0233 vstat |= vmask2;
0234 if ((vstat & vmask3) != 0 || vstat == 0)
0235 vstat |= vmask3;
0236
0237
0238
0239 static const int trk_bld[64] = {
0240 0b1111, 0b0111, 0b0111, 0b0111, 0b1011, 0b0011, 0b1110, 0b0011, 0b0111, 0b0111, 0b0111, 0b0111, 0b1011,
0241 0b0011, 0b0011, 0b0011, 0b1011, 0b1101, 0b0011, 0b0011, 0b1011, 0b0011, 0b0011, 0b0011, 0b1011, 0b0011,
0242 0b0011, 0b0011, 0b1011, 0b0011, 0b0011, 0b0011, 0b1101, 0b1101, 0b1110, 0b0101, 0b1110, 0b1001, 0b1110,
0243 0b0110, 0b0101, 0b0101, 0b0101, 0b0101, 0b1001, 0b1001, 0b0110, 0b0110, 0b1101, 0b1101, 0b0101, 0b0101,
0244 0b1001, 0b1001, 0b1010, 0b1100, 0b0101, 0b0101, 0b0101, 0b0101, 0b1001, 0b1001, 0b1010, 0b0000};
0245
0246 if (not bugAmbigThetaWin_) {
0247
0248
0249 unsigned dth_bad = 0b111111;
0250 if (best_dtheta_valid_arr_1.at(0)) {
0251 dth_bad &= (~(1 << 5));
0252 }
0253 if (best_dtheta_valid_arr_1.at(1)) {
0254 dth_bad &= (~(1 << 2));
0255 }
0256 if (best_dtheta_valid_arr_1.at(2)) {
0257 dth_bad &= (~(1 << 1));
0258 }
0259 if (best_dtheta_valid_arr_1.at(3)) {
0260 dth_bad &= (~(1 << 4));
0261 }
0262 if (best_dtheta_valid_arr_1.at(4)) {
0263 dth_bad &= (~(1 << 0));
0264 }
0265 if (best_dtheta_valid_arr_1.at(5)) {
0266 dth_bad &= (~(1 << 3));
0267 }
0268 emtf_assert(dth_bad < 64);
0269
0270
0271 vstat = trk_bld[dth_bad];
0272 }
0273
0274
0275 for (int istation = 0; istation < emtf::NUM_STATIONS; ++istation) {
0276 if ((vstat & (1 << istation)) == 0) {
0277 st_conv_hits.at(istation).clear();
0278 }
0279 }
0280
0281
0282 int phi_fp = 0;
0283 int theta_fp = 0;
0284 int best_pair = -1;
0285
0286 if ((vstat & (1 << 1)) != 0) {
0287 if (best_dtheta_valid_arr.at(0))
0288 best_pair = 0;
0289 else if (best_dtheta_valid_arr.at(3))
0290 best_pair = 3;
0291 else if (best_dtheta_valid_arr.at(4))
0292 best_pair = 4;
0293 } else if ((vstat & (1 << 2)) != 0) {
0294 if (best_dtheta_valid_arr.at(1))
0295 best_pair = 1;
0296 else if (best_dtheta_valid_arr.at(5))
0297 best_pair = 5;
0298 } else if ((vstat & (1 << 3)) != 0) {
0299 if (best_dtheta_valid_arr.at(2))
0300 best_pair = 2;
0301 }
0302
0303
0304
0305
0306
0307
0308
0309
0310
0311
0312
0313
0314
0315
0316
0317
0318
0319
0320
0321
0322
0323
0324
0325
0326
0327
0328
0329
0330
0331 if (best_pair != -1) {
0332 phi_fp = best_phi_arr.at(best_pair);
0333 theta_fp = best_theta_arr.at(best_pair);
0334 emtf_assert(theta_fp != 0);
0335
0336
0337
0338
0339 struct {
0340 typedef EMTFHit value_type;
0341 bool operator()(const value_type& lhs, const value_type& rhs) const {
0342 return std::abs(lhs.Theta_fp() - theta) < std::abs(rhs.Theta_fp() - theta);
0343 }
0344 int theta;
0345 } less_dtheta_cmp;
0346 less_dtheta_cmp.theta = theta_fp;
0347
0348 for (int istation = 0; istation < emtf::NUM_STATIONS; ++istation) {
0349 std::stable_sort(st_conv_hits.at(istation).begin(), st_conv_hits.at(istation).end(), less_dtheta_cmp);
0350 if (st_conv_hits.at(istation).size() > 1)
0351 st_conv_hits.at(istation).resize(1);
0352 }
0353 }
0354
0355
0356
0357 int old_rank = (track.Rank() << 1);
0358 int rank = ((((old_rank >> 6) & 1) << 6) |
0359 (((old_rank >> 4) & 1) << 4) |
0360 (((old_rank >> 2) & 1) << 2) |
0361 (((vstat >> 0) & 1) << 5) |
0362 (((vstat >> 1) & 1) << 3) |
0363 (((vstat >> 2) & 1) << 1) |
0364 (((vstat >> 3) & 1) << 0)
0365 );
0366
0367 int mode = ((((vstat >> 0) & 1) << 3) |
0368 (((vstat >> 1) & 1) << 2) |
0369 (((vstat >> 2) & 1) << 1) |
0370 (((vstat >> 3) & 1) << 0)
0371 );
0372
0373 int mode_inv = vstat;
0374
0375
0376 if (vstat == 0b0001 || vstat == 0b0010 || vstat == 0b0100 || vstat == 0b1000 || vstat == 0)
0377 rank = 0;
0378
0379
0380 auto isFront = [](int station, int ring, int chamber, int subsystem) {
0381
0382
0383
0384
0385
0386
0387
0388
0389 if (subsystem == L1TMuon::kGEM)
0390 return true;
0391
0392 bool result = false;
0393 bool isOverlapping = !(station == 1 && ring == 3);
0394
0395 if (isOverlapping) {
0396 bool isEven = (chamber % 2 == 0);
0397
0398
0399 result = (station < 3) ? isEven : !isEven;
0400 }
0401 return result;
0402 };
0403
0404
0405 EMTFPtLUT ptlut_data = {};
0406 for (int i = 0; i < emtf::NUM_STATION_PAIRS; ++i) {
0407 ptlut_data.delta_ph[i] = best_dphi_arr.at(i);
0408 ptlut_data.sign_ph[i] = best_dphi_sign_arr.at(i);
0409 ptlut_data.delta_th[i] = best_dtheta_arr.at(i);
0410 ptlut_data.sign_th[i] = best_dtheta_sign_arr.at(i);
0411 }
0412
0413 for (int i = 0; i < emtf::NUM_STATIONS; ++i) {
0414 const auto& v = st_conv_hits.at(i);
0415 ptlut_data.cpattern[i] = v.empty() ? 0 : v.front().Pattern();
0416 ptlut_data.csign[i] = v.empty() ? 0 : v.front().Bend();
0417 ptlut_data.slope[i] = v.empty() ? 0 : v.front().Slope();
0418 ptlut_data.fr[i] =
0419 v.empty() ? 0 : isFront(v.front().Station(), v.front().Ring(), v.front().Chamber(), v.front().Subsystem());
0420 if (i == 0)
0421 ptlut_data.st1_ring2 =
0422 v.empty() ? 0 : (v.front().Station() == 1 && (v.front().Ring() == 2 || v.front().Ring() == 3));
0423 }
0424
0425 for (int i = 0; i < emtf::NUM_STATIONS + 1; ++i) {
0426 ptlut_data.bt_vi[i] = 0;
0427 ptlut_data.bt_hi[i] = 0;
0428 ptlut_data.bt_ci[i] = 0;
0429 ptlut_data.bt_si[i] = 0;
0430 }
0431
0432 for (int i = 0; i < emtf::NUM_STATIONS; ++i) {
0433 const auto& v = st_conv_hits.at(i);
0434 if (!v.empty()) {
0435 int bt_station = v.front().BT_station();
0436 emtf_assert(0 <= bt_station && bt_station <= 4);
0437
0438 int bt_segment = v.front().BT_segment();
0439 ptlut_data.bt_vi[bt_station] = 1;
0440 ptlut_data.bt_hi[bt_station] = (bt_segment >> 5) & 0x3;
0441 ptlut_data.bt_ci[bt_station] = (bt_segment >> 1) & 0xf;
0442 ptlut_data.bt_si[bt_station] = (bt_segment >> 0) & 0x1;
0443 }
0444 }
0445
0446
0447
0448
0449 track.set_rank(rank);
0450 track.set_mode(mode);
0451 track.set_mode_inv(mode_inv);
0452 track.set_phi_fp(phi_fp);
0453 track.set_theta_fp(theta_fp);
0454 track.set_PtLUT(ptlut_data);
0455
0456 track.set_phi_loc(emtf::calc_phi_loc_deg(phi_fp));
0457 track.set_phi_glob(emtf::calc_phi_glob_deg(track.Phi_loc(), track.Sector()));
0458 track.set_theta(emtf::calc_theta_deg_from_int(theta_fp));
0459 track.set_eta(emtf::calc_eta_from_theta_deg(track.Theta(), track.Endcap()));
0460
0461
0462 track.clear_Hits();
0463
0464 EMTFHitCollection tmp_hits = track.Hits();
0465 flatten_container(st_conv_hits, tmp_hits);
0466 track.set_Hits(tmp_hits);
0467 }
0468
0469 void AngleCalculation::calculate_bx(EMTFTrack& track) const {
0470 const int delayBX = bxWindow_ - 1;
0471 emtf_assert(delayBX >= 0);
0472 std::vector<int> counter(delayBX + 1, 0);
0473
0474 for (const auto& conv_hit : track.Hits()) {
0475 for (int i = delayBX; i >= 0; i--) {
0476 if (conv_hit.BX() <= bx_ - i)
0477 counter.at(i) += 1;
0478 }
0479 }
0480
0481 int first_bx = bx_ - delayBX;
0482 int second_bx = 99;
0483 for (int i = delayBX; i >= 0; i--) {
0484 if (counter.at(i) >= 2) {
0485 second_bx = bx_ - i;
0486 break;
0487 }
0488 }
0489 emtf_assert(second_bx != 99);
0490
0491
0492
0493
0494 track.set_first_bx(first_bx);
0495 track.set_second_bx(second_bx);
0496 }
0497
0498 void AngleCalculation::erase_tracks(EMTFTrackCollection& tracks) const {
0499
0500
0501 struct {
0502 typedef EMTFTrack value_type;
0503 bool operator()(const value_type& x) const { return (x.Rank() == 0); }
0504 } rank_zero_pred;
0505
0506
0507 struct {
0508 typedef EMTFTrack value_type;
0509 bool operator()(const value_type& x) const {
0510 return (x.NumHits() == 2 && x.Hits().at(0).BX() != x.Hits().at(1).BX());
0511 }
0512 } two_station_mistime;
0513
0514 tracks.erase(std::remove_if(tracks.begin(), tracks.end(), rank_zero_pred), tracks.end());
0515
0516 if (twoStationSameBX_) {
0517 tracks.erase(std::remove_if(tracks.begin(), tracks.end(), two_station_mistime), tracks.end());
0518 }
0519
0520 for (const auto& track : tracks) {
0521 emtf_assert(!track.Hits().empty());
0522 emtf_assert(track.Hits().size() <= emtf::NUM_STATIONS);
0523 }
0524 }