|
|
|||
Warning, file /acts/Core/include/Acts/EventData/TrackProxy.hpp was not indexed or was modified since last indexation (in which case cross-reference links may be missing, inaccurate or erroneous).
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/TrackParametrization.hpp" 0012 #include "Acts/EventData/MultiTrajectory.hpp" 0013 #include "Acts/EventData/MultiTrajectoryBackendConcept.hpp" 0014 #include "Acts/EventData/ParticleHypothesis.hpp" 0015 #include "Acts/EventData/TrackContainerBackendConcept.hpp" 0016 #include "Acts/EventData/TrackParameters.hpp" 0017 #include "Acts/EventData/TrackProxyConcept.hpp" 0018 #include "Acts/EventData/TrackStatePropMask.hpp" 0019 #include "Acts/Utilities/HashedString.hpp" 0020 #include "Acts/Utilities/TypeTraits.hpp" 0021 #include "Acts/Utilities/UnitVectors.hpp" 0022 0023 #include <iterator> 0024 0025 namespace Acts { 0026 0027 template <TrackContainerBackend track_container_t, 0028 CommonMultiTrajectoryBackend traj_t, 0029 template <typename> class holder_t> 0030 class TrackContainer; 0031 0032 /// Proxy class representing a single track. 0033 /// This class provides a **view** into an associated @ref TrackContainer, and 0034 /// has **reference semantics**. You can think of it as a pointer to a vector 0035 /// of tracks, which exposes an object-oriented interface for accessing the 0036 /// track properties. 0037 /// 0038 /// @tparam track_container_t the container backend 0039 /// @tparam trajectory_t the track state container backend 0040 /// @tparam holder_t ownership management class for the backend 0041 /// @tparam read_only true if this track container is not mutable 0042 template <typename track_container_t, typename trajectory_t, 0043 template <typename> class holder_t, bool read_only = true> 0044 class TrackProxy { 0045 public: 0046 /// Indicates whether this track proxy is read-only or if it can be modified 0047 static constexpr bool ReadOnly = read_only; 0048 0049 /// The track container backend given as a template parameter 0050 using Container = track_container_t; 0051 0052 /// The track state container backend given as a template parameter 0053 using Trajectory = trajectory_t; 0054 0055 /// The index type of the track container 0056 using IndexType = typename Container::IndexType; 0057 0058 /// Sentinel value that indicates an invalid index 0059 static constexpr IndexType kInvalid = Container::kInvalid; 0060 0061 /// Alias for the mutable version of this track proxy, with the same backends 0062 using MutableTrackProxy = 0063 TrackProxy<track_container_t, trajectory_t, holder_t, false>; 0064 0065 /// Alias for the const version of this track proxy, with the same backends 0066 using ConstTrackProxy = 0067 TrackProxy<track_container_t, trajectory_t, holder_t, true>; 0068 0069 /// Alias for an associated const track proxy, with the same backends 0070 using ConstProxyType = ConstTrackProxy; 0071 0072 /// Alias for an associated mutable track state proxy, with the same backends 0073 using TrackStateProxy = typename Trajectory::TrackStateProxy; 0074 0075 /// Alias for an associated const track state proxy, with the same backends 0076 using ConstTrackStateProxy = typename Trajectory::ConstTrackStateProxy; 0077 0078 /// Map-type for a bound parameter vector. This has reference semantics, i.e. 0079 /// points at a matrix by an internal pointer. 0080 using Parameters = 0081 typename detail_lt::FixedSizeTypes<eBoundSize, false>::CoefficientsMap; 0082 0083 /// Same as @ref Parameters, but with const semantics 0084 using ConstParameters = 0085 typename detail_lt::FixedSizeTypes<eBoundSize, true>::CoefficientsMap; 0086 0087 /// Map-type for a bound covariance. This has reference semantics, i.e. 0088 /// points at a matrix by an internal pointer. 0089 using Covariance = 0090 typename detail_lt::FixedSizeTypes<eBoundSize, false>::CovarianceMap; 0091 0092 /// Same as @ref Covariance, but with const semantics 0093 using ConstCovariance = 0094 typename detail_lt::FixedSizeTypes<eBoundSize, true>::CovarianceMap; 0095 0096 #ifndef DOXYGEN 0097 friend TrackContainer<Container, Trajectory, holder_t>; 0098 friend MutableTrackProxy; 0099 friend ConstTrackProxy; 0100 // Track proxies are friends, not food! 0101 template <typename A, typename B, template <typename> class H, bool R> 0102 friend class TrackProxy; 0103 #endif 0104 0105 /// @anchor track_proxy_construct 0106 /// @name Constructors and assignment operator 0107 /// 0108 /// Public constructors and assignment operators for @c TrackProxy only 0109 /// allow construction from another @c TrackProxy. You should generally 0110 /// not have to construct @c TrackProxy manually. 0111 /// 0112 /// @{ 0113 0114 /// Copy constructor: const to const or mutable to mutable 0115 /// @param other the other track proxy 0116 TrackProxy(const TrackProxy& other) = default; 0117 0118 /// Copy assignment operator: const to const or mutable to mutable 0119 /// @param other the other track proxy 0120 /// @return reference to this track proxy 0121 TrackProxy& operator=(const TrackProxy& other) = default; 0122 0123 /// Constructor from mutable track proxy 0124 /// @note Only available if the track proxy is read-only 0125 /// @param other the other track proxy 0126 explicit TrackProxy(const MutableTrackProxy& other) 0127 requires ReadOnly 0128 : m_container{other.m_container}, m_index{other.m_index} {} 0129 0130 /// Copy assignment operator from mutable track proxy 0131 /// @note Only available if the track proxy is read-only 0132 /// @param other the other track proxy 0133 /// @return reference to this track proxy 0134 TrackProxy& operator=(const MutableTrackProxy& other) 0135 requires ReadOnly 0136 { 0137 m_container = other.m_container; 0138 m_index = other.m_index; 0139 return *this; 0140 } 0141 0142 /// @} 0143 0144 /// Equality operator with another track proxy 0145 /// Checks the container identity and the track index 0146 /// @param other Other track proxy to compare with 0147 /// @return True if the track proxies refer to the same track 0148 bool operator==(const TrackProxy& other) const { 0149 return &(*m_container) == &(*other.m_container) && m_index == other.m_index; 0150 } 0151 0152 /// @anchor track_proxy_props 0153 /// @name TrackProxy properties 0154 /// Methods that give access to the properties of a track represented by 0155 /// @c TrackProxy. 0156 /// 0157 /// Many of these methods come in a @c const and a non-@c const version. The 0158 /// non-@c const version is only available if you have an instance of 0159 /// @c TrackProxy that does not have the @c read_only template parameter set to 0160 /// @c true, even if you hold it as an lvalue. 0161 /// 0162 /// @{ 0163 0164 /// Get the tip index, i.e. the entry point into the track state container 0165 /// @return the tip index by value 0166 IndexType tipIndex() const { 0167 return component<IndexType>(hashString("tipIndex")); 0168 } 0169 0170 /// Index of the stem, i.e. the innermost track state of the track. 0171 /// This might be invalid, signifying that the track state is not 0172 /// forward-linked. 0173 /// @return the stem index 0174 IndexType stemIndex() const { 0175 return component<IndexType>(hashString("stemIndex")); 0176 } 0177 0178 /// Get a mutable reference to the tip index, i.e. the entry point into the 0179 /// track container 0180 /// @note Only available if the track proxy is not read-only 0181 /// @return mutable reference to the tip index 0182 IndexType& tipIndex() 0183 requires(!ReadOnly) 0184 { 0185 return component<IndexType>(hashString("tipIndex")); 0186 } 0187 0188 /// Index of the stem, i.e. the innermost track state of the track. 0189 /// This might be invalid, signifying that the track state is not 0190 /// forward-linked. 0191 /// @note Only available if the track proxy is not read-only 0192 /// @return mutable reference to the stem index 0193 IndexType& stemIndex() 0194 requires(!ReadOnly) 0195 { 0196 return component<IndexType>(hashString("stemIndex")); 0197 } 0198 0199 /// Get the reference surface of the track (e.g. the perigee) 0200 /// @return the reference surface 0201 const Surface& referenceSurface() const { 0202 return *m_container->container().referenceSurface_impl(m_index); 0203 } 0204 0205 // NOLINTBEGIN(performance-unnecessary-value-param) 0206 // looks like a false-positive. clang-tidy believes `srf` is not movable. 0207 /// Set a new reference surface for this track 0208 /// @param srf The surface to set 0209 void setReferenceSurface(std::shared_ptr<const Surface> srf) 0210 requires(!ReadOnly) 0211 { 0212 m_container->container().setReferenceSurface_impl(m_index, std::move(srf)); 0213 } 0214 // NOLINTEND(performance-unnecessary-value-param) 0215 0216 /// Returns whether the track has a reference surface or not 0217 /// @return whether a surface exists or not 0218 bool hasReferenceSurface() const { 0219 // @TODO: This could be more efficient 0220 return m_container->container().referenceSurface_impl(m_index) != nullptr; 0221 } 0222 0223 /// Get the parameters of the track at the reference surface (e.g. perigee). 0224 /// Const version 0225 /// @return Proxy vector for the parameters 0226 ConstParameters parameters() const { 0227 return m_container->parameters(m_index); 0228 } 0229 0230 /// Get the covariance of the track at the reference surface (e.g. perigee). 0231 /// Const version 0232 /// @return Proxy matrix for the covariance 0233 ConstCovariance covariance() const { 0234 return m_container->covariance(m_index); 0235 } 0236 0237 /// Get the parameters of the track at the reference surface (e.g. perigee). 0238 /// Mutable version 0239 /// @note Only available if the track proxy is not read-only 0240 /// @return Proxy vector for the parameters 0241 Parameters parameters() 0242 requires(!ReadOnly) 0243 { 0244 return m_container->parameters(m_index); 0245 } 0246 0247 /// Get the covariance of the track at the reference surface (e.g. perigee). 0248 /// Mutable version 0249 /// @note Only available if the track proxy is not read-only 0250 /// @return Proxy matrix for the covariance 0251 Covariance covariance() 0252 requires(!ReadOnly) 0253 { 0254 return m_container->covariance(m_index); 0255 } 0256 0257 /// Access the theta parameter of the track at the reference surface 0258 /// @return The theta parameter 0259 double theta() const { return parameters()[eBoundTheta]; } 0260 0261 /// Access the phi parameter of the track at the reference surface 0262 /// @return The phi parameter 0263 double phi() const { return parameters()[eBoundPhi]; } 0264 0265 /// Access the loc0 parameter of the track at the reference surface 0266 /// @return The loc0 parameter 0267 double loc0() const { return parameters()[eBoundLoc0]; } 0268 0269 /// Access the loc1 parameter of the track at the reference surface 0270 /// @return The loc1 parameter 0271 double loc1() const { return parameters()[eBoundLoc1]; } 0272 0273 /// Access the time parameter of the track at the reference surface 0274 /// @return The time parameter 0275 double time() const { return parameters()[eBoundTime]; } 0276 0277 /// Access the q/p (curvature) parameter of the track at the reference surface 0278 /// @return The q/p parameter 0279 double qOverP() const { return parameters()[eBoundQOverP]; } 0280 0281 /// Get the particle hypothesis 0282 /// @return the particle hypothesis 0283 ParticleHypothesis particleHypothesis() const { 0284 return m_container->container().particleHypothesis_impl(m_index); 0285 } 0286 0287 /// Set a new particle hypothesis for this track 0288 /// @note Only available if the track proxy is not read-only 0289 /// @param particleHypothesis The particle hypothesis to set 0290 void setParticleHypothesis(const ParticleHypothesis& particleHypothesis) 0291 requires(!ReadOnly) 0292 { 0293 m_container->container().setParticleHypothesis_impl(m_index, 0294 particleHypothesis); 0295 } 0296 0297 /// Get the charge of the tack 0298 /// @note this depends on the charge hypothesis 0299 /// @return The absolute track momentum 0300 double charge() const { return particleHypothesis().qFromQOP(qOverP()); } 0301 0302 /// Get the absolute momentum of the tack 0303 /// @return The absolute track momentum 0304 double absoluteMomentum() const { 0305 return particleHypothesis().extractMomentum(qOverP()); 0306 } 0307 0308 /// Get the transverse momentum of the track 0309 /// @return The track transverse momentum value 0310 double transverseMomentum() const { 0311 return std::sin(theta()) * absoluteMomentum(); 0312 } 0313 0314 /// Get a unit vector along the track direction at the reference surface 0315 /// @return The direction unit vector 0316 Vector3 direction() const { 0317 return makeDirectionFromPhiTheta(phi(), theta()); 0318 } 0319 0320 /// Get the global momentum vector 0321 /// @return the global momentum vector 0322 Vector3 momentum() const { return absoluteMomentum() * direction(); } 0323 0324 /// Get the four-momentum vector: (px, py, pz, e) 0325 /// @return the four-momentum vector 0326 Vector4 fourMomentum() const { 0327 Vector4 p4 = Vector4::Zero(); 0328 0329 Vector3 p3 = momentum(); 0330 p4[eMom0] = p3[eMom0]; 0331 p4[eMom1] = p3[eMom1]; 0332 p4[eMom2] = p3[eMom2]; 0333 0334 float m = particleHypothesis().mass(); 0335 p4[eEnergy] = std::sqrt(m * m + p3.squaredNorm()); 0336 0337 return p4; 0338 } 0339 0340 /// Return the number of track states associated to this track 0341 /// @note This is calculated by iterating over the track states which is 0342 /// somewhat expensive. Consider caching this value if you need It 0343 /// more than once. 0344 /// @return The number of track states 0345 unsigned int nTrackStates() const { 0346 // @TODO: This should probably be cached, distance is expensive 0347 // without random access 0348 if (tipIndex() == kInvalid) { 0349 // no tip index -> no track states 0350 return 0; 0351 } 0352 auto tsRange = trackStatesReversed(); 0353 return std::distance(tsRange.begin(), tsRange.end()); 0354 } 0355 0356 /// Return a mutable reference to the number of measurements for the track. 0357 /// Mutable version 0358 /// @note Only available if the track proxy is not read-only 0359 /// @return The number of measurements 0360 unsigned int& nMeasurements() 0361 requires(!ReadOnly) 0362 { 0363 return component<unsigned int, hashString("nMeasurements")>(); 0364 } 0365 0366 /// Return the number of measurements for the track. Const version 0367 /// @return The number of measurements 0368 unsigned int nMeasurements() const { 0369 return component<unsigned int, hashString("nMeasurements")>(); 0370 } 0371 0372 /// Return a mutable reference to the number of holes for the track. 0373 /// Mutable version 0374 /// @note Only available if the track proxy is not read-only 0375 /// @return The number of holes 0376 unsigned int& nHoles() 0377 requires(!ReadOnly) 0378 { 0379 return component<unsigned int, hashString("nHoles")>(); 0380 } 0381 0382 /// Return the number of measurements for the track. Const version 0383 /// @return The number of measurements 0384 unsigned int nHoles() const { 0385 return component<unsigned int, hashString("nHoles")>(); 0386 } 0387 0388 /// Return a mutable reference to the number of outliers for the track. 0389 /// Mutable version 0390 /// @note Only available if the track proxy is not read-only 0391 /// @return The number of outliers 0392 unsigned int& nOutliers() 0393 requires(!ReadOnly) 0394 { 0395 return component<unsigned int, hashString("nOutliers")>(); 0396 } 0397 0398 /// Return the number of outliers for the track. Const version 0399 /// @return The number of outliers 0400 unsigned int nOutliers() const { 0401 return component<unsigned int, hashString("nOutliers")>(); 0402 } 0403 0404 /// Return a mutable reference to the number of shared hits for the track. 0405 /// Mutable version 0406 /// @note Only available if the track proxy is not read-only 0407 /// @return The number of shared hits 0408 unsigned int& nSharedHits() 0409 requires(!ReadOnly) 0410 { 0411 return component<unsigned int, hashString("nSharedHits")>(); 0412 } 0413 0414 /// Return the number of shared hits for the track. Const version 0415 /// @return The number of shared hits 0416 unsigned int nSharedHits() const { 0417 return component<unsigned int, hashString("nSharedHits")>(); 0418 } 0419 0420 /// Return a mutable reference to the chi squared 0421 /// Mutable version 0422 /// @note Only available if the track proxy is not read-only 0423 /// @return The chi squared 0424 float& chi2() 0425 requires(!ReadOnly) 0426 { 0427 return component<float, hashString("chi2")>(); 0428 } 0429 0430 /// Return the chi squared for the track. Const version 0431 /// @return The chi squared 0432 float chi2() const { return component<float, hashString("chi2")>(); } 0433 0434 /// Return a mutable reference to the number of degrees of freedom for the 0435 /// track. Mutable version 0436 /// @note Only available if the track proxy is not read-only 0437 /// @return The number of degrees of freedom 0438 unsigned int& nDoF() 0439 requires(!ReadOnly) 0440 { 0441 return component<unsigned int, hashString("ndf")>(); 0442 } 0443 0444 /// Return the number of degrees of freedom for the track. Const version 0445 /// @return The number of degrees of freedom 0446 unsigned int nDoF() const { 0447 return component<unsigned int, hashString("ndf")>(); 0448 } 0449 0450 /// Return the index of this track in the track container 0451 /// @note This is separate from the tip index 0452 /// @return the track index 0453 IndexType index() const { return m_index; } 0454 0455 /// @} 0456 0457 /// @anchor track_proxy_track_states 0458 /// @name TrackProxy track state access 0459 /// Methods that give access to the track states of a track represented by @c TrackProxy. 0460 /// @{ 0461 0462 /// Return a const track state proxy to the outermost track state 0463 /// @return The outermost track state proxy 0464 ConstTrackStateProxy outermostTrackState() const { 0465 return m_container->trackStateContainer().getTrackState(tipIndex()); 0466 } 0467 0468 /// Return a mutable track state proxy to the outermost track state 0469 /// @return The outermost track state proxy 0470 TrackStateProxy outermostTrackState() 0471 requires(!ReadOnly) 0472 { 0473 return m_container->trackStateContainer().getTrackState(tipIndex()); 0474 } 0475 0476 /// Return a const track state proxy to the innermost track state 0477 /// @note This is only available, if the track is forward linked 0478 /// @return The innermost track state proxy 0479 auto innermostTrackState() const { 0480 using proxy_t = decltype(m_container->trackStateContainer().getTrackState( 0481 std::declval<IndexType>())); 0482 0483 IndexType stem = component<IndexType, hashString("stemIndex")>(); 0484 if (stem == kInvalid) { 0485 return std::optional<proxy_t>{}; 0486 } else { 0487 return std::optional<proxy_t>{ 0488 m_container->trackStateContainer().getTrackState(stem)}; 0489 } 0490 } 0491 0492 /// Return a mutable track state proxy to the innermost track state 0493 /// @note This is only available, if the track is forward linked 0494 /// @note Only available if the track proxy is not read-only 0495 /// @return The innermost track state proxy 0496 auto innermostTrackState() 0497 requires(!ReadOnly) 0498 { 0499 using proxy_t = decltype(m_container->trackStateContainer().getTrackState( 0500 std::declval<IndexType>())); 0501 0502 IndexType stem = component<IndexType>(hashString("stemIndex")); 0503 if (stem == kInvalid) { 0504 return std::optional<proxy_t>{}; 0505 } else { 0506 return std::optional<proxy_t>{ 0507 m_container->trackStateContainer().getTrackState(stem)}; 0508 } 0509 } 0510 0511 /// Get a range over the track states of this track. Return value is 0512 /// compatible with range based for loop. Const version 0513 /// @note This range is from the outside inwards! 0514 /// @return Track state range to iterate over 0515 auto trackStatesReversed() const { 0516 return m_container->reverseTrackStateRange(m_index); 0517 } 0518 0519 /// Get a range over the track states of this track. Return value is 0520 /// compatible with range based for loop. Mutable version 0521 /// @note Only available if the track proxy is not read-only 0522 /// @note This range is from the outside inwards! 0523 /// @return Track state range to iterate over 0524 auto trackStatesReversed() 0525 requires(!ReadOnly) 0526 { 0527 return m_container->reverseTrackStateRange(m_index); 0528 } 0529 0530 /// Get a range over the track states of this track. Return value is 0531 /// compatible with range based for loop. This overload returns a const-only 0532 /// track state range, which means you cannot modify the track states obtained 0533 /// in the iteration. 0534 /// @note This range is from the inside out! 0535 /// @warning This access direction is only possible if the track states are 0536 /// **forward-linked**. 0537 /// @return Track state range to iterate over 0538 auto trackStates() const { 0539 return m_container->forwardTrackStateRange(m_index); 0540 } 0541 0542 /// Get a range over the track states of this track. 0543 /// Return value is compatible with range based for loop. 0544 /// This overload returns a mutable track state range, which means you 0545 /// can modify the track states obtained in the iteration. 0546 /// @note Only available if the track proxy is not read-only 0547 /// @note This range is from the inside out! 0548 /// @warning This access direction is only possible if the track states are 0549 /// **forward-linked**. 0550 /// @return Track state range to iterate over 0551 auto trackStates() 0552 requires(!ReadOnly) 0553 { 0554 return m_container->forwardTrackStateRange(m_index); 0555 } 0556 0557 /// @} 0558 0559 /// @anchor track_proxy_track_state_manipulation 0560 /// @name TrackProxy track state manipulation 0561 /// Methods that manipulate the track states of a track represented by @c 0562 /// TrackProxy. 0563 /// 0564 /// **Copy Methods Overview:** 0565 /// 0566 /// Three main copy methods are available with different behaviors: 0567 /// - @c copyFrom(): Deep copy including all track states (creates new track 0568 /// states) 0569 /// - @c copyFromWithoutStates(): Copy only track properties, invalidate 0570 /// track state indices 0571 /// - @c copyFromShallow(): Shallow copy sharing the same track states (copy 0572 /// indices only) 0573 /// 0574 /// Choose based on your needs: 0575 /// - Use @c copyFrom() for independent track copies with separate track 0576 /// states 0577 /// - Use @c copyFromWithoutStates() to update track metadata without 0578 /// affecting trajectories 0579 /// - Use @c copyFromShallow() for lightweight copies when track states can 0580 /// be shared 0581 /// @{ 0582 0583 /// Forward connect a track. 0584 /// This means setting indices from the inside out on all track states. 0585 /// @note Only available if the track proxy is not read-only 0586 void linkForward() 0587 requires(!ReadOnly) 0588 { 0589 IndexType last = kInvalid; 0590 for (auto ts : trackStatesReversed()) { 0591 ts.template component<IndexType>(hashString("next")) = last; 0592 last = ts.index(); 0593 } 0594 stemIndex() = last; 0595 } 0596 0597 /// Append a track state to this track. 0598 /// This will modify the tip index to point at the newly created track state, 0599 /// which will be directly after the previous track state at tip index. 0600 /// @note Only available if the track proxy is not read-only 0601 /// @param mask The allocation prop mask for the new track state 0602 /// @return The newly added track state 0603 auto appendTrackState(TrackStatePropMask mask = TrackStatePropMask::All) 0604 requires(!ReadOnly) 0605 { 0606 auto& tsc = m_container->trackStateContainer(); 0607 auto ts = tsc.makeTrackState(mask, tipIndex()); 0608 tipIndex() = ts.index(); 0609 return ts; 0610 } 0611 0612 /// Copy the content of another track proxy into this one 0613 /// @note Only available if the track proxy is not read-only 0614 /// @tparam track_proxy_t the other track proxy's type 0615 /// @param other The track proxy 0616 /// @param copyTrackStates Copy the track state sequence from @p other 0617 template <TrackProxyConcept track_proxy_t> 0618 [[deprecated( 0619 "copyFrom() with copyTrackStates == false is deprecated, use " 0620 "copyFromWithoutStates()")]] 0621 void copyFrom(const track_proxy_t& other, bool copyTrackStates) 0622 requires(!ReadOnly) 0623 { 0624 if (copyTrackStates) { 0625 copyFrom(other); 0626 } else { 0627 // Emulate previous behavior 0628 auto prevTipIndex = tipIndex(); 0629 auto prevStemIndex = stemIndex(); 0630 copyFromWithoutStates(other); 0631 // Reset to previous values 0632 tipIndex() = prevTipIndex; 0633 stemIndex() = prevStemIndex; 0634 } 0635 } 0636 0637 /// Create a complete deep copy of another track, including all track states. 0638 /// This creates new track states in the destination trajectory and copies 0639 /// all data from the source track states. The track state sequence order 0640 /// is preserved. 0641 /// 0642 /// **Implementation details:** 0643 /// - Track states are initially copied in reversed order for efficiency 0644 /// - The track state links are then updated using reverseTrackStates() 0645 /// - As a consequence, the resulting track is forward-linked 0646 /// 0647 /// **What gets copied:** 0648 /// - All track-level properties (parameters, covariance, particle hypothesis, 0649 /// etc.) 0650 /// - Reference surface (shared pointer is copied) 0651 /// - Track summary data (nMeasurements, nHoles, chi2, etc.) 0652 /// - All dynamic track columns 0653 /// - Complete sequence of track states with all their data 0654 /// - All dynamic track state columns 0655 /// 0656 /// **Result:** 0657 /// - The destination track will have newly created track states 0658 /// - tipIndex() and stemIndex() will point to the new track states 0659 /// - Track state indices will be different from the source 0660 /// - All track state data will be identical to the source 0661 /// - The track will be forward-linked (stemIndex() will be valid) 0662 /// 0663 /// @note Only available if the track proxy is not read-only 0664 /// @note Both track containers must have compatible dynamic columns 0665 /// @tparam track_proxy_t the other track proxy's type 0666 /// @param other The source track proxy to copy from 0667 template <TrackProxyConcept track_proxy_t> 0668 void copyFrom(const track_proxy_t& other) 0669 requires(!ReadOnly) 0670 { 0671 copyFromWithoutStates(other); 0672 0673 // append track states (cheap), but they're in the wrong order 0674 for (const auto& srcTrackState : other.trackStatesReversed()) { 0675 auto destTrackState = appendTrackState(srcTrackState.getMask()); 0676 destTrackState.copyFrom(srcTrackState, Acts::TrackStatePropMask::All, 0677 true); 0678 } 0679 0680 // reverse using standard linked list reversal algorithm 0681 reverseTrackStates(); 0682 } 0683 0684 /// Copy track-level properties from another track, but not the track states. 0685 /// This copies all track metadata and properties but leaves the track state 0686 /// sequence unchanged. Useful when you want to copy track properties to an 0687 /// existing track that may already have track states. 0688 /// 0689 /// **What gets copied:** 0690 /// - Track parameters at reference surface 0691 /// - Covariance matrix at reference surface 0692 /// - Particle hypothesis 0693 /// - Reference surface (shared pointer is copied) 0694 /// - Track summary data (nMeasurements, nHoles, nOutliers, nSharedHits, chi2, 0695 /// nDoF) 0696 /// - All dynamic track columns 0697 /// 0698 /// **What does NOT get copied:** 0699 /// - Track states (existing track states remain unchanged in the container) 0700 /// 0701 /// **Result:** 0702 /// - All track-level properties are updated to match the source 0703 /// - tipIndex() and stemIndex() are set to kInvalid (track states become 0704 /// inaccessible) 0705 /// - Existing track states remain in the container but are no longer linked 0706 /// to this track 0707 /// - nTrackStates() will return 0 due to invalid indices 0708 /// 0709 /// @note Only available if the track proxy is not read-only 0710 /// @note Both track containers must have compatible dynamic columns 0711 /// @tparam track_proxy_t the other track proxy's type 0712 /// @param other The source track proxy to copy properties from 0713 template <TrackProxyConcept track_proxy_t> 0714 void copyFromWithoutStates(const track_proxy_t& other) 0715 requires(!ReadOnly) 0716 { 0717 setParticleHypothesis(other.particleHypothesis()); 0718 0719 if (other.hasReferenceSurface()) { 0720 setReferenceSurface(other.referenceSurface().getSharedPtr()); 0721 parameters() = other.parameters(); 0722 covariance() = other.covariance(); 0723 } else { 0724 setReferenceSurface(nullptr); 0725 } 0726 0727 nMeasurements() = other.nMeasurements(); 0728 nHoles() = other.nHoles(); 0729 nOutliers() = other.nOutliers(); 0730 nSharedHits() = other.nSharedHits(); 0731 chi2() = other.chi2(); 0732 nDoF() = other.nDoF(); 0733 0734 m_container->copyDynamicFrom(m_index, other.m_container->container(), 0735 other.m_index); 0736 0737 tipIndex() = kInvalid; 0738 stemIndex() = kInvalid; 0739 } 0740 0741 /// Creates a *shallow copy* of the track. Track states are not copied, but 0742 /// the resulting track points at the same track states as the original. 0743 /// @note Only available if the track proxy is not read-only 0744 /// @return A shallow copy of this track proxy 0745 [[deprecated("Use copyFromShallow() instead")]] 0746 TrackProxy shallowCopy() 0747 requires(!ReadOnly) 0748 { 0749 auto t = container().makeTrack(); 0750 t.copyFromShallow(*this); 0751 return t; 0752 } 0753 0754 /// Create a shallow copy from another track, sharing the same track states. 0755 /// This copies all track-level properties and makes the destination track 0756 /// point to the same track state sequence as the source. The track states 0757 /// themselves are not duplicated - both tracks will reference the same 0758 /// track state objects in memory. 0759 /// 0760 /// **What gets copied:** 0761 /// - All track-level properties (parameters, covariance, particle hypothesis, 0762 /// etc.) 0763 /// - Reference surface (shared pointer is copied) 0764 /// - Track summary data (nMeasurements, nHoles, chi2, etc.) 0765 /// - All dynamic track columns 0766 /// - tipIndex() and stemIndex() (track state linking information) 0767 /// 0768 /// **What gets shared (not duplicated):** 0769 /// - Track states (both tracks reference the same track state objects) 0770 /// 0771 /// **Result:** 0772 /// - The destination track will have the same nTrackStates() as the source 0773 /// - Both tracks will iterate over the same track state sequence 0774 /// - Modifications to track states will be visible in both tracks 0775 /// - Track state indices will be identical between tracks 0776 /// - The destination track will have a different track index than the source 0777 /// 0778 /// @warning Modifying track states through either track will affect both tracks 0779 /// since they share the same track state objects 0780 /// @warning It is the user's responsibility to ensure that the tip and stem 0781 /// indices from the source track are valid in the destination 0782 /// track's track state container. No validation is performed - 0783 /// invalid indices will lead to undefined behavior when accessing 0784 /// track states 0785 /// @note Only available if the track proxy is not read-only 0786 /// @note Both track containers must have compatible dynamic columns 0787 /// @tparam track_proxy_t the other track proxy's type 0788 /// @param other The source track proxy to create a shallow copy from 0789 template <TrackProxyConcept track_proxy_t> 0790 void copyFromShallow(const track_proxy_t& other) 0791 requires(!ReadOnly) 0792 { 0793 copyFromWithoutStates(other); 0794 tipIndex() = other.tipIndex(); 0795 stemIndex() = other.stemIndex(); 0796 } 0797 0798 /// Reverse the ordering of track states for this track 0799 /// Afterwards, the previous endpoint of the track state sequence will be 0800 /// the "innermost" track state 0801 /// @note Only available if the track proxy is not read-only 0802 /// @note This is dangerous with branching track state sequences, as it will break them 0803 /// @note This also automatically forward-links the track! 0804 /// @param invertJacobians Whether to invert the Jacobians of the track states 0805 void reverseTrackStates(bool invertJacobians = false) 0806 requires(!ReadOnly) 0807 { 0808 IndexType current = tipIndex(); 0809 IndexType next = kInvalid; 0810 IndexType prev = kInvalid; 0811 0812 stemIndex() = tipIndex(); 0813 0814 // @TODO: Maybe refactor to not need this variable if invertJacobians == false 0815 BoundMatrix nextJacobian; 0816 0817 while (current != kInvalid) { 0818 auto ts = m_container->trackStateContainer().getTrackState(current); 0819 prev = ts.previous(); 0820 ts.template component<IndexType>(hashString("next")) = prev; 0821 ts.previous() = next; 0822 if (invertJacobians) { 0823 if (next != kInvalid) { 0824 BoundMatrix curJacobian = ts.jacobian(); 0825 ts.jacobian() = nextJacobian.inverse(); 0826 nextJacobian = curJacobian; 0827 } else { 0828 nextJacobian = ts.jacobian(); 0829 ts.jacobian().setZero(); 0830 } 0831 } 0832 next = current; 0833 tipIndex() = current; 0834 current = prev; 0835 } 0836 } 0837 0838 /// @} 0839 0840 /// @anchor track_proxy_generic_component 0841 /// @name TrackProxy generic component access 0842 /// Methods that give access to generic components of a track represented by 0843 /// @c TrackProxy. Internally, a compile-time hash of the component name is 0844 /// used to identify which component is being requested. Most of the named 0845 /// methods in @ref track_proxy_props "TrackProxy properties" use these 0846 /// methods to retrieve the actual data. 0847 /// 0848 /// A number of overloads exist, where you can either supply the 0849 /// @ref HashedString @c key as a template parameter or a runtime argument. The 0850 /// former has the advantage of being guaranteed to be evaluated at 0851 /// compile-time. 0852 /// 0853 /// @{ 0854 0855 /// Retrieve a mutable reference to a component 0856 /// @tparam T The type of the component to access 0857 /// @tparam key String key for the component to access 0858 /// @return Mutable reference to the component given by @p key 0859 template <typename T, HashedString key> 0860 constexpr T& component() 0861 requires(!ReadOnly) 0862 { 0863 return m_container->template component<T, key>(m_index); 0864 } 0865 0866 /// Retrieve a mutable reference to a component 0867 /// @tparam T The type of the component to access 0868 /// @param key String key for the component to access 0869 /// @return Mutable reference to the component given by @p key 0870 template <typename T> 0871 constexpr T& component(HashedString key) 0872 requires(!ReadOnly) 0873 { 0874 return m_container->template component<T>(key, m_index); 0875 } 0876 0877 /// Retrieve a mutable reference to a component 0878 /// @tparam T The type of the component to access 0879 /// @param key String key for the component to access 0880 /// @note This might hash the @p key at runtime instead of compile-time 0881 /// @return Mutable reference to the component given by @p key 0882 template <typename T> 0883 constexpr T& component(std::string_view key) 0884 requires(!ReadOnly) 0885 { 0886 return m_container->template component<T>(hashStringDynamic(key), m_index); 0887 } 0888 0889 /// Retrieve a const reference to a component 0890 /// @tparam T The type of the component to access 0891 /// @tparam key String key for the component to access 0892 /// @return Const reference to the component given by @p key 0893 template <typename T, HashedString key> 0894 constexpr const T& component() const { 0895 return m_container->template component<T, key>(m_index); 0896 } 0897 0898 /// Retrieve a const reference to a component 0899 /// @tparam T The type of the component to access 0900 /// @param key String key for the component to access 0901 /// @return Const reference to the component given by @p key 0902 template <typename T> 0903 constexpr const T& component(HashedString key) const { 0904 return m_container->template component<T>(key, m_index); 0905 } 0906 0907 /// Retrieve a const reference to a component 0908 /// @tparam T The type of the component to access 0909 /// @param key String key for the component to access 0910 /// @note This might hash the @p key at runtime instead of compile-time 0911 /// @return Const reference to the component given by @p key 0912 template <typename T> 0913 constexpr const T& component(std::string_view key) const { 0914 return m_container->template component<T>(hashStringDynamic(key), m_index); 0915 } 0916 0917 /// @} 0918 0919 /// Return the track parameters at the reference surface 0920 /// @note The parameters are created on the fly 0921 /// @return the track parameters 0922 BoundTrackParameters createParametersAtReference() const { 0923 return BoundTrackParameters(referenceSurface().getSharedPtr(), parameters(), 0924 covariance(), particleHypothesis()); 0925 } 0926 0927 /// Convert a track state into track parameters 0928 /// @note The parameters are created on the fly 0929 /// @param trackState Track state to convert to parameters 0930 /// @return the track parameters 0931 BoundTrackParameters createParametersFromState( 0932 const ConstTrackStateProxy& trackState) const { 0933 return BoundTrackParameters(trackState.referenceSurface().getSharedPtr(), 0934 trackState.parameters(), 0935 trackState.covariance(), particleHypothesis()); 0936 } 0937 0938 /// Return a reference to the track container backend, mutable version. 0939 /// @note Only available if the track proxy is not read-only 0940 /// @return reference to the track container backend 0941 auto& container() 0942 requires(!ReadOnly) 0943 { 0944 return *m_container; 0945 } 0946 0947 /// Return a reference to the track container backend, const version. 0948 /// @return reference to the track container backend 0949 const auto& container() const { return *m_container; } 0950 0951 private: 0952 TrackProxy( 0953 const_if_t<ReadOnly, TrackContainer<Container, Trajectory, holder_t>>& 0954 container, 0955 IndexType itrack) 0956 : m_container{&container}, m_index{itrack} {} 0957 0958 detail_lt::TransitiveConstPointer< 0959 const_if_t<ReadOnly, TrackContainer<Container, Trajectory, holder_t>>> 0960 m_container; 0961 IndexType m_index; 0962 }; 0963 0964 } // namespace Acts
| [ Source navigation ] | [ Diff markup ] | [ Identifier search ] | [ general search ] |
|
This page was automatically generated by the 2.3.7 LXR engine. The LXR team |
|