File indexing completed on 2026-05-27 07:24:12
0001
0002
0003
0004
0005
0006
0007
0008
0009 #pragma once
0010
0011
0012 #include "detray/builders/surface_factory_interface.hpp"
0013 #include "detray/core/detail/data_context.hpp"
0014 #include "detray/definitions/algebra.hpp"
0015 #include "detray/definitions/detail/qualifiers.hpp"
0016 #include "detray/definitions/geometry.hpp"
0017 #include "detray/definitions/indexing.hpp"
0018 #include "detray/geometry/mask.hpp"
0019 #include "detray/geometry/shapes/rectangle2D.hpp"
0020 #include "detray/material/material.hpp"
0021 #include "detray/utils/logging.hpp"
0022
0023
0024 #include <cassert>
0025 #include <iostream>
0026 #include <limits>
0027
0028 namespace detray {
0029
0030
0031
0032
0033 template <concepts::scalar scalar_t>
0034 struct endcap_generator_config {
0035
0036 int m_side{-1};
0037
0038 scalar_t m_center_z{0.f};
0039
0040 scalar_t m_inner_radius{25.f * unit<scalar_t>::mm};
0041
0042 scalar_t m_outer_radius{180.f * unit<scalar_t>::mm};
0043
0044 std::vector<std::vector<scalar_t>> m_mask_values{
0045 {3.f * unit<scalar_t>::mm, 9.5f * unit<scalar_t>::mm,
0046 39.f * unit<scalar_t>::mm},
0047 {6.f * unit<scalar_t>::mm, 10.f * unit<scalar_t>::mm,
0048 39.f * unit<scalar_t>::mm}};
0049
0050 scalar_t m_ring_stagger{2.f * unit<scalar_t>::mm};
0051
0052 std::vector<scalar_t> m_phi_stagger = {4.f * unit<scalar_t>::mm,
0053 4.f * unit<scalar_t>::mm};
0054
0055 std::vector<scalar_t> m_phi_sub_stagger = {0.5f * unit<scalar_t>::mm,
0056 0.5f * unit<scalar_t>::mm};
0057
0058 std::vector<scalar_t> m_tilt = {0.f, 0.f};
0059
0060 std::vector<unsigned int> m_binning = {40u, 68u};
0061
0062
0063
0064 constexpr endcap_generator_config &side(const int s) {
0065 m_side = s;
0066 return *this;
0067 }
0068 constexpr endcap_generator_config ¢er(const scalar_t z) {
0069 m_center_z = std::abs(z);
0070 return *this;
0071 }
0072 constexpr endcap_generator_config &inner_radius(const scalar_t r) {
0073 m_inner_radius = r;
0074 return *this;
0075 }
0076 constexpr endcap_generator_config &outer_radius(const scalar_t r) {
0077 m_outer_radius = r;
0078 return *this;
0079 }
0080 endcap_generator_config &module_bounds(
0081 const std::vector<std::vector<scalar_t>> &bnds) {
0082 m_mask_values.clear();
0083 m_mask_values.push_back(bnds.at(0));
0084 m_mask_values.push_back(bnds.at(1));
0085 return *this;
0086 }
0087 constexpr endcap_generator_config &ring_stagger(const scalar_t rs) {
0088 m_ring_stagger = rs;
0089 return *this;
0090 }
0091 constexpr endcap_generator_config &phi_stagger(
0092 const std::vector<scalar_t> &ps) {
0093 m_phi_stagger = ps;
0094 return *this;
0095 }
0096 constexpr endcap_generator_config &phi_sub_stagger(
0097 const std::vector<scalar_t> &ps) {
0098 m_phi_sub_stagger = ps;
0099 return *this;
0100 }
0101 constexpr endcap_generator_config &module_tilt(
0102 const std::vector<scalar_t> &t) {
0103 m_tilt = t;
0104 return *this;
0105 }
0106 constexpr endcap_generator_config &binning(
0107 const std::vector<unsigned int> &b) {
0108 m_binning = b;
0109 return *this;
0110 }
0111
0112
0113
0114
0115 constexpr int side() const { return m_side; }
0116 constexpr scalar_t center() const { return m_center_z; }
0117 constexpr scalar_t inner_radius() const { return m_inner_radius; }
0118 constexpr scalar_t outer_radius() const { return m_outer_radius; }
0119 const auto &module_bounds() const { return m_mask_values; }
0120 constexpr scalar_t ring_stagger() const { return m_ring_stagger; }
0121 constexpr const auto &phi_stagger() const { return m_phi_stagger; }
0122 constexpr const auto &phi_sub_stagger() const { return m_phi_sub_stagger; }
0123 constexpr const auto &module_tilt() const { return m_tilt; }
0124 constexpr const auto &binning() const { return m_binning; }
0125
0126 };
0127
0128
0129
0130
0131 template <typename detector_t, typename mask_shape_t = trapezoid2D>
0132 class endcap_generator final : public surface_factory_interface<detector_t> {
0133 using algebra_t = typename detector_t::algebra_type;
0134 using scalar_t = dscalar<algebra_t>;
0135 using transform3_t = dtransform3D<algebra_t>;
0136 using point3_t = dpoint3D<algebra_t>;
0137 using vector3_t = dvector3D<algebra_t>;
0138
0139 public:
0140
0141 DETRAY_HOST
0142 explicit endcap_generator(const endcap_generator_config<scalar_t> &cfg)
0143 : m_cfg{cfg} {}
0144
0145
0146 DETRAY_HOST
0147 auto size() const -> dindex override {
0148 dindex n_modules{0u};
0149 for (const unsigned int n_bins : m_cfg.binning()) {
0150 n_modules += n_bins;
0151 }
0152 return n_modules;
0153 }
0154
0155
0156
0157 DETRAY_HOST
0158 void clear() override { };
0159 DETRAY_HOST
0160 void push_back(
0161 surface_data<detector_t> && ) override { }
0162 DETRAY_HOST
0163 auto push_back(std::vector<surface_data<detector_t>> && )
0164 -> void override { }
0165
0166
0167
0168
0169
0170
0171
0172
0173
0174 DETRAY_HOST
0175 auto operator()(typename detector_t::volume_type &volume,
0176 typename detector_t::surface_lookup_container &surfaces,
0177 typename detector_t::transform_container &transforms,
0178 typename detector_t::mask_container &masks,
0179 typename detector_t::geometry_context ctx = {})
0180 -> dindex_range override {
0181 DETRAY_VERBOSE_HOST("Generate silicon tracker endcap modules...");
0182 DETRAY_VERBOSE_HOST("-> Generate " << size() << " surfaces");
0183
0184 using surface_t = typename detector_t::surface_type;
0185 using nav_link_t = typename surface_t::navigation_link;
0186 using mask_link_t = typename surface_t::mask_link;
0187 using material_link_t = typename surface_t::material_link;
0188
0189 auto volume_idx{volume.index()};
0190 const dindex surfaces_offset{static_cast<dindex>(surfaces.size())};
0191 constexpr auto invalid_src_link{detail::invalid_value<std::uint64_t>()};
0192
0193
0194 constexpr auto mask_id{
0195 types::id<typename detector_t::masks, mask<mask_shape_t, algebra_t>>};
0196
0197 constexpr auto no_material{surface_t::material_id::e_none};
0198
0199 const auto mask_volume_link{static_cast<nav_link_t>(volume_idx)};
0200
0201
0202 std::vector<scalar_t> radii{};
0203
0204 scalar_t delta_r{m_cfg.outer_radius() - m_cfg.inner_radius()};
0205
0206
0207
0208
0209
0210
0211 if (m_cfg.binning().size() == 1u) {
0212 radii.push_back(0.5f * (m_cfg.inner_radius() + m_cfg.outer_radius()));
0213 } else {
0214
0215 scalar_t tot_length{0.f};
0216 for (const auto &bounds : m_cfg.module_bounds()) {
0217
0218
0219 tot_length += 2.f * bounds.at(trapezoid2D::e_half_length_2) +
0220 1.f * unit<scalar_t>::mm;
0221 }
0222
0223
0224 scalar_t r_overlap{(tot_length - delta_r) /
0225 static_cast<scalar_t>(m_cfg.module_bounds().size())};
0226
0227
0228 scalar_t prev_r{m_cfg.inner_radius() + r_overlap};
0229 scalar_t prev_hl{0.f};
0230 scalar_t prev_ol{0.f};
0231
0232 for (const auto &bounds : m_cfg.module_bounds()) {
0233 const scalar_t mod_hlength{bounds.at(trapezoid2D::e_half_length_2)};
0234
0235 radii.push_back(prev_r + prev_hl - prev_ol + mod_hlength);
0236 prev_r = radii.back();
0237 prev_ol = 2.f * r_overlap;
0238 prev_hl = mod_hlength;
0239 }
0240 }
0241
0242
0243
0244
0245 for (unsigned int ir = 0u; ir < radii.size(); ++ir) {
0246
0247
0248 const scalar_t center{m_cfg.center()};
0249 const scalar_t staggered_center{
0250 (ir % 2u != 0 ? center + 0.5f * m_cfg.ring_stagger()
0251 : center - 0.5f * m_cfg.ring_stagger())};
0252 const scalar_t rz{radii.size() == 1u ? center : staggered_center};
0253
0254
0255 const scalar_t sub_stagger{m_cfg.phi_sub_stagger().size() > 1u
0256 ? m_cfg.phi_sub_stagger().at(ir)
0257 : 0.f};
0258
0259 std::vector<point3_t> module_positions =
0260 module_positions_ring(rz, radii.at(ir), m_cfg.phi_stagger().at(ir),
0261 sub_stagger, m_cfg.binning().at(ir));
0262
0263
0264 for (const point3_t &mod_position : module_positions) {
0265
0266 mask_link_t mask_link{mask_id, {masks.template size<mask_id>(), 1u}};
0267 material_link_t material_link{no_material, dindex_invalid};
0268
0269
0270 surfaces.push_back({transforms.size(ctx), mask_link, material_link,
0271 volume_idx, surface_id::e_sensitive},
0272 invalid_src_link);
0273
0274
0275
0276
0277 point3_t mod_center{mod_position};
0278 mod_center[2] *= static_cast<scalar_t>(m_cfg.side());
0279
0280 const scalar_t mod_phi{vector::phi(mod_position)};
0281 const vector3_t mod_loc_y{math::cos(mod_phi), math::sin(mod_phi),
0282 static_cast<scalar_t>(0)};
0283
0284 const vector3_t mod_loc_z{static_cast<scalar_t>(0),
0285 static_cast<scalar_t>(0),
0286 static_cast<scalar_t>(m_cfg.side())};
0287 const vector3_t mod_loc_x{vector::cross(mod_loc_y, mod_loc_z)};
0288
0289
0290 transforms.emplace_back(ctx, mod_center, mod_loc_z, mod_loc_x);
0291 }
0292
0293
0294 std::vector<scalar_t> mask_values{m_cfg.module_bounds().at(ir)};
0295
0296
0297 if constexpr (std::is_same_v<mask_shape_t, trapezoid2D>) {
0298 const scalar_t div{
0299 1.f / (2.f * mask_values.at(trapezoid2D::e_half_length_2))};
0300
0301 mask_values.insert(mask_values.begin() + trapezoid2D::e_divisor, div);
0302 }
0303
0304 masks.template emplace_back<mask_id>(empty_context{}, mask_values,
0305 mask_volume_link);
0306 }
0307
0308 return {surfaces_offset, static_cast<dindex>(surfaces.size())};
0309 }
0310
0311 private:
0312
0313
0314
0315
0316
0317
0318
0319
0320
0321 auto module_positions_ring(const scalar_t z, const scalar_t radius,
0322 const scalar_t phi_stagger,
0323 const scalar_t phi_sub_stagger,
0324 const unsigned int n_phi_bins) const {
0325
0326 std::vector<point3_t> mod_positions;
0327 mod_positions.reserve(n_phi_bins);
0328
0329
0330 const scalar_t phi_step{2.f * constant<scalar_t>::pi /
0331 static_cast<scalar_t>(n_phi_bins)};
0332 const scalar_t min_phi{-constant<scalar_t>::pi + 0.5f * phi_step};
0333
0334 for (unsigned int iphi = 0u; iphi < n_phi_bins; ++iphi) {
0335
0336 scalar_t rzs{0.f};
0337
0338
0339
0340 if (phi_sub_stagger != 0.f && (n_phi_bins % 4u == 0u)) {
0341
0342 if ((iphi % 4u) == 0u) {
0343 rzs = phi_sub_stagger;
0344 } else if (((iphi + 1u) % 4u) == 0u) {
0345 rzs = -phi_sub_stagger;
0346 }
0347 }
0348
0349 const scalar_t phi{min_phi + static_cast<scalar_t>(iphi) * phi_step};
0350
0351 const scalar_t rz{iphi % 2u != 0u ? z - 0.5f * phi_stagger
0352 : z + 0.5f * phi_stagger};
0353 mod_positions.push_back(
0354 {radius * math::cos(phi), radius * math::sin(phi), rz + rzs});
0355 }
0356 return mod_positions;
0357 }
0358
0359
0360 endcap_generator_config<scalar_t> m_cfg{};
0361 };
0362
0363 }