Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 12:20:53

0001 #include "L1Trigger/L1TMuonBarrel/interface/L1TMuonBarrelKalmanSectorProcessor.h"
0002 
0003 L1TMuonBarrelKalmanSectorProcessor::L1TMuonBarrelKalmanSectorProcessor(const edm::ParameterSet& iConfig, int sector)
0004     : verbose_(iConfig.getParameter<int>("verbose")), sector_(sector) {
0005   std::vector<int> wheels = iConfig.getParameter<std::vector<int> >("wheelsToProcess");
0006   for (const auto wheel : wheels)
0007     regions_.push_back(
0008         L1TMuonBarrelKalmanRegionModule(iConfig.getParameter<edm::ParameterSet>("regionSettings"), wheel, sector));
0009 }
0010 
0011 L1TMuonBarrelKalmanSectorProcessor::~L1TMuonBarrelKalmanSectorProcessor() {}
0012 
0013 L1MuKBMTrackCollection L1TMuonBarrelKalmanSectorProcessor::process(L1TMuonBarrelKalmanAlgo* trackMaker,
0014                                                                    const L1MuKBMTCombinedStubRefVector& stubsAll,
0015                                                                    int bx) {
0016   L1MuKBMTrackCollection tracksM2;
0017   L1MuKBMTrackCollection tracksM1;
0018   L1MuKBMTrackCollection tracks0;
0019   L1MuKBMTrackCollection tracksP1;
0020   L1MuKBMTrackCollection tracksP2;
0021 
0022   for (auto& region : regions_) {
0023     L1MuKBMTrackCollection tmp = region.process(trackMaker, stubsAll, bx);
0024     if (region.wheel() == -2)
0025       tracksM2.insert(tracksM2.end(), tmp.begin(), tmp.end());
0026     if (region.wheel() == -1)
0027       tracksM1.insert(tracksM1.end(), tmp.begin(), tmp.end());
0028     if (region.wheel() == 0)
0029       tracks0.insert(tracks0.end(), tmp.begin(), tmp.end());
0030     if (region.wheel() == 1)
0031       tracksP1.insert(tracksP1.end(), tmp.begin(), tmp.end());
0032     if (region.wheel() == 2)
0033       tracksP2.insert(tracksP2.end(), tmp.begin(), tmp.end());
0034   }
0035 
0036   L1MuKBMTrackCollection out = wedgeSort(tracksM2, tracksM1, tracks0, tracksP1, tracksP2);
0037   if (verbose_ == 1)
0038     verbose(trackMaker, out);
0039 
0040   return out;
0041 }
0042 
0043 L1TMuonBarrelKalmanSectorProcessor::bmtf_out L1TMuonBarrelKalmanSectorProcessor::makeWord(
0044     L1TMuonBarrelKalmanAlgo* trackMaker, const L1MuKBMTrackCollection& tracks) {
0045   L1TMuonBarrelKalmanSectorProcessor::bmtf_out out;
0046   out.pt_1 = 0;
0047   out.qual_1 = 0;
0048   out.eta_1 = 0;
0049   out.HF_1 = 0;
0050   out.phi_1 = 0;
0051   out.bx0_1 = 0;
0052   out.charge_1 = 0;
0053   out.chargeValid_1 = 0;
0054   out.dxy_1 = 0;
0055   out.addr1_1 = 3;
0056   out.addr2_1 = 15;
0057   out.addr3_1 = 15;
0058   out.addr4_1 = 15;
0059   out.reserved_1 = 0;
0060   out.wheel_1 = 0;
0061   out.ptSTA_1 = 0;
0062   out.SE_1 = 0;
0063 
0064   out.pt_2 = 0;
0065   out.qual_2 = 0;
0066   out.eta_2 = 0;
0067   out.HF_2 = 0;
0068   out.phi_2 = 0;
0069   out.bx0_2 = 0;
0070   out.charge_2 = 0;
0071   out.chargeValid_2 = 0;
0072   out.dxy_2 = 0;
0073   out.addr1_2 = 3;
0074   out.addr2_2 = 15;
0075   out.addr3_2 = 15;
0076   out.addr4_2 = 15;
0077   out.reserved_2 = 0;
0078   out.wheel_2 = 0;
0079   out.ptSTA_2 = 0;
0080   out.SE_2 = 0;
0081 
0082   out.pt_3 = 0;
0083   out.qual_3 = 0;
0084   out.eta_3 = 0;
0085   out.HF_3 = 0;
0086   out.phi_3 = 0;
0087   out.bx0_3 = 0;
0088   out.charge_3 = 0;
0089   out.chargeValid_3 = 0;
0090   out.dxy_3 = 0;
0091   out.addr1_3 = 3;
0092   out.addr2_3 = 15;
0093   out.addr3_3 = 15;
0094   out.addr4_3 = 15;
0095   out.reserved_3 = 0;
0096   out.wheel_3 = 0;
0097   out.ptSTA_3 = 0;
0098   out.SE_3 = 0;
0099 
0100   if (!tracks.empty()) {
0101     l1t::RegionalMuonCand mu = trackMaker->convertToBMTF(tracks[0]);
0102     out.pt_1 = mu.hwPt();
0103     out.qual_1 = mu.hwQual();
0104     out.eta_1 = mu.hwEta();
0105     out.HF_1 = mu.hwHF();
0106     out.phi_1 = mu.hwPhi();
0107     out.bx0_1 = 0;
0108     out.charge_1 = mu.hwSign();
0109     out.chargeValid_1 = mu.hwSignValid();
0110     out.dxy_1 = mu.hwDXY();
0111     out.addr1_1 = mu.trackSubAddress(l1t::RegionalMuonCand::kStat1);
0112     out.addr2_1 = mu.trackSubAddress(l1t::RegionalMuonCand::kStat2);
0113     out.addr3_1 = mu.trackSubAddress(l1t::RegionalMuonCand::kStat3);
0114     out.addr4_1 = mu.trackSubAddress(l1t::RegionalMuonCand::kStat4);
0115     out.wheel_1 = mu.trackSubAddress(l1t::RegionalMuonCand::kWheelSide) * (1 << 2) +
0116                   mu.trackSubAddress(l1t::RegionalMuonCand::kWheelNum);
0117     out.ptSTA_1 = mu.hwPtUnconstrained();
0118   }
0119 
0120   if (tracks.size() > 1) {
0121     l1t::RegionalMuonCand mu = trackMaker->convertToBMTF(tracks[1]);
0122     out.pt_2 = mu.hwPt();
0123     out.qual_2 = mu.hwQual();
0124     out.eta_2 = mu.hwEta();
0125     out.HF_2 = mu.hwHF();
0126     out.phi_2 = mu.hwPhi();
0127     out.bx0_2 = 0;
0128     out.charge_2 = mu.hwSign();
0129     out.chargeValid_2 = mu.hwSignValid();
0130     out.dxy_2 = mu.hwDXY();
0131     out.addr1_2 = mu.trackSubAddress(l1t::RegionalMuonCand::kStat1);
0132     out.addr2_2 = mu.trackSubAddress(l1t::RegionalMuonCand::kStat2);
0133     out.addr3_2 = mu.trackSubAddress(l1t::RegionalMuonCand::kStat3);
0134     out.addr4_2 = mu.trackSubAddress(l1t::RegionalMuonCand::kStat4);
0135     out.wheel_2 = mu.trackSubAddress(l1t::RegionalMuonCand::kWheelSide) * (1 << 2) +
0136                   mu.trackSubAddress(l1t::RegionalMuonCand::kWheelNum);
0137 
0138     out.ptSTA_2 = mu.hwPtUnconstrained();
0139   }
0140 
0141   if (tracks.size() > 2) {
0142     l1t::RegionalMuonCand mu = trackMaker->convertToBMTF(tracks[2]);
0143     out.pt_3 = mu.hwPt();
0144     out.qual_3 = mu.hwQual();
0145     out.eta_3 = mu.hwEta();
0146     out.HF_3 = mu.hwHF();
0147     out.phi_3 = mu.hwPhi();
0148     out.bx0_3 = 0;
0149     out.charge_3 = mu.hwSign();
0150     out.chargeValid_3 = mu.hwSignValid();
0151     out.dxy_3 = mu.hwDXY();
0152     out.addr1_3 = mu.trackSubAddress(l1t::RegionalMuonCand::kStat1);
0153     out.addr2_3 = mu.trackSubAddress(l1t::RegionalMuonCand::kStat2);
0154     out.addr3_3 = mu.trackSubAddress(l1t::RegionalMuonCand::kStat3);
0155     out.addr4_3 = mu.trackSubAddress(l1t::RegionalMuonCand::kStat4);
0156     out.wheel_3 = mu.trackSubAddress(l1t::RegionalMuonCand::kWheelSide) * (1 << 2) +
0157                   mu.trackSubAddress(l1t::RegionalMuonCand::kWheelNum);
0158     out.ptSTA_3 = mu.hwPtUnconstrained();
0159   }
0160   return out;
0161 }
0162 
0163 void L1TMuonBarrelKalmanSectorProcessor::verbose(L1TMuonBarrelKalmanAlgo* algo, const L1MuKBMTrackCollection& tracks) {
0164   L1TMuonBarrelKalmanSectorProcessor::bmtf_out out = makeWord(algo, tracks);
0165   std::cout << "O " << sector_ << " " << out.pt_1 << " " << out.qual_1 << " " << out.eta_1 << " " << out.HF_1 << " "
0166             << out.phi_1 << " " << out.charge_1 << " " << out.chargeValid_1 << " " << out.dxy_1 << " " << out.addr1_1
0167             << " " << out.addr2_1 << " " << out.addr3_1 << " " << out.addr4_1 << " " << out.wheel_1 << " "
0168             << out.ptSTA_1 << " " << out.pt_2 << " " << out.qual_2 << " " << out.eta_2 << " " << out.HF_2 << " "
0169             << out.phi_2 << " " << out.charge_2 << " " << out.chargeValid_2 << " " << out.dxy_2 << " " << out.addr1_2
0170             << " " << out.addr2_2 << " " << out.addr3_2 << " " << out.addr4_2 << " " << out.wheel_2 << " "
0171             << out.ptSTA_2 << " " << out.pt_3 << " " << out.qual_3 << " " << out.eta_3 << " " << out.HF_3 << " "
0172             << out.phi_3 << " " << out.charge_3 << " " << out.chargeValid_3 << " " << out.dxy_3 << " " << out.addr1_3
0173             << " " << out.addr2_3 << " " << out.addr3_3 << " " << out.addr4_3 << " " << out.wheel_3 << " "
0174             << out.ptSTA_3 << std::endl;
0175 }
0176 
0177 // L1MuKBMTrackCollection L1TMuonBarrelKalmanSectorProcessor::cleanAndSort(const L1MuKBMTrackCollection& pretracks,uint keep) {
0178 //   L1MuKBMTrackCollection out;
0179 
0180 //   if (verbose_)
0181 //     printf" -----Preselected Kalman Tracks for sector %d----- \n",sector_);
0182 
0183 //   for(const auto& track1 : pretracks) {
0184 //     if (verbose_)
0185 //       printf("Pre Track charge=%d pt=%f eta=%f phi=%f curvature=%d curvature STA =%d stubs=%d bitmask=%d rank=%d chi=%d pts=%f %f\n",track1.charge(),track1.pt(),track1.eta(),track1.phi(),track1.curvatureAtVertex(),track1.curvatureAtMuon(),int(track1.stubs().size()),track1.hitPattern(),track1.rank(),track1.approxChi2(),track1.pt(),track1.ptUnconstrained());
0186 
0187 //     bool keep=true;
0188 //     for(const auto& track2 : pretracks) {
0189 //       if (track1==track2)
0190 //  continue;
0191 //       if (!track1.overlapTrack(track2))
0192 //  continue;
0193 
0194 //       if (track1.rank()<track2.rank())
0195 //  keep=false;
0196 
0197 //       if ((track1.rank()==track2.rank()) && (fabs(track1.stubs()[0]->whNum())<fabs(track2.stubs()[0]->whNum())))
0198 //  keep=false;
0199 //     }
0200 //     if (keep)
0201 //       out.push_back(track1);
0202 //   }
0203 
0204 //   if (verbose_) {
0205 //   printf(" -----Algo Result Kalman Tracks-----\n");
0206 //   for (const auto& track1 :out)
0207 //     printf("Final Kalman Track charge=%d pt=%f eta=%f phi=%f curvature=%d curvature STA =%d stubs=%d chi2=%d pts=%f %f\n",track1.charge(),track1.pt(),track1.eta(),track1.phi(),track1.curvatureAtVertex(),track1.curvatureAtMuon(),int(track1.stubs().size()),track1.approxChi2(),track1.pt(),track1.ptUnconstrained());
0208 //   }
0209 
0210 //   TrackSorter sorter;
0211 //   if (!out.empty())
0212 //     std::sort(out.begin(),out.end(),sorter);
0213 
0214 //   L1MuKBMTrackCollection exported;
0215 //   for (uint i=0;i<out.size();++i)
0216 //     if (i<=keep)
0217 //       exported.push_back(out[i]);
0218 //   return exported;
0219 // }
0220 
0221 L1MuKBMTrackCollection L1TMuonBarrelKalmanSectorProcessor::cleanNeighbor(const L1MuKBMTrackCollection& coll1,
0222                                                                          const L1MuKBMTrackCollection& coll2) {
0223   L1MuKBMTrackCollection out;
0224 
0225   for (const auto& track1 : coll1) {
0226     /*
0227     if (verbose_)
0228       printf(
0229           "Pre Track charge=%d pt=%f eta=%f phi=%f curvature=%d curvature STA =%d stubs=%d bitmask=%d rank=%d chi=%d "
0230           "pts=%f %f\n",
0231           track1.charge(),
0232           track1.pt(),
0233           track1.eta(),
0234           track1.phi(),
0235           track1.curvatureAtVertex(),
0236           track1.curvatureAtMuon(),
0237           int(track1.stubs().size()),
0238           track1.hitPattern(),
0239           track1.rank(),
0240           track1.approxChi2(),
0241           track1.pt(),
0242           track1.ptUnconstrained());
0243     */
0244     bool keep = true;
0245     for (const auto& track2 : coll2) {
0246       if (!track1.overlapTrack(track2))
0247         continue;
0248 
0249       if (track1.rank() < track2.rank())
0250         keep = false;
0251 
0252       if ((track1.rank() == track2.rank()) && (fabs(track1.stubs()[0]->whNum()) < fabs(track2.stubs()[0]->whNum())))
0253         keep = false;
0254     }
0255     if (keep)
0256       out.push_back(track1);
0257     else {
0258       L1MuKBMTrack temp = track1;
0259       temp.setPtEtaPhi(0, 0, 0);
0260       temp.setRank(0);
0261       out.push_back(temp);
0262     }
0263   }
0264 
0265   return out;
0266 }
0267 
0268 L1MuKBMTrackCollection L1TMuonBarrelKalmanSectorProcessor::cleanNeighbors(const L1MuKBMTrackCollection& coll1,
0269                                                                           const L1MuKBMTrackCollection& coll2,
0270                                                                           const L1MuKBMTrackCollection& coll3) {
0271   L1MuKBMTrackCollection out;
0272 
0273   for (const auto& track1 : coll1) {
0274     // if (verbose_)
0275     //   printf("Pre Track charge=%d pt=%f eta=%f phi=%f curvature=%d curvature STA =%d stubs=%d bitmask=%d rank=%d chi=%d pts=%f %f\n",track1.charge(),track1.pt(),track1.eta(),track1.phi(),track1.curvatureAtVertex(),track1.curvatureAtMuon(),int(track1.stubs().size()),track1.hitPattern(),track1.rank(),track1.approxChi2(),track1.pt(),track1.ptUnconstrained());
0276 
0277     bool keep = true;
0278     for (const auto& track2 : coll2) {
0279       if (!track1.overlapTrack(track2))
0280         continue;
0281 
0282       if (track1.rank() < track2.rank())
0283         keep = false;
0284 
0285       if ((track1.rank() == track2.rank()) && (fabs(track1.stubs()[0]->whNum()) < fabs(track2.stubs()[0]->whNum())))
0286         keep = false;
0287     }
0288 
0289     for (const auto& track2 : coll3) {
0290       if (!track1.overlapTrack(track2))
0291         continue;
0292 
0293       if (track1.rank() < track2.rank())
0294         keep = false;
0295 
0296       if ((track1.rank() == track2.rank()) && (fabs(track1.stubs()[0]->whNum()) < fabs(track2.stubs()[0]->whNum())))
0297         keep = false;
0298     }
0299 
0300     if (keep)
0301       out.push_back(track1);
0302     else {
0303       L1MuKBMTrack temp = track1;
0304       temp.setPtEtaPhi(0, 0, 0);
0305       temp.setRank(0);
0306       out.push_back(temp);
0307     }
0308   }
0309 
0310   return out;
0311 }
0312 
0313 void swap(std::map<uint, double>& info, std::map<uint, L1MuKBMTrack>& trackInfo, uint i, uint j) {
0314   double pt1 = info[i];
0315   double pt2 = info[j];
0316 
0317   if (pt2 > pt1) {
0318     info[j] = pt1;
0319     info[i] = pt2;
0320 
0321     if (pt1 != 0.0 && pt2 != 0.0) {
0322       L1MuKBMTrack tmp = trackInfo[i];
0323       trackInfo[i] = trackInfo[j];
0324       trackInfo[j] = tmp;
0325     }
0326     if (pt2 != 0.0 && pt1 == 0.0) {
0327       trackInfo[i] = trackInfo[j];
0328     }
0329   }
0330 }
0331 
0332 L1MuKBMTrackCollection L1TMuonBarrelKalmanSectorProcessor::wedgeSort(const L1MuKBMTrackCollection& preminus2,
0333                                                                      const L1MuKBMTrackCollection& preminus1,
0334                                                                      const L1MuKBMTrackCollection& prezero,
0335                                                                      const L1MuKBMTrackCollection& preone,
0336                                                                      const L1MuKBMTrackCollection& pretwo) {
0337   //first clean
0338 
0339   L1MuKBMTrackCollection minus2 = cleanNeighbor(preminus2, preminus1);
0340   L1MuKBMTrackCollection minus1 = cleanNeighbors(preminus1, preminus2, prezero);
0341   L1MuKBMTrackCollection zero = cleanNeighbors(prezero, preminus1, preone);
0342   L1MuKBMTrackCollection plus1 = cleanNeighbors(preone, prezero, pretwo);
0343   L1MuKBMTrackCollection plus2 = cleanNeighbor(pretwo, preone);
0344 
0345   std::map<uint, L1MuKBMTrack> trackInfo;
0346   std::map<uint, double> ptInfo;
0347 
0348   for (uint i = 0; i < 10; ++i)
0349     ptInfo[i] = 0.0;
0350 
0351   if (minus2.size() > 1) {
0352     ptInfo[0] = minus2[0].pt();
0353     trackInfo[0] = minus2[0];
0354     ptInfo[1] = minus2[1].pt();
0355     trackInfo[1] = minus2[1];
0356   } else if (minus2.size() == 1) {
0357     if (minus2[0].stubs()[0]->tag()) {
0358       ptInfo[1] = minus2[0].pt();
0359       trackInfo[1] = minus2[0];
0360     } else {
0361       ptInfo[0] = minus2[0].pt();
0362       trackInfo[0] = minus2[0];
0363     }
0364   }
0365 
0366   if (minus1.size() > 1) {
0367     ptInfo[2] = minus1[0].pt();
0368     trackInfo[2] = minus1[0];
0369     ptInfo[3] = minus1[1].pt();
0370     trackInfo[3] = minus1[1];
0371   } else if (minus1.size() == 1) {
0372     if (minus1[0].stubs()[0]->tag()) {
0373       ptInfo[3] = minus1[0].pt();
0374       trackInfo[3] = minus1[0];
0375     } else {
0376       ptInfo[2] = minus1[0].pt();
0377       trackInfo[2] = minus1[0];
0378     }
0379   }
0380 
0381   if (zero.size() > 1) {
0382     ptInfo[4] = zero[0].pt();
0383     trackInfo[4] = zero[0];
0384     ptInfo[5] = zero[1].pt();
0385     trackInfo[5] = zero[1];
0386   } else if (zero.size() == 1) {
0387     if (zero[0].stubs()[0]->tag()) {
0388       ptInfo[5] = zero[0].pt();
0389       trackInfo[5] = zero[0];
0390     } else {
0391       ptInfo[4] = zero[0].pt();
0392       trackInfo[4] = zero[0];
0393     }
0394   }
0395 
0396   if (plus1.size() > 1) {
0397     ptInfo[6] = plus1[0].pt();
0398     trackInfo[6] = plus1[0];
0399     ptInfo[7] = plus1[1].pt();
0400     trackInfo[7] = plus1[1];
0401   } else if (plus1.size() == 1) {
0402     if (plus1[0].stubs()[0]->tag()) {
0403       ptInfo[7] = plus1[0].pt();
0404       trackInfo[7] = plus1[0];
0405     } else {
0406       ptInfo[6] = plus1[0].pt();
0407       trackInfo[6] = plus1[0];
0408     }
0409   }
0410 
0411   if (plus2.size() > 1) {
0412     ptInfo[8] = plus2[0].pt();
0413     trackInfo[8] = plus2[0];
0414     ptInfo[9] = plus2[1].pt();
0415     trackInfo[9] = plus2[1];
0416   } else if (plus2.size() == 1) {
0417     if (plus2[0].stubs()[0]->tag()) {
0418       ptInfo[9] = plus2[0].pt();
0419       trackInfo[9] = plus2[0];
0420     } else {
0421       ptInfo[8] = plus2[0].pt();
0422       trackInfo[8] = plus2[0];
0423     }
0424   }
0425 
0426   //Now my glorious partial bitonic != power of two sorter
0427   swap(ptInfo, trackInfo, 0, 5);
0428   swap(ptInfo, trackInfo, 1, 6);
0429   swap(ptInfo, trackInfo, 2, 7);
0430   swap(ptInfo, trackInfo, 3, 8);
0431   swap(ptInfo, trackInfo, 4, 9);
0432 
0433   swap(ptInfo, trackInfo, 0, 3);
0434   swap(ptInfo, trackInfo, 1, 4);
0435   swap(ptInfo, trackInfo, 5, 8);
0436   swap(ptInfo, trackInfo, 6, 9);
0437 
0438   swap(ptInfo, trackInfo, 0, 2);
0439   swap(ptInfo, trackInfo, 3, 6);
0440   swap(ptInfo, trackInfo, 7, 9);
0441 
0442   swap(ptInfo, trackInfo, 0, 1);
0443   swap(ptInfo, trackInfo, 2, 4);
0444   swap(ptInfo, trackInfo, 5, 7);
0445   swap(ptInfo, trackInfo, 8, 9);
0446 
0447   swap(ptInfo, trackInfo, 1, 2);
0448   swap(ptInfo, trackInfo, 3, 5);
0449   swap(ptInfo, trackInfo, 4, 6);
0450   swap(ptInfo, trackInfo, 7, 8);
0451 
0452   swap(ptInfo, trackInfo, 1, 3);
0453   swap(ptInfo, trackInfo, 2, 5);
0454   swap(ptInfo, trackInfo, 4, 7);
0455   swap(ptInfo, trackInfo, 6, 8);
0456 
0457   swap(ptInfo, trackInfo, 2, 3);
0458 
0459   L1MuKBMTrackCollection out;
0460   if (ptInfo[0] != 0.0)
0461     out.push_back(trackInfo[0]);
0462   if (ptInfo[1] != 0.0)
0463     out.push_back(trackInfo[1]);
0464   if (ptInfo[2] != 0.0)
0465     out.push_back(trackInfo[2]);
0466 
0467   return out;
0468 }