GenericDetector/src/ProtoLayerCreator.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 "ActsExamples/GenericDetector/ProtoLayerCreator.hpp"
0010
0011 using Acts::VectorHelpers::phi;
0012
0013 namespace ActsExamples::Generic {
0014
0015 std::vector<ProtoLayerSurfaces> ProtoLayerCreator::centralProtoLayers(
0016     const Acts::GeometryContext& gctx) const {
0017   // create the vector
0018   std::vector<ProtoLayerSurfaces> cpLayers;
0019
0020   // ----------------------- central layers -------------------------
0021   // the central layers
0022   std::size_t numcLayers = m_cfg.centralLayerRadii.size();
0023   if (numcLayers != 0u) {
0024     ACTS_DEBUG("Configured to build " << numcLayers
0025                                       << " active central layers.");
0026     cpLayers.reserve(numcLayers);
0027     // loop through
0028     for (std::size_t icl = 0; icl < numcLayers; ++icl) {
0029       // layer R/Z
0030       double layerR = m_cfg.centralLayerRadii.at(icl);
0031       // some screen output
0032       ACTS_DEBUG("Build layer " << icl << " with target radius = " << layerR);
0033
0034       // prepare the Surface vector
0035       std::vector<std::shared_ptr<Acts::Surface>> sVector;
0036       // assign the current envelope
0037       double layerEnvelopeCoverZ =
0038           !m_cfg.centralLayerEnvelopes.empty()
0039               ? m_cfg.centralLayerEnvelopes.at(icl).second
0040               : 0.;
0041       // module size & tilt
0042       double modulePhiTilt = m_cfg.centralModuleTiltPhi.at(icl);
0043       double moduleHalfX = m_cfg.centralModuleHalfX.at(icl);
0044       double moduleHalfY = m_cfg.centralModuleHalfY.at(icl);
0045       double moduleThickness = m_cfg.centralModuleThickness.at(icl);
0046       // create the shared module
0047       auto moduleBounds =
0048           std::make_shared<Acts::RectangleBounds>(moduleHalfX, moduleHalfY);
0049       std::size_t nCentralModules =
0050           m_cfg.centralModuleBinningSchema.at(icl).first *
0051           m_cfg.centralModuleBinningSchema.at(icl).second;
0052
0053       ACTS_DEBUG("- number of modules "
0054                  << nCentralModules << " ( from "
0055                  << m_cfg.centralModuleBinningSchema.at(icl).first << " x "
0056                  << m_cfg.centralModuleBinningSchema.at(icl).second << " )");
0057
0058       sVector.reserve(nCentralModules);
0059
0060       // prepartation :
0061       // create the Module material from input
0062       std::shared_ptr<const Acts::ISurfaceMaterial> moduleMaterialPtr = nullptr;
0063       if (!m_cfg.centralModuleMaterial.empty()) {
0064         // get the sensor material from configuration
0065         moduleMaterialPtr = m_cfg.centralModuleMaterial.at(icl);
0066       }
0067
0068       // confirm
0069       if (m_cfg.centralModulePositions.at(icl).size() != nCentralModules) {
0070         ACTS_WARNING("Mismatching module numbers, configuration error!");
0071         ACTS_WARNING("- Binning schema suggests : " << nCentralModules);
0072         ACTS_WARNING("- Positions provided are  : "
0073                      << m_cfg.centralModulePositions.at(icl).size());
0074       }
0075       // loop over the position, create the modules
0076       for (auto& moduleCenter : m_cfg.centralModulePositions.at(icl)) {
0077         // create the association transform
0078         double modulePhi = phi(moduleCenter);
0079         // the local z axis is the normal vector
0080         Acts::Vector3 moduleLocalZ(std::cos(modulePhi + modulePhiTilt),
0081                                    std::sin(modulePhi + modulePhiTilt), 0.);
0082         // the local y axis is the global z axis
0083         Acts::Vector3 moduleLocalY(0., 0., 1);
0084         // the local x axis the normal to local y,z
0085         Acts::Vector3 moduleLocalX(-std::sin(modulePhi + modulePhiTilt),
0086                                    std::cos(modulePhi + modulePhiTilt), 0.);
0087         // create the RotationMatrix
0088         Acts::RotationMatrix3 moduleRotation;
0089         moduleRotation.col(0) = moduleLocalX;
0090         moduleRotation.col(1) = moduleLocalY;
0091         moduleRotation.col(2) = moduleLocalZ;
0092         // get the moduleTransform
0093         Acts::Transform3 moduleTransform(Acts::Translation3(moduleCenter) *
0094                                          moduleRotation);
0095         // stereo angle if necessary
0096         if (!m_cfg.centralModuleFrontsideStereo.empty() &&
0097             m_cfg.centralModuleFrontsideStereo.at(icl) != 0.) {
0098           // twist by the stereo angle
0099           double stereo = m_cfg.centralModuleFrontsideStereo.at(icl);
0100           moduleTransform *= Acts::AngleAxis3(-stereo, Acts::Vector3::UnitZ());
0101         }
0102         // create the module
0103         auto moduleElement = m_cfg.detectorElementFactory(
0104             moduleTransform, moduleBounds, moduleThickness, moduleMaterialPtr);
0105         // register the surface
0106         sVector.push_back(moduleElement->surface().getSharedPtr());
0107         // IF double modules exist
0108         // and the backside one (if configured to do so)
0109         if (!m_cfg.centralModuleBacksideGap.empty()) {
0110           // create the module identifier
0111
0112           Acts::Vector3 bsModuleCenter =
0113               moduleCenter +
0114               m_cfg.centralModuleBacksideGap.at(icl) * moduleLocalZ;
0115           Acts::Transform3 bsModuleTransform(
0116               Acts::Translation3(bsModuleCenter) * moduleRotation);
0117           // apply the stereo
0118           if (!m_cfg.centralModuleBacksideStereo.empty()) {
0119             // twist by the stereo angle
0120             double stereoBackSide = m_cfg.centralModuleBacksideStereo.at(icl);
0121             bsModuleTransform *=
0122                 Acts::AngleAxis3(-stereoBackSide, Acts::Vector3::UnitZ());
0123           }
0124           // create the backseide moulde
0125           auto bsModuleElement =
0126               m_cfg.detectorElementFactory(bsModuleTransform, moduleBounds,
0127                                            moduleThickness, moduleMaterialPtr);
0128           // register the backside surface
0129           sVector.push_back(bsModuleElement->surface().getSharedPtr());
0130         }
0131       }
0132
0133       std::size_t phiBins = m_cfg.centralModuleBinningSchema.at(icl).first;
0134       phiBins *= m_cfg.centralLayerBinMultipliers.first;
0135       std::size_t zBins = m_cfg.centralModuleBinningSchema.at(icl).second;
0136       zBins *= m_cfg.centralLayerBinMultipliers.second;
0137       // create the surface array - it will also fill the accessible binmember
0138       // cache if available
0139       Acts::MutableProtoLayer pl(gctx, sVector);
0140       pl.envelope[Acts::AxisDirection::AxisR] = {m_cfg.approachSurfaceEnvelope,
0141                                                  m_cfg.approachSurfaceEnvelope};
0142       pl.envelope[Acts::AxisDirection::AxisZ] = {layerEnvelopeCoverZ,
0143                                                  layerEnvelopeCoverZ};
0144
0145       // Record the proto layer and the surfaces for the later layer building
0146       ProtoLayerSurfaces pls{std::move(pl), sVector, phiBins, zBins};
0147       cpLayers.push_back(std::move(pls));
0148     }
0149   }
0150   return cpLayers;
0151 }
0152
0153 std::vector<ProtoLayerSurfaces> ProtoLayerCreator::negativeProtoLayers(
0154     const Acts::GeometryContext& gctx) const {
0155   return createProtoLayers(gctx, -1);
0156 }
0157
0158 std::vector<ProtoLayerSurfaces> ProtoLayerCreator::positiveProtoLayers(
0159     const Acts::GeometryContext& gctx) const {
0160   return createProtoLayers(gctx, 1);
0161 }
0162
0163 ProtoLayerCreator::ProtoLayerCreator(const ProtoLayerCreator::Config& cfg,
0164                                      std::unique_ptr<const Acts::Logger> log)
0165     : m_cfg(cfg), m_logger(std::move(log)) {
0166   if (!m_cfg.detectorElementFactory) {
0167     throw std::invalid_argument("Detector element factory is not set");
0168   }
0169 }
0170
0171 std::vector<ProtoLayerSurfaces> ProtoLayerCreator::createProtoLayers(
0172     const Acts::GeometryContext& gctx, int side) const {
0173   // Count the current detector modules identifiers
0174   // the return layers
0175   std::vector<ProtoLayerSurfaces> epLayers;
0176   // -------------------------------- endcap type layers
0177   // pos/neg layers
0178   std::size_t numpnLayers = m_cfg.posnegLayerPositionsZ.size();
0179   if (numpnLayers != 0u) {
0180     ACTS_DEBUG("Configured to build 2 * "
0181                << numpnLayers << " passive positive/negative side layers.");
0182     epLayers.reserve(numpnLayers);
0183
0184     /// this is the loop over the layer positions
0185     for (std::size_t ipnl = 0; ipnl < numpnLayers; ++ipnl) {
0186       // some screen output
0187       ACTS_VERBOSE("- building layer "
0188                    << ipnl << " and " << numpnLayers + ipnl << " at z = "
0189                    << side * m_cfg.posnegLayerPositionsZ.at(ipnl));
0190       /// some preparation work
0191       // define the layer envelope
0192       double layerEnvelopeR = m_cfg.posnegLayerEnvelopeR.at(ipnl);
0193       // prepare for the r binning
0194       std::vector<std::shared_ptr<Acts::Surface>> esVector;
0195       // now fill the vectors
0196       std::size_t ipnR = 0;
0197       for (auto& discModulePositions : m_cfg.posnegModulePositions.at(ipnl)) {
0198         ACTS_VERBOSE("- building ring " << ipnR << " for this layer.");
0199         // now prepare all the shared stuff
0200         // (0) module specifications
0201         double moduleThickness = m_cfg.posnegModuleThickness.at(ipnl).at(ipnR);
0202         double moduleMinHalfX = m_cfg.posnegModuleMinHalfX.at(ipnl).at(ipnR);
0203         double moduleMaxHalfX = 0.;
0204         if (m_cfg.posnegModuleMaxHalfX.size() > ipnl &&
0205             m_cfg.posnegModuleMaxHalfX.at(ipnl).size() > ipnR) {
0206           moduleMaxHalfX = m_cfg.posnegModuleMaxHalfX.at(ipnl).at(ipnR);
0207         }
0208         double moduleHalfY = m_cfg.posnegModuleHalfY.at(ipnl).at(ipnR);
0209         // (1) module bounds
0210         // create the bounds
0211         std::shared_ptr<const Acts::PlanarBounds> moduleBounds;
0212         if (moduleMaxHalfX != 0. && moduleMinHalfX != moduleMaxHalfX) {
0213           moduleBounds = std::make_shared<Acts::TrapezoidBounds>(
0214               moduleMinHalfX, moduleMaxHalfX, moduleHalfY);
0215         } else {
0216           moduleBounds = std::make_shared<Acts::RectangleBounds>(moduleMinHalfX,
0217                                                                  moduleHalfY);
0218         }
0219         // (2)) module material
0220         // create the Module material from input
0221         std::shared_ptr<const Acts::ISurfaceMaterial> moduleMaterialPtr =
0222             nullptr;
0223         if (!m_cfg.posnegModuleMaterial.empty()) {
0224           // and create the shared pointer
0225           moduleMaterialPtr = m_cfg.posnegModuleMaterial.at(ipnl).at(ipnR);
0226         }
0227
0228         // low loop over the phi positions and build the stuff
0229         for (auto& ringModulePosition : discModulePositions) {
0230           // the module transform from the position
0231           double modulePhi = phi(ringModulePosition);
0232           // the center position of the modules
0233           Acts::Vector3 moduleCenter(ringModulePosition);
0234           moduleCenter.z() *= side;
0235           // the rotation matrix of the module
0236           Acts::Vector3 moduleLocalY(std::cos(modulePhi), std::sin(modulePhi),
0237                                      0.);
0238           // take different axis to have the same readout direction
0239           Acts::Vector3 moduleLocalZ(0., 0., side * 1.);
0240           Acts::Vector3 moduleLocalX = moduleLocalY.cross(moduleLocalZ);
0241           // local rotation matrices
0242           // create the RotationMatrix - negative side
0243           Acts::RotationMatrix3 moduleRotation;
0244           moduleRotation.col(0) = moduleLocalX;
0245           moduleRotation.col(1) = moduleLocalY;
0246           moduleRotation.col(2) = moduleLocalZ;
0247           // the transforms for the front module
0248           const Acts::Transform3 moduleTransform(
0249               Acts::Translation3(moduleCenter) * moduleRotation);
0250
0251           // create the module
0252           auto moduleElement =
0253               m_cfg.detectorElementFactory(moduleTransform, moduleBounds,
0254                                            moduleThickness, moduleMaterialPtr);
0255           // register the surface
0256           esVector.push_back(moduleElement->surface().getSharedPtr());
0257           // now deal with the potential backside
0258           if (!m_cfg.posnegModuleBacksideGap.empty()) {
0259             // the new centers
0260             moduleCenter =
0261                 moduleCenter +
0262                 m_cfg.posnegModuleBacksideGap.at(ipnl).at(ipnR) * moduleLocalZ;
0263             // the backside transforms
0264             Acts::Transform3 bsModuleTransform(
0265                 Acts::Translation3(moduleCenter) * moduleRotation);
0266             // apply the stereo
0267             if (!m_cfg.posnegModuleBacksideStereo.empty()) {
0268               // twist by the stereo angle
0269               double stereoBackSide =
0270                   m_cfg.posnegModuleBacksideStereo.at(ipnl).at(ipnR);
0271               bsModuleTransform *=
0272                   Acts::AngleAxis3(-stereoBackSide, Acts::Vector3::UnitZ());
0273             }
0274             // everything is set for the next module
0275             auto bsModuleElement = m_cfg.detectorElementFactory(
0276                 bsModuleTransform, moduleBounds, moduleThickness,
0277                 moduleMaterialPtr);
0278             // register the backside surface
0279             esVector.push_back(bsModuleElement->surface().getSharedPtr());
0280           }
0281         }
0282         // counter of rings
0283         ++ipnR;
0284       }
0285       // the binning
0286       std::size_t layerBinsR = m_cfg.posnegModulePhiBins.at(ipnl).size();
0287       // never multiply 1 single r-bin, does not make sense
0288       if (layerBinsR > 1) {
0289         // multiply with the given bin multiplier
0290         layerBinsR *= m_cfg.posnegLayerBinMultipliers.first;
0291       }
0292       std::size_t layerBinsPhi = 0;
0293       // take the minimum phi bins in that layer
0294       for (unsigned int phiBins : m_cfg.posnegModulePhiBins.at(ipnl)) {
0295         layerBinsPhi = phiBins < layerBinsPhi ? phiBins : layerBinsPhi;
0296         layerBinsPhi *= m_cfg.posnegLayerBinMultipliers.second;
0297       }
0298       // create the layers with the surface arrays
0299       Acts::MutableProtoLayer ple(gctx, esVector);
0300       ple.envelope[Acts::AxisDirection::AxisR] = {layerEnvelopeR,
0301                                                   layerEnvelopeR};
0302       ple.envelope[Acts::AxisDirection::AxisZ] = {
0303           m_cfg.approachSurfaceEnvelope, m_cfg.approachSurfaceEnvelope};
0304
0305       // push it into the layer vector
0306       ProtoLayerSurfaces ples{std::move(ple), esVector, layerBinsR,
0307                               layerBinsPhi};
0308       epLayers.push_back(std::move(ples));
0309     }
0310   }
0311   return epLayers;
0312 }
0313
0314 }  // namespace ActsExamples::Generic