Analysis Software
Documentation for sPHENIX simulation software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Acts::LineSurface Class Reference

#include <acts/blob/sPHENIX/Core/include/Acts/Surfaces/LineSurface.hpp>

+ Inheritance diagram for Acts::LineSurface:
+ Collaboration diagram for Acts::LineSurface:

Public Member Functions

 ~LineSurface () override=default
 
 LineSurface ()=delete
 
LineSurfaceoperator= (const LineSurface &other)
 
Vector3 normal (const GeometryContext &gctx, const Vector2 &lposition) const final
 
Vector3 binningPosition (const GeometryContext &gctx, BinningValue bValue) const final
 
RotationMatrix3 referenceFrame (const GeometryContext &gctx, const Vector3 &position, const Vector3 &direction) const final
 
BoundToFreeMatrix boundToFreeJacobian (const GeometryContext &gctx, const BoundVector &boundParams) const final
 
FreeToPathMatrix freeToPathDerivative (const GeometryContext &gctx, const FreeVector &parameters) const final
 
Vector3 localToGlobal (const GeometryContext &gctx, const Vector2 &lposition, const Vector3 &direction) const final
 
Result< Vector2globalToLocal (const GeometryContext &gctx, const Vector3 &position, const Vector3 &direction, double tolerance=s_onSurfaceTolerance) const final
 
SurfaceMultiIntersection intersect (const GeometryContext &gctx, const Vector3 &position, const Vector3 &direction, const BoundaryCheck &bcheck=false, ActsScalar tolerance=s_onSurfaceTolerance) const final
 
double pathCorrection (const GeometryContext &gctx, const Vector3 &position, const Vector3 &direction) const override
 
const SurfaceBoundsbounds () const final
 This method returns the bounds of the surface by reference.
 
std::string name () const override
 Return properly formatted class name for screen output.
 
AlignmentToPathMatrix alignmentToPathDerivative (const GeometryContext &gctx, const FreeVector &parameters) const final
 
ActsMatrix< 2, 3 > localCartesianToBoundLocalDerivative (const GeometryContext &gctx, const Vector3 &position) const final
 
Vector3 lineDirection (const GeometryContext &gctx) const
 
- Public Member Functions inherited from Acts::Surface
virtual ~Surface ()
 
std::shared_ptr< SurfacegetSharedPtr ()
 
std::shared_ptr< const SurfacegetSharedPtr () const
 
Surfaceoperator= (const Surface &other)
 
virtual bool operator== (const Surface &other) const
 
virtual bool operator!= (const Surface &sf) const
 
virtual SurfaceType type () const =0
 Return method for the Surface type to avoid dynamic casts.
 
virtual const Transform3transform (const GeometryContext &gctx) const
 
virtual Vector3 center (const GeometryContext &gctx) const
 
virtual Vector3 normal (const GeometryContext &gctx, const Vector3 &position) const
 
virtual Vector3 normal (const GeometryContext &gctx) const
 
const DetectorElementBaseassociatedDetectorElement () const
 
const LayerassociatedLayer () const
 
void associateLayer (const Layer &lay)
 
const ISurfaceMaterialsurfaceMaterial () const
 
const std::shared_ptr< const
ISurfaceMaterial > & 
surfaceMaterialSharedPtr () const
 
void assignDetectorElement (const DetectorElementBase &detelement)
 
void assignSurfaceMaterial (std::shared_ptr< const ISurfaceMaterial > material)
 
bool isOnSurface (const GeometryContext &gctx, const Vector3 &position, const Vector3 &direction, const BoundaryCheck &bcheck=true) const
 
virtual bool insideBounds (const Vector2 &lposition, const BoundaryCheck &bcheck=true) const
 
virtual FreeToBoundMatrix freeToBoundJacobian (const GeometryContext &gctx, const FreeVector &parameters) const
 
virtual std::ostream & toStream (const GeometryContext &gctx, std::ostream &sl) const
 
std::string toString (const GeometryContext &gctx) const
 
virtual Polyhedron polyhedronRepresentation (const GeometryContext &gctx, size_t lseg) const =0
 
