L1TMuonEndCap/src/PatternRecognition.cc

0001 #include "L1Trigger/L1TMuonEndCap/interface/PatternRecognition.h"
0002 #include "DataFormats/L1TMuon/interface/L1TMuonSubsystems.h"
0003
0004 #include "helper.h"  // to_hex, to_binary
0005
0006 namespace {
0007   const int padding_w_st1 = 15;
0008   const int padding_w_st3 = 7;
0009   const int padding_extra_w_st1 = padding_w_st1 - padding_w_st3;
0010 }  // namespace
0011
0012 void PatternRecognition::configure(int verbose,
0013                                    int endcap,
0014                                    int sector,
0015                                    int bx,
0016                                    int bxWindow,
0017                                    const std::vector<std::string>& pattDefinitions,
0018                                    const std::vector<std::string>& symPattDefinitions,
0019                                    bool useSymPatterns,
0020                                    int maxRoadsPerZone,
0021                                    bool useSecondEarliest) {
0022   verbose_ = verbose;
0023   endcap_ = endcap;
0024   sector_ = sector;
0025   bx_ = bx;
0026
0027   bxWindow_ = bxWindow;
0028   pattDefinitions_ = pattDefinitions;
0029   symPattDefinitions_ = symPattDefinitions;
0030   useSymPatterns_ = useSymPatterns;
0031   maxRoadsPerZone_ = maxRoadsPerZone;
0032   useSecondEarliest_ = useSecondEarliest;
0033
0034   configure_details();
0035 }
0036
0037 void PatternRecognition::configure_details() {
0038   emtf_assert(1 <= bxWindow_ && bxWindow_ <= 3);  // only work for BX windows <= 3
0039
0040   patterns_.clear();
0041
0042   // Parse pattern definitions
0043   if (!useSymPatterns_) {
0044     // Normal patterns
0045     for (const auto& s : pattDefinitions_) {
0046       const std::vector<std::string>& tokens = split_string(s, ',', ':');  // split by comma or colon
0047       emtf_assert(tokens.size() == 9);                                     // want to find 9 numbers
0048
0049       std::vector<std::string>::const_iterator tokens_it = tokens.begin();
0050
0051       // Get the 9 integers
0052       // straightness, hits in ME1, hits in ME2, hits in ME3, hits in ME4
0053       int straightness = std::stoi(*tokens_it++);
0054       int st1_max = std::stoi(*tokens_it++);
0055       int st1_min = std::stoi(*tokens_it++);
0056       int st2_max = std::stoi(*tokens_it++);
0057       int st2_min = std::stoi(*tokens_it++);
0058       int st3_max = std::stoi(*tokens_it++);
0059       int st3_min = std::stoi(*tokens_it++);
0060       int st4_max = std::stoi(*tokens_it++);
0061       int st4_min = std::stoi(*tokens_it++);
0062
0063       // There can only be one zone hit in the key station in the pattern
0064       // and it has to be this magic number
0065       emtf_assert(st2_max == padding_w_st3 && st2_min == padding_w_st3);
0066
0067       // There is extra "padding" in st1 w.r.t st2,3,4
0068       // Add the extra padding to st2,3,4
0069       st2_max += padding_extra_w_st1;
0070       st2_min += padding_extra_w_st1;
0071       st3_max += padding_extra_w_st1;
0072       st3_min += padding_extra_w_st1;
0073       st4_max += padding_extra_w_st1;
0074       st4_min += padding_extra_w_st1;
0075
0076       // Create a pattern
0077       PhiMemoryImage pattern;
0078       pattern.set_straightness(straightness);
0079       int i = 0;
0080
0081       for (i = st1_min; i <= st1_max; i++)
0082         pattern.set_bit(0, i);
0083       for (i = st2_min; i <= st2_max; i++)
0084         pattern.set_bit(1, i);
0085       for (i = st3_min; i <= st3_max; i++)
0086         pattern.set_bit(2, i);
0087       for (i = st4_min; i <= st4_max; i++)
0088         pattern.set_bit(3, i);
0089
0090       // Remove the extra padding
0091       pattern.rotr(padding_extra_w_st1);
0092       patterns_.push_back(pattern);
0093     }
0094     emtf_assert(patterns_.size() == pattDefinitions_.size());
0095
0096   } else {
0097     // Symmetrical patterns
0098     for (const auto& s : symPattDefinitions_) {
0099       const std::vector<std::string>& tokens = split_string(s, ',', ':');  // split by comma or colon
0100       emtf_assert(tokens.size() == 17);                                    // want to find 17 numbers
0101
0102       std::vector<std::string>::const_iterator tokens_it = tokens.begin();
0103
0104       // Get the 17 integers
0105       // straightness, hits in ME1, hits in ME2, hits in ME3, hits in ME4
0106       int straightness = std::stoi(*tokens_it++);
0107       int st1_max1 = std::stoi(*tokens_it++);
0108       int st1_min1 = std::stoi(*tokens_it++);
0109       int st1_max2 = std::stoi(*tokens_it++);
0110       int st1_min2 = std::stoi(*tokens_it++);
0111       int st2_max1 = std::stoi(*tokens_it++);
0112       int st2_min1 = std::stoi(*tokens_it++);
0113       int st2_max2 = std::stoi(*tokens_it++);
0114       int st2_min2 = std::stoi(*tokens_it++);
0115       int st3_max1 = std::stoi(*tokens_it++);
0116       int st3_min1 = std::stoi(*tokens_it++);
0117       int st3_max2 = std::stoi(*tokens_it++);
0118       int st3_min2 = std::stoi(*tokens_it++);
0119       int st4_max1 = std::stoi(*tokens_it++);
0120       int st4_min1 = std::stoi(*tokens_it++);
0121       int st4_max2 = std::stoi(*tokens_it++);
0122       int st4_min2 = std::stoi(*tokens_it++);
0123
0124       // There can only be one zone hit in the key station in the pattern
0125       // and it has to be this magic number
0126       emtf_assert(st2_max1 == padding_w_st3 && st2_min1 == padding_w_st3);
0127       emtf_assert(st2_max2 == padding_w_st3 && st2_min2 == padding_w_st3);
0128
0129       // There is extra "padding" in st1 w.r.t st2,3,4
0130       // Add the extra padding to st2,3,4
0131       st2_max1 += padding_extra_w_st1;
0132       st2_min1 += padding_extra_w_st1;
0133       st2_max2 += padding_extra_w_st1;
0134       st2_min2 += padding_extra_w_st1;
0135       st3_max1 += padding_extra_w_st1;
0136       st3_min1 += padding_extra_w_st1;
0137       st3_max2 += padding_extra_w_st1;
0138       st3_min2 += padding_extra_w_st1;
0139       st4_max1 += padding_extra_w_st1;
0140       st4_min1 += padding_extra_w_st1;
0141       st4_max2 += padding_extra_w_st1;
0142       st4_min2 += padding_extra_w_st1;
0143
0144       // Create a pattern
0145       PhiMemoryImage pattern;
0146       pattern.set_straightness(straightness);
0147       int i = 0;
0148
0149       for (i = st1_min1; i <= st1_max1; i++)
0150         pattern.set_bit(0, i);
0151       for (i = st1_min2; i <= st1_max2; i++)
0152         pattern.set_bit(0, i);
0153       for (i = st2_min1; i <= st2_max1; i++)
0154         pattern.set_bit(1, i);
0155       for (i = st2_min2; i <= st2_max2; i++)
0156         pattern.set_bit(1, i);
0157       for (i = st3_min1; i <= st3_max1; i++)
0158         pattern.set_bit(2, i);
0159       for (i = st3_min2; i <= st3_max2; i++)
0160         pattern.set_bit(2, i);
0161       for (i = st4_min1; i <= st4_max1; i++)
0162         pattern.set_bit(3, i);
0163       for (i = st4_min2; i <= st4_max2; i++)
0164         pattern.set_bit(3, i);
0165
0166       // Remove the extra padding
0167       pattern.rotr(padding_extra_w_st1);
0168       patterns_.push_back(pattern);
0169     }
0170     emtf_assert(patterns_.size() == symPattDefinitions_.size());
0171   }
0172
0173   if (verbose_ > 2) {  // debug
0174     for (const auto& pattern : patterns_) {
0175       std::cout << "Pattern straightness: " << pattern.get_straightness() << " image: " << std::endl;
0176       std::cout << pattern << std::endl;
0177     }
0178   }
0179 }
0180
0181 void PatternRecognition::process(const std::deque<EMTFHitCollection>& extended_conv_hits,
0182                                  std::map<pattern_ref_t, int>& patt_lifetime_map,
0183                                  emtf::zone_array<EMTFRoadCollection>& zone_roads) const {
0184   // Exit if no hits
0185   int num_conv_hits = 0;
0186   for (const auto& conv_hits : extended_conv_hits)
0187     num_conv_hits += conv_hits.size();
0188   bool early_exit = (num_conv_hits == 0) && (patt_lifetime_map.empty());
0189
0190   if (early_exit)
0191     return;
0192
0193   if (verbose_ > 0) {  // debug
0194     for (const auto& conv_hits : extended_conv_hits) {
0195       for (const auto& conv_hit : conv_hits) {
0196         if (conv_hit.Subsystem() == L1TMuon::kCSC) {
0197           std::cout << "CSC hit st: " << conv_hit.PC_station() << " ch: " << conv_hit.PC_chamber()
0198                     << " ph: " << conv_hit.Phi_fp() << " th: " << conv_hit.Theta_fp() << " strip: " << conv_hit.Strip()
0199                     << " wire: " << conv_hit.Wire() << " cpat: " << conv_hit.Pattern()
0200                     << " zone_hit: " << conv_hit.Zone_hit() << " zone_code: " << conv_hit.Zone_code()
0201                     << " bx: " << conv_hit.BX() << std::endl;
0202         }
0203       }
0204     }
0205
0206     for (const auto& conv_hits : extended_conv_hits) {
0207       for (const auto& conv_hit : conv_hits) {
0208         if (conv_hit.Subsystem() == L1TMuon::kRPC) {
0209           std::cout << "RPC hit st: " << conv_hit.PC_station() << " ch: " << conv_hit.PC_chamber()
0210                     << " ph>>2: " << (conv_hit.Phi_fp() >> 2) << " th>>2: " << (conv_hit.Theta_fp() >> 2)
0211                     << " strip: " << conv_hit.Strip() << " roll: " << conv_hit.Roll() << " cpat: " << conv_hit.Pattern()
0212                     << " bx: " << conv_hit.BX() << std::endl;
0213         }
0214       }
0215     }
0216
0217     for (const auto& conv_hits : extended_conv_hits) {
0218       for (const auto& conv_hit : conv_hits) {
0219         if (conv_hit.Subsystem() == L1TMuon::kGEM) {
0220           std::cout << "GEM hit st: " << conv_hit.PC_station() << " ch: " << conv_hit.PC_chamber()
0221                     << " ph: " << conv_hit.Phi_fp() << " th: " << conv_hit.Theta_fp() << " strip: " << conv_hit.Strip()
0222                     << " roll: " << conv_hit.Roll() << " cpat: " << conv_hit.Pattern() << " bx: " << conv_hit.BX()
0223                     << std::endl;
0224         }
0225       }
0226     }
0227   }  // end debug
0228
0229   // Perform pattern recognition in each zone
0230   emtf::zone_array<PhiMemoryImage> zone_images;
0231
0232   for (int izone = 0; izone < emtf::NUM_ZONES; ++izone) {
0233     // Skip the zone if no hits and no patterns
0234     if (is_zone_empty(izone + 1, extended_conv_hits, patt_lifetime_map))
0235       continue;
0236
0237     // Make zone images
0238     make_zone_image(izone + 1, extended_conv_hits, zone_images.at(izone));
0239
0240     // Detect patterns
0241     process_single_zone(izone + 1, zone_images.at(izone), patt_lifetime_map, zone_roads.at(izone));
0242   }
0243
0244   if (verbose_ > 2) {  // debug
0245     for (int izone = emtf::NUM_ZONES; izone >= 1; --izone) {
0246       std::cout << "zone: " << izone << std::endl;
0247       std::cout << zone_images.at(izone - 1) << std::endl;
0248     }
0249     //for (const auto& kv : patt_lifetime_map) {
0250     //  std::cout << "zone: " << kv.first.at(0) << " izhit: " << kv.first.at(1) << " ipatt: " << kv.first.at(2) << " lifetime: " << kv.second << std::endl;
0251     //}
0252   }
0253
0254   if (verbose_ > 0) {  // debug
0255     for (const auto& roads : zone_roads) {
0256       for (const auto& road : reversed(roads)) {
0257         std::cout << "pattern: z: " << road.Zone() - 1 << " ph: " << road.Key_zhit()
0258                   << " q: " << to_hex(road.Quality_code()) << " ly: " << to_binary(road.Layer_code(), 3)
0259                   << " str: " << to_binary(road.Straightness(), 3) << " bx: " << road.BX() << std::endl;
0260       }
0261     }
0262   }
0263
0264   // Sort patterns and select best three patterns in each zone
0265   for (int izone = 0; izone < emtf::NUM_ZONES; ++izone) {
0266     sort_single_zone(zone_roads.at(izone));
0267   }
0268 }
0269
0270 bool PatternRecognition::is_zone_empty(int zone,
0271                                        const std::deque<EMTFHitCollection>& extended_conv_hits,
0272                                        const std::map<pattern_ref_t, int>& patt_lifetime_map) const {
0273   int izone = zone - 1;
0274   int num_conv_hits = 0;
0275   int num_patts = 0;
0276
0277   std::deque<EMTFHitCollection>::const_iterator ext_conv_hits_it = extended_conv_hits.begin();
0278   std::deque<EMTFHitCollection>::const_iterator ext_conv_hits_end = extended_conv_hits.end();
0279
0280   for (; ext_conv_hits_it != ext_conv_hits_end; ++ext_conv_hits_it) {
0281     EMTFHitCollection::const_iterator conv_hits_it = ext_conv_hits_it->begin();
0282     EMTFHitCollection::const_iterator conv_hits_end = ext_conv_hits_it->end();
0283
0284     for (; conv_hits_it != conv_hits_end; ++conv_hits_it) {
0285       emtf_assert(conv_hits_it->PC_segment() <=
0286                   4);  // With 2 unique LCTs per chamber, 4 possible strip/wire combinations
0287
0288       if (conv_hits_it->Subsystem() == L1TMuon::kRPC)
0289         continue;  // Don't use RPC hits for pattern formation
0290
0291       if (conv_hits_it->Subsystem() == L1TMuon::kGEM)
0292         continue;  // Don't use GEM hits for pattern formation
0293
0294       if (conv_hits_it->Zone_code() & (1 << izone)) {  // hit belongs to this zone
0295         num_conv_hits += 1;
0296       }
0297     }  // end loop over conv_hits
0298   }    // end loop over extended_conv_hits
0299
0300   std::map<pattern_ref_t, int>::const_iterator patt_lifetime_map_it = patt_lifetime_map.begin();
0301   std::map<pattern_ref_t, int>::const_iterator patt_lifetime_map_end = patt_lifetime_map.end();
0302
0303   for (; patt_lifetime_map_it != patt_lifetime_map_end; ++patt_lifetime_map_it) {
0304     if (patt_lifetime_map_it->first.at(0) == zone) {
0305       num_patts += 1;
0306     }
0307   }  // end loop over patt_lifetime_map
0308
0309   return (num_conv_hits == 0) && (num_patts == 0);
0310 }
0311
0312 void PatternRecognition::make_zone_image(int zone,
0313                                          const std::deque<EMTFHitCollection>& extended_conv_hits,
0314                                          PhiMemoryImage& image) const {
0315   int izone = zone - 1;
0316
0317   std::deque<EMTFHitCollection>::const_iterator ext_conv_hits_it = extended_conv_hits.begin();
0318   std::deque<EMTFHitCollection>::const_iterator ext_conv_hits_end = extended_conv_hits.end();
0319
0320   for (; ext_conv_hits_it != ext_conv_hits_end; ++ext_conv_hits_it) {
0321     EMTFHitCollection::const_iterator conv_hits_it = ext_conv_hits_it->begin();
0322     EMTFHitCollection::const_iterator conv_hits_end = ext_conv_hits_it->end();
0323
0324     for (; conv_hits_it != conv_hits_end; ++conv_hits_it) {
0325       if (conv_hits_it->Subsystem() == L1TMuon::kRPC)
0326         continue;  // Don't use RPC hits for pattern formation
0327
0328       if (conv_hits_it->Subsystem() == L1TMuon::kGEM)
0329         continue;  // Don't use GEM hits for pattern formation
0330
0331       if (conv_hits_it->Zone_code() & (1 << izone)) {  // hit belongs to this zone
0332         unsigned int layer = conv_hits_it->Station() - 1;
0333         unsigned int bit = conv_hits_it->Zone_hit();
0334         image.set_bit(layer, bit);
0335       }
0336     }  // end loop over conv_hits
0337   }    // end loop over extended_conv_hits
0338 }
0339
0340 void PatternRecognition::process_single_zone(int zone,
0341                                              PhiMemoryImage cloned_image,
0342                                              std::map<pattern_ref_t, int>& patt_lifetime_map,
0343                                              EMTFRoadCollection& roads) const {
0344   roads.clear();
0345
0346   const int drift_time = bxWindow_ - 1;
0347   const int npatterns = patterns_.size();
0348
0349   // The zone hit image is rotated/shifted before comparing with patterns
0350   // First, rotate left/shift up to the zone hit in the key station
0351   cloned_image.rotl(padding_w_st3);
0352
0353   for (int izhit = 0; izhit < emtf::NUM_ZONE_HITS; ++izhit) {
0354     // The zone hit image is rotated/shift before comparing with patterns
0355     // For every zone hit, rotate right/shift down by one
0356     if (izhit > 0)
0357       cloned_image.rotr(1);
0358
0359     int max_quality_code = -1;
0360     EMTFRoad tmp_road;
0361
0362     // Compare with patterns
0363     for (int ipatt = 0; ipatt < npatterns; ++ipatt) {
0364       const PhiMemoryImage& patt = patterns_.at(ipatt);
0365       const pattern_ref_t patt_ref = {{zone, izhit, ipatt}};  // due to GCC bug, use {{}} instead of {}
0366       int straightness = patt.get_straightness();
0367
0368       bool is_lifetime_up = false;
0369
0370       int layer_code = patt.op_and(cloned_image);  // kind of like AND operator
0371       bool more_than_one = (layer_code != 0) && (layer_code != 1) && (layer_code != 2) && (layer_code != 4);
0372       bool more_than_zero = (layer_code != 0);
0373
0374       if (more_than_zero) {
0375         // Insert this pattern
0376         auto ins = patt_lifetime_map.insert({patt_ref, 0});
0377
0378         if (!useSecondEarliest_) {
0379           // Use earliest
0380           auto patt_ins = ins.first;  // iterator of patt_lifetime_map pointing to this pattern
0381           bool patt_exists = !ins.second;
0382
0383           if (patt_exists) {                       // if exists
0384             if (patt_ins->second == drift_time) {  // is lifetime up?
0385               is_lifetime_up = true;
0386             }
0387           }
0388           patt_ins->second += 1;  // bx starts counting at any hit in the pattern, even single
0389
0390         } else {
0391           // Use 2nd earliest
0392           auto patt_ins = ins.first;  // iterator of patt_lifetime_map pointing to this pattern
0393           int bx_shifter = patt_ins->second;
0394           int bx2 = bool(bx_shifter & (1 << 2));
0395           int bx1 = bool(bx_shifter & (1 << 1));
0396           int bx0 = bool(bx_shifter & (1 << 0));
0397
0398           // is lifetime up?
0399           if (drift_time == 2 && bx2 == 0 && bx1 == 1) {
0400             is_lifetime_up = true;
0401           } else if (drift_time == 1 && bx1 == 0 && bx0 == 1) {
0402             is_lifetime_up = true;
0403           } else if (drift_time == 0) {
0404             is_lifetime_up = true;
0405           }
0406
0407           bx2 = bx1;
0408           bx1 = bx0;
0409           bx0 = more_than_zero;  // put 1 in shifter when one layer is hit
0410           bx_shifter = (bx2 << 2) | (bx1 << 1) | bx0;
0411           patt_ins->second = bx_shifter;
0412         }
0413
0414       } else {
0415         // Zero hit
0416         patt_lifetime_map.erase(patt_ref);  // erase if exists
0417       }
0418
0419       // If lifetime is up, and not single-layer hit patterns (stations 3&4 considered
0420       // as a single layer), find quality of this pattern
0421       if (is_lifetime_up && more_than_one) {
0422         // This quality code scheme is giving almost-equal priority to
0423         // more stations and better straightness
0424         // Station 1 has higher weight, station 2 lower, stations 3&4 lowest
0425         int quality_code =
0426             ((((straightness >> 2) & 1) << 5) | (((straightness >> 1) & 1) << 3) | (((straightness >> 0) & 1) << 1) |
0427              (((layer_code >> 2) & 1) << 4) | (((layer_code >> 1) & 1) << 2) | (((layer_code >> 0) & 1) << 0));
0428
0429         // Create a road (fired pattern)
0430         EMTFRoad road;
0431         road.set_endcap((endcap_ == 1) ? 1 : -1);
0432         road.set_sector(sector_);
0433         road.set_sector_idx((endcap_ == 1) ? sector_ - 1 : sector_ + 5);
0434         road.set_bx(bx_ - drift_time);
0435
0436         road.set_zone(patt_ref.at(0));
0437         road.set_key_zhit(patt_ref.at(1));
0438         road.set_pattern(patt_ref.at(2));
0439
0440         road.set_straightness(straightness);
0441         road.set_layer_code(layer_code);
0442         road.set_quality_code(quality_code);
0443
0444         // Find max quality code in a given key_zhit
0445         if (max_quality_code < quality_code) {
0446           max_quality_code = quality_code;
0447           tmp_road = road;
0448         }
0449       }  // end if is_lifetime_up
0450
0451     }  // end loop over patterns
0452
0453     // Output road
0454     if (max_quality_code != -1) {
0455       roads.push_back(tmp_road);
0456     }
0457
0458   }  // end loop over zone hits
0459
0460   // Ghost cancellation logic by considering neighbor patterns
0461
0462   if (!roads.empty()) {
0463     std::array<int, emtf::NUM_ZONE_HITS> quality_codes;
0464     quality_codes.fill(0);
0465
0466     EMTFRoadCollection::iterator roads_it = roads.begin();
0467     EMTFRoadCollection::iterator roads_end = roads.end();
0468
0469     for (; roads_it != roads_end; ++roads_it) {
0470       quality_codes.at(roads_it->Key_zhit()) = roads_it->Quality_code();
0471     }
0472
0473     roads_it = roads.begin();
0474     roads_end = roads.end();
0475
0476     for (; roads_it != roads_end; ++roads_it) {
0477       int izhit = roads_it->Key_zhit();
0478
0479       // Center quality is the current one
0480       int qc = quality_codes.at(izhit);
0481
0482       // Left and right qualities are the neighbors
0483       // Protect against the right end and left end special cases
0484       int ql = (izhit == emtf::NUM_ZONE_HITS - 1) ? 0 : quality_codes.at(izhit + 1);
0485       int qr = (izhit == 0) ? 0 : quality_codes.at(izhit - 1);
0486
0487       // Cancellation conditions
0488       if (qc <= ql || qc < qr) {        // this pattern is lower quality than neighbors
0489         roads_it->set_quality_code(0);  // cancel
0490       }
0491     }
0492   }
0493
0494   // Erase roads with quality_code == 0
0495   // using erase-remove idiom
0496   struct {
0497     typedef EMTFRoad value_type;
0498     bool operator()(const value_type& x) const { return (x.Quality_code() == 0); }
0499   } quality_code_zero_pred;
0500
0501   roads.erase(std::remove_if(roads.begin(), roads.end(), quality_code_zero_pred), roads.end());
0502 }
0503
0504 void PatternRecognition::sort_single_zone(EMTFRoadCollection& roads) const {
0505   // First, order by key_zhit (highest to lowest)
0506   struct {
0507     typedef EMTFRoad value_type;
0508     bool operator()(const value_type& lhs, const value_type& rhs) const { return lhs.Key_zhit() > rhs.Key_zhit(); }
0509   } greater_zhit_cmp;
0510
0511   std::sort(roads.begin(), roads.end(), greater_zhit_cmp);
0512
0513   // Second, sort by quality_code (highest to lowest), but preserving the original order if qualities are equal
0514   struct {
0515     typedef EMTFRoad value_type;
0516     bool operator()(const value_type& lhs, const value_type& rhs) const {
0517       return lhs.Quality_code() > rhs.Quality_code();
0518     }
0519   } greater_quality_cmp;
0520
0521   std::stable_sort(roads.begin(), roads.end(), greater_quality_cmp);
0522
0523   // Finally, select 3 best
0524   const size_t n = maxRoadsPerZone_;
0525   if (roads.size() > n) {
0526     roads.erase(roads.begin() + n, roads.end());
0527   }
0528   emtf_assert(roads.size() <= n);
0529
0530   // Assign the winner variable
0531   for (unsigned iroad = 0; iroad < roads.size(); ++iroad) {
0532     roads.at(iroad).set_winner(iroad);
0533   }
0534 }