diff --git a/InnerDetector/InDetDetDescr/PixelGeoModelXml/src/PixelGmxInterface.cxx b/InnerDetector/InDetDetDescr/PixelGeoModelXml/src/PixelGmxInterface.cxx index 66fd2d83aea81453c326b234341428d79ba8126d..a70922526c3dc504052c54f63c19b2ef0a27293c 100644 --- a/InnerDetector/InDetDetDescr/PixelGeoModelXml/src/PixelGmxInterface.cxx +++ b/InnerDetector/InDetDetDescr/PixelGeoModelXml/src/PixelGmxInterface.cxx @@ -212,7 +212,13 @@ void PixelGmxInterface::addSensor(std::string typeName, // // Create the detector element and add to the DetectorManager // - const InDetDD::SiDetectorDesign *design = m_detectorManager->getDesign(m_geometryMap[typeName]); + auto it = m_geometryMap.find(typeName); + if(it == m_geometryMap.end()) { + ATH_MSG_ERROR("addSensor: Error: Readout sensor type " << typeName << " not found."); + throw std::runtime_error("readout sensor type " + typeName + " not found."); + } + const InDetDD::SiDetectorDesign *design = m_detectorManager->getDesign(it->second); + ATH_MSG_VERBOSE("Adding sensor with design: " << typeName << " " << design); if (design == nullptr) { ATH_MSG_ERROR("addSensor: Error: Readout sensor type " << typeName << " not found."); throw std::runtime_error("readout sensor type " + typeName + " not found."); diff --git a/InnerDetector/InDetDetDescr/SCT_ReadoutGeometry/SCT_ReadoutGeometry/StripStereoAnnulusDesign.h b/InnerDetector/InDetDetDescr/SCT_ReadoutGeometry/SCT_ReadoutGeometry/StripStereoAnnulusDesign.h index 708a3b79154aec5d3bbca85f9c8191861b06b1f0..74fc5aa57087181c57a63fcaf5b8f4bb5f7a1617 100644 --- a/InnerDetector/InDetDetDescr/SCT_ReadoutGeometry/SCT_ReadoutGeometry/StripStereoAnnulusDesign.h +++ b/InnerDetector/InDetDetDescr/SCT_ReadoutGeometry/SCT_ReadoutGeometry/StripStereoAnnulusDesign.h @@ -204,6 +204,11 @@ StripStereoAnnulusDesign(const SiDetectorDesign::Axis &stripDirection, double pitch(const SiCellId &cellId) const; double stripLength(const SiCellId &cellId) const; + double minR() const; + double maxR() const; + double phiWidth() const; + double stereo() const; + // Give upper and lower boundaries, and length, of dead area double deadAreaUpperBoundary() const; double deadAreaLowerBoundary() const; @@ -301,6 +306,24 @@ inline SiReadoutCellId StripStereoAnnulusDesign::readoutIdOfCell(const SiCellId return SiReadoutCellId(strip, row); } + +inline double StripStereoAnnulusDesign::minR() const { + return m_stripStartRadius[0]; +} + +inline double StripStereoAnnulusDesign::maxR() const { + return m_stripEndRadius.back(); +} + + +inline double StripStereoAnnulusDesign::phiWidth() const { + return m_nStrips[0] * m_pitch[0]; +} + +inline double StripStereoAnnulusDesign::stereo() const { + return m_stereo; +} + inline int StripStereoAnnulusDesign::row(int stripId1Dim) const { //This method is scheduled for deletion diff --git a/InnerDetector/InDetDetDescr/StripGeoModelXml/src/StripGmxInterface.cxx b/InnerDetector/InDetDetDescr/StripGeoModelXml/src/StripGmxInterface.cxx index adcb499e4f010e27390ae5b6acbef034f7ae2c46..36cf8541f08339abf7cb86c6b551aa63238b9a06 100644 --- a/InnerDetector/InDetDetDescr/StripGeoModelXml/src/StripGmxInterface.cxx +++ b/InnerDetector/InDetDetDescr/StripGeoModelXml/src/StripGmxInterface.cxx @@ -479,11 +479,12 @@ void StripGmxInterface::addSplitSensor(std::string typeName, splitTypeName += "_" + std::to_string(updatedIndex["side"]); } - const InDetDD::SiDetectorDesign *design = m_geometryMap[splitTypeName]; - if (design == nullptr) { + auto it = m_geometryMap.find(splitTypeName); + if(it == m_geometryMap.end() || it->second == nullptr) { ATH_MSG_ERROR("addSplitSensor: Error: Readout sensor type " << typeName << " not found."); throw std::runtime_error("readout sensor type " + typeName + " not found."); } + const InDetDD::SiDetectorDesign *design = it->second; m_detectorManager->addDetectorElement(new InDetDD::SiDetectorElement(id, design, fpv, m_commonItems)); @@ -540,11 +541,12 @@ void StripGmxInterface::addSensor(std::string typeName, // // Create the detector element and add to the DetectorManager // - const InDetDD::SiDetectorDesign *design = m_geometryMap[typeName]; - if (design == nullptr) { + auto it = m_geometryMap.find(typeName); + if(it == m_geometryMap.end() || it->second == nullptr) { ATH_MSG_ERROR("addSensor: Error: Readout sensor type " << typeName << " not found."); throw std::runtime_error("readout sensor type " + typeName + " not found."); } + const InDetDD::SiDetectorDesign *design = it->second; m_detectorManager->addDetectorElement(new InDetDD::SiDetectorElement(id, design, fpv, m_commonItems)); diff --git a/Tracking/Acts/ActsGeometry/ActsGeometry/ActsDetectorElement.h b/Tracking/Acts/ActsGeometry/ActsGeometry/ActsDetectorElement.h index 7984804aa18380d8754e4e68e75baee77cdef7e6..0c83470b123c1077df6ee74d48c3dd4221d1748f 100644 --- a/Tracking/Acts/ActsGeometry/ActsGeometry/ActsDetectorElement.h +++ b/Tracking/Acts/ActsGeometry/ActsGeometry/ActsDetectorElement.h @@ -5,30 +5,28 @@ #ifndef ACTSGEOMETRY_ACTSDETECTORELEMENT_H #define ACTSGEOMETRY_ACTSDETECTORELEMENT_H -//Amg Eigen plugin includes +// Amg Eigen plugin includes #include "EventPrimitives/EventPrimitives.h" #include "GeoPrimitives/GeoPrimitives.h" // ATHENA INCLUDES -#include "TRT_ReadoutGeometry/TRT_BaseElement.h" -#include "InDetReadoutGeometry/SiDetectorElement.h" #include "InDetIdentifier/PixelID.h" #include "InDetIdentifier/SCT_ID.h" +#include "InDetReadoutGeometry/SiDetectorElement.h" +#include "TRT_ReadoutGeometry/TRT_BaseElement.h" +#include "TrkSurfaces/AnnulusBounds.h" // ACTS #include "Acts/Geometry/DetectorElementBase.hpp" #include "Acts/Geometry/GeometryContext.hpp" // STL -#include <mutex> #include <iostream> - -// BOOST -#include <boost/variant.hpp> +#include <mutex> namespace Acts { - class SurfaceBounds; +class SurfaceBounds; } class ActsTrackingGeometrySvc; @@ -38,96 +36,53 @@ class IdentityHelper; /// @class ActsDetectorElement /// -class ActsDetectorElement : public Acts::DetectorElementBase -{ - using DetElemVariant = boost::variant<const InDetDD::SiDetectorElement*, const InDetDD::TRT_BaseElement*>; +class ActsDetectorElement : public Acts::DetectorElementBase { public: - enum class Subdetector { Pixel, SCT, TRT, Size }; - - ActsDetectorElement(const InDetDD::SiDetectorElement* detElem); + ActsDetectorElement(const InDetDD::SiDetectorElement &detElem); /// Constructor for a straw surface. /// @param transform Transform to the straw system - ActsDetectorElement(const Acts::Transform3& trf, - const InDetDD::TRT_BaseElement* detElem, - const Identifier& id // we need explicit ID here b/c of straws - ); + ActsDetectorElement( + const Acts::Transform3 &trf, const InDetDD::TRT_BaseElement &detElem, + const Identifier &id // we need explicit ID here b/c of straws + ); /// Destructor virtual ~ActsDetectorElement() {} /// Identifier - Identifier - identify() const; + Identifier identify() const; /// Return local to global transform associated with this identifier - - void - storeTransform(ActsAlignmentStore* gas) const; + + void storeTransform(ActsAlignmentStore *gas) const; virtual const Acts::Transform3 & transform(const Acts::GeometryContext &gctx) const final override; - /// Return surface associated with this identifier, which should come from the - virtual const Acts::Surface& - surface() const final override; + virtual const Acts::Surface &surface() const final override; - /// Return a shared pointer on the ATLAS surface associated with this identifier, - const Trk::Surface& - atlasSurface() const; + /// Return a shared pointer on the ATLAS surface associated with this + /// identifier, + const Trk::Surface &atlasSurface() const; /// Returns the thickness of the module - virtual double - thickness() const final override; - - // optimize this, it doesn't ever change so can be saved on construction - Subdetector getSubdetector() const { - if(m_detElement.which() == 1) { - return Subdetector::TRT; - } - - // this should be Pixel or SCT - const auto& de = *boost::get<const InDetDD::SiDetectorElement*>(m_detElement); - if(de.isPixel()) { - return Subdetector::Pixel; - } - else { - return Subdetector::SCT; - } - } + virtual double thickness() const final override; IdentityHelper identityHelper() const; - -private: - /// Returns default transform. For TRT this is static and set in constructor. - /// For silicon detectors it is calulated from GM, and stored. Thus the method + /// For silicon detectors it is calulated from GM, and stored. Thus the method /// is not const. The store is mutexed. - const Acts::Transform3& - getDefaultTransformMutexed() const; - - struct IdVisitor : public boost::static_visitor<Identifier> - { - explicit IdVisitor(const Identifier& id) : m_explicitIdentifier(id) {} - - Identifier operator()(const InDetDD::SiDetectorElement* detElem) const - { - return detElem->identify(); // easy, det element has identifier - } - - Identifier operator()(const InDetDD::TRT_BaseElement*) const - { - // we got the identifier in constructrion, because it identifies - // the STRAW - return m_explicitIdentifier; - } + const Acts::Transform3 &getDefaultTransformMutexed() const; - Identifier m_explicitIdentifier; - }; + /// Returns the underllying GeoModel detectorelement that this one + /// is based on. + const GeoVDetectorElement *upstreamDetectorElement() const; +private: /// Detector element as variant - DetElemVariant m_detElement; + const GeoVDetectorElement *m_detElement; /// Boundaries of the detector element std::shared_ptr<const Acts::SurfaceBounds> m_bounds; /// Thickness of this detector element @@ -139,14 +94,11 @@ private: // this is pretty much only used single threaded, so // the mutex does not hurt mutable std::mutex m_cacheMutex; - mutable std::shared_ptr<const Acts::Transform3> m_defTransform; - - Identifier m_explicitIdentifier; - - // this is threadsafe! - //mutable Gaudi::Hive::ContextSpecificData<Acts::Transform3> m_ctxSpecificTransform; + mutable std::optional<Acts::Transform3> m_defTransform; + Acts::Transform3 m_extraTransform{Acts::Transform3::Identity()}; + Identifier m_explicitIdentifier; }; #endif diff --git a/Tracking/Acts/ActsGeometry/ActsGeometry/ActsExtrapolationTool.h b/Tracking/Acts/ActsGeometry/ActsGeometry/ActsExtrapolationTool.h index 1f3d35a14c51a5caa7fb4860f8a843e75041c060..05233d031e2bb3f7a2510cceba5b6e379b9fa0e6 100644 --- a/Tracking/Acts/ActsGeometry/ActsGeometry/ActsExtrapolationTool.h +++ b/Tracking/Acts/ActsGeometry/ActsGeometry/ActsExtrapolationTool.h @@ -6,6 +6,7 @@ #define ACTSGEOMETRY_ACTSEXTRAPOLATIONTOOL_H // ATHENA +#include "GeoPrimitives/GeoPrimitives.h" #include "AthenaBaseComps/AthAlgTool.h" #include "GaudiKernel/IInterface.h" #include "GaudiKernel/ServiceHandle.h" diff --git a/Tracking/Acts/ActsGeometry/ActsGeometry/ActsLayerBuilder.h b/Tracking/Acts/ActsGeometry/ActsGeometry/ActsLayerBuilder.h index 9e1de50e238a7e10337317877cb5a63bbecb79c8..ab14d3d2babbcf663c310d86ca239731dad3a7ad 100644 --- a/Tracking/Acts/ActsGeometry/ActsGeometry/ActsLayerBuilder.h +++ b/Tracking/Acts/ActsGeometry/ActsGeometry/ActsLayerBuilder.h @@ -31,6 +31,10 @@ class LayerCreator; class ActsLayerBuilder : public Acts::ILayerBuilder { public: + enum class Mode { + Undefined, Pixel, SCT, TRT, ITkPixelInner, ITkPixelOuter, ITkStrip, + }; + using ElementVector = std::vector<std::shared_ptr<const ActsDetectorElement>>; @@ -41,8 +45,7 @@ public: { /// string based identification std::string configurationName = "undefined"; - ActsDetectorElement::Subdetector subdetector - = ActsDetectorElement::Subdetector::Pixel; + Mode mode = Mode::Undefined; const InDetDD::SiDetectorManager* mng; std::shared_ptr<const Acts::LayerCreator> layerCreator = nullptr; /// the binning type of the contained surfaces in phi @@ -58,6 +61,8 @@ public: std::pair<size_t, size_t> endcapMaterialBins = {20, 5}; std::pair<size_t, size_t> barrelMaterialBins = {10, 10}; + + bool objDebugOutput = false; }; /// Constructor @@ -65,12 +70,7 @@ public: /// @param logger the local logging instance ActsLayerBuilder(const Config& cfg, std::unique_ptr<const Acts::Logger> logger - = Acts::getDefaultLogger("GMLayBldr", Acts::Logging::INFO)) - : m_logger(std::move(logger)) - { - // std::cout << "GMLB construct" << std::endl; - m_cfg = cfg; - } + = Acts::getDefaultLogger("GMLayBldr", Acts::Logging::INFO)); /// Destructor ~ActsLayerBuilder() {} @@ -84,6 +84,7 @@ public: const Acts::LayerVector positiveLayers(const Acts::GeometryContext& gctx) const override; + /// Name identification // const std::string& // identification() const final; @@ -139,7 +140,13 @@ private: // @param layers is goint to be filled // @param type is the indication which ones to build -1 | 0 | 1 void - buildLayers(const Acts::GeometryContext& gctx, Acts::LayerVector& layersOutput, int type = 0); + buildBarrel(const Acts::GeometryContext& gctx, Acts::LayerVector& layersOutput); + + void + buildEndcap(const Acts::GeometryContext& gctx, Acts::LayerVector& layersOutput, int type = 0); + }; +std::ostream& operator<<(std::ostream& os, const ActsLayerBuilder::Mode& mode); + #endif diff --git a/Tracking/Acts/ActsGeometry/ActsGeometry/ActsTrackingGeometrySvc.h b/Tracking/Acts/ActsGeometry/ActsGeometry/ActsTrackingGeometrySvc.h index b9c5d4b4ef55606e337ec1b542edcb1ffafad5fa..f7d6eacf2573b1012a2e526feef1947fe42b677a 100644 --- a/Tracking/Acts/ActsGeometry/ActsGeometry/ActsTrackingGeometrySvc.h +++ b/Tracking/Acts/ActsGeometry/ActsGeometry/ActsTrackingGeometrySvc.h @@ -13,6 +13,7 @@ // PACKAGE #include "ActsGeometryInterfaces/IActsTrackingGeometrySvc.h" #include "ActsGeometryInterfaces/IActsTrackingVolumeBuilder.h" +#include "ActsGeometry/ActsLayerBuilder.h" // STL #include <map> @@ -60,8 +61,12 @@ public: getNominalAlignmentStore() const override; private: + ActsLayerBuilder::Config + makeLayerBuilderConfig(const InDetDD::InDetDetectorManager* manager); + std::shared_ptr<const Acts::ILayerBuilder> - makeLayerBuilder(const InDetDD::InDetDetectorManager* manager); + makeStrawLayerBuilder(const InDetDD::InDetDetectorManager* manager); + std::shared_ptr<Acts::TrackingVolume> makeSCTTRTAssembly(const Acts::GeometryContext& gctx, const Acts::ILayerBuilder& sct_lb, @@ -72,6 +77,8 @@ private: const InDetDD::SiDetectorManager* p_pixelManager; const InDetDD::SiDetectorManager* p_SCTManager; const InDetDD::TRT_DetectorManager* p_TRTManager; + const InDetDD::SiDetectorManager* p_ITkPixelManager; + const InDetDD::SiDetectorManager* p_ITkStripManager; std::shared_ptr<std::vector<std::shared_ptr<const ActsDetectorElement>>> m_elementStore; std::shared_ptr<const Acts::TrackingGeometry> m_trackingGeometry; @@ -81,6 +88,7 @@ private: std::unique_ptr<const ActsAlignmentStore> m_nominalAlignmentStore{nullptr}; Gaudi::Property<bool> m_useMaterialMap{this, "UseMaterialMap", false, ""}; + Gaudi::Property<bool> m_objDebugOutput{this, "ObjDebugOutput", false, ""}; Gaudi::Property<std::string> m_materialMapInputFile{this, "MaterialMapInputFile", "", ""}; Gaudi::Property<std::vector<size_t>> m_barrelMaterialBins{this, "BarrelMaterialBins", {10, 10}}; Gaudi::Property<std::vector<size_t>> m_endcapMaterialBins{this, "EndcapMaterialBins", {5, 20}}; diff --git a/Tracking/Acts/ActsGeometry/CMakeLists.txt b/Tracking/Acts/ActsGeometry/CMakeLists.txt index 1b3a337f357a914532b41e1a4506ec184fc82211..d0963d06a818bb574239a165da6b5a77defde39f 100755 --- a/Tracking/Acts/ActsGeometry/CMakeLists.txt +++ b/Tracking/Acts/ActsGeometry/CMakeLists.txt @@ -41,6 +41,7 @@ atlas_add_library( ActsGeometryLib PRIVATE_LINK_LIBRARIES StoreGateLib TRT_ReadoutGeometry + SCT_ReadoutGeometry TrkGeometry TrkSurfaces TrkFitterUtils diff --git a/Tracking/Acts/ActsGeometry/python/ActsGeometryConfig.py b/Tracking/Acts/ActsGeometry/python/ActsGeometryConfig.py index 447da6ed4a5679243f51b51ef54109b370b76050..655bc88e0c79a550f3602979d2a81f8d20996987 100644 --- a/Tracking/Acts/ActsGeometry/python/ActsGeometryConfig.py +++ b/Tracking/Acts/ActsGeometry/python/ActsGeometryConfig.py @@ -3,7 +3,7 @@ from AthenaConfiguration.ComponentAccumulator import ComponentAccumulator from AthenaConfiguration.ComponentFactory import CompFactory -def ActsTrackingGeometrySvcCfg(configFlags, name = "ActsTrackingGeometrySvc" ) : +def ActsTrackingGeometrySvcCfg(configFlags, name = "ActsTrackingGeometrySvc", **kwargs) : result = ComponentAccumulator() Acts_ActsTrackingGeometrySvc = CompFactory.ActsTrackingGeometrySvc @@ -16,13 +16,18 @@ def ActsTrackingGeometrySvcCfg(configFlags, name = "ActsTrackingGeometrySvc" ) : subDetectors += ["TRT"] if configFlags.Detector.GeometryCalo: subDetectors += ["Calo"] - # need to configure calo geometry, otherwise we get a crash from LArGeoAlgsNV.LArGMConfig import LArGMCfg result.merge(LArGMCfg(configFlags)) from TileGeoModel.TileGMConfig import TileGMCfg result.merge(TileGMCfg(configFlags)) + if configFlags.Detector.GeometryITkPixel: + subDetectors += ["ITkPixel"] + if configFlags.Detector.GeometryITkStrip: + subDetectors += ["ITkStrip"] + + idSub = [sd in subDetectors for sd in ("Pixel", "SCT", "TRT")] if any(idSub): # ANY of the ID subdetectors are on => we require GM sources @@ -35,10 +40,7 @@ def ActsTrackingGeometrySvcCfg(configFlags, name = "ActsTrackingGeometrySvc" ) : from AthenaCommon.Logging import log log.warning("ConfigFlags indicate %s should be built. Not all ID subdetectors are set, but I'll set all of them up to capture the extra setup happening here.", ", ".join(subDetectors)) - actsTrackingGeometrySvc = Acts_ActsTrackingGeometrySvc(name, BuildSubDetectors=subDetectors) - - - actsTrackingGeometrySvc = Acts_ActsTrackingGeometrySvc(name, BuildSubDetectors=subDetectors) + actsTrackingGeometrySvc = Acts_ActsTrackingGeometrySvc(name, BuildSubDetectors=subDetectors, **kwargs) if configFlags.TrackingGeometry.MaterialSource == "Input": actsTrackingGeometrySvc.UseMaterialMap = True @@ -200,3 +202,21 @@ def ActsObjWriterToolCfg(name= "ActsObjWriterTool", **kwargs) : result.addPublicTool(ActsObjWriterTool, primary=True) return result + + +def ActsExtrapolationAlgCfg(configFlags, name = "ActsExtrapolationAlg", **kwargs): + result = ComponentAccumulator() + + if "ExtrapolationTool" not in kwargs: + extrapTool = ActsExtrapolationToolCfg(configFlags) + kwargs["ExtrapolationTool"] = extrapTool.getPrimary() + result.merge(extrapTool) + + propStepWriterSvc = ActsPropStepRootWriterSvcCfg(configFlags) + result.merge(propStepWriterSvc) + + ActsExtrapolationAlg = CompFactory.ActsExtrapolationAlg + alg = ActsExtrapolationAlg(name, **kwargs) + result.addEventAlgo(alg) + + return result \ No newline at end of file diff --git a/Tracking/Acts/ActsGeometry/share/ActsExtrapolationAlg.py b/Tracking/Acts/ActsGeometry/share/ActsExtrapolationAlg.py index 5df264b388bacd22c114319f3c928c0a916cb6a2..a54938751eb511a2d61a10e1d3664c1fde947776 100644 --- a/Tracking/Acts/ActsGeometry/share/ActsExtrapolationAlg.py +++ b/Tracking/Acts/ActsGeometry/share/ActsExtrapolationAlg.py @@ -10,25 +10,7 @@ Copyright (C) 2002-2021 CERN for the benefit of the ATLAS collaboration from AthenaConfiguration.ComponentAccumulator import ComponentAccumulator from AthenaConfiguration.ComponentFactory import CompFactory -from ActsGeometry.ActsGeometryConfig import ActsExtrapolationToolCfg -from ActsGeometry.ActsGeometryConfig import ActsAlignmentCondAlgCfg, ActsPropStepRootWriterSvcCfg - -def ActsExtrapolationAlgCfg(configFlags, name = "ActsExtrapolationAlg", **kwargs): - result = ComponentAccumulator() - - if "ExtrapolationTool" not in kwargs: - extrapTool = ActsExtrapolationToolCfg(configFlags) - kwargs["ExtrapolationTool"] = extrapTool.getPrimary() - result.merge(extrapTool) - - propStepWriterSvc = ActsPropStepRootWriterSvcCfg(configFlags) - result.merge(propStepWriterSvc) - - ActsExtrapolationAlg = CompFactory.ActsExtrapolationAlg - alg = ActsExtrapolationAlg(name, **kwargs) - result.addEventAlgo(alg) - - return result +from ActsGeometry.ActsGeometryConfig import ActsAlignmentCondAlgCfg, ActsPropStepRootWriterSvcCfg, ActsExtrapolationAlgCfg, ActsExtrapolationToolCfg if "__main__" == __name__: from AthenaCommon.Configurable import Configurable diff --git a/Tracking/Acts/ActsGeometry/share/ITkTest.py b/Tracking/Acts/ActsGeometry/share/ITkTest.py new file mode 100755 index 0000000000000000000000000000000000000000..3ae2a077ac4e058640fa388dbfa0ea217270bcfb --- /dev/null +++ b/Tracking/Acts/ActsGeometry/share/ITkTest.py @@ -0,0 +1,74 @@ +#!/usr/bin/env python +"""Run ACTS geometry construction for ITk + +Copyright (C) 2002-2021 CERN for the benefit of the ATLAS collaboration +""" +from AthenaCommon.Configurable import Configurable +from AthenaConfiguration.AllConfigFlags import ConfigFlags + +# Set up logging and new style config +Configurable.configurableRun3Behavior = True + +ConfigFlags.GeoModel.useLocalGeometry = True + +detectors = [ + "ITkPixel", + "ITkStrip" +] + +from AthenaConfiguration.DetectorConfigFlags import setupDetectorsFromList +setupDetectorsFromList(ConfigFlags, detectors, toggle_geometry=True) + +from AthenaConfiguration.Enums import ProductionStep +ConfigFlags.Common.ProductionStep = ProductionStep.Simulation +ConfigFlags.GeoModel.AtlasVersion = "ATLAS-P2-ITK-24-00-00" +ConfigFlags.IOVDb.GlobalTag = "OFLCOND-MC16-SDR-14" +ConfigFlags.GeoModel.Align.Dynamic = False +ConfigFlags.Input.Files = ['/cvmfs/atlas-nightlies.cern.ch/repo/data/data-art/SimCoreTests/valid1.410000.PowhegPythiaEvtGen_P2012_ttbar_hdamp172p5_nonallhad.evgen.EVNT.e4993.EVNT.08166201._000012.pool.root.1'] + +ConfigFlags.Concurrency.NumThreads = 64 +ConfigFlags.Concurrency.NumConcurrentEvents = 64 + +ConfigFlags.lock() + +from AthenaConfiguration.MainServicesConfig import MainServicesCfg +acc = MainServicesCfg(ConfigFlags) +from AthenaConfiguration.ComponentFactory import CompFactory + + +from PixelGeoModelXml.ITkPixelGeoModelConfig import ITkPixelGeometryCfg +itkPixel = ITkPixelGeometryCfg(ConfigFlags) +acc.merge(itkPixel) + +from StripGeoModelXml.ITkStripGeoModelConfig import ITkStripGeometryCfg +itkStrip = ITkStripGeometryCfg(ConfigFlags) +acc.merge(itkStrip) + +from BeamPipeGeoModel.BeamPipeGMConfig import BeamPipeGeometryCfg +acc.merge(BeamPipeGeometryCfg(ConfigFlags)) + +from ActsGeometry.ActsGeometryConfig import ActsAlignmentCondAlgCfg, ActsPropStepRootWriterSvcCfg, ActsExtrapolationAlgCfg, ActsExtrapolationToolCfg, ActsTrackingGeometrySvcCfg + +from ActsGeometry.ActsGeometryConfig import NominalAlignmentCondAlgCfg +acc.merge(NominalAlignmentCondAlgCfg(ConfigFlags)) + +from AthenaCommon.Constants import INFO, VERBOSE +tgSvc = ActsTrackingGeometrySvcCfg(ConfigFlags, + OutputLevel=VERBOSE, + ObjDebugOutput=True) +acc.merge(tgSvc) + +alg = ActsExtrapolationAlgCfg(ConfigFlags, + OutputLevel=INFO, + NParticlesPerEvent = int(100), + WritePropStep = True, + EtaRange = [-5, 5], + PtRange = [20, 100]) + +acc.merge(alg) + +acc.printConfig() + + + +acc.run(200) diff --git a/Tracking/Acts/ActsGeometry/src/ActsATLASConverterTool.cxx b/Tracking/Acts/ActsGeometry/src/ActsATLASConverterTool.cxx index 43dd14611d7eba5d72a64fc0bbe652130f9f9c29..34451e009c5a70d00a00b82abd21ec4a329b6e7d 100644 --- a/Tracking/Acts/ActsGeometry/src/ActsATLASConverterTool.cxx +++ b/Tracking/Acts/ActsGeometry/src/ActsATLASConverterTool.cxx @@ -56,7 +56,8 @@ ActsATLASConverterTool::initialize() { const auto* actsElement = dynamic_cast<const ActsDetectorElement*>(surface->associatedDetectorElement()); // Conversion from Acts to ATLAS surface impossible for the TRT so the TRT surface are not sure in this map - if(actsElement && actsElement->getSubdetector() != ActsDetectorElement::Subdetector::TRT + bool isTRT = (dynamic_cast<const InDetDD::TRT_BaseElement*>(actsElement->upstreamDetectorElement()) != nullptr); + if(actsElement && !isTRT && m_actsSurfaceMap.find(actsElement->identify()) == m_actsSurfaceMap.end()){ m_actsSurfaceMap.insert(std::pair<Identifier, const Acts::Surface &>(actsElement->identify(), *surface)); } diff --git a/Tracking/Acts/ActsGeometry/src/ActsDetectorElement.cxx b/Tracking/Acts/ActsGeometry/src/ActsDetectorElement.cxx index ece8cdc081646fc8aadbfc4834b303e46b98de57..cdd6ed9e10eeeabdbd6db021d0deead7d847e1e2 100644 --- a/Tracking/Acts/ActsGeometry/src/ActsDetectorElement.cxx +++ b/Tracking/Acts/ActsGeometry/src/ActsDetectorElement.cxx @@ -5,28 +5,36 @@ #include "ActsGeometry/ActsDetectorElement.h" // ATHENA -#include "TRT_ReadoutGeometry/TRT_EndcapElement.h" -#include "TRT_ReadoutGeometry/TRT_BarrelElement.h" #include "ActsInterop/IdentityHelper.h" +#include "GeoPrimitives/CLHEPtoEigenConverter.h" +#include "GeoPrimitives/GeoPrimitives.h" +#include "SCT_ReadoutGeometry/StripBoxDesign.h" +#include "SCT_ReadoutGeometry/StripStereoAnnulusDesign.h" +#include "TRT_ReadoutGeometry/TRT_BarrelElement.h" +#include "TRT_ReadoutGeometry/TRT_EndcapElement.h" +#include "TrkSurfaces/AnnulusBounds.h" #include "TrkSurfaces/RectangleBounds.h" #include "TrkSurfaces/Surface.h" #include "TrkSurfaces/SurfaceBounds.h" #include "TrkSurfaces/TrapezoidBounds.h" -#include "GeoPrimitives/CLHEPtoEigenConverter.h" // PACKAGE -#include "ActsGeometry/ActsTrackingGeometrySvc.h" #include "ActsGeometry/ActsAlignmentStore.h" #include "ActsGeometry/ActsGeometryContext.h" +#include "ActsGeometry/ActsTrackingGeometrySvc.h" // ACTS -#include "Acts/Surfaces/StrawSurface.hpp" +#include "Acts/Definitions/Units.hpp" +#include "Acts/Geometry/GeometryContext.hpp" +#include "Acts/Surfaces/AnnulusBounds.hpp" +#include "Acts/Surfaces/DiscSurface.hpp" #include "Acts/Surfaces/LineBounds.hpp" #include "Acts/Surfaces/PlaneSurface.hpp" #include "Acts/Surfaces/RectangleBounds.hpp" +#include "Acts/Surfaces/StrawSurface.hpp" #include "Acts/Surfaces/TrapezoidBounds.hpp" -#include "Acts/Geometry/GeometryContext.hpp" -#include "Acts/Definitions/Units.hpp" +#include "Acts/Visualization/ObjVisualization3D.hpp" +#include "Acts/Visualization/PlyVisualization3D.hpp" // STL #include <mutex> @@ -35,237 +43,274 @@ #include <boost/variant.hpp> #include <boost/variant/get.hpp> - -using Acts::Transform3; using Acts::Surface; +using Acts::Transform3; using namespace Acts::UnitLiterals; constexpr double length_unit = 1_mm; ActsDetectorElement::ActsDetectorElement( - const InDetDD::SiDetectorElement* detElem) -{ - m_detElement = detElem; + const InDetDD::SiDetectorElement &detElem) { + m_detElement = &detElem; - //auto center = detElem->center(); - auto boundsType = detElem->bounds().type(); + auto boundsType = detElem.bounds().type(); // thickness - m_thickness = detElem->thickness(); + m_thickness = detElem.thickness(); if (boundsType == Trk::SurfaceBounds::Rectangle) { - const InDetDD::SiDetectorDesign &design = detElem->design(); - double hlX = design.width()/2. * length_unit; - double hlY = design.length()/2. * length_unit; + const InDetDD::SiDetectorDesign &design = detElem.design(); + double hlX = design.width() / 2. * length_unit; + double hlY = design.length() / 2. * length_unit; - auto rectangleBounds = std::make_shared<const Acts::RectangleBounds>( - hlX, hlY); + auto rectangleBounds = + std::make_shared<const Acts::RectangleBounds>(hlX, hlY); m_bounds = rectangleBounds; - m_surface = Acts::Surface::makeShared<Acts::PlaneSurface>(rectangleBounds, *this); + if (const auto *bd = dynamic_cast<const InDetDD::StripBoxDesign *>(&design); + bd != nullptr) { + // extra shift for split row modules + m_extraTransform = Acts::Transform3::Identity() * + Acts::AngleAxis3{-M_PI / 2., Acts::Vector3::UnitZ()} * + bd->moduleShift() * + Acts::AngleAxis3{M_PI / 2., Acts::Vector3::UnitZ()}; + } - } else if (boundsType == Trk::SurfaceBounds::Trapezoid) { + m_surface = + Acts::Surface::makeShared<Acts::PlaneSurface>(rectangleBounds, *this); - std::string shapeString = detElem->getMaterialGeom()->getLogVol()->getShape()->type(); - //std::cout << __FUNCTION__ << "tapezoid, GeoLogVol -> shape says: " << shapeString << std::endl; + } else if (boundsType == Trk::SurfaceBounds::Trapezoid) { - const InDetDD::SiDetectorDesign &design = detElem->design(); + const InDetDD::SiDetectorDesign &design = detElem.design(); - double minHlX = design.minWidth()/2. * length_unit; - double maxHlX = design.maxWidth()/2. * length_unit; - double hlY = design.length()/2. * length_unit; + double minHlX = design.minWidth() / 2. * length_unit; + double maxHlX = design.maxWidth() / 2. * length_unit; + double hlY = design.length() / 2. * length_unit; - auto trapezoidBounds = std::make_shared<const Acts::TrapezoidBounds>( - minHlX, maxHlX, hlY); + auto trapezoidBounds = + std::make_shared<const Acts::TrapezoidBounds>(minHlX, maxHlX, hlY); m_bounds = trapezoidBounds; - m_surface = Acts::Surface::makeShared<Acts::PlaneSurface>(trapezoidBounds, *this); + m_surface = + Acts::Surface::makeShared<Acts::PlaneSurface>(trapezoidBounds, *this); + + } else if (boundsType == Trk::SurfaceBounds::Annulus) { + + const InDetDD::SiDetectorDesign &design = detElem.design(); + const auto *annulus = + dynamic_cast<const InDetDD::StripStereoAnnulusDesign *>(&design); + if (annulus == nullptr) { + throw std::domain_error("ActsDetectorElement got inconsistent surface"); + } + + double phi = annulus->phiWidth(); + double phiS = annulus->stereo(); + double R = annulus->waferCentreR(); + double maxR = annulus->maxR(); + double minR = annulus->minR(); + + double phiAvg = + 0; // phiAvg is the bounds-internal local rotation. We don't want one + + // // phi is the total opening angle, set up symmetric phi bounds + double phiMax = phi / 2.; + double phiMin = -phiMax; + + // // need to rotate pi/2 to reproduce ABXY orientation, phiS so that phi=0 + // is center and symmetric + double phiShift = M_PI / 2. - phiS; + + Amg::Vector2D originStripXYRotated(R * (1 - std::cos(phiS)), + R * std::sin(-phiS)); + + auto annulusBounds = std::make_shared<Acts::AnnulusBounds>( + minR, maxR, phiMin, phiMax, originStripXYRotated, phiAvg); + m_bounds = annulusBounds; + + m_surface = + Acts::Surface::makeShared<Acts::DiscSurface>(annulusBounds, *this); + + Amg::Vector2D origin2D = annulusBounds->moduleOrigin(); + Amg::Translation3D transl(Amg::Vector3D(origin2D.x(), origin2D.y(), 0)); + Amg::Rotation3D rot(Amg::AngleAxis3D(-phiShift, Amg::Vector3D::UnitZ())); + Amg::Transform3D originTrf; + originTrf = transl * rot; + m_extraTransform = originTrf.inverse(); } else { - throw std::domain_error("ActsDetectorElement does not support this surface type"); + std::cout << boundsType << std::endl; + throw std::domain_error( + "ActsDetectorElement does not support this surface type"); } } ActsDetectorElement::ActsDetectorElement( - const Acts::Transform3& trf, - const InDetDD::TRT_BaseElement* detElem, - const Identifier& id) -{ - m_detElement = detElem; - m_defTransform = std::make_shared<const Acts::Transform3>(trf); + const Acts::Transform3 &trf, const InDetDD::TRT_BaseElement &detElem, + const Identifier &id) { + m_detElement = &detElem; + m_defTransform = trf; m_explicitIdentifier = id; // we know this is a straw - double length = detElem->strawLength()*0.5 * length_unit; + double length = detElem.strawLength() * 0.5 * length_unit; // we need to find the radius - auto ecElem = dynamic_cast<const InDetDD::TRT_EndcapElement*>(detElem); - auto brlElem = dynamic_cast<const InDetDD::TRT_BarrelElement*>(detElem); + auto ecElem = dynamic_cast<const InDetDD::TRT_EndcapElement *>(&detElem); + auto brlElem = dynamic_cast<const InDetDD::TRT_BarrelElement *>(&detElem); double innerTubeRadius{0.}; if (ecElem) { innerTubeRadius = ecElem->getDescriptor()->innerTubeRadius() * length_unit; - } - else { - if (brlElem){ - innerTubeRadius = brlElem->getDescriptor()->innerTubeRadius() * length_unit; + } else { + if (brlElem) { + innerTubeRadius = + brlElem->getDescriptor()->innerTubeRadius() * length_unit; } else { - throw std::runtime_error("Cannot get tube radius for element in ActsDetectorElement c'tor"); + throw std::runtime_error( + "Cannot get tube radius for element in ActsDetectorElement c'tor"); } } - auto lineBounds = std::make_shared<const Acts::LineBounds>(innerTubeRadius, length); + auto lineBounds = + std::make_shared<const Acts::LineBounds>(innerTubeRadius, length); m_bounds = lineBounds; m_surface = Acts::Surface::makeShared<Acts::StrawSurface>(lineBounds, *this); } -IdentityHelper -ActsDetectorElement::identityHelper() const -{ - size_t which = m_detElement.which(); - if (which == 0) { - return IdentityHelper(boost::get<const InDetDD::SiDetectorElement*>(m_detElement)); +IdentityHelper ActsDetectorElement::identityHelper() const { + if (const auto *detElem = + dynamic_cast<const InDetDD::SiDetectorElement *>(m_detElement); + detElem != nullptr) { + return IdentityHelper(detElem); } else { throw std::domain_error("Cannot get IdentityHelper for TRT element"); } } -const Acts::Transform3& -ActsDetectorElement::transform(const Acts::GeometryContext& anygctx) const -{ +const Acts::Transform3 & +ActsDetectorElement::transform(const Acts::GeometryContext &anygctx) const { // any cast to known context type - const ActsGeometryContext* gctx = anygctx.get<const ActsGeometryContext*>(); + const ActsGeometryContext *gctx = anygctx.get<const ActsGeometryContext *>(); - // This is needed for initial geometry construction. At that point, we don't have a - // consistent view of the geometry yet, and thus we can't populate an alignment store - // at that time. + // This is needed for initial geometry construction. At that point, we don't + // have a consistent view of the geometry yet, and thus we can't populate an + // alignment store at that time. if (gctx->construction) { // this should only happen at initialize (1 thread, but mutex anyway) return getDefaultTransformMutexed(); } // unpack the alignment store from the context - const ActsAlignmentStore* alignmentStore = gctx->alignmentStore; + const ActsAlignmentStore *alignmentStore = gctx->alignmentStore; // no GAS, is this initialization? assert(alignmentStore != nullptr); // get the correct cached transform // units should be fine here since we converted at construction - const Acts::Transform3* cachedTrf = alignmentStore->getTransform(this); + const Acts::Transform3 *cachedTrf = alignmentStore->getTransform(this); assert(cachedTrf != nullptr); return *cachedTrf; } -void -ActsDetectorElement::storeTransform(ActsAlignmentStore* gas) const -{ - struct get_transform : public boost::static_visitor<Acts::Transform3> - { - get_transform(ActsAlignmentStore* gas2, const Acts::Transform3* trtTrf) - : m_store(gas2), m_trtTrf(trtTrf) {} - - Acts::Transform3 operator()(const InDetDD::SiDetectorElement* detElem) const - { - Amg::Transform3D g2l - = detElem->getMaterialGeom()->getAbsoluteTransform(m_store) - * Amg::CLHEPTransformToEigen(detElem->recoToHitTransform()); - - // need to make sure translation has correct units - g2l.translation() *= length_unit; - - return g2l; - } - - Acts::Transform3 operator()(const InDetDD::TRT_BaseElement*) const - { - return *m_trtTrf; - } - - ActsAlignmentStore* m_store; - const Acts::Transform3* m_trtTrf; - }; - - Acts::Transform3 trf - = boost::apply_visitor(get_transform(gas, m_defTransform.get()), m_detElement); +void ActsDetectorElement::storeTransform(ActsAlignmentStore *gas) const { + Amg::Transform3D trf; + if (const auto *detElem = + dynamic_cast<const InDetDD::SiDetectorElement *>(m_detElement); + detElem != nullptr) { + Amg::Transform3D l2g = + detElem->getMaterialGeom()->getAbsoluteTransform(gas) * + Amg::CLHEPTransformToEigen(detElem->recoToHitTransform()); + + // need to make sure translation has correct units + l2g.translation() *= 1.0 / CLHEP::mm * length_unit; + + l2g = l2g * m_extraTransform; + + trf = l2g; + } else if (const auto *detElem = + dynamic_cast<const InDetDD::TRT_BaseElement *>(m_detElement); + detElem != nullptr) { + // So far: NO ALIGNMENT for the ACTS TRT version. Default transform set in + // constructor, should be safe to access without mutex. + trf = *m_defTransform; + } else { + throw std::runtime_error{"Unknown detector element type"}; + } gas->setTransform(this, trf); if (gas->getTransform(this) == nullptr) { - throw std::runtime_error("Detector element was unable to store transform in GAS"); + throw std::runtime_error( + "Detector element was unable to store transform in GAS"); } - } -const Acts::Transform3& -ActsDetectorElement::getDefaultTransformMutexed() const -{ - struct get_default_transform : public boost::static_visitor<Acts::Transform3> - { - get_default_transform(const Acts::Transform3* trtTrf) : m_trtTrf(trtTrf) {} - - Acts::Transform3 operator()(const InDetDD::SiDetectorElement* detElem) const - { - Amg::Transform3D g2l - = detElem->getMaterialGeom()->getDefAbsoluteTransform() - * Amg::CLHEPTransformToEigen(detElem->recoToHitTransform()); - - // need to make sure translation has correct units - g2l.translation() *= length_unit; - - return g2l; - } - - Acts::Transform3 operator()(const InDetDD::TRT_BaseElement*) const - { - return *m_trtTrf; - } - - const Acts::Transform3* m_trtTrf; - }; - +const Acts::Transform3 & +ActsDetectorElement::getDefaultTransformMutexed() const { std::lock_guard<std::mutex> guard(m_cacheMutex); if (m_defTransform) { return *m_defTransform; } // transform not yet set - m_defTransform - = std::make_shared<const Acts::Transform3>( - boost::apply_visitor(get_default_transform(m_defTransform.get()), m_detElement)); + if (const auto *detElem = + dynamic_cast<const InDetDD::SiDetectorElement *>(m_detElement); + detElem != nullptr) { + Amg::Transform3D l2g = + detElem->getMaterialGeom()->getDefAbsoluteTransform() * + Amg::CLHEPTransformToEigen(detElem->recoToHitTransform()); + + l2g.translation() *= 1.0 / CLHEP::mm * length_unit; + + l2g = l2g * m_extraTransform; + + m_defTransform = l2g; + } else if (const auto *detElem = + dynamic_cast<const InDetDD::TRT_BaseElement *>(m_detElement); + detElem != nullptr) { + throw std::logic_error{ + "TRT transform should have been set in the constructor"}; + } else { + throw std::runtime_error{"Unknown detector element type"}; + } return *m_defTransform; } -const Acts::Surface& -ActsDetectorElement::surface() const -{ +const Acts::Surface &ActsDetectorElement::surface() const { return (*m_surface); } -const Trk::Surface& -ActsDetectorElement::atlasSurface() const -{ - size_t which = m_detElement.which(); - if (which == 0) { - return (boost::get<const InDetDD::SiDetectorElement *>(m_detElement))->surface(); - } - else { +const Trk::Surface &ActsDetectorElement::atlasSurface() const { + if (const auto *detElem = + dynamic_cast<const InDetDD::SiDetectorElement *>(m_detElement); + detElem != nullptr) { + return detElem->surface(); + } else { throw std::domain_error("Cannot get surface for TRT element"); } } -double -ActsDetectorElement::thickness() const -{ - return m_thickness; -} +double ActsDetectorElement::thickness() const { return m_thickness; } -Identifier -ActsDetectorElement::identify() const -{ - return boost::apply_visitor(IdVisitor(m_explicitIdentifier), m_detElement); +Identifier ActsDetectorElement::identify() const { + if (const auto *detElem = + dynamic_cast<const InDetDD::SiDetectorElement *>(m_detElement); + detElem != nullptr) { + return detElem->identify(); + } else if (dynamic_cast<const InDetDD::TRT_BaseElement *>(m_detElement) != + nullptr) { + return m_explicitIdentifier; + } else { + throw std::runtime_error{"Unknown detector element type"}; + } } + +const GeoVDetectorElement * +ActsDetectorElement::upstreamDetectorElement() const { + return m_detElement; +} \ No newline at end of file diff --git a/Tracking/Acts/ActsGeometry/src/ActsLayerBuilder.cxx b/Tracking/Acts/ActsGeometry/src/ActsLayerBuilder.cxx index ae761115a4b12142bec53d5f3e0512449f528c3e..0379f600df9fc81f90d19b9804f9e64a78b02b3a 100644 --- a/Tracking/Acts/ActsGeometry/src/ActsLayerBuilder.cxx +++ b/Tracking/Acts/ActsGeometry/src/ActsLayerBuilder.cxx @@ -11,325 +11,562 @@ #include "ActsInterop/IdentityHelper.h" // ACTS +#include "Acts/Definitions/Algebra.hpp" +#include "Acts/Definitions/Units.hpp" +#include "Acts/Geometry/ApproachDescriptor.hpp" +#include "Acts/Geometry/GenericApproachDescriptor.hpp" +#include "Acts/Geometry/GeometryContext.hpp" +#include "Acts/Geometry/LayerCreator.hpp" +#include "Acts/Geometry/ProtoLayer.hpp" #include "Acts/Material/ProtoSurfaceMaterial.hpp" #include "Acts/Surfaces/CylinderSurface.hpp" #include "Acts/Surfaces/DiscSurface.hpp" -#include "Acts/Geometry/GenericApproachDescriptor.hpp" -#include "Acts/Geometry/ApproachDescriptor.hpp" -#include "Acts/Geometry/ProtoLayer.hpp" -#include "Acts/Geometry/LayerCreator.hpp" -#include "Acts/Definitions/Algebra.hpp" -#include "Acts/Geometry/GeometryContext.hpp" -#include "Acts/Definitions/Units.hpp" #include "Acts/Utilities/BinningType.hpp" +#include "Acts/Surfaces/AnnulusBounds.hpp" + +#include "Acts/Visualization/GeometryView3D.hpp" +#include "Acts/Visualization/ObjVisualization3D.hpp" +#include <iterator> +#include <unordered_map> using Acts::Surface; using Acts::Transform3; using Acts::Translation3; - using namespace Acts::UnitLiterals; -const Acts::LayerVector -ActsLayerBuilder::negativeLayers(const Acts::GeometryContext& gctx) const +ActsLayerBuilder::ActsLayerBuilder(const Config& cfg, + std::unique_ptr<const Acts::Logger> logger) + : m_logger(std::move(logger)) { + m_cfg = cfg; + if(m_cfg.mode == Mode::Undefined) { + throw std::invalid_argument("ALB mode is undefined"); + } +} + +const Acts::LayerVector +ActsLayerBuilder::negativeLayers(const Acts::GeometryContext &gctx) const { ACTS_VERBOSE("Building negative layers"); // @todo Remove this hack once the m_elementStore mess is sorted out - auto mutableThis = const_cast<ActsLayerBuilder*>(this); + auto mutableThis = const_cast<ActsLayerBuilder *>(this); Acts::LayerVector nVector; - mutableThis->buildLayers(gctx, nVector, -1); + mutableThis->buildEndcap(gctx, nVector, -1); return nVector; } const Acts::LayerVector -ActsLayerBuilder::centralLayers(const Acts::GeometryContext& gctx) const -{ +ActsLayerBuilder::centralLayers(const Acts::GeometryContext &gctx) const { ACTS_VERBOSE("Building central layers"); // @todo Remove this hack once the m_elementStore mess is sorted out - auto mutableThis = const_cast<ActsLayerBuilder*>(this); + auto mutableThis = const_cast<ActsLayerBuilder *>(this); Acts::LayerVector cVector; - mutableThis->buildLayers(gctx, cVector, 0); + mutableThis->buildBarrel(gctx, cVector); return cVector; } const Acts::LayerVector -ActsLayerBuilder::positiveLayers(const Acts::GeometryContext& gctx) const -{ +ActsLayerBuilder::positiveLayers(const Acts::GeometryContext &gctx) const { ACTS_VERBOSE("Building positive layers"); // @todo Remove this hack once the m_elementStore mess is sorted out - auto mutableThis = const_cast<ActsLayerBuilder*>(this); + auto mutableThis = const_cast<ActsLayerBuilder *>(this); Acts::LayerVector pVector; - mutableThis->buildLayers(gctx, pVector, 1); + mutableThis->buildEndcap(gctx, pVector, 1); return pVector; - } std::vector<std::shared_ptr<const ActsDetectorElement>> -ActsLayerBuilder::getDetectorElements() const -{ - auto siDetMng = static_cast<const InDetDD::SiDetectorManager*>(m_cfg.mng); +ActsLayerBuilder::getDetectorElements() const { + ACTS_VERBOSE("Retrieving detector elements from detector manager"); + if (m_cfg.mng == nullptr) { + ACTS_ERROR("Manager is null"); + throw std::runtime_error{"Detector manager is null"}; + } + auto siDetMng = static_cast<const InDetDD::SiDetectorManager *>(m_cfg.mng); + ACTS_VERBOSE("Detector manager has " + << std::distance(siDetMng->getDetectorElementBegin(), + siDetMng->getDetectorElementEnd()) + << " elements"); std::vector<std::shared_ptr<const ActsDetectorElement>> elements; InDetDD::SiDetectorElementCollection::const_iterator iter; for (iter = siDetMng->getDetectorElementBegin(); - iter != siDetMng->getDetectorElementEnd(); - ++iter) { - const InDetDD::SiDetectorElement* siDetElement = *iter; + iter != siDetMng->getDetectorElementEnd(); ++iter) { + const InDetDD::SiDetectorElement *siDetElement = + dynamic_cast<InDetDD::SiDetectorElement *>(*iter); + if (siDetElement == nullptr) { + ACTS_ERROR("Detector element was nullptr"); + throw std::runtime_error{"Corrupt detector element collection"}; + } elements.push_back( - std::make_shared<const ActsDetectorElement>( - siDetElement)); + std::make_shared<const ActsDetectorElement>(*siDetElement)); } + ACTS_VERBOSE("Retrieved " << elements.size() << " elements"); return elements; } -void -ActsLayerBuilder::buildLayers(const Acts::GeometryContext& gctx, - Acts::LayerVector& layersOutput, int type) -{ - - std::vector<std::shared_ptr<const ActsDetectorElement>> elements = getDetectorElements(); +void ActsLayerBuilder::buildBarrel(const Acts::GeometryContext &gctx, + Acts::LayerVector &layersOutput) { + ACTS_VERBOSE("Build layers: BARREL"); + std::vector<std::shared_ptr<const ActsDetectorElement>> elements = + getDetectorElements(); - std::map<std::pair<int, int>, std::vector<std::shared_ptr<const Surface>>> layers; + std::map<int, std::vector<std::shared_ptr<const Acts::Surface>>> layers{}; for (const auto &element : elements) { IdentityHelper id = element->identityHelper(); + if (id.bec() != 0) + continue; - // wrong subdetector - // want barrel but not barrel - if (type == 0 && id.bec() != 0) continue; - // want endcap but is barrel - if (type != 0 && id.bec() == 0) continue; - // want endcap, but wrong side - if (type != 0 && type * id.bec() < 0) continue; + if(m_cfg.mode == Mode::ITkPixelInner) { + if(id.layer_disk() >= 2) continue; + } + if(m_cfg.mode == Mode::ITkPixelOuter) { + if(id.layer_disk() <= 1) continue; + } m_cfg.elementStore->push_back(element); - // const PixelID* idhlp = dynamic_cast<const - // PixelID*>(siDetElement->getIdHelper()); - // int layer_disk = idhlp->layer_disk(siDetElement->identify()); - // int eta_module = idhlp->eta_module(siDetElement->identify()); - - // std::cout << "SELECTED: "; - // std::cout << element->bec() << " "; - // std::cout << "(" << element->layer_disk() << " " << layer_disk << ") "; - // std::cout << "(" << element->eta_module() << " " << eta_module << ") "; - // std::cout << element->phi_module() << " "; - // std::cout << "Z = " << element->surface().center().z() << std::endl; - int elementLayer; elementLayer = id.layer_disk(); - std::pair<int, int> layerKey(elementLayer, id.bec()); + layers[elementLayer].push_back(element->surface().getSharedPtr()); + } - if (layers.count(layerKey) == 0) { - // layer not set yet - layers.insert(std::make_pair(layerKey, std::vector<std::shared_ptr<const Surface>>())); - } + // // @TODO: Exclude BCM, should happen via separate detector manager - // push into correct layer - layers.at(layerKey).push_back(element->surface().getSharedPtr()); + if (m_cfg.objDebugOutput) { + std::string lab; + for (auto &[key, layerSurfaces] : layers) { - } + std::stringstream name; - // @TODO: Maybe exclude BLM (BCM?) - ACTS_VERBOSE("Found " << layers.size() << " layers"); + + name << "obj/" << m_cfg.mode << "_brl_" << std::setfill('0') << std::setw(2) + << std::to_string(key) + "_bec" + << ".obj"; + + std::ofstream ofs{name.str()}; + Acts::ObjVisualization3D vis{}; + Acts::ViewConfig vc = Acts::s_viewSensitive; + vc.nSegments = 200; + for (const auto &surface : layerSurfaces) { + Acts::GeometryView3D::drawSurface(vis, *surface, gctx, + Acts::Transform3::Identity(), vc); + } + + vis.write(ofs); + } + } + + ACTS_VERBOSE("Found " << layers.size() << " barrel layers"); // pre loop for more concise output. only do if we actually print it - if(logger().doPrint(Acts::Logging::VERBOSE)) { + if (logger().doPrint(Acts::Logging::VERBOSE)) { size_t n = 1; - for (const auto& layerPair : layers) { - const std::vector<std::shared_ptr<const Surface>>& layerSurfaces = layerPair.second; - auto key = layerPair.first; - Acts::ProtoLayer pl(gctx, layerSurfaces); - ACTS_VERBOSE("Layer #" << n << " with layerKey: (" - << key.first << ", " << key.second << ")"); - if (type == 0) { // BARREL - ACTS_VERBOSE(" -> at rMin / rMax: " << pl.min(Acts::binR) << " / " << pl.max(Acts::binR)); - ACTS_VERBOSE(" -> at zMin / zMax: " << pl.min(Acts::binZ) << " / " << pl.max(Acts::binZ)); - } - else { - ACTS_VERBOSE(" -> at zMin / zMax: " << pl.min(Acts::binZ) << " / " << pl.max(Acts::binZ)); - ACTS_VERBOSE(" -> at rMin / rMax: " << pl.min(Acts::binR) << " / " << pl.max(Acts::binR)); - } + for (const auto &[key, surfaces] : layers) { + Acts::ProtoLayer pl(gctx, surfaces); + ACTS_VERBOSE("Layer #" << n << " with layerKey: (" << key << ")"); + ACTS_VERBOSE(" -> at rMin / rMax: " << pl.min(Acts::binR) << " / " + << pl.max(Acts::binR)); + ACTS_VERBOSE(" -> at zMin / zMax: " << pl.min(Acts::binZ) << " / " + << pl.max(Acts::binZ)); n++; } } - for (const auto& layerPair : layers) { + for (const auto &[key, surfaces] : layers) { std::unique_ptr<Acts::ApproachDescriptor> approachDescriptor = nullptr; std::shared_ptr<const Acts::ProtoSurfaceMaterial> materialProxy = nullptr; - // use ref here, copy later - const std::vector<std::shared_ptr<const Surface>>& layerSurfaces = layerPair.second; + // layers and extent are determined, build actual layer + Acts::ProtoLayer pl(gctx, surfaces); + pl.envelope[Acts::binR] = std::make_pair(2_mm, 2_mm); + pl.envelope[Acts::binZ] = std::make_pair(2_mm, 2_mm); + + double layerZ = pl.medium(Acts::binZ, true); + double layerHalfZ = 0.5 * pl.range(Acts::binZ); + + Acts::Transform3 transform(Translation3(0., 0., -layerZ)); + // set up approach descriptor + + std::shared_ptr<Acts::CylinderSurface> innerBoundary = + Acts::Surface::makeShared<Acts::CylinderSurface>( + transform, pl.min(Acts::binR), layerHalfZ); + + std::shared_ptr<Acts::CylinderSurface> outerBoundary = + Acts::Surface::makeShared<Acts::CylinderSurface>( + transform, pl.max(Acts::binR), layerHalfZ); + + std::shared_ptr<Acts::CylinderSurface> centralSurface = + Acts::Surface::makeShared<Acts::CylinderSurface>( + transform, (pl.min(Acts::binR) + pl.max(Acts::binR)) / 2., + layerHalfZ); + + size_t binsPhi = m_cfg.barrelMaterialBins.first; + size_t binsZ = m_cfg.barrelMaterialBins.second; + + Acts::BinUtility materialBinUtil(binsPhi, -M_PI, M_PI, Acts::closed, + Acts::binPhi); + materialBinUtil += Acts::BinUtility(binsZ, -layerHalfZ, layerHalfZ, + Acts::open, Acts::binZ, transform); + + materialProxy = + std::make_shared<const Acts::ProtoSurfaceMaterial>(materialBinUtil); + + ACTS_VERBOSE("[L] Layer is marked to carry support material on Surface ( " + "inner=0 / center=1 / outer=2 ) : " + << "inner"); + ACTS_VERBOSE("with binning: [" << binsPhi << ", " << binsZ << "]"); + + ACTS_VERBOSE("Created ApproachSurfaces for cylinder layer at:"); + ACTS_VERBOSE(" - inner: R=" << pl.min(Acts::binR)); + ACTS_VERBOSE(" - central: R=" << (pl.min(Acts::binR) + pl.max(Acts::binR)) / + 2.); + ACTS_VERBOSE(" - outer: R=" << pl.max(Acts::binR)); + + // set material on inner + innerBoundary->assignSurfaceMaterial(materialProxy); + + std::vector<std::shared_ptr<const Acts::Surface>> aSurfaces; + aSurfaces.push_back(std::move(innerBoundary)); + aSurfaces.push_back(std::move(centralSurface)); + aSurfaces.push_back(std::move(outerBoundary)); + + approachDescriptor = + std::make_unique<Acts::GenericApproachDescriptor>(std::move(aSurfaces)); + + std::shared_ptr<Acts::Layer> layer; + if(m_cfg.mode == Mode::ITkStrip) { + // @TODO: This needs to be fixed! + size_t nBinsPhi = 1; + size_t nBinsZ = 1; + layer = m_cfg.layerCreator->cylinderLayer( + gctx, surfaces, nBinsPhi, nBinsZ, pl, transform, + std::move(approachDescriptor)); + } + else { + layer = m_cfg.layerCreator->cylinderLayer( + gctx, surfaces, Acts::equidistant, Acts::equidistant, pl, transform, + std::move(approachDescriptor)); + } - if (type == 0) { // BARREL - // layers and extent are determined, build actual layer - Acts::ProtoLayer pl(gctx, layerSurfaces); - pl.envelope[Acts::binR] = std::make_pair(0_mm, 0_mm); - pl.envelope[Acts::binZ] = std::make_pair(20_mm, 20_mm); + layersOutput.push_back(layer); + } +} + +void ActsLayerBuilder::buildEndcap(const Acts::GeometryContext &gctx, + Acts::LayerVector &layersOutput, int type) { + + ACTS_VERBOSE("Build layers: " << (type < 0 ? "NEGATIVE" : "POSITIVE") + << " ENDCAP"); + std::vector<std::shared_ptr<const ActsDetectorElement>> elements = + getDetectorElements(); + std::map<std::tuple<int, int, int>, std::vector<const Acts::Surface *>> + initialLayers{}; + + for (const auto &element : elements) { + + IdentityHelper id = element->identityHelper(); + if (id.bec() * type <= 0) { + continue; + } + + if(m_cfg.mode == Mode::ITkPixelInner) { + if(id.layer_disk() >= 3) continue; + } + if(m_cfg.mode == Mode::ITkPixelOuter) { + if(id.layer_disk() <= 2) continue; + } + + m_cfg.elementStore->push_back(element); + + std::tuple<int, int, int> key{id.bec(), id.layer_disk(), id.eta_module()}; + + initialLayers[key].push_back(&element->surface()); + } + + ACTS_VERBOSE("Found " << initialLayers.size() << " " + << (type < 0 ? "NEGATIVE" : "POSITIVE") + << " ENDCAP inital layers"); + + std::vector<Acts::ProtoLayer> protoLayers; + protoLayers.reserve(initialLayers.size()); - double binPosZ = 0.5 * (pl.min(Acts::binZ) + pl.max(Acts::binZ)); - double envZShift = 0.5 * (-pl.envelope[Acts::binZ].first + pl.envelope[Acts::binZ].second); - double layerZ = binPosZ + envZShift; - double layerHalfZ - = std::abs(pl.max(Acts::binZ) + pl.envelope[Acts::binZ].second - layerZ); + for (const auto &[key, surfaces] : initialLayers) { + auto &pl = protoLayers.emplace_back(gctx, surfaces); + pl.envelope[Acts::binR] = std::make_pair(2_mm, 2_mm); + pl.envelope[Acts::binZ] = std::make_pair(2_mm, 2_mm); + } + + // sort proto layers by their medium z position + std::sort(protoLayers.begin(), protoLayers.end(), + [type](const Acts::ProtoLayer &a, const Acts::ProtoLayer &b) { + double midA = (a.min(Acts::binZ) + a.max(Acts::binZ)) / 2.0; + double midB = (b.min(Acts::binZ) + b.max(Acts::binZ)) / 2.0; + if (type < 0) { + return midA < midB; + } else { + return midA > midB; + } + }); + + auto plPrintZ = [](const auto &pl) -> std::string { + std::stringstream ss; + double zMid = (pl.min(Acts::binZ) + pl.max(Acts::binZ)) / 2.0; + ss << " < " << pl.min(Acts::binZ) << " | " << zMid << " | " + << pl.max(Acts::binZ) << " > "; + return ss.str(); + }; + if (logger().doPrint(Acts::Logging::VERBOSE)) { + for (const auto &pl : protoLayers) { + + ACTS_VERBOSE(" -> at < zMin | zMid | zMax >: " << plPrintZ(pl)); + + ACTS_VERBOSE(" -> at rMin / rMax: " << pl.min(Acts::binR) << " / " + << pl.max(Acts::binR)); + } + } - Acts::Transform3 transform(Translation3(0., 0., -layerZ)); - // set up approach descriptor + std::vector<Acts::ProtoLayer> mergedProtoLayers; + mergedProtoLayers.push_back(protoLayers.front()); + + std::vector<const Acts::Surface *> surfaces; + for (size_t i = 1; i < protoLayers.size(); i++) { + auto &pl = protoLayers[i]; + auto &pl_prev = mergedProtoLayers.back(); + // do they intersect? + ACTS_VERBOSE("Compare: " << plPrintZ(pl_prev) << " and " << plPrintZ(pl)); + bool overlap = (pl.min(Acts::binZ) <= pl_prev.max(Acts::binZ) && + pl.max(Acts::binZ) >= pl_prev.min(Acts::binZ)); + ACTS_VERBOSE(" -> overlap? " << (overlap ? "yes" : "no")); + if (overlap) { + ACTS_VERBOSE(" ===> merging"); + surfaces.clear(); + surfaces.reserve(pl.surfaces().size() + pl_prev.surfaces().size()); + surfaces.insert(surfaces.end(), pl.surfaces().begin(), + pl.surfaces().end()); + surfaces.insert(surfaces.end(), pl_prev.surfaces().begin(), + pl_prev.surfaces().end()); + mergedProtoLayers.pop_back(); + auto &new_pl = mergedProtoLayers.emplace_back(gctx, std::move(surfaces)); + new_pl.envelope[Acts::binR] = pl.envelope[Acts::binR]; + new_pl.envelope[Acts::binZ] = pl.envelope[Acts::binZ]; + } else { + mergedProtoLayers.push_back(std::move(pl)); + } + } - std::shared_ptr<Acts::CylinderSurface> innerBoundary - = Acts::Surface::makeShared<Acts::CylinderSurface>(transform, pl.min(Acts::binR), layerHalfZ); + ACTS_VERBOSE("" << mergedProtoLayers.size() << " " + << (type < 0 ? "NEGATIVE" : "POSITIVE") + << " ENDCAP layers remain after merging"); - std::shared_ptr<Acts::CylinderSurface> outerBoundary - = Acts::Surface::makeShared<Acts::CylinderSurface>(transform, pl.max(Acts::binR), layerHalfZ); + for (size_t i = 0; i < mergedProtoLayers.size(); i++) { - std::shared_ptr<Acts::CylinderSurface> centralSurface - = Acts::Surface::makeShared<Acts::CylinderSurface>(transform, (pl.min(Acts::binR) + pl.max(Acts::binR))/2., layerHalfZ); + std::stringstream ss; + ss << "obj/" << m_cfg.mode << "_" << (type < 0 ? "neg" : "pos") << "_disk_" << std::setfill('0') + << std::setw(2) << i << ".obj"; - size_t binsPhi = m_cfg.barrelMaterialBins.first; - size_t binsZ = m_cfg.barrelMaterialBins.second; + std::ofstream ofs{ss.str()}; + Acts::ObjVisualization3D vis{}; + Acts::ViewConfig vc = Acts::s_viewSensitive; + vc.nSegments = 200; + for (const auto &surface : mergedProtoLayers[i].surfaces()) { + Acts::GeometryView3D::drawSurface(vis, *surface, gctx, + Acts::Transform3::Identity(), vc); - Acts::BinUtility materialBinUtil( - binsPhi, -M_PI, M_PI, Acts::closed, Acts::binPhi); - materialBinUtil += Acts::BinUtility( - binsZ, -layerHalfZ, layerHalfZ, Acts::open, Acts::binZ, transform); - - materialProxy - = std::make_shared<const Acts::ProtoSurfaceMaterial>(materialBinUtil); - - ACTS_VERBOSE("[L] Layer is marked to carry support material on Surface ( " - "inner=0 / center=1 / outer=2 ) : " << "inner"); - ACTS_VERBOSE("with binning: [" << binsPhi << ", " << binsZ << "]"); - - ACTS_VERBOSE("Created ApproachSurfaces for cylinder layer at:"); - ACTS_VERBOSE(" - inner: R=" << pl.min(Acts::binR)); - ACTS_VERBOSE(" - central: R=" << (pl.min(Acts::binR) + pl.max(Acts::binR))/2.); - ACTS_VERBOSE(" - outer: R=" << pl.max(Acts::binR)); - - // set material on inner - // @TODO: make this configurable somehow - innerBoundary->assignSurfaceMaterial(materialProxy); - - std::vector<std::shared_ptr<const Acts::Surface>> aSurfaces; - aSurfaces.push_back(std::move(innerBoundary)); - aSurfaces.push_back(std::move(centralSurface)); - aSurfaces.push_back(std::move(outerBoundary)); - - approachDescriptor - = std::make_unique<Acts::GenericApproachDescriptor>(std::move(aSurfaces)); - - auto layer = m_cfg.layerCreator->cylinderLayer(gctx, - layerSurfaces, - Acts::equidistant, - Acts::equidistant, - pl, - transform, - std::move(approachDescriptor)); - - layersOutput.push_back(layer); - } else { // ENDCAP - Acts::ProtoLayer pl(gctx, layerSurfaces); - pl.envelope[Acts::binR] = std::make_pair(0_mm, 0_mm); - pl.envelope[Acts::binZ] = std::make_pair(10_mm, 10_mm); - - // copied from layercreator - double layerZ - = 0.5 * (pl.min(Acts::binZ) - pl.envelope[Acts::binZ].first + pl.max(Acts::binZ) - + pl.envelope[Acts::binZ].second); - double layerThickness = (pl.max(Acts::binZ) - pl.min(Acts::binZ)) - + pl.envelope[Acts::binZ].first + pl.envelope[Acts::binZ].second; - - double layerZInner = layerZ - layerThickness/2.; - double layerZOuter = layerZ + layerThickness/2.; - - if (std::abs(layerZInner) > std::abs(layerZOuter)) std::swap(layerZInner, layerZOuter); - - Acts::Transform3 transformNominal(Translation3(0., 0., layerZ)); - Acts::Transform3 transformInner(Translation3(0., 0., layerZInner)); - Acts::Transform3 transformOuter(Translation3(0., 0., layerZOuter)); - - std::shared_ptr<Acts::DiscSurface> innerBoundary = Acts::Surface::makeShared<Acts::DiscSurface>(transformInner, pl.min(Acts::binR), pl.max(Acts::binR)); - - std::shared_ptr<Acts::DiscSurface> nominalSurface = Acts::Surface::makeShared<Acts::DiscSurface>(transformNominal, pl.min(Acts::binR), pl.max(Acts::binR)); - - std::shared_ptr<Acts::DiscSurface> outerBoundary = Acts::Surface::makeShared<Acts::DiscSurface>(transformOuter, pl.min(Acts::binR), pl.max(Acts::binR)); - - size_t matBinsPhi = m_cfg.endcapMaterialBins.first; - size_t matBinsR = m_cfg.endcapMaterialBins.second; - - Acts::BinUtility materialBinUtil( - matBinsPhi, -M_PI, M_PI, Acts::closed, Acts::binPhi); - materialBinUtil += Acts::BinUtility( - matBinsR, pl.min(Acts::binR), pl.max(Acts::binR), Acts::open, Acts::binR, transformNominal); - - materialProxy = std::make_shared<const Acts::ProtoSurfaceMaterial>(materialBinUtil); - - ACTS_VERBOSE("[L] Layer is marked to carry support material on Surface ( " - "inner=0 / center=1 / outer=2 ) : " - << "inner"); - ACTS_VERBOSE("with binning: [" << matBinsPhi << ", " << matBinsR << "]"); - - ACTS_VERBOSE("Created ApproachSurfaces for disc layer at:"); - ACTS_VERBOSE(" - inner: Z=" << layerZInner); - ACTS_VERBOSE(" - central: Z=" << layerZ); - ACTS_VERBOSE(" - outer: Z=" << layerZOuter); - - // set material on inner - // @TODO: make this configurable somehow - innerBoundary->assignSurfaceMaterial(materialProxy); - - int nModPhi = std::numeric_limits<int>::max(); - int nModR = 0; - - // want to figure out bins in phi - for (const auto &srf : layerSurfaces) - { - auto elm = dynamic_cast<const ActsDetectorElement *>(srf->associatedDetectorElement()); - if (elm) - { - auto id = elm->identityHelper(); - int phi_mod_max = id.phi_module_max(); - int eta_mod_max = id.eta_module_max(); - nModPhi = std::min(nModPhi, phi_mod_max + 1); - nModR = eta_mod_max + 1; - } + } + + vis.write(ofs); + } + + + + std::vector<std::shared_ptr<const Surface>> ownedSurfaces; + for (const auto &pl : mergedProtoLayers) { + + std::unique_ptr<Acts::ApproachDescriptor> approachDescriptor = nullptr; + std::shared_ptr<const Acts::ProtoSurfaceMaterial> materialProxy = nullptr; + + double layerZ = pl.medium(Acts::binZ); + double layerHalfZ = 0.5 * pl.range(Acts::binZ); + double layerThickness = pl.range(Acts::binZ); + + double layerZInner = layerZ - layerHalfZ; + double layerZOuter = layerZ + layerHalfZ; + + if (std::abs(layerZInner) > std::abs(layerZOuter)) + std::swap(layerZInner, layerZOuter); + + std::vector<std::shared_ptr<const Acts::Surface>> aSurfaces; + + Acts::Transform3 transformNominal(Translation3(0., 0., layerZ)); + Acts::Transform3 transformInner(Translation3(0., 0., layerZInner)); + Acts::Transform3 transformOuter(Translation3(0., 0., layerZOuter)); + + std::shared_ptr<Acts::DiscSurface> innerBoundary = + Acts::Surface::makeShared<Acts::DiscSurface>( + transformInner, pl.min(Acts::binR), pl.max(Acts::binR)); + aSurfaces.push_back(innerBoundary); + + std::shared_ptr<Acts::DiscSurface> nominalSurface = + Acts::Surface::makeShared<Acts::DiscSurface>( + transformNominal, pl.min(Acts::binR), pl.max(Acts::binR)); + aSurfaces.push_back(nominalSurface); + + std::shared_ptr<Acts::DiscSurface> outerBoundary = + Acts::Surface::makeShared<Acts::DiscSurface>( + transformOuter, pl.min(Acts::binR), pl.max(Acts::binR)); + aSurfaces.push_back(outerBoundary); + + if(layerThickness > 2_mm) { + ACTS_VERBOSE("Wide disc layer ("<< layerThickness << ") => adding cylinder like approach surfaces"); + Acts::Transform3 trf{Translation3{0, 0, layerZ}}; + auto cylinderInner = + Acts::Surface::makeShared<Acts::CylinderSurface>( + trf, pl.min(Acts::binR), layerHalfZ); + aSurfaces.push_back(cylinderInner); + + auto cylinderOuter = + Acts::Surface::makeShared<Acts::CylinderSurface>( + trf, pl.max(Acts::binR), layerHalfZ); + aSurfaces.push_back(cylinderOuter); + } + + + size_t matBinsPhi = m_cfg.endcapMaterialBins.first; + size_t matBinsR = m_cfg.endcapMaterialBins.second; + + Acts::BinUtility materialBinUtil(matBinsPhi, -M_PI, M_PI, Acts::closed, + Acts::binPhi); + materialBinUtil += + Acts::BinUtility(matBinsR, pl.min(Acts::binR), pl.max(Acts::binR), + Acts::open, Acts::binR, transformNominal); + + materialProxy = + std::make_shared<const Acts::ProtoSurfaceMaterial>(materialBinUtil); + + ACTS_VERBOSE("[L] Layer is marked to carry support material on Surface ( " + "inner=0 / center=1 / outer=2 ) : " + << "inner"); + ACTS_VERBOSE("with binning: [" << matBinsPhi << ", " << matBinsR << "]"); + + ACTS_VERBOSE("Created ApproachSurfaces for disc layer at:"); + ACTS_VERBOSE(" - inner: Z=" << layerZInner); + ACTS_VERBOSE(" - central: Z=" << layerZ); + ACTS_VERBOSE(" - outer: Z=" << layerZOuter); + + // set material on inner + innerBoundary->assignSurfaceMaterial(materialProxy); + + + // std::set<int> ringIds; + bool isITk = m_cfg.mode == Mode::ITkPixelInner || + m_cfg.mode == Mode::ITkPixelOuter || + m_cfg.mode == Mode::ITkStrip; + + std::map<int, std::set<int>> phiModuleByRing; + + // want to figure out bins in phi + for (const auto &srf : pl.surfaces()) { + auto elm = dynamic_cast<const ActsDetectorElement *>( + srf->associatedDetectorElement()); + if (elm) { + auto id = elm->identityHelper(); + int ring_number; + if(m_cfg.mode == Mode::ITkStrip || !isITk) { + ring_number = id.eta_module(); + } + else { + ring_number = id.layer_disk(); + } + phiModuleByRing[ring_number].insert(id.phi_module()); } + } + size_t nModPhi = std::numeric_limits<size_t>::max(); + for(const auto& [ring, phiModules] : phiModuleByRing) { + nModPhi = std::min(nModPhi, phiModules.size()); + } - ACTS_VERBOSE("Identifier reports: " << nModPhi << " is lowest for " << nModR << " r-rings"); + size_t nModR = phiModuleByRing.size(); - // the modules in the innermost r-rings are exactly shifted by one half module width - // since it's the same number of modules, this gives binning trouble. Reduce bins - // by half: about 2 module pairs should be in each bin. This should be fine. + ACTS_VERBOSE("Identifier reports: " << nModPhi << " is lowest for " << nModR + << " r-rings"); + + size_t nBinsPhi = nModPhi; + size_t nBinsR = nModR; + if(!isITk) { + // In the ID, the modules in the innermost r-rings are exactly shifted by + // one half module width since it's the same number of modules, this gives + // binning trouble. Reduce bins by half: about 2 module pairs should be in + // each bin. This should be fine. // @TODO: evaluate - size_t nBinsPhi = nModPhi/2.; - size_t nBinsR = nModR; + nBinsPhi /= 2.0; + } + if(m_cfg.mode == Mode::ITkStrip) { + // @FIXME: ACTS surface binning for disc surfaces is currently not working + // correctly, see: + // https://github.com/acts-project/acts/issues/851. + // Use a single bin as a workaround. + nBinsR = 1; + nBinsPhi = 1; + } - ACTS_VERBOSE("Creating r x phi binned layer with " << nBinsR << " x " << nBinsPhi << " bins"); + ACTS_VERBOSE("Creating r x phi binned layer with " << nBinsR << " x " + << nBinsPhi << " bins"); - std::vector<std::shared_ptr<const Acts::Surface>> aSurfaces; - aSurfaces.push_back(std::move(innerBoundary)); - aSurfaces.push_back(std::move(nominalSurface)); - aSurfaces.push_back(std::move(outerBoundary)); - approachDescriptor - = std::make_unique<Acts::GenericApproachDescriptor>(aSurfaces); + approachDescriptor = + std::make_unique<Acts::GenericApproachDescriptor>(aSurfaces); - auto layer = m_cfg.layerCreator->discLayer(gctx, - layerSurfaces, - nBinsR, - nBinsPhi, - pl, - transformNominal, - std::move(approachDescriptor)); + // get ownership of pl surfaces + ownedSurfaces.clear(); + ownedSurfaces.reserve(pl.surfaces().size()); + std::transform(pl.surfaces().begin(), pl.surfaces().end(), + std::back_inserter(ownedSurfaces), + [](const auto &s) { return s->getSharedPtr(); }); - layersOutput.push_back(layer); - } + auto layer = m_cfg.layerCreator->discLayer(gctx, ownedSurfaces, nBinsR, + nBinsPhi, pl, transformNominal, + std::move(approachDescriptor)); + + layersOutput.push_back(layer); } } + +std::ostream& operator<<(std::ostream& os, const ActsLayerBuilder::Mode& mode) { + + using Mode = ActsLayerBuilder::Mode; + switch(mode) { + case Mode::Undefined: + os << "Undefined"; + break; + case Mode::Pixel: + os << "Pixel"; + break; + case Mode::SCT: + os << "SCT"; + break; + case Mode::TRT: + os << "TRT"; + break; + case Mode::ITkPixelInner: + os << "ITkPixelInner"; + break; + case Mode::ITkPixelOuter: + os << "ITkPixelOuter"; + break; + case Mode::ITkStrip: + os << "ITkStrip"; + break; + } + + return os; +} \ No newline at end of file diff --git a/Tracking/Acts/ActsGeometry/src/ActsStrawLayerBuilder.cxx b/Tracking/Acts/ActsGeometry/src/ActsStrawLayerBuilder.cxx index 408a07f71395278721e27664c9dbbf1d7831f15a..43dc951b753ab0240a8f11cc262d1136e2142e95 100644 --- a/Tracking/Acts/ActsGeometry/src/ActsStrawLayerBuilder.cxx +++ b/Tracking/Acts/ActsGeometry/src/ActsStrawLayerBuilder.cxx @@ -120,7 +120,7 @@ ActsStrawLayerBuilder::centralLayers(const Acts::GeometryContext& gctx) istraw); auto elem = std::make_shared<const ActsDetectorElement>( - trf, brlElem, straw_id); + trf, *brlElem, straw_id); m_cfg.elementStore->push_back(elem); @@ -221,7 +221,7 @@ ActsStrawLayerBuilder::endcapLayers(const Acts::GeometryContext& gctx, int side) auto elem = std::make_shared<const ActsDetectorElement>( - trf, ecElem, straw_id); + trf, *ecElem, straw_id); m_cfg.elementStore->push_back(elem); diff --git a/Tracking/Acts/ActsGeometry/src/ActsTrackingGeometrySvc.cxx b/Tracking/Acts/ActsGeometry/src/ActsTrackingGeometrySvc.cxx index 3a8eb1fff9f17bb430508ca75774b3fb84b12069..fee64a5d2b1ba7cd4e7fb9eefb422ff396327d00 100644 --- a/Tracking/Acts/ActsGeometry/src/ActsTrackingGeometrySvc.cxx +++ b/Tracking/Acts/ActsGeometry/src/ActsTrackingGeometrySvc.cxx @@ -5,99 +5,112 @@ #include "ActsGeometry/ActsTrackingGeometrySvc.h" // ATHENA -#include "InDetReadoutGeometry/SiDetectorManager.h" -#include "TRT_ReadoutGeometry/TRT_DetectorManager.h" -#include "StoreGate/StoreGateSvc.h" #include "GaudiKernel/EventContext.h" +#include "GeoPrimitives/GeoPrimitives.h" #include "InDetIdentifier/TRT_ID.h" +#include "InDetReadoutGeometry/SiDetectorManager.h" +#include "StoreGate/StoreGateSvc.h" +#include "TRT_ReadoutGeometry/TRT_DetectorManager.h" // ACTS -#include "Acts/Geometry/TrackingGeometry.hpp" -#include "Acts/Geometry/TrackingVolume.hpp" +#include "Acts/ActsVersion.hpp" +#include "Acts/Geometry/CylinderVolumeBounds.hpp" +#include "Acts/Geometry/CylinderVolumeBuilder.hpp" +#include "Acts/Geometry/CylinderVolumeHelper.hpp" #include "Acts/Geometry/ITrackingVolumeBuilder.hpp" #include "Acts/Geometry/LayerArrayCreator.hpp" -#include "Acts/Geometry/SurfaceArrayCreator.hpp" #include "Acts/Geometry/LayerCreator.hpp" -#include "Acts/Geometry/TrackingVolumeArrayCreator.hpp" -#include "Acts/Geometry/CylinderVolumeHelper.hpp" +#include "Acts/Geometry/SurfaceArrayCreator.hpp" +#include "Acts/Geometry/TrackingGeometry.hpp" #include "Acts/Geometry/TrackingGeometryBuilder.hpp" -#include "Acts/Geometry/CylinderVolumeBuilder.hpp" -#include "Acts/Geometry/CylinderVolumeBounds.hpp" -#include <Acts/Plugins/Json/MaterialMapJsonConverter.hpp> -#include <Acts/Plugins/Json/JsonMaterialDecorator.hpp> +#include "Acts/Geometry/TrackingVolume.hpp" +#include "Acts/Geometry/TrackingVolumeArrayCreator.hpp" #include "Acts/Utilities/Logger.hpp" -#include "Acts/ActsVersion.hpp" +#include "Acts/Definitions/Units.hpp" +#include <Acts/Plugins/Json/JsonMaterialDecorator.hpp> +#include <Acts/Plugins/Json/MaterialMapJsonConverter.hpp> // PACKAGE +#include "ActsGeometry/ActsAlignmentStore.h" +#include "ActsGeometry/ActsDetectorElement.h" +#include "ActsGeometry/ActsGeometryContext.h" #include "ActsGeometry/ActsLayerBuilder.h" #include "ActsGeometry/ActsStrawLayerBuilder.h" -#include "ActsGeometry/ActsDetectorElement.h" #include "ActsInterop/IdentityHelper.h" #include "ActsInterop/Logger.h" -#include "ActsGeometry/ActsAlignmentStore.h" -#include "ActsGeometry/ActsGeometryContext.h" +using namespace Acts::UnitLiterals; -ActsTrackingGeometrySvc::ActsTrackingGeometrySvc(const std::string& name, ISvcLocator* svc) - : base_class(name,svc), - m_detStore("StoreGateSvc/DetectorStore", name) -{ - m_elementStore = std::make_shared<std::vector<std::shared_ptr<const ActsDetectorElement>>>(); +ActsTrackingGeometrySvc::ActsTrackingGeometrySvc(const std::string &name, + ISvcLocator *svc) + : base_class(name, svc), m_detStore("StoreGateSvc/DetectorStore", name) { + m_elementStore = std::make_shared< + std::vector<std::shared_ptr<const ActsDetectorElement>>>(); } -StatusCode -ActsTrackingGeometrySvc::initialize() -{ +StatusCode ActsTrackingGeometrySvc::initialize() { ATH_MSG_INFO(name() << " is initializing"); - // FIXME: ActsCaloTrackingVolumeBuilder holds ReadHandle to CaloDetDescrManager. - // Hopefully this service is never called before that object is available. + // FIXME: ActsCaloTrackingVolumeBuilder holds ReadHandle to + // CaloDetDescrManager. Hopefully this service is never called before that + // object is available. m_autoRetrieveTools = false; m_checkToolDeps = false; - ATH_MSG_INFO("Acts version is: v" << Acts::VersionMajor << "." - << Acts::VersionMinor << "." - << Acts::VersionPatch - << " [" << Acts::CommitHash << "]"); - - ATH_CHECK ( m_detStore->retrieve(p_pixelManager, "Pixel") ); - ATH_CHECK ( m_detStore->retrieve(p_SCTManager, "SCT") ); - ATH_CHECK ( m_detStore->retrieve(p_TRTManager, "TRT") ); + ATH_MSG_INFO("ACTS version is: v" + << Acts::VersionMajor << "." << Acts::VersionMinor << "." + << Acts::VersionPatch << " [" << Acts::CommitHash << "]"); - if (m_detStore->retrieve(m_TRT_idHelper, "TRT_ID").isFailure()) { - msg(MSG::ERROR) << "Could not retrieve TRT ID Helper" << endmsg; + // load which subdetectors to build from property + std::set<std::string> buildSubdet(m_buildSubdetectors.begin(), + m_buildSubdetectors.end()); + ATH_MSG_INFO("Configured to build " << buildSubdet.size() + << " subdetectors:"); + for (const auto &s : buildSubdet) { + ATH_MSG_INFO(" - " << s); } + ATH_MSG_DEBUG("Loading detector manager(s)"); + if (buildSubdet.find("Pixel") != buildSubdet.end()) { + ATH_CHECK(m_detStore->retrieve(p_pixelManager, "Pixel")); + } + if (buildSubdet.find("SCT") != buildSubdet.end()) { + ATH_CHECK(m_detStore->retrieve(p_SCTManager, "SCT")); + } + if (buildSubdet.find("TRT") != buildSubdet.end()) { + ATH_CHECK(m_detStore->retrieve(p_TRTManager, "TRT")); + ATH_CHECK(m_detStore->retrieve(m_TRT_idHelper, "TRT_ID")); + } + if (buildSubdet.find("ITkPixel") != buildSubdet.end()) { + ATH_CHECK(m_detStore->retrieve(p_ITkPixelManager, "ITkPixel")); + } + if (buildSubdet.find("ITkStrip") != buildSubdet.end()) { + ATH_CHECK(m_detStore->retrieve(p_ITkStripManager, "ITkStrip")); + } + ATH_MSG_DEBUG("Setting up ACTS geometry helpers"); Acts::LayerArrayCreator::Config lacCfg; auto layerArrayCreator = std::make_shared<const Acts::LayerArrayCreator>( lacCfg, makeActsAthenaLogger(this, "LayArrCrtr", "ActsTGSvc")); Acts::TrackingVolumeArrayCreator::Config tvcCfg; - auto trackingVolumeArrayCreator - = std::make_shared<const Acts::TrackingVolumeArrayCreator>( + auto trackingVolumeArrayCreator = + std::make_shared<const Acts::TrackingVolumeArrayCreator>( tvcCfg, makeActsAthenaLogger(this, "TrkVolArrCrtr", "ActsTGSvc")); Acts::CylinderVolumeHelper::Config cvhConfig; - cvhConfig.layerArrayCreator = layerArrayCreator; + cvhConfig.layerArrayCreator = layerArrayCreator; cvhConfig.trackingVolumeArrayCreator = trackingVolumeArrayCreator; - auto cylinderVolumeHelper - = std::make_shared<const Acts::CylinderVolumeHelper>( - cvhConfig, makeActsAthenaLogger(this, "CylVolHlpr", "ActsTGSvc")); - - // load which subdetectors to build from property - std::set<std::string> buildSubdet(m_buildSubdetectors.begin(), m_buildSubdetectors.end()); - ATH_MSG_INFO("Configured to build:"); - for(const auto& s : buildSubdet) { - ATH_MSG_INFO(" - " << s); - } + auto cylinderVolumeHelper = + std::make_shared<const Acts::CylinderVolumeHelper>( + cvhConfig, makeActsAthenaLogger(this, "CylVolHlpr", "ActsTGSvc")); Acts::TrackingGeometryBuilder::Config tgbConfig; - tgbConfig.trackingVolumeHelper = cylinderVolumeHelper; + tgbConfig.trackingVolumeHelper = cylinderVolumeHelper; - if(m_useMaterialMap){ + if (m_useMaterialMap) { std::shared_ptr<const Acts::IMaterialDecorator> matDeco = nullptr; std::string matFile = m_materialMapInputFile; if (matFile.find(".json") != std::string::npos) { @@ -108,119 +121,214 @@ ActsTrackingGeometrySvc::initialize() jsonGeoConvConfig, m_materialMapInputFile, Acts::Logging::INFO); } tgbConfig.materialDecorator = matDeco; - } + } + + // default geometry context, this is nominal + ActsGeometryContext constructionContext; + constructionContext.construction = true; try { // PIXEL - if(buildSubdet.count("Pixel") > 0) { - tgbConfig.trackingVolumeBuilders.push_back([&]( - const auto & gctx, const auto& inner, const auto&) { - auto lb = makeLayerBuilder(p_pixelManager); - Acts::CylinderVolumeBuilder::Config cvbConfig; - cvbConfig.layerEnvelopeR = {5, 5}; - cvbConfig.layerEnvelopeZ = 2; - cvbConfig.trackingVolumeHelper = cylinderVolumeHelper; - cvbConfig.volumeSignature = 0; - cvbConfig.volumeName = "Pixel"; - cvbConfig.layerBuilder = lb; - cvbConfig.buildToRadiusZero = true; - - Acts::CylinderVolumeBuilder cvb(cvbConfig, - makeActsAthenaLogger(this, "CylVolBldr", "ActsTGSvc")); - - return cvb.trackingVolume(gctx, inner); + if (buildSubdet.count("Pixel") > 0) { + tgbConfig.trackingVolumeBuilders.push_back([&](const auto &gctx, + const auto &inner, + const auto &) { + auto cfg = makeLayerBuilderConfig(p_pixelManager); + cfg.mode = ActsLayerBuilder::Mode::Pixel; + auto lb = std::make_shared<ActsLayerBuilder>( + cfg, makeActsAthenaLogger(this, "PixelGMSLayBldr", "ActsTGSvc")); + Acts::CylinderVolumeBuilder::Config cvbConfig; + cvbConfig.layerEnvelopeR = {5_mm, 5_mm}; + cvbConfig.layerEnvelopeZ = 1_mm; + cvbConfig.trackingVolumeHelper = cylinderVolumeHelper; + cvbConfig.volumeSignature = 0; + cvbConfig.volumeName = "Pixel"; + cvbConfig.layerBuilder = lb; + cvbConfig.buildToRadiusZero = true; + + Acts::CylinderVolumeBuilder cvb( + cvbConfig, makeActsAthenaLogger(this, "CylVolBldr", "ActsTGSvc")); + + return cvb.trackingVolume(gctx, inner); }); } - + + // ITK PIXEL + if (buildSubdet.count("ITkPixel") > 0) { + tgbConfig.trackingVolumeBuilders.push_back( + [&](const auto &gctx, const auto &inner, const auto &) { + auto cfg = makeLayerBuilderConfig(p_ITkPixelManager); + cfg.mode = ActsLayerBuilder::Mode::ITkPixelInner; + cfg.objDebugOutput = m_objDebugOutput; + auto lb = std::make_shared<ActsLayerBuilder>( + cfg, makeActsAthenaLogger(this, "ITkPxInLb", "ActsTGSvc")); + + Acts::CylinderVolumeBuilder::Config cvbConfig; + cvbConfig.layerEnvelopeR = {5_mm, 5_mm}; + cvbConfig.layerEnvelopeZ = 1_mm; + cvbConfig.trackingVolumeHelper = cylinderVolumeHelper; + cvbConfig.volumeSignature = 0; + cvbConfig.volumeName = "ITkPixelInner"; + cvbConfig.layerBuilder = lb; + cvbConfig.buildToRadiusZero = true; + + Acts::CylinderVolumeBuilder cvb( + cvbConfig, + makeActsAthenaLogger(this, "CylVolBldr", "ActsTGSvc")); + + return cvb.trackingVolume(gctx, inner); + }); + + tgbConfig.trackingVolumeBuilders.push_back( + [&](const auto &gctx, const auto &inner, const auto &) { + auto cfg = makeLayerBuilderConfig(p_ITkPixelManager); + cfg.mode = ActsLayerBuilder::Mode::ITkPixelOuter; + cfg.objDebugOutput = m_objDebugOutput; + auto lb = std::make_shared<ActsLayerBuilder>( + cfg, makeActsAthenaLogger(this, "ITkPxOtLb", "ActsTGSvc")); + + Acts::CylinderVolumeBuilder::Config cvbConfig; + cvbConfig.layerEnvelopeR = {5_mm, 5_mm}; + cvbConfig.layerEnvelopeZ = 1_mm; + cvbConfig.trackingVolumeHelper = cylinderVolumeHelper; + cvbConfig.volumeSignature = 0; + cvbConfig.volumeName = "ITkPixelOuter"; + cvbConfig.layerBuilder = lb; + cvbConfig.buildToRadiusZero = false; + + Acts::CylinderVolumeBuilder cvb( + cvbConfig, + makeActsAthenaLogger(this, "CylVolBldr", "ActsTGSvc")); + + return cvb.trackingVolume(gctx, inner); + }); + } + + // ITK STRIP + if (buildSubdet.count("ITkStrip") > 0) { + tgbConfig.trackingVolumeBuilders.push_back( + [&](const auto &gctx, const auto &inner, const auto &) { + auto cfg = makeLayerBuilderConfig(p_ITkStripManager); + cfg.mode = ActsLayerBuilder::Mode::ITkStrip; + cfg.objDebugOutput = m_objDebugOutput; + auto lb = std::make_shared<ActsLayerBuilder>( + cfg, makeActsAthenaLogger(this, "ITkStripLB", "ActsTGSvc")); + + Acts::CylinderVolumeBuilder::Config cvbConfig; + cvbConfig.layerEnvelopeR = {5_mm, 5_mm}; + cvbConfig.layerEnvelopeZ = 1_mm; + cvbConfig.trackingVolumeHelper = cylinderVolumeHelper; + cvbConfig.volumeSignature = 0; + cvbConfig.volumeName = "ITkStrip"; + cvbConfig.layerBuilder = lb; + cvbConfig.buildToRadiusZero = buildSubdet.count("ITkPixel") == 0; + + Acts::CylinderVolumeBuilder cvb( + cvbConfig, + makeActsAthenaLogger(this, "CylVolBldr", "ActsTGSvc")); + + return cvb.trackingVolume(gctx, inner); + }); + } + bool buildSCT = buildSubdet.count("SCT") > 0; bool buildTRT = buildSubdet.count("TRT") > 0; - if(buildSCT && buildTRT) { + if (buildSCT && buildTRT) { // building both we need to take care - tgbConfig.trackingVolumeBuilders.push_back([&]( - const auto& gctx, const auto& inner, const auto&) { - auto sct_lb = makeLayerBuilder(p_SCTManager); - auto trt_lb = makeLayerBuilder(p_TRTManager); + tgbConfig.trackingVolumeBuilders.push_back( + [&](const auto &gctx, const auto &inner, const auto &) { + auto cfg = makeLayerBuilderConfig(p_SCTManager); + cfg.mode = ActsLayerBuilder::Mode::SCT; + auto sct_lb = std::make_shared<ActsLayerBuilder>( + cfg, makeActsAthenaLogger(this, "SCTGMSLayBldr", "ActsTGSvc")); - return makeSCTTRTAssembly(gctx, *sct_lb, *trt_lb, *cylinderVolumeHelper, inner); - }); + auto trt_lb = makeStrawLayerBuilder(p_TRTManager); + return makeSCTTRTAssembly(gctx, *sct_lb, *trt_lb, + *cylinderVolumeHelper, inner); + }); } else if (buildSCT) { - tgbConfig.trackingVolumeBuilders.push_back([&]( - const auto& gctx, const auto& inner, const auto&) { - auto lb = makeLayerBuilder(p_SCTManager); - - Acts::CylinderVolumeBuilder::Config cvbConfig; - cvbConfig.layerEnvelopeR = {5, 5}; - cvbConfig.layerEnvelopeZ = 2; - cvbConfig.trackingVolumeHelper = cylinderVolumeHelper; - cvbConfig.volumeSignature = 0; - cvbConfig.volumeName = "SCT"; - cvbConfig.layerBuilder = lb; - cvbConfig.buildToRadiusZero = false; - - Acts::CylinderVolumeBuilder cvb(cvbConfig, - makeActsAthenaLogger(this, "SCTCylVolBldr", "ActsTGSvc")); - - return cvb.trackingVolume(gctx, inner); - }); + tgbConfig.trackingVolumeBuilders.push_back( + [&](const auto &gctx, const auto &inner, const auto &) { + auto lbCfg = makeLayerBuilderConfig(p_SCTManager); + lbCfg.mode = ActsLayerBuilder::Mode::SCT; + auto lb = std::make_shared<ActsLayerBuilder>( + lbCfg, + makeActsAthenaLogger(this, "SCTGMSLayBldr", "ActsTGSvc")); + + Acts::CylinderVolumeBuilder::Config cvbConfig; + cvbConfig.layerEnvelopeR = {5_mm, 5_mm}; + cvbConfig.layerEnvelopeZ = 2_mm; + cvbConfig.trackingVolumeHelper = cylinderVolumeHelper; + cvbConfig.volumeSignature = 0; + cvbConfig.volumeName = "SCT"; + cvbConfig.layerBuilder = lb; + cvbConfig.buildToRadiusZero = false; + + Acts::CylinderVolumeBuilder cvb( + cvbConfig, + makeActsAthenaLogger(this, "SCTCylVolBldr", "ActsTGSvc")); + + return cvb.trackingVolume(gctx, inner); + }); } else if (buildTRT) { - tgbConfig.trackingVolumeBuilders.push_back([&]( - const auto& gctx, const auto& inner, const auto&) { - auto lb = makeLayerBuilder(p_TRTManager); - Acts::CylinderVolumeBuilder::Config cvbConfig; - cvbConfig.layerEnvelopeR = {5, 5}; - cvbConfig.layerEnvelopeZ = 2; - cvbConfig.trackingVolumeHelper = cylinderVolumeHelper; - cvbConfig.volumeSignature = 0; - cvbConfig.volumeName = "TRT"; - cvbConfig.layerBuilder = lb; - cvbConfig.buildToRadiusZero = false; - - Acts::CylinderVolumeBuilder cvb(cvbConfig, - makeActsAthenaLogger(this, "TRTCylVolBldr", "ActsTGSvc")); - - return cvb.trackingVolume(gctx, inner); - }); + tgbConfig.trackingVolumeBuilders.push_back( + [&](const auto &gctx, const auto &inner, const auto &) { + auto lb = makeStrawLayerBuilder(p_TRTManager); + Acts::CylinderVolumeBuilder::Config cvbConfig; + cvbConfig.layerEnvelopeR = {5_mm, 5_mm}; + cvbConfig.layerEnvelopeZ = 2_mm; + cvbConfig.trackingVolumeHelper = cylinderVolumeHelper; + cvbConfig.volumeSignature = 0; + cvbConfig.volumeName = "TRT"; + cvbConfig.layerBuilder = lb; + cvbConfig.buildToRadiusZero = false; + + Acts::CylinderVolumeBuilder cvb( + cvbConfig, + makeActsAthenaLogger(this, "TRTCylVolBldr", "ActsTGSvc")); + + return cvb.trackingVolume(gctx, inner); + }); } // Calo - if(buildSubdet.count("Calo") > 0) { - tgbConfig.trackingVolumeBuilders.push_back([&]( - const auto& gctx, const auto& inner, const auto&) { - return m_caloVolumeBuilder->trackingVolume(gctx, inner, nullptr); - }); + if (buildSubdet.count("Calo") > 0) { + tgbConfig.trackingVolumeBuilders.push_back( + [&](const auto &gctx, const auto &inner, const auto &) { + return m_caloVolumeBuilder->trackingVolume(gctx, inner, nullptr); + }); } - } - catch (const std::exception& e) { + } catch (const std::exception &e) { ATH_MSG_ERROR("Encountered error when building Acts tracking geometry"); ATH_MSG_ERROR(e.what()); return StatusCode::FAILURE; } - - auto trackingGeometryBuilder - = std::make_shared<const Acts::TrackingGeometryBuilder>(tgbConfig, - makeActsAthenaLogger(this, "TrkGeomBldr", "ActsTGSvc")); - - - // default geometry context, this is nominal - ActsGeometryContext constructionContext; - constructionContext.construction = true; + auto trackingGeometryBuilder = + std::make_shared<const Acts::TrackingGeometryBuilder>( + tgbConfig, makeActsAthenaLogger(this, "TrkGeomBldr", "ActsTGSvc")); ATH_MSG_VERBOSE("Begin building process"); - m_trackingGeometry = trackingGeometryBuilder - ->trackingGeometry(constructionContext.context()); + m_trackingGeometry = + trackingGeometryBuilder->trackingGeometry(constructionContext.context()); ATH_MSG_VERBOSE("Building process completed"); + if (!m_trackingGeometry) { + ATH_MSG_ERROR("No ACTS tracking geometry was built. Cannot proceeed"); + return StatusCode::FAILURE; + } + ATH_MSG_VERBOSE("Building nominal alignment store"); - ActsAlignmentStore* nominalAlignmentStore = new ActsAlignmentStore(); + ActsAlignmentStore *nominalAlignmentStore = new ActsAlignmentStore(); populateAlignmentStore(nominalAlignmentStore); // manage ownership - m_nominalAlignmentStore = std::unique_ptr<const ActsAlignmentStore>(nominalAlignmentStore); + m_nominalAlignmentStore = + std::unique_ptr<const ActsAlignmentStore>(nominalAlignmentStore); ATH_MSG_INFO("Acts TrackingGeometry construction completed"); @@ -234,159 +342,153 @@ ActsTrackingGeometrySvc::trackingGeometry() { return m_trackingGeometry; } -std::shared_ptr<const Acts::ILayerBuilder> -ActsTrackingGeometrySvc::makeLayerBuilder(const InDetDD::InDetDetectorManager* manager) -{ - std::string managerName = manager->getName(); - - std::shared_ptr<const Acts::ILayerBuilder> gmLayerBuilder; - if (manager->getName() == "TRT") { - auto matcher = [](const Acts::GeometryContext& /*gctx*/, Acts::BinningValue /*bValue*/, const Acts::Surface* /*aS*/, - const Acts::Surface* /*bS*/) -> bool { - return false; - }; - - Acts::SurfaceArrayCreator::Config sacCfg; - sacCfg.surfaceMatcher = matcher; - sacCfg.doPhiBinningOptimization = false; - - auto surfaceArrayCreator = std::make_shared<Acts::SurfaceArrayCreator>( - sacCfg, - makeActsAthenaLogger(this, managerName+"SrfArrCrtr", "ActsTGSvc")); - Acts::LayerCreator::Config lcCfg; - lcCfg.surfaceArrayCreator = surfaceArrayCreator; - auto layerCreator = std::make_shared<Acts::LayerCreator>( - lcCfg, makeActsAthenaLogger(this, managerName+"LayCrtr", "ActsTGSvc")); - - ActsStrawLayerBuilder::Config cfg; - cfg.mng = static_cast<const InDetDD::TRT_DetectorManager*>(manager); - cfg.elementStore = m_elementStore; - cfg.layerCreator = layerCreator; - cfg.idHelper = m_TRT_idHelper; - gmLayerBuilder = std::make_shared<const ActsStrawLayerBuilder>(cfg, - makeActsAthenaLogger(this, managerName+"GMSLayBldr", "ActsTGSvc")); - - //gmLayerBuilder->centralLayers(); - //gmLayerBuilder->negativeLayers(); - //gmLayerBuilder->positiveLayers(); - } - else { - auto matcher = [](const Acts::GeometryContext& /*gctx*/, Acts::BinningValue bValue, - const Acts::Surface* aS, const Acts::Surface* bS) -> bool { - - auto a = dynamic_cast<const ActsDetectorElement*>(aS->associatedDetectorElement()); - auto b = dynamic_cast<const ActsDetectorElement*>(bS->associatedDetectorElement()); - if ((not a) or (not b)){ - throw std::runtime_error("Cast of surface associated element to ActsDetectorElement failed in ActsTrackingGeometrySvc::makeVolumeBuilder"); - } - - IdentityHelper idA = a->identityHelper(); - IdentityHelper idB = b->identityHelper(); - - // check if same bec - // can't be same if not - if(idA.bec() != idB.bec()) return false; - - if (bValue == Acts::binPhi) { - //std::cout << idA.phi_module() << " <-> " << idB.phi_module() << std::endl; - return idA.phi_module() == idB.phi_module(); - } - - if (bValue == Acts::binZ) { - return (idA.eta_module() == idB.eta_module()) - && (idA.layer_disk() == idB.layer_disk()) - && (idA.bec() == idB.bec()); - } - - if (bValue == Acts::binR) { - return (idA.eta_module() == idB.eta_module()) - && (idA.layer_disk() == idB.layer_disk()) - && (idB.bec() == idA.bec()); - } - - return false; - }; +std::shared_ptr<const Acts::ILayerBuilder> +ActsTrackingGeometrySvc::makeStrawLayerBuilder( + const InDetDD::InDetDetectorManager *manager) { - Acts::SurfaceArrayCreator::Config sacCfg; - sacCfg.surfaceMatcher = matcher; + std::string managerName = manager->getName(); + auto matcher = [](const Acts::GeometryContext & /*gctx*/, + Acts::BinningValue /*bValue*/, const Acts::Surface * /*aS*/, + const Acts::Surface * + /*bS*/) -> bool { return false; }; + + Acts::SurfaceArrayCreator::Config sacCfg; + sacCfg.surfaceMatcher = matcher; + sacCfg.doPhiBinningOptimization = false; + + auto surfaceArrayCreator = std::make_shared<Acts::SurfaceArrayCreator>( + sacCfg, + makeActsAthenaLogger(this, managerName + "SrfArrCrtr", "ActsTGSvc")); + Acts::LayerCreator::Config lcCfg; + lcCfg.surfaceArrayCreator = surfaceArrayCreator; + auto layerCreator = std::make_shared<Acts::LayerCreator>( + lcCfg, makeActsAthenaLogger(this, managerName + "LayCrtr", "ActsTGSvc")); + + ActsStrawLayerBuilder::Config cfg; + cfg.mng = static_cast<const InDetDD::TRT_DetectorManager *>(manager); + cfg.elementStore = m_elementStore; + cfg.layerCreator = layerCreator; + cfg.idHelper = m_TRT_idHelper; + return std::make_shared<const ActsStrawLayerBuilder>( + cfg, makeActsAthenaLogger(this, managerName + "GMSLayBldr", "ActsTGSvc")); +} - auto surfaceArrayCreator = std::make_shared<Acts::SurfaceArrayCreator>( - sacCfg, - makeActsAthenaLogger(this, managerName+"SrfArrCrtr", "ActsTGSvc")); - Acts::LayerCreator::Config lcCfg; - lcCfg.surfaceArrayCreator = surfaceArrayCreator; - auto layerCreator = std::make_shared<Acts::LayerCreator>( - lcCfg, makeActsAthenaLogger(this, managerName+"LayCrtr", "ActsTGSvc")); +ActsLayerBuilder::Config ActsTrackingGeometrySvc::makeLayerBuilderConfig( + const InDetDD::InDetDetectorManager *manager) { + std::string managerName = manager->getName(); + std::shared_ptr<const Acts::ILayerBuilder> gmLayerBuilder; + auto matcher = [](const Acts::GeometryContext & /*gctx*/, + Acts::BinningValue bValue, const Acts::Surface *aS, + const Acts::Surface *bS) -> bool { + auto a = dynamic_cast<const ActsDetectorElement *>( + aS->associatedDetectorElement()); + auto b = dynamic_cast<const ActsDetectorElement *>( + bS->associatedDetectorElement()); + if ((not a) or (not b)) { + throw std::runtime_error( + "Cast of surface associated element to ActsDetectorElement failed " + "in ActsTrackingGeometrySvc::makeVolumeBuilder"); + } + IdentityHelper idA = a->identityHelper(); + IdentityHelper idB = b->identityHelper(); - ActsLayerBuilder::Config cfg; + // check if same bec + // can't be same if not + if (idA.bec() != idB.bec()) + return false; - if(managerName == "Pixel") { - cfg.subdetector = ActsDetectorElement::Subdetector::Pixel; - } - else { - cfg.subdetector = ActsDetectorElement::Subdetector::SCT; + if (bValue == Acts::binPhi) { + // std::cout << idA.phi_module() << " <-> " << idB.phi_module() << + // std::endl; + return idA.phi_module() == idB.phi_module(); } - // set bins from configuration - if (m_barrelMaterialBins.size() != 2) { - throw std::invalid_argument("Number of barrel material bin counts != 2"); + if (bValue == Acts::binZ) { + return (idA.eta_module() == idB.eta_module()) && + (idA.layer_disk() == idB.layer_disk()) && (idA.bec() == idB.bec()); } - std::vector<size_t> brlBins(m_barrelMaterialBins); - cfg.barrelMaterialBins = {brlBins.at(0), - brlBins.at(1)}; - if (m_endcapMaterialBins.size() != 2) { - throw std::invalid_argument("Number of endcap material bin counts != 2"); + if (bValue == Acts::binR) { + return (idA.eta_module() == idB.eta_module()) && + (idA.layer_disk() == idB.layer_disk()) && (idB.bec() == idA.bec()); } - std::vector<size_t> ecBins(m_endcapMaterialBins); - cfg.endcapMaterialBins = {ecBins.at(0), - ecBins.at(1)}; - cfg.mng = static_cast<const InDetDD::SiDetectorManager*>(manager); - // use class member element store - cfg.elementStore = m_elementStore; - cfg.layerCreator = layerCreator; + return false; + }; + + Acts::SurfaceArrayCreator::Config sacCfg; + sacCfg.surfaceMatcher = matcher; - gmLayerBuilder = std::make_shared<const ActsLayerBuilder>(cfg, - makeActsAthenaLogger(this, managerName+"GMLayBldr", "ActsTGSvc")); + auto surfaceArrayCreator = std::make_shared<Acts::SurfaceArrayCreator>( + sacCfg, + makeActsAthenaLogger(this, managerName + "SrfArrCrtr", "ActsTGSvc")); + Acts::LayerCreator::Config lcCfg; + lcCfg.surfaceArrayCreator = surfaceArrayCreator; + auto layerCreator = std::make_shared<Acts::LayerCreator>( + lcCfg, makeActsAthenaLogger(this, managerName + "LayCrtr", "ActsTGSvc")); + + ActsLayerBuilder::Config cfg; + + // set bins from configuration + if (m_barrelMaterialBins.size() != 2) { + throw std::invalid_argument("Number of barrel material bin counts != 2"); } + std::vector<size_t> brlBins(m_barrelMaterialBins); + cfg.barrelMaterialBins = {brlBins.at(0), brlBins.at(1)}; + if (m_endcapMaterialBins.size() != 2) { + throw std::invalid_argument("Number of endcap material bin counts != 2"); + } + std::vector<size_t> ecBins(m_endcapMaterialBins); + cfg.endcapMaterialBins = {ecBins.at(0), ecBins.at(1)}; - return gmLayerBuilder; + cfg.mng = static_cast<const InDetDD::SiDetectorManager *>(manager); + // use class member element store + cfg.elementStore = m_elementStore; + cfg.layerCreator = layerCreator; + // gmLayerBuilder = std::make_shared<const ActsLayerBuilder>( + // cfg, makeActsAthenaLogger(this, managerName + "GMLayBldr", + // "ActsTGSvc")); + // return gmLayerBuilder; + return cfg; } - std::shared_ptr<Acts::TrackingVolume> -ActsTrackingGeometrySvc::makeSCTTRTAssembly(const Acts::GeometryContext& gctx, - const Acts::ILayerBuilder& sct_lb, const Acts::ILayerBuilder& trt_lb, - const Acts::CylinderVolumeHelper& cvh, const std::shared_ptr<const Acts::TrackingVolume>& pixel) { +ActsTrackingGeometrySvc::makeSCTTRTAssembly( + const Acts::GeometryContext &gctx, const Acts::ILayerBuilder &sct_lb, + const Acts::ILayerBuilder &trt_lb, const Acts::CylinderVolumeHelper &cvh, + const std::shared_ptr<const Acts::TrackingVolume> &pixel) { Acts::CylinderVolumeBuilder::Config cvbCfg; - Acts::CylinderVolumeBuilder cvb(cvbCfg, makeActsAthenaLogger(this, "SCTTRTCVB", "ActsTGSvc")); + Acts::CylinderVolumeBuilder cvb( + cvbCfg, makeActsAthenaLogger(this, "SCTTRTCVB", "ActsTGSvc")); ATH_MSG_VERBOSE("Making SCT negative layers: "); - Acts::VolumeConfig sctNegEC = cvb.analyzeContent(gctx, sct_lb.negativeLayers(gctx), {}); + Acts::VolumeConfig sctNegEC = + cvb.analyzeContent(gctx, sct_lb.negativeLayers(gctx), {}); ATH_MSG_VERBOSE("Making SCT positive layers: "); - Acts::VolumeConfig sctPosEC = cvb.analyzeContent(gctx, sct_lb.positiveLayers(gctx), {}); + Acts::VolumeConfig sctPosEC = + cvb.analyzeContent(gctx, sct_lb.positiveLayers(gctx), {}); ATH_MSG_VERBOSE("Making SCT central layers: "); - Acts::VolumeConfig sctBrl = cvb.analyzeContent(gctx, sct_lb.centralLayers(gctx), {}); - + Acts::VolumeConfig sctBrl = + cvb.analyzeContent(gctx, sct_lb.centralLayers(gctx), {}); ATH_MSG_VERBOSE("Making TRT negative layers: "); - Acts::VolumeConfig trtNegEC = cvb.analyzeContent(gctx, trt_lb.negativeLayers(gctx), {}); + Acts::VolumeConfig trtNegEC = + cvb.analyzeContent(gctx, trt_lb.negativeLayers(gctx), {}); ATH_MSG_VERBOSE("Making TRT positive layers: "); - Acts::VolumeConfig trtPosEC = cvb.analyzeContent(gctx, trt_lb.positiveLayers(gctx), {}); + Acts::VolumeConfig trtPosEC = + cvb.analyzeContent(gctx, trt_lb.positiveLayers(gctx), {}); ATH_MSG_VERBOSE("Making TRT central layers: "); - Acts::VolumeConfig trtBrl = cvb.analyzeContent(gctx, trt_lb.centralLayers(gctx), {}); - + Acts::VolumeConfig trtBrl = + cvb.analyzeContent(gctx, trt_lb.centralLayers(gctx), {}); - // synchronize trt - + double absZMinEC = std::min(std::abs(trtNegEC.zMax), std::abs(trtPosEC.zMin)); double absZMaxEC = std::max(std::abs(trtNegEC.zMin), std::abs(trtPosEC.zMax)); @@ -394,36 +496,46 @@ ActsTrackingGeometrySvc::makeSCTTRTAssembly(const Acts::GeometryContext& gctx, trtNegEC.zMax = -absZMinEC; trtPosEC.zMin = absZMinEC; trtPosEC.zMax = absZMaxEC; - + using CVBBV = Acts::CylinderVolumeBounds::BoundValues; // if pixel is present, shrink SCT volumes in R - if(pixel) { + if (pixel) { ATH_MSG_VERBOSE("Shrinking SCT to fit around Pixel"); - auto pixelBounds = dynamic_cast<const Acts::CylinderVolumeBounds*>(&pixel->volumeBounds()); + auto pixelBounds = dynamic_cast<const Acts::CylinderVolumeBounds *>( + &pixel->volumeBounds()); sctNegEC.rMin = pixelBounds->get(CVBBV::eMaxR); sctPosEC.rMin = pixelBounds->get(CVBBV::eMaxR); sctBrl.rMin = pixelBounds->get(CVBBV::eMaxR); - } - else { + } else { ATH_MSG_VERBOSE("Pixel is not configured, not wrapping"); } - ATH_MSG_VERBOSE("SCT Volume Configuration:"); - ATH_MSG_VERBOSE("- SCT::NegativeEndcap: " << sctNegEC.layers.size() << " layers, " << sctNegEC.toString()); - ATH_MSG_VERBOSE("- SCT::Barrel: " << sctBrl.layers.size() << " layers, " << sctBrl.toString()); - ATH_MSG_VERBOSE("- SCT::PositiveEncap: " << sctPosEC.layers.size() << " layers, " << sctPosEC.toString()); + ATH_MSG_VERBOSE("- SCT::NegativeEndcap: " << sctNegEC.layers.size() + << " layers, " + << sctNegEC.toString()); + ATH_MSG_VERBOSE("- SCT::Barrel: " << sctBrl.layers.size() << " layers, " + << sctBrl.toString()); + ATH_MSG_VERBOSE("- SCT::PositiveEncap: " << sctPosEC.layers.size() + << " layers, " + << sctPosEC.toString()); ATH_MSG_VERBOSE("TRT Volume Configuration:"); - ATH_MSG_VERBOSE("- TRT::NegativeEndcap: " << trtNegEC.layers.size() << " layers, " << trtNegEC.toString()); - ATH_MSG_VERBOSE("- TRT::Barrel: " << trtBrl.layers.size() << " layers, " << trtBrl.toString()); - ATH_MSG_VERBOSE("- TRT::PositiveEncap: " << trtPosEC.layers.size() << " layers, " << trtPosEC.toString()); - - // harmonize SCT BRL <-> EC, normally the CVB does this, but we're skipping that - sctBrl.zMax = (sctBrl.zMax + sctPosEC.zMin)/2.; + ATH_MSG_VERBOSE("- TRT::NegativeEndcap: " << trtNegEC.layers.size() + << " layers, " + << trtNegEC.toString()); + ATH_MSG_VERBOSE("- TRT::Barrel: " << trtBrl.layers.size() << " layers, " + << trtBrl.toString()); + ATH_MSG_VERBOSE("- TRT::PositiveEncap: " << trtPosEC.layers.size() + << " layers, " + << trtPosEC.toString()); + + // harmonize SCT BRL <-> EC, normally the CVB does this, but we're skipping + // that + sctBrl.zMax = (sctBrl.zMax + sctPosEC.zMin) / 2.; sctBrl.zMin = -sctBrl.zMax; - + // and now harmonize everything // inflate TRT Barrel to match SCT trtBrl.zMin = sctBrl.zMin; @@ -450,97 +562,106 @@ ActsTrackingGeometrySvc::makeSCTTRTAssembly(const Acts::GeometryContext& gctx, ATH_MSG_VERBOSE("Dimensions after synchronization between SCT and TRT"); ATH_MSG_VERBOSE("SCT Volume Configuration:"); - ATH_MSG_VERBOSE("- SCT::NegativeEndcap: " << sctNegEC.layers.size() << " layers, " << sctNegEC.toString()); - ATH_MSG_VERBOSE("- SCT::Barrel: " << sctBrl.layers.size() << " layers, " << sctBrl.toString()); - ATH_MSG_VERBOSE("- SCT::PositiveEncap: " << sctPosEC.layers.size() << " layers, " << sctPosEC.toString()); + ATH_MSG_VERBOSE("- SCT::NegativeEndcap: " << sctNegEC.layers.size() + << " layers, " + << sctNegEC.toString()); + ATH_MSG_VERBOSE("- SCT::Barrel: " << sctBrl.layers.size() << " layers, " + << sctBrl.toString()); + ATH_MSG_VERBOSE("- SCT::PositiveEncap: " << sctPosEC.layers.size() + << " layers, " + << sctPosEC.toString()); ATH_MSG_VERBOSE("TRT Volume Configuration:"); - ATH_MSG_VERBOSE("- TRT::NegativeEndcap: " << trtNegEC.layers.size() << " layers, " << trtNegEC.toString()); - ATH_MSG_VERBOSE("- TRT::Barrel: " << trtBrl.layers.size() << " layers, " << trtBrl.toString()); - ATH_MSG_VERBOSE("- TRT::PositiveEncap: " << trtPosEC.layers.size() << " layers, " << trtPosEC.toString()); - - - auto makeTVol = [&](const auto& vConf, const auto& name) { - return cvh.createTrackingVolume(gctx, - vConf.layers, - {}, + ATH_MSG_VERBOSE("- TRT::NegativeEndcap: " << trtNegEC.layers.size() + << " layers, " + << trtNegEC.toString()); + ATH_MSG_VERBOSE("- TRT::Barrel: " << trtBrl.layers.size() << " layers, " + << trtBrl.toString()); + ATH_MSG_VERBOSE("- TRT::PositiveEncap: " << trtPosEC.layers.size() + << " layers, " + << trtPosEC.toString()); + + auto makeTVol = [&](const auto &vConf, const auto &name) { + return cvh.createTrackingVolume(gctx, vConf.layers, {}, nullptr, // no material - vConf.rMin, - vConf.rMax, - vConf.zMin, - vConf.zMax, - name - ); + vConf.rMin, vConf.rMax, vConf.zMin, + vConf.zMax, name); }; // now turn them into actual TrackingVolumes auto tvSctNegEC = makeTVol(sctNegEC, "SCT::NegativeEndcap"); - auto tvSctBrl = makeTVol(sctBrl, "SCT::Barrel"); + auto tvSctBrl = makeTVol(sctBrl, "SCT::Barrel"); auto tvSctPosEC = makeTVol(sctPosEC, "SCT::PositiveEndcap"); auto tvTrtNegEC = makeTVol(trtNegEC, "TRT::NegativeEndcap"); - auto tvTrtBrl = makeTVol(trtBrl, "TRT::Barrel"); + auto tvTrtBrl = makeTVol(trtBrl, "TRT::Barrel"); auto tvTrtPosEC = makeTVol(trtPosEC, "TRT::PositiveEndcap"); // combine the endcaps and the barrels, respetively - auto negEC = cvh.createContainerTrackingVolume(gctx, {tvSctNegEC, tvTrtNegEC}); - auto posEC = cvh.createContainerTrackingVolume(gctx, {tvSctPosEC, tvTrtPosEC}); + auto negEC = + cvh.createContainerTrackingVolume(gctx, {tvSctNegEC, tvTrtNegEC}); + auto posEC = + cvh.createContainerTrackingVolume(gctx, {tvSctPosEC, tvTrtPosEC}); auto barrel = cvh.createContainerTrackingVolume(gctx, {tvSctBrl, tvTrtBrl}); // and now combine all of those into one container for the assembly - auto container = cvh.createContainerTrackingVolume(gctx, {negEC, barrel, posEC}); - - // if pixel is present, add positive and negative gap volumes so we can wrap it all - if(pixel) { - auto containerBounds = dynamic_cast<const Acts::CylinderVolumeBounds*>(&container->volumeBounds()); - auto pixelBounds = dynamic_cast<const Acts::CylinderVolumeBounds*>(&pixel->volumeBounds()); - std::vector<std::shared_ptr<Acts::TrackingVolume>> noVolumes; - - auto posGap = cvh.createGapTrackingVolume(gctx, - noVolumes, - nullptr, // no material, - pixelBounds->get(CVBBV::eMinR), - pixelBounds->get(CVBBV::eMaxR), - pixelBounds->get(CVBBV::eHalfLengthZ), - containerBounds->get(CVBBV::eHalfLengthZ), - 0, // material layers, - true, // cylinder - "Pixel::PositiveGap"); - auto negGap = cvh.createGapTrackingVolume(gctx, - noVolumes, - nullptr, // no material, - pixelBounds->get(CVBBV::eMinR), - pixelBounds->get(CVBBV::eMaxR), - -containerBounds->get(CVBBV::eHalfLengthZ), - -pixelBounds->get(CVBBV::eHalfLengthZ), - 0, // material layers, - true, // cylinder - "Pixel::NegativeGap"); - - auto pixelContainer = cvh.createContainerTrackingVolume(gctx, {negGap, pixel, posGap}); + auto container = + cvh.createContainerTrackingVolume(gctx, {negEC, barrel, posEC}); + + // if pixel is present, add positive and negative gap volumes so we can wrap + // it all + if (pixel) { + auto containerBounds = dynamic_cast<const Acts::CylinderVolumeBounds *>( + &container->volumeBounds()); + auto pixelBounds = dynamic_cast<const Acts::CylinderVolumeBounds *>( + &pixel->volumeBounds()); + std::vector<std::shared_ptr<Acts::TrackingVolume>> noVolumes; + + auto posGap = cvh.createGapTrackingVolume( + gctx, noVolumes, + nullptr, // no material, + pixelBounds->get(CVBBV::eMinR), pixelBounds->get(CVBBV::eMaxR), + pixelBounds->get(CVBBV::eHalfLengthZ), + containerBounds->get(CVBBV::eHalfLengthZ), + 0, // material layers, + true, // cylinder + "Pixel::PositiveGap"); + auto negGap = cvh.createGapTrackingVolume( + gctx, noVolumes, + nullptr, // no material, + pixelBounds->get(CVBBV::eMinR), pixelBounds->get(CVBBV::eMaxR), + -containerBounds->get(CVBBV::eHalfLengthZ), + -pixelBounds->get(CVBBV::eHalfLengthZ), + 0, // material layers, + true, // cylinder + "Pixel::NegativeGap"); + + auto pixelContainer = + cvh.createContainerTrackingVolume(gctx, {negGap, pixel, posGap}); // and now create one container that contains Pixel+SCT+TRT - container = cvh.createContainerTrackingVolume(gctx, {pixelContainer, container}); + container = + cvh.createContainerTrackingVolume(gctx, {pixelContainer, container}); } - return container; } -void -ActsTrackingGeometrySvc::populateAlignmentStore(ActsAlignmentStore *store) const -{ - // populate the alignment store with all detector elements - m_trackingGeometry->visitSurfaces( - [store](const Acts::Surface* srf) { - const Acts::DetectorElementBase* detElem = srf->associatedDetectorElement(); - const auto* gmde = dynamic_cast<const ActsDetectorElement*>(detElem); +void ActsTrackingGeometrySvc::populateAlignmentStore( + ActsAlignmentStore *store) const { + ATH_MSG_DEBUG("Populate the alignment store with all detector elements"); + size_t nElements = 0; + m_trackingGeometry->visitSurfaces([store, + &nElements](const Acts::Surface *srf) { + const Acts::DetectorElementBase *detElem = srf->associatedDetectorElement(); + const auto *gmde = dynamic_cast<const ActsDetectorElement *>(detElem); gmde->storeTransform(store); + nElements++; }); + ATH_MSG_DEBUG("Populated with " << nElements << " elements"); } -const ActsAlignmentStore* -ActsTrackingGeometrySvc::getNominalAlignmentStore() const -{ +const ActsAlignmentStore * +ActsTrackingGeometrySvc::getNominalAlignmentStore() const { return m_nominalAlignmentStore.get(); } diff --git a/Tracking/Acts/ActsGeometry/src/ActsWriteTrackingGeometryTransforms.cxx b/Tracking/Acts/ActsGeometry/src/ActsWriteTrackingGeometryTransforms.cxx index 5dd1dce5c8300a3b1de7da7957efb6ae89dc01fb..ae94ba7e1cab2a3ba698a5cd442c10811ca1a121 100644 --- a/Tracking/Acts/ActsGeometry/src/ActsWriteTrackingGeometryTransforms.cxx +++ b/Tracking/Acts/ActsGeometry/src/ActsWriteTrackingGeometryTransforms.cxx @@ -55,10 +55,15 @@ StatusCode ActsWriteTrackingGeometryTransforms::execute() { const Acts::DetectorElementBase *detElem = srf->associatedDetectorElement(); const auto *gmde = static_cast<const ActsDetectorElement *>(detElem); - if(gmde->getSubdetector() == ActsDetectorElement::Subdetector::TRT) { + if(dynamic_cast<const InDetDD::TRT_BaseElement*>(gmde->upstreamDetectorElement()) != nullptr) { return; } + const auto side = dynamic_cast<const InDetDD::SiDetectorElement*>(gmde->upstreamDetectorElement()); + if(side == nullptr) { + throw std::runtime_error{"Not TRT but not Si either"}; // this shouldn't happen + } + gid geoID = srf->geometryId(); os << geoID.volume() << ";"; os << geoID.boundary() << ";"; @@ -68,10 +73,10 @@ StatusCode ActsWriteTrackingGeometryTransforms::execute() { os << ctx.eventID().event_number() << ";"; - if(gmde->getSubdetector() == ActsDetectorElement::Subdetector::Pixel) { + if(side->isPixel()) { os << 0; } - else if(gmde->getSubdetector() == ActsDetectorElement::Subdetector::SCT) { + else if(side->isSCT()) { os << 1; } diff --git a/Tracking/Acts/ActsTrkFitting/src/ActsKalmanFitter.cxx b/Tracking/Acts/ActsTrkFitting/src/ActsKalmanFitter.cxx index bc46108a10aca0b349a419cad516f9dd9f8cfb9c..c038f57fe88477a2887b848140e621d0c2583c2a 100755 --- a/Tracking/Acts/ActsTrkFitting/src/ActsKalmanFitter.cxx +++ b/Tracking/Acts/ActsTrkFitting/src/ActsKalmanFitter.cxx @@ -429,7 +429,9 @@ ActsKalmanFitter::makeTrack(const EventContext& ctx, Acts::GeometryContext& tgCo auto flag = state.typeFlags(); if (state.referenceSurface().associatedDetectorElement() != nullptr) { const auto* actsElement = dynamic_cast<const ActsDetectorElement*>(state.referenceSurface().associatedDetectorElement()); - if (actsElement != nullptr && actsElement->getSubdetector() != ActsDetectorElement::Subdetector::TRT){ + if (actsElement != nullptr + && dynamic_cast<const InDetDD::TRT_BaseElement*>(actsElement->upstreamDetectorElement()) == nullptr) { + const auto* detElem = dynamic_cast<const InDetDD::SiDetectorElement*>(actsElement->upstreamDetectorElement()); // We need to determine the type of state std::bitset<Trk::TrackStateOnSurface::NumberOfTrackStateOnSurfaceTypes> typePattern; const Trk::TrackParameters *parm; @@ -444,10 +446,10 @@ ActsKalmanFitter::makeTrack(const EventContext& ctx, Acts::GeometryContext& tgCo // Check if this is a hole, a dead sensors or a state outside the sensor boundary if(boundaryCheck == Trk::BoundaryCheckResult::DeadElement){ - if (actsElement->getSubdetector() == ActsDetectorElement::Subdetector::Pixel){ + if (detElem->isPixel()) { ++numberOfDeadPixel; } - else if (actsElement->getSubdetector() == ActsDetectorElement::Subdetector::SCT){ + else if (detElem->isSCT()) { ++numberOfDeadSCT; } // Dead sensors states are not stored