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(cos(modulePhi + modulePhiTilt),
0081                                    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(-sin(modulePhi + modulePhiTilt),
0086                                    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(cos(modulePhi), sin(modulePhi), 0.);
0237           // take different axis to have the same readout direction
0238           Acts::Vector3 moduleLocalZ(0., 0., side * 1.);
0239           Acts::Vector3 moduleLocalX = moduleLocalY.cross(moduleLocalZ);
0240           // local rotation matrices
0241           // create the RotationMatrix - negative side
0242           Acts::RotationMatrix3 moduleRotation;
0243           moduleRotation.col(0) = moduleLocalX;
0244           moduleRotation.col(1) = moduleLocalY;
0245           moduleRotation.col(2) = moduleLocalZ;
0246           // the transforms for the front module
0247           const Acts::Transform3 moduleTransform(
0248               Acts::Translation3(moduleCenter) * moduleRotation);
0249
0250           // create the module
0251           auto moduleElement =
0252               m_cfg.detectorElementFactory(moduleTransform, moduleBounds,
0253                                            moduleThickness, moduleMaterialPtr);
0254           // register the surface
0255           esVector.push_back(moduleElement->surface().getSharedPtr());
0256           // now deal with the potential backside
0257           if (!m_cfg.posnegModuleBacksideGap.empty()) {
0258             // the new centers
0259             moduleCenter =
0260                 moduleCenter +
0261                 m_cfg.posnegModuleBacksideGap.at(ipnl).at(ipnR) * moduleLocalZ;
0262             // the backside transforms
0263             Acts::Transform3 bsModuleTransform(
0264                 Acts::Translation3(moduleCenter) * moduleRotation);
0265             // apply the stereo
0266             if (!m_cfg.posnegModuleBacksideStereo.empty()) {
0267               // twist by the stereo angle
0268               double stereoBackSide =
0269                   m_cfg.posnegModuleBacksideStereo.at(ipnl).at(ipnR);
0270               bsModuleTransform *=
0271                   Acts::AngleAxis3(-stereoBackSide, Acts::Vector3::UnitZ());
0272             }
0273             // everything is set for the next module
0274             auto bsModuleElement = m_cfg.detectorElementFactory(
0275                 bsModuleTransform, moduleBounds, moduleThickness,
0276                 moduleMaterialPtr);
0277             // register the backside surface
0278             esVector.push_back(bsModuleElement->surface().getSharedPtr());
0279           }
0280         }
0281         // counter of rings
0282         ++ipnR;
0283       }
0284       // the binning
0285       std::size_t layerBinsR = m_cfg.posnegModulePhiBins.at(ipnl).size();
0286       // never multiply 1 single r-bin, does not make sense
0287       if (layerBinsR > 1) {
0288         // multiply with the given bin multiplier
0289         layerBinsR *= m_cfg.posnegLayerBinMultipliers.first;
0290       }
0291       std::size_t layerBinsPhi = 0;
0292       // take the minimum phi bins in that layer
0293       for (unsigned int phiBins : m_cfg.posnegModulePhiBins.at(ipnl)) {
0294         layerBinsPhi = phiBins < layerBinsPhi ? phiBins : layerBinsPhi;
0295         layerBinsPhi *= m_cfg.posnegLayerBinMultipliers.second;
0296       }
0297       // create the layers with the surface arrays
0298       Acts::MutableProtoLayer ple(gctx, esVector);
0299       ple.envelope[Acts::AxisDirection::AxisR] = {layerEnvelopeR,
0300                                                   layerEnvelopeR};
0301       ple.envelope[Acts::AxisDirection::AxisZ] = {
0302           m_cfg.approachSurfaceEnvelope, m_cfg.approachSurfaceEnvelope};
0303
0304       // push it into the layer vector
0305       ProtoLayerSurfaces ples{std::move(ple), esVector, layerBinsR,
0306                               layerBinsPhi};
0307       epLayers.push_back(std::move(ples));
0308     }
0309   }
0310   return epLayers;
0311 }
0312
0313 }  // namespace ActsExamples::Generic