AlignmentToBoundMatrix alignmentToBoundDerivative (const GeometryContext &gctx, const FreeVector &parameters, const FreeVector &pathDerivative) const
 
- Public Member Functions inherited from Acts::GeometryObject
 GeometryObject ()=default
 Defaulted constructor.
 
 GeometryObject (const GeometryObject &)=default
 Defaulted copy constructor.
 
 GeometryObject (const GeometryIdentifier &geometryId)
 
GeometryObjectoperator= (const GeometryObject &geometryId)
 
const GeometryIdentifiergeometryId () const
 
virtual double binningPositionValue (const GeometryContext &gctx, BinningValue bValue) const
 
void assignGeometryId (const GeometryIdentifier &geometryId)
 

Protected Member Functions

 LineSurface (const Transform3 &transform, double radius, double halez)
 
 LineSurface (const Transform3 &transform, std::shared_ptr< const LineBounds > lbounds=nullptr)
 
 LineSurface (std::shared_ptr< const LineBounds > lbounds, const DetectorElementBase &detelement)
 
 LineSurface (const LineSurface &other)
 
 LineSurface (const GeometryContext &gctx, const LineSurface &other, const Transform3 &shift)
 
- Protected Member Functions inherited from Acts::Surface
 Surface (const Transform3 &transform=Transform3::Identity())
 
 Surface (const Surface &other)
 
 Surface (const DetectorElementBase &detelement)
 
 Surface (const GeometryContext &gctx, const Surface &other, const Transform3 &shift)
 

Protected Attributes

std::shared_ptr< const LineBoundsm_bounds
 bounds (shared)
 
- Protected Attributes inherited from Acts::Surface
Transform3 m_transform = Transform3::Identity()
 
const DetectorElementBasem_associatedDetElement {nullptr}
 Pointer to the a DetectorElementBase.
 
const Layerm_associatedLayer {nullptr}
 
const TrackingVolumem_associatedTrackingVolume {nullptr}
 
std::shared_ptr< const
ISurfaceMaterial
m_surfaceMaterial
 Possibility to attach a material description.
 
- Protected Attributes inherited from Acts::GeometryObject
GeometryIdentifier m_geometryId
 

Private Member Functions

bool globalToLocalPlain (const GeometryContext &gctx, const Vector3 &position, const Vector3 &direction, Vector2 &lposition) const
 

Private Attributes

friend Surface
 

Additional Inherited Members

- Public Types inherited from Acts::Surface
enum  SurfaceType {
  Cone = 0, Cylinder = 1, Disc = 2, Perigee = 3,
  Plane = 4, Straw = 5, Curvilinear = 6, Other = 7
}
 
- Static Public Member Functions inherited from Acts::Surface
template<class T , typename... Args>
static std::shared_ptr< TmakeShared (Args &&...args)
 
- Static Public Attributes inherited from Acts::Surface
static std::array< std::string,
SurfaceType::Other > 
s_surfaceTypeNames
 Helper strings for screen output.
 

Detailed Description

Base class for a linear surfaces in the TrackingGeometry to describe dirft tube, straw like detectors or the Perigee It inherits from Surface.

Note
It leaves the type() method virtual, so it can not be instantiated
LineSurface.png

Definition at line 41 of file LineSurface.hpp.

View newest version in sPHENIX GitHub at line 41 of file LineSurface.hpp

Constructor & Destructor Documentation

Acts::LineSurface::LineSurface ( const Transform3 transform,
double  radius,
double  halez 
)
protected

Constructor from Transform3 and bounds

Parameters
transformThe transform that positions the surface in the global frame
radiusThe straw radius
halezThe half length in z

Definition at line 33 of file LineSurface.cpp.

View newest version in sPHENIX GitHub at line 33 of file LineSurface.cpp

Acts::LineSurface::LineSurface ( const Transform3 transform,
std::shared_ptr< const LineBounds lbounds = nullptr 
)
protected

Constructor from Transform3 and a shared bounds object

Parameters
transformThe transform that positions the surface in the global frame
lboundsThe bounds describing the straw dimensions, can be optionally nullptr

Definition at line 39 of file LineSurface.cpp.

