Back to home page

Project CMSSW displayed by LXR

 
 

    


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;                        // bit width of ph, full precision
0007   const int bw_th = 7;                          // bit width of th
0008   const int invalid_dtheta = (1 << bw_th) - 1;  // = 127
0009   const int invalid_dphi = (1 << bw_fph) - 1;   // = 8191
0010 }  // namespace
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);  // pass by reference
0038 
0039     EMTFTrackCollection::iterator tracks_it = tracks.begin();
0040     EMTFTrackCollection::iterator tracks_end = tracks.end();
0041 
0042     // Calculate deltas
0043     for (; tracks_it != tracks_end; ++tracks_it) {
0044       calculate_angles(*tracks_it, izone);
0045     }
0046 
0047     // Erase tracks with rank = 0
0048     // Erase hits that are not selected as the best phi and theta in each station
0049     // Erase two-station tracks with hits in different BX (2018)
0050     erase_tracks(tracks);
0051 
0052     tracks_it = tracks.begin();
0053     tracks_end = tracks.end();
0054 
0055     // Calculate bx
0056     // (in the firmware, this happens during best track selection.)
0057     for (; tracks_it != tracks_end; ++tracks_it) {
0058       calculate_bx(*tracks_it);
0059     }
0060   }  // end loop over zones
0061 
0062   if (verbose_ > 0) {  // debug
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   // Group track hits by station
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);  // ambiguity in theta is max 4
0092     } else {
0093       emtf_assert(st_conv_hits.at(istation).size() <= 2);  // ambiguity in theta is max 2
0094     }
0095   }
0096   emtf_assert(st_conv_hits.size() == emtf::NUM_STATIONS);
0097 
0098   // Best theta deltas and phi deltas
0099   // from 0 to 5: dtheta12, dtheta13, dtheta14, dtheta23, dtheta24, dtheta34
0100   std::array<int, emtf::NUM_STATION_PAIRS>
0101       best_dtheta_arr;  // Best of up to 4 dTheta values per pair of stations (with duplicate thetas)
0102   std::array<int, emtf::NUM_STATION_PAIRS> best_dtheta_sign_arr;
0103   std::array<int, emtf::NUM_STATION_PAIRS>
0104       best_dphi_arr;  // Not really "best" - there is only one dPhi value per pair of stations
0105   std::array<int, emtf::NUM_STATION_PAIRS> best_dphi_sign_arr;
0106 
0107   // Best angles
0108   // from 0 to 5: ME2,      ME3,      ME4,      ME2,      ME2,      ME3
0109   //              dtheta12, dtheta13, dtheta14, dtheta23, dtheta24, dtheta34
0110   std::array<int, emtf::NUM_STATION_PAIRS> best_theta_arr;
0111   std::array<int, emtf::NUM_STATION_PAIRS> best_phi_arr;
0112 
0113   // Keep track of which pair is valid
0114   std::array<bool, emtf::NUM_STATION_PAIRS> best_dtheta_valid_arr;
0115   // std::array<bool, emtf::NUM_STATION_PAIRS> best_has_rpc_arr;  // Not used currently
0116 
0117   // Initialize
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   // best_has_rpc_arr     .fill(false);
0126 
0127   auto abs_diff = [](int a, int b) { return std::abs(a - b); };
0128 
0129   // Calculate angles
0130   int ipair = 0;
0131 
0132   for (int ist1 = 0; ist1 + 1 < emtf::NUM_STATIONS; ++ist1) {       // station A
0133     for (int ist2 = ist1 + 1; ist2 < emtf::NUM_STATIONS; ++ist2) {  // station B
0134       const EMTFHitCollection& conv_hitsA = st_conv_hits.at(ist1);
0135       const EMTFHitCollection& conv_hitsB = st_conv_hits.at(ist2);
0136 
0137       // More than 1 hit per station when hit has ambigous theta
0138       for (const auto& conv_hitA : conv_hitsA) {
0139         for (const auto& conv_hitB : conv_hitsB) {
0140           // bool has_rpc = (conv_hitA.Subsystem() == TriggerPrimitive::kRPC || conv_hitB.Subsystem() == TriggerPrimitive::kRPC);
0141 
0142           // Calculate theta deltas
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);  // sign
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             // best_has_rpc_arr.at(ipair) = has_rpc;  // FW doesn't currently check whether a segment is CSC or RPC
0155 
0156             // first 3 pairs, use station B
0157             // last 3 pairs, use station A
0158             best_theta_arr.at(ipair) = (ipair < 3) ? thB : thA;
0159           }
0160 
0161           // Calculate phi deltas
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             // first 3 pairs, use station B
0172             // last 3 pairs, use station A
0173             best_phi_arr.at(ipair) = (ipair < 3) ? phB : phA;
0174           }
0175         }  // end loop over conv_hits in station B
0176       }  // end loop over conv_hits in station A
0177 
0178       ++ipair;
0179     }  // end loop over station B
0180   }  // end loop over station A
0181   emtf_assert(ipair == emtf::NUM_STATION_PAIRS);
0182 
0183   // Apply cuts on dtheta
0184 
0185   // There is a possible bug in FW. After a dtheta pair fails the theta window
0186   // cut, the valid flag of the pair is not updated. Later on, theta from
0187   // this pair is used to assign the precise theta of the track.
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)  // Tighter theta window for Zone 0 (Ring 1), where there are no RPCs
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   // Find valid segments
0200   // vmask contains valid station mask = {ME4,ME3,ME2,ME1}. "0b" prefix for binary.
0201   int vmask1 = 0, vmask2 = 0, vmask3 = 0;
0202 
0203   if (best_dtheta_valid_arr_1.at(0)) {
0204     vmask1 |= 0b0011;  // 12
0205   }
0206   if (best_dtheta_valid_arr_1.at(1)) {
0207     vmask1 |= 0b0101;  // 13
0208   }
0209   if (best_dtheta_valid_arr_1.at(2)) {
0210     vmask1 |= 0b1001;  // 14
0211   }
0212   if (best_dtheta_valid_arr_1.at(3)) {
0213     vmask2 |= 0b0110;  // 23
0214   }
0215   if (best_dtheta_valid_arr_1.at(4)) {
0216     vmask2 |= 0b1010;  // 24
0217   }
0218   if (best_dtheta_valid_arr_1.at(5)) {
0219     vmask3 |= 0b1100;  // 34
0220   }
0221 
0222   // merge station masks only if they share bits
0223   // Station 1 hits pass if any dTheta1X values pass
0224   // Station 2 hits pass if any dTheta2X values pass, *EXCEPT* the following cases:
0225   //           Only {13, 24} pass, only {13, 24, 34} pass,
0226   //           Only {14, 23} pass, only {14, 23, 34} pass.
0227   // Station 3 hits pass if any dTheta3X values pass, *EXCEPT* the following cases:
0228   //           Only {12, 34} pass, only {14, 23} pass.
0229   // Station 4 hits pass if any dTheta4X values pass, *EXCEPT* the following cases:
0230   //           Only {12, 34} pass, only {13, 24} pass.
0231   int vstat = vmask1;  // valid stations based on th coordinates
0232   if ((vstat & vmask2) != 0 || vstat == 0)
0233     vstat |= vmask2;
0234   if ((vstat & vmask3) != 0 || vstat == 0)
0235     vstat |= vmask3;
0236 
0237   // Truth table to remove ambiguity in passing the dTheta window cut when there are
0238   // two LCTs in the same station with the same phi value, but different theta values
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_) {  // Fixed at the beginning of 2018
0247     // construct bad delta word
0248     // dth_bad = {12,23,34,13,14,24}
0249     unsigned dth_bad = 0b111111;  // "1" is bad. if valid, change to "0" (good)
0250     if (best_dtheta_valid_arr_1.at(0)) {
0251       dth_bad &= (~(1 << 5));  // 12
0252     }
0253     if (best_dtheta_valid_arr_1.at(1)) {
0254       dth_bad &= (~(1 << 2));  // 13
0255     }
0256     if (best_dtheta_valid_arr_1.at(2)) {
0257       dth_bad &= (~(1 << 1));  // 14
0258     }
0259     if (best_dtheta_valid_arr_1.at(3)) {
0260       dth_bad &= (~(1 << 4));  // 23
0261     }
0262     if (best_dtheta_valid_arr_1.at(4)) {
0263       dth_bad &= (~(1 << 0));  // 24
0264     }
0265     if (best_dtheta_valid_arr_1.at(5)) {
0266       dth_bad &= (~(1 << 3));  // 34
0267     }
0268     emtf_assert(dth_bad < 64);
0269 
0270     // extract station mask from LUT
0271     vstat = trk_bld[dth_bad];
0272   }
0273 
0274   // remove valid flag for station if hit does not pass the dTheta mask
0275   for (int istation = 0; istation < emtf::NUM_STATIONS; ++istation) {
0276     if ((vstat & (1 << istation)) == 0) {  // station bit not set
0277       st_conv_hits.at(istation).clear();
0278     }
0279   }
0280 
0281   // assign precise phi and theta for the track
0282   int phi_fp = 0;
0283   int theta_fp = 0;
0284   int best_pair = -1;
0285 
0286   if ((vstat & (1 << 1)) != 0) {      // ME2 present
0287     if (best_dtheta_valid_arr.at(0))  // 12
0288       best_pair = 0;
0289     else if (best_dtheta_valid_arr.at(3))  // 23
0290       best_pair = 3;
0291     else if (best_dtheta_valid_arr.at(4))  // 24
0292       best_pair = 4;
0293   } else if ((vstat & (1 << 2)) != 0) {  // ME3 present
0294     if (best_dtheta_valid_arr.at(1))     // 13
0295       best_pair = 1;
0296     else if (best_dtheta_valid_arr.at(5))  // 34
0297       best_pair = 5;
0298   } else if ((vstat & (1 << 3)) != 0) {  // ME4 present
0299     if (best_dtheta_valid_arr.at(2))     // 14
0300       best_pair = 2;
0301   }
0302 
0303   // // Possible logic preferring CSC LCTs for the track theta and phi assignment
0304   // if ((vstat & (1<<1)) != 0) {            // ME2 present
0305   //   if (!best_has_rpc_arr.at(0) && best_dtheta_valid_arr.at(0))      // 12, no RPC
0306   //     best_pair = 0;
0307   //   else if (!best_has_rpc_arr.at(3) && best_dtheta_valid_arr.at(3)) // 23, no RPC
0308   //     best_pair = 3;
0309   //   else if (!best_has_rpc_arr.at(4) && best_dtheta_valid_arr.at(4)) // 24, no RPC
0310   //     best_pair = 4;
0311   //   else if (best_dtheta_valid_arr.at(0)) // 12, has RPC
0312   //     best_pair = 0;
0313   //   else if (best_dtheta_valid_arr.at(3)) // 23, has RPC
0314   //     best_pair = 3;
0315   //   else if (best_dtheta_valid_arr.at(4)) // 24, has RPC
0316   //     best_pair = 4;
0317   // } else if ((vstat & (1<<2)) != 0) {     // ME3 present
0318   //   if (!best_has_rpc_arr.at(1) && best_dtheta_valid_arr.at(1))      // 13, no RPC
0319   //     best_pair = 1;
0320   //   else if (!best_has_rpc_arr.at(5) && best_dtheta_valid_arr.at(5)) // 34, no RPC
0321   //     best_pair = 5;
0322   //   else if (best_dtheta_valid_arr.at(1)) // 13, has RPC
0323   //     best_pair = 1;
0324   //   else if (best_dtheta_valid_arr.at(5)) // 34, has RPC
0325   //     best_pair = 5;
0326   // } else if ((vstat & (1<<3)) != 0) {     // ME4 present
0327   //   if (best_dtheta_valid_arr.at(2))      // 14
0328   //     best_pair = 2;
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     // In firmware, the track is associated to LCTs by the segment number, which
0337     // identifies the best strip, but does not resolve the ambiguity in theta.
0338     // In emulator, this additional logic also resolves the ambiguity in theta.
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;  // capture
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);  // just the minimum in dtheta
0352     }
0353   }
0354 
0355   // update rank taking into account available stations after theta deltas
0356   // keep straightness as it was
0357   int old_rank = (track.Rank() << 1);  // output rank is one bit longer than input rank, to accomodate ME4 separately
0358   int rank = ((((old_rank >> 6) & 1) << 6) |  // straightness
0359               (((old_rank >> 4) & 1) << 4) |  // straightness
0360               (((old_rank >> 2) & 1) << 2) |  // straightness
0361               (((vstat >> 0) & 1) << 5) |     // ME1
0362               (((vstat >> 1) & 1) << 3) |     // ME2
0363               (((vstat >> 2) & 1) << 1) |     // ME3
0364               (((vstat >> 3) & 1) << 0)       // ME4
0365   );
0366 
0367   int mode = ((((vstat >> 0) & 1) << 3) |  // ME1
0368               (((vstat >> 1) & 1) << 2) |  // ME2
0369               (((vstat >> 2) & 1) << 1) |  // ME3
0370               (((vstat >> 3) & 1) << 0)    // ME4
0371   );
0372 
0373   int mode_inv = vstat;
0374 
0375   // if less than 2 segments, kill rank
0376   if (vstat == 0b0001 || vstat == 0b0010 || vstat == 0b0100 || vstat == 0b1000 || vstat == 0)
0377     rank = 0;
0378 
0379   // From RecoMuon/DetLayers/src/MuonCSCDetLayerGeometryBuilder.cc
0380   auto isFront = [](int station, int ring, int chamber, int subsystem) {
0381     // // RPCs are behind the CSCs in stations 1, 3, and 4; in front in 2
0382     // if (subsystem == TriggerPrimitive::kRPC)
0383     //   return (station == 2);
0384 
0385     // In EMTF firmware, RPC hits are treated as if they came from the corresponding
0386     // CSC chamber, so the FR bit assignment is the same as for CSCs - AWB 06.06.17
0387 
0388     // GEMs are in front of the CSCs
0389     if (subsystem == L1TMuon::kGEM)
0390       return true;
0391 
0392     bool result = false;
0393     bool isOverlapping = !(station == 1 && ring == 3);
0394     // not overlapping means back
0395     if (isOverlapping) {
0396       bool isEven = (chamber % 2 == 0);
0397       // odd chambers are bolted to the iron, which faces
0398       // forward in 1&2, backward in 3&4, so...
0399       result = (station < 3) ? isEven : !isEven;
0400     }
0401     return result;
0402   };
0403 
0404   // Fill ptlut_data
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();  // Automatically set to 0 for RPCs
0416     ptlut_data.csign[i] = v.empty() ? 0 : v.front().Bend();        // Automatically set to 0 for RPCs
0417     ptlut_data.slope[i] = v.empty() ? 0 : v.front().Slope();       // Automatically set to 0 for RPCs
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) {  // 'bt' arrays use 5-station convention
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   // Output
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   // Only keep the best segments
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;  // Count stubs delayed by i BX or more
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) {  // If 2 or more stubs are delayed by i BX or more
0485       second_bx = bx_ - i;     // if i == delayBX, analyze immediately
0486       break;
0487     }
0488   }
0489   emtf_assert(second_bx != 99);
0490 
0491   // ___________________________________________________________________________
0492   // Output
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   // Erase tracks with rank == 0
0500   // using erase-remove idiom
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   // Erase two-station tracks with hits in different BX
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_) {  // Modified at the beginning of 2018
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 }