CommonHelpers/src/CylindricalTrackingGeometry.cpp

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 "ActsTests/CommonHelpers/CylindricalTrackingGeometry.hpp"
0010
0011 #include "Acts/Geometry/Blueprint.hpp"
0012 #include "Acts/Geometry/BlueprintOptions.hpp"
0013 #include "Acts/Geometry/ContainerBlueprintNode.hpp"
0014 #include "Acts/Geometry/CylinderVolumeBounds.hpp"
0015 #include "Acts/Geometry/CylinderVolumeBuilder.hpp"
0016 #include "Acts/Geometry/CylinderVolumeHelper.hpp"
0017 #include "Acts/Geometry/LayerBlueprintNode.hpp"
0018 #include "Acts/Geometry/MaterialDesignatorBlueprintNode.hpp"
0019 #include "Acts/Geometry/PassiveLayerBuilder.hpp"
0020 #include "Acts/Geometry/SurfaceArrayCreator.hpp"
0021 #include "Acts/Geometry/TrackingVolumeArrayCreator.hpp"
0022 #include "Acts/Material/BinnedSurfaceMaterial.hpp"
0023 #include "Acts/Material/HomogeneousSurfaceMaterial.hpp"
0024 #include "Acts/Navigation/SurfaceArrayNavigationPolicy.hpp"
0025 #include "Acts/Navigation/TryAllNavigationPolicy.hpp"
0026 #include "Acts/Surfaces/TrapezoidBounds.hpp"
0027 #include "Acts/Utilities/AxisDefinitions.hpp"
0028 #include "ActsTests/CommonHelpers/PredefinedMaterials.hpp"
0029
0030 using namespace Acts;
0031
0032 namespace ActsTests {
0033
0034 std::vector<Surface*> CylindricalTrackingGeometry::surfacesRing(
0035     DetectorStore& detStore, double moduleHalfXminY, double moduleHalfXmaxY,
0036     double moduleHalfY, double moduleThickness, double moduleTilt,
0037     double ringRadius, double ringZ, double zStagger, int nPhi) {
0038   std::vector<Surface*> layerSurfaces;
0039
0040   // Module material from input
0041   MaterialSlab moduleMaterial(makeSilicon(), moduleThickness);
0042
0043   // Create a new surface material
0044   std::shared_ptr<const ISurfaceMaterial> moduleMaterialPtr =
0045       std::shared_ptr<const ISurfaceMaterial>(
0046           new HomogeneousSurfaceMaterial(moduleMaterial));
0047
0048   // The rectangle/trapezoid bounds for all modules
0049   std::shared_ptr<PlanarBounds> mBounds = nullptr;
0050   if (moduleHalfXminY == moduleHalfXmaxY) {
0051     mBounds = std::make_shared<RectangleBounds>(moduleHalfXminY, moduleHalfY);
0052   } else {
0053     mBounds = std::make_shared<TrapezoidBounds>(moduleHalfXminY,
0054                                                 moduleHalfXmaxY, moduleHalfY);
0055   }
0056
0057   double phiStep = 2 * std::numbers::pi / nPhi;
0058
0059   for (int im = 0; im < nPhi; ++im) {
0060     // Get the moduleTransform
0061     double phi = -std::numbers::pi + im * phiStep;
0062     auto mModuleTransform = Transform3(
0063         Translation3(ringRadius * std::cos(phi), ringRadius * std::sin(phi),
0064                      ringZ + (im % 2) * zStagger) *
0065         AngleAxis3(phi - std::numbers::pi / 2., Vector3::UnitZ()) *
0066         AngleAxis3(moduleTilt, Vector3::UnitY()));
0067
0068     // Create the detector element
0069     auto detElement = std::make_unique<DetectorElementStub>(
0070         mModuleTransform, mBounds, moduleThickness, moduleMaterialPtr);
0071
0072     layerSurfaces.push_back(&detElement->surface());
0073     detStore.push_back(std::move(detElement));
0074   }
0075
0076   return layerSurfaces;
0077 }
0078
0079 std::vector<Surface*> CylindricalTrackingGeometry::surfacesCylinder(
0080     DetectorStore& detStore, double moduleHalfX, double moduleHalfY,
0081     double moduleThickness, double moduleTiltPhi, double layerRadius,
0082     double radialStagger, double longitudinalOverlap,
0083     const std::pair<int, int>& binningSchema) {
0084   std::vector<Surface*> layerSurfaces;
0085
0086   // Module material from input
0087   MaterialSlab moduleMaterial(makeSilicon(), moduleThickness);
0088
0089   // Create a new surface material
0090   std::shared_ptr<const ISurfaceMaterial> moduleMaterialPtr =
0091       std::shared_ptr<const ISurfaceMaterial>(
0092           new HomogeneousSurfaceMaterial(moduleMaterial));
0093
0094   // The rectangle bounds for all modules
0095   auto mBounds = std::make_shared<RectangleBounds>(moduleHalfX, moduleHalfY);
0096
0097   // Create the module centers
0098   auto moduleCenters =
0099       modulePositionsCylinder(layerRadius, radialStagger, moduleHalfY,
0100                               longitudinalOverlap, binningSchema);
0101
0102   for (auto& mCenter : moduleCenters) {
0103     // The association transform
0104     double modulePhi = VectorHelpers::phi(mCenter);
0105     // Local z axis is the normal vector
0106     Vector3 moduleLocalZ(std::cos(modulePhi + moduleTiltPhi),
0107                          std::sin(modulePhi + moduleTiltPhi), 0.);
0108     // Local y axis is the global z axis
0109     Vector3 moduleLocalY(0., 0., 1);
0110     // Local x axis the normal to local y,z
0111     Vector3 moduleLocalX(-std::sin(modulePhi + moduleTiltPhi),
0112                          std::cos(modulePhi + moduleTiltPhi), 0.);
0113     // Create the RotationMatrix
0114     RotationMatrix3 moduleRotation;
0115     moduleRotation.col(0) = moduleLocalX;
0116     moduleRotation.col(1) = moduleLocalY;
0117     moduleRotation.col(2) = moduleLocalZ;
0118     // Get the moduleTransform
0119     auto mModuleTransform = Transform3(Translation3(mCenter) * moduleRotation);
0120     // Create the detector element
0121     auto detElement = std::make_unique<DetectorElementStub>(
0122         mModuleTransform, mBounds, moduleThickness, moduleMaterialPtr);
0123
0124     layerSurfaces.push_back(&detElement->surface());
0125     detStore.push_back(std::move(detElement));
0126   }
0127   return layerSurfaces;
0128 }
0129
0130 /// Helper method for cylinder layer
0131 /// create the positions for module surfaces on a cylinder
0132 std::vector<Vector3> CylindricalTrackingGeometry::modulePositionsCylinder(
0133     double radius, double zStagger, double moduleHalfLength, double lOverlap,
0134     const std::pair<int, int>& binningSchema) {
0135   int nPhiBins = binningSchema.first;
0136   int nZbins = binningSchema.second;
0137   // prepare the return value
0138   std::vector<Vector3> mPositions;
0139   mPositions.reserve(nPhiBins * nZbins);
0140   // prep work
0141   double phiStep = 2 * std::numbers::pi / (nPhiBins);
0142   double minPhi = -std::numbers::pi + phiStep / 2.;
0143   double zStart = -0.5 * (nZbins - 1) * (2 * moduleHalfLength - lOverlap);
0144   double zStep = 2 * std::abs(zStart) / (nZbins - 1);
0145   // loop over the bins
0146   for (std::size_t zBin = 0; zBin < static_cast<std::size_t>(nZbins); ++zBin) {
0147     // prepare z and r
0148     double moduleZ = zStart + zBin * zStep;
0149     double moduleR =
0150         (zBin % 2) != 0u ? radius - 0.5 * zStagger : radius + 0.5 * zStagger;
0151     for (std::size_t phiBin = 0; phiBin < static_cast<std::size_t>(nPhiBins);
0152          ++phiBin) {
0153       // calculate the current phi value
0154       double modulePhi = minPhi + phiBin * phiStep;
0155       mPositions.push_back(Vector3(moduleR * std::cos(modulePhi),
0156                                    moduleR * std::sin(modulePhi), moduleZ));
0157     }
0158   }
0159   return mPositions;
0160 }
0161
0162 std::shared_ptr<TrackingGeometry> CylindricalTrackingGeometry::buildGen1(
0163     const Logger& /*logger*/) {
0164   using namespace UnitLiterals;
0165
0166   Logging::Level surfaceLLevel = Logging::INFO;
0167   Logging::Level layerLLevel = Logging::INFO;
0168   Logging::Level volumeLLevel = Logging::INFO;
0169
0170   // configure surface array creator
0171   auto surfaceArrayCreator = std::make_shared<const SurfaceArrayCreator>(
0172       getDefaultLogger("SurfaceArrayCreator", surfaceLLevel));
0173   // configure the layer creator that uses the surface array creator
0174   LayerCreator::Config lcConfig;
0175   lcConfig.surfaceArrayCreator = surfaceArrayCreator;
0176   auto layerCreator = std::make_shared<const LayerCreator>(
0177       lcConfig, getDefaultLogger("LayerCreator", layerLLevel));
0178   // configure the layer array creator
0179   LayerArrayCreator::Config lacConfig;
0180   auto layerArrayCreator = std::make_shared<const LayerArrayCreator>(
0181       lacConfig, getDefaultLogger("LayerArrayCreator", layerLLevel));
0182
0183   // tracking volume array creator
0184   TrackingVolumeArrayCreator::Config tvacConfig;
0185   auto tVolumeArrayCreator = std::make_shared<const TrackingVolumeArrayCreator>(
0186       tvacConfig, getDefaultLogger("TrackingVolumeArrayCreator", volumeLLevel));
0187   // configure the cylinder volume helper
0188   CylinderVolumeHelper::Config cvhConfig;
0189   cvhConfig.layerArrayCreator = layerArrayCreator;
0190   cvhConfig.trackingVolumeArrayCreator = tVolumeArrayCreator;
0191   auto cylinderVolumeHelper = std::make_shared<const CylinderVolumeHelper>(
0192       cvhConfig, getDefaultLogger("CylinderVolumeHelper", volumeLLevel));
0193
0194   // ----------------- build a beam pipe -----------------------------------
0195   MaterialSlab beamPipeMaterial(makeBeryllium(), 0.8_mm);
0196   PassiveLayerBuilder::Config bplConfig;
0197   bplConfig.layerIdentification = "BeamPipe";
0198   bplConfig.centralLayerRadii = std::vector<double>(1, kBeamPipeRadius);
0199   bplConfig.centralLayerHalflengthZ =
0200       std::vector<double>(1, kBeamPipeHalfLengthZ);
0201   bplConfig.centralLayerThickness = std::vector<double>(1, kBeamPipeThickness);
0202   bplConfig.centralLayerMaterial = {
0203       std::make_shared<const HomogeneousSurfaceMaterial>(beamPipeMaterial)};
0204   auto beamPipeBuilder = std::make_shared<const PassiveLayerBuilder>(
0205       bplConfig, getDefaultLogger("BeamPipeLayerBuilder", layerLLevel));
0206   // create the volume for the beam pipe
0207   CylinderVolumeBuilder::Config bpvConfig;
0208   bpvConfig.trackingVolumeHelper = cylinderVolumeHelper;
0209   bpvConfig.volumeName = "BeamPipe";
0210   bpvConfig.layerBuilder = beamPipeBuilder;
0211   bpvConfig.layerEnvelopeR = {1_mm, 1_mm};
0212   bpvConfig.buildToRadiusZero = true;
0213   auto beamPipeVolumeBuilder = std::make_shared<const CylinderVolumeBuilder>(
0214       bpvConfig, getDefaultLogger("BeamPipeVolumeBuilder", volumeLLevel));
0215
0216   // create the bounds and the volume
0217   auto beamPipeBounds =
0218       std::make_shared<const CylinderVolumeBounds>(0., 25., 1100.);
0219   auto beamPipeVolume = beamPipeVolumeBuilder->trackingVolume(
0220       geoContext, nullptr, beamPipeBounds);
0221
0222   //-------------------------------------------------------------------------------------
0223   // some prep work for the material
0224   // Layer material properties - thickness, X0, L0, A, Z, Rho
0225   MaterialSlab lProperties(makeSilicon(), 1.5_mm);
0226
0227   auto layerMaterialPtr =
0228       std::make_shared<HomogeneousSurfaceMaterial>(lProperties);
0229
0230   std::vector<LayerPtr> pLayers;
0231
0232   for (std::size_t ilp = 0; ilp < kLayerRadii.size(); ++ilp) {
0233     std::vector<Surface*> layerSurfacesMutable =
0234         surfacesCylinder(detectorStore, kModuleHalfX[ilp], kModuleHalfY[ilp],
0235                          kModuleThickness[ilp], kModuleTiltPhi[ilp],
0236                          kLayerRadii[ilp], 2_mm, 5_mm, kLayerBinning[ilp]);
0237
0238     std::vector<const Surface*> layerSurfaces;
0239     std::ranges::transform(layerSurfacesMutable,
0240                            std::back_inserter(layerSurfaces),
0241                            [](Surface* s) { return s; });
0242
0243     // Make a shared version out of it
0244     std::vector<std::shared_ptr<const Surface>> layerSurfacePtrs;
0245     layerSurfacePtrs.reserve(layerSurfaces.size());
0246     for (auto& sf : layerSurfaces) {
0247       layerSurfacePtrs.push_back(sf->getSharedPtr());
0248     }
0249
0250     // create the layer and store it
0251     ProtoLayer protoLayer(geoContext, layerSurfaces);
0252     protoLayer.envelope[AxisDirection::AxisR] = {0.5, 0.5};
0253     auto pLayer = layerCreator->cylinderLayer(
0254         geoContext, std::move(layerSurfacePtrs), kLayerBinning[ilp].first,
0255         kLayerBinning[ilp].second, protoLayer);
0256     auto approachSurfaces = pLayer->approachDescriptor()->containedSurfaces();
0257     auto mutableOuterSurface = const_cast<Surface*>(approachSurfaces.at(1));
0258     mutableOuterSurface->assignSurfaceMaterial(layerMaterialPtr);
0259     /// now push back the layer
0260     pLayers.push_back(pLayer);
0261
0262   }  // loop over layers
0263
0264   // layer array
0265   auto pLayerArray = layerArrayCreator->layerArray(
0266       geoContext, pLayers, 25., 300., arbitrary, AxisDirection::AxisR);
0267   auto pVolumeBounds = std::make_shared<CylinderVolumeBounds>(25., 300., 1100.);
0268   // create the Tracking volume
0269   auto pVolume = std::make_shared<TrackingVolume>(
0270       Transform3::Identity(), pVolumeBounds, nullptr, std::move(pLayerArray),
0271       nullptr, MutableTrackingVolumeVector{}, "Pixel::Barrel");
0272
0273   // The combined volume
0274   auto detectorVolume = cylinderVolumeHelper->createContainerTrackingVolume(
0275       geoContext, {beamPipeVolume, pVolume});
0276
0277   // create and return the geometry
0278   return std::make_shared<TrackingGeometry>(detectorVolume);
0279 }
0280
0281 std::shared_ptr<TrackingGeometry> CylindricalTrackingGeometry::buildGen3(
0282     const Logger& logger) {
0283   using namespace Experimental;
0284   using namespace UnitLiterals;
0285   using enum CylinderVolumeBounds::Face;
0286   using enum AxisDirection;
0287   using LayerType = LayerBlueprintNode::LayerType;
0288
0289   const MaterialSlab lProperties(makeSilicon(), 1.5_mm);
0290
0291   // Create a binned material in 2 bins - irregularly in z, 2 bins in phi
0292   std::vector<float> binEdges = {// empirical bin edges. these are not checked!
0293                                  -476.5, 0, 476.5};
0294
0295   BinUtility binUtility(2u, -std::numbers::pi, std::numbers::pi,
0296                         BinningOption::closed, AxisPhi);
0297
0298   binUtility += Acts::BinUtility(binEdges, BinningOption::open, AxisZ);
0299
0300   std::vector<MaterialSlab> materialSlabs0 = {lProperties, lProperties};
0301   std::vector<MaterialSlab> materialSlabs1 = {lProperties, lProperties};
0302
0303   auto binnedMaterial = std::make_shared<BinnedSurfaceMaterial>(
0304       binUtility, std::vector{materialSlabs0, materialSlabs1}, 0.,
0305       MappingType::Default);
0306
0307   Blueprint::Config cfg;
0308   cfg.envelope = ExtentEnvelope{{
0309       .z = {20_mm, 20_mm},
0310       .r = {0_mm, 20_mm},
0311   }};
0312   Blueprint root{cfg};
0313
0314   root.addCylinderContainer("Detector", AxisR, [&](auto& detector) {
0315     auto beampipeBounds =
0316         std::make_unique<CylinderVolumeBounds>(0_mm, kBeamPipeRadius, 100_mm);
0317     auto beampipe = std::make_unique<TrackingVolume>(
0318         Transform3::Identity(), std::move(beampipeBounds), "Beampipe");
0319
0320     detector.addMaterial("BeampipeMaterial", [&](auto& bpMat) {
0321       MaterialSlab beamPipeMaterial(makeBeryllium(), kBeamPipeThickness);
0322       bpMat.configureFace(OuterCylinder, binnedMaterial);
0323       bpMat.addStaticVolume(std::move(beampipe));
0324     });
0325
0326     auto layerMaterialPtr =
0327         std::make_shared<HomogeneousSurfaceMaterial>(lProperties);
0328
0329     for (std::size_t ilp = 0; ilp < kLayerRadii.size(); ++ilp) {
0330       std::vector<Surface*> layerSurfaces =
0331           surfacesCylinder(detectorStore, kModuleHalfX[ilp], kModuleHalfY[ilp],
0332                            kModuleThickness[ilp], kModuleTiltPhi[ilp],
0333                            kLayerRadii[ilp], 2_mm, 5_mm, kLayerBinning[ilp]);
0334
0335       // Make a shared version out of it
0336       std::vector<std::shared_ptr<Surface>> layerSurfacePtrs;
0337       layerSurfacePtrs.reserve(layerSurfaces.size());
0338       for (auto& sf : layerSurfaces) {
0339         layerSurfacePtrs.push_back(sf->getSharedPtr());
0340       }
0341
0342       // create the layer and store it
0343       MutableProtoLayer protoLayer(geoContext, layerSurfaces);
0344       protoLayer.envelope[AxisR] = {0.5, 0.5};
0345
0346       std::string layerName = "L" + std::to_string(ilp);
0347       detector.addMaterial(layerName + "_Mat", [&](auto& mat) {
0348         mat.configureFace(OuterCylinder, layerMaterialPtr);
0349         mat.addLayer(layerName, [&](auto& layer) {
0350           layer.setProtoLayer(protoLayer);
0351           layer.setLayerType(LayerType::Cylinder);
0352           layer.setEnvelope(
0353               ExtentEnvelope{{.z = {5_mm, 5_mm}, .r = {5_mm, 5_mm}}});
0354
0355           using SrfArrayNavPol = SurfaceArrayNavigationPolicy;
0356
0357           layer.setNavigationPolicyFactory(
0358               NavigationPolicyFactory{}
0359                   .add<TryAllNavigationPolicy>(
0360                       TryAllNavigationPolicy::Config{.sensitives = false})
0361                   .add<SrfArrayNavPol>(SrfArrayNavPol::Config{
0362                       .layerType = SrfArrayNavPol::LayerType::Cylinder,
0363                       .bins = kLayerBinning[ilp]})
0364                   .asUniquePtr());
0365         });
0366       });
0367
0368     }  // loop over layers
0369   });
0370
0371   BlueprintOptions opts;
0372   return root.construct(opts, geoContext, logger);
0373 }
0374
0375 std::shared_ptr<TrackingGeometry> CylindricalTrackingGeometry::operator()(
0376     const Logger& logger) {
0377   if (gen3) {
0378     return buildGen3(logger);
0379   } else {
0380     return buildGen1(logger);
0381   }
0382 }
0383 }  // namespace ActsTests