View newest version in sPHENIX GitHub at line 39 of file LineSurface.cpp

Acts::LineSurface::LineSurface ( std::shared_ptr< const LineBounds lbounds,
const DetectorElementBase detelement 
)
protected

Constructor from DetectorElementBase : Element proxy

Parameters
lboundsThe bounds describing the straw dimensions
detelementfor which this surface is (at least) one representation

Definition at line 43 of file LineSurface.cpp.

View newest version in sPHENIX GitHub at line 43 of file LineSurface.cpp

References m_bounds, and throw_assert.

Acts::LineSurface::LineSurface ( const LineSurface other)
protected

Copy constructor

Parameters
otherThe source surface for copying

Definition at line 49 of file LineSurface.cpp.

View newest version in sPHENIX GitHub at line 49 of file LineSurface.cpp

Acts::LineSurface::LineSurface ( const GeometryContext gctx,
const LineSurface other,
const Transform3 shift 
)
protected

Copy constructor - with shift

Parameters
gctxThe current geometry context object, e.g. alignment
otheris the source cone surface
shiftis the additional transform applied after copying

Definition at line 52 of file LineSurface.cpp.

View newest version in sPHENIX GitHub at line 52 of file LineSurface.cpp

Acts::LineSurface::~LineSurface ( )
overridedefault
Acts::LineSurface::LineSurface ( )
delete

Member Function Documentation

Acts::AlignmentToPathMatrix Acts::LineSurface::alignmentToPathDerivative ( const GeometryContext gctx,
const FreeVector parameters 
) const
finalvirtual

Calculate the derivative of path length at the geometry constraint or point-of-closest-approach w.r.t. alignment parameters of the surface (i.e. local frame origin in global 3D Cartesian coordinates and its rotation represented with extrinsic Euler angles)

Parameters
gctxThe current geometry context object, e.g. alignment
parametersis the free parameters
Returns
Derivative of path length w.r.t. the alignment parameters

Reimplemented from Acts::Surface.

Definition at line 280 of file LineSurface.cpp.

View newest version in sPHENIX GitHub at line 280 of file LineSurface.cpp

References dz, Acts::eAlignmentCenter0, Acts::eAlignmentRotation0, Acts::eFreeDir0, Acts::eFreePos0, norm, position, Acts::detail::rotationToLocalAxesDerivative(), and Acts::Test::transform.

Referenced by Acts::Test::BOOST_AUTO_TEST_CASE().

+ Here is the call graph for this function:

+ Here is the caller graph for this function:

Acts::Vector3 Acts::LineSurface::binningPosition ( const GeometryContext gctx,
BinningValue  bValue 
) const
finalvirtual

The binning position is the position calculated for a certain binning type

Parameters
gctxThe current geometry context object, e.g. alignment
bValueis the binning type to be used
Returns
position that can beused for this binning

Implements Acts::GeometryObject.

Definition at line 130 of file LineSurface.cpp.

View newest version in sPHENIX GitHub at line 130 of file LineSurface.cpp

Referenced by Acts::Test::BOOST_AUTO_TEST_CASE().

+ Here is the caller graph for this function:

const Acts::SurfaceBounds & Acts::LineSurface::bounds ( ) const
finalvirtual

This method returns the bounds of the surface by reference.

Implements Acts::Surface.

Definition at line 141 of file LineSurface.cpp.

View newest version in sPHENIX GitHub at line 141 of file LineSurface.cpp

References Acts::s_noBounds.

Referenced by Acts::Test::BOOST_AUTO_TEST_CASE().

+ Here is the caller graph for this function:

Acts::BoundToFreeMatrix Acts::LineSurface::boundToFreeJacobian ( const GeometryContext gctx,
const BoundVector boundParams 
) const
finalvirtual

Calculate the jacobian from local to global which the surface knows best, hence the calculation is done here.

Parameters
gctxThe current geometry context object, e.g. alignment
boundParamsis the bound parameters vector
Returns
Jacobian from local to global

Reimplemented from Acts::Surface.

Definition at line 197 of file LineSurface.cpp.

