From c3ef0b5f9791f323406ae4547f1c96b8d11ebb6c Mon Sep 17 00:00:00 2001 From: Federico Meloni Date: Wed, 3 Jun 2026 14:57:04 +0200 Subject: [PATCH 01/10] initial attempt at extrapolating to calo --- .../include/k4ActsTracking/Helpers.hxx | 38 +++++ .../include/k4ActsTracking/IActsGeoSvc.h | 18 +++ k4ActsTracking/src/components/ActsGeoSvc.cpp | 140 ++++++++++++++++++ k4ActsTracking/src/components/ActsGeoSvc.h | 12 ++ .../src/components/CKFTrackingAlg.cpp | 62 ++++++++ k4ActsTracking/src/components/Helpers.cxx | 66 +++++++++ 6 files changed, 336 insertions(+) diff --git a/k4ActsTracking/include/k4ActsTracking/Helpers.hxx b/k4ActsTracking/include/k4ActsTracking/Helpers.hxx index 6b391b75..1fe302c6 100644 --- a/k4ActsTracking/include/k4ActsTracking/Helpers.hxx +++ b/k4ActsTracking/include/k4ActsTracking/Helpers.hxx @@ -38,16 +38,31 @@ #include #include #include +#include +#include #include +#include +#include +#include #include #include #include "Acts/EventData/ParticleHypothesis.hpp" // ACTSTracking +#include "k4ActsTracking/IActsGeoSvc.h" #include "k4ActsTracking/SourceLink.hxx" +// Standard +#include + namespace ACTSTracking { + /// Propagator used to extrapolate fitted tracks out to the calorimeter face. + /// It uses a VoidNavigator (no tracking geometry), so it can reach target + /// surfaces that lie outside the tracking-geometry world volume. Material + /// between the tracker and the calorimeter is not accounted for (field-only). + using CaloFacePropagator = Acts::Propagator, Acts::VoidNavigator>; + using TrackResult = Acts::TrackContainer::TrackProxy; @@ -105,4 +120,27 @@ namespace ACTSTracking { */ Acts::ParticleHypothesis convertParticle(const edm4hep::MCParticle mcParticle); + //! Extrapolate track parameters to the calorimeter face. + /** + * Selects, among the calorimeter-face surfaces, the one the track reaches first + * (smallest positive path length with a valid, bounds-checked intersection) and + * propagates the given parameters to it. The barrel is a ring of planar + * surfaces (one per polygon side) and the endcaps are flat discs; intersecting + * all candidates naturally handles both the barrel-sector choice and the + * barrel/endcap transition. + * + * \param propagator Field-only propagator (VoidNavigator) + * \param start Track parameters to extrapolate from (e.g. at the last hit) + * \param surfaces Calorimeter-face surfaces from IActsGeoSvc + * \param gctx Geometry context + * \param mctx Magnetic-field context + * + * \return Bound track parameters at the calorimeter face, or std::nullopt if no + * surface is reached or the propagation fails. + */ + std::optional extrapolateToCaloFace( + const CaloFacePropagator& propagator, const Acts::BoundTrackParameters& start, + const IActsGeoSvc::CaloFaceSurfaces& surfaces, const Acts::GeometryContext& gctx, + const Acts::MagneticFieldContext& mctx); + } // namespace ACTSTracking diff --git a/k4ActsTracking/include/k4ActsTracking/IActsGeoSvc.h b/k4ActsTracking/include/k4ActsTracking/IActsGeoSvc.h index 63a44626..0baa0cb6 100644 --- a/k4ActsTracking/include/k4ActsTracking/IActsGeoSvc.h +++ b/k4ActsTracking/include/k4ActsTracking/IActsGeoSvc.h @@ -27,6 +27,7 @@ #include #include #include +#include namespace dd4hep { namespace rec { @@ -44,6 +45,22 @@ class GAUDI_API IActsGeoSvc : virtual public IService { public: using CellIDSurfaceMap = std::unordered_map; + /// Surfaces approximating the inner face of the electromagnetic calorimeter, + /// derived from the DD4hep geometry. These live outside the tracking-geometry + /// world volume and are meant as target surfaces for a geometry-free + /// extrapolation (see ACTSTracking::extrapolateToCaloFace). + /// + /// The barrel is a regular polygon, so its inner face is modelled as one + /// planar surface per polygon side rather than a cylinder. The endcaps are + /// flat discs at constant z. + struct CaloFaceSurfaces { + std::vector> barrelFaces; ///< one plane per polygon side + std::shared_ptr endcapPos; ///< disc at +z, may be null + std::shared_ptr endcapNeg; ///< disc at -z, may be null + + bool empty() const { return barrelFaces.empty() && !endcapPos && !endcapNeg; } + }; + public: DeclareInterfaceID(IActsGeoSvc, 1, 0); @@ -51,6 +68,7 @@ class GAUDI_API IActsGeoSvc : virtual public IService { virtual std::shared_ptr magneticField() const = 0; virtual const CellIDSurfaceMap& cellIdToSurfaceMap() const = 0; virtual std::string cellIDEncodingString() const = 0; + virtual const CaloFaceSurfaces& caloFaceSurfaces() const = 0; virtual ~IActsGeoSvc() = default; }; diff --git a/k4ActsTracking/src/components/ActsGeoSvc.cpp b/k4ActsTracking/src/components/ActsGeoSvc.cpp index 97806db3..778777b7 100644 --- a/k4ActsTracking/src/components/ActsGeoSvc.cpp +++ b/k4ActsTracking/src/components/ActsGeoSvc.cpp @@ -37,6 +37,9 @@ #include #include #include +#include +#include +#include #include #include #include @@ -45,7 +48,9 @@ #include #include +#include #include +#include #include @@ -54,6 +59,8 @@ #include #include +#include +#include template <> struct fmt::formatter : fmt::ostream_formatter {}; @@ -159,5 +166,138 @@ StatusCode ActsGeoSvc::initialize() { vis.write(m_objDumpFileName.value()); } + if (m_buildCaloSurfaces.value()) { + buildCaloFaceSurfaces(); + } + return StatusCode::SUCCESS; } + +void ActsGeoSvc::buildCaloFaceSurfaces() { + using dd4hep::DetType; + using dd4hep::rec::LayeredCalorimeterData; + namespace UC = Acts::UnitConstants; + + // LayeredCalorimeterData::extent[] is stored in DD4hep native length units + // (cm), so convert to ACTS units with the same scale the blueprint uses. + const double lengthScale = Acts::UnitConstants::cm / dd4hep::cm; + + const auto* dd4hepDet = m_geoSvc->getDetector(); + + // Locate the electromagnetic calorimeter sub-detectors purely via DetType + // flags, so that this works across detector concepts without hard-coding + // DetElement names. The dimensions themselves come from the standard DDRec + // LayeredCalorimeterData extension (the same source Pandora uses). + const auto ecalBarrel = + dd4hepDet->detectors(DetType::CALORIMETER | DetType::ELECTROMAGNETIC | DetType::BARREL, DetType::FORWARD); + const auto ecalEndcap = + dd4hepDet->detectors(DetType::CALORIMETER | DetType::ELECTROMAGNETIC | DetType::ENDCAP, 0); + + // Barrel circumradius (corner radius), needed both for the barrel faces and to + // guarantee the endcap discs reach far enough to avoid hermeticity gaps at the + // barrel/endcap junction. + double barrelCircumradius = 0.0; + double barrelHalfZ = 0.0; + + // --- Barrel: one planar surface per polygon side -------------------------- + if (ecalBarrel.empty()) { + warning() << "No electromagnetic barrel calorimeter found via DetType flags; " + "no barrel calo-face surfaces will be built." + << endmsg; + } else { + const auto* caloData = ecalBarrel.front().extension(false); + if (caloData == nullptr) { + warning() << "ECAL barrel DetElement has no LayeredCalorimeterData extension; " + "skipping barrel calo-face surfaces." + << endmsg; + } else { + const int nSides = caloData->inner_symmetry > 0 ? caloData->inner_symmetry : 0; + // For a BarrelLayout the calorimeter is centred on z = 0 and extent[] is + // {rmin, rmax, zmin=0, zmax=half_length}, so the half-length is extent[3]. + const double apothem = caloData->extent[0] * lengthScale; // perpendicular distance to inner face + barrelHalfZ = caloData->extent[3] * lengthScale; // barrel half-length + const double phi0 = caloData->inner_phi0; // azimuth of first inner-face normal + + if (nSides < 3) { + warning() << fmt::format("ECAL barrel has unusable inner_symmetry={}; skipping barrel faces.", + caloData->inner_symmetry) + << endmsg; + } else { + const double dPhi = 2 * std::numbers::pi / nSides; + barrelCircumradius = apothem / std::cos(std::numbers::pi / nSides); + const double halfWidth = apothem * std::tan(std::numbers::pi / nSides); + + m_caloFaceSurfaces.barrelFaces.reserve(nSides); + for (int i = 0; i < nSides; ++i) { + const double phi = phi0 + i * dPhi; + const double cphi = std::cos(phi); + const double sphi = std::sin(phi); + Acts::Vector3 normal{cphi, sphi, 0}; // local z (surface normal, radial) + Acts::Vector3 localX{-sphi, cphi, 0}; // tangential + Acts::Vector3 localY{0, 0, 1}; // along global z + Acts::Vector3 center = apothem * normal; // barrel centred on z = 0 + + Acts::Transform3 transform = Acts::Transform3::Identity(); + transform.linear().col(0) = localX; + transform.linear().col(1) = localY; + transform.linear().col(2) = normal; + transform.translation() = center; + + auto bounds = std::make_shared(halfWidth, barrelHalfZ); + m_caloFaceSurfaces.barrelFaces.push_back(Acts::Surface::makeShared(transform, bounds)); + } + info() << fmt::format( + "Built ECAL barrel calo face: {} planar faces, apothem={:.1f} mm, circumradius={:.1f} mm, " + "halfZ={:.1f} mm, phi0={:.4f}", + nSides, apothem / UC::mm, barrelCircumradius / UC::mm, barrelHalfZ / UC::mm, phi0) + << endmsg; + } + } + } + + // --- Endcaps: flat discs, widened to stay hermetic with the barrel -------- + if (ecalEndcap.empty()) { + warning() << "No electromagnetic endcap calorimeter found via DetType flags; " + "no endcap calo-face surfaces will be built." + << endmsg; + } else { + const auto* caloData = ecalEndcap.front().extension(false); + if (caloData == nullptr) { + warning() << "ECAL endcap DetElement has no LayeredCalorimeterData extension; " + "skipping endcap calo-face surfaces." + << endmsg; + } else { + const double rMin = caloData->extent[0] * lengthScale; + double rMax = caloData->extent[1] * lengthScale; + const double zEndcap = caloData->extent[2] * lengthScale; // inner-face z (zmin) + + // Guarantee the disc reaches at least the barrel corner radius so there is + // no gap at the barrel/endcap junction. + if (barrelCircumradius > rMax) { + warning() << fmt::format( + "ECAL endcap rMax ({:.1f} mm) is smaller than the barrel circumradius ({:.1f} mm); " + "extending the endcap disc to close the hermeticity gap.", + rMax / UC::mm, barrelCircumradius / UC::mm) + << endmsg; + rMax = barrelCircumradius; + } + + Acts::Transform3 tPos = Acts::Transform3::Identity(); + tPos.translation() = Acts::Vector3{0, 0, zEndcap}; + Acts::Transform3 tNeg = Acts::Transform3::Identity(); + tNeg.translation() = Acts::Vector3{0, 0, -zEndcap}; + + m_caloFaceSurfaces.endcapPos = Acts::Surface::makeShared(tPos, rMin, rMax); + m_caloFaceSurfaces.endcapNeg = Acts::Surface::makeShared(tNeg, rMin, rMax); + + info() << fmt::format("Built ECAL endcap calo faces: discs at z=+/-{:.1f} mm, rMin={:.1f} mm, rMax={:.1f} mm", + zEndcap / UC::mm, rMin / UC::mm, rMax / UC::mm) + << endmsg; + } + } + + info() << fmt::format("Calo-face surfaces built: {} barrel faces, {} endcap discs", + m_caloFaceSurfaces.barrelFaces.size(), + (m_caloFaceSurfaces.endcapPos ? 1 : 0) + (m_caloFaceSurfaces.endcapNeg ? 1 : 0)) + << endmsg; +} diff --git a/k4ActsTracking/src/components/ActsGeoSvc.h b/k4ActsTracking/src/components/ActsGeoSvc.h index fcc62d46..26bf40a5 100644 --- a/k4ActsTracking/src/components/ActsGeoSvc.h +++ b/k4ActsTracking/src/components/ActsGeoSvc.h @@ -53,6 +53,8 @@ class ActsGeoSvc : public extends { std::shared_ptr magneticField() const override; + const CaloFaceSurfaces& caloFaceSurfaces() const override { return m_caloFaceSurfaces; } + ActsGeoSvc(const std::string& name, ISvcLocator* svcLoc); ~ActsGeoSvc() = default; @@ -66,6 +68,9 @@ class ActsGeoSvc : public extends { Gaudi::Property m_encodingStringConstant{ this, "EncodingStringVariable", "GlobalTrackerReadoutID", "Name of the DD4hep constant holding the CellID encoding string."}; + Gaudi::Property m_buildCaloSurfaces{ + this, "BuildCaloSurfaces", true, + "Whether to build the ECAL inner-face surfaces (for track extrapolation to the calorimeter face)."}; const CellIDSurfaceMap& cellIdToSurfaceMap() const override { return m_cellIDToSurface; } std::string cellIDEncodingString() const override { return m_cellIDEncodingString; } @@ -75,12 +80,19 @@ class ActsGeoSvc : public extends { using BlueprintPopulationFunc = void(const std::string&, Acts::Experimental::Blueprint&, BlueprintBuilder&); + /// Build the ECAL inner-face surfaces (m_caloFaceSurfaces) from the DD4hep + /// geometry. Surfaces are located via DetType flags and dimensioned from the + /// dd4hep::rec::LayeredCalorimeterData extension. Missing ECAL sub-detectors + /// or extensions are skipped with a warning rather than treated as an error. + void buildCaloFaceSurfaces(); + SmartIF m_geoSvc; std::shared_ptr m_trackingGeo{nullptr}; std::shared_ptr m_magneticField{nullptr}; std::unordered_map m_cellIDToSurface{}; std::unordered_map m_bluePrintPopulationFuncs{}; std::string m_cellIDEncodingString{}; + CaloFaceSurfaces m_caloFaceSurfaces{}; }; inline std::shared_ptr ActsGeoSvc::trackingGeometry() const { return m_trackingGeo; } diff --git a/k4ActsTracking/src/components/CKFTrackingAlg.cpp b/k4ActsTracking/src/components/CKFTrackingAlg.cpp index 55b765ff..98a08d4b 100644 --- a/k4ActsTracking/src/components/CKFTrackingAlg.cpp +++ b/k4ActsTracking/src/components/CKFTrackingAlg.cpp @@ -58,6 +58,7 @@ #include #include #include +#include #include #include #include @@ -143,6 +144,9 @@ struct CKFTrackingAlg final ///@{ Gaudi::Property m_runCKF{this, "RunCKF", true, "Run tracking using CKF. False means stop at seeding."}; Gaudi::Property m_propagateBackward{this, "PropagateBackward", false, "Extrapolates tracks towards beamline."}; + Gaudi::Property m_extrapolateToCalo{ + this, "ExtrapolateToCalo", true, + "Extrapolate fitted tracks to the calorimeter face and add an AtCalorimeter track state."}; ///@} /// @name Seed-finding configuration @@ -219,6 +223,11 @@ struct CKFTrackingAlg final // once in initialize() and reused (read-only) across events and threads. std::optional m_trackFinder{}; + // Field-only propagator (no tracking geometry) used to extrapolate fitted + // tracks out to the calorimeter face, which lies outside the tracking + // geometry. Built once in initialize() and reused read-only across threads. + std::optional m_caloPropagator{}; + k4ActsTracking::CellIDSelector m_seedSelector{}; mutable std::mutex m_seedMutex{}; @@ -264,6 +273,20 @@ StatusCode CKFTrackingAlg::initialize() { Propagator propagator(std::move(stepper), std::move(navigator)); m_trackFinder.emplace(std::move(propagator)); + // The calorimeter face lies outside the tracking geometry, so extrapolation + // there uses a geometry-free (VoidNavigator) propagator. Only build it when + // requested and when the geometry service actually provides calo surfaces. + if (m_extrapolateToCalo) { + if (m_actsGeoSvc->caloFaceSurfaces().empty()) { + warning() << "ExtrapolateToCalo requested but ActsGeoSvc provides no calorimeter-face surfaces; " + "no AtCalorimeter track states will be produced." + << endmsg; + } else { + Acts::EigenStepper<> caloStepper(m_actsGeoSvc->magneticField()); + m_caloPropagator.emplace(std::move(caloStepper), Acts::VoidNavigator{}); + } + } + return StatusCode::SUCCESS; } @@ -640,6 +663,45 @@ StatusCode CKFTrackingAlg::tracking(const std::vectormagneticField(), magCache); + + // Extrapolate the fitted track to the calorimeter face and add an + // AtCalorimeter track state for Pandora / ParticleFlow. Starts from the + // outermost smoothed state (closest to the calorimeter). + if (m_caloPropagator) { + std::optional startParams; + for (const auto& state : trackTip.trackStatesReversed()) { + if (state.hasSmoothed()) { + startParams.emplace(state.referenceSurface().getSharedPtr(), state.smoothed(), + state.smoothedCovariance(), trackTip.particleHypothesis()); + break; + } + } + + if (startParams) { + const Acts::MagneticFieldContext magCtx{}; + auto caloParams = ACTSTracking::extrapolateToCaloFace(*m_caloPropagator, *startParams, + m_actsGeoSvc->caloFaceSurfaces(), geoCtx, magCtx); + if (caloParams) { + const Acts::Vector3 caloPos = caloParams->position(geoCtx); + auto fieldRes = m_actsGeoSvc->magneticField()->getField(caloPos, magCache); + const double Bz = fieldRes.ok() ? (*fieldRes)[2] / Acts::UnitConstants::T : 0.0; + auto caloState = ACTSTracking::ACTS2edm4hep_trackState(edm4hep::TrackState::AtCalorimeter, *caloParams, Bz); + // The calo-face parameters are local to the target surface, so the + // edm4hep D0/Z0 from the generic conversion are not meaningful here. + // Express the state at the impact point instead: set the reference + // point to the global calo-face position (D0 = Z0 = 0 there). + caloState.referencePoint = edm4hep::Vector3f(caloPos.x(), caloPos.y(), caloPos.z()); + caloState.D0 = 0; + caloState.Z0 = 0; + track.addToTrackStates(caloState); + } else { + debug() << "Extrapolation to the calorimeter face did not reach a surface; " + "no AtCalorimeter state added for this track." + << endmsg; + } + } + } + { std::lock_guard lock{m_trackMutex}; trackCollection.push_back(track); diff --git a/k4ActsTracking/src/components/Helpers.cxx b/k4ActsTracking/src/components/Helpers.cxx index 16acd5de..e675f460 100644 --- a/k4ActsTracking/src/components/Helpers.cxx +++ b/k4ActsTracking/src/components/Helpers.cxx @@ -31,6 +31,12 @@ // ACTS #include #include +#include +#include +#include +#include + +#include // ACTSTracking #include "config.h.in" @@ -234,4 +240,64 @@ namespace ACTSTracking { return Acts::ParticleHypothesis{pdg, mass, charge_type}; } + std::optional extrapolateToCaloFace( + const CaloFacePropagator& propagator, const Acts::BoundTrackParameters& start, + const IActsGeoSvc::CaloFaceSurfaces& surfaces, const Acts::GeometryContext& gctx, + const Acts::MagneticFieldContext& mctx) { + if (surfaces.empty()) { + return std::nullopt; + } + + const Acts::Vector3 position = start.position(gctx); + const Acts::Vector3 direction = start.direction(); + + // Allow a small slack at face edges / the barrel-endcap seam so tracks + // crossing right at a boundary are not lost. + constexpr double tolerance = 1.0 * Acts::UnitConstants::mm; + const Acts::BoundaryTolerance boundaryTolerance = Acts::BoundaryTolerance::AbsoluteEuclidean(tolerance); + + // Build the candidate list: every barrel face plus the endcap disc on the + // side the track is heading towards. + std::vector candidates; + candidates.reserve(surfaces.barrelFaces.size() + 1); + for (const auto& face : surfaces.barrelFaces) { + candidates.push_back(face.get()); + } + const auto& endcap = (direction.z() >= 0) ? surfaces.endcapPos : surfaces.endcapNeg; + if (endcap) { + candidates.push_back(endcap.get()); + } + + // Pick the surface reached first along the track direction. + const Acts::Surface* target = nullptr; + double bestPath = std::numeric_limits::max(); + constexpr double minPath = 1e-3; // ignore intersections essentially at the start point + for (const Acts::Surface* surface : candidates) { + const auto multiIntersection = surface->intersect(gctx, position, direction, boundaryTolerance); + for (const auto& intersection : multiIntersection) { + if (!intersection.isValid()) { + continue; + } + const double path = intersection.pathLength(); + if (path > minPath && path < bestPath) { + bestPath = path; + target = surface; + } + } + } + + if (target == nullptr) { + return std::nullopt; + } + + Acts::PropagatorPlainOptions options{gctx, mctx}; + options.maxSteps = 10000; + + auto result = propagator.propagateToSurface(start, *target, options); + if (!result.ok()) { + return std::nullopt; + } + return result.value(); + } + } // namespace ACTSTracking From 13a7682eaf2c96e950f80ba5be40266603d31ad8 Mon Sep 17 00:00:00 2001 From: Federico Meloni Date: Wed, 3 Jun 2026 15:04:28 +0200 Subject: [PATCH 02/10] Apply clang-format to calo-extrapolation sources Fixes the pre-commit clang-format CI failure on PR #64. Co-Authored-By: Claude Opus 4.8 --- .../include/k4ActsTracking/Helpers.hxx | 9 ++++--- k4ActsTracking/src/components/ActsGeoSvc.cpp | 25 +++++++++---------- .../src/components/CKFTrackingAlg.cpp | 17 +++++++------ k4ActsTracking/src/components/Helpers.cxx | 19 +++++++------- 4 files changed, 36 insertions(+), 34 deletions(-) diff --git a/k4ActsTracking/include/k4ActsTracking/Helpers.hxx b/k4ActsTracking/include/k4ActsTracking/Helpers.hxx index 1fe302c6..24386b6e 100644 --- a/k4ActsTracking/include/k4ActsTracking/Helpers.hxx +++ b/k4ActsTracking/include/k4ActsTracking/Helpers.hxx @@ -138,9 +138,10 @@ namespace ACTSTracking { * \return Bound track parameters at the calorimeter face, or std::nullopt if no * surface is reached or the propagation fails. */ - std::optional extrapolateToCaloFace( - const CaloFacePropagator& propagator, const Acts::BoundTrackParameters& start, - const IActsGeoSvc::CaloFaceSurfaces& surfaces, const Acts::GeometryContext& gctx, - const Acts::MagneticFieldContext& mctx); + std::optional extrapolateToCaloFace(const CaloFacePropagator& propagator, + const Acts::BoundTrackParameters& start, + const IActsGeoSvc::CaloFaceSurfaces& surfaces, + const Acts::GeometryContext& gctx, + const Acts::MagneticFieldContext& mctx); } // namespace ACTSTracking diff --git a/k4ActsTracking/src/components/ActsGeoSvc.cpp b/k4ActsTracking/src/components/ActsGeoSvc.cpp index 778777b7..dea36c26 100644 --- a/k4ActsTracking/src/components/ActsGeoSvc.cpp +++ b/k4ActsTracking/src/components/ActsGeoSvc.cpp @@ -190,8 +190,7 @@ void ActsGeoSvc::buildCaloFaceSurfaces() { // LayeredCalorimeterData extension (the same source Pandora uses). const auto ecalBarrel = dd4hepDet->detectors(DetType::CALORIMETER | DetType::ELECTROMAGNETIC | DetType::BARREL, DetType::FORWARD); - const auto ecalEndcap = - dd4hepDet->detectors(DetType::CALORIMETER | DetType::ELECTROMAGNETIC | DetType::ENDCAP, 0); + const auto ecalEndcap = dd4hepDet->detectors(DetType::CALORIMETER | DetType::ELECTROMAGNETIC | DetType::ENDCAP, 0); // Barrel circumradius (corner radius), needed both for the barrel faces and to // guarantee the endcap discs reach far enough to avoid hermeticity gaps at the @@ -229,19 +228,19 @@ void ActsGeoSvc::buildCaloFaceSurfaces() { m_caloFaceSurfaces.barrelFaces.reserve(nSides); for (int i = 0; i < nSides; ++i) { - const double phi = phi0 + i * dPhi; - const double cphi = std::cos(phi); - const double sphi = std::sin(phi); - Acts::Vector3 normal{cphi, sphi, 0}; // local z (surface normal, radial) - Acts::Vector3 localX{-sphi, cphi, 0}; // tangential - Acts::Vector3 localY{0, 0, 1}; // along global z - Acts::Vector3 center = apothem * normal; // barrel centred on z = 0 + const double phi = phi0 + i * dPhi; + const double cphi = std::cos(phi); + const double sphi = std::sin(phi); + Acts::Vector3 normal{cphi, sphi, 0}; // local z (surface normal, radial) + Acts::Vector3 localX{-sphi, cphi, 0}; // tangential + Acts::Vector3 localY{0, 0, 1}; // along global z + Acts::Vector3 center = apothem * normal; // barrel centred on z = 0 Acts::Transform3 transform = Acts::Transform3::Identity(); - transform.linear().col(0) = localX; - transform.linear().col(1) = localY; - transform.linear().col(2) = normal; - transform.translation() = center; + transform.linear().col(0) = localX; + transform.linear().col(1) = localY; + transform.linear().col(2) = normal; + transform.translation() = center; auto bounds = std::make_shared(halfWidth, barrelHalfZ); m_caloFaceSurfaces.barrelFaces.push_back(Acts::Surface::makeShared(transform, bounds)); diff --git a/k4ActsTracking/src/components/CKFTrackingAlg.cpp b/k4ActsTracking/src/components/CKFTrackingAlg.cpp index 98a08d4b..3396d2b4 100644 --- a/k4ActsTracking/src/components/CKFTrackingAlg.cpp +++ b/k4ActsTracking/src/components/CKFTrackingAlg.cpp @@ -189,15 +189,15 @@ struct CKFTrackingAlg final /// @name Track-fit initial error estimates ///@{ Gaudi::Property m_initialTrackError_pos{this, "InitialTrackError_Pos", 10 * Acts::UnitConstants::um, - "Initial track error for local position."}; + "Initial track error for local position."}; Gaudi::Property m_initialTrackError_phi{this, "InitialTrackError_Phi", 1 * Acts::UnitConstants::degree, - "Initial track error for phi."}; + "Initial track error for phi."}; Gaudi::Property m_initialTrackError_relP{this, "InitialTrackError_RelP", 0.25, - "Initial track error for momentum (relative)."}; + "Initial track error for momentum (relative)."}; Gaudi::Property m_initialTrackError_lambda{this, "InitialTrackError_Lambda", 1 * Acts::UnitConstants::degree, - "Initial track error for lambda."}; + "Initial track error for lambda."}; Gaudi::Property m_initialTrackError_time{this, "InitialTrackError_Time", 100 * Acts::UnitConstants::ns, - "Initial track error for time."}; + "Initial track error for time."}; Gaudi::Property m_CKF_chi2CutOff{this, "CKF_Chi2CutOff", 15, "Maximum local chi2 contribution."}; Gaudi::Property m_CKF_numMeasurementsCutOff{this, "CKF_NumMeasurementsCutOff", 10, "Maximum measurements on a single surface."}; @@ -671,8 +671,8 @@ StatusCode CKFTrackingAlg::tracking(const std::vector startParams; for (const auto& state : trackTip.trackStatesReversed()) { if (state.hasSmoothed()) { - startParams.emplace(state.referenceSurface().getSharedPtr(), state.smoothed(), - state.smoothedCovariance(), trackTip.particleHypothesis()); + startParams.emplace(state.referenceSurface().getSharedPtr(), state.smoothed(), state.smoothedCovariance(), + trackTip.particleHypothesis()); break; } } @@ -685,7 +685,8 @@ StatusCode CKFTrackingAlg::tracking(const std::vectorposition(geoCtx); auto fieldRes = m_actsGeoSvc->magneticField()->getField(caloPos, magCache); const double Bz = fieldRes.ok() ? (*fieldRes)[2] / Acts::UnitConstants::T : 0.0; - auto caloState = ACTSTracking::ACTS2edm4hep_trackState(edm4hep::TrackState::AtCalorimeter, *caloParams, Bz); + auto caloState = + ACTSTracking::ACTS2edm4hep_trackState(edm4hep::TrackState::AtCalorimeter, *caloParams, Bz); // The calo-face parameters are local to the target surface, so the // edm4hep D0/Z0 from the generic conversion are not meaningful here. // Express the state at the impact point instead: set the reference diff --git a/k4ActsTracking/src/components/Helpers.cxx b/k4ActsTracking/src/components/Helpers.cxx index e675f460..d3e49c4e 100644 --- a/k4ActsTracking/src/components/Helpers.cxx +++ b/k4ActsTracking/src/components/Helpers.cxx @@ -240,10 +240,11 @@ namespace ACTSTracking { return Acts::ParticleHypothesis{pdg, mass, charge_type}; } - std::optional extrapolateToCaloFace( - const CaloFacePropagator& propagator, const Acts::BoundTrackParameters& start, - const IActsGeoSvc::CaloFaceSurfaces& surfaces, const Acts::GeometryContext& gctx, - const Acts::MagneticFieldContext& mctx) { + std::optional extrapolateToCaloFace(const CaloFacePropagator& propagator, + const Acts::BoundTrackParameters& start, + const IActsGeoSvc::CaloFaceSurfaces& surfaces, + const Acts::GeometryContext& gctx, + const Acts::MagneticFieldContext& mctx) { if (surfaces.empty()) { return std::nullopt; } @@ -253,8 +254,8 @@ namespace ACTSTracking { // Allow a small slack at face edges / the barrel-endcap seam so tracks // crossing right at a boundary are not lost. - constexpr double tolerance = 1.0 * Acts::UnitConstants::mm; - const Acts::BoundaryTolerance boundaryTolerance = Acts::BoundaryTolerance::AbsoluteEuclidean(tolerance); + constexpr double tolerance = 1.0 * Acts::UnitConstants::mm; + const Acts::BoundaryTolerance boundaryTolerance = Acts::BoundaryTolerance::AbsoluteEuclidean(tolerance); // Build the candidate list: every barrel face plus the endcap disc on the // side the track is heading towards. @@ -269,9 +270,9 @@ namespace ACTSTracking { } // Pick the surface reached first along the track direction. - const Acts::Surface* target = nullptr; - double bestPath = std::numeric_limits::max(); - constexpr double minPath = 1e-3; // ignore intersections essentially at the start point + const Acts::Surface* target = nullptr; + double bestPath = std::numeric_limits::max(); + constexpr double minPath = 1e-3; // ignore intersections essentially at the start point for (const Acts::Surface* surface : candidates) { const auto multiIntersection = surface->intersect(gctx, position, direction, boundaryTolerance); for (const auto& intersection : multiIntersection) { From 374d36029ca357fe3498e1ab1e1a48bb5c37d3b9 Mon Sep 17 00:00:00 2001 From: Federico Meloni Date: Wed, 3 Jun 2026 15:10:45 +0200 Subject: [PATCH 03/10] Restore initial-track-error alignment to match CI clang-format The earlier clang-format pass used a newer clang-format that realigned the pre-existing m_initialTrackError_* continuation lines; CI's version disagrees. Restore the original alignment. Co-Authored-By: Claude Opus 4.8 --- k4ActsTracking/src/components/CKFTrackingAlg.cpp | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/k4ActsTracking/src/components/CKFTrackingAlg.cpp b/k4ActsTracking/src/components/CKFTrackingAlg.cpp index 3396d2b4..d93e74ce 100644 --- a/k4ActsTracking/src/components/CKFTrackingAlg.cpp +++ b/k4ActsTracking/src/components/CKFTrackingAlg.cpp @@ -189,15 +189,15 @@ struct CKFTrackingAlg final /// @name Track-fit initial error estimates ///@{ Gaudi::Property m_initialTrackError_pos{this, "InitialTrackError_Pos", 10 * Acts::UnitConstants::um, - "Initial track error for local position."}; + "Initial track error for local position."}; Gaudi::Property m_initialTrackError_phi{this, "InitialTrackError_Phi", 1 * Acts::UnitConstants::degree, - "Initial track error for phi."}; + "Initial track error for phi."}; Gaudi::Property m_initialTrackError_relP{this, "InitialTrackError_RelP", 0.25, - "Initial track error for momentum (relative)."}; + "Initial track error for momentum (relative)."}; Gaudi::Property m_initialTrackError_lambda{this, "InitialTrackError_Lambda", 1 * Acts::UnitConstants::degree, - "Initial track error for lambda."}; + "Initial track error for lambda."}; Gaudi::Property m_initialTrackError_time{this, "InitialTrackError_Time", 100 * Acts::UnitConstants::ns, - "Initial track error for time."}; + "Initial track error for time."}; Gaudi::Property m_CKF_chi2CutOff{this, "CKF_Chi2CutOff", 15, "Maximum local chi2 contribution."}; Gaudi::Property m_CKF_numMeasurementsCutOff{this, "CKF_NumMeasurementsCutOff", 10, "Maximum measurements on a single surface."}; From 734f846633d0a35ba8438e6d84184146c8a4f434 Mon Sep 17 00:00:00 2001 From: Federico Meloni Date: Wed, 3 Jun 2026 21:09:38 +0200 Subject: [PATCH 04/10] get track parameters at IP for both new and old --- .../ACTSSeededCKFTrackingAlg.hxx | 6 ++- .../components/ACTSSeededCKFTrackingAlg.cxx | 32 ++++++++++++++-- .../src/components/CKFTrackingAlg.cpp | 38 +++++++++++++++++++ 3 files changed, 70 insertions(+), 6 deletions(-) diff --git a/k4ActsTracking/include/k4ActsTracking/ACTSSeededCKFTrackingAlg.hxx b/k4ActsTracking/include/k4ActsTracking/ACTSSeededCKFTrackingAlg.hxx index 1651b900..d08203c3 100644 --- a/k4ActsTracking/include/k4ActsTracking/ACTSSeededCKFTrackingAlg.hxx +++ b/k4ActsTracking/include/k4ActsTracking/ACTSSeededCKFTrackingAlg.hxx @@ -60,6 +60,7 @@ namespace Acts { template class SeedFinder; class SeedFinderOptions; + class PerigeeSurface; } // namespace Acts /** @@ -116,8 +117,9 @@ public: Acts::MagneticFieldProvider::Cache& magCache) const; StatusCode tracking(const std::vector& paramseeds, const CKF& trackFinder, - const TrackFinderOptions& ckfOptions, Acts::MagneticFieldProvider::Cache& magCache, - edm4hep::TrackCollection& trackCollection) const; + const TrackFinderOptions& ckfOptions, const Propagator& extrapPropagator, + const Acts::PerigeeSurface& perigeeSurface, Propagator::Options<>& extrapOptions, + Acts::MagneticFieldProvider::Cache& magCache, edm4hep::TrackCollection& trackCollection) const; protected: /** diff --git a/k4ActsTracking/src/components/ACTSSeededCKFTrackingAlg.cxx b/k4ActsTracking/src/components/ACTSSeededCKFTrackingAlg.cxx index 27f87a93..06fc6d45 100644 --- a/k4ActsTracking/src/components/ACTSSeededCKFTrackingAlg.cxx +++ b/k4ActsTracking/src/components/ACTSSeededCKFTrackingAlg.cxx @@ -36,6 +36,7 @@ #include #include #include +#include // TBB #include @@ -294,6 +295,12 @@ std::tuple ACTSSeededCKFTrac Propagator propagator(std::move(stepper), std::move(navigator)); CKF trackFinder(std::move(propagator)); + // For extrapolating the fitted track back to the IP (perigee surface) + Stepper extrapStepper(magneticField()); + Navigator extrapNavigator(navigatorCfg); + Propagator extrapPropagator(std::move(extrapStepper), std::move(extrapNavigator)); + auto perigeeSurface = Acts::Surface::makeShared(Acts::Vector3{0., 0., 0.}); + // Set the options Acts::MeasurementSelector::Config measurementSelectorCfg = { {Acts::GeometryIdentifier(), {{}, {m_CKF_chi2CutOff}, {(std::size_t)(m_CKF_numMeasurementsCutOff)}}}}; @@ -304,9 +311,11 @@ std::tuple ACTSSeededCKFTrac pOptions.direction = Acts::Direction::Backward(); } - // Construct a perigee surface as the target surface - std::shared_ptr perigeeSurface = - Acts::Surface::makeShared(Acts::Vector3{0., 0., 0.}); + Propagator::Options<> extrapOptions{geometryContext(), magneticFieldContext()}; + extrapOptions.maxSteps = 10000; + if (m_propagateBackward) { + extrapOptions.direction = Acts::Direction::Backward(); + } Acts::GainMatrixUpdater kfUpdater; @@ -365,7 +374,9 @@ std::tuple ACTSSeededCKFTrac if (!m_runCKF) continue; - if (!tracking(paramseeds, trackFinder, ckfOptions, magCache, trackCollection).isSuccess()) { + if (!tracking(paramseeds, trackFinder, ckfOptions, extrapPropagator, *perigeeSurface, extrapOptions, magCache, + trackCollection) + .isSuccess()) { warning() << "Tracking failed for this event" << endmsg; } } @@ -389,6 +400,9 @@ std::tuple ACTSSeededCKFTrac // CKF tracking, StatusCode ACTSSeededCKFTrackingAlg::tracking(const std::vector& paramseeds, const CKF& trackFinder, const TrackFinderOptions& ckfOptions, + const Propagator& extrapPropagator, + const Acts::PerigeeSurface& perigeeSurface, + Propagator::Options<>& extrapOptions, Acts::MagneticFieldProvider::Cache& magCache, edm4hep::TrackCollection& trackCollection) const { // Initialize track finder @@ -416,6 +430,16 @@ StatusCode ACTSSeededCKFTrackingAlg::tracking(const std::vector #include #include +#include // TBB #include @@ -228,6 +229,13 @@ struct CKFTrackingAlg final // geometry. Built once in initialize() and reused read-only across threads. std::optional m_caloPropagator{}; + // Propagator (with tracking-geometry navigator) used to extrapolate the + // fitted track back to the perigee surface at the IP, so the track parameters + // (in particular D0/Z0) are expressed there. Built once in initialize() and + // reused read-only across threads. + std::optional m_extrapPropagator{}; + std::shared_ptr m_perigeeSurface{}; + k4ActsTracking::CellIDSelector m_seedSelector{}; mutable std::mutex m_seedMutex{}; @@ -273,6 +281,15 @@ StatusCode CKFTrackingAlg::initialize() { Propagator propagator(std::move(stepper), std::move(navigator)); m_trackFinder.emplace(std::move(propagator)); + // Build the propagator used to extrapolate fitted tracks back to the IP + // perigee surface. It shares the tracking-geometry navigator configuration + // with the CKF so it can navigate through the detector to the beam line. + Stepper extrapStepper(m_actsGeoSvc->magneticField()); + Navigator extrapNavigator(navigatorCfg); + Propagator extrapPropagator(std::move(extrapStepper), std::move(extrapNavigator)); + m_extrapPropagator.emplace(std::move(extrapPropagator)); + m_perigeeSurface = Acts::Surface::makeShared(Acts::Vector3{0., 0., 0.}); + // The calorimeter face lies outside the tracking geometry, so extrapolation // there uses a geometry-free (VoidNavigator) propagator. Only build it when // requested and when the geometry service actually provides calo surfaces. @@ -662,6 +679,27 @@ StatusCode CKFTrackingAlg::tracking(const std::vector extrapOptions{geoCtx, magCtx}; + extrapOptions.maxSteps = 10000; + if (m_propagateBackward) { + extrapOptions.direction = Acts::Direction::Backward(); + } + auto extrapResult = + Acts::extrapolateTrackToReferenceSurface(trackTip, *m_perigeeSurface, *m_extrapPropagator, extrapOptions, + Acts::TrackExtrapolationStrategy::firstOrLast); + if (!extrapResult.ok()) { + warning() << "Track extrapolation to perigee failed: " << extrapResult.error() << endmsg; + continue; + } + } + auto track = ACTSTracking::ACTS2edm4hep_track(trackTip, m_actsGeoSvc->magneticField(), magCache); // Extrapolate the fitted track to the calorimeter face and add an From d8ff6e75b1b8883ab18953c93262d71ffa831c73 Mon Sep 17 00:00:00 2001 From: Federico Meloni Date: Tue, 9 Jun 2026 15:31:54 +0200 Subject: [PATCH 05/10] Fix onLayer callback return type for acts main compatibility acts tightened the OnLayerReturnsNode concept in BlueprintBuilder.hpp to require the callback result to be exactly std::shared_ptr (std::same_as), so unsetXYCoG returning the derived std::shared_ptr no longer satisfies it and ElementLayerAssembler::onLayer fails to resolve. Return the base BlueprintNodePtr instead; `return layer;` upcasts implicitly. This satisfies the current concept and remains valid if acts-project/acts#5563 (which relaxes same_as -> convertible_to) is merged. Co-Authored-By: Claude Opus 4.8 --- .../src/components/DD4hepBlueprintConstruction.cpp | 9 +++++++-- 1 file changed, 7 insertions(+), 2 deletions(-) diff --git a/k4ActsTracking/src/components/DD4hepBlueprintConstruction.cpp b/k4ActsTracking/src/components/DD4hepBlueprintConstruction.cpp index 30865ac9..8277086e 100644 --- a/k4ActsTracking/src/components/DD4hepBlueprintConstruction.cpp +++ b/k4ActsTracking/src/components/DD4hepBlueprintConstruction.cpp @@ -36,6 +36,7 @@ #include #include +using Acts::Experimental::BlueprintNode; using Acts::Experimental::ContainerBlueprintNode; using Acts::Experimental::CylinderContainerBlueprintNode; using Acts::Experimental::LayerBlueprintNode; @@ -232,8 +233,12 @@ namespace Blueprints { /// the center of gravity for auto-sizing. This is useful for cases where the /// detectors have has an odd number of modules, which shifts them off the /// z-axis with the default sizing - std::shared_ptr unsetXYCoG(const std::optional&, - std::shared_ptr layer) { + // Return type is the base BlueprintNodePtr (not the derived LayerBlueprintNode) + // so the callback satisfies acts' OnLayerReturnsNode concept, which requires + // the result to be exactly std::shared_ptr. `return layer;` + // upcasts implicitly. + std::shared_ptr unsetXYCoG(const std::optional&, + std::shared_ptr layer) { layer->setUseCenterOfGravity(false, false, true); return layer; } From c9a143323d5d89606d11994e9c5ce419f1a5e7d2 Mon Sep 17 00:00:00 2001 From: Federico Meloni Date: Sun, 21 Jun 2026 10:20:09 +0200 Subject: [PATCH 06/10] fix merge --- k4ActsTracking/src/components/ACTSSeededCKFTrackingAlg.cxx | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/k4ActsTracking/src/components/ACTSSeededCKFTrackingAlg.cxx b/k4ActsTracking/src/components/ACTSSeededCKFTrackingAlg.cxx index 8fe24d33..72d00309 100644 --- a/k4ActsTracking/src/components/ACTSSeededCKFTrackingAlg.cxx +++ b/k4ActsTracking/src/components/ACTSSeededCKFTrackingAlg.cxx @@ -480,7 +480,9 @@ std::tuple ACTSSeededCKFTrac if (!m_runCKF) return; - if (!tracking(paramseeds, trackFinder, ckfOptions, localMagCache, trackCollection).isSuccess()) { + if (!tracking(paramseeds, trackFinder, ckfOptions, extrapPropagator, *perigeeSurface, extrapOptions, localMagCache, + trackCollection) + .isSuccess()) { warning() << "Tracking failed for this event" << endmsg; } }; // parallelSeedingAndTracking From 87a48f110e01cf2a906e7da2bf80f696ce464566 Mon Sep 17 00:00:00 2001 From: Federico Meloni Date: Mon, 22 Jun 2026 12:40:42 +0200 Subject: [PATCH 07/10] attempt to include calo in tracking geo and use standard navigation --- .../include/k4ActsTracking/Helpers.hxx | 64 +++-- .../include/k4ActsTracking/IActsGeoSvc.h | 31 ++- k4ActsTracking/src/components/ActsGeoSvc.cpp | 231 +++++++++++------- k4ActsTracking/src/components/ActsGeoSvc.h | 8 +- .../src/components/CKFTrackingAlg.cpp | 101 ++++++-- .../DD4hepBlueprintConstruction.cpp | 145 +++++++++-- .../components/DD4hepBlueprintConstruction.h | 12 +- k4ActsTracking/src/components/Helpers.cxx | 100 ++++---- 8 files changed, 483 insertions(+), 209 deletions(-) diff --git a/k4ActsTracking/include/k4ActsTracking/Helpers.hxx b/k4ActsTracking/include/k4ActsTracking/Helpers.hxx index bfc831aa..d6ee5748 100644 --- a/k4ActsTracking/include/k4ActsTracking/Helpers.hxx +++ b/k4ActsTracking/include/k4ActsTracking/Helpers.hxx @@ -43,9 +43,10 @@ #include #include #include +#include #include +#include #include -#include #include #include @@ -55,14 +56,15 @@ // Standard #include +#include namespace ACTSTracking { - /// Propagator used to extrapolate fitted tracks out to the calorimeter face. - /// It uses a VoidNavigator (no tracking geometry), so it can reach target - /// surfaces that lie outside the tracking-geometry world volume. Material - /// between the tracker and the calorimeter is not accounted for (field-only). - using CaloFacePropagator = Acts::Propagator, Acts::VoidNavigator>; + /// Geometry-aware propagator used to extrapolate fitted tracks out to the + /// calorimeter face. It navigates the tracking geometry (which now includes + /// the passive calo volumes), so material along the way is accounted for and + /// the actual curved trajectory selects which calo surface is hit. + using CaloFacePropagator = Acts::Propagator, Acts::Navigator>; using TrackResult = Acts::TrackContainer::TrackProxy; @@ -121,28 +123,42 @@ namespace ACTSTracking { */ Acts::ParticleHypothesis convertParticle(const edm4hep::MCParticle mcParticle); + //! Outcome of a calorimeter-face extrapolation. + enum class CaloExtrapolationStatus { + Ok, ///< reached a calo-face surface + NoSurfaces, ///< no calo-face surfaces are configured + NotReached, ///< propagation finished without reaching a calo-face surface + PropagationError, ///< the propagation itself failed + }; + + //! Result of a calorimeter-face extrapolation. + struct CaloExtrapolationResult { + std::optional params{}; + CaloExtrapolationStatus status{CaloExtrapolationStatus::NotReached}; + }; + //! Extrapolate track parameters to the calorimeter face. /** - * Selects, among the calorimeter-face surfaces, the one the track reaches first - * (smallest positive path length with a valid, bounds-checked intersection) and - * propagates the given parameters to it. The barrel is a ring of planar - * surfaces (one per polygon side) and the endcaps are flat discs; intersecting - * all candidates naturally handles both the barrel-sector choice and the - * barrel/endcap transition. + * Propagates the given parameters through the tracking geometry (which contains + * the passive calo volumes) and terminates as soon as the actual trajectory + * reaches one of the calorimeter-face surfaces, identified by their geometry + * ids. Because the propagation follows the real curved trajectory through the + * geometry, no straight-line pre-selection of the target surface is needed and + * material handling integrates naturally once the calo volumes carry material. * - * \param propagator Field-only propagator (VoidNavigator) - * \param start Track parameters to extrapolate from (e.g. at the last hit) - * \param surfaces Calorimeter-face surfaces from IActsGeoSvc - * \param gctx Geometry context - * \param mctx Magnetic-field context + * \param propagator Geometry-aware propagator + * \param start Track parameters to extrapolate from (e.g. at the last hit) + * \param caloSurfaceGeoIds Geometry ids of the calo-face surfaces (from IActsGeoSvc) + * \param gctx Geometry context + * \param mctx Magnetic-field context * - * \return Bound track parameters at the calorimeter face, or std::nullopt if no - * surface is reached or the propagation fails. + * \return Bound track parameters at the calorimeter face together with a status + * describing why the extrapolation succeeded or failed. */ - std::optional extrapolateToCaloFace(const CaloFacePropagator& propagator, - const Acts::BoundTrackParameters& start, - const IActsGeoSvc::CaloFaceSurfaces& surfaces, - const Acts::GeometryContext& gctx, - const Acts::MagneticFieldContext& mctx); + CaloExtrapolationResult extrapolateToCaloFace(const CaloFacePropagator& propagator, + const Acts::BoundTrackParameters& start, + const std::vector& caloSurfaceGeoIds, + const Acts::GeometryContext& gctx, + const Acts::MagneticFieldContext& mctx); } // namespace ACTSTracking diff --git a/k4ActsTracking/include/k4ActsTracking/IActsGeoSvc.h b/k4ActsTracking/include/k4ActsTracking/IActsGeoSvc.h index 0baa0cb6..d07f0eba 100644 --- a/k4ActsTracking/include/k4ActsTracking/IActsGeoSvc.h +++ b/k4ActsTracking/include/k4ActsTracking/IActsGeoSvc.h @@ -24,6 +24,8 @@ #include +#include + #include #include #include @@ -46,17 +48,27 @@ class GAUDI_API IActsGeoSvc : virtual public IService { using CellIDSurfaceMap = std::unordered_map; /// Surfaces approximating the inner face of the electromagnetic calorimeter, - /// derived from the DD4hep geometry. These live outside the tracking-geometry - /// world volume and are meant as target surfaces for a geometry-free - /// extrapolation (see ACTSTracking::extrapolateToCaloFace). + /// derived from the DD4hep geometry. They are inserted into the tracking + /// geometry as passive surfaces of dedicated calo volumes, so that the + /// standard geometry-aware propagation can navigate to them (see + /// ACTSTracking::extrapolateToCaloFace). /// /// The barrel is a regular polygon, so its inner face is modelled as one /// planar surface per polygon side rather than a cylinder. The endcaps are - /// flat discs at constant z. + /// flat discs at constant z. The struct also carries the bounding cylinder + /// dimensions used to build the enclosing TrackingVolumes (ACTS units). struct CaloFaceSurfaces { - std::vector> barrelFaces; ///< one plane per polygon side - std::shared_ptr endcapPos; ///< disc at +z, may be null - std::shared_ptr endcapNeg; ///< disc at -z, may be null + std::vector> barrelFaces; ///< one plane per polygon side + std::shared_ptr endcapPos; ///< disc at +z, may be null + std::shared_ptr endcapNeg; ///< disc at -z, may be null + + // Bounding-cylinder dimensions for the enclosing TrackingVolumes. + double barrelRMin = 0; ///< inner radius of the barrel calo volume (~apothem) + double barrelRMax = 0; ///< outer radius of the barrel calo volume (~circumradius) + double barrelHalfZ = 0; ///< barrel half-length along z + double endcapRMin = 0; ///< inner radius of the endcap disc volumes + double endcapRMax = 0; ///< outer radius of the endcap disc volumes + double endcapZ = 0; ///< |z| of the endcap inner face bool empty() const { return barrelFaces.empty() && !endcapPos && !endcapNeg; } }; @@ -70,6 +82,11 @@ class GAUDI_API IActsGeoSvc : virtual public IService { virtual std::string cellIDEncodingString() const = 0; virtual const CaloFaceSurfaces& caloFaceSurfaces() const = 0; + /// Geometry identifiers of the calorimeter-face surfaces, valid after the + /// tracking geometry has been constructed. Used by the extrapolation to + /// recognise when the propagation has reached the calo face. + virtual const std::vector& caloSurfaceGeoIds() const = 0; + virtual ~IActsGeoSvc() = default; }; diff --git a/k4ActsTracking/src/components/ActsGeoSvc.cpp b/k4ActsTracking/src/components/ActsGeoSvc.cpp index dea36c26..7c98ea35 100644 --- a/k4ActsTracking/src/components/ActsGeoSvc.cpp +++ b/k4ActsTracking/src/components/ActsGeoSvc.cpp @@ -106,6 +106,13 @@ StatusCode ActsGeoSvc::initialize() { {.dd4hepDetector = dd4hepDet, .lengthScale = Acts::UnitConstants::cm / dd4hep::cm, .gctx = gctxt}, gaudiLogger->cloneWithSuffix("|BlpBld")}; + // The calo-face surfaces are inserted into the tracking geometry as passive + // surfaces of dedicated calo volumes, so they must be built before the + // blueprint is populated and constructed. + if (m_buildCaloSurfaces.value()) { + buildCaloFaceSurfaces(); + } + using Acts::Experimental::Blueprint; using Acts::Experimental::BlueprintOptions; using namespace Acts::UnitLiterals; @@ -120,7 +127,7 @@ StatusCode ActsGeoSvc::initialize() { debug() << fmt::format("Getting Blueprint construction function for detector: {}", detName) << endmsg; if (const auto it = m_bluePrintPopulationFuncs.find(detName); it != m_bluePrintPopulationFuncs.end()) { auto bluePrintFunc = it->second; - bluePrintFunc(detName, root, builder); + bluePrintFunc(detName, root, builder, m_caloFaceSurfaces); } else { error() << fmt::format("Cannot find a Blueprint construction function for detector: {}", detName) << endmsg; return StatusCode::FAILURE; @@ -133,9 +140,15 @@ StatusCode ActsGeoSvc::initialize() { std::size_t nSurfaces = 0; m_trackingGeo->visitSurfaces([&](const Acts::Surface* surface) { + // Skip surfaces that are not backed by a DD4hep detector element, such as + // the passive calorimeter-face surfaces inserted by buildCaloFaceSurfaces. + const auto* actsDetElemPtr = + dynamic_cast(surface->surfacePlacement()); + if (actsDetElemPtr == nullptr) { + return; + } nSurfaces++; - const auto& actsDetElem = dynamic_cast(*surface->surfacePlacement()); - const auto& detElem = actsDetElem.sourceElement(); + const auto& detElem = actsDetElemPtr->sourceElement(); verbose() << fmt::format("Adding Acts surface {} pointing to dd4hep DetElement {}", surface->geometryId(), detElem.volumeID()) << endmsg; @@ -166,8 +179,22 @@ StatusCode ActsGeoSvc::initialize() { vis.write(m_objDumpFileName.value()); } + // The calo-face surfaces are the very ones inserted into the calo volumes, so + // after construction they carry their assigned geometry identifiers. Collect + // them for the extrapolation aborter to recognise the calo face. + m_caloSurfaceGeoIds.clear(); if (m_buildCaloSurfaces.value()) { - buildCaloFaceSurfaces(); + auto collect = [&](const std::shared_ptr& surface) { + if (surface) { + m_caloSurfaceGeoIds.push_back(surface->geometryId()); + } + }; + for (const auto& face : m_caloFaceSurfaces.barrelFaces) { + collect(face); + } + collect(m_caloFaceSurfaces.endcapPos); + collect(m_caloFaceSurfaces.endcapNeg); + info() << fmt::format("Collected {} calorimeter-face surface geometry ids.", m_caloSurfaceGeoIds.size()) << endmsg; } return StatusCode::SUCCESS; @@ -192,107 +219,143 @@ void ActsGeoSvc::buildCaloFaceSurfaces() { dd4hepDet->detectors(DetType::CALORIMETER | DetType::ELECTROMAGNETIC | DetType::BARREL, DetType::FORWARD); const auto ecalEndcap = dd4hepDet->detectors(DetType::CALORIMETER | DetType::ELECTROMAGNETIC | DetType::ENDCAP, 0); - // Barrel circumradius (corner radius), needed both for the barrel faces and to - // guarantee the endcap discs reach far enough to avoid hermeticity gaps at the - // barrel/endcap junction. + // --- Pass 1: extract the barrel and endcap dimensions -------------------- + // Surfaces are created in pass 2, after a corner-gap correction that needs + // both the barrel half-length and the endcap inner-face z. + bool haveBarrel = false; + int nSides = 0; + double apothem = 0.0; + double barrelHalfZRaw = 0.0; + double phi0 = 0.0; double barrelCircumradius = 0.0; - double barrelHalfZ = 0.0; + double halfWidth = 0.0; - // --- Barrel: one planar surface per polygon side -------------------------- if (ecalBarrel.empty()) { warning() << "No electromagnetic barrel calorimeter found via DetType flags; " "no barrel calo-face surfaces will be built." << endmsg; + } else if (const auto* caloData = ecalBarrel.front().extension(false); + caloData == nullptr) { + warning() << "ECAL barrel DetElement has no LayeredCalorimeterData extension; " + "skipping barrel calo-face surfaces." + << endmsg; } else { - const auto* caloData = ecalBarrel.front().extension(false); - if (caloData == nullptr) { - warning() << "ECAL barrel DetElement has no LayeredCalorimeterData extension; " - "skipping barrel calo-face surfaces." + nSides = caloData->inner_symmetry > 0 ? caloData->inner_symmetry : 0; + // For a BarrelLayout the calorimeter is centred on z = 0 and extent[] is + // {rmin, rmax, zmin=0, zmax=half_length}, so the half-length is extent[3]. + apothem = caloData->extent[0] * lengthScale; // perpendicular distance to inner face + barrelHalfZRaw = caloData->extent[3] * lengthScale; // barrel half-length + phi0 = caloData->inner_phi0; // azimuth of first inner-face normal + if (nSides < 3) { + warning() << fmt::format("ECAL barrel has unusable inner_symmetry={}; skipping barrel faces.", + caloData->inner_symmetry) << endmsg; } else { - const int nSides = caloData->inner_symmetry > 0 ? caloData->inner_symmetry : 0; - // For a BarrelLayout the calorimeter is centred on z = 0 and extent[] is - // {rmin, rmax, zmin=0, zmax=half_length}, so the half-length is extent[3]. - const double apothem = caloData->extent[0] * lengthScale; // perpendicular distance to inner face - barrelHalfZ = caloData->extent[3] * lengthScale; // barrel half-length - const double phi0 = caloData->inner_phi0; // azimuth of first inner-face normal - - if (nSides < 3) { - warning() << fmt::format("ECAL barrel has unusable inner_symmetry={}; skipping barrel faces.", - caloData->inner_symmetry) - << endmsg; - } else { - const double dPhi = 2 * std::numbers::pi / nSides; - barrelCircumradius = apothem / std::cos(std::numbers::pi / nSides); - const double halfWidth = apothem * std::tan(std::numbers::pi / nSides); - - m_caloFaceSurfaces.barrelFaces.reserve(nSides); - for (int i = 0; i < nSides; ++i) { - const double phi = phi0 + i * dPhi; - const double cphi = std::cos(phi); - const double sphi = std::sin(phi); - Acts::Vector3 normal{cphi, sphi, 0}; // local z (surface normal, radial) - Acts::Vector3 localX{-sphi, cphi, 0}; // tangential - Acts::Vector3 localY{0, 0, 1}; // along global z - Acts::Vector3 center = apothem * normal; // barrel centred on z = 0 - - Acts::Transform3 transform = Acts::Transform3::Identity(); - transform.linear().col(0) = localX; - transform.linear().col(1) = localY; - transform.linear().col(2) = normal; - transform.translation() = center; - - auto bounds = std::make_shared(halfWidth, barrelHalfZ); - m_caloFaceSurfaces.barrelFaces.push_back(Acts::Surface::makeShared(transform, bounds)); - } - info() << fmt::format( - "Built ECAL barrel calo face: {} planar faces, apothem={:.1f} mm, circumradius={:.1f} mm, " - "halfZ={:.1f} mm, phi0={:.4f}", - nSides, apothem / UC::mm, barrelCircumradius / UC::mm, barrelHalfZ / UC::mm, phi0) - << endmsg; - } + barrelCircumradius = apothem / std::cos(std::numbers::pi / nSides); + halfWidth = apothem * std::tan(std::numbers::pi / nSides); + haveBarrel = true; } } - // --- Endcaps: flat discs, widened to stay hermetic with the barrel -------- + bool haveEndcap = false; + double rMin = 0.0; + double rMax = 0.0; + double zEndcap = 0.0; + if (ecalEndcap.empty()) { warning() << "No electromagnetic endcap calorimeter found via DetType flags; " "no endcap calo-face surfaces will be built." << endmsg; + } else if (const auto* caloData = ecalEndcap.front().extension(false); + caloData == nullptr) { + warning() << "ECAL endcap DetElement has no LayeredCalorimeterData extension; " + "skipping endcap calo-face surfaces." + << endmsg; } else { - const auto* caloData = ecalEndcap.front().extension(false); - if (caloData == nullptr) { - warning() << "ECAL endcap DetElement has no LayeredCalorimeterData extension; " - "skipping endcap calo-face surfaces." + rMin = caloData->extent[0] * lengthScale; + rMax = caloData->extent[1] * lengthScale; + zEndcap = caloData->extent[2] * lengthScale; // inner-face z (zmin) + // Guarantee the disc reaches at least the barrel corner radius so there is + // no gap at the barrel/endcap junction. + if (barrelCircumradius > rMax) { + warning() << fmt::format( + "ECAL endcap rMax ({:.1f} mm) is smaller than the barrel circumradius ({:.1f} mm); " + "extending the endcap disc to close the hermeticity gap.", + rMax / UC::mm, barrelCircumradius / UC::mm) << endmsg; - } else { - const double rMin = caloData->extent[0] * lengthScale; - double rMax = caloData->extent[1] * lengthScale; - const double zEndcap = caloData->extent[2] * lengthScale; // inner-face z (zmin) - - // Guarantee the disc reaches at least the barrel corner radius so there is - // no gap at the barrel/endcap junction. - if (barrelCircumradius > rMax) { - warning() << fmt::format( - "ECAL endcap rMax ({:.1f} mm) is smaller than the barrel circumradius ({:.1f} mm); " - "extending the endcap disc to close the hermeticity gap.", - rMax / UC::mm, barrelCircumradius / UC::mm) - << endmsg; - rMax = barrelCircumradius; - } - - Acts::Transform3 tPos = Acts::Transform3::Identity(); - tPos.translation() = Acts::Vector3{0, 0, zEndcap}; - Acts::Transform3 tNeg = Acts::Transform3::Identity(); - tNeg.translation() = Acts::Vector3{0, 0, -zEndcap}; + rMax = barrelCircumradius; + } + haveEndcap = true; + } - m_caloFaceSurfaces.endcapPos = Acts::Surface::makeShared(tPos, rMin, rMax); - m_caloFaceSurfaces.endcapNeg = Acts::Surface::makeShared(tNeg, rMin, rMax); + // The ECAL barrel inner face and the endcap inner face meet at a hermetic + // corner, so the barrel half-length and the endcap z are almost equal. When + // the calo surfaces are placed in the tracking geometry, the barrel volume + // and the endcap volume are stacked along z and must not overlap. Shorten the + // barrel face slightly so the endcap disc sits clearly beyond the barrel + // (z-)extent; tracks crossing the trimmed corner strip are still caught by the + // endcap disc. + constexpr double cornerGap = 5.0 * UC::mm; + double barrelHalfZ = barrelHalfZRaw; + if (haveBarrel && haveEndcap && (zEndcap - barrelHalfZRaw) < cornerGap) { + barrelHalfZ = std::max(0.0, zEndcap - cornerGap); + info() << fmt::format( + "Shortening ECAL barrel face half-length from {:.1f} mm to {:.1f} mm to clear the endcap disc " + "at z={:.1f} mm (corner gap {:.1f} mm).", + barrelHalfZRaw / UC::mm, barrelHalfZ / UC::mm, zEndcap / UC::mm, cornerGap / UC::mm) + << endmsg; + } - info() << fmt::format("Built ECAL endcap calo faces: discs at z=+/-{:.1f} mm, rMin={:.1f} mm, rMax={:.1f} mm", - zEndcap / UC::mm, rMin / UC::mm, rMax / UC::mm) - << endmsg; + // --- Pass 2: build the surfaces ------------------------------------------ + if (haveBarrel) { + const double dPhi = 2 * std::numbers::pi / nSides; + m_caloFaceSurfaces.barrelFaces.reserve(nSides); + for (int i = 0; i < nSides; ++i) { + const double phi = phi0 + i * dPhi; + const double cphi = std::cos(phi); + const double sphi = std::sin(phi); + Acts::Vector3 normal{cphi, sphi, 0}; // local z (surface normal, radial) + Acts::Vector3 localX{-sphi, cphi, 0}; // tangential + Acts::Vector3 localY{0, 0, 1}; // along global z + Acts::Vector3 center = apothem * normal; // barrel centred on z = 0 + + Acts::Transform3 transform = Acts::Transform3::Identity(); + transform.linear().col(0) = localX; + transform.linear().col(1) = localY; + transform.linear().col(2) = normal; + transform.translation() = center; + + auto bounds = std::make_shared(halfWidth, barrelHalfZ); + m_caloFaceSurfaces.barrelFaces.push_back(Acts::Surface::makeShared(transform, bounds)); } + // Bounding cylinder for the barrel calo volume: from the inner face + // (apothem) out to the polygon corners (circumradius), trimmed half-length. + m_caloFaceSurfaces.barrelRMin = apothem; + m_caloFaceSurfaces.barrelRMax = barrelCircumradius; + m_caloFaceSurfaces.barrelHalfZ = barrelHalfZ; + info() << fmt::format( + "Built ECAL barrel calo face: {} planar faces, apothem={:.1f} mm, circumradius={:.1f} mm, " + "halfZ={:.1f} mm, phi0={:.4f}", + nSides, apothem / UC::mm, barrelCircumradius / UC::mm, barrelHalfZ / UC::mm, phi0) + << endmsg; + } + + if (haveEndcap) { + Acts::Transform3 tPos = Acts::Transform3::Identity(); + tPos.translation() = Acts::Vector3{0, 0, zEndcap}; + Acts::Transform3 tNeg = Acts::Transform3::Identity(); + tNeg.translation() = Acts::Vector3{0, 0, -zEndcap}; + + m_caloFaceSurfaces.endcapPos = Acts::Surface::makeShared(tPos, rMin, rMax); + m_caloFaceSurfaces.endcapNeg = Acts::Surface::makeShared(tNeg, rMin, rMax); + + m_caloFaceSurfaces.endcapRMin = rMin; + m_caloFaceSurfaces.endcapRMax = rMax; + m_caloFaceSurfaces.endcapZ = zEndcap; + + info() << fmt::format("Built ECAL endcap calo faces: discs at z=+/-{:.1f} mm, rMin={:.1f} mm, rMax={:.1f} mm", + zEndcap / UC::mm, rMin / UC::mm, rMax / UC::mm) + << endmsg; } info() << fmt::format("Calo-face surfaces built: {} barrel faces, {} endcap discs", diff --git a/k4ActsTracking/src/components/ActsGeoSvc.h b/k4ActsTracking/src/components/ActsGeoSvc.h index 26bf40a5..a7bd4224 100644 --- a/k4ActsTracking/src/components/ActsGeoSvc.h +++ b/k4ActsTracking/src/components/ActsGeoSvc.h @@ -55,6 +55,8 @@ class ActsGeoSvc : public extends { const CaloFaceSurfaces& caloFaceSurfaces() const override { return m_caloFaceSurfaces; } + const std::vector& caloSurfaceGeoIds() const override { return m_caloSurfaceGeoIds; } + ActsGeoSvc(const std::string& name, ISvcLocator* svcLoc); ~ActsGeoSvc() = default; @@ -78,12 +80,15 @@ class ActsGeoSvc : public extends { private: using BlueprintBuilder = ActsPlugins::DD4hep::BlueprintBuilder; - using BlueprintPopulationFunc = void(const std::string&, Acts::Experimental::Blueprint&, BlueprintBuilder&); + using BlueprintPopulationFunc = + void(const std::string&, Acts::Experimental::Blueprint&, BlueprintBuilder&, const CaloFaceSurfaces&); /// Build the ECAL inner-face surfaces (m_caloFaceSurfaces) from the DD4hep /// geometry. Surfaces are located via DetType flags and dimensioned from the /// dd4hep::rec::LayeredCalorimeterData extension. Missing ECAL sub-detectors /// or extensions are skipped with a warning rather than treated as an error. + /// Must run before the blueprint is constructed, since the surfaces are + /// inserted into the tracking geometry as passive calo volumes. void buildCaloFaceSurfaces(); SmartIF m_geoSvc; @@ -93,6 +98,7 @@ class ActsGeoSvc : public extends { std::unordered_map m_bluePrintPopulationFuncs{}; std::string m_cellIDEncodingString{}; CaloFaceSurfaces m_caloFaceSurfaces{}; + std::vector m_caloSurfaceGeoIds{}; }; inline std::shared_ptr ActsGeoSvc::trackingGeometry() const { return m_trackingGeo; } diff --git a/k4ActsTracking/src/components/CKFTrackingAlg.cpp b/k4ActsTracking/src/components/CKFTrackingAlg.cpp index 87b02aef..6a7d030f 100644 --- a/k4ActsTracking/src/components/CKFTrackingAlg.cpp +++ b/k4ActsTracking/src/components/CKFTrackingAlg.cpp @@ -52,6 +52,7 @@ #include #include #include +#include #include #include #include @@ -59,7 +60,6 @@ #include #include #include -#include #include #include #include @@ -84,6 +84,7 @@ #include // Standard +#include #include #include #include @@ -122,6 +123,8 @@ struct CKFTrackingAlg final StatusCode initialize() override; + StatusCode finalize() override; + std::tuple operator()( const edm4hep::TrackerHitPlaneCollection& trackerHitCollection, const edm4hep::TrackerHitSimTrackerHitLinkCollection& trackerHitRelations) const override; @@ -329,9 +332,11 @@ struct CKFTrackingAlg final // once in initialize() and reused (read-only) across events and threads. std::optional m_trackFinder{}; - // Field-only propagator (no tracking geometry) used to extrapolate fitted - // tracks out to the calorimeter face, which lies outside the tracking - // geometry. Built once in initialize() and reused read-only across threads. + // Geometry-aware propagator used to extrapolate fitted tracks out to the + // calorimeter face. The calo inner-face surfaces are part of the tracking + // geometry (passive surfaces of dedicated calo volumes), so the propagation + // follows the real trajectory through the detector and terminates on the + // surface actually reached. Built once in initialize(), reused across threads. std::optional m_caloPropagator{}; // Propagator (with tracking-geometry navigator) used to extrapolate the @@ -343,6 +348,15 @@ struct CKFTrackingAlg final k4ActsTracking::CellIDSelector m_seedSelector{}; + // Calorimeter-face extrapolation monitoring. operator() is const and runs on + // many threads, so the counters are mutable and atomic. Summarised in + // finalize() to report the rate at which the extrapolation fails. + mutable std::atomic m_caloAttempts{0}; ///< tracks with a usable start state + mutable std::atomic m_caloNoStartState{0}; ///< tracks without a measured smoothed state + mutable std::atomic m_caloNotReached{0}; ///< propagation did not reach a calo face + mutable std::atomic m_caloPropFailed{0}; ///< propagation itself failed + mutable std::atomic m_caloOk{0}; ///< reached a calo face + mutable std::mutex m_seedMutex{}; mutable std::mutex m_trackMutex{}; }; @@ -395,23 +409,49 @@ StatusCode CKFTrackingAlg::initialize() { m_extrapPropagator.emplace(std::move(extrapPropagator)); m_perigeeSurface = Acts::Surface::makeShared(Acts::Vector3{0., 0., 0.}); - // The calorimeter face lies outside the tracking geometry, so extrapolation - // there uses a geometry-free (VoidNavigator) propagator. Only build it when + // The calorimeter inner-face surfaces are part of the tracking geometry, so + // extrapolation there uses a geometry-aware propagator. The calo surfaces are + // passive, so the navigator must resolve passive surfaces. Only build it when // requested and when the geometry service actually provides calo surfaces. if (m_extrapolateToCalo) { - if (m_actsGeoSvc->caloFaceSurfaces().empty()) { + if (m_actsGeoSvc->caloSurfaceGeoIds().empty()) { warning() << "ExtrapolateToCalo requested but ActsGeoSvc provides no calorimeter-face surfaces; " "no AtCalorimeter track states will be produced." << endmsg; } else { - Acts::EigenStepper<> caloStepper(m_actsGeoSvc->magneticField()); - m_caloPropagator.emplace(std::move(caloStepper), Acts::VoidNavigator{}); + Navigator::Config caloNavigatorCfg{m_actsGeoSvc->trackingGeometry()}; + caloNavigatorCfg.resolvePassive = true; + caloNavigatorCfg.resolveMaterial = true; + caloNavigatorCfg.resolveSensitive = true; + + Stepper caloStepper(m_actsGeoSvc->magneticField()); + Navigator caloNavigator(caloNavigatorCfg); + m_caloPropagator.emplace(std::move(caloStepper), std::move(caloNavigator)); } } return StatusCode::SUCCESS; } +StatusCode CKFTrackingAlg::finalize() { + if (m_extrapolateToCalo) { + const std::size_t attempts = m_caloAttempts.load(); + const std::size_t ok = m_caloOk.load(); + const std::size_t notReached = m_caloNotReached.load(); + const std::size_t propFailed = m_caloPropFailed.load(); + const std::size_t noStart = m_caloNoStartState.load(); + const std::size_t failed = notReached + propFailed; + const double failRate = attempts > 0 ? static_cast(failed) / static_cast(attempts) : 0.0; + + info() << fmt::format( + "Calorimeter-face extrapolation summary: {} attempts, {} reached the face, {} failed " + "({:.2f}%: {} not reached, {} propagation errors); {} tracks had no measured smoothed start state.", + attempts, ok, failed, 100.0 * failRate, notReached, propFailed, noStart) + << endmsg; + } + return StatusCode::SUCCESS; +} + std::tuple CKFTrackingAlg::operator()( const edm4hep::TrackerHitPlaneCollection& trackerHitCollection, const edm4hep::TrackerHitSimTrackerHitLinkCollection& /*trackerHitRelations*/) const { @@ -946,28 +986,40 @@ StatusCode CKFTrackingAlg::tracking(const std::vectormagneticField(), magCache); // Extrapolate the fitted track to the calorimeter face and add an - // AtCalorimeter track state for Pandora / ParticleFlow. Starts from the - // outermost smoothed state (closest to the calorimeter). + // AtCalorimeter track state for Pandora / ParticleFlow. Start from the + // outermost smoothed state that carries a real measurement (closest to + // the calorimeter), so the extrapolation begins from a genuine fitted + // hit rather than a hole, outlier or material-only state. if (m_caloPropagator) { std::optional startParams; for (const auto& state : trackTip.trackStatesReversed()) { - if (state.hasSmoothed()) { + const auto flags = state.typeFlags(); + if (state.hasSmoothed() && flags.test(Acts::TrackStateFlag::HasMeasurement) && + !flags.test(Acts::TrackStateFlag::IsOutlier)) { startParams.emplace(state.referenceSurface().getSharedPtr(), state.smoothed(), state.smoothedCovariance(), trackTip.particleHypothesis()); break; } } - if (startParams) { + if (!startParams) { + ++m_caloNoStartState; + debug() << "No measured smoothed state available; no AtCalorimeter state added for this track." << endmsg; + } else { + ++m_caloAttempts; const Acts::MagneticFieldContext magCtx{}; - auto caloParams = ACTSTracking::extrapolateToCaloFace(*m_caloPropagator, *startParams, - m_actsGeoSvc->caloFaceSurfaces(), geoCtx, magCtx); - if (caloParams) { - const Acts::Vector3 caloPos = caloParams->position(geoCtx); + const auto caloResult = ACTSTracking::extrapolateToCaloFace( + *m_caloPropagator, *startParams, m_actsGeoSvc->caloSurfaceGeoIds(), geoCtx, magCtx); + + using ACTSTracking::CaloExtrapolationStatus; + switch (caloResult.status) { + case CaloExtrapolationStatus::Ok: { + ++m_caloOk; + const Acts::Vector3 caloPos = caloResult.params->position(geoCtx); auto fieldRes = m_actsGeoSvc->magneticField()->getField(caloPos, magCache); const double Bz = fieldRes.ok() ? (*fieldRes)[2] / Acts::UnitConstants::T : 0.0; auto caloState = - ACTSTracking::ACTS2edm4hep_trackState(edm4hep::TrackState::AtCalorimeter, *caloParams, Bz); + ACTSTracking::ACTS2edm4hep_trackState(edm4hep::TrackState::AtCalorimeter, *caloResult.params, Bz); // The calo-face parameters are local to the target surface, so the // edm4hep D0/Z0 from the generic conversion are not meaningful here. // Express the state at the impact point instead: set the reference @@ -976,10 +1028,21 @@ StatusCode CKFTrackingAlg::tracking(const std::vector #include #include #include #include #include +#include #include #include +#include +#include +#include #include +#include +#include #include #include #include #include +#include #include #include #include @@ -573,28 +581,117 @@ namespace Blueprints { return innerTracker; } + /// Navigation policy factory for the passive calo volumes. Each calo volume + /// holds only a handful of explicitly added surfaces (the polygon barrel + /// faces, or a single endcap disc). A TryAll policy (portals plus all passive + /// surfaces) is the simplest robust choice here: with so few surfaces there is + /// no benefit to a binned SurfaceArray, and TryAll needs no binning + /// configuration that could be mis-set. + std::shared_ptr makeCaloNavigationPolicyFactory() { + return std::make_shared( + Acts::NavigationPolicyFactory{}.add(Acts::TryAllNavigationPolicy::Config{})); + } + + /// Add the calorimeter barrel as a passive static volume to @p parent (the + /// radial container around the tracker). The volume is a cylinder enclosing + /// the regular-polygon inner face, with one planar surface per polygon side. + /// + /// A barrel of planar surfaces would normally be a Cylinder layer (a + /// cylindrical layer volume whose modules are binned into a SurfaceArray, as + /// the tracker barrels are). Here the face is just a few polygon planes, so a + /// hand-built static volume navigated with a TryAll policy is simpler and + /// avoids picking a SurfaceArray binning for a non-cylindrical polygon; it is + /// not a workaround for any missing layer type. + void addCaloBarrel(BlueprintNode& parent, const IActsGeoSvc::CaloFaceSurfaces& calo) { + constexpr double pad = 1_mm; + auto bounds = + std::make_shared(std::max(0.0, calo.barrelRMin - pad), calo.barrelRMax + pad, + calo.barrelHalfZ + pad); + auto vol = std::make_unique(Acts::Transform3::Identity(), std::move(bounds), "CaloBarrel"); + for (const auto& face : calo.barrelFaces) { + vol->addSurface(face); + } + parent.addStaticVolume(std::move(vol)).setNavigationPolicyFactory(makeCaloNavigationPolicyFactory()); + } + + /// Add one calorimeter endcap disc as a passive static volume to @p parent + /// (the top-level z container). The volume abuts the central region in z and + /// extends out beyond the disc face; it shares the central radial extent so + /// it stacks cleanly along z. + void addCaloEndcap(BlueprintNode& parent, const IActsGeoSvc::CaloFaceSurfaces& calo, bool positive) { + const auto& disc = positive ? calo.endcapPos : calo.endcapNeg; + if (!disc) { + return; + } + constexpr double pad = 1_mm; + // The central region (which holds the calo barrel) reaches barrelHalfZ plus + // the barrel volume's z-padding. Start the endcap just beyond that so the + // top-level z-stack sees a small gap rather than an overlap. The endcap + // disc (at endcapZ, which the surface builder keeps clear of the barrel) + // then sits comfortably inside the volume. + const double zInner = calo.barrelHalfZ + 2 * pad; + const double zOuter = calo.endcapZ + 10_mm; // beyond the disc face + const double halfZ = std::max(5_mm, (zOuter - zInner) / 2.0); + const double zc = (zInner + zOuter) / 2.0; + + // Span the full radius (0 .. barrel circumradius) so the volume shares the + // central radial extent and the z-stack does not need radial gap shells. + auto bounds = std::make_shared(0.0, calo.barrelRMax + pad, halfZ); + + Acts::Transform3 transform = Acts::Transform3::Identity(); + transform.translation() = Acts::Vector3{0, 0, positive ? zc : -zc}; + auto vol = std::make_unique(transform, std::move(bounds), + positive ? "CaloEndcapPos" : "CaloEndcapNeg"); + vol->addSurface(disc); + parent.addStaticVolume(std::move(vol)).setNavigationPolicyFactory(makeCaloNavigationPolicyFactory()); + } + } // namespace Blueprints namespace MuColl { namespace MAIA_v0 { void populateBlueprint(const std::string& detName, Acts::Experimental::Blueprint& root, - ActsPlugins::DD4hep::BlueprintBuilder& builder) { - auto& outer = root.addCylinderContainer(detName, AxisR); - Blueprints::addCylindricalBeampipe(outer); - - // NOTE: Need to set rather small padding here for the R-direction, because - // the innermost two layers are a double layer for which the cylindrical - // volumes are overlapping otherwise - auto vertexBarrel = Blueprints::makeGroupedBarrel(builder, "VertexBarrel", std::regex{"layer_\\d"}, "ZYX", - Blueprints::kTightBarrelEnvelope); - auto vertex = Blueprints::attachEndcaps(builder, std::move(vertexBarrel), Blueprints::DoubleBarrelLayerVertexSpec, - "Vertex"); - - auto innerTracker = Blueprints::makeNestedInnerTracker(builder, std::move(vertex)); - outer.addChild(innerTracker); + ActsPlugins::DD4hep::BlueprintBuilder& builder, const IActsGeoSvc::CaloFaceSurfaces& calo) { + // Build the tracker detectors as radial children of the supplied + // container. + auto buildTrackers = [&](ContainerBlueprintNode& outer) { + Blueprints::addCylindricalBeampipe(outer); + + // NOTE: Need to set rather small padding here for the R-direction, + // because the innermost two layers are a double layer for which the + // cylindrical volumes are overlapping otherwise + auto vertexBarrel = Blueprints::makeGroupedBarrel(builder, "VertexBarrel", std::regex{"layer_\\d"}, "ZYX", + Blueprints::kTightBarrelEnvelope); + auto vertex = Blueprints::attachEndcaps(builder, std::move(vertexBarrel), + Blueprints::DoubleBarrelLayerVertexSpec, "Vertex"); + + auto innerTracker = Blueprints::makeNestedInnerTracker(builder, std::move(vertex)); + outer.addChild(innerTracker); + + auto outerTracker = Blueprints::makeRegularTracker(builder, Blueprints::OuterTrackerSpec, "OuterTracker"); + outer.addChild(outerTracker); + }; + + if (calo.empty()) { + // No calorimeter face: keep the original purely-radial layout. + auto& outer = root.addCylinderContainer(detName, AxisR); + buildTrackers(outer); + return; + } - auto outerTracker = Blueprints::makeRegularTracker(builder, Blueprints::OuterTrackerSpec, "OuterTracker"); - outer.addChild(outerTracker); + // The calorimeter wraps the tracker: its endcaps reach to small radius at + // large |z| where the tracker does not extend. This is expressed as a + // top-level z-stack [calo -endcap | central (tracker + calo barrel) | + // calo +endcap], with the calo barrel as the outermost radial child of + // the central region. + auto& world = root.addCylinderContainer(detName, AxisZ); + Blueprints::addCaloEndcap(world, calo, /*positive=*/false); + auto& central = world.addCylinderContainer(detName + "Central", AxisR); + buildTrackers(central); + if (!calo.barrelFaces.empty()) { + Blueprints::addCaloBarrel(central, calo); + } + Blueprints::addCaloEndcap(world, calo, /*positive=*/true); } } // namespace MAIA_v0 } // namespace MuColl @@ -602,7 +699,9 @@ namespace MuColl { namespace FCCee { namespace ILD_FCCee_v01 { void populateBlueprint(const std::string& detName, Acts::Experimental::Blueprint& root, - ActsPlugins::DD4hep::BlueprintBuilder& builder) { + ActsPlugins::DD4hep::BlueprintBuilder& builder, + [[maybe_unused]] const IActsGeoSvc::CaloFaceSurfaces& calo) { + // TODO: integrate the calo face (see MAIA_v0); deferred until validated. auto& outer = root.addCylinderContainer(detName, AxisR); Blueprints::addCylindricalBeampipe(outer); @@ -631,7 +730,9 @@ namespace FCCee { namespace ILD_FCCee_v02 { void populateBlueprint(const std::string& detName, Acts::Experimental::Blueprint& root, - ActsPlugins::DD4hep::BlueprintBuilder& builder) { + ActsPlugins::DD4hep::BlueprintBuilder& builder, + [[maybe_unused]] const IActsGeoSvc::CaloFaceSurfaces& calo) { + // TODO: integrate the calo face (see MAIA_v0); deferred until validated. auto& outer = root.addCylinderContainer(detName, AxisR); Blueprints::addCylindricalBeampipe(outer); @@ -650,7 +751,9 @@ namespace FCCee { namespace CLD_o2_v07 { void populateBlueprint(const std::string& detName, Acts::Experimental::Blueprint& root, - ActsPlugins::DD4hep::BlueprintBuilder& builder) { + ActsPlugins::DD4hep::BlueprintBuilder& builder, + [[maybe_unused]] const IActsGeoSvc::CaloFaceSurfaces& calo) { + // TODO: integrate the calo face (see MAIA_v0); deferred until validated. auto& outer = root.addCylinderContainer(detName, AxisR); Blueprints::addCylindricalBeampipe(outer); auto vtxBarrel = Blueprints::makeBarrel(builder, Blueprints::UngroupedDoubleBarrelLayerVertexSpec, @@ -672,7 +775,9 @@ namespace FCCee { namespace LUXE { namespace LUXE_v0 { void populateBlueprint(const std::string& detName, Acts::Experimental::Blueprint& root, - ActsPlugins::DD4hep::BlueprintBuilder& builder) { + ActsPlugins::DD4hep::BlueprintBuilder& builder, + [[maybe_unused]] const IActsGeoSvc::CaloFaceSurfaces& calo) { + // No electromagnetic calorimeter face integration for LUXE. auto& tracker = root.addCuboidContainer(detName, AxisZ); auto envelope = Acts::ExtentEnvelope{}.set(AxisZ, {0.4_mm, 0.4_mm}).set(AxisX, {0.4_mm, 0.4_mm}).set(AxisY, {0.4_mm, 0.4_mm}); diff --git a/k4ActsTracking/src/components/DD4hepBlueprintConstruction.h b/k4ActsTracking/src/components/DD4hepBlueprintConstruction.h index e180a68a..bec6f37a 100644 --- a/k4ActsTracking/src/components/DD4hepBlueprintConstruction.h +++ b/k4ActsTracking/src/components/DD4hepBlueprintConstruction.h @@ -19,6 +19,8 @@ #ifndef K4ACTSTRACKING_DD4HEPBLUEPRINTCONSTRUCTION_H #define K4ACTSTRACKING_DD4HEPBLUEPRINTCONSTRUCTION_H +#include "k4ActsTracking/IActsGeoSvc.h" + #include #include @@ -32,24 +34,24 @@ namespace Acts::Experimental { namespace MuColl { namespace MAIA_v0 { void populateBlueprint(const std::string& detName, Acts::Experimental::Blueprint& root, - ActsPlugins::DD4hep::BlueprintBuilder& builder); + ActsPlugins::DD4hep::BlueprintBuilder& builder, const IActsGeoSvc::CaloFaceSurfaces& calo); } } // namespace MuColl namespace FCCee { namespace ILD_FCCee_v01 { void populateBlueprint(const std::string& detName, Acts::Experimental::Blueprint& root, - ActsPlugins::DD4hep::BlueprintBuilder& builder); + ActsPlugins::DD4hep::BlueprintBuilder& builder, const IActsGeoSvc::CaloFaceSurfaces& calo); } namespace ILD_FCCee_v02 { void populateBlueprint(const std::string& detName, Acts::Experimental::Blueprint& root, - ActsPlugins::DD4hep::BlueprintBuilder& builder); + ActsPlugins::DD4hep::BlueprintBuilder& builder, const IActsGeoSvc::CaloFaceSurfaces& calo); } namespace CLD_o2_v07 { void populateBlueprint(const std::string& detName, Acts::Experimental::Blueprint& root, - ActsPlugins::DD4hep::BlueprintBuilder& builder); + ActsPlugins::DD4hep::BlueprintBuilder& builder, const IActsGeoSvc::CaloFaceSurfaces& calo); } } // namespace FCCee @@ -57,7 +59,7 @@ namespace FCCee { namespace LUXE { namespace LUXE_v0 { void populateBlueprint(const std::string& detName, Acts::Experimental::Blueprint& root, - ActsPlugins::DD4hep::BlueprintBuilder& builder); + ActsPlugins::DD4hep::BlueprintBuilder& builder, const IActsGeoSvc::CaloFaceSurfaces& calo); } } // namespace LUXE diff --git a/k4ActsTracking/src/components/Helpers.cxx b/k4ActsTracking/src/components/Helpers.cxx index d3e49c4e..3438aea6 100644 --- a/k4ActsTracking/src/components/Helpers.cxx +++ b/k4ActsTracking/src/components/Helpers.cxx @@ -31,11 +31,14 @@ // ACTS #include #include +#include #include #include #include #include +#include +#include #include // ACTSTracking @@ -240,65 +243,64 @@ namespace ACTSTracking { return Acts::ParticleHypothesis{pdg, mass, charge_type}; } - std::optional extrapolateToCaloFace(const CaloFacePropagator& propagator, - const Acts::BoundTrackParameters& start, - const IActsGeoSvc::CaloFaceSurfaces& surfaces, - const Acts::GeometryContext& gctx, - const Acts::MagneticFieldContext& mctx) { - if (surfaces.empty()) { - return std::nullopt; + namespace { + /// Abort condition for the calorimeter-face extrapolation: terminate the + /// propagation as soon as the navigator's current surface is one of the + /// calorimeter-face surfaces. Works with the geometry navigator, which sets + /// the current surface as it visits the calo volumes' passive surfaces. + struct CaloSurfaceReached { + const std::vector* caloIds = nullptr; + + template + bool checkAbort(propagator_state_t& state, const stepper_t& /*stepper*/, const navigator_t& navigator, + const Acts::Logger& /*logger*/) const { + if (caloIds == nullptr) { + return false; + } + const Acts::Surface* current = navigator.currentSurface(state.navigation); + if (current == nullptr) { + return false; + } + return std::find(caloIds->begin(), caloIds->end(), current->geometryId()) != caloIds->end(); + } + }; + } // namespace + + CaloExtrapolationResult extrapolateToCaloFace(const CaloFacePropagator& propagator, + const Acts::BoundTrackParameters& start, + const std::vector& caloSurfaceGeoIds, + const Acts::GeometryContext& gctx, + const Acts::MagneticFieldContext& mctx) { + if (caloSurfaceGeoIds.empty()) { + return {std::nullopt, CaloExtrapolationStatus::NoSurfaces}; } - const Acts::Vector3 position = start.position(gctx); - const Acts::Vector3 direction = start.direction(); + using ActorList = Acts::ActorList; + using Options = CaloFacePropagator::Options; - // Allow a small slack at face edges / the barrel-endcap seam so tracks - // crossing right at a boundary are not lost. - constexpr double tolerance = 1.0 * Acts::UnitConstants::mm; - const Acts::BoundaryTolerance boundaryTolerance = Acts::BoundaryTolerance::AbsoluteEuclidean(tolerance); + Options options{gctx, mctx}; + options.maxSteps = 10000; + options.actorList.get().caloIds = &caloSurfaceGeoIds; - // Build the candidate list: every barrel face plus the endcap disc on the - // side the track is heading towards. - std::vector candidates; - candidates.reserve(surfaces.barrelFaces.size() + 1); - for (const auto& face : surfaces.barrelFaces) { - candidates.push_back(face.get()); - } - const auto& endcap = (direction.z() >= 0) ? surfaces.endcapPos : surfaces.endcapNeg; - if (endcap) { - candidates.push_back(endcap.get()); + auto result = propagator.propagate(start, options); + if (!result.ok()) { + return {std::nullopt, CaloExtrapolationStatus::PropagationError}; } - // Pick the surface reached first along the track direction. - const Acts::Surface* target = nullptr; - double bestPath = std::numeric_limits::max(); - constexpr double minPath = 1e-3; // ignore intersections essentially at the start point - for (const Acts::Surface* surface : candidates) { - const auto multiIntersection = surface->intersect(gctx, position, direction, boundaryTolerance); - for (const auto& intersection : multiIntersection) { - if (!intersection.isValid()) { - continue; - } - const double path = intersection.pathLength(); - if (path > minPath && path < bestPath) { - bestPath = path; - target = surface; - } - } + const auto& output = result.value(); + if (!output.endParameters.has_value()) { + return {std::nullopt, CaloExtrapolationStatus::NotReached}; } - if (target == nullptr) { - return std::nullopt; + // The propagation may also terminate at the world boundary; only treat it as + // a success if it actually finished on a calo-face surface. + const auto& endParams = output.endParameters.value(); + const auto endId = endParams.referenceSurface().geometryId(); + if (std::find(caloSurfaceGeoIds.begin(), caloSurfaceGeoIds.end(), endId) == caloSurfaceGeoIds.end()) { + return {std::nullopt, CaloExtrapolationStatus::NotReached}; } - Acts::PropagatorPlainOptions options{gctx, mctx}; - options.maxSteps = 10000; - - auto result = propagator.propagateToSurface(start, *target, options); - if (!result.ok()) { - return std::nullopt; - } - return result.value(); + return {endParams, CaloExtrapolationStatus::Ok}; } } // namespace ACTSTracking From 2c3d590e118427fb537d366a38c76a88d3773c06 Mon Sep 17 00:00:00 2001 From: Federico Meloni Date: Mon, 22 Jun 2026 13:09:03 +0200 Subject: [PATCH 08/10] clang-formats --- .../include/k4ActsTracking/Helpers.hxx | 2 +- k4ActsTracking/src/components/ActsGeoSvc.cpp | 21 +++--- k4ActsTracking/src/components/ActsGeoSvc.h | 4 +- .../src/components/CKFTrackingAlg.cpp | 74 +++++++++---------- .../DD4hepBlueprintConstruction.cpp | 21 +++--- k4ActsTracking/src/components/Helpers.cxx | 2 +- 6 files changed, 60 insertions(+), 64 deletions(-) diff --git a/k4ActsTracking/include/k4ActsTracking/Helpers.hxx b/k4ActsTracking/include/k4ActsTracking/Helpers.hxx index d6ee5748..ef321de8 100644 --- a/k4ActsTracking/include/k4ActsTracking/Helpers.hxx +++ b/k4ActsTracking/include/k4ActsTracking/Helpers.hxx @@ -41,9 +41,9 @@ #include #include #include +#include #include #include -#include #include #include #include diff --git a/k4ActsTracking/src/components/ActsGeoSvc.cpp b/k4ActsTracking/src/components/ActsGeoSvc.cpp index 7c98ea35..4bb38cf4 100644 --- a/k4ActsTracking/src/components/ActsGeoSvc.cpp +++ b/k4ActsTracking/src/components/ActsGeoSvc.cpp @@ -142,8 +142,7 @@ StatusCode ActsGeoSvc::initialize() { m_trackingGeo->visitSurfaces([&](const Acts::Surface* surface) { // Skip surfaces that are not backed by a DD4hep detector element, such as // the passive calorimeter-face surfaces inserted by buildCaloFaceSurfaces. - const auto* actsDetElemPtr = - dynamic_cast(surface->surfacePlacement()); + const auto* actsDetElemPtr = dynamic_cast(surface->surfacePlacement()); if (actsDetElemPtr == nullptr) { return; } @@ -222,20 +221,19 @@ void ActsGeoSvc::buildCaloFaceSurfaces() { // --- Pass 1: extract the barrel and endcap dimensions -------------------- // Surfaces are created in pass 2, after a corner-gap correction that needs // both the barrel half-length and the endcap inner-face z. - bool haveBarrel = false; - int nSides = 0; - double apothem = 0.0; - double barrelHalfZRaw = 0.0; - double phi0 = 0.0; + bool haveBarrel = false; + int nSides = 0; + double apothem = 0.0; + double barrelHalfZRaw = 0.0; + double phi0 = 0.0; double barrelCircumradius = 0.0; - double halfWidth = 0.0; + double halfWidth = 0.0; if (ecalBarrel.empty()) { warning() << "No electromagnetic barrel calorimeter found via DetType flags; " "no barrel calo-face surfaces will be built." << endmsg; - } else if (const auto* caloData = ecalBarrel.front().extension(false); - caloData == nullptr) { + } else if (const auto* caloData = ecalBarrel.front().extension(false); caloData == nullptr) { warning() << "ECAL barrel DetElement has no LayeredCalorimeterData extension; " "skipping barrel calo-face surfaces." << endmsg; @@ -266,8 +264,7 @@ void ActsGeoSvc::buildCaloFaceSurfaces() { warning() << "No electromagnetic endcap calorimeter found via DetType flags; " "no endcap calo-face surfaces will be built." << endmsg; - } else if (const auto* caloData = ecalEndcap.front().extension(false); - caloData == nullptr) { + } else if (const auto* caloData = ecalEndcap.front().extension(false); caloData == nullptr) { warning() << "ECAL endcap DetElement has no LayeredCalorimeterData extension; " "skipping endcap calo-face surfaces." << endmsg; diff --git a/k4ActsTracking/src/components/ActsGeoSvc.h b/k4ActsTracking/src/components/ActsGeoSvc.h index a7bd4224..0e4fb588 100644 --- a/k4ActsTracking/src/components/ActsGeoSvc.h +++ b/k4ActsTracking/src/components/ActsGeoSvc.h @@ -80,8 +80,8 @@ class ActsGeoSvc : public extends { private: using BlueprintBuilder = ActsPlugins::DD4hep::BlueprintBuilder; - using BlueprintPopulationFunc = - void(const std::string&, Acts::Experimental::Blueprint&, BlueprintBuilder&, const CaloFaceSurfaces&); + using BlueprintPopulationFunc = void(const std::string&, Acts::Experimental::Blueprint&, BlueprintBuilder&, + const CaloFaceSurfaces&); /// Build the ECAL inner-face surfaces (m_caloFaceSurfaces) from the DD4hep /// geometry. Surfaces are located via DetType flags and dimensioned from the diff --git a/k4ActsTracking/src/components/CKFTrackingAlg.cpp b/k4ActsTracking/src/components/CKFTrackingAlg.cpp index 6a7d030f..97ca2246 100644 --- a/k4ActsTracking/src/components/CKFTrackingAlg.cpp +++ b/k4ActsTracking/src/components/CKFTrackingAlg.cpp @@ -51,8 +51,8 @@ #include #include #include -#include #include +#include #include #include #include @@ -351,11 +351,11 @@ struct CKFTrackingAlg final // Calorimeter-face extrapolation monitoring. operator() is const and runs on // many threads, so the counters are mutable and atomic. Summarised in // finalize() to report the rate at which the extrapolation fails. - mutable std::atomic m_caloAttempts{0}; ///< tracks with a usable start state - mutable std::atomic m_caloNoStartState{0}; ///< tracks without a measured smoothed state - mutable std::atomic m_caloNotReached{0}; ///< propagation did not reach a calo face - mutable std::atomic m_caloPropFailed{0}; ///< propagation itself failed - mutable std::atomic m_caloOk{0}; ///< reached a calo face + mutable std::atomic m_caloAttempts{0}; ///< tracks with a usable start state + mutable std::atomic m_caloNoStartState{0}; ///< tracks without a measured smoothed state + mutable std::atomic m_caloNotReached{0}; ///< propagation did not reach a calo face + mutable std::atomic m_caloPropFailed{0}; ///< propagation itself failed + mutable std::atomic m_caloOk{0}; ///< reached a calo face mutable std::mutex m_seedMutex{}; mutable std::mutex m_trackMutex{}; @@ -1008,41 +1008,41 @@ StatusCode CKFTrackingAlg::tracking(const std::vectorcaloSurfaceGeoIds(), geoCtx, magCtx); using ACTSTracking::CaloExtrapolationStatus; switch (caloResult.status) { - case CaloExtrapolationStatus::Ok: { - ++m_caloOk; - const Acts::Vector3 caloPos = caloResult.params->position(geoCtx); - auto fieldRes = m_actsGeoSvc->magneticField()->getField(caloPos, magCache); - const double Bz = fieldRes.ok() ? (*fieldRes)[2] / Acts::UnitConstants::T : 0.0; - auto caloState = - ACTSTracking::ACTS2edm4hep_trackState(edm4hep::TrackState::AtCalorimeter, *caloResult.params, Bz); - // The calo-face parameters are local to the target surface, so the - // edm4hep D0/Z0 from the generic conversion are not meaningful here. - // Express the state at the impact point instead: set the reference - // point to the global calo-face position (D0 = Z0 = 0 there). - caloState.referencePoint = edm4hep::Vector3f(caloPos.x(), caloPos.y(), caloPos.z()); - caloState.D0 = 0; - caloState.Z0 = 0; - track.addToTrackStates(caloState); - break; - } - case CaloExtrapolationStatus::NotReached: - case CaloExtrapolationStatus::NoSurfaces: - ++m_caloNotReached; - debug() << "Extrapolation to the calorimeter face did not reach a surface; " - "no AtCalorimeter state added for this track." - << endmsg; - break; - case CaloExtrapolationStatus::PropagationError: - ++m_caloPropFailed; - debug() << "Extrapolation to the calorimeter face failed during propagation; " - "no AtCalorimeter state added for this track." - << endmsg; - break; + case CaloExtrapolationStatus::Ok: { + ++m_caloOk; + const Acts::Vector3 caloPos = caloResult.params->position(geoCtx); + auto fieldRes = m_actsGeoSvc->magneticField()->getField(caloPos, magCache); + const double Bz = fieldRes.ok() ? (*fieldRes)[2] / Acts::UnitConstants::T : 0.0; + auto caloState = + ACTSTracking::ACTS2edm4hep_trackState(edm4hep::TrackState::AtCalorimeter, *caloResult.params, Bz); + // The calo-face parameters are local to the target surface, so the + // edm4hep D0/Z0 from the generic conversion are not meaningful here. + // Express the state at the impact point instead: set the reference + // point to the global calo-face position (D0 = Z0 = 0 there). + caloState.referencePoint = edm4hep::Vector3f(caloPos.x(), caloPos.y(), caloPos.z()); + caloState.D0 = 0; + caloState.Z0 = 0; + track.addToTrackStates(caloState); + break; + } + case CaloExtrapolationStatus::NotReached: + case CaloExtrapolationStatus::NoSurfaces: + ++m_caloNotReached; + debug() << "Extrapolation to the calorimeter face did not reach a surface; " + "no AtCalorimeter state added for this track." + << endmsg; + break; + case CaloExtrapolationStatus::PropagationError: + ++m_caloPropFailed; + debug() << "Extrapolation to the calorimeter face failed during propagation; " + "no AtCalorimeter state added for this track." + << endmsg; + break; } } } diff --git a/k4ActsTracking/src/components/DD4hepBlueprintConstruction.cpp b/k4ActsTracking/src/components/DD4hepBlueprintConstruction.cpp index 68fa2c20..a2f38645 100644 --- a/k4ActsTracking/src/components/DD4hepBlueprintConstruction.cpp +++ b/k4ActsTracking/src/components/DD4hepBlueprintConstruction.cpp @@ -603,10 +603,9 @@ namespace Blueprints { /// avoids picking a SurfaceArray binning for a non-cylindrical polygon; it is /// not a workaround for any missing layer type. void addCaloBarrel(BlueprintNode& parent, const IActsGeoSvc::CaloFaceSurfaces& calo) { - constexpr double pad = 1_mm; - auto bounds = - std::make_shared(std::max(0.0, calo.barrelRMin - pad), calo.barrelRMax + pad, - calo.barrelHalfZ + pad); + constexpr double pad = 1_mm; + auto bounds = std::make_shared(std::max(0.0, calo.barrelRMin - pad), + calo.barrelRMax + pad, calo.barrelHalfZ + pad); auto vol = std::make_unique(Acts::Transform3::Identity(), std::move(bounds), "CaloBarrel"); for (const auto& face : calo.barrelFaces) { vol->addSurface(face); @@ -641,7 +640,7 @@ namespace Blueprints { Acts::Transform3 transform = Acts::Transform3::Identity(); transform.translation() = Acts::Vector3{0, 0, positive ? zc : -zc}; auto vol = std::make_unique(transform, std::move(bounds), - positive ? "CaloEndcapPos" : "CaloEndcapNeg"); + positive ? "CaloEndcapPos" : "CaloEndcapNeg"); vol->addSurface(disc); parent.addStaticVolume(std::move(vol)).setNavigationPolicyFactory(makeCaloNavigationPolicyFactory()); } @@ -662,8 +661,8 @@ namespace MuColl { // cylindrical volumes are overlapping otherwise auto vertexBarrel = Blueprints::makeGroupedBarrel(builder, "VertexBarrel", std::regex{"layer_\\d"}, "ZYX", Blueprints::kTightBarrelEnvelope); - auto vertex = Blueprints::attachEndcaps(builder, std::move(vertexBarrel), - Blueprints::DoubleBarrelLayerVertexSpec, "Vertex"); + auto vertex = Blueprints::attachEndcaps(builder, std::move(vertexBarrel), + Blueprints::DoubleBarrelLayerVertexSpec, "Vertex"); auto innerTracker = Blueprints::makeNestedInnerTracker(builder, std::move(vertex)); outer.addChild(innerTracker); @@ -699,7 +698,7 @@ namespace MuColl { namespace FCCee { namespace ILD_FCCee_v01 { void populateBlueprint(const std::string& detName, Acts::Experimental::Blueprint& root, - ActsPlugins::DD4hep::BlueprintBuilder& builder, + ActsPlugins::DD4hep::BlueprintBuilder& builder, [[maybe_unused]] const IActsGeoSvc::CaloFaceSurfaces& calo) { // TODO: integrate the calo face (see MAIA_v0); deferred until validated. auto& outer = root.addCylinderContainer(detName, AxisR); @@ -730,7 +729,7 @@ namespace FCCee { namespace ILD_FCCee_v02 { void populateBlueprint(const std::string& detName, Acts::Experimental::Blueprint& root, - ActsPlugins::DD4hep::BlueprintBuilder& builder, + ActsPlugins::DD4hep::BlueprintBuilder& builder, [[maybe_unused]] const IActsGeoSvc::CaloFaceSurfaces& calo) { // TODO: integrate the calo face (see MAIA_v0); deferred until validated. auto& outer = root.addCylinderContainer(detName, AxisR); @@ -751,7 +750,7 @@ namespace FCCee { namespace CLD_o2_v07 { void populateBlueprint(const std::string& detName, Acts::Experimental::Blueprint& root, - ActsPlugins::DD4hep::BlueprintBuilder& builder, + ActsPlugins::DD4hep::BlueprintBuilder& builder, [[maybe_unused]] const IActsGeoSvc::CaloFaceSurfaces& calo) { // TODO: integrate the calo face (see MAIA_v0); deferred until validated. auto& outer = root.addCylinderContainer(detName, AxisR); @@ -775,7 +774,7 @@ namespace FCCee { namespace LUXE { namespace LUXE_v0 { void populateBlueprint(const std::string& detName, Acts::Experimental::Blueprint& root, - ActsPlugins::DD4hep::BlueprintBuilder& builder, + ActsPlugins::DD4hep::BlueprintBuilder& builder, [[maybe_unused]] const IActsGeoSvc::CaloFaceSurfaces& calo) { // No electromagnetic calorimeter face integration for LUXE. auto& tracker = root.addCuboidContainer(detName, AxisZ); diff --git a/k4ActsTracking/src/components/Helpers.cxx b/k4ActsTracking/src/components/Helpers.cxx index 3438aea6..f47e35ff 100644 --- a/k4ActsTracking/src/components/Helpers.cxx +++ b/k4ActsTracking/src/components/Helpers.cxx @@ -279,7 +279,7 @@ namespace ACTSTracking { using Options = CaloFacePropagator::Options; Options options{gctx, mctx}; - options.maxSteps = 10000; + options.maxSteps = 10000; options.actorList.get().caloIds = &caloSurfaceGeoIds; auto result = propagator.propagate(start, options); From 02d80d27f24b7f86f7d0f7fc875b9d1717e64cc6 Mon Sep 17 00:00:00 2001 From: Federico Meloni Date: Mon, 22 Jun 2026 17:43:39 +0200 Subject: [PATCH 09/10] add calo for CLD too --- .../DD4hepBlueprintConstruction.cpp | 50 +++++++++++++------ 1 file changed, 35 insertions(+), 15 deletions(-) diff --git a/k4ActsTracking/src/components/DD4hepBlueprintConstruction.cpp b/k4ActsTracking/src/components/DD4hepBlueprintConstruction.cpp index a2f38645..b52e35f3 100644 --- a/k4ActsTracking/src/components/DD4hepBlueprintConstruction.cpp +++ b/k4ActsTracking/src/components/DD4hepBlueprintConstruction.cpp @@ -750,23 +750,43 @@ namespace FCCee { namespace CLD_o2_v07 { void populateBlueprint(const std::string& detName, Acts::Experimental::Blueprint& root, - ActsPlugins::DD4hep::BlueprintBuilder& builder, - [[maybe_unused]] const IActsGeoSvc::CaloFaceSurfaces& calo) { - // TODO: integrate the calo face (see MAIA_v0); deferred until validated. - auto& outer = root.addCylinderContainer(detName, AxisR); - Blueprints::addCylindricalBeampipe(outer); - auto vtxBarrel = Blueprints::makeBarrel(builder, Blueprints::UngroupedDoubleBarrelLayerVertexSpec, - Blueprints::kBarrelEnvelope, Blueprints::doubleLayerKey); - auto vertex = Blueprints::attachEndcaps(builder, std::move(vtxBarrel), - Blueprints::UngroupedDoubleBarrelLayerVertexSpec, "Vertex"); + ActsPlugins::DD4hep::BlueprintBuilder& builder, const IActsGeoSvc::CaloFaceSurfaces& calo) { + // Build the tracker detectors as radial children of the supplied + // container. + auto buildTrackers = [&](ContainerBlueprintNode& outer) { + Blueprints::addCylindricalBeampipe(outer); + auto vtxBarrel = Blueprints::makeBarrel(builder, Blueprints::UngroupedDoubleBarrelLayerVertexSpec, + Blueprints::kBarrelEnvelope, Blueprints::doubleLayerKey); + auto vertex = Blueprints::attachEndcaps(builder, std::move(vtxBarrel), + Blueprints::UngroupedDoubleBarrelLayerVertexSpec, "Vertex"); - auto innerTracker = - Blueprints::makeNestedInnerTracker(builder, std::move(vertex), Blueprints::UngroupedNestedInnerTrackerSpec); - outer.addChild(innerTracker); + auto innerTracker = + Blueprints::makeNestedInnerTracker(builder, std::move(vertex), Blueprints::UngroupedNestedInnerTrackerSpec); + outer.addChild(innerTracker); - auto outerTracker = - Blueprints::makeRegularTracker(builder, Blueprints::UngroupedOuterTrackerSpec, "OuterTracker"); - outer.addChild(outerTracker); + auto outerTracker = + Blueprints::makeRegularTracker(builder, Blueprints::UngroupedOuterTrackerSpec, "OuterTracker"); + outer.addChild(outerTracker); + }; + + if (calo.empty()) { + // No calorimeter face: keep the original purely-radial layout. + auto& outer = root.addCylinderContainer(detName, AxisR); + buildTrackers(outer); + return; + } + + // The calorimeter wraps the tracker, so the top level is a z-stack + // [calo -endcap | central (tracker + calo barrel) | calo +endcap] (see + // MAIA_v0 for details). + auto& world = root.addCylinderContainer(detName, AxisZ); + Blueprints::addCaloEndcap(world, calo, /*positive=*/false); + auto& central = world.addCylinderContainer(detName + "Central", AxisR); + buildTrackers(central); + if (!calo.barrelFaces.empty()) { + Blueprints::addCaloBarrel(central, calo); + } + Blueprints::addCaloEndcap(world, calo, /*positive=*/true); } } // namespace CLD_o2_v07 } // namespace FCCee From a15349acf32018461f81c29f567d110f77d42f94 Mon Sep 17 00:00:00 2001 From: Federico Meloni Date: Mon, 22 Jun 2026 20:08:38 +0200 Subject: [PATCH 10/10] add calo for ILD too --- .../DD4hepBlueprintConstruction.cpp | 116 ++++++++++++------ 1 file changed, 77 insertions(+), 39 deletions(-) diff --git a/k4ActsTracking/src/components/DD4hepBlueprintConstruction.cpp b/k4ActsTracking/src/components/DD4hepBlueprintConstruction.cpp index b52e35f3..1c014d08 100644 --- a/k4ActsTracking/src/components/DD4hepBlueprintConstruction.cpp +++ b/k4ActsTracking/src/components/DD4hepBlueprintConstruction.cpp @@ -698,53 +698,91 @@ namespace MuColl { namespace FCCee { namespace ILD_FCCee_v01 { void populateBlueprint(const std::string& detName, Acts::Experimental::Blueprint& root, - ActsPlugins::DD4hep::BlueprintBuilder& builder, - [[maybe_unused]] const IActsGeoSvc::CaloFaceSurfaces& calo) { - // TODO: integrate the calo face (see MAIA_v0); deferred until validated. - auto& outer = root.addCylinderContainer(detName, AxisR); - - Blueprints::addCylindricalBeampipe(outer); - - auto vtxBarrel = Blueprints::makeBarrel(builder, Blueprints::UngroupedDoubleBarrelLayerVertexSpec, - Blueprints::kBarrelEnvelope, Blueprints::doubleLayerKey); - auto vertex = Blueprints::attachEndcaps(builder, std::move(vtxBarrel), - Blueprints::UngroupedDoubleBarrelLayerVertexSpec, "Vertex"); - - auto innerTrackerBarrel = Blueprints::makeBarrel(builder, Blueprints::UngroupedInnerTrackerSpec); - innerTrackerBarrel->addChild(vertex); - - auto innerTrackerEndcap = Blueprints::attachEndcaps(builder, std::move(innerTrackerBarrel), - Blueprints::UngroupedInnerTrackerSpec, "InnerTrackerEndcap"); - outer.addChild(innerTrackerEndcap); - - // TODO: this is not yet properly working only part of the SET show up in - // the exporte .obj geometry. This usually indicates some issues with the - // AxisDirection, but that would mean that there are different - // AxisDirections in play for the SET geometry - auto set = - Blueprints::makeBarrel(builder, Blueprints::SETSpec, Blueprints::kBarrelEnvelope, Blueprints::doubleLayerKey); - outer.addChild(set); + ActsPlugins::DD4hep::BlueprintBuilder& builder, const IActsGeoSvc::CaloFaceSurfaces& calo) { + // Build the tracker detectors as radial children of the supplied + // container. + auto buildTrackers = [&](ContainerBlueprintNode& outer) { + Blueprints::addCylindricalBeampipe(outer); + + auto vtxBarrel = Blueprints::makeBarrel(builder, Blueprints::UngroupedDoubleBarrelLayerVertexSpec, + Blueprints::kBarrelEnvelope, Blueprints::doubleLayerKey); + auto vertex = Blueprints::attachEndcaps(builder, std::move(vtxBarrel), + Blueprints::UngroupedDoubleBarrelLayerVertexSpec, "Vertex"); + + auto innerTrackerBarrel = Blueprints::makeBarrel(builder, Blueprints::UngroupedInnerTrackerSpec); + innerTrackerBarrel->addChild(vertex); + + auto innerTrackerEndcap = Blueprints::attachEndcaps( + builder, std::move(innerTrackerBarrel), Blueprints::UngroupedInnerTrackerSpec, "InnerTrackerEndcap"); + outer.addChild(innerTrackerEndcap); + + // TODO: this is not yet properly working only part of the SET show up in + // the exporte .obj geometry. This usually indicates some issues with the + // AxisDirection, but that would mean that there are different + // AxisDirections in play for the SET geometry + auto set = Blueprints::makeBarrel(builder, Blueprints::SETSpec, Blueprints::kBarrelEnvelope, + Blueprints::doubleLayerKey); + outer.addChild(set); + }; + + if (calo.empty()) { + // No calorimeter face: keep the original purely-radial layout. + auto& outer = root.addCylinderContainer(detName, AxisR); + buildTrackers(outer); + return; + } + + // The calorimeter wraps the tracker, so the top level is a z-stack + // [calo -endcap | central (tracker + calo barrel) | calo +endcap] (see + // MAIA_v0 for details). + auto& world = root.addCylinderContainer(detName, AxisZ); + Blueprints::addCaloEndcap(world, calo, /*positive=*/false); + auto& central = world.addCylinderContainer(detName + "Central", AxisR); + buildTrackers(central); + if (!calo.barrelFaces.empty()) { + Blueprints::addCaloBarrel(central, calo); + } + Blueprints::addCaloEndcap(world, calo, /*positive=*/true); } } // namespace ILD_FCCee_v01 namespace ILD_FCCee_v02 { void populateBlueprint(const std::string& detName, Acts::Experimental::Blueprint& root, - ActsPlugins::DD4hep::BlueprintBuilder& builder, - [[maybe_unused]] const IActsGeoSvc::CaloFaceSurfaces& calo) { - // TODO: integrate the calo face (see MAIA_v0); deferred until validated. - auto& outer = root.addCylinderContainer(detName, AxisR); + ActsPlugins::DD4hep::BlueprintBuilder& builder, const IActsGeoSvc::CaloFaceSurfaces& calo) { + // Build the tracker detectors as radial children of the supplied + // container. + auto buildTrackers = [&](ContainerBlueprintNode& outer) { + Blueprints::addCylindricalBeampipe(outer); + auto vtxBarrel = Blueprints::makeBarrel(builder, Blueprints::UngroupedDoubleBarrelLayerVertexSpec, + Blueprints::kBarrelEnvelope, Blueprints::doubleLayerKey); + auto vertex = Blueprints::attachEndcaps(builder, std::move(vtxBarrel), + Blueprints::UngroupedDoubleBarrelLayerVertexSpec, "Vertex"); + + auto innerTracker = + Blueprints::makeNestedInnerTracker(builder, std::move(vertex), Blueprints::UngroupedNestedInnerTrackerSpec); + outer.addChild(innerTracker); - Blueprints::addCylindricalBeampipe(outer); - auto vtxBarrel = Blueprints::makeBarrel(builder, Blueprints::UngroupedDoubleBarrelLayerVertexSpec, - Blueprints::kBarrelEnvelope, Blueprints::doubleLayerKey); - auto vertex = Blueprints::attachEndcaps(builder, std::move(vtxBarrel), - Blueprints::UngroupedDoubleBarrelLayerVertexSpec, "Vertex"); + // TODO: Add SET (see V01 for caveats) + }; - auto innerTracker = - Blueprints::makeNestedInnerTracker(builder, std::move(vertex), Blueprints::UngroupedNestedInnerTrackerSpec); - outer.addChild(innerTracker); + if (calo.empty()) { + // No calorimeter face: keep the original purely-radial layout. + auto& outer = root.addCylinderContainer(detName, AxisR); + buildTrackers(outer); + return; + } - // TODO: Add SET (see V01 for caveats) + // The calorimeter wraps the tracker, so the top level is a z-stack + // [calo -endcap | central (tracker + calo barrel) | calo +endcap] (see + // MAIA_v0 for details). + auto& world = root.addCylinderContainer(detName, AxisZ); + Blueprints::addCaloEndcap(world, calo, /*positive=*/false); + auto& central = world.addCylinderContainer(detName + "Central", AxisR); + buildTrackers(central); + if (!calo.barrelFaces.empty()) { + Blueprints::addCaloBarrel(central, calo); + } + Blueprints::addCaloEndcap(world, calo, /*positive=*/true); } } // namespace ILD_FCCee_v02