File indexing completed on 2023-11-04 00:02:30
0001 #include "RecoTracker/MkFitCore/interface/HitStructures.h"
0002
0003 #include "RecoTracker/MkFitCore/interface/IterationConfig.h"
0004 #include "Matriplex/Memory.h"
0005
0006 #include "Debug.h"
0007
0008 namespace mkfit {
0009
0010 void LayerOfHits::Initializator::setup(float qmin, float qmax, float dq) {
0011 assert(qmax > qmin);
0012 float extent = qmax - qmin;
0013 m_nq = std::ceil(extent / dq);
0014 float extra = 0.5f * (m_nq * dq - extent);
0015 m_qmin = qmin - extra;
0016 m_qmax = qmax + extra;
0017 }
0018
0019 LayerOfHits::Initializator::Initializator(const LayerInfo &li, float qmin, float qmax, unsigned int nq)
0020 : m_linfo(li), m_qmin(qmin), m_qmax(qmax), m_nq(nq) {}
0021
0022 LayerOfHits::Initializator::Initializator(const LayerInfo &li, float qmin, float qmax, float dq) : m_linfo(li) {
0023 setup(qmin, qmax, dq);
0024 }
0025
0026 LayerOfHits::Initializator::Initializator(const LayerInfo &li) : m_linfo(li) {
0027 if (li.is_barrel())
0028 setup(li.zmin(), li.zmax(), li.q_bin());
0029 else
0030 setup(li.rin(), li.rout(), li.q_bin());
0031 }
0032
0033 LayerOfHits::LayerOfHits(const LayerOfHits::Initializator &i)
0034 : m_ax_phi(-Const::PI, Const::PI),
0035 m_ax_eta(i.m_qmin, i.m_qmax, i.m_nq),
0036 m_binnor(m_ax_phi, m_ax_eta, true, false)
0037 {
0038 m_layer_info = &i.m_linfo;
0039 m_is_barrel = m_layer_info->is_barrel();
0040
0041 m_dead_bins.resize(m_ax_eta.size_of_N() * m_ax_phi.size_of_N());
0042 }
0043
0044 LayerOfHits::~LayerOfHits() {
0045 #ifdef COPY_SORTED_HITS
0046 free_hits();
0047 #endif
0048 }
0049
0050 #ifdef COPY_SORTED_HITS
0051 void LayerOfHits::alloc_hits(int size) {
0052 m_hits = (Hit *)Matriplex::aligned_alloc64(sizeof(Hit) * size);
0053 m_capacity = size;
0054 for (int ihit = 0; ihit < m_capacity; ihit++) {
0055 m_hits[ihit] = Hit();
0056 }
0057 }
0058
0059 void LayerOfHits::free_hits() { std::free(m_hits); }
0060 #endif
0061
0062
0063
0064 void LayerOfHits::reset() {
0065 m_hit_infos.clear();
0066 m_ext_idcs.clear();
0067 m_min_ext_idx = std::numeric_limits<unsigned int>::max();
0068 m_max_ext_idx = std::numeric_limits<unsigned int>::min();
0069 m_n_hits = 0;
0070 m_binnor.reset_contents();
0071 }
0072
0073
0074
0075 void LayerOfHits::suckInHits(const HitVec &hitv) {
0076 m_ext_hits = &hitv;
0077 m_n_hits = hitv.size();
0078
0079 m_binnor.begin_registration(m_n_hits);
0080
0081 #ifdef COPY_SORTED_HITS
0082 if (m_capacity < m_n_hits) {
0083 free_hits();
0084 alloc_hits(m_n_hits);
0085 }
0086 #endif
0087
0088 std::vector<HitInfo> hinfos;
0089 if (Config::usePhiQArrays) {
0090 hinfos.reserve(m_n_hits);
0091 m_hit_infos.reserve(m_n_hits);
0092 }
0093
0094 for (unsigned int i = 0; i < m_n_hits; ++i) {
0095 const Hit &h = hitv[i];
0096
0097 float phi = h.phi();
0098 float q = m_is_barrel ? h.z() : h.r();
0099
0100 m_binnor.register_entry_safe(phi, q);
0101
0102 if (Config::usePhiQArrays) {
0103 const float sqrt3 = std::sqrt(3);
0104 float half_length, qbar;
0105 if (m_is_barrel) {
0106 half_length = sqrt3 * std::sqrt(h.ezz());
0107 qbar = h.r();
0108 } else {
0109 half_length = sqrt3 * std::sqrt(h.exx() + h.eyy());
0110 qbar = h.z();
0111 }
0112 hinfos.emplace_back(HitInfo({phi, q, half_length, qbar}));
0113 }
0114 }
0115
0116 m_binnor.finalize_registration();
0117
0118 for (unsigned int i = 0; i < m_n_hits; ++i) {
0119 unsigned int j = m_binnor.m_ranks[i];
0120 #ifdef COPY_SORTED_HITS
0121 memcpy(&m_hits[i], &hitv[j], sizeof(Hit));
0122 #endif
0123 if (Config::usePhiQArrays) {
0124 m_hit_infos.emplace_back(hinfos[j]);
0125 }
0126 }
0127 }
0128
0129
0130
0131 void LayerOfHits::suckInDeads(const DeadVec &deadv) {
0132 m_dead_bins.assign(m_dead_bins.size(), false);
0133
0134 for (const auto &d : deadv) {
0135 bin_index_t q_bin_1 = qBinChecked(d.q1);
0136 bin_index_t q_bin_2 = qBinChecked(d.q2) + 1;
0137 bin_index_t phi_bin_1 = phiBin(d.phi1);
0138 bin_index_t phi_bin_2 = phiMaskApply(phiBin(d.phi2) + 1);
0139
0140 for (bin_index_t q_bin = q_bin_1; q_bin != q_bin_2; q_bin++) {
0141 const unsigned int qoff = q_bin * m_ax_phi.size_of_N();
0142 for (bin_index_t pb = phi_bin_1; pb != phi_bin_2; pb = phiMaskApply(pb + 1)) {
0143 m_dead_bins[qoff + pb] = true;
0144 }
0145 }
0146 }
0147 }
0148
0149
0150
0151 void LayerOfHits::beginRegistrationOfHits(const HitVec &hitv) {
0152 m_ext_hits = &hitv;
0153 m_n_hits = 0;
0154
0155 m_binnor.begin_registration(128);
0156 }
0157
0158 void LayerOfHits::registerHit(unsigned int idx) {
0159 const Hit &h = (*m_ext_hits)[idx];
0160
0161 m_ext_idcs.push_back(idx);
0162 m_min_ext_idx = std::min(m_min_ext_idx, idx);
0163 m_max_ext_idx = std::max(m_max_ext_idx, idx);
0164
0165 float phi = h.phi();
0166 float q = m_is_barrel ? h.z() : h.r();
0167
0168 m_binnor.register_entry_safe(phi, q);
0169
0170 if (Config::usePhiQArrays) {
0171 const float sqrt3 = std::sqrt(3);
0172 float half_length, qbar;
0173 if (m_is_barrel) {
0174 half_length = sqrt3 * std::sqrt(h.ezz());
0175 qbar = h.r();
0176 } else {
0177 half_length = sqrt3 * std::sqrt(h.exx() + h.eyy());
0178 qbar = h.z();
0179 }
0180 m_hit_infos.emplace_back(HitInfo({phi, q, half_length, qbar}));
0181 }
0182 }
0183
0184 void LayerOfHits::endRegistrationOfHits(bool build_original_to_internal_map) {
0185 m_n_hits = m_ext_idcs.size();
0186 if (m_n_hits == 0)
0187 return;
0188
0189 m_binnor.finalize_registration();
0190
0191
0192
0193 #ifdef COPY_SORTED_HITS
0194 if (m_capacity < m_n_hits) {
0195 free_hits();
0196 alloc_hits(m_n_hits);
0197 }
0198 #endif
0199
0200 std::vector<HitInfo> hinfos;
0201 if (Config::usePhiQArrays) {
0202 hinfos.swap(m_hit_infos);
0203 m_hit_infos.reserve(m_n_hits);
0204 }
0205
0206 for (unsigned int i = 0; i < m_n_hits; ++i) {
0207 unsigned int j = m_binnor.m_ranks[i];
0208 unsigned int k = m_ext_idcs[j];
0209
0210 #ifdef COPY_SORTED_HITS
0211 memcpy(&m_hits[i], &hitv[k], sizeof(Hit));
0212 #endif
0213
0214 if (Config::usePhiQArrays) {
0215 m_hit_infos.emplace_back(hinfos[j]);
0216 }
0217
0218
0219 m_binnor.m_ranks[i] = k;
0220 }
0221
0222 if (build_original_to_internal_map) {
0223 if (m_max_ext_idx - m_min_ext_idx + 1 > 8 * m_n_hits) {
0224
0225
0226
0227
0228
0229
0230
0231 printf(
0232 "LayerOfHits::endRegistrationOfHits() original_to_internal index map vector is largish: m_n_hits=%d, "
0233 "map_vector_size=%d\n",
0234 m_n_hits,
0235 m_max_ext_idx - m_min_ext_idx + 1);
0236 }
0237
0238 m_ext_idcs.resize(m_max_ext_idx - m_min_ext_idx + 1);
0239 for (unsigned int i = 0; i < m_n_hits; ++i) {
0240 m_ext_idcs[m_hit_ranks[i] - m_min_ext_idx] = i;
0241 }
0242 }
0243
0244
0245
0246
0247
0248 }
0249
0250 void LayerOfHits::printBins() {
0251 for (bin_index_t qb = 0; qb <= m_ax_eta.m_last_N_bin; ++qb) {
0252 printf("%c bin %d\n", is_barrel() ? 'Z' : 'R', qb);
0253 for (bin_index_t pb = 0; pb <= m_ax_phi.m_last_N_bin; ++pb) {
0254 if (pb % 8 == 0)
0255 printf(" Phi %4d: ", pb);
0256 auto content = m_binnor.get_content(pb, qb);
0257 printf("%5d,%4d %s", content.first, content.count, ((pb + 1) % 8 == 0) ? "\n" : "");
0258 }
0259 }
0260 }
0261
0262
0263
0264
0265
0266 EventOfHits::EventOfHits(const TrackerInfo &trk_inf) : m_n_layers(trk_inf.n_layers()) {
0267 m_layers_of_hits.reserve(trk_inf.n_layers());
0268 for (int ii = 0; ii < trk_inf.n_layers(); ++ii) {
0269 m_layers_of_hits.emplace_back(LayerOfHits::Initializator(trk_inf.layer(ii)));
0270 }
0271 }
0272
0273 }