View newest version in sPHENIX GitHub at line 197 of file LineSurface.cpp

References Acts::VectorHelpers::cross(), Acts::eBoundLoc0, Acts::eBoundPhi, Acts::eBoundQOverP, Acts::eBoundTheta, Acts::eBoundTime, Acts::eFreeDir0, Acts::eFreeDir1, Acts::eFreeDir2, Acts::eFreePos0, Acts::eFreeQOverP, Acts::eFreeTime, position, and Acts::detail::transformBoundToFreeParameters().

+ Here is the call graph for this function:

Acts::FreeToPathMatrix Acts::LineSurface::freeToPathDerivative ( const GeometryContext gctx,
const FreeVector parameters 
) const
finalvirtual

Calculate the derivative of path length at the geometry constraint or point-of-closest-approach w.r.t. free parameters

Parameters
gctxThe current geometry context object, e.g. alignment
parametersis the free parameters
Returns
Derivative of path length w.r.t. free parameters

Reimplemented from Acts::Surface.

Definition at line 250 of file LineSurface.cpp.

View newest version in sPHENIX GitHub at line 250 of file LineSurface.cpp

References dz, Acts::eFreeDir0, Acts::eFreePos0, norm, and position.

Acts::Result< Acts::Vector2 > Acts::LineSurface::globalToLocal ( const GeometryContext gctx,
const Vector3 position,
const Vector3 direction,
double  tolerance = s_onSurfaceTolerance 
) const
finalvirtual

Specified for LineSurface: global to local method without dynamic memory allocation.

This method is the true global -> local transformation. It makes use of globalToLocal and indicates the sign of the Acts::eBoundLoc0 by the given momentum direction.

The calculation of the sign of the radius (or $ d_0 $) can be done as follows: May $ \vec d = \vec m - \vec c $ denote the difference between the center of the line and the global position of the measurement/predicted state. Then, $ \vec d $ lies in the so-called measurement plane. The latter is determined by the two orthogonal vectors $ \vec{\texttt{measY}} = \vec{e}_z $ and $ \vec{\texttt{measX}} = \vec{\texttt{measY}} \times \frac{\vec{p}}{|\vec{p}|} $.

The sign of the radius (or $ d_{0} $ ) is then defined by the projection of $ \vec{d} $ on $ \vec{measX} $:
$ sign = -sign(\vec{d} \cdot \vec{measX}) $

SignOfDriftCircleD0.gif
Parameters
gctxThe current geometry context object, e.g. alignment
positionglobal 3D position - considered to be on surface but not inside bounds (check is done)
directionglobal 3D momentum direction (optionally ignored)
tolerance(unused)
Returns
A Result<Vector2>, which is set to !ok() if the position is not the point of closest approach to the line surface.

Implements Acts::Surface.

Definition at line 79 of file LineSurface.cpp.

View newest version in sPHENIX GitHub at line 79 of file LineSurface.cpp

References Acts::VectorHelpers::perp(), and Acts::Test::transform.

Referenced by Acts::Test::BOOST_AUTO_TEST_CASE().

+ Here is the call graph for this function:

+ Here is the caller graph for this function:

bool Acts::LineSurface::globalToLocalPlain ( const GeometryContext gctx,
const Vector3 position,
const Vector3 direction,
Vector2 lposition 
) const
private

helper function to apply the globalToLocal with out transform

Parameters
gctxThe current geometry context object, e.g. alignment
positionis the global position
directionis the momentum direction
lpositionis the local position to be filled
Acts::SurfaceMultiIntersection Acts::LineSurface::intersect ( const GeometryContext gctx,
const Vector3 position,
const Vector3 direction,
const BoundaryCheck bcheck = false,
ActsScalar  tolerance = s_onSurfaceTolerance 
) const
finalvirtual

Calculate the straight-line intersection with the line surface.

Mathematical motivation:

Given two lines in parameteric form:

$ \vec l_{a}(u) = \vec m_a + u \cdot \vec e_{a} $

$ \vec l_{b}(\mu) = \vec m_b + \mu \cdot \vec e_{b} $

The vector between any two points on the two lines is given by:

