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