File indexing completed on 2024-10-17 22:59:00
0001 #include "RecoTracker/MkFitCore/interface/cms_common_macros.h"
0002 #include "RecoTracker/MkFitCMS/interface/MkStdSeqs.h"
0003
0004 #include "RecoTracker/MkFitCore/interface/HitStructures.h"
0005 #include "RecoTracker/MkFitCore/interface/IterationConfig.h"
0006 #include "RecoTracker/MkFitCore/interface/MkJob.h"
0007 #include "RecoTracker/MkFitCore/interface/TrackStructures.h"
0008
0009 #include "RecoTracker/MkFitCore/interface/binnor.h"
0010
0011 namespace mkfit {
0012
0013 namespace StdSeq {
0014
0015
0016
0017
0018
0019 void loadDeads(EventOfHits &eoh, const std::vector<DeadVec> &deadvectors) {
0020 for (size_t il = 0; il < deadvectors.size(); il++) {
0021 eoh.suckInDeads(int(il), deadvectors[il]);
0022 }
0023 }
0024
0025
0026
0027
0028
0029 void cmssw_LoadHits_Begin(EventOfHits &eoh, const std::vector<const HitVec *> &orig_hitvectors) {
0030 eoh.reset();
0031 for (int i = 0; i < eoh.nLayers(); ++i) {
0032 auto &&l = eoh[i];
0033 l.beginRegistrationOfHits(*orig_hitvectors[l.is_pixel() ? 0 : 1]);
0034 }
0035 }
0036
0037
0038
0039
0040
0041
0042
0043
0044
0045 void cmssw_LoadHits_End(EventOfHits &eoh) {
0046 for (int i = 0; i < eoh.nLayers(); ++i) {
0047 auto &&l = eoh[i];
0048 l.endRegistrationOfHits(false);
0049 }
0050 }
0051
0052
0053
0054
0055
0056 void cmssw_Map_TrackHitIndices(const EventOfHits &eoh, TrackVec &seeds) {
0057 for (auto &&track : seeds) {
0058 for (int i = 0; i < track.nTotalHits(); ++i) {
0059 const int hitidx = track.getHitIdx(i);
0060 const int hitlyr = track.getHitLyr(i);
0061 if (hitidx >= 0) {
0062 const auto &loh = eoh[hitlyr];
0063 track.setHitIdx(i, loh.getHitIndexFromOriginal(hitidx));
0064 }
0065 }
0066 }
0067 }
0068
0069 void cmssw_ReMap_TrackHitIndices(const EventOfHits &eoh, TrackVec &out_tracks) {
0070 for (auto &&track : out_tracks) {
0071 for (int i = 0; i < track.nTotalHits(); ++i) {
0072 const int hitidx = track.getHitIdx(i);
0073 const int hitlyr = track.getHitLyr(i);
0074 if (hitidx >= 0) {
0075 const auto &loh = eoh[hitlyr];
0076 track.setHitIdx(i, loh.getOriginalHitIndex(hitidx));
0077 }
0078 }
0079 }
0080 }
0081
0082
0083
0084
0085 int clean_cms_seedtracks_iter(TrackVec &seeds, const IterationConfig &itrcfg, const BeamSpot &bspot) {
0086 using Algo = TrackBase::TrackAlgorithm;
0087
0088 const float etamax_brl = Config::c_etamax_brl;
0089 const float dpt_common = Config::c_dpt_common;
0090
0091 const float dzmax_bh = itrcfg.sc_dzmax_bh;
0092 const float drmax_bh = itrcfg.sc_drmax_bh;
0093 const float dzmax_eh = itrcfg.sc_dzmax_eh;
0094 const float drmax_eh = itrcfg.sc_drmax_eh;
0095 const float dzmax_bl = itrcfg.sc_dzmax_bl;
0096 const float drmax_bl = itrcfg.sc_drmax_bl;
0097 const float dzmax_el = itrcfg.sc_dzmax_el;
0098 const float drmax_el = itrcfg.sc_drmax_el;
0099
0100 const float ptmin_hpt = itrcfg.sc_ptthr_hpt;
0101
0102 const float dzmax2_inv_bh = 1.f / (dzmax_bh * dzmax_bh);
0103 const float drmax2_inv_bh = 1.f / (drmax_bh * drmax_bh);
0104 const float dzmax2_inv_eh = 1.f / (dzmax_eh * dzmax_eh);
0105 const float drmax2_inv_eh = 1.f / (drmax_eh * drmax_eh);
0106 const float dzmax2_inv_bl = 1.f / (dzmax_bl * dzmax_bl);
0107 const float drmax2_inv_bl = 1.f / (drmax_bl * drmax_bl);
0108 const float dzmax2_inv_el = 1.f / (dzmax_el * dzmax_el);
0109 const float drmax2_inv_el = 1.f / (drmax_el * drmax_el);
0110
0111
0112
0113
0114 const bool merge_hits = true;
0115 const float ptmax_merge_lowPtQuad = 0.2;
0116 const float etamin_merge_lowPtQuad = 1.5;
0117
0118 if (seeds.empty())
0119 return 0;
0120
0121 const int ns = seeds.size();
0122 #ifdef DEBUG
0123 std::cout << "before seed cleaning " << seeds.size() << std::endl;
0124 #endif
0125 TrackVec cleanSeedTracks;
0126 cleanSeedTracks.reserve(ns);
0127 std::vector<bool> writetrack(ns, true);
0128
0129 const float invR1GeV = 1.f / Config::track1GeVradius;
0130
0131 std::vector<int> nHits(ns);
0132 std::vector<int> charge(ns);
0133 std::vector<float> oldPhi(ns);
0134 std::vector<float> pos2(ns);
0135 std::vector<float> eta(ns);
0136 std::vector<float> ctheta(ns);
0137 std::vector<float> invptq(ns);
0138 std::vector<float> pt(ns);
0139 std::vector<float> x(ns);
0140 std::vector<float> y(ns);
0141 std::vector<float> z(ns);
0142 std::vector<float> d0(ns);
0143 int i1, i2;
0144
0145 axis_pow2_u1<float, unsigned short, 16, 8> ax_phi(-Const::PI, Const::PI);
0146 axis<float, unsigned short, 8, 8> ax_eta(-3.0, 3.0, 30u);
0147 binnor<unsigned int, decltype(ax_phi), decltype(ax_eta), 24, 8> phi_eta_binnor(ax_phi, ax_eta);
0148
0149 phi_eta_binnor.begin_registration(ns);
0150
0151 for (int ts = 0; ts < ns; ts++) {
0152 const Track &tk = seeds[ts];
0153 nHits[ts] = tk.nFoundHits();
0154 charge[ts] = tk.charge();
0155 oldPhi[ts] = tk.momPhi();
0156 pos2[ts] = std::pow(tk.x(), 2) + std::pow(tk.y(), 2);
0157 eta[ts] = tk.momEta();
0158 ctheta[ts] = 1.f / std::tan(tk.theta());
0159 invptq[ts] = tk.charge() * tk.invpT();
0160 pt[ts] = tk.pT();
0161 x[ts] = tk.x();
0162 y[ts] = tk.y();
0163 z[ts] = tk.z();
0164 d0[ts] = tk.d0BeamSpot(bspot.x, bspot.y);
0165
0166 phi_eta_binnor.register_entry_safe(oldPhi[ts], eta[ts]);
0167
0168 }
0169
0170 phi_eta_binnor.finalize_registration();
0171
0172 for (int sorted_ts = 0; sorted_ts < ns; sorted_ts++) {
0173 int ts = phi_eta_binnor.m_ranks[sorted_ts];
0174
0175 if (not writetrack[ts])
0176 continue;
0177
0178 const float oldPhi1 = oldPhi[ts];
0179 const float pos2_first = pos2[ts];
0180 const float eta1 = eta[ts];
0181 const float pt1 = pt[ts];
0182 const float invptq_first = invptq[ts];
0183
0184
0185 int n_ovlp_hits_added = 0;
0186
0187 auto phi_rng = ax_phi.from_R_rdr_to_N_bins(oldPhi[ts], 0.08f);
0188 auto eta_rng = ax_eta.from_R_rdr_to_N_bins(eta[ts], .1f);
0189
0190 for (auto i_phi = phi_rng.begin; i_phi != phi_rng.end; i_phi = ax_phi.next_N_bin(i_phi)) {
0191 for (auto i_eta = eta_rng.begin; i_eta != eta_rng.end; i_eta = ax_eta.next_N_bin(i_eta)) {
0192 const auto cbin = phi_eta_binnor.get_content(i_phi, i_eta);
0193 for (auto i = cbin.first; i < cbin.end(); ++i) {
0194 int tss = phi_eta_binnor.m_ranks[i];
0195
0196 if (not writetrack[ts])
0197 continue;
0198 if (not writetrack[tss])
0199 continue;
0200 if (tss == ts)
0201 continue;
0202
0203 const float pt2 = pt[tss];
0204
0205
0206 if (charge[tss] != charge[ts])
0207 continue;
0208
0209 const float thisDPt = std::abs(pt2 - pt1);
0210
0211 if (thisDPt > dpt_common * pt1)
0212 continue;
0213
0214 const float eta2 = eta[tss];
0215 const float deta2 = std::pow(eta1 - eta2, 2);
0216
0217 const float oldPhi2 = oldPhi[tss];
0218
0219 const float pos2_second = pos2[tss];
0220 const float thisDXYSign05 = pos2_second > pos2_first ? -0.5f : 0.5f;
0221
0222 const float thisDXY = thisDXYSign05 * sqrt(std::pow(x[ts] - x[tss], 2) + std::pow(y[ts] - y[tss], 2));
0223
0224 const float invptq_second = invptq[tss];
0225
0226 const float newPhi1 = oldPhi1 - thisDXY * invR1GeV * invptq_first;
0227 const float newPhi2 = oldPhi2 + thisDXY * invR1GeV * invptq_second;
0228
0229 const float dphi = cdist(std::abs(newPhi1 - newPhi2));
0230
0231 const float dr2 = deta2 + dphi * dphi;
0232
0233 const float thisDZ = z[ts] - z[tss] - thisDXY * (ctheta[ts] + ctheta[tss]);
0234 const float dz2 = thisDZ * thisDZ;
0235
0236
0237
0238 bool overlapping = false;
0239 if (std::abs(eta1) < etamax_brl) {
0240 if (pt1 > ptmin_hpt) {
0241 if (dz2 * dzmax2_inv_bh + dr2 * drmax2_inv_bh < 1.0f)
0242 overlapping = true;
0243 } else {
0244 if (dz2 * dzmax2_inv_bl + dr2 * drmax2_inv_bl < 1.0f)
0245 overlapping = true;
0246 }
0247 } else {
0248 if (pt1 > ptmin_hpt) {
0249 if (dz2 * dzmax2_inv_eh + dr2 * drmax2_inv_eh < 1.0f)
0250 overlapping = true;
0251 } else {
0252 if (dz2 * dzmax2_inv_el + dr2 * drmax2_inv_el < 1.0f)
0253 overlapping = true;
0254 }
0255 }
0256
0257 if (overlapping) {
0258
0259 i1 = ts;
0260 i2 = tss;
0261 if (d0[tss] > d0[ts])
0262 writetrack[tss] = false;
0263 else {
0264 writetrack[ts] = false;
0265 i2 = ts;
0266 i1 = tss;
0267 }
0268
0269
0270
0271 Track &tk = seeds[i1];
0272 if (merge_hits && tk.nTotalHits() < Track::Status::kMaxSeedHits &&
0273 (Algo(itrcfg.m_track_algorithm) != Algo::lowPtQuadStep ||
0274 (pt1 < ptmax_merge_lowPtQuad && std::abs(eta1) > etamin_merge_lowPtQuad))) {
0275 const Track &tk2 = seeds[i2];
0276
0277 float fakeChi2 = 0.0;
0278
0279 for (int j = 0; j < tk2.nTotalHits(); ++j) {
0280 int hitidx = tk2.getHitIdx(j);
0281 int hitlyr = tk2.getHitLyr(j);
0282 if (hitidx >= 0) {
0283 bool unique = true;
0284 for (int i = 0; i < tk.nTotalHits(); ++i) {
0285 if ((hitidx == tk.getHitIdx(i)) && (hitlyr == tk.getHitLyr(i))) {
0286 unique = false;
0287 break;
0288 }
0289 }
0290 if (unique) {
0291 tk.addHitIdx(tk2.getHitIdx(j), tk2.getHitLyr(j), fakeChi2);
0292 ++n_ovlp_hits_added;
0293 if (tk.nTotalHits() >= Track::Status::kMaxSeedHits)
0294 break;
0295 }
0296 }
0297 }
0298 }
0299 if (n_ovlp_hits_added > 0)
0300 tk.sortHitsByLayer();
0301 }
0302 }
0303 }
0304 }
0305
0306 if (writetrack[ts]) {
0307 cleanSeedTracks.emplace_back(seeds[ts]);
0308 }
0309 }
0310
0311 seeds.swap(cleanSeedTracks);
0312
0313 #ifdef DEBUG
0314 {
0315 const int ns2 = seeds.size();
0316 printf("Number of CMS seeds before %d --> after %d cleaning\n", ns, ns2);
0317
0318 for (int it = 0; it < ns2; it++) {
0319 const Track &ss = seeds[it];
0320 printf(" %3i q=%+i pT=%7.3f eta=% 7.3f nHits=%i label=% i\n",
0321 it,
0322 ss.charge(),
0323 ss.pT(),
0324 ss.momEta(),
0325 ss.nFoundHits(),
0326 ss.label());
0327 }
0328 }
0329 #endif
0330
0331 #ifdef DEBUG
0332 std::cout << "AFTER seed cleaning " << seeds.size() << std::endl;
0333 #endif
0334
0335 return seeds.size();
0336 }
0337
0338 namespace {
0339 CMS_SA_ALLOW struct register_seed_cleaners {
0340 register_seed_cleaners() {
0341 IterationConfig::register_seed_cleaner("phase1:default", clean_cms_seedtracks_iter);
0342 }
0343 } rsc_instance;
0344 }
0345
0346
0347
0348
0349
0350 void remove_duplicates(TrackVec &tracks) {
0351 tracks.erase(std::remove_if(tracks.begin(), tracks.end(), [](auto track) { return track.getDuplicateValue(); }),
0352 tracks.end());
0353 }
0354
0355 void clean_duplicates(TrackVec &tracks, const IterationConfig &) {
0356 const auto ntracks = tracks.size();
0357 float eta1, phi1, pt1, deta, dphi, dr2;
0358
0359 if (ntracks == 0) {
0360 return;
0361 }
0362 for (auto itrack = 0U; itrack < ntracks - 1; itrack++) {
0363 auto &track = tracks[itrack];
0364 eta1 = track.momEta();
0365 phi1 = track.momPhi();
0366 pt1 = track.pT();
0367 for (auto jtrack = itrack + 1; jtrack < ntracks; jtrack++) {
0368 auto &track2 = tracks[jtrack];
0369 if (track.label() == track2.label())
0370 continue;
0371 if (track.algoint() != track2.algoint())
0372 continue;
0373
0374 deta = std::abs(track2.momEta() - eta1);
0375 if (deta > Config::maxdEta)
0376 continue;
0377
0378 dphi = std::abs(squashPhiMinimal(phi1 - track2.momPhi()));
0379 if (dphi > Config::maxdPhi)
0380 continue;
0381
0382 float maxdR = Config::maxdR;
0383 float maxdRSquared = maxdR * maxdR;
0384 if (std::abs(eta1) > 2.5f)
0385 maxdRSquared *= 16.0f;
0386 else if (std::abs(eta1) > 1.44f)
0387 maxdRSquared *= 9.0f;
0388 dr2 = dphi * dphi + deta * deta;
0389 if (dr2 < maxdRSquared) {
0390
0391 if (track.score() > track2.score())
0392 track2.setDuplicateValue(true);
0393 else
0394 track.setDuplicateValue(true);
0395 continue;
0396 } else {
0397 if (pt1 == 0)
0398 continue;
0399 if (track2.pT() == 0)
0400 continue;
0401
0402 if (std::abs((1 / track2.pT()) - (1 / pt1)) < Config::maxdPt) {
0403 if (Config::useHitsForDuplicates) {
0404 float numHitsShared = 0;
0405 for (int ihit2 = 0; ihit2 < track2.nTotalHits(); ihit2++) {
0406 const int hitidx2 = track2.getHitIdx(ihit2);
0407 const int hitlyr2 = track2.getHitLyr(ihit2);
0408 if (hitidx2 >= 0) {
0409 auto const it = std::find_if(track.beginHitsOnTrack(),
0410 track.endHitsOnTrack(),
0411 [&hitidx2, &hitlyr2](const HitOnTrack &element) {
0412 return (element.index == hitidx2 && element.layer == hitlyr2);
0413 });
0414 if (it != track.endHitsOnTrack())
0415 numHitsShared++;
0416 }
0417 }
0418
0419 float fracHitsShared = numHitsShared / std::min(track.nFoundHits(), track2.nFoundHits());
0420
0421 if (fracHitsShared < Config::minFracHitsShared)
0422 continue;
0423 }
0424
0425 if (track.score() > track2.score())
0426 track2.setDuplicateValue(true);
0427 else
0428 track.setDuplicateValue(true);
0429 }
0430 }
0431 }
0432 }
0433
0434 remove_duplicates(tracks);
0435 }
0436
0437
0438
0439
0440
0441 void clean_duplicates_sharedhits(TrackVec &tracks, const IterationConfig &itconf) {
0442 const float fraction = itconf.dc_fracSharedHits;
0443 const auto ntracks = tracks.size();
0444
0445 std::vector<float> ctheta(ntracks);
0446 std::multimap<int, int> hitMap;
0447
0448 for (auto itrack = 0U; itrack < ntracks; itrack++) {
0449 auto &trk = tracks[itrack];
0450 ctheta[itrack] = 1.f / std::tan(trk.theta());
0451 for (int i = 0; i < trk.nTotalHits(); ++i) {
0452 if (trk.getHitIdx(i) < 0)
0453 continue;
0454 int a = trk.getHitLyr(i);
0455 int b = trk.getHitIdx(i);
0456 hitMap.insert(std::make_pair(b * 1000 + a, i > 0 ? itrack : -itrack));
0457 }
0458 }
0459
0460 for (auto itrack = 0U; itrack < ntracks; itrack++) {
0461 auto &trk = tracks[itrack];
0462 auto phi1 = trk.momPhi();
0463 auto ctheta1 = ctheta[itrack];
0464
0465 std::map<int, int> sharingMap;
0466 for (int i = 0; i < trk.nTotalHits(); ++i) {
0467 if (trk.getHitIdx(i) < 0)
0468 continue;
0469 int a = trk.getHitLyr(i);
0470 int b = trk.getHitIdx(i);
0471 auto range = hitMap.equal_range(b * 1000 + a);
0472 for (auto it = range.first; it != range.second; ++it) {
0473 if (std::abs(it->second) >= (int)itrack)
0474 continue;
0475 if (i == 0 && it->second < 0)
0476 continue;
0477 sharingMap[std::abs(it->second)]++;
0478 }
0479 }
0480
0481 for (const auto &elem : sharingMap) {
0482 auto &track2 = tracks[elem.first];
0483
0484
0485 auto dctheta = std::abs(ctheta[elem.first] - ctheta1);
0486 if (dctheta > 1.)
0487 continue;
0488
0489 auto dphi = std::abs(squashPhiMinimal(phi1 - track2.momPhi()));
0490 if (dphi > 1.)
0491 continue;
0492
0493 if (elem.second >= std::min(trk.nFoundHits(), track2.nFoundHits()) * fraction) {
0494 if (trk.score() > track2.score())
0495 track2.setDuplicateValue(true);
0496 else
0497 trk.setDuplicateValue(true);
0498 }
0499 }
0500 }
0501
0502 remove_duplicates(tracks);
0503 }
0504
0505 void clean_duplicates_sharedhits_pixelseed(TrackVec &tracks, const IterationConfig &itconf) {
0506 const float fraction = itconf.dc_fracSharedHits;
0507 const float drth_central = itconf.dc_drth_central;
0508 const float drth_obarrel = itconf.dc_drth_obarrel;
0509 const float drth_forward = itconf.dc_drth_forward;
0510 const auto ntracks = tracks.size();
0511
0512 std::vector<float> ctheta(ntracks);
0513 for (auto itrack = 0U; itrack < ntracks; itrack++) {
0514 auto &trk = tracks[itrack];
0515 ctheta[itrack] = 1.f / std::tan(trk.theta());
0516 }
0517
0518 float phi1, invpt1, dctheta, ctheta1, dphi, dr2;
0519 for (auto itrack = 0U; itrack < ntracks; itrack++) {
0520 auto &trk = tracks[itrack];
0521 phi1 = trk.momPhi();
0522 invpt1 = trk.invpT();
0523 ctheta1 = ctheta[itrack];
0524 for (auto jtrack = itrack + 1; jtrack < ntracks; jtrack++) {
0525 auto &track2 = tracks[jtrack];
0526 if (trk.label() == track2.label())
0527 continue;
0528
0529 dctheta = std::abs(ctheta[jtrack] - ctheta1);
0530
0531 if (dctheta > Config::maxdcth)
0532 continue;
0533
0534 dphi = std::abs(squashPhiMinimal(phi1 - track2.momPhi()));
0535
0536 if (dphi > Config::maxdphi)
0537 continue;
0538
0539 float maxdRSquared = drth_central * drth_central;
0540 if (std::abs(ctheta1) > Config::maxcth_fw)
0541 maxdRSquared = drth_forward * drth_forward;
0542 else if (std::abs(ctheta1) > Config::maxcth_ob)
0543 maxdRSquared = drth_obarrel * drth_obarrel;
0544 dr2 = dphi * dphi + dctheta * dctheta;
0545 if (dr2 < maxdRSquared) {
0546
0547 if (trk.score() > track2.score())
0548 track2.setDuplicateValue(true);
0549 else
0550 trk.setDuplicateValue(true);
0551 continue;
0552 }
0553
0554 if (std::abs(track2.invpT() - invpt1) > Config::maxd1pt)
0555 continue;
0556
0557 auto sharedCount = 0;
0558 auto sharedFirst = 0;
0559 const auto minFoundHits = std::min(trk.nFoundHits(), track2.nFoundHits());
0560
0561 for (int i = 0; i < trk.nTotalHits(); ++i) {
0562 if (trk.getHitIdx(i) < 0)
0563 continue;
0564 const int a = trk.getHitLyr(i);
0565 const int b = trk.getHitIdx(i);
0566 for (int j = 0; j < track2.nTotalHits(); ++j) {
0567 if (track2.getHitIdx(j) < 0)
0568 continue;
0569 const int c = track2.getHitLyr(j);
0570 const int d = track2.getHitIdx(j);
0571
0572
0573 if (a == c && b == d)
0574 sharedCount += 1;
0575 if (j == 0 && i == 0 && a == c && b == d)
0576 sharedFirst += 1;
0577
0578 if ((sharedCount - sharedFirst) >= ((minFoundHits - sharedFirst) * fraction))
0579 continue;
0580 }
0581 if ((sharedCount - sharedFirst) >= ((minFoundHits - sharedFirst) * fraction))
0582 continue;
0583 }
0584
0585
0586 if ((sharedCount - sharedFirst) >= ((minFoundHits - sharedFirst) * fraction)) {
0587 if (trk.score() > track2.score())
0588 track2.setDuplicateValue(true);
0589 else
0590 trk.setDuplicateValue(true);
0591 }
0592 }
0593 }
0594
0595 remove_duplicates(tracks);
0596 }
0597
0598 namespace {
0599 CMS_SA_ALLOW struct register_duplicate_cleaners {
0600 register_duplicate_cleaners() {
0601 IterationConfig::register_duplicate_cleaner("phase1:clean_duplicates", clean_duplicates);
0602 IterationConfig::register_duplicate_cleaner("phase1:clean_duplicates_sharedhits",
0603 clean_duplicates_sharedhits);
0604 IterationConfig::register_duplicate_cleaner("phase1:clean_duplicates_sharedhits_pixelseed",
0605 clean_duplicates_sharedhits_pixelseed);
0606 }
0607 } rdc_instance;
0608 }
0609
0610
0611
0612
0613
0614
0615
0616 template <class TRACK>
0617 bool qfilter_n_hits(const TRACK &t, const MkJob &j) {
0618 int seedHits = t.getNSeedHits();
0619 int seedReduction = (seedHits <= 5) ? 2 : 3;
0620 return t.nFoundHits() - seedReduction >= j.params_cur().minHitsQF;
0621 }
0622
0623
0624 template <class TRACK>
0625 bool qfilter_n_hits_pixseed(const TRACK &t, const MkJob &j) {
0626 return t.nFoundHits() >= j.params_cur().minHitsQF;
0627 }
0628
0629
0630
0631 template <class TRACK>
0632 bool qfilter_n_layers(const TRACK &t, const MkJob &j) {
0633 const BeamSpot &bspot = j.m_beam_spot;
0634 const TrackerInfo &trk_inf = j.m_trk_info;
0635 int enhits = t.nHitsByTypeEncoded(trk_inf);
0636 int npixhits = t.nPixelDecoded(enhits);
0637 int enlyrs = t.nLayersByTypeEncoded(trk_inf);
0638 int npixlyrs = t.nPixelDecoded(enlyrs);
0639 int nmatlyrs = t.nTotMatchDecoded(enlyrs);
0640 int llyr = t.getLastFoundHitLyr();
0641 int lplyr = t.getLastFoundPixelHitLyr();
0642 float invpt = t.invpT();
0643
0644
0645 float invptmin = 1.43;
0646 float d0BS = t.d0BeamSpot(bspot.x, bspot.y);
0647 float d0_max = 0.1;
0648
0649
0650 bool endsInsidePix = (llyr == 2 || llyr == 18 || llyr == 45);
0651
0652 bool lastInsidePix = ((0 <= lplyr && lplyr < 3) || (18 <= lplyr && lplyr < 20) || (45 <= lplyr && lplyr < 47));
0653
0654 return !(((npixhits <= 3 || npixlyrs <= 3) && endsInsidePix &&
0655 (invpt < invptmin || (invpt >= invptmin && std::abs(d0BS) > d0_max))) ||
0656 ((npixlyrs <= 3 && nmatlyrs <= 6) && lastInsidePix && llyr != lplyr && std::abs(d0BS) > d0_max));
0657 }
0658
0659
0660
0661 template <class TRACK>
0662 bool qfilter_pixelLessFwd(const TRACK &t, const MkJob &j) {
0663 const BeamSpot &bspot = j.m_beam_spot;
0664 const TrackerInfo &tk_info = j.m_trk_info;
0665 float d0BS = t.d0BeamSpot(bspot.x, bspot.y);
0666 float d0_max = 0.05;
0667
0668 int encoded;
0669 encoded = t.nLayersByTypeEncoded(tk_info);
0670 int nLyrs = t.nTotMatchDecoded(encoded);
0671 encoded = t.nHitsByTypeEncoded(tk_info);
0672 int nHits = t.nTotMatchDecoded(encoded);
0673
0674
0675 int seedReduction = (t.getNSeedHits() <= 5) ? 2 : 3;
0676
0677
0678 float invpt = t.invpT();
0679 float invptmin = 1.11;
0680
0681 float thetasym = std::abs(t.theta() - Const::PIOver2);
0682 float thetasymmin = 1.11;
0683
0684
0685 return (((t.nFoundHits() - seedReduction >= 4 && invpt < invptmin) ||
0686 (t.nFoundHits() - seedReduction >= 3 && invpt > invptmin && thetasym <= thetasymmin) ||
0687 (t.nFoundHits() - seedReduction >= 4 && invpt > invptmin && thetasym > thetasymmin)) &&
0688 !((nLyrs <= 4 || nHits <= 4) && std::abs(d0BS) > d0_max && invpt < invptmin));
0689 }
0690
0691
0692
0693 template <class TRACK>
0694 bool qfilter_pixelLessBkwd(const TRACK &t, const MkJob &j) {
0695 const BeamSpot &bspot = j.m_beam_spot;
0696 const TrackerInfo &tk_info = j.m_trk_info;
0697 float d0BS = t.d0BeamSpot(bspot.x, bspot.y);
0698 float d0_max = 0.1;
0699
0700 int encoded;
0701 encoded = t.nLayersByTypeEncoded(tk_info);
0702 int nLyrs = t.nTotMatchDecoded(encoded);
0703 encoded = t.nHitsByTypeEncoded(tk_info);
0704 int nHits = t.nTotMatchDecoded(encoded);
0705
0706
0707 float invpt = t.invpT();
0708 float invptmin = 1.11;
0709
0710 float thetasym = std::abs(t.theta() - Const::PIOver2);
0711 float thetasymmin_l = 0.80;
0712 float thetasymmin_h = 1.11;
0713
0714
0715 return !(
0716 ((nLyrs <= 3 || nHits <= 3)) ||
0717 ((nLyrs <= 4 || nHits <= 4) && (invpt < invptmin || (thetasym > thetasymmin_l && std::abs(d0BS) > d0_max))) ||
0718 ((nLyrs <= 5 || nHits <= 5) && (invpt > invptmin && thetasym > thetasymmin_h && std::abs(d0BS) > d0_max)));
0719 }
0720
0721 namespace {
0722 CMS_SA_ALLOW struct register_quality_filters {
0723 register_quality_filters() {
0724 IterationConfig::register_candidate_filter("phase1:qfilter_n_hits", qfilter_n_hits<TrackCand>);
0725 IterationConfig::register_candidate_filter("phase1:qfilter_n_hits_pixseed",
0726 qfilter_n_hits_pixseed<TrackCand>);
0727 IterationConfig::register_candidate_filter("phase1:qfilter_n_layers", qfilter_n_layers<TrackCand>);
0728 IterationConfig::register_candidate_filter("phase1:qfilter_pixelLessFwd", qfilter_pixelLessFwd<TrackCand>);
0729 IterationConfig::register_candidate_filter("phase1:qfilter_pixelLessBkwd", qfilter_pixelLessBkwd<TrackCand>);
0730 }
0731 } rqf_instance;
0732 }
0733
0734
0735
0736
0737
0738 float trackScoreDefault(const int nfoundhits,
0739 const int ntailholes,
0740 const int noverlaphits,
0741 const int nmisshits,
0742 const float chi2,
0743 const float pt,
0744 const bool inFindCandidates) {
0745 float maxBonus = 8.0;
0746 float bonus = Config::validHitSlope_ * nfoundhits + Config::validHitBonus_;
0747 float penalty = Config::missingHitPenalty_;
0748 float tailPenalty = Config::tailMissingHitPenalty_;
0749 float overlapBonus = Config::overlapHitBonus_;
0750 if (pt < 0.9) {
0751 penalty *= inFindCandidates ? 1.7f : 1.5f;
0752 bonus = std::min(bonus * (inFindCandidates ? 0.9f : 1.0f), maxBonus);
0753 }
0754 float score =
0755 bonus * nfoundhits + overlapBonus * noverlaphits - penalty * nmisshits - tailPenalty * ntailholes - chi2;
0756 return score;
0757 }
0758
0759 namespace {
0760 CMS_SA_ALLOW struct register_track_scorers {
0761 register_track_scorers() {
0762 IterationConfig::register_track_scorer("default", trackScoreDefault);
0763 IterationConfig::register_track_scorer("phase1:default", trackScoreDefault);
0764 }
0765 } rts_instance;
0766 }
0767
0768 }
0769 }