$ \vec s(u, \mu) = \vec l_{b} - l_{a} = \vec m_{ab} + \mu \cdot \vec e_{b} - u \cdot \vec e_{a} $,

where $ \vec m_{ab} = \vec m_{b} - \vec m_{a} $.

$ \vec s(u_0, \mu_0) $ denotes the vector between the two closest points

$ \vec l_{a,0} = l_{a}(u_0) $ and $ \vec l_{b,0} = l_{b}(\mu_0) $

and is perpendicular to both, $ \vec e_{a} $ and $ \vec e_{b} $.

This results in a system of two linear equations:

  • (i) $ 0 = \vec s(u_0, \mu_0) \cdot \vec e_a = \vec m_{ab} \cdot \vec e_a + \mu_0 \vec e_a \cdot \vec e_b - u_0 $
  • (ii) $ 0 = \vec s(u_0, \mu_0) \cdot \vec e_b = \vec m_{ab} \cdot \vec e_b + \mu_0 - u_0 \vec e_b \cdot \vec e_a $

Solving (i) and (ii) for $ u $ and $ \mu_0 $ yields:

  • $ u_0 = \frac{(\vec m_{ab} \cdot \vec e_a)-(\vec m_{ab} \cdot \vec e_b)(\vec e_a \cdot \vec e_b)}{1-(\vec e_a \cdot \vec e_b)^2} $
  • $ \mu_0 = - \frac{(\vec m_{ab} \cdot \vec e_b)-(\vec m_{ab} \cdot \vec e_a)(\vec e_a \cdot \vec e_b)}{1-(\vec e_a \cdot \vec e_b)^2} $

The function checks if $ u_0 \simeq 0$ to check if the current position is at the point of closest approach, i.e. the intersection point, in which case it will return an onSurace intersection result. Otherwise, the path length from position to the point of closest approach ( $ u_0 $) is returned in a reachable intersection.

Parameters
gctxThe current geometry context object, e.g. alignment
positionThe global position as a starting point
directionThe global direction at the starting point
Note
expected to be normalized
Parameters
bcheckThe boundary check directive for the estimate
tolerancethe tolerance used for the intersection
Returns
is the intersection object

Implements Acts::Surface.

Definition at line 148 of file LineSurface.cpp.

View newest version in sPHENIX GitHub at line 148 of file LineSurface.cpp

References Acts::LineBounds::eHalfLengthZ, Acts::LineBounds::eR, Acts::Intersection< 3 >::invalid(), norm, position, utils::status, Acts::Test::tolerance, Acts::Test::transform, and physmon_ckf_tracking::u.

Referenced by Acts::Test::BOOST_AUTO_TEST_CASE().

+ Here is the call graph for this function:

+ Here is the caller graph for this function:

Acts::Vector3 Acts::LineSurface::lineDirection ( const GeometryContext gctx) const

Definition at line 323 of file LineSurface.cpp.

View newest version in sPHENIX GitHub at line 323 of file LineSurface.cpp

References Acts::Test::transform.

Acts::ActsMatrix< 2, 3 > Acts::LineSurface::localCartesianToBoundLocalDerivative ( const GeometryContext gctx,
const Vector3 position 
) const
finalvirtual

Calculate the derivative of bound track parameters local position w.r.t. position in local 3D Cartesian coordinates

Parameters
gctxThe current geometry context object, e.g. alignment
positionThe position of the parameters in global
Returns
Derivative of bound local position w.r.t. position in local 3D cartesian coordinates

Implements Acts::Surface.

Definition at line 311 of file LineSurface.cpp.

View newest version in sPHENIX GitHub at line 311 of file LineSurface.cpp

References Acts::VectorHelpers::phi(), position, and Acts::Test::transform.

Referenced by Acts::Test::BOOST_AUTO_TEST_CASE().

+ Here is the call graph for this function:

+ Here is the caller graph for this function:

Acts::Vector3 Acts::LineSurface::localToGlobal ( const GeometryContext gctx,
const Vector2 lposition,
const Vector3 direction 
) const
finalvirtual

Local to global transformation

