diff --git a/.github/workflows/checks.yml b/.github/workflows/checks.yml index edd31d051a4..e84044776a8 100644 --- a/.github/workflows/checks.yml +++ b/.github/workflows/checks.yml @@ -145,7 +145,7 @@ jobs: python-version: '3.12' - name: Install dependencies run: > - pip install -r CI/requirements_fpe_masks.txt + pip install -r CI/fpe_masks/requirements.txt - name: Check run: > CI/check_fpe_masks.py --token ${{ secrets.GITHUB_TOKEN }} diff --git a/.github/workflows/milestone.yml b/.github/workflows/milestone.yml new file mode 100644 index 00000000000..de5ee80356c --- /dev/null +++ b/.github/workflows/milestone.yml @@ -0,0 +1,58 @@ +name: Check PR milestone + +on: + pull_request_target: + types: [milestoned, demilestoned, opened, reopened] + branches: + - main + +jobs: + check_milestone: + if: ${{ github.event.issue.pull_request || github.event.pull_request }} + + runs-on: ubuntu-latest + + steps: + - uses: actions/checkout@v4 + + - name: "Check for milestone on PR" + uses: actions/github-script@v7 + with: + script: | + let milestone = context.payload.pull_request.milestone; + + if(context.payload.action === 'opened' || context.payload.action === 'reopened') { + if(milestone !== null) { + core.notice(`Milestone is ${milestone.title}`); + } else { + const milestones = await github.rest.issues.listMilestones({ + owner: context.repo.owner, + repo: context.repo.repo, + state: "open" + }); + for (const default_milestone of milestones.data) { + if (default_milestone.title === "next") { + core.notice(`No milestone set, setting default milestone: ${default_milestone.title}`); + + await github.rest.issues.update({ + owner: context.repo.owner, + repo: context.repo.repo, + issue_number: context.issue.number, + milestone: default_milestone.number + }); + return; + } + } + + core.warning("Could not find default milestone named 'next'"); + + } + + } + else { + if(milestone !== null) { + core.notice(`Milestone is ${milestone.title}`); + } else { + core.setFailed("No milestone: Please add a version milestone"); + } + } diff --git a/.github/workflows/update-pip-requirements.yml b/.github/workflows/update-pip-requirements.yml new file mode 100755 index 00000000000..444ade7b287 --- /dev/null +++ b/.github/workflows/update-pip-requirements.yml @@ -0,0 +1,76 @@ +name: Update Pip Requirements + +on: + workflow_dispatch: # Allow running on-demand + schedule: + # Runs every Sunday at 1:23 UTC + - cron: '23 1 * * 0' + +jobs: + update-pip-requirements: + # This action checks all monitored (specified in `folder_list` below) requirements.in + # files and generates a new requirements.txt file with pip-compile. If any + # requirements changed, a PR is created with the changes. + runs-on: ubuntu-latest + env: + # This branch will receive updates each time the workflow runs + # It doesn't matter if it's deleted when merged, it'll be re-created + BRANCH_NAME: auto-dependency-upgrades + steps: + - uses: actions/checkout@v4 + + - name: Set up Python + uses: actions/setup-python@v5 + with: + python-version: 3.12 + + - name: Install pip-tools + run: pip install pip-tools + + - name: Compile all requirements.txt + run: | + # Update this list after adding/removing requirements-files + folder_list=( + CI/clang_tidy + CI/fpe_masks + docs + Examples/Python/tests/ + Examples/Scripts + ) + for folder in "${folder_list[@]}"; do + pip-compile "${folder}/requirements.in" > "${folder}/requirements.txt" + done + + - name: Detect changes + id: changes + run: + # This output boolean tells us if the dependencies have actually changed + echo "count=$(git status --porcelain=v1 2>/dev/null | wc -l)" >> $GITHUB_OUTPUT + + - name: Commit & push changes + # Only push if changes exist + if: steps.changes.outputs.count > 0 + run: | + git config user.name github-actions + git config user.email github-actions@github.com + git add . + git commit -m "Weekly Update: Regenerate requirements.txt" + git push -f origin ${{ github.ref_name }}:$BRANCH_NAME + env: + GITHUB_TOKEN: ${{ secrets.GITHUB_TOKEN }} + + - name: Open pull request if needed + if: steps.changes.outputs.count > 0 + env: + GITHUB_TOKEN: ${{ secrets.GITHUB_TOKEN }} + # Only open a PR if the branch is not attached to an existing one + run: | + PR=$(gh pr list --head $BRANCH_NAME --json number -q '.[0].number') + if [ -z $PR ]; then + gh pr create \ + --head $BRANCH_NAME \ + --title "chore: automated python requirements upgrades" \ + --body "Full log: https://github.com/${{ github.repository }}/actions/runs/${{ github.run_id }}" + else + echo "Pull request already exists, won't create a new one." + fi diff --git a/CI/clang_tidy/requirements.txt b/CI/clang_tidy/requirements.txt index 6bd13b68d95..d0d78f5904e 100644 --- a/CI/clang_tidy/requirements.txt +++ b/CI/clang_tidy/requirements.txt @@ -4,41 +4,41 @@ # # pip-compile CI/clang_tidy/requirements.in # -annotated-types==0.6.0 +annotated-types==0.7.0 # via pydantic appdirs==1.4.4 # via fs -codereport==0.3.2 +codereport==0.4.0 # via -r CI/clang_tidy/requirements.in fs==2.4.16 # via codereport -jinja2==3.1.2 +jinja2==3.1.4 # via codereport markdown-it-py==3.0.0 # via rich -markupsafe==2.1.3 +markupsafe==3.0.2 # via jinja2 mdurl==0.1.2 # via markdown-it-py -pydantic==2.5.2 +pydantic==2.9.2 # via -r CI/clang_tidy/requirements.in -pydantic-core==2.14.5 +pydantic-core==2.23.4 # via pydantic -pygments==2.17.2 +pygments==2.18.0 # via # codereport # rich python-slugify==6.1.2 # via codereport -pyyaml==6.0.1 +pyyaml==6.0.2 # via -r CI/clang_tidy/requirements.in -rich==13.7.0 +rich==13.9.4 # via -r CI/clang_tidy/requirements.in six==1.16.0 # via fs text-unidecode==1.3 # via python-slugify -typing-extensions==4.8.0 +typing-extensions==4.12.2 # via # pydantic # pydantic-core diff --git a/CI/requirements_fpe_masks.in b/CI/fpe_masks/requirements.in similarity index 100% rename from CI/requirements_fpe_masks.in rename to CI/fpe_masks/requirements.in diff --git a/CI/fpe_masks/requirements.txt b/CI/fpe_masks/requirements.txt new file mode 100644 index 00000000000..01f945e3d74 --- /dev/null +++ b/CI/fpe_masks/requirements.txt @@ -0,0 +1,62 @@ +# +# This file is autogenerated by pip-compile with Python 3.12 +# by the following command: +# +# pip-compile CI/fpe_masks/requirements.in +# +aiohappyeyeballs==2.4.3 + # via aiohttp +aiohttp==3.11.6 + # via -r CI/fpe_masks/requirements.in +aiosignal==1.3.1 + # via aiohttp +attrs==24.2.0 + # via aiohttp +cffi==1.17.1 + # via cryptography +click==8.1.7 + # via typer +cryptography==43.0.3 + # via pyjwt +frozenlist==1.5.0 + # via + # aiohttp + # aiosignal +gidgethub==5.3.0 + # via -r CI/fpe_masks/requirements.in +idna==3.10 + # via yarl +markdown-it-py==3.0.0 + # via rich +mdurl==0.1.2 + # via markdown-it-py +multidict==6.1.0 + # via + # aiohttp + # yarl +propcache==0.2.0 + # via + # aiohttp + # yarl +pycparser==2.22 + # via cffi +pygments==2.18.0 + # via rich +pyjwt[crypto]==2.10.0 + # via + # gidgethub + # pyjwt +rich==13.9.4 + # via + # -r CI/fpe_masks/requirements.in + # typer +shellingham==1.5.4 + # via typer +typer==0.13.1 + # via -r CI/fpe_masks/requirements.in +typing-extensions==4.12.2 + # via typer +uritemplate==4.1.1 + # via gidgethub +yarl==1.17.2 + # via aiohttp diff --git a/CI/requirements_fpe_masks.txt b/CI/requirements_fpe_masks.txt deleted file mode 100644 index 07e86bd58d5..00000000000 --- a/CI/requirements_fpe_masks.txt +++ /dev/null @@ -1,50 +0,0 @@ -# -# This file is autogenerated by pip-compile with Python 3.12 -# by the following command: -# -# pip-compile CI/requirements_fpe_masks.in -# -aiohttp==3.9.1 - # via -r CI/requirements_fpe_masks.in -aiosignal==1.3.1 - # via aiohttp -attrs==23.1.0 - # via aiohttp -cffi==1.15.1 - # via cryptography -click==8.1.4 - # via typer -cryptography==41.0.1 - # via pyjwt -frozenlist==1.4.0 - # via - # aiohttp - # aiosignal -gidgethub==5.3.0 - # via -r CI/requirements_fpe_masks.in -idna==3.4 - # via yarl -markdown-it-py==3.0.0 - # via rich -mdurl==0.1.2 - # via markdown-it-py -multidict==6.0.4 - # via - # aiohttp - # yarl -pycparser==2.21 - # via cffi -pygments==2.15.1 - # via rich -pyjwt[crypto]==2.7.0 - # via gidgethub -rich==13.4.2 - # via -r CI/requirements_fpe_masks.in -typer==0.9.0 - # via -r CI/requirements_fpe_masks.in -typing-extensions==4.7.1 - # via typer -uritemplate==4.1.1 - # via gidgethub -yarl==1.9.2 - # via aiohttp diff --git a/Core/include/Acts/EventData/detail/TrackParametersUtils.hpp b/Core/include/Acts/EventData/detail/TrackParametersUtils.hpp new file mode 100644 index 00000000000..2efcef0e763 --- /dev/null +++ b/Core/include/Acts/EventData/detail/TrackParametersUtils.hpp @@ -0,0 +1,49 @@ +// This file is part of the ACTS project. +// +// Copyright (C) 2016 CERN for the benefit of the ACTS project +// +// This Source Code Form is subject to the terms of the Mozilla Public +// License, v. 2.0. If a copy of the MPL was not distributed with this +// file, You can obtain one at https://mozilla.org/MPL/2.0/. + +#pragma once + +#include "Acts/EventData/GenericBoundTrackParameters.hpp" +#include "Acts/EventData/TrackParametersConcept.hpp" + +namespace Acts::detail { + +/// @brief Shorthand for Bound or Free track parameters +template +concept isBoundOrFreeTrackParams = + Acts::FreeTrackParametersConcept || + Acts::BoundTrackParametersConcept; + +/// @brief Shorthand for GenericBoundTrackParameters +template +concept isGenericBoundTrackParams = + std::same_as>; + +/// @brief Concept that restricts the type of the +/// accumulation grid cell +template +concept TrackParamsGrid = requires { + typename grid_t::value_type::first_type; + typename grid_t::value_type::second_type; + + requires isBoundOrFreeTrackParams< + typename grid_t::value_type::first_type::element_type>; + requires isBoundOrFreeTrackParams< + typename grid_t::value_type::second_type::element_type>; + + requires requires(typename grid_t::value_type val) { + { + val.first + } -> std::same_as< + std::shared_ptr&>; + { val.second } -> std::same_as; + }; +}; + +} // namespace Acts::detail diff --git a/Core/include/Acts/Seeding/EstimateTrackParamsFromSeed.hpp b/Core/include/Acts/Seeding/EstimateTrackParamsFromSeed.hpp index c59a55ad533..e4618296b50 100644 --- a/Core/include/Acts/Seeding/EstimateTrackParamsFromSeed.hpp +++ b/Core/include/Acts/Seeding/EstimateTrackParamsFromSeed.hpp @@ -8,17 +8,14 @@ #pragma once +#include "Acts/Definitions/Algebra.hpp" #include "Acts/Definitions/TrackParametrization.hpp" #include "Acts/Definitions/Units.hpp" #include "Acts/Geometry/GeometryContext.hpp" #include "Acts/Surfaces/Surface.hpp" -#include "Acts/Utilities/Logger.hpp" -#include "Acts/Utilities/MathHelpers.hpp" #include "Acts/Utilities/Zip.hpp" #include -#include -#include #include #include @@ -27,17 +24,16 @@ namespace Acts { /// Estimate the full track parameters from three space points /// /// This method is based on the conformal map transformation. It estimates the -/// full free track parameters, i.e. (x, y, z, t, dx, dy, dz, q/p) at the -/// bottom space point. The bottom space is assumed to be the first element -/// in the range defined by the iterators. The magnetic field (which might be -/// along any direction) is also necessary for the momentum estimation. +/// full free track parameters, i.e. (x, y, z, t, dx, dy, dz, q/p) at the bottom +/// space point. The bottom space is assumed to be the first element in the +/// range defined by the iterators. The magnetic field (which might be along any +/// direction) is also necessary for the momentum estimation. /// /// This is a purely spatial estimation, i.e. the time parameter will be set to /// 0. /// -/// It resembles the method used in ATLAS for the track parameters -/// estimated from seed, i.e. the function InDet::SiTrackMaker_xk::getAtaPlane -/// here: +/// It resembles the method used in ATLAS for the track parameters estimated +/// from seed, i.e. the function InDet::SiTrackMaker_xk::getAtaPlane here: /// https://acode-browser.usatlas.bnl.gov/lxr/source/athena/InnerDetector/InDetRecTools/SiTrackMakerTool_xk/src/SiTrackMaker_xk.cxx /// /// @tparam spacepoint_iterator_t The type of space point iterator @@ -55,14 +51,13 @@ FreeVector estimateTrackParamsFromSeed(const Vector3& sp0, const Vector3& sp1, /// Estimate the full track parameters from three space points /// /// This method is based on the conformal map transformation. It estimates the -/// full free track parameters, i.e. (x, y, z, t, dx, dy, dz, q/p) at the -/// bottom space point. The bottom space is assumed to be the first element -/// in the range defined by the iterators. The magnetic field (which might be -/// along any direction) is also necessary for the momentum estimation. +/// full free track parameters, i.e. (x, y, z, t, dx, dy, dz, q/p) at the bottom +/// space point. The bottom space is assumed to be the first element in the +/// range defined by the iterators. The magnetic field (which might be along any +/// direction) is also necessary for the momentum estimation. /// -/// It resembles the method used in ATLAS for the track parameters -/// estimated from seed, i.e. the function InDet::SiTrackMaker_xk::getAtaPlane -/// here: +/// It resembles the method used in ATLAS for the track parameters estimated +/// from seed, i.e. the function InDet::SiTrackMaker_xk::getAtaPlane here: /// https://acode-browser.usatlas.bnl.gov/lxr/source/athena/InnerDetector/InDetRecTools/SiTrackMakerTool_xk/src/SiTrackMaker_xk.cxx /// /// @tparam spacepoint_iterator_t The type of space point iterator @@ -107,150 +102,53 @@ FreeVector estimateTrackParamsFromSeed(spacepoint_range_t spRange, /// This method is based on the conformal map transformation. It estimates the /// full bound track parameters, i.e. (loc0, loc1, phi, theta, q/p, t) at the /// bottom space point. The bottom space is assumed to be the first element -/// in the range defined by the iterators. It must lie on the surface -/// provided for the representation of the bound track parameters. The magnetic -/// field (which might be along any direction) is also necessary for the -/// momentum estimation. +/// in the range defined by the iterators. It must lie on the surface provided +/// for the representation of the bound track parameters. The magnetic field +/// (which might be along any direction) is also necessary for the momentum +/// estimation. /// -/// It resembles the method used in ATLAS for the track parameters -/// estimated from seed, i.e. the function InDet::SiTrackMaker_xk::getAtaPlane -/// here: +/// It resembles the method used in ATLAS for the track parameters estimated +/// from seed, i.e. the function InDet::SiTrackMaker_xk::getAtaPlane here: /// https://acode-browser.usatlas.bnl.gov/lxr/source/athena/InnerDetector/InDetRecTools/SiTrackMakerTool_xk/src/SiTrackMaker_xk.cxx /// /// @tparam spacepoint_iterator_t The type of space point iterator /// /// @param gctx is the geometry context -/// @param spBegin is the begin iterator for the space points -/// @param spEnd is the end iterator for the space points +/// @param spRange is the range of space points /// @param surface is the surface of the bottom space point. The estimated bound /// track parameters will be represented also at this surface /// @param bField is the magnetic field vector -/// @param logger A logger instance /// -/// @return optional bound parameters -template -std::optional estimateTrackParamsFromSeed( - const GeometryContext& gctx, spacepoint_iterator_t spBegin, - spacepoint_iterator_t spEnd, const Surface& surface, const Vector3& bField, - const Acts::Logger& logger = getDummyLogger()) { - // Check the number of provided space points - std::size_t numSP = std::distance(spBegin, spEnd); - if (numSP != 3) { - ACTS_ERROR("There should be exactly three space points provided."); - return std::nullopt; - } - - // The global positions of the bottom, middle and space points - std::array spGlobalPositions = {Vector3::Zero(), Vector3::Zero(), - Vector3::Zero()}; - std::array, 3> spGlobalTimes = { - std::nullopt, std::nullopt, std::nullopt}; - // The first, second and third space point are assumed to be bottom, middle - // and top space point, respectively - for (std::size_t isp = 0; isp < 3; ++isp) { - spacepoint_iterator_t it = std::next(spBegin, isp); - if (*it == nullptr) { - ACTS_ERROR("Empty space point found. This should not happen."); - return std::nullopt; - } - const auto& sp = *it; - spGlobalPositions[isp] = Vector3(sp->x(), sp->y(), sp->z()); - spGlobalTimes[isp] = sp->t(); - } - - // Define a new coordinate frame with its origin at the bottom space point, z - // axis long the magnetic field direction and y axis perpendicular to vector - // from the bottom to middle space point. Hence, the projection of the middle - // space point on the transverse plane will be located at the x axis of the - // new frame. - Vector3 relVec = spGlobalPositions[1] - spGlobalPositions[0]; - Vector3 newZAxis = bField.normalized(); - Vector3 newYAxis = newZAxis.cross(relVec).normalized(); - Vector3 newXAxis = newYAxis.cross(newZAxis); - RotationMatrix3 rotation; - rotation.col(0) = newXAxis; - rotation.col(1) = newYAxis; - rotation.col(2) = newZAxis; - // The center of the new frame is at the bottom space point - Translation3 trans(spGlobalPositions[0]); - // The transform which constructs the new frame - Transform3 transform(trans * rotation); - - // The coordinate of the middle and top space point in the new frame - Vector3 local1 = transform.inverse() * spGlobalPositions[1]; - Vector3 local2 = transform.inverse() * spGlobalPositions[2]; - - // In the new frame the bottom sp is at the origin, while the middle - // sp in along the x axis. As such, the x-coordinate of the circle is - // at: x-middle / 2. - // The y coordinate can be found by using the straight line passing - // between the mid point between the middle and top sp and perpendicular to - // the line connecting them - Vector2 circleCenter; - circleCenter(0) = 0.5 * local1(0); +/// @return bound parameters +template +Result estimateTrackParamsFromSeed(const GeometryContext& gctx, + spacepoint_range_t spRange, + const Surface& surface, + const Vector3& bField) { + FreeVector freeParams = estimateTrackParamsFromSeed(spRange, bField); - ActsScalar deltaX21 = local2(0) - local1(0); - ActsScalar sumX21 = local2(0) + local1(0); - // straight line connecting the two points - // y = a * x + c (we don't care about c right now) - // we simply need the slope - // we compute 1./a since this is what we need for the following computation - ActsScalar ia = deltaX21 / local2(1); - // Perpendicular line is then y = -1/a *x + b - // we can evaluate b given we know a already by imposing - // the line passes through P = (0.5 * (x2 + x1), 0.5 * y2) - ActsScalar b = 0.5 * (local2(1) + ia * sumX21); - circleCenter(1) = -ia * circleCenter(0) + b; - // Radius is a signed distance between circleCenter and first sp, which is at - // (0, 0) in the new frame. Sign depends on the slope a (positive vs negative) - int sign = ia > 0 ? -1 : 1; - const ActsScalar R = circleCenter.norm(); - ActsScalar invTanTheta = - local2.z() / (2 * R * std::asin(local2.head<2>().norm() / (2 * R))); - // The momentum direction in the new frame (the center of the circle has the - // coordinate (-1.*A/(2*B), 1./(2*B))) - ActsScalar A = -circleCenter(0) / circleCenter(1); - Vector3 transDirection(1., A, fastHypot(1, A) * invTanTheta); - // Transform it back to the original frame - Vector3 direction = rotation * transDirection.normalized(); + const auto* sp0 = *spRange.begin(); + Vector3 origin = Vector3(sp0->x(), sp0->y(), sp0->z()); + Vector3 direction = freeParams.segment<3>(eFreeDir0); - // Initialize the bound parameters vector BoundVector params = BoundVector::Zero(); - - // The estimated phi and theta params[eBoundPhi] = VectorHelpers::phi(direction); params[eBoundTheta] = VectorHelpers::theta(direction); + params[eBoundQOverP] = freeParams[eFreeQOverP]; // Transform the bottom space point to local coordinates of the provided // surface - auto lpResult = surface.globalToLocal(gctx, spGlobalPositions[0], direction); + auto lpResult = surface.globalToLocal(gctx, origin, direction); if (!lpResult.ok()) { - ACTS_ERROR( - "Global to local transformation did not succeed. Please make sure the " - "bottom space point lies on the provided surface."); - return std::nullopt; + return Result::failure(lpResult.error()); } Vector2 bottomLocalPos = lpResult.value(); // The estimated loc0 and loc1 params[eBoundLoc0] = bottomLocalPos.x(); params[eBoundLoc1] = bottomLocalPos.y(); - params[eBoundTime] = spGlobalTimes[0].value_or(0.); - - ActsScalar bFieldStrength = bField.norm(); - // The estimated q/pt in [GeV/c]^-1 (note that the pt is the projection of - // momentum on the transverse plane of the new frame) - ActsScalar qOverPt = sign / (bFieldStrength * R); - // The estimated q/p in [GeV/c]^-1 - params[eBoundQOverP] = qOverPt / fastHypot(1., invTanTheta); + params[eBoundTime] = sp0->t().value_or(0); - if (params.hasNaN()) { - ACTS_ERROR( - "The NaN value exists at the estimated track parameters from seed with" - << "\nbottom sp: " << spGlobalPositions[0] << "\nmiddle sp: " - << spGlobalPositions[1] << "\ntop sp: " << spGlobalPositions[2]); - return std::nullopt; - } - return params; + return Result::success(params); } /// Configuration for the estimation of the covariance matrix of the track diff --git a/Core/include/Acts/TrackFinding/TrackParamsLookupAccumulator.hpp b/Core/include/Acts/TrackFinding/TrackParamsLookupAccumulator.hpp new file mode 100644 index 00000000000..a1e1173335d --- /dev/null +++ b/Core/include/Acts/TrackFinding/TrackParamsLookupAccumulator.hpp @@ -0,0 +1,179 @@ +// This file is part of the ACTS project. +// +// Copyright (C) 2016 CERN for the benefit of the ACTS project +// +// This Source Code Form is subject to the terms of the Mozilla Public +// License, v. 2.0. If a copy of the MPL was not distributed with this +// file, You can obtain one at https://mozilla.org/MPL/2.0/. + +#pragma once + +#include "Acts/EventData/detail/TrackParametersUtils.hpp" +#include "Acts/Geometry/GeometryContext.hpp" + +#include +#include +#include +#include +#include + +namespace Acts { + +/// @brief Class to accumulate and average track lookup tables +/// +/// @tparam Grid type for track parameters accumulation +/// +/// This class is used to accumulate track parameters in +/// reference layer grids and average them to create a lookup +/// table for track parameter estimation in seeding +/// +/// @note Geometry context is left to be handled by the user +/// outside of accumulation +template +class TrackParamsLookupAccumulator { + public: + using LookupGrid = grid_t; + using TrackParameters = typename std::pointer_traits< + typename grid_t::value_type::first_type>::element_type; + + /// @brief Constructor + explicit TrackParamsLookupAccumulator(grid_t grid) + : m_grid(std::move(grid)) {} + + /// @brief Add track parameters to the accumulator + /// + /// @param ipTrackParameters Track parameters at the IP + /// @param refTrackParameters Track parameters at the reference layer + /// @param position Local position of the track hit on the reference layer + void addTrack(const TrackParameters& ipTrackParameters, + const TrackParameters& refTrackParameters, + const Vector2& position) { + std::lock_guard lock(m_gridMutex); + + auto bin = m_grid.localBinsFromPosition(position); + + if (m_countGrid[bin] == 0) { + m_grid.atLocalBins(bin).first = + std::make_shared(ipTrackParameters); + m_grid.atLocalBins(bin).second = + std::make_shared(refTrackParameters); + + m_countGrid.at(bin)++; + return; + } + + *m_grid.atLocalBins(bin).first = + addTrackParameters(*m_grid.atLocalBins(bin).first, ipTrackParameters); + *m_grid.atLocalBins(bin).second = + addTrackParameters(*m_grid.atLocalBins(bin).second, refTrackParameters); + m_countGrid.at(bin)++; + } + + /// @brief Finalize the lookup table + /// + /// @return Grid with the bin track parameters averaged + LookupGrid finalizeLookup() { + auto meanTrack = [&](const TrackParameters& track, std::size_t count) { + if constexpr (detail::isGenericBoundTrackParams) { + Acts::GeometryContext gctx; + + auto res = TrackParameters::create( + track.referenceSurface().getSharedPtr(), gctx, + track.fourPosition(gctx) / count, track.momentum().normalized(), + count * track.charge() / track.momentum().norm(), + track.covariance(), track.particleHypothesis()); + + if (!res.ok()) { + throw std::invalid_argument("Bound track grid finalization failed"); + } + return res.value(); + } else { + return TrackParameters(track.fourPosition() / count, + track.momentum().normalized(), + count * track.charge() / track.momentum().norm(), + track.covariance(), track.particleHypothesis()); + } + }; + + for (auto [bin, count] : m_countGrid) { + if (count == 0) { + continue; + } + *m_grid.atLocalBins(bin).first = + meanTrack(*m_grid.atLocalBins(bin).first, count); + *m_grid.atLocalBins(bin).second = + meanTrack(*m_grid.atLocalBins(bin).second, count); + } + + return m_grid; + } + + private: + /// @brief Add two track parameters + /// + /// @param a First track parameter in the sum + /// @param b Second track parameter in the sum + /// + /// @return Sum of track parameters a + b + /// + /// @note Covariances of the track parameters + /// are not added and instead assumed to be + /// generated by the same random process for + /// both a and b, making its averaging redundant + TrackParameters addTrackParameters(const TrackParameters& a, + const TrackParameters& b) { + if (a.particleHypothesis() != b.particleHypothesis()) { + throw std::invalid_argument( + "Cannot accumulate track parameters with different particle " + "hypotheses"); + } + if (a.charge() != b.charge()) { + throw std::invalid_argument( + "Cannot accumulate track parameters with different charges"); + } + if constexpr (detail::isGenericBoundTrackParams) { + if (a.referenceSurface() != b.referenceSurface()) { + throw std::invalid_argument( + "Cannot accumulate bound track parameters with different reference " + "surfaces"); + } + } + + Acts::Vector3 momentum = a.momentum() + b.momentum(); + + // Assume track parameters being i.i.d. + if constexpr (detail::isGenericBoundTrackParams) { + Acts::GeometryContext gctx; + + Acts::Vector4 fourPosition = a.fourPosition(gctx) + b.fourPosition(gctx); + + auto res = TrackParameters::create( + a.referenceSurface().getSharedPtr(), gctx, fourPosition, + momentum.normalized(), a.charge() / momentum.norm(), a.covariance(), + a.particleHypothesis()); + + if (!res.ok()) { + throw std::runtime_error("Invalid bound track parameters"); + } + return res.value(); + } else { + Acts::Vector4 fourPosition = a.fourPosition() + b.fourPosition(); + return TrackParameters(fourPosition, momentum.normalized(), + a.charge() / momentum.norm(), a.covariance(), + a.particleHypothesis()); + } + } + + /// Grids to accumulate IP and reference + /// layer track parameters + LookupGrid m_grid; + + /// Mutex for protecting grid access + std::mutex m_gridMutex; + + /// Map to keep the accumulation count + /// in the occupied grid bins + std::map, std::size_t> m_countGrid; +}; + +} // namespace Acts diff --git a/Core/include/Acts/TrackFitting/GlobalChiSquareFitter.hpp b/Core/include/Acts/TrackFitting/GlobalChiSquareFitter.hpp index 8442071a23f..a7b28d5e8f3 100644 --- a/Core/include/Acts/TrackFitting/GlobalChiSquareFitter.hpp +++ b/Core/include/Acts/TrackFitting/GlobalChiSquareFitter.hpp @@ -1408,9 +1408,9 @@ class Gx2Fitter { if ((gx2fOptions.relChi2changeCutOff != 0) && (nUpdate > 0) && (std::abs(extendedSystem.chi2() / oldChi2sum - 1) < gx2fOptions.relChi2changeCutOff)) { - ACTS_INFO("Abort with relChi2changeCutOff after " - << nUpdate + 1 << "/" << gx2fOptions.nUpdateMax - << " iterations."); + ACTS_DEBUG("Abort with relChi2changeCutOff after " + << nUpdate + 1 << "/" << gx2fOptions.nUpdateMax + << " iterations."); updateGx2fCovarianceParams(fullCovariancePredicted, extendedSystem); break; } diff --git a/Core/src/Seeding/EstimateTrackParamsFromSeed.cpp b/Core/src/Seeding/EstimateTrackParamsFromSeed.cpp index 055691bbb2f..b071d5dd7e6 100644 --- a/Core/src/Seeding/EstimateTrackParamsFromSeed.cpp +++ b/Core/src/Seeding/EstimateTrackParamsFromSeed.cpp @@ -9,6 +9,7 @@ #include "Acts/Seeding/EstimateTrackParamsFromSeed.hpp" #include "Acts/Definitions/TrackParametrization.hpp" +#include "Acts/Utilities/MathHelpers.hpp" #include diff --git a/Examples/Algorithms/TrackFinding/CMakeLists.txt b/Examples/Algorithms/TrackFinding/CMakeLists.txt index 38781494068..33b924e35aa 100644 --- a/Examples/Algorithms/TrackFinding/CMakeLists.txt +++ b/Examples/Algorithms/TrackFinding/CMakeLists.txt @@ -10,6 +10,7 @@ add_library( src/TrackParamsEstimationAlgorithm.cpp src/MuonHoughSeeder.cpp src/GbtsSeedingAlgorithm.cpp + src/TrackParamsLookupEstimation.cpp ) target_include_directories( diff --git a/Examples/Algorithms/TrackFinding/include/ActsExamples/TrackFinding/ITrackParamsLookupReader.hpp b/Examples/Algorithms/TrackFinding/include/ActsExamples/TrackFinding/ITrackParamsLookupReader.hpp new file mode 100644 index 00000000000..5777f90ef18 --- /dev/null +++ b/Examples/Algorithms/TrackFinding/include/ActsExamples/TrackFinding/ITrackParamsLookupReader.hpp @@ -0,0 +1,27 @@ +// This file is part of the ACTS project. +// +// Copyright (C) 2016 CERN for the benefit of the ACTS project +// +// This Source Code Form is subject to the terms of the Mozilla Public +// License, v. 2.0. If a copy of the MPL was not distributed with this +// file, You can obtain one at https://mozilla.org/MPL/2.0/. + +#pragma once + +#include "ActsExamples/TrackFinding/TrackParamsLookupTable.hpp" + +namespace ActsExamples { + +/// @brief Interface for reading track parameter lookup tables +class ITrackParamsLookupReader { + public: + /// Virtual Destructor + virtual ~ITrackParamsLookupReader() = default; + + /// Reader method + /// + /// @param path the path to the file to read + virtual TrackParamsLookup readLookup(const std::string& path) const = 0; +}; + +} // namespace ActsExamples diff --git a/Examples/Algorithms/TrackFinding/include/ActsExamples/TrackFinding/ITrackParamsLookupWriter.hpp b/Examples/Algorithms/TrackFinding/include/ActsExamples/TrackFinding/ITrackParamsLookupWriter.hpp new file mode 100644 index 00000000000..9cddadae062 --- /dev/null +++ b/Examples/Algorithms/TrackFinding/include/ActsExamples/TrackFinding/ITrackParamsLookupWriter.hpp @@ -0,0 +1,27 @@ +// This file is part of the ACTS project. +// +// Copyright (C) 2016 CERN for the benefit of the ACTS project +// +// This Source Code Form is subject to the terms of the Mozilla Public +// License, v. 2.0. If a copy of the MPL was not distributed with this +// file, You can obtain one at https://mozilla.org/MPL/2.0/. + +#pragma once + +#include "ActsExamples/TrackFinding/TrackParamsLookupTable.hpp" + +namespace ActsExamples { + +/// @brief Interface for writing track parameter lookup tables +class ITrackParamsLookupWriter { + public: + /// Virtual Destructor + virtual ~ITrackParamsLookupWriter() = default; + + /// Writer method + /// + /// @param lookup track lookup to write + virtual void writeLookup(const TrackParamsLookup& lookup) const = 0; +}; + +} // namespace ActsExamples diff --git a/Examples/Algorithms/TrackFinding/include/ActsExamples/TrackFinding/TrackParamsLookupEstimation.hpp b/Examples/Algorithms/TrackFinding/include/ActsExamples/TrackFinding/TrackParamsLookupEstimation.hpp new file mode 100644 index 00000000000..e50991cd944 --- /dev/null +++ b/Examples/Algorithms/TrackFinding/include/ActsExamples/TrackFinding/TrackParamsLookupEstimation.hpp @@ -0,0 +1,78 @@ +// This file is part of the ACTS project. +// +// Copyright (C) 2016 CERN for the benefit of the ACTS project +// +// This Source Code Form is subject to the terms of the Mozilla Public +// License, v. 2.0. If a copy of the MPL was not distributed with this +// file, You can obtain one at https://mozilla.org/MPL/2.0/. + +#pragma once + +#include "Acts/TrackFinding/TrackParamsLookupAccumulator.hpp" +#include "ActsExamples/EventData/SimHit.hpp" +#include "ActsExamples/EventData/SimParticle.hpp" +#include "ActsExamples/Framework/DataHandle.hpp" +#include "ActsExamples/Framework/IAlgorithm.hpp" +#include "ActsExamples/TrackFinding/ITrackParamsLookupWriter.hpp" + +#include + +namespace ActsExamples { + +/// @brief Algorithm to estimate track parameters lookup tables +/// +/// This algorithm is used to estimate track parameters lookup tables +/// for track parameter estimation in seeding. The algorithm imposes +/// grids onto the reference tracking layers and accumulates track +/// parameters in the grid bins. The track parameters are then averaged +/// to create a lookup table for track parameter estimation in seeding. +class TrackParamsLookupEstimation : public IAlgorithm { + public: + using TrackParamsLookupAccumulator = + Acts::TrackParamsLookupAccumulator; + + /// @brief Nested configuration struct + struct Config { + /// Reference tracking layers + std::unordered_map + refLayers; + /// Binning of the grid to be emposed + /// onto the reference layers + std::pair bins; + /// Input SimHit container + std::string inputHits = "InputHits"; + /// Input SimParticle container + std::string inputParticles = "InputParticles"; + /// Track lookup writers + std::vector> + trackLookupGridWriters{}; + }; + + /// @brief Constructor + TrackParamsLookupEstimation(const Config& config, Acts::Logging::Level level); + + /// @brief The execute method + ProcessCode execute(const AlgorithmContext& ctx) const override; + + ProcessCode finalize() override; + + /// Get readonly access to the config parameters + const Config& config() const { return m_cfg; } + + private: + /// Configuration + Config m_cfg; + + /// Input data handles + ReadDataHandle m_inputParticles{this, + "InputSimParticles"}; + + ReadDataHandle m_inputSimHits{this, "InputSimHits"}; + + /// Accumulators for the track parameters + std::unordered_map> + m_accumulators; +}; + +} // namespace ActsExamples diff --git a/Examples/Algorithms/TrackFinding/include/ActsExamples/TrackFinding/TrackParamsLookupTable.hpp b/Examples/Algorithms/TrackFinding/include/ActsExamples/TrackFinding/TrackParamsLookupTable.hpp new file mode 100644 index 00000000000..5b9c1fabd83 --- /dev/null +++ b/Examples/Algorithms/TrackFinding/include/ActsExamples/TrackFinding/TrackParamsLookupTable.hpp @@ -0,0 +1,44 @@ +// This file is part of the ACTS project. +// +// Copyright (C) 2016 CERN for the benefit of the ACTS project +// +// This Source Code Form is subject to the terms of the Mozilla Public +// License, v. 2.0. If a copy of the MPL was not distributed with this +// file, You can obtain one at https://mozilla.org/MPL/2.0/. + +#pragma once + +#include "Acts/EventData/TrackParameters.hpp" +#include "Acts/Utilities/Grid.hpp" +#include "Acts/Utilities/GridAxisGenerators.hpp" + +#include +#include + +namespace ActsExamples { + +using TrackParamsLookupPair = + std::pair, + std::shared_ptr>; + +/// @brief Track parameters lookup table axis used +/// in the track estimation algorithm +using TrackParamsLookupAxis = + Acts::Axis; + +/// @brief Track parameters lookup table axis generator +/// used in the track estimation algorithm +using TrackParamsLookupAxisGen = Acts::GridAxisGenerators::EqOpenEqOpen; + +/// @brief Lookup grid for track parameters estimation +/// in a given layer +using TrackParamsLookupGrid = + Acts::Grid; + +/// @brief Lookup table for track parameters estimation +/// in the track estimation algorithm +using TrackParamsLookup = + std::unordered_map; + +} // namespace ActsExamples diff --git a/Examples/Algorithms/TrackFinding/src/TrackParamsEstimationAlgorithm.cpp b/Examples/Algorithms/TrackFinding/src/TrackParamsEstimationAlgorithm.cpp index dfb0269a56d..f2711ce1d1a 100644 --- a/Examples/Algorithms/TrackFinding/src/TrackParamsEstimationAlgorithm.cpp +++ b/Examples/Algorithms/TrackFinding/src/TrackParamsEstimationAlgorithm.cpp @@ -21,12 +21,10 @@ #include "ActsExamples/EventData/Track.hpp" #include "ActsExamples/Framework/AlgorithmContext.hpp" -#include #include #include #include #include -#include #include #include @@ -119,16 +117,14 @@ ProcessCode TrackParamsEstimationAlgorithm::execute( } // Estimate the track parameters from seed - auto optParams = Acts::estimateTrackParamsFromSeed( - ctx.geoContext, seed.sp().begin(), seed.sp().end(), *surface, field, - logger()); - if (!optParams.has_value()) { - ACTS_WARNING("Estimation of track parameters for seed " << iseed - << " failed."); + const auto paramsResult = Acts::estimateTrackParamsFromSeed( + ctx.geoContext, seed.sp(), *surface, field); + if (!paramsResult.ok()) { + ACTS_WARNING("Skip track because param estimation failed " + << paramsResult.error()); continue; } - - const auto& params = optParams.value(); + const auto& params = *paramsResult; Acts::EstimateTrackParamCovarianceConfig config{ .initialSigmas = diff --git a/Examples/Algorithms/TrackFinding/src/TrackParamsLookupEstimation.cpp b/Examples/Algorithms/TrackFinding/src/TrackParamsLookupEstimation.cpp new file mode 100644 index 00000000000..4e4b283ddf9 --- /dev/null +++ b/Examples/Algorithms/TrackFinding/src/TrackParamsLookupEstimation.cpp @@ -0,0 +1,104 @@ +// This file is part of the ACTS project. +// +// Copyright (C) 2016 CERN for the benefit of the ACTS project +// +// This Source Code Form is subject to the terms of the Mozilla Public +// License, v. 2.0. If a copy of the MPL was not distributed with this +// file, You can obtain one at https://mozilla.org/MPL/2.0/. + +#include "ActsExamples/TrackFinding/TrackParamsLookupEstimation.hpp" + +#include "Acts/Surfaces/RectangleBounds.hpp" +#include "ActsExamples/Framework/ProcessCode.hpp" + +ActsExamples::TrackParamsLookupEstimation::TrackParamsLookupEstimation( + const Config& config, Acts::Logging::Level level) + : IAlgorithm("TrackParamsLookupEstimation", level), m_cfg(config) { + // Iterate over the reference layers and create + // track parameter accumulators + for (const auto& [geoId, refSurface] : m_cfg.refLayers) { + // Get bounds to construct the accumulator grid + auto bounds = + dynamic_cast(&refSurface->bounds()); + + if (bounds == nullptr) { + throw std::invalid_argument("Only rectangle bounds supported"); + } + if (refSurface->type() != Acts::Surface::SurfaceType::Plane) { + throw std::invalid_argument("Only plane surfaces supported"); + } + + // Initialize the accumulator grid + auto halfX = bounds->halfLengthX(); + auto halfY = bounds->halfLengthY(); + + TrackParamsLookupAxisGen axisGen{ + {-halfX, halfX}, m_cfg.bins.first, {-halfY, halfY}, m_cfg.bins.second}; + + // Each reference layer has its own accumulator + m_accumulators[geoId] = std::make_unique( + TrackParamsLookupGrid(axisGen())); + } + + m_inputParticles.initialize(m_cfg.inputParticles); + m_inputSimHits.initialize(m_cfg.inputHits); +} + +ActsExamples::ProcessCode +ActsExamples::TrackParamsLookupEstimation::finalize() { + // Finiliaze the lookup tables and write them + ActsExamples::TrackParamsLookup lookup; + for (auto& [id, acc] : m_accumulators) { + lookup.insert({id, acc->finalizeLookup()}); + } + for (const auto& writer : m_cfg.trackLookupGridWriters) { + writer->writeLookup(lookup); + } + + return ActsExamples::ProcessCode::SUCCESS; +}; + +ActsExamples::ProcessCode ActsExamples::TrackParamsLookupEstimation::execute( + const ActsExamples::AlgorithmContext& ctx) const { + // Get the particles and hits + const auto& particles = m_inputParticles(ctx); + const auto& hits = m_inputSimHits(ctx); + + // Iterate over the reference layer hits and + // accumulate the track parameters + for (const auto& [geoId, refSurface] : m_cfg.refLayers) { + // Get reference layer hits + auto refLayerHits = hits.equal_range(geoId); + + for (auto hit = refLayerHits.first; hit != refLayerHits.second; ++hit) { + // Get the corresponding particle + const auto& id = hit->particleId(); + const auto& particle = particles.find(id); + + if (particle == particles.end()) { + throw std::invalid_argument("Particle not found"); + } + + // Hit stores the reference layer parameters + auto refLayerPars = Acts::CurvilinearTrackParameters( + hit->fourPosition(), hit->direction(), particle->qOverP(), + std::nullopt, particle->hypothesis()); + + // Particle stores the IP parameters + auto ipPars = Acts::CurvilinearTrackParameters( + particle->fourPosition(), particle->direction(), particle->qOverP(), + std::nullopt, particle->hypothesis()); + + // Get the local position of the hit + auto localPos = refSurface + ->globalToLocal(ctx.geoContext, hit->position(), + Acts::Vector3{0, 1, 0}) + .value(); + + // Add the track parameters to the accumulator grid + m_accumulators.at(geoId)->addTrack(ipPars, refLayerPars, localPos); + } + } + + return ActsExamples::ProcessCode::SUCCESS; +} diff --git a/Examples/Algorithms/TrackFindingExaTrkX/src/PrototracksToParameters.cpp b/Examples/Algorithms/TrackFindingExaTrkX/src/PrototracksToParameters.cpp index 9ffe49ae460..53fe1cdbdad 100644 --- a/Examples/Algorithms/TrackFindingExaTrkX/src/PrototracksToParameters.cpp +++ b/Examples/Algorithms/TrackFindingExaTrkX/src/PrototracksToParameters.cpp @@ -167,17 +167,17 @@ ProcessCode PrototracksToParameters::execute( continue; } - auto pars = Acts::estimateTrackParamsFromSeed( - ctx.geoContext, seed.sp().begin(), seed.sp().end(), surface, field); - - if (not pars) { + auto parsResult = Acts::estimateTrackParamsFromSeed( + ctx.geoContext, seed.sp(), surface, field); + if (!parsResult.ok()) { ACTS_WARNING("Skip track because of bad params"); } + const auto &pars = *parsResult; seededTracks.push_back(track); seeds.emplace_back(std::move(seed)); parameters.push_back(Acts::BoundTrackParameters( - surface.getSharedPtr(), *pars, m_covariance, m_cfg.particleHypothesis)); + surface.getSharedPtr(), pars, m_covariance, m_cfg.particleHypothesis)); } if (skippedTracks > 0) { diff --git a/Examples/Io/Json/include/ActsExamples/Io/Json/JsonTrackParamsLookupReader.hpp b/Examples/Io/Json/include/ActsExamples/Io/Json/JsonTrackParamsLookupReader.hpp new file mode 100644 index 00000000000..94abea8b903 --- /dev/null +++ b/Examples/Io/Json/include/ActsExamples/Io/Json/JsonTrackParamsLookupReader.hpp @@ -0,0 +1,101 @@ +// This file is part of the ACTS project. +// +// Copyright (C) 2016 CERN for the benefit of the ACTS project +// +// This Source Code Form is subject to the terms of the Mozilla Public +// License, v. 2.0. If a copy of the MPL was not distributed with this +// file, You can obtain one at https://mozilla.org/MPL/2.0/. + +#pragma once + +#include "Acts/Plugins/Json/GridJsonConverter.hpp" +#include "Acts/Surfaces/RectangleBounds.hpp" +#include "ActsExamples/TrackFinding/ITrackParamsLookupReader.hpp" + +#include + +#include + +namespace ActsExamples { + +/// @brief Json reader for track parameter lookup tables +/// +/// This reader is used to read track parameter lookup tables +/// from a json file to be later used in track parameter estimation +/// for seeding +class JsonTrackParamsLookupReader final : public ITrackParamsLookupReader { + public: + /// @brief Nested configuration struct + struct Config { + /// Reference tracking layers + std::unordered_map + refLayers; + /// Binning of the grid to be emposed + /// onto the reference layers + std::pair bins; + }; + + explicit JsonTrackParamsLookupReader(const Config& config) : m_cfg(config) {}; + + ~JsonTrackParamsLookupReader() override = default; + + /// @brief Read the lookup from a json file + /// + /// @param path path to the json file + /// + /// @return lookup table for track parameter estimation + TrackParamsLookup readLookup(const std::string& path) const override { + // Read the json file + std::ifstream ifj(path); + nlohmann::json jLookup; + ifj >> jLookup; + + TrackParamsLookup lookup; + // Iterate over the json and deserialize the grids + for (const auto& jGrid : jLookup) { + Acts::GeometryIdentifier id(jGrid["geo_id"]); + + if (!m_cfg.refLayers.contains(id)) { + throw std::invalid_argument("Geometry identifier not found"); + } + + const auto* refSurface = m_cfg.refLayers.at(id); + + // Get bounds to construct the lookup grid + auto bounds = + dynamic_cast(&refSurface->bounds()); + + if (bounds == nullptr) { + throw std::invalid_argument("Only rectangle bounds supported"); + } + + // Axis is not deserilizable, so we need to recreate it + auto halfX = bounds->halfLengthX(); + auto halfY = bounds->halfLengthY(); + + TrackParamsLookupAxisGen axisGen{{-halfX, halfX}, + m_cfg.bins.first, + {-halfY, halfY}, + m_cfg.bins.second}; + + // Deserialize the grid + TrackParamsLookupGrid grid = + Acts::GridJsonConverter::fromJson( + jGrid["grid"], axisGen); + + lookup.try_emplace(id, std::move(grid)); + } + + return lookup; + }; + + /// Readonly access to the config + const Config& config() const { return m_cfg; } + + private: + /// The config of the writer + Config m_cfg; +}; + +} // namespace ActsExamples diff --git a/Examples/Io/Json/include/ActsExamples/Io/Json/JsonTrackParamsLookupWriter.hpp b/Examples/Io/Json/include/ActsExamples/Io/Json/JsonTrackParamsLookupWriter.hpp new file mode 100644 index 00000000000..63a2c083618 --- /dev/null +++ b/Examples/Io/Json/include/ActsExamples/Io/Json/JsonTrackParamsLookupWriter.hpp @@ -0,0 +1,69 @@ +// This file is part of the ACTS project. +// +// Copyright (C) 2016 CERN for the benefit of the ACTS project +// +// This Source Code Form is subject to the terms of the Mozilla Public +// License, v. 2.0. If a copy of the MPL was not distributed with this +// file, You can obtain one at https://mozilla.org/MPL/2.0/. + +#pragma once + +#include "Acts/Plugins/Json/GridJsonConverter.hpp" +#include "ActsExamples/TrackFinding/ITrackParamsLookupWriter.hpp" + +#include + +#include + +namespace ActsExamples { + +/// @brief Json writer for track parameter lookup tables +/// +/// This writer is used to write track parameter lookup tables +/// to a json file to be later used in track parameter estimation +/// for seeding +class JsonTrackParamsLookupWriter final : public ITrackParamsLookupWriter { + public: + /// @brief Nested configuration struct + struct Config { + /// Output file name + std::string path; + }; + + /// Constructor + /// + /// @param config The configuration struct of the writer + explicit JsonTrackParamsLookupWriter(const Config& config) : m_cfg(config) {}; + + /// Virtual destructor + ~JsonTrackParamsLookupWriter() override = default; + + /// Write out track parameters lookup table + /// + /// @param lookup The lookup to write + void writeLookup(const TrackParamsLookup& lookup) const override { + nlohmann::json jLookup; + + // Iterate over the lookup and serialize the grids + for (const auto& [id, grid] : lookup) { + nlohmann::json jGrid; + jGrid["geo_id"] = id.value(); + jGrid["grid"] = Acts::GridJsonConverter::toJson(grid); + + jLookup.push_back(jGrid); + } + + // Write the json file + std::ofstream ofj(m_cfg.path, std::ios::out); + ofj << std::setw(4) << jLookup << std::endl; + }; + + /// Readonly access to the config + const Config& config() const { return m_cfg; } + + private: + /// The config of the writer + Config m_cfg; +}; + +} // namespace ActsExamples diff --git a/Examples/Python/src/Json.cpp b/Examples/Python/src/Json.cpp index be367612ac2..12baaa9f6d3 100644 --- a/Examples/Python/src/Json.cpp +++ b/Examples/Python/src/Json.cpp @@ -19,6 +19,8 @@ #include "ActsExamples/Io/Json/JsonMaterialWriter.hpp" #include "ActsExamples/Io/Json/JsonSurfacesReader.hpp" #include "ActsExamples/Io/Json/JsonSurfacesWriter.hpp" +#include "ActsExamples/Io/Json/JsonTrackParamsLookupReader.hpp" +#include "ActsExamples/Io/Json/JsonTrackParamsLookupWriter.hpp" #include #include @@ -37,6 +39,11 @@ class IMaterialDecorator; namespace ActsExamples { class IMaterialWriter; class IWriter; + +namespace Experimental { +class ITrackParamsLookupWriter; +} // namespace Experimental + } // namespace ActsExamples namespace py = pybind11; @@ -111,6 +118,50 @@ void addJson(Context& ctx) { ACTS_PYTHON_STRUCT_END(); } + { + using IWriter = ActsExamples::ITrackParamsLookupWriter; + using Writer = ActsExamples::JsonTrackParamsLookupWriter; + using Config = Writer::Config; + + auto cls = py::class_>( + mex, "JsonTrackParamsLookupWriter") + .def(py::init(), py::arg("config")) + .def("writeLookup", &Writer::writeLookup) + .def_property_readonly("config", &Writer::config); + + auto c = py::class_(cls, "Config") + .def(py::init<>()) + .def(py::init(), py::arg("path")); + + ACTS_PYTHON_STRUCT_BEGIN(c, Config); + ACTS_PYTHON_MEMBER(path); + ACTS_PYTHON_STRUCT_END(); + } + + { + using IReader = ActsExamples::ITrackParamsLookupReader; + using Reader = ActsExamples::JsonTrackParamsLookupReader; + using Config = Reader::Config; + + auto cls = py::class_>( + mex, "JsonTrackParamsLookupReader") + .def(py::init(), py::arg("config")) + .def("readLookup", &Reader::readLookup) + .def_property_readonly("config", &Reader::config); + + auto c = py::class_(cls, "Config") + .def(py::init<>()) + .def(py::init, + std::pair>(), + py::arg("refLayers"), py::arg("bins")); + + ACTS_PYTHON_STRUCT_BEGIN(c, Config); + ACTS_PYTHON_MEMBER(refLayers); + ACTS_PYTHON_MEMBER(bins); + ACTS_PYTHON_STRUCT_END(); + } + { auto cls = py::class_ #include @@ -279,6 +281,14 @@ void addOutput(Context& ctx) { py::class_>( mex, "IMaterialWriter"); + py::class_>( + mex, "ITrackParamsLookupWriter"); + + py::class_>( + mex, "ITrackParamsLookupReader"); + { using Writer = ActsExamples::RootMaterialWriter; auto w = py::class_>( diff --git a/Examples/Python/src/TrackFinding.cpp b/Examples/Python/src/TrackFinding.cpp index ad3ad364ce7..54eb7e76641 100644 --- a/Examples/Python/src/TrackFinding.cpp +++ b/Examples/Python/src/TrackFinding.cpp @@ -28,6 +28,7 @@ #include "ActsExamples/TrackFinding/SpacePointMaker.hpp" #include "ActsExamples/TrackFinding/TrackFindingAlgorithm.hpp" #include "ActsExamples/TrackFinding/TrackParamsEstimationAlgorithm.hpp" +#include "ActsExamples/TrackFinding/TrackParamsLookupEstimation.hpp" #include #include @@ -293,6 +294,11 @@ void addTrackFinding(Context& ctx) { magneticField, bFieldMin, initialSigmas, initialSigmaPtRel, initialVarInflation, noTimeVarInflation, particleHypothesis); + ACTS_PYTHON_DECLARE_ALGORITHM(ActsExamples::TrackParamsLookupEstimation, mex, + "TrackParamsLookupEstimation", refLayers, bins, + inputHits, inputParticles, + trackLookupGridWriters); + { using Alg = ActsExamples::TrackFindingAlgorithm; using Config = Alg::Config; diff --git a/Examples/Python/tests/requirements.txt b/Examples/Python/tests/requirements.txt index c4067c4adf6..535d3458274 100644 --- a/Examples/Python/tests/requirements.txt +++ b/Examples/Python/tests/requirements.txt @@ -4,42 +4,46 @@ # # pip-compile Examples/Python/tests/requirements.in # -awkward==2.6.1 +awkward==2.7.0 # via # -r Examples/Python/tests/requirements.in # uproot -awkward-cpp==29 +awkward-cpp==41 # via awkward -execnet==2.0.2 +cramjam==2.9.0 + # via uproot +execnet==2.1.1 # via pytest-xdist -fsspec==2024.2.0 +fsspec==2024.10.0 # via # awkward # uproot iniconfig==2.0.0 # via pytest -numpy==1.26.4 +numpy==2.1.3 # via # awkward # awkward-cpp # uproot -packaging==23.2 +packaging==24.2 # via # awkward # pytest # uproot -pluggy==1.4.0 +pluggy==1.5.0 # via pytest -pytest==8.0.0 +pytest==8.3.3 # via # -r Examples/Python/tests/requirements.in # pytest-check # pytest-xdist -pytest-check==2.3.1 +pytest-check==2.4.1 # via -r Examples/Python/tests/requirements.in -pytest-xdist==3.5.0 +pytest-xdist==3.6.1 # via -r Examples/Python/tests/requirements.in -pyyaml==6.0.1 +pyyaml==6.0.2 # via -r Examples/Python/tests/requirements.in -uproot==5.2.2 +uproot==5.5.0 # via -r Examples/Python/tests/requirements.in +xxhash==3.5.0 + # via uproot diff --git a/Examples/Python/tests/root_file_hashes.txt b/Examples/Python/tests/root_file_hashes.txt index 45337b2f762..4674ce1ce2e 100644 --- a/Examples/Python/tests/root_file_hashes.txt +++ b/Examples/Python/tests/root_file_hashes.txt @@ -41,11 +41,11 @@ test_ckf_tracks_example[generic-truth_estimated]__tracksummary_ckf.root: 8e0116c test_ckf_tracks_example[generic-truth_estimated]__performance_seeding.root: 1facb05c066221f6361b61f015cdf0918e94d9f3fce2269ec7b6a4dffeb2bc7e test_ckf_tracks_example[generic-truth_smeared]__trackstates_ckf.root: 82a6744980553e6274df78eea15f0dec22676b1c04e14afc3828bff9bbf5e1b1 test_ckf_tracks_example[generic-truth_smeared]__tracksummary_ckf.root: 06d6ae1d05cb611b19df3c59531997c9b0108f5ef6027d76c4827bd2d9edb921 -test_ckf_tracks_example[odd-full_seeding]__trackstates_ckf.root: 17c48c5a61b1a5495d91336cdf06f9c24e50d81349c1f31d7c70ffff5810a376 -test_ckf_tracks_example[odd-full_seeding]__tracksummary_ckf.root: b5805e54030ab8ac80a8c0a764700c65433dc659783fc8ff3b2c96e512a1d045 +test_ckf_tracks_example[odd-full_seeding]__trackstates_ckf.root: 53a96c225ed78bd7da1491ae3ff4173f0e6c8bffe27e3eab1ab7adfe55426080 +test_ckf_tracks_example[odd-full_seeding]__tracksummary_ckf.root: 2a98d8ec8fae97e18f4661580b71885f558d7222f94bca5edfbf5cdb595021f7 test_ckf_tracks_example[odd-full_seeding]__performance_seeding_trees.root: 43c58577aafe07645e5660c4f43904efadf91d8cda45c5c04c248bbe0f59814f -test_ckf_tracks_example[odd-truth_estimated]__trackstates_ckf.root: 86be5a086d2a87dfde9320bb880bd0788d733ea9727cb5ee6dc0282ec4be39f4 -test_ckf_tracks_example[odd-truth_estimated]__tracksummary_ckf.root: ffce6a73f16986cb3f0386d4a8c1e0ff6f0b4130b9bb12d1af0eb905d000e3e9 +test_ckf_tracks_example[odd-truth_estimated]__trackstates_ckf.root: 7ec5f4b24fa69bfc9d812a6f224d74d36df2d0a42aef791ae0116f309d6f1433 +test_ckf_tracks_example[odd-truth_estimated]__tracksummary_ckf.root: 1eaae038ced2cc5c757480ca42eab60cdaff14d812c34a807a841267d6bfa110 test_ckf_tracks_example[odd-truth_estimated]__performance_seeding.root: 1a36b7017e59f1c08602ef3c2cb0483c51df248f112e3780c66594110719c575 test_ckf_tracks_example[odd-truth_smeared]__trackstates_ckf.root: 35a65e15a6f479f628a96f56ee78e1ac371d71a686ee0c974944d681499fe6bd test_ckf_tracks_example[odd-truth_smeared]__tracksummary_ckf.root: 3e257de624674fa9a19dcc72598c78c29a52633821acaa56dc2aa39a1395f1b5 diff --git a/Examples/Scripts/Python/full_chain_test.py b/Examples/Scripts/Python/full_chain_test.py new file mode 100755 index 00000000000..f0010bd498e --- /dev/null +++ b/Examples/Scripts/Python/full_chain_test.py @@ -0,0 +1,805 @@ +#!/usr/bin/env python3 + +import sys, os, argparse, pathlib + + +def parse_args(): + from acts.examples.reconstruction import SeedingAlgorithm + + parser = argparse.ArgumentParser( + description=""" +Script to test the full chain ACTS simulation and reconstruction. + +This script is provided for interactive developer testing only. +It is not intended (and not supported) for end user use, automated testing, +and certainly should never be called in production. The Python API is the +proper way to access the ActsExamples from scripts. The other Examples/Scripts +are much better examples of how to do that. physmon in the CI is the proper +way to do automated integration tests. This script is only for the case of +interactive testing with one-off configuration specified by command-line options. +""" + ) + parser.add_argument( + "-G", + "--generic-detector", + action="store_true", + help="Use generic detector geometry and config", + ) + parser.add_argument( + "--odd", + default=True, + action=argparse.BooleanOptionalAction, + help="Use Open Data Detector geometry and config (default unless overridden by -G or -A). Requires ACTS_BUILD_ODD.", + ) + parser.add_argument( + "-A", + "--itk", + action="store_true", + help="Use ATLAS ITk geometry and config. Requires acts-itk/ in current directory.", + ) + parser.add_argument( + "-g", + "--geant4", + action="store_true", + help="Use Geant4 instead of Fatras for detector simulation", + ) + parser.add_argument( + "--edm4hep", + type=pathlib.Path, + help="Use edm4hep inputs", + ) + parser.add_argument( + "-b", + "--bf-constant", + action="store_true", + help="Use constant 2T B-field also for ITk; and don't include material map", + ) + parser.add_argument( + "-j", + "--threads", + type=int, + default=-1, + help="Number of parallel threads, negative for automatic (default).", + ) + parser.add_argument( + "-o", + "--output-dir", + "--output", + default=None, + type=pathlib.Path, + help="Directory to write outputs to", + ) + parser.add_argument( + "-O", + "--output-detail", + action="count", + default=0, + help="fewer output files. Use -OO for more output files. Use -OOO to disable all output.", + ) + parser.add_argument( + "-c", + "--output-csv", + action="count", + default=0, + help="Use CSV output instead of ROOT. Specify -cc to output both.", + ) + parser.add_argument( + "-n", + "--events", + type=int, + default=100, + help="The number of events to process (default=%(default)d).", + ) + parser.add_argument( + "-s", + "--skip", + type=int, + default=0, + help="Number of events to skip (default=%(default)d)", + ) + # Many of the following option names were inherited from the old examples binaries and full_chain_odd.py. + # To maintain compatibility, both option names are supported. + parser.add_argument( + "-N", + "--gen-nparticles", + "--gun-particles", + type=int, + default=4, + help="Number of generated particles per vertex from the particle gun (default=%(default)d).", + ) + parser.add_argument( + "-M", + "--gen-nvertices", + "--gun-multiplicity", + "--ttbar-pu", + type=int, + default=200, + help="Number of vertices per event (multiplicity) from the particle gun; or number of pileup events (default=%(default)d)", + ) + parser.add_argument( + "-t", + "--ttbar-pu200", + "--ttbar", + action="store_true", + help="Generate ttbar + mu=200 pile-up using Pythia8", + ) + parser.add_argument( + "-p", + "--gen-pt-range", + "--gun-pt-range", + default="1:10", + help="transverse momentum (pT) range (min:max) of the particle gun in GeV (default=%(default)s)", + ) + parser.add_argument( + "--gen-eta-range", + "--gun-eta-range", + help="Eta range (min:max) of the particle gun (default -2:2 (Generic), -3:3 (ODD), -4:4 (ITk))", + ) + parser.add_argument( + "--gen-cos-theta", + action="store_true", + help="Sample eta as cos(theta) and not uniform", + ) + parser.add_argument( + "-r", + "--random-seed", + type=int, + default=42, + help="Random number seed (default=%(default)d)", + ) + parser.add_argument( + "-F", + "--disable-fpemon", + action="store_true", + help="sets ACTS_SEQUENCER_DISABLE_FPEMON=1", + ) + parser.add_argument( + "-l", + "--loglevel", + type=int, + default=2, + help="The output log level. Please set the wished number (0 = VERBOSE, 1 = DEBUG, 2 = INFO (default), 3 = WARNING, 4 = ERROR, 5 = FATAL).", + ) + parser.add_argument( + "-d", + "--dump-args-calls", + action="store_true", + help="Show pybind function call details", + ) + parser.add_argument( + "--digi-config", + type=pathlib.Path, + help="Digitization configuration file", + ) + parser.add_argument( + "--material-config", + type=pathlib.Path, + help="Material map configuration file", + ) + parser.add_argument( + "-S", + "--seeding-algorithm", + action=EnumAction, + enum=SeedingAlgorithm, + default=SeedingAlgorithm.Default, + help="Select the seeding algorithm to use", + ) + parser.add_argument( + "--ckf", + default=True, + action=argparse.BooleanOptionalAction, + help="Switch CKF on/off", + ) + parser.add_argument( + "--reco", + default=True, + action=argparse.BooleanOptionalAction, + help="Switch reco on/off", + ) + parser.add_argument( + "--vertexing", + default=True, + action=argparse.BooleanOptionalAction, + help="Switch vertexing on/off", + ) + parser.add_argument( + "--MLSeedFilter", + action="store_true", + help="Use the ML seed filter to select seed after the seeding step", + ) + parser.add_argument( + "--ambi-solver", + type=str, + choices=["greedy", "scoring", "ML", "none"], + default="greedy", + help="Set which ambiguity solver to use (default=%(default)s)", + ) + parser.add_argument( + "--ambi-config", + type=pathlib.Path, + default=pathlib.Path.cwd() / "ambi_config.json", + help="Set the configuration file for the Score Based ambiguity resolution (default=%(default)s)", + ) + return parser.parse_args() + + +def full_chain(args): + import acts + + # keep these in memory after we return the sequence + global detector, trackingGeometry, decorators, field, rnd + global logger + + if args.disable_fpemon: + os.environ["ACTS_SEQUENCER_DISABLE_FPEMON"] = "1" + + if args.dump_args_calls: + acts.examples.dump_args_calls(locals()) + + logger = acts.logging.getLogger("full_chain_test") + + nDetArgs = [args.generic_detector, args.odd, args.itk].count(True) + if nDetArgs == 0: + args.generic_detector = True + elif nDetArgs == 2: + args.odd = False + nDetArgs = [args.generic_detector, args.odd, args.itk].count(True) + if nDetArgs != 1: + logger.fatal("require exactly one of: --generic-detector --odd --itk") + sys.exit(2) + if args.generic_detector: + detname = "gen" + elif args.itk: + detname = "itk" + elif args.odd: + detname = "odd" + + u = acts.UnitConstants + + if args.output_detail == 3: + outputDirLess = None + elif args.output_dir is None: + outputDirLess = pathlib.Path.cwd() / f"{detname}_output" + else: + outputDirLess = args.output_dir + + outputDir = None if args.output_detail == 1 else outputDirLess + outputDirMore = None if args.output_detail in (0, 1) else outputDirLess + + outputDirRoot = outputDir if args.output_csv != 1 else None + outputDirLessRoot = outputDirLess if args.output_csv != 1 else None + outputDirMoreRoot = outputDirMore if args.output_csv != 1 else None + outputDirCsv = outputDir if args.output_csv != 0 else None + outputDirLessCsv = outputDirLess if args.output_csv != 0 else None + outputDirMoreCsv = outputDirMore if args.output_csv != 0 else None + + # fmt: off + if args.generic_detector: + etaRange = (-2.0, 2.0) + ptMin = 0.5 * u.GeV + rhoMax = 24.0 * u.mm + geo_dir = pathlib.Path(acts.__file__).resolve().parent.parent.parent.parent.parent + if args.loglevel <= 2: + logger.info(f"Load Generic Detector from {geo_dir}") + if args.digi_config is None: + args.digi_config = geo_dir / "Examples/Algorithms/Digitization/share/default-smearing-config-generic.json" + seedingConfigFile = geo_dir / "Examples/Algorithms/TrackFinding/share/geoSelection-genericDetector.json" + args.bf_constant = True + detector, trackingGeometry, decorators = acts.examples.GenericDetector.create() + elif args.odd: + import acts.examples.odd + etaRange = (-3.0, 3.0) + ptMin = 1.0 * u.GeV + rhoMax = 24.0 * u.mm + beamTime = 1.0 * u.ns + geo_dir = acts.examples.odd.getOpenDataDetectorDirectory() + if args.loglevel <= 2: + logger.info(f"Load Open Data Detector from {geo_dir.resolve()}") + if args.digi_config is None: + args.digi_config = geo_dir / "config/odd-digi-smearing-config.json" + seedingConfigFile = geo_dir / "config/odd-seeding-config.json" + if args.material_config is None: + args.material_config = geo_dir / "data/odd-material-maps.root" + args.bf_constant = True + detector, trackingGeometry, decorators = acts.examples.odd.getOpenDataDetector( + odd_dir=geo_dir, + mdecorator=acts.IMaterialDecorator.fromFile(args.material_config), + ) + elif args.itk: + import acts.examples.itk as itk + etaRange = (-4.0, 4.0) + ptMin = 1.0 * u.GeV + rhoMax = 28.0 * u.mm + beamTime = 5.0 * u.ns + geo_dir = pathlib.Path("acts-itk") + if args.loglevel <= 2: + logger.info(f"Load ATLAS ITk from {geo_dir.resolve()}") + if args.digi_config is None: + args.digi_config = geo_dir / "itk-hgtd/itk-smearing-config.json" + seedingConfigFile = geo_dir / "itk-hgtd/geoSelection-ITk.json" + # args.material_config defaulted in itk.buildITkGeometry: geo_dir / "itk-hgtd/material-maps-ITk-HGTD.json" + bFieldFile = geo_dir / "bfield/ATLAS-BField-xyz.root" + detector, trackingGeometry, decorators = itk.buildITkGeometry( + geo_dir, + customMaterialFile=args.material_config, + material=not args.bf_constant, + logLevel=acts.logging.Level(args.loglevel), + ) + # fmt: on + + if args.bf_constant: + field = acts.ConstantBField(acts.Vector3(0.0, 0.0, 2.0 * u.T)) + else: + logger.info("Create magnetic field map from %s" % str(bFieldFile)) + field = acts.examples.MagneticFieldMapXyz(str(bFieldFile)) + rnd = acts.examples.RandomNumbers(seed=42) + + from acts.examples.simulation import ( + MomentumConfig, + EtaConfig, + PhiConfig, + ParticleConfig, + ParticleSelectorConfig, + addDigitization, + addParticleSelection, + ) + + s = acts.examples.Sequencer( + events=args.events, + skip=args.skip, + numThreads=args.threads if not (args.geant4 and args.threads == -1) else 1, + logLevel=acts.logging.Level(args.loglevel), + outputDir="" if outputDirLess is None else str(outputDirLess), + ) + + # is this needed? + for d in decorators: + s.addContextDecorator(d) + + preSelectParticles = ( + ParticleSelectorConfig( + rho=(0.0 * u.mm, rhoMax), + absZ=(0.0 * u.mm, 1.0 * u.m), + eta=etaRange, + pt=(150 * u.MeV, None), + ) + if args.edm4hep or args.geant4 or args.ttbar_pu200 + else ParticleSelectorConfig() + ) + + postSelectParticles = ParticleSelectorConfig( + pt=(ptMin, None), + eta=etaRange if not args.generic_detector else (None, None), + measurements=(9, None), + removeNeutral=True, + ) + + if args.edm4hep: + import acts.examples.edm4hep + + edm4hepReader = acts.examples.edm4hep.EDM4hepReader( + inputPath=str(args.edm4hep), + inputSimHits=[ + "PixelBarrelReadout", + "PixelEndcapReadout", + "ShortStripBarrelReadout", + "ShortStripEndcapReadout", + "LongStripBarrelReadout", + "LongStripEndcapReadout", + ], + outputParticlesGenerator="particles_input", + outputParticlesInitial="particles_initial", + outputParticlesFinal="particles_final", + outputSimHits="simhits", + graphvizOutput="graphviz", + dd4hepDetector=detector, + trackingGeometry=trackingGeometry, + sortSimHitsInTime=True, + level=acts.logging.INFO, + ) + s.addReader(edm4hepReader) + s.addWhiteboardAlias("particles", edm4hepReader.config.outputParticlesGenerator) + + addParticleSelection( + s, + config=preSelectParticles, + inputParticles="particles", + outputParticles="particles_selected", + ) + + else: + + if not args.ttbar_pu200: + from acts.examples.simulation import addParticleGun + + addParticleGun( + s, + MomentumConfig( + *strToRange(args.gen_pt_range, "--gen-pt-range", u.GeV), + transverse=True, + ), + EtaConfig( + *( + strToRange(args.gen_eta_range, "--gen-eta-range") + if args.gen_eta_range + else etaRange + ), + uniform=( + not args.gen_cos_theta + if args.gen_cos_theta or not args.odd + else None + ), + ), + PhiConfig(0.0, 360.0 * u.degree) if not args.itk else PhiConfig(), + ParticleConfig( + args.gen_nparticles, acts.PdgParticle.eMuon, randomizeCharge=True + ), + vtxGen=( + acts.examples.GaussianVertexGenerator( + mean=acts.Vector4(0, 0, 0, 0), + stddev=acts.Vector4( + 0.0125 * u.mm, 0.0125 * u.mm, 55.5 * u.mm, 1.0 * u.ns + ), + ) + if args.odd + else None + ), + multiplicity=args.gen_nvertices, + rnd=rnd, + outputDirRoot=outputDirMoreRoot, + outputDirCsv=outputDirMoreCsv, + ) + else: + from acts.examples.simulation import addPythia8 + + addPythia8( + s, + hardProcess=["Top:qqbar2ttbar=on"], + npileup=args.gen_nvertices, + vtxGen=acts.examples.GaussianVertexGenerator( + stddev=acts.Vector4( + 0.0125 * u.mm, 0.0125 * u.mm, 55.5 * u.mm, 5.0 * u.ns + ), + mean=acts.Vector4(0, 0, 0, 0), + ), + rnd=rnd, + outputDirRoot=outputDirRoot, + outputDirCsv=outputDirCsv, + ) + + if not args.geant4: + from acts.examples.simulation import addFatras + + addFatras( + s, + trackingGeometry, + field, + rnd=rnd, + preSelectParticles=preSelectParticles, + postSelectParticles=postSelectParticles, + outputDirRoot=outputDirRoot, + outputDirCsv=outputDirCsv, + ) + else: + if s.config.numThreads != 1: + logger.fatal( + f"Geant 4 simulation does not support multi-threading (threads={s.config.numThreads})" + ) + sys.exit(2) + + from acts.examples.simulation import addGeant4 + + # Pythia can sometime simulate particles outside the world volume, a cut on the Z of the track help mitigate this effect + # Older version of G4 might not work, this as has been tested on version `geant4-11-00-patch-03` + # For more detail see issue #1578 + addGeant4( + s, + detector, + trackingGeometry, + field, + rnd=rnd, + preSelectParticles=preSelectParticles, + postSelectParticles=postSelectParticles, + killVolume=trackingGeometry.highestTrackingVolume, + killAfterTime=25 * u.ns, + outputDirRoot=outputDirRoot, + outputDirCsv=outputDirCsv, + ) + + addDigitization( + s, + trackingGeometry, + field, + digiConfigFile=args.digi_config, + rnd=rnd, + outputDirRoot=outputDirRoot, + outputDirCsv=outputDirCsv, + ) + + if not args.reco: + return s + + from acts.examples.reconstruction import ( + addSeeding, + ParticleSmearingSigmas, + addCKFTracks, + CkfConfig, + SeedingAlgorithm, + TrackSelectorConfig, + addAmbiguityResolution, + AmbiguityResolutionConfig, + addVertexFitting, + VertexFinder, + ) + + if args.itk and args.seeding_algorithm == SeedingAlgorithm.Default: + seedingAlgConfig = itk.itkSeedingAlgConfig( + itk.InputSpacePointsType.PixelSpacePoints + ) + else: + seedingAlgConfig = [] + + addSeeding( + s, + trackingGeometry, + field, + *seedingAlgConfig, + seedingAlgorithm=args.seeding_algorithm, + **( + dict( + particleSmearingSigmas=ParticleSmearingSigmas(ptRel=0.01), + rnd=rnd, + ) + if args.seeding_algorithm == SeedingAlgorithm.TruthSmeared + else {} + ), + initialSigmas=[ + 1 * u.mm, + 1 * u.mm, + 1 * u.degree, + 1 * u.degree, + 0.1 * u.e / u.GeV, + 1 * u.ns, + ], + initialSigmaPtRel=0.1, + initialVarInflation=[1.0] * 6, + geoSelectionConfigFile=seedingConfigFile, + outputDirRoot=outputDirLessRoot, + outputDirCsv=outputDirLessCsv, + ) + + if args.MLSeedFilter: + from acts.examples.reconstruction import ( + addSeedFilterML, + SeedFilterMLDBScanConfig, + ) + + addSeedFilterML( + s, + SeedFilterMLDBScanConfig( + epsilonDBScan=0.03, minPointsDBScan=2, minSeedScore=0.1 + ), + onnxModelFile=str( + geo_dir + / "Examples/Scripts/Python/MLAmbiguityResolution/seedDuplicateClassifier.onnx" + ), + outputDirRoot=outputDirLessRoot, + outputDirCsv=outputDirLessCsv, + ) + + if not args.ckf: + return s + + if args.seeding_algorithm != SeedingAlgorithm.TruthSmeared: + ckfConfig = CkfConfig( + seedDeduplication=True, + stayOnSeed=True, + ) + else: + ckfConfig = CkfConfig() + + if not args.itk: + trackSelectorConfig = TrackSelectorConfig( + pt=(ptMin if args.ttbar_pu200 else 0.0, None), + absEta=(None, 3.0), + loc0=(-4.0 * u.mm, 4.0 * u.mm), + nMeasurementsMin=7, + maxHoles=2, + maxOutliers=2, + ) + ckfConfig = ckfConfig._replace( + chi2CutOffMeasurement=15.0, + chi2CutOffOutlier=25.0, + numMeasurementsCutOff=10, + ) + else: + # fmt: off + trackSelectorConfig = ( + TrackSelectorConfig(absEta=(None, 2.0), pt=(0.9 * u.GeV, None), nMeasurementsMin=9, maxHoles=2, maxOutliers=2, maxSharedHits=2), + TrackSelectorConfig(absEta=(None, 2.6), pt=(0.4 * u.GeV, None), nMeasurementsMin=8, maxHoles=2, maxOutliers=2, maxSharedHits=2), + TrackSelectorConfig(absEta=(None, 4.0), pt=(0.4 * u.GeV, None), nMeasurementsMin=7, maxHoles=2, maxOutliers=2, maxSharedHits=2), + ) + # fmt: on + + if args.odd: + ckfConfig = ckfConfig._replace( + pixelVolumes=[16, 17, 18], + stripVolumes=[23, 24, 25], + maxPixelHoles=1, + maxStripHoles=2, + constrainToVolumes=[ + 2, # beam pipe + 32, + 4, # beam pip gap + 16, + 17, + 18, # pixel + 20, # PST + 23, + 24, + 25, # short strip + 26, + 8, # long strip gap + 28, + 29, + 30, # long strip + ], + ) + elif args.itk: + ckfConfig = ckfConfig._replace( + # ITk volumes from Noemi's plot + pixelVolumes=[8, 9, 10, 13, 14, 15, 16, 18, 19, 20], + stripVolumes=[22, 23, 24], + maxPixelHoles=1, + maxStripHoles=2, + ) + + if args.output_detail == 1: + writeDetail = dict(writeTrackSummary=False) + elif args.output_detail == 2: + writeDetail = dict(writeTrackStates=True) + else: + writeDetail = {} + + if args.odd and args.output_detail != 1: + writeCovMat = dict(writeCovMat=True) + else: + writeCovMat = {} + + addCKFTracks( + s, + trackingGeometry, + field, + trackSelectorConfig=trackSelectorConfig, + ckfConfig=ckfConfig, + **writeDetail, + **writeCovMat, + outputDirRoot=outputDirLessRoot, + outputDirCsv=outputDirLessCsv, + ) + + if args.ambi_solver == "ML": + + from acts.examples.reconstruction import ( + addAmbiguityResolutionML, + AmbiguityResolutionMLConfig, + ) + + addAmbiguityResolutionML( + s, + AmbiguityResolutionMLConfig( + maximumSharedHits=3, maximumIterations=1000000, nMeasurementsMin=7 + ), + onnxModelFile=str( + geo_dir + / "Examples/Scripts/Python/MLAmbiguityResolution/duplicateClassifier.onnx" + ), + outputDirRoot=outputDirLessRoot, + outputDirCsv=outputDirLessCsv, + ) + + elif args.ambi_solver == "scoring": + + from acts.examples.reconstruction import ( + addScoreBasedAmbiguityResolution, + ScoreBasedAmbiguityResolutionConfig, + ) + import math + + addScoreBasedAmbiguityResolution( + s, + ScoreBasedAmbiguityResolutionConfig( + minScore=0, + minScoreSharedTracks=1, + maxShared=2, + maxSharedTracksPerMeasurement=2, + pTMax=1400, + pTMin=0.5, + phiMax=math.pi, + phiMin=-math.pi, + etaMax=4, + etaMin=-4, + useAmbiguityFunction=False, + ), + ambiVolumeFile=args.ambi_config, + **writeCovMat, + outputDirRoot=outputDirLessRoot, + outputDirCsv=outputDirLessCsv, + ) + + elif args.ambi_solver == "greedy": + + addAmbiguityResolution( + s, + AmbiguityResolutionConfig( + maximumSharedHits=3, + maximumIterations=10000 if args.itk else 1000000, + nMeasurementsMin=6 if args.itk else 7, + ), + **writeDetail, + **writeCovMat, + outputDirRoot=outputDirLessRoot, + outputDirCsv=outputDirLessCsv, + ) + + if args.vertexing: + addVertexFitting( + s, + field, + vertexFinder=VertexFinder.AMVF, + outputDirRoot=outputDirLessRoot, + ) + + return s + + +def strToRange(s: str, optName: str, unit: float = 1.0): + global logger + try: + range = [float(e) * unit if e != "" else None for e in s.split(":")] + except ValueError: + range = [] + if len(range) == 1: + range.append(range[0]) # 100 -> 100:100 + if len(range) != 2: + logger.fatal(f"bad option value: {optName} {s}") + sys.exit(2) + return range + + +# Graciously taken from https://stackoverflow.com/a/60750535/4280680 (via seeding.py) +class EnumAction(argparse.Action): + """ + Argparse action for handling Enums + """ + + def __init__(self, **kwargs): + import enum + + # Pop off the type value + enum_type = kwargs.pop("enum", None) + + # Ensure an Enum subclass is provided + if enum_type is None: + raise ValueError("type must be assigned an Enum when using EnumAction") + if not issubclass(enum_type, enum.Enum): + raise TypeError("type must be an Enum when using EnumAction") + + # Generate choices from the Enum + kwargs.setdefault("choices", tuple(e.name for e in enum_type)) + + super(EnumAction, self).__init__(**kwargs) + + self._enum = enum_type + + def __call__(self, parser, namespace, values, option_string=None): + for e in self._enum: + if e.name == values: + setattr(namespace, self.dest, e) + break + else: + raise ValueError("%s is not a validly enumerated algorithm." % values) + + +# main program: parse arguments, setup sequence, and run the full chain +full_chain(parse_args()).run() diff --git a/Examples/Scripts/Python/telescope_track_params_lookup_generation.py b/Examples/Scripts/Python/telescope_track_params_lookup_generation.py new file mode 100644 index 00000000000..ecdffc20ec3 --- /dev/null +++ b/Examples/Scripts/Python/telescope_track_params_lookup_generation.py @@ -0,0 +1,111 @@ +#!/usr/bin/env python3 + +import argparse + +import acts +import acts.examples +from acts.examples.simulation import ( + addParticleGun, + addFatras, + MomentumConfig, + EtaConfig, + PhiConfig, + ParticleConfig, +) + +u = acts.UnitConstants + + +def estimateLookup(trackingGeometry, numEvents, outputPath): + + # Set up the dipole magnetic field + field = acts.ConstantBField(acts.Vector3(50 * u.T, 0, 0)) + + # Fatras simulation of muons + rnd = acts.examples.RandomNumbers(seed=42) + + s = acts.examples.Sequencer( + events=numEvents, numThreads=1, logLevel=acts.logging.INFO + ) + + vertexGen = acts.examples.GaussianVertexGenerator( + stddev=acts.Vector4(0, 0, 0, 0), mean=acts.Vector4(0, 9, 0, 0) + ) + + addParticleGun( + s=s, + etaConfig=EtaConfig(10.0, 10.0), + phiConfig=PhiConfig(0, 0), + momentumConfig=MomentumConfig(0.5 * u.GeV, 10 * u.GeV), + particleConfig=ParticleConfig(1, acts.PdgParticle.eMuon, False), + multiplicity=1, + rnd=rnd, + vtxGen=vertexGen, + ) + + addFatras( + s, + trackingGeometry, + field, + inputParticles="particles_input", + outputSimHits="sim_hits", + rnd=rnd, + preSelectParticles=None, + ) + + # Set up the track lookup grid writer + jsonWriterConfig = acts.examples.JsonTrackParamsLookupWriter.Config(path=outputPath) + jsonWriter = acts.examples.JsonTrackParamsLookupWriter(jsonWriterConfig) + + # Set up the track estimation algorithm + surfaces = list(trackingGeometry.geoIdSurfaceMap().values()) + refSurface = surfaces[0] + refGeometryId = refSurface.geometryId() + + trackEstConfig = acts.examples.TrackParamsLookupEstimation.Config( + refLayers={refGeometryId: refSurface}, + bins=(1, 1000), + inputHits="sim_hits", + inputParticles="particles_input", + trackLookupGridWriters=[jsonWriter], + ) + trackEstAlg = acts.examples.TrackParamsLookupEstimation( + trackEstConfig, acts.logging.INFO + ) + + s.addAlgorithm(trackEstAlg) + + s.run() + + +if __name__ == "__main__": + p = argparse.ArgumentParser() + + p.add_argument( + "-n", + "--events", + type=int, + default=100000, + help="Number of events for lookup estimation", + ) + p.add_argument( + "-o", + "--output", + type=str, + default="lookup.json", + help="Output lookup file name", + ) + + args = p.parse_args() + + # Initialize the geometry + detector, trackingGeometry, decorators = acts.examples.TelescopeDetector.create( + bounds=[4, 10], + positions=[30, 60, 90], + stereos=[0, 0, 0], + binValue=2, + surfaceType=0, + ) + + # Estimate the lookup + estimateLookup(trackingGeometry, args.events, args.output) diff --git a/Examples/Scripts/requirements.txt b/Examples/Scripts/requirements.txt index 371122cba63..05f1e49374f 100644 --- a/Examples/Scripts/requirements.txt +++ b/Examples/Scripts/requirements.txt @@ -4,45 +4,51 @@ # # pip-compile Examples/Scripts/requirements.in # -annotated-types==0.6.0 +annotated-types==0.7.0 # via pydantic -awkward==2.5.0 +awkward==2.7.0 # via # -r Examples/Scripts/requirements.in # uproot -awkward-cpp==26 +awkward-cpp==41 # via awkward -boost-histogram==1.4.0 +boost-histogram==1.5.0 # via hist click==8.1.7 # via # histoprint # typer -contourpy==1.2.0 +contourpy==1.3.1 # via matplotlib +cramjam==2.9.0 + # via uproot cycler==0.12.1 # via matplotlib -fonttools==4.46.0 +fonttools==4.55.0 # via matplotlib -hist==2.7.2 +fsspec==2024.10.0 + # via + # awkward + # uproot +hist==2.8.0 # via -r Examples/Scripts/requirements.in -histoprint==2.4.0 +histoprint==2.5.0 # via hist -kiwisolver==1.4.5 +kiwisolver==1.4.7 # via matplotlib markdown-it-py==3.0.0 # via rich -matplotlib==3.8.2 +matplotlib==3.9.2 # via # -r Examples/Scripts/requirements.in # mplhep mdurl==0.1.2 # via markdown-it-py -mplhep==0.3.31 +mplhep==0.3.55 # via -r Examples/Scripts/requirements.in -mplhep-data==0.0.3 +mplhep-data==0.0.4 # via mplhep -numpy==1.26.2 +numpy==2.1.3 # via # awkward # awkward-cpp @@ -56,50 +62,56 @@ numpy==1.26.2 # scipy # uhi # uproot -packaging==23.2 +packaging==24.2 # via # awkward # matplotlib # mplhep # uproot -pandas==2.1.3 +pandas==2.2.3 # via -r Examples/Scripts/requirements.in -pillow==10.1.0 +pillow==11.0.0 # via matplotlib -pydantic==2.5.2 +pydantic==2.9.2 # via -r Examples/Scripts/requirements.in -pydantic-core==2.14.5 +pydantic-core==2.23.4 # via pydantic -pygments==2.17.2 +pygments==2.18.0 # via rich -pyparsing==3.1.1 +pyparsing==3.2.0 # via matplotlib -python-dateutil==2.8.2 +python-dateutil==2.9.0.post0 # via # matplotlib # pandas -pytz==2023.3.post1 +pytz==2024.2 # via pandas -pyyaml==6.0.1 - # via -r Examples/Scripts/requirements.in -rich==13.7.0 +pyyaml==6.0.2 # via -r Examples/Scripts/requirements.in -scipy==1.11.4 +rich==13.9.4 + # via + # -r Examples/Scripts/requirements.in + # typer +scipy==1.14.1 # via -r Examples/Scripts/requirements.in +shellingham==1.5.4 + # via typer six==1.16.0 # via python-dateutil -typer==0.9.0 +typer==0.13.1 # via -r Examples/Scripts/requirements.in -typing-extensions==4.8.0 +typing-extensions==4.12.2 # via # pydantic # pydantic-core # typer -tzdata==2023.3 +tzdata==2024.2 # via pandas -uhi==0.4.0 +uhi==0.5.0 # via # histoprint # mplhep -uproot==5.1.2 +uproot==5.5.0 # via -r Examples/Scripts/requirements.in +xxhash==3.5.0 + # via uproot diff --git a/Plugins/Json/include/Acts/Plugins/Json/TrackParametersJsonConverter.hpp b/Plugins/Json/include/Acts/Plugins/Json/TrackParametersJsonConverter.hpp index ebf7d5c6054..d9670858566 100644 --- a/Plugins/Json/include/Acts/Plugins/Json/TrackParametersJsonConverter.hpp +++ b/Plugins/Json/include/Acts/Plugins/Json/TrackParametersJsonConverter.hpp @@ -8,28 +8,13 @@ #pragma once -#include "Acts/EventData/TrackParameters.hpp" -#include "Acts/Plugins/Json/ActsJson.hpp" +#include "Acts/Definitions/PdgParticle.hpp" +#include "Acts/EventData/GenericBoundTrackParameters.hpp" +#include "Acts/EventData/detail/TrackParametersUtils.hpp" #include "Acts/Plugins/Json/SurfaceJsonConverter.hpp" #include -namespace { - -// Alias to bound adl_serializer specialization -// only to track parameters -template -concept TrackParameters = Acts::FreeTrackParametersConcept || - Acts::BoundTrackParametersConcept; - -// Shorthand for bound track parameters -template -concept IsGenericBound = - std::same_as>; - -} // namespace - namespace Acts { NLOHMANN_JSON_SERIALIZE_ENUM(Acts::PdgParticle, @@ -67,7 +52,7 @@ namespace nlohmann { /// convention is followed. /// /// @tparam parameters_t The track parameters type -template +template struct adl_serializer { /// Covariance matrix type attached to the parameters using CovarianceMatrix = typename parameters_t::CovarianceMatrix; @@ -101,7 +86,7 @@ struct adl_serializer { // Bound track parameters have // reference surface attached // and position takes a geometry context - if constexpr (IsGenericBound) { + if constexpr (Acts::detail::isGenericBoundTrackParams) { Acts::GeometryContext gctx; j["position"] = t.fourPosition(gctx); @@ -152,7 +137,7 @@ struct adl_serializer { // reference surface attached // and constructor is hidden // behind a factory method - if constexpr (IsGenericBound) { + if constexpr (Acts::detail::isGenericBoundTrackParams) { Acts::GeometryContext gctx; auto referenceSurface = Acts::SurfaceJsonConverter::fromJson(j.at("referenceSurface")); @@ -178,7 +163,7 @@ struct adl_serializer { /// convention is followed. /// /// @tparam parameters_t The track parameters type -template +template struct adl_serializer> { using CovarianceMatrix = typename parameters_t::CovarianceMatrix; static void to_json(nlohmann::json& j, @@ -202,7 +187,7 @@ struct adl_serializer> { /// convention is followed. /// /// @tparam parameters_t The track parameters type -template +template struct adl_serializer> { using CovarianceMatrix = typename parameters_t::CovarianceMatrix; static void to_json(nlohmann::json& j, diff --git a/Tests/UnitTests/Core/Seeding/EstimateTrackParamsFromSeedTest.cpp b/Tests/UnitTests/Core/Seeding/EstimateTrackParamsFromSeedTest.cpp index 9aab6c1c70b..45dbbd9d0bb 100644 --- a/Tests/UnitTests/Core/Seeding/EstimateTrackParamsFromSeedTest.cpp +++ b/Tests/UnitTests/Core/Seeding/EstimateTrackParamsFromSeedTest.cpp @@ -174,11 +174,10 @@ BOOST_AUTO_TEST_CASE(trackparameters_estimation_test) { BOOST_CHECK(!estFreeParams.hasNaN()); // Test the bound track parameters estimator - auto fullParamsOpt = estimateTrackParamsFromSeed( - geoCtx, spacePointPtrs.begin(), spacePointPtrs.end(), - *bottomSurface, bField, *logger); - BOOST_REQUIRE(fullParamsOpt.has_value()); - const auto& estFullParams = fullParamsOpt.value(); + auto estFullParamsResult = estimateTrackParamsFromSeed( + geoCtx, spacePointPtrs, *bottomSurface, bField); + BOOST_CHECK(estFullParamsResult.ok()); + const auto& estFullParams = estFullParamsResult.value(); BOOST_TEST_INFO( "The estimated full track parameters at the bottom space point: " "\n" diff --git a/Tests/UnitTests/Core/TrackFinding/CMakeLists.txt b/Tests/UnitTests/Core/TrackFinding/CMakeLists.txt index 523200299ac..bf0c520d672 100644 --- a/Tests/UnitTests/Core/TrackFinding/CMakeLists.txt +++ b/Tests/UnitTests/Core/TrackFinding/CMakeLists.txt @@ -1,2 +1,3 @@ add_unittest(CombinatorialKalmanFilter CombinatorialKalmanFilterTests.cpp) add_unittest(TrackSelector TrackSelectorTests.cpp) +add_unittest(TrackParamsLookupAccumulator TrackParamsLookupAccumulatorTests.cpp) diff --git a/Tests/UnitTests/Core/TrackFinding/TrackParamsLookupAccumulatorTests.cpp b/Tests/UnitTests/Core/TrackFinding/TrackParamsLookupAccumulatorTests.cpp new file mode 100644 index 00000000000..10e12747fce --- /dev/null +++ b/Tests/UnitTests/Core/TrackFinding/TrackParamsLookupAccumulatorTests.cpp @@ -0,0 +1,221 @@ +// This file is part of the ACTS project. +// +// Copyright (C) 2016 CERN for the benefit of the ACTS project +// +// This Source Code Form is subject to the terms of the Mozilla Public +// License, v. 2.0. If a copy of the MPL was not distributed with this +// file, You can obtain one at https://mozilla.org/MPL/2.0/. + +#include +#include +#include + +#include "Acts/Definitions/Algebra.hpp" +#include "Acts/EventData/ParticleHypothesis.hpp" +#include "Acts/EventData/TrackParameters.hpp" +#include "Acts/Geometry/GeometryContext.hpp" +#include "Acts/Surfaces/PlaneSurface.hpp" +#include "Acts/Surfaces/RectangleBounds.hpp" +#include "Acts/Surfaces/Surface.hpp" +#include "Acts/Tests/CommonHelpers/FloatComparisons.hpp" +#include "Acts/TrackFinding/TrackParamsLookupAccumulator.hpp" +#include "Acts/Utilities/AxisFwd.hpp" +#include "Acts/Utilities/Grid.hpp" +#include "Acts/Utilities/GridAxisGenerators.hpp" + +#include +#include +#include +#include +#include +#include + +BOOST_AUTO_TEST_SUITE(TrackParamsLookupAccumulator) + +Acts::GeometryContext gctx; + +using Axis = + Acts::Axis; +using AxisGen = Acts::GridAxisGenerators::EqOpenEqOpen; + +using CellBound = std::pair, + std::shared_ptr>; + +using GridBound = Acts::Grid; +using AccBound = Acts::TrackParamsLookupAccumulator; + +using CellCurvilinear = + std::pair, + std::shared_ptr>; + +using GridCurvilinear = Acts::Grid; +using AccCurvilinear = Acts::TrackParamsLookupAccumulator; + +using CellFree = std::pair, + std::shared_ptr>; + +using GridFree = Acts::Grid; +using AccFree = Acts::TrackParamsLookupAccumulator; + +AxisGen axisGen{{-1, 1}, 2, {-1, 1}, 2}; + +BOOST_AUTO_TEST_CASE(Exceptions) { + // Instantiate grid + GridBound grid(axisGen()); + AccBound acc(grid); + + // Create a reference surface for bound parameters + auto transform = Acts::Transform3::Identity(); + auto bounds1 = std::make_shared(1, 1); + auto bounds2 = std::make_shared(2, 2); + + auto surf1 = + Acts::Surface::makeShared(transform, bounds1); + + auto surf2 = + Acts::Surface::makeShared(transform, bounds2); + + // Create parameters to accumulate + Acts::Vector4 pos{1, 2, 0, 4}; + Acts::Vector3 dir{1, 0, 0}; + Acts::ActsScalar P = 1; + + auto hypothesis1 = Acts::ParticleHypothesis::electron(); + auto hypothesis2 = Acts::ParticleHypothesis::muon(); + + auto pars1 = Acts::BoundTrackParameters::create(surf1, gctx, pos, dir, 1. / P, + std::nullopt, hypothesis1) + .value(); + + auto pars2 = Acts::BoundTrackParameters::create(surf2, gctx, pos, dir, 1. / P, + std::nullopt, hypothesis1) + .value(); + + auto pars3 = Acts::BoundTrackParameters::create(surf1, gctx, pos, dir, 1. / P, + std::nullopt, hypothesis2) + .value(); + + auto pars4 = Acts::BoundTrackParameters::create( + surf1, gctx, pos, dir, -1. / P, std::nullopt, hypothesis2) + .value(); + + // Get the point of the grid + auto bin = grid.localBinsFromGlobalBin(2); + auto center = grid.binCenter(bin); + Acts::Vector2 loc{center.at(0), center.at(1)}; + + // Fill in grid + acc.addTrack(pars1, pars1, loc); + + // Different reference surfaces + BOOST_CHECK_THROW(acc.addTrack(pars2, pars2, loc), std::invalid_argument); + + // Different particle hypotheses + BOOST_CHECK_THROW(acc.addTrack(pars3, pars3, loc), std::invalid_argument); + + // Different charges + BOOST_CHECK_THROW(acc.addTrack(pars4, pars4, loc), std::invalid_argument); +} + +BOOST_AUTO_TEST_CASE(Accumulation) { + // Instantiate grids + GridBound gridBound(axisGen()); + AccBound accBound(gridBound); + + GridCurvilinear gridCurvilinear(axisGen()); + AccCurvilinear accCurvilinear(gridCurvilinear); + + GridFree gridFree(axisGen()); + AccFree accFree(gridFree); + + // Create a reference surface for bound parameters + auto transform = Acts::Transform3::Identity(); + auto bounds = std::make_shared(1, 1); + auto surf = Acts::Surface::makeShared(transform, bounds); + + auto hypothesis = Acts::ParticleHypothesis::electron(); + + std::vector avgPoss; + std::vector avgMoms; + Acts::Vector4 pos{1, 2, 0, 4}; + for (std::size_t i = 0; i < gridBound.size(); i++) { + // Create parameters to accumulate + std::array fourPositions = {pos * (i + 1), pos * (i + 2), + pos * (i + 3), pos * (i + 4)}; + + std::array thetas = { + std::numbers::pi / (i + 1), std::numbers::pi / (i + 2), + std::numbers::pi / (i + 3), std::numbers::pi / (i + 4)}; + + std::array phis = { + 2 * std::numbers::pi / (i + 1), 2 * std::numbers::pi / (i + 2), + 2 * std::numbers::pi / (i + 3), 2 * std::numbers::pi / (i + 4)}; + + Acts::ActsScalar P = 1.5 * (i + 1); + + // Get the point of the grid + auto bin = gridBound.localBinsFromGlobalBin(i); + auto center = gridBound.binCenter(bin); + Acts::Vector2 loc{center.at(0), center.at(1)}; + + // Accumulate + Acts::Vector4 avgPos{0, 0, 0, 0}; + Acts::Vector3 avgMom{0, 0, 0}; + for (std::size_t j = 0; j < 4; j++) { + Acts::Vector3 direction{std::sin(thetas.at(j)) * std::cos(phis.at(j)), + std::sin(thetas.at(j)) * std::sin(phis.at(j)), + std::cos(thetas.at(j))}; + + avgPos += fourPositions.at(j); + avgMom += P * direction; + + // Fill in each grid + auto parsBound = Acts::BoundTrackParameters::create( + surf, gctx, fourPositions.at(j), direction, 1. / P, + std::nullopt, hypothesis) + .value(); + + auto parsCurvilinear = Acts::CurvilinearTrackParameters( + fourPositions.at(j), direction, 1. / P, std::nullopt, hypothesis); + + auto parsFree = Acts::FreeTrackParameters( + fourPositions.at(j), direction, 1. / P, std::nullopt, hypothesis); + + accBound.addTrack(parsBound, parsBound, loc); + accCurvilinear.addTrack(parsCurvilinear, parsCurvilinear, loc); + accFree.addTrack(parsFree, parsFree, loc); + } + avgPoss.push_back(avgPos / fourPositions.size()); + avgMoms.push_back(avgMom / fourPositions.size()); + } + + // Finalize and compare + GridBound avgGridBound = accBound.finalizeLookup(); + GridCurvilinear avgGridCurvilinear = accCurvilinear.finalizeLookup(); + GridFree avgGridFree = accFree.finalizeLookup(); + for (std::size_t i = 0; i < avgGridBound.size(); i++) { + auto [ipBound, refBound] = avgGridBound.at(i); + auto [ipCurvilinear, refCurvilinear] = avgGridCurvilinear.at(i); + auto [ipFree, refFree] = avgGridFree.at(i); + + Acts::Vector4 avgPos = avgPoss.at(i); + + Acts::Vector3 avgMom = avgMoms.at(i); + Acts::Vector3 avgDir = avgMom.normalized(); + Acts::ActsScalar avgP = avgMom.norm(); + + CHECK_CLOSE_ABS(ipBound->fourPosition(gctx), avgPos, 1e-3); + CHECK_CLOSE_ABS(ipBound->direction(), avgDir, 1e-3); + CHECK_CLOSE_ABS(ipBound->absoluteMomentum(), avgP, 1e-3); + + CHECK_CLOSE_ABS(ipCurvilinear->fourPosition(), avgPos, 1e-3); + CHECK_CLOSE_ABS(ipCurvilinear->direction(), avgDir, 1e-3); + CHECK_CLOSE_ABS(ipCurvilinear->absoluteMomentum(), avgP, 1e-3); + + CHECK_CLOSE_ABS(ipFree->fourPosition(), avgPos, 1e-3); + CHECK_CLOSE_ABS(ipFree->direction(), avgDir, 1e-3); + CHECK_CLOSE_ABS(ipFree->absoluteMomentum(), avgP, 1e-3); + } +} + +BOOST_AUTO_TEST_SUITE_END() diff --git a/docs/requirements.in b/docs/requirements.in index cd5b9dff585..1633d176167 100644 --- a/docs/requirements.in +++ b/docs/requirements.in @@ -1,6 +1,6 @@ breathe myst-parser -sphinx +sphinx==7.3.7 sphinx_rtd_theme docutils rich diff --git a/docs/requirements.txt b/docs/requirements.txt index dd60e5e5724..0ac5d41d6b6 100644 --- a/docs/requirements.txt +++ b/docs/requirements.txt @@ -4,52 +4,54 @@ # # pip-compile docs/requirements.in # -aiohttp==3.9.1 +aiohappyeyeballs==2.4.3 + # via aiohttp +aiohttp==3.11.5 # via -r docs/requirements.in aiosignal==1.3.1 # via aiohttp -alabaster==0.7.13 +alabaster==0.7.16 # via sphinx -annotated-types==0.6.0 +annotated-types==0.7.0 # via pydantic -attrs==23.1.0 +attrs==24.2.0 # via aiohttp -babel==2.13.1 +babel==2.16.0 # via sphinx breathe==4.35.0 # via -r docs/requirements.in -certifi==2023.11.17 +certifi==2024.8.30 # via requests -cffi==1.16.0 +cffi==1.17.1 # via cryptography -charset-normalizer==3.3.2 +charset-normalizer==3.4.0 # via requests click==8.1.7 # via typer -cryptography==41.0.7 +cryptography==43.0.3 # via pyjwt -docutils==0.20.1 +docutils==0.21.2 # via # -r docs/requirements.in # breathe # myst-parser # sphinx # sphinx-rtd-theme -frozenlist==1.4.0 +frozenlist==1.5.0 # via # aiohttp # aiosignal -fsspec==2023.12.0 +fsspec==2024.10.0 # via -r docs/requirements.in gidgethub==5.3.0 # via -r docs/requirements.in -idna==3.6 +idna==3.10 # via # requests # yarl imagesize==1.4.1 # via sphinx -jinja2==3.1.2 +jinja2==3.1.4 # via # -r docs/requirements.in # myst-parser @@ -59,85 +61,85 @@ markdown-it-py==3.0.0 # mdit-py-plugins # myst-parser # rich -markupsafe==2.1.3 +markupsafe==3.0.2 # via jinja2 -mdit-py-plugins==0.4.0 +mdit-py-plugins==0.4.2 # via myst-parser mdurl==0.1.2 # via markdown-it-py -multidict==6.0.4 +multidict==6.1.0 # via # aiohttp # yarl -myst-parser==2.0.0 +myst-parser==4.0.0 # via -r docs/requirements.in -packaging==23.2 +packaging==24.2 # via sphinx -pycparser==2.21 +propcache==0.2.0 + # via + # aiohttp + # yarl +pycparser==2.22 # via cffi -pydantic==2.5.2 +pydantic==2.9.2 # via -r docs/requirements.in -pydantic-core==2.14.5 +pydantic-core==2.23.4 # via pydantic -pygments==2.17.2 +pygments==2.18.0 # via # rich # sphinx -pyjwt[crypto]==2.8.0 +pyjwt[crypto]==2.10.0 # via gidgethub -python-dotenv==1.0.0 +python-dotenv==1.0.1 # via -r docs/requirements.in -pyyaml==6.0.1 +pyyaml==6.0.2 # via myst-parser -requests==2.31.0 +requests==2.32.3 # via sphinx -rich==13.7.0 - # via -r docs/requirements.in +rich==13.9.4 + # via + # -r docs/requirements.in + # typer +shellingham==1.5.4 + # via typer snowballstemmer==2.2.0 # via sphinx -sphinx==7.2.6 +sphinx==7.3.7 # via # -r docs/requirements.in # breathe # myst-parser # sphinx-rtd-theme - # sphinxcontrib-applehelp - # sphinxcontrib-devhelp - # sphinxcontrib-htmlhelp # sphinxcontrib-jquery - # sphinxcontrib-qthelp - # sphinxcontrib-serializinghtml -sphinx-rtd-theme==2.0.0 +sphinx-rtd-theme==3.0.2 # via -r docs/requirements.in -sphinxcontrib-applehelp==1.0.7 +sphinxcontrib-applehelp==2.0.0 # via sphinx -sphinxcontrib-devhelp==1.0.5 +sphinxcontrib-devhelp==2.0.0 # via sphinx -sphinxcontrib-htmlhelp==2.0.4 +sphinxcontrib-htmlhelp==2.1.0 # via sphinx sphinxcontrib-jquery==4.1 # via sphinx-rtd-theme sphinxcontrib-jsmath==1.0.1 # via sphinx -sphinxcontrib-qthelp==1.0.6 +sphinxcontrib-qthelp==2.0.0 # via sphinx -sphinxcontrib-serializinghtml==1.1.9 +sphinxcontrib-serializinghtml==2.0.0 # via sphinx toml==0.10.2 # via -r docs/requirements.in -typer==0.9.0 +typer==0.13.1 # via -r docs/requirements.in -typing-extensions==4.8.0 +typing-extensions==4.12.2 # via # pydantic # pydantic-core # typer uritemplate==4.1.1 # via gidgethub -urllib3==2.1.0 +urllib3==2.2.3 # via requests -yarl==1.9.3 +yarl==1.17.2 # via aiohttp - -# The following packages are considered to be unsafe in a requirements file: -# setuptools