Back to home page

EIC code displayed by LXR
 
 

    

File indexing completed on 2025-01-18 09:11:17 

0001 // This file is part of the ACTS project.
0002 //
0003 // Copyright (C) 2016 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 https://mozilla.org/MPL/2.0/.
0008 
0009 #include "Acts/Detector/detail/BlueprintHelper.hpp"
0010 
0011 #include "Acts/Definitions/Algebra.hpp"
0012 #include "Acts/Definitions/Tolerance.hpp"
0013 #include "Acts/Geometry/VolumeBounds.hpp"
0014 
0015 #include <algorithm>
0016 #include <array>
0017 
0018 namespace {
0019 
0020 std::array<Acts::Vector3, 2u> endPointsXYZ(
0021     const Acts::Experimental::Blueprint::Node& node, Acts::AxisDirection bVal) {
0022   unsigned int bIdx = 0;
0023   switch (bVal) {
0024     case Acts::AxisDirection::AxisX:
0025       bIdx = 0;
0026       break;
0027     case Acts::AxisDirection::AxisY:
0028       bIdx = 1;
0029       break;
0030     case Acts::AxisDirection::AxisZ:
0031       bIdx = 2;
0032       break;
0033     default:
0034       break;
0035   }
0036   Acts::Vector3 axis = node.transform.rotation().col(bIdx);
0037   auto halfL = node.boundaryValues[bIdx];
0038   Acts::Vector3 center = node.transform.translation();
0039   Acts::Vector3 p0 = center - halfL * axis;
0040   Acts::Vector3 p1 = center + halfL * axis;
0041   return {p0, p1};
0042 }
0043 
0044 }  // namespace
0045 
0046 void Acts::Experimental::detail::BlueprintHelper::sort(Blueprint::Node& node,
0047                                                        bool recursive) {
0048   if (node.children.size() < 2u) {
0049     return;
0050   }
0051   // Sort along x, y, z
0052   if (node.binning.size() == 1) {
0053     auto bVal = node.binning.front();
0054     // x,y,z binning along the axis
0055     if (bVal == AxisDirection::AxisX || bVal == AxisDirection::AxisY ||
0056         bVal == AxisDirection::AxisZ) {
0057       Vector3 nodeCenter = node.transform.translation();
0058       Vector3 nodeSortAxis = node.transform.rotation().col(toUnderlying(bVal));
0059       std::ranges::sort(node.children, {}, [&](const auto& c) {
0060         return (c->transform.translation() - nodeCenter).dot(nodeSortAxis);
0061       });
0062     } else if (bVal == AxisDirection::AxisR &&
0063                node.boundsType == VolumeBounds::eCylinder) {
0064       std::ranges::sort(node.children, {}, [](const auto& c) {
0065         return c->boundaryValues[0] + c->boundaryValues[1];
0066       });
0067     }
0068   }
0069 
0070   // Sort the children
0071   if (recursive) {
0072     for (auto& child : node.children) {
0073       sort(*child, true);
0074     }
0075   }
0076 }
0077 
0078 void Acts::Experimental::detail::BlueprintHelper::fillGaps(
0079     Blueprint::Node& node, bool adjustToParent) {
0080   // Return if this is a leaf node
0081   if (node.isLeaf()) {
0082     return;
0083   }
0084   if (node.boundsType == VolumeBounds::eCylinder && node.binning.size() == 1) {
0085     fillGapsCylindrical(node, adjustToParent);
0086   } else if (node.boundsType == VolumeBounds::eCuboid &&
0087              node.binning.size() == 1) {
0088     // Doesn't look like NOT adjusting to parent
0089     // makes sense. The gaps are not going
0090     // to be filled in non-binned directions
0091     fillGapsCuboidal(node, adjustToParent);
0092   } else {
0093     throw std::runtime_error(
0094         "BlueprintHelper: gap filling is not implemented for "
0095         "this boundary type");
0096   }
0097 }
0098 
0099 void Acts::Experimental::detail::BlueprintHelper::fillGapsCylindrical(
0100     Blueprint::Node& node, bool adjustToParent) {
0101   // Nodes must be sorted
0102   sort(node, false);
0103 
0104   // Container values
0105   auto cInnerR = node.boundaryValues[0];
0106   auto cOuterR = node.boundaryValues[1];
0107   auto cHalfZ = node.boundaryValues[2];
0108 
0109   std::vector<std::unique_ptr<Blueprint::Node>> gaps;
0110   // Only 1D binning implemented for the moment
0111   if (AxisDirection bVal = node.binning.front(); bVal == AxisDirection::AxisZ) {
0112     // adjust inner/outer radius
0113     if (adjustToParent) {
0114       std::for_each(node.children.begin(), node.children.end(),
0115                     [&](auto& child) {
0116                       child->boundaryValues[0] = cInnerR;
0117                       child->boundaryValues[1] = cOuterR;
0118                     });
0119     }
0120     auto [negC, posC] = endPointsXYZ(node, bVal);
0121     // Assume sorted along the local z axis
0122     unsigned int igap = 0;
0123     for (auto& child : node.children) {
0124       auto [neg, pos] = endPointsXYZ(*child, bVal);
0125       double gapSpan = (neg - negC).norm();
0126       if (gapSpan > s_onSurfaceTolerance) {
0127         // Fill a gap node
0128         auto gapName = node.name + "_gap_" + std::to_string(igap);
0129         auto gapTransform = Transform3::Identity();
0130         gapTransform.rotate(node.transform.rotation());
0131         gapTransform.pretranslate(0.5 * (neg + negC));
0132         auto gap = std::make_unique<Blueprint::Node>(
0133             gapName, gapTransform, VolumeBounds::eCylinder,
0134             std::vector<double>{cInnerR, cOuterR, 0.5 * gapSpan});
0135         gaps.push_back(std::move(gap));
0136         ++igap;
0137       }
0138       // Set to new current negative value
0139       negC = pos;
0140     }
0141     // Check if a last one needs to be filled
0142     double gapSpan = (negC - posC).norm();
0143     if (gapSpan > s_onSurfaceTolerance) {
0144       // Fill a gap node
0145       auto gapName = node.name + "_gap_" + std::to_string(igap);
0146       auto gapTransform = Transform3::Identity();
0147       gapTransform.rotate(node.transform.rotation());
0148       gapTransform.pretranslate(0.5 * (negC + posC));
0149       auto gap = std::make_unique<Blueprint::Node>(
0150           gapName, gapTransform, VolumeBounds::eCylinder,
0151           std::vector<double>{cInnerR, cOuterR, 0.5 * gapSpan});
0152       gaps.push_back(std::move(gap));
0153     }
0154 
0155   } else if (bVal == AxisDirection::AxisR) {
0156     // We have binning in R present
0157     if (adjustToParent) {
0158       std::for_each(node.children.begin(), node.children.end(),
0159                     [&](auto& child) {
0160                       child->transform = node.transform;
0161                       child->boundaryValues[2] = cHalfZ;
0162                     });
0163     }
0164     // Fill the gaps in R
0165     unsigned int igap = 0;
0166     double lastR = cInnerR;
0167     for (auto& child : node.children) {
0168       double iR = child->boundaryValues[0];
0169       if (std::abs(iR - lastR) > s_onSurfaceTolerance) {
0170         auto gap = std::make_unique<Blueprint::Node>(
0171             node.name + "_gap_" + std::to_string(igap), node.transform,
0172             VolumeBounds::eCylinder, std::vector<double>{lastR, iR, cHalfZ});
0173         gaps.push_back(std::move(gap));
0174         ++igap;
0175       }
0176       // Set to new outer radius
0177       lastR = child->boundaryValues[1];
0178     }
0179     // Check if a last one needs to be filled
0180     if (std::abs(lastR - cOuterR) > s_onSurfaceTolerance) {
0181       auto gap = std::make_unique<Blueprint::Node>(
0182           node.name + "_gap_" + std::to_string(igap), node.transform,
0183           VolumeBounds::eCylinder, std::vector<double>{lastR, cOuterR, cHalfZ});
0184       gaps.push_back(std::move(gap));
0185     }
0186   } else {
0187     throw std::runtime_error(
0188         "BlueprintHelper: gap filling not implemented for "
0189         "cylinder and this binning type.");
0190   }
0191 
0192   // Insert
0193   for (auto& gap : gaps) {
0194     node.add(std::move(gap));
0195   }
0196 
0197   // Sort again after inserting
0198   sort(node, false);
0199   // Fill the gaps recursively
0200   for (auto& child : node.children) {
0201     fillGaps(*child, adjustToParent);
0202   }
0203 }
0204 
0205 void Acts::Experimental::detail::BlueprintHelper::fillGapsCuboidal(
0206     Blueprint::Node& node, bool adjustToParent) {
0207   // Nodes must be sorted
0208   sort(node, false);
0209 
0210   // Cuboidal detector binnings
0211   std::array<Acts::AxisDirection, 3u> allowedBinVals = {
0212       AxisDirection::AxisX, AxisDirection::AxisY, AxisDirection::AxisZ};
0213 
0214   std::vector<std::unique_ptr<Blueprint::Node>> gaps;
0215   auto binVal = node.binning.front();
0216 
0217   // adjust non-binned directions
0218   if (adjustToParent) {
0219     std::for_each(node.children.begin(), node.children.end(), [&](auto& child) {
0220       for (auto bv : allowedBinVals) {
0221         if (bv != binVal) {
0222           // Both boundary values and translation
0223           // have to be adjusted
0224           child->boundaryValues[toUnderlying(bv)] =
0225               node.boundaryValues[toUnderlying(bv)];
0226           child->transform.translation()[toUnderlying(bv)] =
0227               node.transform.translation()[toUnderlying(bv)];
0228         }
0229       }
0230     });
0231   }
0232   auto [negC, posC] = endPointsXYZ(node, binVal);
0233 
0234   // Assume sorted along the local binned axis
0235   unsigned int igap = 0;
0236   for (auto& child : node.children) {
0237     auto [neg, pos] = endPointsXYZ(*child, binVal);
0238     double gapSpan = (neg - negC).norm();
0239     if (gapSpan > s_onSurfaceTolerance) {
0240       // Fill a gap node
0241       auto gapName = node.name + "_gap_" + std::to_string(igap);
0242       auto gapTransform = Transform3::Identity();
0243       gapTransform.rotate(node.transform.rotation());
0244       gapTransform.pretranslate(0.5 * (neg + negC));
0245       std::vector<double> gapBounds{0, 0, 0};
0246       gapBounds[toUnderlying(binVal)] = 0.5 * gapSpan;
0247       for (auto bv : allowedBinVals) {
0248         if (bv != binVal) {
0249           gapBounds[toUnderlying(bv)] = node.boundaryValues[toUnderlying(bv)];
0250         }
0251       }
0252       auto gap = std::make_unique<Blueprint::Node>(
0253           gapName, gapTransform, VolumeBounds::eCuboid, gapBounds);
0254       gaps.push_back(std::move(gap));
0255       ++igap;
0256     }
0257     // Set to new current negative value
0258     negC = pos;
0259   }
0260   // Check if a last one needs to be filled
0261   double gapSpan = (negC - posC).norm();
0262   if (gapSpan > s_onSurfaceTolerance) {
0263     // Fill a gap node
0264     auto gapName = node.name + "_gap_" + std::to_string(igap);
0265     auto gapTransform = Transform3::Identity();
0266     gapTransform.rotate(node.transform.rotation());
0267     gapTransform.pretranslate(0.5 * (negC + posC));
0268     std::vector<double> gapBounds{0, 0, 0};
0269     gapBounds[toUnderlying(binVal)] = 0.5 * gapSpan;
0270     for (auto bv : allowedBinVals) {
0271       if (bv != binVal) {
0272         gapBounds[toUnderlying(bv)] = node.boundaryValues[toUnderlying(bv)];
0273       }
0274     }
0275     auto gap = std::make_unique<Blueprint::Node>(
0276         gapName, gapTransform, VolumeBounds::eCuboid, gapBounds);
0277     gaps.push_back(std::move(gap));
0278   }
0279 
0280   // Insert
0281   for (auto& gap : gaps) {
0282     node.add(std::move(gap));
0283   }
0284 
0285   // Sort again after inserting
0286   sort(node, false);
0287   // Fill the gaps recursively
0288   for (auto& child : node.children) {
0289     fillGaps(*child, adjustToParent);
0290   }
0291 }
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.3.7 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!