ActsExamples/EventData/AverageSimHits.hpp

0001 // This file is part of the ACTS project.
0002 //
0003 // Copyright (C) 2016 CERN for the benefit of the ACTS project
0004 //
0005 // This Source Code Form is subject to the terms of the Mozilla Public
0006 // License, v. 2.0. If a copy of the MPL was not distributed with this
0007 // file, You can obtain one at https://mozilla.org/MPL/2.0/.
0008
0009 #pragma once
0010
0011 #include "Acts/Definitions/Algebra.hpp"
0012 #include "Acts/Surfaces/Surface.hpp"
0013 #include "Acts/Utilities/Logger.hpp"
0014 #include "ActsExamples/EventData/Index.hpp"
0015 #include "ActsExamples/EventData/SimHit.hpp"
0016 #include "ActsExamples/Utilities/Range.hpp"
0017
0018 #include <tuple>
0019
0020 namespace ActsExamples {
0021
0022 /// A range within a hit-simhits map.
0023 using HitSimHitsRange = Range<IndexMultimap<Index>::const_iterator>;
0024
0025 /// Create (average) truth representation for selected simulated hits.
0026 ///
0027 /// @param gCtx The geometry context for this
0028 /// @param surface The reference surface of the measurement
0029 /// @param simHits The simulated hits container
0030 /// @param hitSimHitsRange Selection of simulated hits from the container
0031 /// @return a local position, a 4D global position, a direction
0032 ///
0033 /// If more than one simulated hit is selected, the average truth information is
0034 /// returned.
0035 inline std::tuple<Acts::Vector2, Acts::Vector4, Acts::Vector3> averageSimHits(
0036     const Acts::GeometryContext& gCtx, const Acts::Surface& surface,
0037     const SimHitContainer& simHits, const HitSimHitsRange& hitSimHitsRange,
0038     const Acts::Logger& logger) {
0039   Acts::Vector2 avgLocal = Acts::Vector2::Zero();
0040   Acts::Vector4 avgPos4 = Acts::Vector4::Zero();
0041   Acts::Vector3 avgDir = Acts::Vector3::Zero();
0042
0043   std::size_t n = 0u;
0044   for (auto [_, simHitIdx] : hitSimHitsRange) {
0045     n += 1u;
0046
0047     // we assume that the indices are within valid ranges so we do not need to
0048     // check their validity again.
0049     const auto& simHit = *simHits.nth(simHitIdx);
0050
0051     // We use the thickness of the surface as tolerance, because Geant4
0052     // treats the Surfaces as volumes and thus it is not ensured, that each hit
0053     // lies exactly on the Acts::Surface
0054     const double tolerance = surface.isSensitive() ? surface.thickness()
0055                                                    : Acts::s_onSurfaceTolerance;
0056
0057     // transforming first to local positions and average that ensures that the
0058     // averaged position is still on the surface. the averaged global position
0059     // might not be on the surface anymore.
0060     auto result = surface.globalToLocal(gCtx, simHit.position(),
0061                                         simHit.direction(), tolerance);
0062     if (result.ok()) {
0063       avgLocal += result.value();
0064     } else {
0065       ACTS_WARNING("While averaging simhit, hit "
0066                    << simHitIdx << " is not on the corresponding surface "
0067                    << surface.geometryId() << "; use [0,0] as local position");
0068     }
0069     // global position should already be at the intersection. no need to perform
0070     // an additional intersection call.
0071     avgPos4 += simHit.fourPosition();
0072     avgDir += simHit.direction();
0073   }
0074
0075   // only need to average if there are at least two inputs
0076   if (2u <= n) {
0077     double scale = 1.0 / n;
0078     avgLocal *= scale;
0079     avgPos4 *= scale;
0080     avgDir.normalize();
0081   }
0082
0083   return {avgLocal, avgPos4, avgDir};
0084 }
0085
0086 }  // namespace ActsExamples