add ctest to check if the edges of all surfaces are inside layer boundaries

Tracking/Acts/FaserActsGeometry/CMakeLists.txt

@@ -56,6 +56,8 @@ atlas_add_component( FaserActsGeometry
                          src/FaserActsAlignmentCondAlg.cxx
                          src/NominalAlignmentCondAlg.cxx
                          src/FaserActsVolumeMappingTool.cxx
+                         src/FaserActsGeometryBoundaryTestAlg.h
+                         src/FaserActsGeometryBoundaryTestAlg.cxx
                          src/components/*.cxx
                          PUBLIC_HEADERS FaserActsGeometry
                          INCLUDE_DIRS ${CLHEP_INCLUDE_DIRS} ${EIGEN_INCLUDE_DIRS} ${BOOST_INCLUDE_DIRS}
@@ -84,3 +86,7 @@ atlas_install_headers( FaserActsGeometry )
 atlas_install_python_modules( python/*.py )
 atlas_install_scripts( test/*.py )
 
+atlas_add_test( FaserActsGeometryBoundary_test
+                SCRIPT python ${CMAKE_CURRENT_SOURCE_DIR}/test/FaserActsGeometryBoundary_test.py
+                PROPERTIES WORKING_DIRECTORY ${CMAKE_BINARY_DIR}
+                PROPERTIES TIMEOUT 300 )

Tracking/Acts/FaserActsGeometry/src/FaserActsGeometryBoundaryTestAlg.cxx

+/*
+  Copyright (C) 2002-2023 CERN for the benefit of the ATLAS and FASER collaborations
+*/
+
+
+#include "FaserActsGeometryBoundaryTestAlg.h"
+#include "FaserActsGeometry/FaserActsGeometryContext.h"
+#include "Acts/Surfaces/Surface.hpp"
+#include "Acts/Surfaces/PlaneSurface.hpp"
+#include "Acts/Surfaces/RectangleBounds.hpp"
+#include "Acts/Geometry/TrackingGeometry.hpp"
+#include "Acts/Geometry/DetectorElementBase.hpp"
+#include "Acts/Geometry/TrackingVolume.hpp"
+
+
+FaserActsGeometryBoundaryTestAlg::FaserActsGeometryBoundaryTestAlg(const std::string &name, ISvcLocator *pSvcLocator)
+    : AthReentrantAlgorithm(name, pSvcLocator) {}
+
+
+StatusCode FaserActsGeometryBoundaryTestAlg::initialize() {
+  ATH_CHECK(m_trackingGeometryTool.retrieve());
+  return StatusCode::SUCCESS;
+}
+
+
+StatusCode FaserActsGeometryBoundaryTestAlg::execute(const EventContext &ctx) const {
+
+  std::shared_ptr<const Acts::TrackingGeometry> trackingGeometry
+      = m_trackingGeometryTool->trackingGeometry();
+
+  const Acts::GeometryContext gctx =
+      m_trackingGeometryTool->getGeometryContext(ctx).context();
+
+  // loop over all tracking layer and check if contained detector elements are insie boundaries
+  // return StatusCode::Failure so that test fails, if the edge of any detector element is outside boundaries
+  const Acts::TrackingVolume *trackingVolume = trackingGeometry->highestTrackingVolume();
+  for (auto volume : trackingVolume->confinedVolumes()->arrayObjects()) {
+    for (const auto& layer : volume->confinedLayers()->arrayObjects()) {
+      if (layer->layerType() == Acts::LayerType::active) {
+        // get inner and outer boundaries from approach descriptor:
+        // - innerBoundary is approach suface with minimum z position
+        // - outerBoundary is approach suface with maximum z position
+        Acts::SurfaceVector approachSurfaces = layer->approachDescriptor()->containedSurfaces();
+        const Acts::Surface *innerApproachSurface = *std::min_element(
+            approachSurfaces.begin(), approachSurfaces.end(),
+            [&gctx](const auto &lhs, const auto &rhs) {
+                return lhs->center(gctx).z() < rhs->center(gctx).z();
+            }
+        );
+        const Acts::Surface *outerApproachSurface = *std::max_element(
+            approachSurfaces.begin(), approachSurfaces.end(),
+            [&gctx](const auto &lhs, const auto &rhs) {
+                return lhs->center(gctx).z() < rhs->center(gctx).z();
+            }
+        );
+        double zInnerBoundary = innerApproachSurface->center(gctx).z();
+        double zOuterBoundary = outerApproachSurface->center(gctx).z();
+
+        // loop over surface array and check if all edges are between inner and outer boundary
+        for (const Acts::Surface *surface : layer->surfaceArray()->surfaces()) {
+          auto planeSurface = dynamic_cast<const Acts::PlaneSurface *>(surface);
+          // make sure surface has a associated detector element (there are other active detector elements 
+          // e.g. containing material or at the tracking volume boundaries)
+          if (surface->associatedDetectorElement() != nullptr) {
+            auto bounds = dynamic_cast<const Acts::RectangleBounds*>(&surface->bounds());
+            const Acts::Vector2 min = bounds->min();
+            const Acts::Vector2 max = bounds->max();
+            // create dummpy momentum vector: local to global transformation requires  momentum vector,
+            // which is ignored for PlaneSurface
+            Acts::Vector3 dummyMomentum {1., 1., 1.};
+            // get global position at all edges of the surface
+            std::vector<Acts::Vector3> edges {};
+            edges.push_back(planeSurface->localToGlobal(gctx, {min.x(), min.y()}, dummyMomentum));
+            edges.push_back(planeSurface->localToGlobal(gctx, {min.x(), max.y()}, dummyMomentum));
+            edges.push_back(planeSurface->localToGlobal(gctx, {max.x(), min.y()}, dummyMomentum));
+            edges.push_back(planeSurface->localToGlobal(gctx, {max.x(), max.y()}, dummyMomentum));
+            for (const Acts::Vector3 &edgePosition : edges) {
+              if ((edgePosition.z() < zInnerBoundary) || (edgePosition.z() > zOuterBoundary)) {
+                std::cout << "?? surface outside boundaries\n";
+                std::cout << "inner Boundary: " << zInnerBoundary << std::endl;
+                std::cout << "outer Boundary: " << zOuterBoundary << std::endl;
+                std::cout << "edge: " << edgePosition.x() << ", " << edgePosition.y() << ", " << edgePosition.z() << std::endl;
+                return StatusCode::FAILURE;
+              }
+            }
+          }
+        }
+      }
+    }
+  }
+
+  return StatusCode::SUCCESS;
+}
+
+StatusCode FaserActsGeometryBoundaryTestAlg::finalize() {
+  return StatusCode::SUCCESS;
+}

Tracking/Acts/FaserActsGeometry/src/FaserActsGeometryBoundaryTestAlg.h

+/*
+  Copyright (C) 2002-2023 CERN for the benefit of the ATLAS and FASER collaborations
+*/
+
+#ifndef FASERACTSGEOMETRY_FASERACTSGEOMETRYBOUNDARYTESTALG_H
+#define FASERACTSGEOMETRY_FASERACTSGEOMETRYBOUNDARYTESTALG_H
+
+#include "AthenaBaseComps/AthReentrantAlgorithm.h"
+#include "FaserActsGeometryInterfaces/IFaserActsTrackingGeometryTool.h"
+
+
+class FaserActsGeometryBoundaryTestAlg : public AthReentrantAlgorithm {
+public:
+  FaserActsGeometryBoundaryTestAlg(const std::string &name, ISvcLocator *pSvcLocator);
+  virtual StatusCode initialize() override final;
+  virtual StatusCode execute(const EventContext &ctx) const override final;
+  virtual StatusCode finalize() override final;
+
+private:
+  ToolHandle<IFaserActsTrackingGeometryTool> m_trackingGeometryTool {
+    this, "TrackingGeometryTool", "FaserActsTrackingGeometryTool"};
+};
+
+
+#endif // FASERACTSGEOMETRY_FASERACTSGEOMETRYBOUNDARYTESTALG_H
\ No newline at end of file

Tracking/Acts/FaserActsGeometry/src/components/FaserActsGeometry_entries.cxx

@@ -16,6 +16,7 @@
 #include "../FaserActsMaterialMapping.h"
 #include "../FaserActsSurfaceMappingTool.h"
 #include "../FaserActsMaterialTrackWriterSvc.h"
+#include "../FaserActsGeometryBoundaryTestAlg.h"
 
 DECLARE_COMPONENT( FaserActsTrackingGeometrySvc )
 DECLARE_COMPONENT( FaserActsTrackingGeometryTool )
@@ -30,3 +31,4 @@ DECLARE_COMPONENT( FaserActsMaterialJsonWriterTool )
 DECLARE_COMPONENT( FaserActsMaterialTrackWriterSvc )
 DECLARE_COMPONENT( FaserActsMaterialMapping )
 DECLARE_COMPONENT( FaserActsSurfaceMappingTool )
+DECLARE_COMPONENT( FaserActsGeometryBoundaryTestAlg )

Tracking/Acts/FaserActsGeometry/test/FaserActsGeometryBoundary_test.py

+#!/usr/bin/env python
+"""
+Copyright (C) 2002-2023 CERN for the benefit of the ATLAS and FASER collaboration
+"""
+
+
+def FaserActsGeometryBoundaryTestCfg(flags, name="FaserActsGeometryBoundaryTestAlg", **kwargs):
+
+    from FaserSCT_GeoModel.FaserSCT_GeoModelConfig import FaserSCT_GeometryCfg
+    from MagFieldServices.MagFieldServicesConfig import MagneticFieldSvcCfg
+    from FaserActsGeometry.ActsGeometryConfig import ActsTrackingGeometryToolCfg
+    from AthenaConfiguration.ComponentFactory import CompFactory
+
+    acc = FaserSCT_GeometryCfg(flags)
+    acc.merge(MagneticFieldSvcCfg(flags))
+    result, actsTrackingGeometryTool = ActsTrackingGeometryToolCfg(flags)
+    test_alg = CompFactory.FaserActsGeometryBoundaryTestAlg
+    test_alg.TrackingGeometryTool = actsTrackingGeometryTool
+    acc.merge(result)
+    acc.addEventAlgo(test_alg(name, **kwargs))
+    return acc
+
+
+if __name__ == "__main__":
+    
+    import sys
+    from AthenaCommon.Configurable import Configurable
+    from CalypsoConfiguration.AllConfigFlags import ConfigFlags
+    from CalypsoConfiguration.MainServicesConfig import MainServicesCfg
+    from FaserGeoModel.FaserGeoModelConfig import FaserGeometryCfg
+
+    Configurable.configurableRun3Behavior = True
+    ConfigFlags.Input.isMC = False
+    ConfigFlags.GeoModel.Align.Dynamic = False
+    ConfigFlags.IOVDb.DatabaseInstance = "CONDBR3"
+    ConfigFlags.IOVDb.GlobalTag = "OFLCOND-FASER-04"
+    ConfigFlags.GeoModel.FaserVersion = "FASERNU-03"
+    ConfigFlags.Detector.GeometryFaserSCT = True
+    ConfigFlags.lock()
+
+    acc = MainServicesCfg(ConfigFlags)
+    acc.merge(FaserGeometryCfg(ConfigFlags))
+    acc.merge(FaserActsGeometryBoundaryTestCfg(ConfigFlags))
+
+    replicaSvc = acc.getService("DBReplicaSvc")
+    replicaSvc.COOLSQLiteVetoPattern = ""
+    replicaSvc.UseCOOLSQLite = True
+    replicaSvc.UseCOOLFrontier = False
+    replicaSvc.UseGeomSQLite = True
+
+    sys.exit(int(acc.run(maxEvents=1).isFailure()))