Acts/Seeding/GbtsDataStorage.hpp

0001 // This file is part of the Acts project.
0002 //
0003 // Copyright (C) 2021 CERN for the benefit of the Acts project
0004 //
0005 // This Source Code Form is subject to the terms of the Mozilla Public
0006 // License, v. 2.0. If a copy of the MPL was not distributed with this
0007 // file, You can obtain one at http://mozilla.org/MPL/2.0/.
0008
0009 #pragma once
0010
0011 // TODO: update to C++17 style
0012 #include "Acts/Definitions/Algebra.hpp"
0013 #include "Acts/Seeding/GbtsGeometry.hpp"
0014
0015 #include <algorithm>
0016 #include <map>
0017 #include <vector>
0018
0019 namespace Acts {
0020
0021 constexpr std::size_t MAX_SEG_PER_NODE = 1000;  // was 30
0022 constexpr std::size_t N_SEG_CONNS = 6;          // was 6
0023
0024 // new sp struct
0025 template <typename space_point_t>
0026 struct GbtsSP {
0027   const space_point_t *SP;  // want inside to have pointer
0028   int gbtsID;
0029   int combined_ID;
0030   GbtsSP(const space_point_t *sp, int id, int combined_id)
0031       : SP(sp), gbtsID(id), combined_ID{combined_id} {
0032     if (SP->sourceLinks().size() == 1) {  // pixels have 1 SL
0033       m_isPixel = true;
0034     } else {
0035       m_isPixel = false;
0036     }
0037     m_phi = std::atan(SP->x() / SP->y());
0038   };
0039   bool isPixel() const { return m_isPixel; }
0040   bool isSCT() const { return !m_isPixel; }
0041   float phi() const { return m_phi; }
0042   bool m_isPixel;
0043   float m_phi;
0044 };
0045
0046 template <typename space_point_t>
0047 class GbtsNode {
0048  public:
0049   struct CompareByPhi {
0050     bool operator()(const GbtsNode<space_point_t> *n1,
0051                     const GbtsNode<space_point_t> *n2) {
0052       return (n1->m_spGbts.phi() < n2->m_spGbts.phi());
0053     }
0054   };
0055
0056   GbtsNode(const GbtsSP<space_point_t> &spGbts, float minT = -100.0,
0057            float maxT = 100.0)
0058       : m_spGbts(spGbts), m_minCutOnTau(minT), m_maxCutOnTau(maxT) {}
0059
0060   inline void addIn(int i) {
0061     if (m_in.size() < MAX_SEG_PER_NODE) {
0062       m_in.push_back(i);
0063     }
0064   }
0065
0066   inline void addOut(int i) {
0067     if (m_out.size() < MAX_SEG_PER_NODE) {
0068       m_out.push_back(i);
0069     }
0070   }
0071
0072   inline bool isConnector() const {
0073     if (m_in.empty() || m_out.empty()) {
0074       return false;
0075     }
0076     return true;
0077   }
0078
0079   inline bool isFull() const {
0080     if (m_in.size() == MAX_SEG_PER_NODE && m_out.size() == MAX_SEG_PER_NODE) {
0081       return true;
0082     } else {
0083       return false;
0084     }
0085   }
0086
0087   const GbtsSP<space_point_t> &m_spGbts;
0088
0089   std::vector<unsigned int> m_in;  // indices of the edges in the edge storage
0090   std::vector<unsigned int> m_out;
0091   float m_minCutOnTau, m_maxCutOnTau;
0092 };
0093
0094 template <typename space_point_t>
0095 class GbtsEtaBin {
0096  public:
0097   GbtsEtaBin() { m_vn.clear(); }
0098
0099   ~GbtsEtaBin() {
0100     for (typename std::vector<GbtsNode<space_point_t> *>::iterator it =
0101              m_vn.begin();
0102          it != m_vn.end(); ++it) {
0103       delete (*it);
0104     }
0105   }
0106
0107   void sortByPhi() {
0108     std::sort(m_vn.begin(), m_vn.end(),
0109               typename Acts::GbtsNode<space_point_t>::CompareByPhi());
0110   }
0111
0112   bool empty() const { return m_vn.empty(); }
0113
0114   void generatePhiIndexing(float dphi) {
0115     for (unsigned int nIdx = 0; nIdx < m_vn.size(); nIdx++) {
0116       GbtsNode<space_point_t> *pN = m_vn.at(nIdx);
0117       // float phi = pN->m_sp.phi();
0118       // float phi = (std::atan(pN->m_sp.x() / pN->m_sp.y()));
0119       float phi = pN->m_spGbts.phi();
0120       if (phi <= M_PI - dphi) {
0121         continue;
0122       }
0123
0124       m_vPhiNodes.push_back(
0125           std::pair<float, unsigned int>(phi - 2 * M_PI, nIdx));
0126     }
0127
0128     for (unsigned int nIdx = 0; nIdx < m_vn.size(); nIdx++) {
0129       GbtsNode<space_point_t> *pN = m_vn.at(nIdx);
0130       float phi = pN->m_spGbts.phi();
0131       m_vPhiNodes.push_back(std::pair<float, unsigned int>(phi, nIdx));
0132     }
0133
0134     for (unsigned int nIdx = 0; nIdx < m_vn.size(); nIdx++) {
0135       GbtsNode<space_point_t> *pN = m_vn.at(nIdx);
0136       float phi = pN->m_spGbts.phi();
0137       if (phi >= -M_PI + dphi) {
0138         break;
0139       }
0140       m_vPhiNodes.push_back(
0141           std::pair<float, unsigned int>(phi + 2 * M_PI, nIdx));
0142     }
0143   }
0144
0145   std::vector<GbtsNode<space_point_t> *> m_vn;
0146   // TODO change to
0147   // std::vector<std::unique_ptr<GbtsNode<space_point_t>>> m_vn;
0148   std::vector<std::pair<float, unsigned int>> m_vPhiNodes;
0149 };
0150
0151 template <typename space_point_t>
0152 class GbtsDataStorage {
0153  public:
0154   GbtsDataStorage(const GbtsGeometry<space_point_t> &g) : m_geo(g) {
0155     m_etaBins.reserve(g.num_bins());
0156     for (int k = 0; k < g.num_bins(); k++) {
0157       m_etaBins.emplace_back(GbtsEtaBin<space_point_t>());
0158     }
0159   }
0160
0161   int addSpacePoint(const GbtsSP<space_point_t> &sp, bool useClusterWidth) {
0162     const GbtsLayer<space_point_t> *pL =
0163         m_geo.getGbtsLayerByKey(sp.combined_ID);
0164
0165     if (pL == nullptr) {
0166       return -1;
0167     }
0168
0169     int binIndex = pL->getEtaBin(sp.SP->z(), sp.SP->r());
0170
0171     if (binIndex == -1) {
0172       return -2;
0173     }
0174
0175     bool isBarrel = (pL->m_layer.m_type == 0);
0176
0177     if (isBarrel) {
0178       float min_tau = -100.0;
0179       float max_tau = 100.0;
0180       // can't do this bit yet as dont have cluster width
0181       if (useClusterWidth) {
0182         float cluster_width = 1;  // temporary while cluster width not available
0183         min_tau = 6.7 * (cluster_width - 0.2);
0184         max_tau =
0185             1.6 + 0.15 / (cluster_width + 0.2) + 6.1 * (cluster_width - 0.2);
0186       }
0187
0188       m_etaBins.at(binIndex).m_vn.push_back(new GbtsNode<space_point_t>(
0189           sp, min_tau, max_tau));  // adding ftf member to nodes
0190     } else {
0191       if (useClusterWidth) {
0192         float cluster_width = 1;  // temporary while cluster width not available
0193         if (cluster_width > 0.2) {
0194           return -3;
0195         }
0196       }
0197       m_etaBins.at(binIndex).m_vn.push_back(new GbtsNode<space_point_t>(sp));
0198     }
0199
0200     return 0;
0201   }
0202
0203   // for safety to prevent passing as copy
0204   GbtsDataStorage(const GbtsDataStorage &) = delete;
0205   GbtsDataStorage &operator=(const GbtsDataStorage &) = delete;
0206
0207   unsigned int numberOfNodes() const {
0208     unsigned int n = 0;
0209
0210     for (auto &b : m_etaBins) {
0211       n += b.m_vn.size();
0212     }
0213     return n;
0214   }
0215
0216   void getConnectingNodes(std::vector<const GbtsNode<space_point_t> *> &vn) {
0217     vn.clear();
0218     vn.reserve(numberOfNodes());
0219     for (const auto &b : m_etaBins) {
0220       for (typename std::vector<GbtsNode<space_point_t> *>::const_iterator nIt =
0221                b.m_vn.begin();
0222            nIt != b.m_vn.end(); ++nIt) {
0223         if ((*nIt)->m_in.empty()) {
0224           continue;
0225         }
0226         if ((*nIt)->m_out.empty()) {
0227           continue;
0228         }
0229         vn.push_back(*nIt);
0230       }
0231     }
0232   }
0233
0234   void sortByPhi() {
0235     for (auto &b : m_etaBins) {
0236       b.sortByPhi();
0237     }
0238   }
0239
0240   void generatePhiIndexing(float dphi) {
0241     for (auto &b : m_etaBins) {
0242       b.generatePhiIndexing(dphi);
0243     }
0244   }
0245
0246   const GbtsEtaBin<space_point_t> &getEtaBin(int idx) const {
0247     if (idx >= static_cast<int>(m_etaBins.size())) {
0248       idx = idx - 1;
0249     }
0250     return m_etaBins.at(idx);
0251   }
0252
0253  protected:
0254   const GbtsGeometry<space_point_t> &m_geo;
0255
0256   std::vector<GbtsEtaBin<space_point_t>> m_etaBins;
0257 };
0258
0259 template <typename space_point_t>
0260 class GbtsEdge {
0261  public:
0262   struct CompareLevel {
0263    public:
0264     bool operator()(const GbtsEdge *pS1, const GbtsEdge *pS2) {
0265       return pS1->m_level > pS2->m_level;
0266     }
0267   };
0268
0269   GbtsEdge(GbtsNode<space_point_t> *n1, GbtsNode<space_point_t> *n2, float p1,
0270            float p2, float p3, float p4)
0271       : m_n1(n1), m_n2(n2), m_level(1), m_next(1) {
0272     m_p[0] = p1;
0273     m_p[1] = p2;
0274     m_p[2] = p3;
0275     m_p[3] = p4;
0276   }
0277
0278   GbtsEdge() : m_n1(nullptr), m_n2(nullptr), m_level(-1), m_next(-1) {}
0279
0280   // GbtsEdge(const GbtsEdge<space_point_t> &e)
0281   //     : m_n1(e.m_n1), m_n2(e.m_n2) {}
0282
0283   // inline void initialize(GbtsNode<space_point_t> *n1,
0284   //                        GbtsNode<space_point_t> *n2) {
0285   //   m_n1 = n1;
0286   //   m_n2 = n2;
0287   //   m_level = 1;
0288   //   m_next = 1;
0289   //   m_nNei = 0;
0290   // }
0291
0292   GbtsNode<space_point_t> *m_n1{nullptr};
0293   GbtsNode<space_point_t> *m_n2{nullptr};
0294
0295   signed char m_level{}, m_next{};
0296
0297   unsigned char m_nNei{0};
0298   float m_p[4]{};
0299
0300   unsigned int m_vNei[N_SEG_CONNS]{};  // global indices of the connected edges
0301 };
0302
0303 }  // namespace Acts