![]() |
|
|||
File indexing completed on 2025-07-01 07:52:23
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 /// @return True if the track proxies refer to the same track 0147 bool operator==(const TrackProxy& other) const { 0148 return &(*m_container) == &(*other.m_container) && m_index == other.m_index; 0149 } 0150 0151 /// @anchor track_proxy_props 0152 /// @name TrackProxy properties 0153 /// Methods that give access to the properties of a track represented by 0154 /// @c TrackProxy. 0155 /// 0156 /// Many of these methods come in a @c const and a non-@c const version. The 0157 /// non-@c const version is only available if you have an instance of 0158 /// @c TrackProxy that does not have the @c read_only template parameter set to 0159 /// @c true, even if you hold it as an lvalue. 0160 /// 0161 /// @{ 0162 0163 /// Get the tip index, i.e. the entry point into the track state container 0164 /// @return the tip index by value 0165 IndexType tipIndex() const { 0166 return component<IndexType>(hashString("tipIndex")); 0167 } 0168 0169 /// Index of the stem, i.e. the innermost track state of the track. 0170 /// This might be invalid, signifying that the track state is not 0171 /// forward-linked. 0172 /// @return the stem index 0173 IndexType stemIndex() const { 0174 return component<IndexType>(hashString("stemIndex")); 0175 } 0176 0177 /// Get a mutable reference to the tip index, i.e. the entry point into the 0178 /// track container 0179 /// @note Only available if the track proxy is not read-only 0180 /// @return mutable reference to the tip index 0181 IndexType& tipIndex() 0182 requires(!ReadOnly) 0183 { 0184 return component<IndexType>(hashString("tipIndex")); 0185 } 0186 0187 /// Index of the stem, i.e. the innermost track state of the track. 0188 /// This might be invalid, signifying that the track state is not 0189 /// forward-linked. 0190 /// @note Only available if the track proxy is not read-only 0191 /// @return mutable reference to the stem index 0192 IndexType& stemIndex() 0193 requires(!ReadOnly) 0194 { 0195 return component<IndexType>(hashString("stemIndex")); 0196 } 0197 0198 /// Get the reference surface of the track (e.g. the perigee) 0199 /// @return the reference surface 0200 const Surface& referenceSurface() const { 0201 return *m_container->container().referenceSurface_impl(m_index); 0202 } 0203 0204 // NOLINTBEGIN(performance-unnecessary-value-param) 0205 // looks like a false-positive. clang-tidy believes `srf` is not movable. 0206 /// Set a new reference surface for this track 0207 /// @param srf The surface to set 0208 void setReferenceSurface(std::shared_ptr<const Surface> srf) 0209 requires(!ReadOnly) 0210 { 0211 m_container->container().setReferenceSurface_impl(m_index, std::move(srf)); 0212 } 0213 // NOLINTEND(performance-unnecessary-value-param) 0214 0215 /// Returns whether the track has a reference surface or not 0216 /// @return whether a surface exists or not 0217 bool hasReferenceSurface() const { 0218 // @TODO: This could be more efficient 0219 return m_container->container().referenceSurface_impl(m_index) != nullptr; 0220 } 0221 0222 /// Get the parameters of the track at the reference surface (e.g. perigee). 0223 /// Const version 0224 /// @return Proxy vector for the parameters 0225 ConstParameters parameters() const { 0226 return m_container->parameters(m_index); 0227 } 0228 0229 /// Get the covariance of the track at the reference surface (e.g. perigee). 0230 /// Const version 0231 /// @return Proxy matrix for the covariance 0232 ConstCovariance covariance() const { 0233 return m_container->covariance(m_index); 0234 } 0235 0236 /// Get the parameters of the track at the reference surface (e.g. perigee). 0237 /// Mutable version 0238 /// @note Only available if the track proxy is not read-only 0239 /// @return Proxy vector for the parameters 0240 Parameters parameters() 0241 requires(!ReadOnly) 0242 { 0243 return m_container->parameters(m_index); 0244 } 0245 0246 /// Get the covariance of the track at the reference surface (e.g. perigee). 0247 /// Mutable version 0248 /// @note Only available if the track proxy is not read-only 0249 /// @return Proxy matrix for the covariance 0250 Covariance covariance() 0251 requires(!ReadOnly) 0252 { 0253 return m_container->covariance(m_index); 0254 } 0255 0256 /// Access the theta parameter of the track at the reference surface 0257 /// @return The theta parameter 0258 double theta() const { return parameters()[eBoundTheta]; } 0259 0260 /// Access the phi parameter of the track at the reference surface 0261 /// @return The phi parameter 0262 double phi() const { return parameters()[eBoundPhi]; } 0263 0264 /// Access the loc0 parameter of the track at the reference surface 0265 /// @return The loc0 parameter 0266 double loc0() const { return parameters()[eBoundLoc0]; } 0267 0268 /// Access the loc1 parameter of the track at the reference surface 0269 /// @return The loc1 parameter 0270 double loc1() const { return parameters()[eBoundLoc1]; } 0271 0272 /// Access the time parameter of the track at the reference surface 0273 /// @return The time parameter 0274 double time() const { return parameters()[eBoundTime]; } 0275 0276 /// Access the q/p (curvature) parameter of the track at the reference surface 0277 /// @return The q/p parameter 0278 double qOverP() const { return parameters()[eBoundQOverP]; } 0279 0280 /// Get the particle hypothesis 0281 /// @return the particle hypothesis 0282 ParticleHypothesis particleHypothesis() const { 0283 return m_container->container().particleHypothesis_impl(m_index); 0284 } 0285 0286 /// Set a new particle hypothesis for this track 0287 /// @note Only available if the track proxy is not read-only 0288 /// @param particleHypothesis The particle hypothesis to set 0289 void setParticleHypothesis(const ParticleHypothesis& particleHypothesis) 0290 requires(!ReadOnly) 0291 { 0292 m_container->container().setParticleHypothesis_impl(m_index, 0293 particleHypothesis); 0294 } 0295 0296 /// Get the charge of the tack 0297 /// @note this depends on the charge hypothesis 0298 /// @return The absolute track momentum 0299 double charge() const { return particleHypothesis().qFromQOP(qOverP()); } 0300 0301 /// Get the absolute momentum of the tack 0302 /// @return The absolute track momentum 0303 double absoluteMomentum() const { 0304 return particleHypothesis().extractMomentum(qOverP()); 0305 } 0306 0307 /// Get the transverse momentum of the track 0308 /// @return The track transverse momentum value 0309 double transverseMomentum() const { 0310 return std::sin(theta()) * absoluteMomentum(); 0311 } 0312 0313 /// Get a unit vector along the track direction at the reference surface 0314 /// @return The direction unit vector 0315 Vector3 direction() const { 0316 return makeDirectionFromPhiTheta(phi(), theta()); 0317 } 0318 0319 /// Get the global momentum vector 0320 /// @return the global momentum vector 0321 Vector3 momentum() const { return absoluteMomentum() * direction(); } 0322 0323 /// Return the number of track states associated to this track 0324 /// @note This is calculated by iterating over the track states which is 0325 /// somewhat expensive. Consider caching this value if you need It 0326 /// more than once. 0327 /// @return The number of track states 0328 unsigned int nTrackStates() const { 0329 // @TODO: This should probably be cached, distance is expensive 0330 // without random access 0331 if (tipIndex() == kInvalid) { 0332 // no tip index -> no track states 0333 return 0; 0334 } 0335 auto tsRange = trackStatesReversed(); 0336 return std::distance(tsRange.begin(), tsRange.end()); 0337 } 0338 0339 /// Return a mutable reference to the number of measurements for the track. 0340 /// Mutable version 0341 /// @note Only available if the track proxy is not read-only 0342 /// @return The number of measurements 0343 unsigned int& nMeasurements() 0344 requires(!ReadOnly) 0345 { 0346 return component<unsigned int, hashString("nMeasurements")>(); 0347 } 0348 0349 /// Return the number of measurements for the track. Const version 0350 /// @return The number of measurements 0351 unsigned int nMeasurements() const { 0352 return component<unsigned int, hashString("nMeasurements")>(); 0353 } 0354 0355 /// Return a mutable reference to the number of holes for the track. 0356 /// Mutable version 0357 /// @note Only available if the track proxy is not read-only 0358 /// @return The number of holes 0359 unsigned int& nHoles() 0360 requires(!ReadOnly) 0361 { 0362 return component<unsigned int, hashString("nHoles")>(); 0363 } 0364 0365 /// Return the number of measurements for the track. Const version 0366 /// @return The number of measurements 0367 unsigned int nHoles() const { 0368 return component<unsigned int, hashString("nHoles")>(); 0369 } 0370 0371 /// Return a mutable reference to the number of outliers for the track. 0372 /// Mutable version 0373 /// @note Only available if the track proxy is not read-only 0374 /// @return The number of outliers 0375 unsigned int& nOutliers() 0376 requires(!ReadOnly) 0377 { 0378 return component<unsigned int, hashString("nOutliers")>(); 0379 } 0380 0381 /// Return the number of outliers for the track. Const version 0382 /// @return The number of outliers 0383 unsigned int nOutliers() const { 0384 return component<unsigned int, hashString("nOutliers")>(); 0385 } 0386 0387 /// Return a mutable reference to the number of shared hits for the track. 0388 /// Mutable version 0389 /// @note Only available if the track proxy is not read-only 0390 /// @return The number of shared hits 0391 unsigned int& nSharedHits() 0392 requires(!ReadOnly) 0393 { 0394 return component<unsigned int, hashString("nSharedHits")>(); 0395 } 0396 0397 /// Return the number of shared hits for the track. Const version 0398 /// @return The number of shared hits 0399 unsigned int nSharedHits() const { 0400 return component<unsigned int, hashString("nSharedHits")>(); 0401 } 0402 0403 /// Return a mutable reference to the chi squared 0404 /// Mutable version 0405 /// @note Only available if the track proxy is not read-only 0406 /// @return The chi squared 0407 float& chi2() 0408 requires(!ReadOnly) 0409 { 0410 return component<float, hashString("chi2")>(); 0411 } 0412 0413 /// Return the chi squared for the track. Const version 0414 /// @return The chi squared 0415 float chi2() const { return component<float, hashString("chi2")>(); } 0416 0417 /// Return a mutable reference to the number of degrees of freedom for the 0418 /// track. Mutable version 0419 /// @note Only available if the track proxy is not read-only 0420 /// @return The number of degrees of freedom 0421 unsigned int& nDoF() 0422 requires(!ReadOnly) 0423 { 0424 return component<unsigned int, hashString("ndf")>(); 0425 } 0426 0427 /// Return the number of degrees of freedom for the track. Const version 0428 /// @return The number of degrees of freedom 0429 unsigned int nDoF() const { 0430 return component<unsigned int, hashString("ndf")>(); 0431 } 0432 0433 /// Return the index of this track in the track container 0434 /// @note This is separate from the tip index 0435 /// @return the track index 0436 IndexType index() const { return m_index; } 0437 0438 /// @} 0439 0440 /// @anchor track_proxy_track_states 0441 /// @name TrackProxy track state access 0442 /// Methods that give access to the track states of a track represented by @c TrackProxy. 0443 /// @{ 0444 0445 /// Return a const track state proxy to the outermost track state 0446 /// @return The outermost track state proxy 0447 ConstTrackStateProxy outermostTrackState() const { 0448 return m_container->trackStateContainer().getTrackState(tipIndex()); 0449 } 0450 0451 /// Return a mutable track state proxy to the outermost track state 0452 /// @return The outermost track state proxy 0453 TrackStateProxy outermostTrackState() 0454 requires(!ReadOnly) 0455 { 0456 return m_container->trackStateContainer().getTrackState(tipIndex()); 0457 } 0458 0459 /// Return a const track state proxy to the innermost track state 0460 /// @note This is only available, if the track is forward linked 0461 /// @return The innermost track state proxy 0462 auto innermostTrackState() const { 0463 using proxy_t = decltype(m_container->trackStateContainer().getTrackState( 0464 std::declval<IndexType>())); 0465 0466 IndexType stem = component<IndexType, hashString("stemIndex")>(); 0467 if (stem == kInvalid) { 0468 return std::optional<proxy_t>{}; 0469 } else { 0470 return std::optional<proxy_t>{ 0471 m_container->trackStateContainer().getTrackState(stem)}; 0472 } 0473 } 0474 0475 /// Return a mutable track state proxy to the innermost track state 0476 /// @note This is only available, if the track is forward linked 0477 /// @note Only available if the track proxy is not read-only 0478 /// @return The innermost track state proxy 0479 auto innermostTrackState() 0480 requires(!ReadOnly) 0481 { 0482 using proxy_t = decltype(m_container->trackStateContainer().getTrackState( 0483 std::declval<IndexType>())); 0484 0485 IndexType stem = component<IndexType>(hashString("stemIndex")); 0486 if (stem == kInvalid) { 0487 return std::optional<proxy_t>{}; 0488 } else { 0489 return std::optional<proxy_t>{ 0490 m_container->trackStateContainer().getTrackState(stem)}; 0491 } 0492 } 0493 0494 /// Get a range over the track states of this track. Return value is 0495 /// compatible with range based for loop. Const version 0496 /// @note This range is from the outside inwards! 0497 /// @return Track state range to iterate over 0498 auto trackStatesReversed() const { 0499 return m_container->reverseTrackStateRange(m_index); 0500 } 0501 0502 /// Get a range over the track states of this track. Return value is 0503 /// compatible with range based for loop. Mutable version 0504 /// @note Only available if the track proxy is not read-only 0505 /// @note This range is from the outside inwards! 0506 /// @return Track state range to iterate over 0507 auto trackStatesReversed() 0508 requires(!ReadOnly) 0509 { 0510 return m_container->reverseTrackStateRange(m_index); 0511 } 0512 0513 /// Get a range over the track states of this track. Return value is 0514 /// compatible with range based for loop. This overload returns a const-only 0515 /// track state range, which means you cannot modify the track states obtained 0516 /// in the iteration. 0517 /// @note This range is from the inside out! 0518 /// @warning This access direction is only possible if the track states are 0519 /// **forward-linked**. 0520 /// @return Track state range to iterate over 0521 auto trackStates() const { 0522 return m_container->forwardTrackStateRange(m_index); 0523 } 0524 0525 /// Get a range over the track states of this track. 0526 /// Return value is compatible with range based for loop. 0527 /// This overload returns a mutable track state range, which means you 0528 /// can modify the track states obtained in the iteration. 0529 /// @note Only available if the track proxy is not read-only 0530 /// @note This range is from the inside out! 0531 /// @warning This access direction is only possible if the track states are 0532 /// **forward-linked**. 0533 /// @return Track state range to iterate over 0534 auto trackStates() 0535 requires(!ReadOnly) 0536 { 0537 return m_container->forwardTrackStateRange(m_index); 0538 } 0539 0540 /// @} 0541 0542 /// @anchor track_proxy_track_state_manipulation 0543 /// @name TrackProxy track state manipulation 0544 /// Methods that manipulate the track states of a track represented by @c TrackProxy. 0545 /// @{ 0546 0547 /// Forward connect a track. 0548 /// This means setting indices from the inside out on all track states. 0549 /// @note Only available if the track proxy is not read-only 0550 void linkForward() 0551 requires(!ReadOnly) 0552 { 0553 IndexType last = kInvalid; 0554 for (auto ts : trackStatesReversed()) { 0555 ts.template component<IndexType>(hashString("next")) = last; 0556 last = ts.index(); 0557 } 0558 stemIndex() = last; 0559 } 0560 0561 /// Append a track state to this track. 0562 /// This will modify the tip index to point at the newly created track state, 0563 /// which will be directly after the previous track state at tip index. 0564 /// @note Only available if the track proxy is not read-only 0565 /// @param mask The allocation prop mask for the new track state 0566 /// @return The newly added track state 0567 auto appendTrackState(TrackStatePropMask mask = TrackStatePropMask::All) 0568 requires(!ReadOnly) 0569 { 0570 auto& tsc = m_container->trackStateContainer(); 0571 auto ts = tsc.makeTrackState(mask, tipIndex()); 0572 tipIndex() = ts.index(); 0573 return ts; 0574 } 0575 0576 /// Copy the content of another track proxy into this one 0577 /// @note Only available if the track proxy is not read-only 0578 /// @tparam track_proxy_t the other track proxy's type 0579 /// @param other The track proxy 0580 /// @param copyTrackStates Copy the track state sequence from @p other 0581 template <TrackProxyConcept track_proxy_t> 0582 void copyFrom(const track_proxy_t& other, bool copyTrackStates = true) 0583 requires(!ReadOnly) 0584 { 0585 // @TODO: Add constraint on which track proxies are allowed, 0586 // this is only implicit right now 0587 0588 if (copyTrackStates) { 0589 // append track states (cheap), but they're in the wrong order 0590 for (const auto& srcTrackState : other.trackStatesReversed()) { 0591 auto destTrackState = appendTrackState(srcTrackState.getMask()); 0592 destTrackState.copyFrom(srcTrackState, Acts::TrackStatePropMask::All, 0593 true); 0594 } 0595 0596 // reverse using standard linked list reversal algorithm 0597 reverseTrackStates(); 0598 } 0599 0600 setParticleHypothesis(other.particleHypothesis()); 0601 0602 if (other.hasReferenceSurface()) { 0603 setReferenceSurface(other.referenceSurface().getSharedPtr()); 0604 parameters() = other.parameters(); 0605 covariance() = other.covariance(); 0606 } else { 0607 setReferenceSurface(nullptr); 0608 } 0609 0610 nMeasurements() = other.nMeasurements(); 0611 nHoles() = other.nHoles(); 0612 nOutliers() = other.nOutliers(); 0613 nSharedHits() = other.nSharedHits(); 0614 chi2() = other.chi2(); 0615 nDoF() = other.nDoF(); 0616 0617 m_container->copyDynamicFrom(m_index, other.m_container->container(), 0618 other.m_index); 0619 } 0620 0621 /// Creates a *shallow copy* of the track. Track states are not copied, but 0622 /// the resulting track points at the same track states as the original. 0623 /// @note Only available if the track proxy is not read-only 0624 TrackProxy shallowCopy() 0625 requires(!ReadOnly) 0626 { 0627 auto ts = container().makeTrack(); 0628 ts.copyFrom(*this, false); 0629 ts.tipIndex() = tipIndex(); 0630 ts.stemIndex() = stemIndex(); 0631 return ts; 0632 } 0633 0634 /// Reverse the ordering of track states for this track 0635 /// Afterwards, the previous endpoint of the track state sequence will be the 0636 /// "innermost" track state 0637 /// @note Only available if the track proxy is not read-only 0638 /// @note This is dangerous with branching track state sequences, as it will break them 0639 /// @note This also automatically forward-links the track! 0640 /// @param invertJacobians Whether to invert the Jacobians of the track states 0641 void reverseTrackStates(bool invertJacobians = false) 0642 requires(!ReadOnly) 0643 { 0644 IndexType current = tipIndex(); 0645 IndexType next = kInvalid; 0646 IndexType prev = kInvalid; 0647 0648 stemIndex() = tipIndex(); 0649 0650 // @TODO: Maybe refactor to not need this variable if invertJacobians == false 0651 BoundMatrix nextJacobian; 0652 0653 while (current != kInvalid) { 0654 auto ts = m_container->trackStateContainer().getTrackState(current); 0655 prev = ts.previous(); 0656 ts.template component<IndexType>(hashString("next")) = prev; 0657 ts.previous() = next; 0658 if (invertJacobians) { 0659 if (next != kInvalid) { 0660 BoundMatrix curJacobian = ts.jacobian(); 0661 ts.jacobian() = nextJacobian.inverse(); 0662 nextJacobian = curJacobian; 0663 } else { 0664 nextJacobian = ts.jacobian(); 0665 ts.jacobian().setZero(); 0666 } 0667 } 0668 next = current; 0669 tipIndex() = current; 0670 current = prev; 0671 } 0672 } 0673 0674 /// @} 0675 0676 /// @anchor track_proxy_generic_component 0677 /// @name TrackProxy generic component access 0678 /// Methods that give access to generic components of a track represented by 0679 /// @c TrackProxy. Internally, a compile-time hash of the component name is 0680 /// used to identify which component is being requested. Most of the named 0681 /// methods in @ref track_proxy_props "TrackProxy properties" use these 0682 /// methods to retrieve the actual data. 0683 /// 0684 /// A number of overloads exist, where you can either supply the 0685 /// @ref HashedString @c key as a template parameter or a runtime argument. The 0686 /// former has the advantage of being guaranteed to be evaluated at 0687 /// compile-time. 0688 /// 0689 /// @{ 0690 0691 /// Retrieve a mutable reference to a component 0692 /// @tparam T The type of the component to access 0693 /// @tparam key String key for the component to access 0694 /// @return Mutable reference to the component given by @p key 0695 template <typename T, HashedString key> 0696 constexpr T& component() 0697 requires(!ReadOnly) 0698 { 0699 return m_container->template component<T, key>(m_index); 0700 } 0701 0702 /// Retrieve a mutable reference to a component 0703 /// @tparam T The type of the component to access 0704 /// @param key String key for the component to access 0705 /// @return Mutable reference to the component given by @p key 0706 template <typename T> 0707 constexpr T& component(HashedString key) 0708 requires(!ReadOnly) 0709 { 0710 return m_container->template component<T>(key, m_index); 0711 } 0712 0713 /// Retrieve a mutable reference to a component 0714 /// @tparam T The type of the component to access 0715 /// @param key String key for the component to access 0716 /// @note This might hash the @p key at runtime instead of compile-time 0717 /// @return Mutable reference to the component given by @p key 0718 template <typename T> 0719 constexpr T& component(std::string_view key) 0720 requires(!ReadOnly) 0721 { 0722 return m_container->template component<T>(hashStringDynamic(key), m_index); 0723 } 0724 0725 /// Retrieve a const reference to a component 0726 /// @tparam T The type of the component to access 0727 /// @tparam key String key for the component to access 0728 /// @return Const reference to the component given by @p key 0729 template <typename T, HashedString key> 0730 constexpr const T& component() const { 0731 return m_container->template component<T, key>(m_index); 0732 } 0733 0734 /// Retrieve a const reference to a component 0735 /// @tparam T The type of the component to access 0736 /// @param key String key for the component to access 0737 /// @return Const reference to the component given by @p key 0738 template <typename T> 0739 constexpr const T& component(HashedString key) const { 0740 return m_container->template component<T>(key, m_index); 0741 } 0742 0743 /// Retrieve a const reference to a component 0744 /// @tparam T The type of the component to access 0745 /// @param key String key for the component to access 0746 /// @note This might hash the @p key at runtime instead of compile-time 0747 /// @return Const reference to the component given by @p key 0748 template <typename T> 0749 constexpr const T& component(std::string_view key) const { 0750 return m_container->template component<T>(hashStringDynamic(key), m_index); 0751 } 0752 0753 /// @} 0754 0755 /// Return the track parameters at the reference surface 0756 /// @note The parameters are created on the fly 0757 /// @return the track parameters 0758 BoundTrackParameters createParametersAtReference() const { 0759 return BoundTrackParameters(referenceSurface().getSharedPtr(), parameters(), 0760 covariance(), particleHypothesis()); 0761 } 0762 0763 /// Convert a track state into track parameters 0764 /// @note The parameters are created on the fly 0765 /// @return the track parameters 0766 BoundTrackParameters createParametersFromState( 0767 const ConstTrackStateProxy& trackState) const { 0768 return BoundTrackParameters(trackState.referenceSurface().getSharedPtr(), 0769 trackState.parameters(), 0770 trackState.covariance(), particleHypothesis()); 0771 } 0772 0773 /// Return a reference to the track container backend, mutable version. 0774 /// @note Only available if the track proxy is not read-only 0775 /// @return reference to the track container backend 0776 auto& container() 0777 requires(!ReadOnly) 0778 { 0779 return *m_container; 0780 } 0781 0782 /// Return a reference to the track container backend, const version. 0783 /// @return reference to the track container backend 0784 const auto& container() const { return *m_container; } 0785 0786 private: 0787 TrackProxy( 0788 const_if_t<ReadOnly, TrackContainer<Container, Trajectory, holder_t>>& 0789 container, 0790 IndexType itrack) 0791 : m_container{&container}, m_index{itrack} {} 0792 0793 detail_lt::TransitiveConstPointer< 0794 const_if_t<ReadOnly, TrackContainer<Container, Trajectory, holder_t>>> 0795 m_container; 0796 IndexType m_index; 0797 }; 0798 0799 } // 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 |
![]() ![]() |