Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-04-01 23:54:07

0001 #ifndef PHASE2GMT_TRACKMUONMATCHALGO
0002 #define PHASE2GMT_TRACKMUONMATCHALGO
0003 
0004 #include "DataFormats/L1TrackTrigger/interface/TTTrack.h"
0005 #include "DataFormats/L1TrackTrigger/interface/TTTrack_TrackWord.h"
0006 #include "DataFormats/L1TrackTrigger/interface/TTTypes.h"
0007 #include "DataFormats/L1TMuonPhase2/interface/MuonStub.h"
0008 #include "DataFormats/L1TMuonPhase2/interface/TrackerMuon.h"
0009 #include "DataFormats/L1TMuon/interface/RegionalMuonCand.h"
0010 #include "DataFormats/L1TMuon/interface/RegionalMuonCandFwd.h"
0011 #include "DataFormats/L1Trigger/interface/L1TObjComparison.h"
0012 #include "L1Trigger/Phase2L1GMT/interface/TrackConverter.h"
0013 #include "L1Trigger/Phase2L1GMT/interface/PreTrackMatchedMuon.h"
0014 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0015 #include "L1Trigger/Phase2L1GMT/interface/Constants.h"
0016 #include "L1Trigger/Phase2L1GMT/interface/MuonROI.h"
0017 #include "L1Trigger/Phase2L1GMT/interface/ConvertedTTTrack.h"
0018 #include <fstream>
0019 
0020 namespace Phase2L1GMT {
0021 
0022   const unsigned int PHIDIVIDER = 1 << (BITSPHI - BITSSTUBCOORD);
0023   const unsigned int ETADIVIDER = 1 << (BITSETA - BITSSTUBETA);
0024 
0025   typedef struct {
0026     ap_int<BITSSTUBCOORD> coord1;
0027     ap_uint<BITSSIGMACOORD> sigma_coord1;
0028     ap_int<BITSSTUBCOORD> coord2;
0029     ap_uint<BITSSIGMACOORD> sigma_coord2;
0030     ap_int<BITSSTUBETA> eta;
0031     ap_uint<BITSSIGMAETA> sigma_eta1;
0032     ap_uint<BITSSIGMAETA> sigma_eta2;
0033     ap_uint<1> valid;
0034     ap_uint<1> is_barrel;
0035   } propagation_t;
0036 
0037   typedef struct {
0038     ap_uint<BITSMATCHQUALITY - 2> quality;
0039     ap_uint<BITSSTUBID> id;
0040     ap_uint<1> valid;
0041     bool isGlobal;
0042     l1t::RegionalMuonCandRef muRef;
0043     l1t::MuonStubRef stubRef;
0044 
0045   } match_t;
0046 
0047   class TrackMuonMatchAlgorithm {
0048   public:
0049     TrackMuonMatchAlgorithm(const edm::ParameterSet& iConfig) : verbose_(iConfig.getParameter<int>("verbose")) {}
0050 
0051     ~TrackMuonMatchAlgorithm() {}
0052 
0053     std::vector<PreTrackMatchedMuon> processNonant(const std::vector<ConvertedTTTrack>& convertedTracks,
0054                                                    const std::vector<MuonROI>& rois) {
0055       std::vector<PreTrackMatchedMuon> preMuons;
0056       for (const auto& track : convertedTracks) {
0057         PreTrackMatchedMuon mu = processTrack(track, rois);
0058         if (mu.valid() && preMuons.size() < 16)
0059           preMuons.push_back(mu);
0060       }
0061       std::vector<PreTrackMatchedMuon> cleanedMuons = clean(preMuons);
0062       return cleanedMuons;
0063     }
0064 
0065     std::vector<PreTrackMatchedMuon> cleanNeighbor(const std::vector<PreTrackMatchedMuon>& muons,
0066                                                    const std::vector<PreTrackMatchedMuon>& muonsPrevious,
0067                                                    const std::vector<PreTrackMatchedMuon>& muonsNext,
0068                                                    bool equality) {
0069       std::vector<PreTrackMatchedMuon> out;
0070 
0071       if (muons.empty())
0072         return out;
0073 
0074       if (verbose_ == 1) {
0075         printf("-----Cleaning Up Muons in the neighbours\n");
0076         printf("Before:\n");
0077       }
0078 
0079       for (uint i = 0; i < muons.size(); ++i) {
0080         if (verbose_ == 1) {
0081           muons[i].print();
0082         }
0083         ap_uint<5> mask = 0x1f;
0084         for (uint j = 0; j < muonsPrevious.size(); ++j) {
0085           mask = mask & cleanMuon(muons[i], muonsPrevious[j], equality);
0086         }
0087         for (uint j = 0; j < muonsNext.size(); ++j) {
0088           mask = mask & cleanMuon(muons[i], muonsNext[j], equality);
0089         }
0090         if (mask) {
0091           if (verbose_ == 1)
0092             printf("kept\n");
0093           out.push_back(muons[i]);
0094         } else {
0095           if (verbose_ == 1)
0096             printf("discarded\n");
0097         }
0098       }
0099       return out;
0100     }
0101 
0102     std::vector<l1t::TrackerMuon> convert(std::vector<PreTrackMatchedMuon>& muons, uint maximum) {
0103       std::vector<l1t::TrackerMuon> out;
0104       for (const auto& mu : muons) {
0105         if (out.size() == maximum)
0106           break;
0107         l1t::TrackerMuon muon(mu.trkPtr(), mu.charge(), mu.pt(), mu.eta(), mu.phi(), mu.z0(), mu.d0(), mu.quality());
0108         //muon.setMuonRef(mu.muonRef());
0109         for (const auto& stub : mu.stubs())
0110           muon.addStub(stub);
0111         out.push_back(muon);
0112         if (verbose_ == 1) {
0113           printf("Final Muon:");
0114           muon.print();
0115         }
0116       }
0117       return out;
0118     }
0119 
0120     bool outputGT(std::vector<l1t::TrackerMuon>& muons) {
0121       for (auto& mu : muons) {
0122         wordtype word1 = 0;
0123         wordtype word2 = 0;
0124 
0125         int bstart = 0;
0126         bstart = wordconcat<wordtype>(word1, bstart, mu.hwPt(), BITSGTPT);
0127         bstart = wordconcat<wordtype>(word1, bstart, mu.hwPhi(), BITSGTPHI);
0128         bstart = wordconcat<wordtype>(word1, bstart, mu.hwEta(), BITSGTETA);
0129         bstart = wordconcat<wordtype>(word1, bstart, mu.hwZ0(), BITSGTZ0);
0130         bstart = wordconcat<wordtype>(word1, bstart, (mu.hwD0() >> 2), BITSGTD0);
0131 
0132         bstart = 0;
0133         bstart = wordconcat<wordtype>(word2, bstart, mu.hwCharge(), 1);
0134         bstart = wordconcat<wordtype>(word2, bstart, mu.hwQual(), BITSGTQUALITY);
0135         bstart = wordconcat<wordtype>(word2, bstart, mu.hwIso(), BITSGTISO);
0136         bstart = wordconcat<wordtype>(word2, bstart, mu.hwBeta(), BITSMUONBETA);
0137 
0138         std::array<uint64_t, 2> wordout = {{word1, word2}};
0139         mu.setWord(wordout);
0140       }
0141       return true;
0142     }
0143 
0144     std::vector<l1t::TrackerMuon> sort(std::vector<l1t::TrackerMuon>& muons, uint maximum) {
0145       if (muons.size() < 2)
0146         return muons;
0147 
0148       std::sort(muons.begin(), muons.end(), [](l1t::TrackerMuon a, l1t::TrackerMuon b) { return a.hwPt() > b.hwPt(); });
0149       std::vector<l1t::TrackerMuon> out;
0150       for (unsigned int i = 0; i < muons.size(); ++i) {
0151         out.push_back(muons[i]);
0152         if (i == (maximum - 1))
0153           break;
0154       }
0155 
0156       return out;
0157     }
0158 
0159   private:
0160     int verbose_;
0161 
0162     propagation_t propagate(const ConvertedTTTrack& track, uint layer) {
0163       ap_uint<BITSPROPCOORD> prop_coord1 = 0;
0164       ap_uint<BITSPROPCOORD> prop_coord2 = 0;
0165       ap_uint<BITSPROPSIGMACOORD_A> res0_coord1 = 0;
0166       ap_uint<BITSPROPSIGMACOORD_B> res1_coord1 = 0;
0167       ap_uint<BITSPROPSIGMACOORD_A> res0_coord2 = 0;
0168       ap_uint<BITSPROPSIGMACOORD_B> res1_coord2 = 0;
0169       ap_uint<BITSPROPSIGMAETA_A> res0_eta1 = 0;
0170       ap_uint<BITSPROPSIGMAETA_B> res1_eta = 0;
0171       ap_uint<BITSPROPSIGMAETA_A> res0_eta2 = 0;
0172       ap_uint<1> is_barrel = 0;
0173 
0174       uint reducedAbsEta = track.abseta() / 8;
0175 
0176       if (layer == 0) {
0177         prop_coord1 = lt_prop_coord1_0[reducedAbsEta];
0178         prop_coord2 = lt_prop_coord2_0[reducedAbsEta];
0179         res0_coord1 = lt_res0_coord1_0[reducedAbsEta];
0180         res1_coord1 = lt_res1_coord1_0[reducedAbsEta];
0181         res0_coord2 = lt_res0_coord2_0[reducedAbsEta];
0182         res1_coord2 = lt_res1_coord2_0[reducedAbsEta];
0183         res0_eta1 = lt_res0_eta1_0[reducedAbsEta];
0184         res1_eta = lt_res1_eta_0[reducedAbsEta];
0185         res0_eta2 = lt_res0_eta2_0[reducedAbsEta];
0186         is_barrel = reducedAbsEta < barrelLimit0_ ? 1 : 0;
0187       } else if (layer == 1) {
0188         prop_coord1 = lt_prop_coord1_1[reducedAbsEta];
0189         prop_coord2 = lt_prop_coord2_1[reducedAbsEta];
0190         res0_coord1 = lt_res0_coord1_1[reducedAbsEta];
0191         res1_coord1 = lt_res1_coord1_1[reducedAbsEta];
0192         res0_coord2 = lt_res0_coord2_1[reducedAbsEta];
0193         res1_coord2 = lt_res1_coord2_1[reducedAbsEta];
0194         res0_eta1 = lt_res0_eta1_1[reducedAbsEta];
0195         res1_eta = lt_res1_eta_1[reducedAbsEta];
0196         res0_eta2 = lt_res0_eta2_1[reducedAbsEta];
0197         is_barrel = reducedAbsEta < barrelLimit1_ ? 1 : 0;
0198 
0199       } else if (layer == 2) {
0200         prop_coord1 = lt_prop_coord1_2[reducedAbsEta];
0201         prop_coord2 = lt_prop_coord2_2[reducedAbsEta];
0202         res0_coord1 = lt_res0_coord1_2[reducedAbsEta];
0203         res1_coord1 = lt_res1_coord1_2[reducedAbsEta];
0204         res0_coord2 = lt_res0_coord2_2[reducedAbsEta];
0205         res1_coord2 = lt_res1_coord2_2[reducedAbsEta];
0206         res0_eta1 = lt_res0_eta1_2[reducedAbsEta];
0207         res1_eta = lt_res1_eta_2[reducedAbsEta];
0208         res0_eta2 = lt_res0_eta2_2[reducedAbsEta];
0209         is_barrel = reducedAbsEta < barrelLimit2_ ? 1 : 0;
0210 
0211       } else if (layer == 3) {
0212         prop_coord1 = lt_prop_coord1_3[reducedAbsEta];
0213         prop_coord2 = lt_prop_coord2_3[reducedAbsEta];
0214         res0_coord1 = lt_res0_coord1_3[reducedAbsEta];
0215         res1_coord1 = lt_res1_coord1_3[reducedAbsEta];
0216         res0_coord2 = lt_res0_coord2_3[reducedAbsEta];
0217         res1_coord2 = lt_res1_coord2_3[reducedAbsEta];
0218         res0_eta1 = lt_res0_eta1_3[reducedAbsEta];
0219         res1_eta = lt_res1_eta_3[reducedAbsEta];
0220         res0_eta2 = lt_res0_eta2_3[reducedAbsEta];
0221         is_barrel = reducedAbsEta < barrelLimit3_ ? 1 : 0;
0222 
0223       } else if (layer == 4) {
0224         prop_coord1 = lt_prop_coord1_4[reducedAbsEta];
0225         prop_coord2 = lt_prop_coord2_4[reducedAbsEta];
0226         res0_coord1 = lt_res0_coord1_4[reducedAbsEta];
0227         res1_coord1 = lt_res1_coord1_4[reducedAbsEta];
0228         res0_coord2 = lt_res0_coord2_4[reducedAbsEta];
0229         res1_coord2 = lt_res1_coord2_4[reducedAbsEta];
0230         res0_eta1 = lt_res0_eta1_4[reducedAbsEta];
0231         res1_eta = lt_res1_eta_4[reducedAbsEta];
0232         res0_eta2 = lt_res0_eta2_4[reducedAbsEta];
0233         is_barrel = 0;
0234       }
0235 
0236       propagation_t out;
0237       ap_int<BITSTTCURV> curvature = track.curvature();
0238       ap_int<BITSPHI> phi = track.phi();
0239       ap_int<BITSPROPCOORD + BITSTTCURV> c1kFull = prop_coord1 * curvature;
0240       ap_int<BITSPROPCOORD + BITSTTCURV - 10> c1k = (c1kFull) / 1024;
0241       ap_int<BITSPHI> coord1 = phi - c1k;
0242 
0243       out.coord1 = coord1 / PHIDIVIDER;
0244 
0245       ap_int<BITSPROPCOORD + BITSTTCURV> c2kFull = prop_coord2 * curvature;
0246 
0247       ap_int<BITSPROPCOORD + BITSTTCURV - 10> c2k = (c2kFull) / 1024;
0248       if (is_barrel)
0249         out.coord2 = -c2k / PHIDIVIDER;
0250       else
0251         out.coord2 = (phi - c2k) / PHIDIVIDER;
0252 
0253       ap_int<BITSETA> eta = track.eta();
0254       out.eta = eta / ETADIVIDER;
0255 
0256       ap_uint<2 * BITSTTCURV - 2> curvature2All = curvature * curvature;
0257       ap_uint<BITSTTCURV2> curvature2 = curvature2All / 2;
0258 
0259       //Remember to change emulator with new k2
0260       ap_uint<BITSPROPSIGMACOORD_B + BITSTTCURV2> rescoord1k = (res1_coord1 * curvature2) >> 23;
0261       ap_ufixed<BITSSIGMACOORD, BITSSIGMACOORD, AP_TRN_ZERO, AP_SAT_SYM> sigma_coord1 = res0_coord1 + rescoord1k;
0262       out.sigma_coord1 = ap_uint<BITSSIGMACOORD>(sigma_coord1);
0263 
0264       ap_uint<BITSPROPSIGMACOORD_B + BITSTTCURV2> rescoord2k = (res1_coord2 * curvature2) >> 23;
0265       ap_ufixed<BITSSIGMACOORD, BITSSIGMACOORD, AP_TRN_ZERO, AP_SAT_SYM> sigma_coord2 = res0_coord2 + rescoord2k;
0266       out.sigma_coord2 = ap_uint<BITSSIGMACOORD>(sigma_coord2);
0267 
0268       ap_uint<BITSPROPSIGMAETA_B + BITSTTCURV2> resetak = (res1_eta * curvature2) >> 23;
0269       ap_ufixed<BITSSIGMAETA, BITSSIGMAETA, AP_TRN_ZERO, AP_SAT_SYM> sigma_eta1 = res0_eta1 + resetak;
0270       out.sigma_eta1 = ap_uint<BITSSIGMAETA>(sigma_eta1);
0271       ap_ufixed<BITSSIGMAETA, BITSSIGMAETA, AP_TRN_ZERO, AP_SAT_SYM> sigma_eta2 = res0_eta2 + resetak;
0272       out.sigma_eta2 = ap_uint<BITSSIGMAETA>(sigma_eta2);
0273       out.valid = 1;
0274       out.is_barrel = is_barrel;
0275 
0276       if (verbose_ == 1)
0277 
0278         printf("Propagating to layer %d:is barrel=%d  coords=%d+-%d , %d +-%d etas = %d +- %d +-%d\n",
0279                int(layer),
0280                out.is_barrel.to_int(),
0281                out.coord1.to_int(),
0282                out.sigma_coord1.to_int(),
0283                out.coord2.to_int(),
0284                out.sigma_coord2.to_int(),
0285                out.eta.to_int(),
0286                out.sigma_eta1.to_int(),
0287                out.sigma_eta2.to_int());
0288 
0289       return out;
0290     }
0291 
0292     ap_uint<BITSSIGMAETA + 1> deltaEta(const ap_int<BITSSTUBETA>& eta1, const ap_int<BITSSTUBETA>& eta2) {
0293       ap_fixed<BITSSIGMAETA + 2, BITSSIGMAETA + 2, AP_TRN_ZERO, AP_SAT_SYM> dEta = eta1 - eta2;
0294       if (dEta < 0)
0295         return ap_uint<BITSSIGMAETA + 1>(-dEta);
0296       else
0297         return ap_uint<BITSSIGMAETA + 1>(dEta);
0298     }
0299 
0300     ap_uint<BITSSIGMACOORD + 1> deltaCoord(const ap_int<BITSSTUBCOORD>& phi1, const ap_int<BITSSTUBCOORD>& phi2) {
0301       ap_int<BITSSTUBCOORD> dPhiRoll = phi1 - phi2;
0302       ap_ufixed<BITSSIGMACOORD + 1, BITSSIGMACOORD + 1, AP_TRN_ZERO, AP_SAT_SYM> dPhi;
0303       if (dPhiRoll < 0)
0304         dPhi = ap_ufixed<BITSSIGMACOORD + 1, BITSSIGMACOORD + 1, AP_TRN_ZERO, AP_SAT_SYM>(-dPhiRoll);
0305       else
0306         dPhi = ap_ufixed<BITSSIGMACOORD + 1, BITSSIGMACOORD + 1, AP_TRN_ZERO, AP_SAT_SYM>(dPhiRoll);
0307 
0308       return ap_uint<BITSSIGMACOORD + 1>(dPhi);
0309     }
0310 
0311     match_t match(const propagation_t prop, const l1t::MuonStubRef& stub) {
0312       if (verbose_ == 1) {
0313         printf("Matching to ");
0314         stub->print();
0315       }
0316       //Matching of Coord1
0317       ap_uint<1> coord1Matched;
0318       ap_uint<BITSSIGMACOORD + 1> deltaCoord1 = deltaCoord(prop.coord1, stub->coord1());
0319       if (deltaCoord1 <= prop.sigma_coord1 && (stub->quality() & 0x1)) {
0320         coord1Matched = 1;
0321       } else {
0322         coord1Matched = 0;
0323       }
0324       if (verbose_ == 1)
0325         printf("Coord1 matched=%d delta=%d res=%d\n",
0326                coord1Matched.to_int(),
0327                deltaCoord1.to_int(),
0328                prop.sigma_coord1.to_int());
0329 
0330       //Matching of Coord2
0331       ap_uint<1> coord2Matched;
0332       ap_uint<BITSSIGMACOORD + 1> deltaCoord2 = deltaCoord(prop.coord2, stub->coord2());
0333       if (deltaCoord2 <= prop.sigma_coord2 && (stub->quality() & 0x2)) {
0334         coord2Matched = 1;
0335       } else {
0336         coord2Matched = 0;
0337       }
0338       if (verbose_ == 1)
0339         printf("Coord2 matched=%d delta=%d res=%d\n",
0340                coord2Matched.to_int(),
0341                deltaCoord2.to_int(),
0342                prop.sigma_coord2.to_int());
0343 
0344       //Matching of Eta1
0345 
0346       ap_uint<1> eta1Matched;
0347 
0348       //if we have really bad quality[Barrel no eta]
0349       //increase the resolution
0350       ap_ufixed<BITSSIGMAETA, BITSSIGMAETA, AP_TRN_ZERO, AP_SAT_SYM> prop_sigma_eta1;
0351       if (stub->etaQuality() == 0)
0352         prop_sigma_eta1 = prop.sigma_eta1 + 6;
0353       else
0354         prop_sigma_eta1 = prop.sigma_eta1;
0355 
0356       ap_uint<BITSSIGMAETA + 1> deltaEta1 = deltaEta(prop.eta, stub->eta1());
0357       if (deltaEta1 <= prop_sigma_eta1 && (stub->etaQuality() == 0 || (stub->etaQuality() & 0x1)))
0358         eta1Matched = 1;
0359       else
0360         eta1Matched = 0;
0361 
0362       if (verbose_ == 1)
0363         printf("eta1 matched=%d delta=%d res=%d\n", eta1Matched.to_int(), deltaEta1.to_int(), prop_sigma_eta1.to_int());
0364 
0365       //Matching of Eta2
0366 
0367       ap_uint<1> eta2Matched;
0368 
0369       ap_uint<BITSSIGMAETA + 1> deltaEta2 = deltaEta(prop.eta, stub->eta2());
0370       if (deltaEta2 <= prop.sigma_eta2 && (stub->etaQuality() & 0x2))
0371         eta2Matched = 1;
0372       else
0373         eta2Matched = 0;
0374       match_t out;
0375       out.id = stub->id();
0376 
0377       if (verbose_ == 1)
0378         printf("eta2 matched=%d delta=%d res=%d\n", eta2Matched.to_int(), deltaEta2.to_int(), prop.sigma_eta2.to_int());
0379 
0380       //if barrel, coord1 has to always be matched, coord2 maybe and eta1 is needed if etaQ=0 or then the one that depends on eta quality
0381       if (prop.is_barrel) {
0382         out.valid = (coord1Matched == 1 && (eta1Matched == 1 || eta2Matched == 1));
0383         if (out.valid == 0) {
0384           out.quality = 0;
0385         } else {
0386           out.quality = 32 - deltaCoord1;
0387           if (coord2Matched == 1)
0388             out.quality += 32 - deltaCoord2;
0389         }
0390       }
0391       //if endcap each coordinate is independent except the case where phiQuality=1 and etaQuality==3
0392       else {
0393         bool match1 = (coord1Matched == 1 && eta1Matched == 1);
0394         bool match2 = (coord2Matched == 1 && eta2Matched == 1);
0395         bool match3 =
0396             (coord1Matched == 1 && (eta1Matched || eta2Matched) && stub->etaQuality() == 3 && stub->quality() == 1);
0397         out.valid = match1 || match2 || match3;
0398         if (out.valid == 0)
0399           out.quality = 0;
0400         else {
0401           out.quality = 0;
0402           if (match1 || match3)
0403             out.quality += 32 - deltaCoord1;
0404           if (match2)
0405             out.quality += 32 - deltaCoord2;
0406         }
0407       }
0408       if (verbose_ == 1)
0409         printf("GlobalMatchQuality = %d\n", out.quality.to_int());
0410       out.stubRef = stub;
0411       return out;
0412     }
0413 
0414     match_t propagateAndMatch(const ConvertedTTTrack& track, const l1t::MuonStubRef& stub) {
0415       propagation_t prop = propagate(track, stub->tfLayer());
0416       return match(prop, stub);
0417     }
0418 
0419     match_t getBest(const std::vector<match_t> matches) {
0420       match_t best = matches[0];
0421       for (const auto& m : matches) {
0422         if (m.quality > best.quality)
0423           best = m;
0424       }
0425 
0426       return best;
0427     }
0428 
0429     PreTrackMatchedMuon processTrack(const ConvertedTTTrack& track, const std::vector<MuonROI>& rois) {
0430       std::vector<match_t> matchInfo0;
0431       std::vector<match_t> matchInfo1;
0432       std::vector<match_t> matchInfo2;
0433       std::vector<match_t> matchInfo3;
0434       std::vector<match_t> matchInfo4;
0435 
0436       if (verbose_ == 1 && !rois.empty()) {
0437         printf("-----------processing new track----------");
0438         track.print();
0439       }
0440       for (const auto& roi : rois) {
0441         if (verbose_ == 1) {
0442           printf("New ROI with %d stubs \n", int(roi.stubs().size()));
0443         }
0444         for (const auto& stub : roi.stubs()) {
0445           match_t m = propagateAndMatch(track, stub);
0446           if (m.valid == 1) {
0447             if (roi.isGlobalMuon() && roi.muonRef().isNonnull()) {
0448               m.isGlobal = true;
0449               m.muRef = roi.muonRef();
0450             }
0451 
0452             if (stub->tfLayer() == 0)
0453               matchInfo0.push_back(m);
0454             else if (stub->tfLayer() == 1)
0455               matchInfo1.push_back(m);
0456             else if (stub->tfLayer() == 2)
0457               matchInfo2.push_back(m);
0458             else if (stub->tfLayer() == 3)
0459               matchInfo3.push_back(m);
0460             else if (stub->tfLayer() == 4)
0461               matchInfo4.push_back(m);
0462           }
0463         }
0464       }
0465 
0466       ap_ufixed<6, 6, AP_TRN_ZERO, AP_SAT_SYM> ptPenalty = ap_ufixed<6, 6, AP_TRN_ZERO, AP_SAT_SYM>(track.pt() / 32);
0467 
0468       ap_uint<BITSMATCHQUALITY> quality = 0;
0469       PreTrackMatchedMuon muon(track.charge(), track.pt(), track.eta(), track.phi(), track.z0(), track.d0());
0470 
0471       if (!matchInfo0.empty()) {
0472         match_t b = getBest(matchInfo0);
0473         if (b.valid) {
0474           muon.addStub(b.stubRef);
0475           if (b.isGlobal)
0476             muon.addMuonRef(b.muRef);
0477           quality += b.quality;
0478         }
0479       }
0480       if (!matchInfo1.empty()) {
0481         match_t b = getBest(matchInfo1);
0482         if (b.valid) {
0483           muon.addStub(b.stubRef);
0484           if (b.isGlobal)
0485             muon.addMuonRef(b.muRef);
0486           quality += b.quality;
0487         }
0488       }
0489       if (!matchInfo2.empty()) {
0490         match_t b = getBest(matchInfo2);
0491         if (b.valid) {
0492           muon.addStub(b.stubRef);
0493           if (b.isGlobal)
0494             muon.addMuonRef(b.muRef);
0495           quality += b.quality;
0496         }
0497       }
0498       if (!matchInfo3.empty()) {
0499         match_t b = getBest(matchInfo3);
0500         if (b.valid) {
0501           muon.addStub(b.stubRef);
0502           if (b.isGlobal)
0503             muon.addMuonRef(b.muRef);
0504           quality += b.quality;
0505         }
0506       }
0507       if (!matchInfo4.empty()) {
0508         match_t b = getBest(matchInfo4);
0509         if (b.valid) {
0510           muon.addStub(b.stubRef);
0511           if (b.isGlobal)
0512             muon.addMuonRef(b.muRef);
0513           quality += b.quality;
0514         }
0515       }
0516 
0517       muon.setOfflineQuantities(track.offline_pt(), track.offline_eta(), track.offline_phi());
0518       muon.setTrkPtr(track.trkPtr());
0519 
0520       ap_uint<8> etaAddr = muon.eta() < 0 ? ap_uint<8>(-muon.eta() / 256) : ap_uint<8>((muon.eta()) / 256);
0521       ap_uint<8> ptAddr = muon.pt() > 4095 ? ap_uint<8>(15) : ap_uint<8>(muon.pt() / 256);
0522       ap_uint<8> addr = ptAddr | (etaAddr << 4);
0523       ap_uint<8> qualityCut = lt_tpsID[addr];
0524 
0525       if (quality >= qualityCut) {
0526         muon.setValid(true);
0527         muon.setQuality(quality + ptPenalty);
0528       } else {
0529         muon.setValid(false);
0530         muon.setQuality(0);
0531         muon.resetGlobal();
0532       }
0533       if (verbose_ == 1)
0534         muon.print();
0535 
0536       if (verbose_ == 1 && !rois.empty()) {  //patterns for HLS
0537 
0538         printf("TPS %d", track.trkPtr()->phiSector());
0539         track.printWord();
0540 
0541         for (uint i = 0; i < 16; ++i) {
0542           if (rois.size() > i) {
0543             rois[i].printROILine();
0544           } else {
0545             printf("%08x", 0);
0546             printf("%016lx", 0x1ff000000000000);
0547             printf("%016lx", 0x1ff000000000000);
0548             printf("%016lx", 0x1ff000000000000);
0549             printf("%016lx", 0x1ff000000000000);
0550             printf("%016lx", 0x1ff000000000000);
0551           }
0552         }
0553         muon.printWord();
0554         printf("\n");
0555       }
0556       return muon;
0557     }
0558 
0559     ap_uint<5> cleanMuon(const PreTrackMatchedMuon& mu, const PreTrackMatchedMuon& other, bool eq) {
0560       ap_uint<5> valid = 0;
0561       ap_uint<5> overlap = 0;
0562       if (mu.stubID0() != 511) {
0563         valid = valid | 0x1;
0564         if (mu.stubID0() == other.stubID0())
0565           overlap = overlap | 0x1;
0566       }
0567       if (mu.stubID1() != 511) {
0568         valid = valid | 0x2;
0569         if (mu.stubID1() == other.stubID1())
0570           overlap = overlap | 0x2;
0571       }
0572       if (mu.stubID2() != 511) {
0573         valid = valid | 0x4;
0574         if (mu.stubID2() == other.stubID2())
0575           overlap = overlap | 0x4;
0576       }
0577       if (mu.stubID3() != 511) {
0578         valid = valid | 0x8;
0579         if (mu.stubID3() == other.stubID3())
0580           overlap = overlap | 0x8;
0581       }
0582       if (mu.stubID4() != 511) {
0583         valid = valid | 0x10;
0584         if (mu.stubID4() == other.stubID4())
0585           overlap = overlap | 0x10;
0586       }
0587 
0588       if (((mu.quality() < other.quality()) && (!eq)) || ((mu.quality() <= other.quality()) && (eq)))
0589         return valid & (~overlap);
0590       else
0591         return valid;
0592     }
0593 
0594     std::vector<PreTrackMatchedMuon> clean(std::vector<PreTrackMatchedMuon>& muons) {
0595       std::vector<PreTrackMatchedMuon> out;
0596       if (muons.empty())
0597         return out;
0598       if (verbose_ == 1) {
0599         printf("-----Cleaning Up Muons in the same Nonant\n");
0600         printf("Before:\n");
0601       }
0602       for (uint i = 0; i < muons.size(); ++i) {
0603         if (verbose_ == 1)
0604           muons[i].print();
0605 
0606         ap_uint<5> mask = 0x1f;
0607         for (uint j = 0; j < muons.size(); ++j) {
0608           if (i == j)
0609             continue;
0610           mask = mask & cleanMuon(muons[i], muons[j], false);
0611         }
0612         if (mask) {
0613           if (verbose_ == 1)
0614             printf("kept\n");
0615           out.push_back(muons[i]);
0616         } else {
0617           if (verbose_ == 1)
0618             printf("discarded\n");
0619         }
0620       }
0621       return out;
0622     }
0623   };
0624 }  // namespace Phase2L1GMT
0625 
0626 #endif