Note
for line surfaces the momentum direction is used in order to interpret the drift radius
Parameters
gctxThe current geometry context object, e.g. alignment
lpositionis the local position to be transformed
directionis the global momentum direction (used to sign the closest approach)
Returns
global position by value

Implements Acts::Surface.

Definition at line 65 of file LineSurface.cpp.

View newest version in sPHENIX GitHub at line 65 of file LineSurface.cpp

References Acts::eBoundLoc0, Acts::eBoundLoc1, and Acts::Test::transform.

Referenced by Acts::Test::BOOST_AUTO_TEST_CASE().

+ Here is the caller graph for this function:

std::string Acts::LineSurface::name ( ) const
overridevirtual

Return properly formatted class name for screen output.

Implements Acts::Surface.

Reimplemented in Acts::StrawSurface, and Acts::PerigeeSurface.

Definition at line 105 of file LineSurface.cpp.

View newest version in sPHENIX GitHub at line 105 of file LineSurface.cpp

Referenced by Acts::Test::BOOST_AUTO_TEST_CASE().

+ Here is the caller graph for this function:

Acts::Vector3 Acts::LineSurface::normal ( const GeometryContext gctx,
const Vector2 lposition 
) const
finalvirtual

The normal vector is undefined if we do not know the momentum.

Parameters
gctxThe current geometry context object, e.g. alignment
lpositionis the local position is ignored
Returns
a zero vector

Implements Acts::Surface.

Definition at line 135 of file LineSurface.cpp.

View newest version in sPHENIX GitHub at line 135 of file LineSurface.cpp

Referenced by Acts::Test::BOOST_AUTO_TEST_CASE().

+ Here is the caller graph for this function:

Acts::LineSurface & Acts::LineSurface::operator= ( const LineSurface other)

Assignment operator

Parameters
otheris the source surface dor copying

Definition at line 57 of file LineSurface.cpp.

View newest version in sPHENIX GitHub at line 57 of file LineSurface.cpp

References m_bounds, and Acts::Surface::operator=().

Referenced by Acts::PerigeeSurface::operator=(), and Acts::StrawSurface::operator=().

+ Here is the call graph for this function:

+ Here is the caller graph for this function:

double Acts::LineSurface::pathCorrection ( const GeometryContext gctx,
const Vector3 position,
const Vector3 direction 
) const
overridevirtual

the pathCorrection for derived classes with thickness is by definition 1 for LineSurfaces

Note
input parameters are ignored
there's no material associated to the line surface

Implements Acts::Surface.

Definition at line 124 of file LineSurface.cpp.

View newest version in sPHENIX GitHub at line 124 of file LineSurface.cpp

Referenced by Acts::Test::BOOST_AUTO_TEST_CASE().

+ Here is the caller graph for this function:

Acts::RotationMatrix3 Acts::LineSurface::referenceFrame ( const GeometryContext gctx,
const Vector3 position,
const Vector3 direction 
) const
finalvirtual

Return the measurement frame - this is needed for alignment, in particular

for StraightLine and Perigee Surface

  • the default implementation is the RotationMatrix3 of the transform
Parameters
gctxThe current geometry context object, e.g. alignment
positionis the global position where the measurement frame is constructed
directionis the momentum direction used for the measurement frame construction
Returns
is a rotation matrix that indicates the measurement frame

Reimplemented from Acts::Surface.

Definition at line 109 of file LineSurface.cpp.

View newest version in sPHENIX GitHub at line 109 of file LineSurface.cpp

Referenced by Acts::Test::BOOST_AUTO_TEST_CASE().

+ Here is the caller graph for this function:

Member Data Documentation

std::shared_ptr<const LineBounds> Acts::LineSurface::m_bounds
protected

bounds (shared)

Definition at line 298 of file LineSurface.hpp.

View newest version in sPHENIX GitHub at line 298 of file LineSurface.hpp

Referenced by LineSurface(), Acts::StrawSurface::operator=(), and operator=().

friend Acts::LineSurface::Surface
private

Definition at line 43 of file LineSurface.hpp.

View newest version in sPHENIX GitHub at line 43 of file LineSurface.hpp


The documentation for this class was generated from the following files: