ATLAS Offline Software
Loading...
Searching...
No Matches
PlaneSurface Class Reference

Class for a planaer rectangular or trapezoidal surface in the ATLAS detector. More...

#include <PlaneSurface.h>

Inheritance diagram for PlaneSurface:
Collaboration diagram for PlaneSurface:

Public Member Functions

 PlaneSurface ()
 Default Constructor - needed for persistency.
 PlaneSurface (const PlaneSurface &psf)=default
 Copy Constructor.
 PlaneSurface (PlaneSurface &&psf) noexcept=default
 Move Constructor.
 PlaneSurface (const PlaneSurface &psf, const Amg::Transform3D &transf)
 Copy Constructor with shift.
 PlaneSurface (const Amg::Vector3D &position, const CurvilinearUVT &curvUVT)
 Dedicated Constructor with CurvilinearUVT class.
 PlaneSurface (const TrkDetElementBase &detelement, const Amg::Transform3D &transf)
 Constructor from TrkDetElementBase.
 PlaneSurface (const TrkDetElementBase &detelement)
 Constructor from TrkDetElementBase.
 PlaneSurface (const TrkDetElementBase &detelement, const Identifier &id, const Amg::Transform3D &transf)
 Constructor from TrkDetElementBase and Identifier in case one element holds more surfaces.
 PlaneSurface (const TrkDetElementBase &detelement, const Identifier &id)
 Constructor from TrkDetElementBase and Identifier in case one element holds more surfaces.
 PlaneSurface (const Amg::Transform3D &htrans)
 Constructor for planar Surface without Bounds , reference.
 PlaneSurface (const Amg::Transform3D &htrans, double halephi, double haleta)
 Constructor for Rectangular Planes.
 PlaneSurface (const Amg::Transform3D &htrans, double minhalephi, double maxhalephi, double haleta)
 Constructor for Trapezoidal Planes.
 PlaneSurface (const Amg::Transform3D &htrans, std::shared_ptr< const Trk::SurfaceBounds > sbounds)
 Constructor for Planes with shared object.
PlaneSurfaceoperator= (const PlaneSurface &psf)=default
 Assignment operator.
PlaneSurfaceoperator= (PlaneSurface &&psf) noexcept=default
 Move assignment operator.
virtual ~PlaneSurface ()
 Destructor.
virtual bool operator== (const Surface &sf) const override
 Equality operator.
bool operator== (const PlaneSurface &cf) const
virtual PlaneSurfaceclone () const override
 Virtual constructor.
virtual constexpr SurfaceType type () const override final
 Return the surface type.
virtual Surface::ChargedTrackParametersUniquePtr createUniqueTrackParameters (double l1, double l2, double phi, double theta, double qop, std::optional< AmgSymMatrix(5)> cov=std::nullopt) const override final
 Use the Surface as a ParametersBase constructor, from local parameters - charged.
virtual Surface::ChargedTrackParametersUniquePtr createUniqueTrackParameters (const Amg::Vector3D &position, const Amg::Vector3D &momentum, double charge, std::optional< AmgSymMatrix(5)> cov=std::nullopt) const override final
 Use the Surface as a ParametersBase constructor, from global parameters - charged.
virtual NeutralTrackParametersUniquePtr createUniqueNeutralParameters (double l1, double l2, double phi, double theta, double oop, std::optional< AmgSymMatrix(5)> cov=std::nullopt) const override final
 Use the Surface as a ParametersBase constructor, from local parameters - neutral.
virtual NeutralTrackParametersUniquePtr createUniqueNeutralParameters (const Amg::Vector3D &position, const Amg::Vector3D &momentum, double charge=0., std::optional< AmgSymMatrix(5)> cov=std::nullopt) const override final
 Use the Surface as a ParametersBase constructor, from global parameters.
template<int DIM, class T>
std::unique_ptr< ParametersT< DIM, T, PlaneSurface > > createUniqueParameters (double l1, double l2, double phi, double theta, double qop, std::optional< AmgSymMatrix(DIM)> cov=std::nullopt) const
 Use the Surface as a ParametersBase constructor, from local parameters.
template<int DIM, class T>
std::unique_ptr< ParametersT< DIM, T, PlaneSurface > > createUniqueParameters (const Amg::Vector3D &position, const Amg::Vector3D &momentum, double charge, std::optional< AmgSymMatrix(DIM)> cov=std::nullopt) const
 Use the Surface as a ParametersBase constructor, from global parameters.
virtual const SurfaceBounds & bounds () const override final
 This method returns the bounds by reference, static NoBounds in case of no boundaries.
virtual bool insideBounds (const Amg::Vector2D &locpos, double tol1=0., double tol2=0.) const override
 This method calls the inside() method of the Bounds.
virtual bool insideBoundsCheck (const Amg::Vector2D &locpos, const BoundaryCheck &bchk) const override final
virtual bool isOnSurface (const Amg::Vector3D &glopo, const BoundaryCheck &bchk=true, double tol1=0., double tol2=0.) const override final
 This method returns true if the GlobalPosition is on the Surface for both, within or without check of whether the local position is inside boundaries or not.
virtual void localToGlobal (const Amg::Vector2D &locp, const Amg::Vector3D &mom, Amg::Vector3D &glob) const override final
 Specified for PlaneSurface: LocalToGlobal method without dynamic memory allocation.
virtual bool globalToLocal (const Amg::Vector3D &glob, const Amg::Vector3D &mom, Amg::Vector2D &loc) const override final
 Specified for PlaneSurface: GlobalToLocal method without dynamic memory allocation - boolean checks if on surface.
void localToGlobalDirection (const Trk::LocalDirection &locdir, Amg::Vector3D &globdir) const
 This method transforms a local direction wrt the plane to a global direction.
void globalToLocalDirection (const Amg::Vector3D &glodir, Trk::LocalDirection &locdir) const
 This method transforms the global direction to a local direction wrt the plane.
virtual Intersection straightLineIntersection (const Amg::Vector3D &pos, const Amg::Vector3D &dir, bool forceDir, Trk::BoundaryCheck bchk) const override final
 fast straight line intersection schema - standard: provides closest intersection and (signed) path length forceDir is to provide the closest forward solution
virtual DistanceSolution straightLineDistanceEstimate (const Amg::Vector3D &pos, const Amg::Vector3D &dir) const override final
 fast straight line distance evaluation to Surface
virtual DistanceSolution straightLineDistanceEstimate (const Amg::Vector3D &pos, const Amg::Vector3D &dir, bool Bound) const override final
 fast straight line distance evaluation to Surface - with bound option
virtual std::string name () const override
 Return properly formatted class name for screen output.

Static Public Attributes

static constexpr SurfaceType staticType = SurfaceType::Plane
 The surface type static constexpr.

Protected Attributes

std::shared_ptr< const SurfaceBounds > m_bounds
 bounds (shared)

Static Protected Attributes

static const NoBounds s_boundless

Friends

template<class SURFACE, class BOUNDS_CNV>
class ::BoundSurfaceCnv_p1
 < data members
template<class SURFACE, class BOUNDS_CNV>
class ::BoundSurfaceCnv_p2

Detailed Description

Class for a planaer rectangular or trapezoidal surface in the ATLAS detector.

It inherits from Surface.

The Trk::PlaneSurface extends the Surface class with the possibility to convert in addition to local to global positions, also local to global direction (vice versa). The definition with of a local direciton with respect to a plane can be found in the dedicated Trk::LocalDirection class of the TrkEventPrimitives package.

Author
Andre.nosp@m.as.S.nosp@m.alzbu.nosp@m.rger.nosp@m.@cern.nosp@m..ch
Christos Anastopoulos (Thread safety and interface cleanup)
Shaun Roe (interface cleanup)

Definition at line 63 of file PlaneSurface.h.

Constructor & Destructor Documentation

◆ PlaneSurface() [1/13]

Trk::PlaneSurface::PlaneSurface ( )

Default Constructor - needed for persistency.

Definition at line 36 of file PlaneSurface.cxx.

37 : Trk::Surface()
38 , m_bounds(nullptr)
39{}
std::shared_ptr< const SurfaceBounds > m_bounds
bounds (shared)

◆ PlaneSurface() [2/13]

Trk::PlaneSurface::PlaneSurface ( const PlaneSurface & psf)
default

Copy Constructor.

◆ PlaneSurface() [3/13]

Trk::PlaneSurface::PlaneSurface ( PlaneSurface && psf)
defaultnoexcept

Move Constructor.

◆ PlaneSurface() [4/13]

Trk::PlaneSurface::PlaneSurface ( const PlaneSurface & psf,
const Amg::Transform3D & transf )

Copy Constructor with shift.

Definition at line 43 of file PlaneSurface.cxx.

44 : Trk::Surface(psf, transf)
45 , m_bounds(psf.m_bounds)
46{}

◆ PlaneSurface() [5/13]

ATH_FLATTEN Trk::PlaneSurface::PlaneSurface ( const Amg::Vector3D & position,
const CurvilinearUVT & curvUVT )

Dedicated Constructor with CurvilinearUVT class.

Definition at line 55 of file PlaneSurface.cxx.

56 : Trk::Surface()
57 , m_bounds(nullptr) // curvilinear surfaces are boundless
58{
59 Amg::Translation3D curvilinearTranslation(position.x(), position.y(), position.z());
60 // create the rotation
61 Amg::RotationMatrix3D curvilinearRotation;
62 curvilinearRotation.col(0) = curvUVT.curvU();
63 curvilinearRotation.col(1) = curvUVT.curvV();
64 curvilinearRotation.col(2) = curvUVT.curvT();
66 transform = curvilinearRotation;
67 transform.pretranslate(position);
68 Trk::Surface::m_transforms = std::make_unique<Transforms>(transform);
69}
const Amg::Transform3D & transform() const
Returns HepGeom::Transform3D by reference.
std::unique_ptr< Transforms > m_transforms
Unique Pointer to the Transforms struct.
Eigen::Matrix< double, 3, 3 > RotationMatrix3D
Eigen::Affine3d Transform3D
Eigen::Translation< double, 3 > Translation3D

◆ PlaneSurface() [6/13]

Trk::PlaneSurface::PlaneSurface ( const TrkDetElementBase & detelement,
const Amg::Transform3D & transf )

Constructor from TrkDetElementBase.

Definition at line 72 of file PlaneSurface.cxx.

73 : Trk::Surface(detelement)
74 , m_bounds(nullptr)
75{
76 Trk::Surface::m_transforms = std::make_unique<Transforms>(transf);
77}

◆ PlaneSurface() [7/13]

Trk::PlaneSurface::PlaneSurface ( const TrkDetElementBase & detelement)

Constructor from TrkDetElementBase.

Definition at line 80 of file PlaneSurface.cxx.

81 : Trk::Surface(detelement)
82 , m_bounds(nullptr)
83{
84 //
85}

◆ PlaneSurface() [8/13]

Trk::PlaneSurface::PlaneSurface ( const TrkDetElementBase & detelement,
const Identifier & id,
const Amg::Transform3D & transf )

Constructor from TrkDetElementBase and Identifier in case one element holds more surfaces.

Definition at line 88 of file PlaneSurface.cxx.

91 : Trk::Surface(detelement, id)
92 , m_bounds(nullptr)
93{
94 Trk::Surface::m_transforms = std::make_unique<Transforms>(transf);
95}

◆ PlaneSurface() [9/13]

Trk::PlaneSurface::PlaneSurface ( const TrkDetElementBase & detelement,
const Identifier & id )

Constructor from TrkDetElementBase and Identifier in case one element holds more surfaces.

Definition at line 98 of file PlaneSurface.cxx.

100 : Trk::Surface(detelement, id)
101 , m_bounds(nullptr)
102{
103 //
104}

◆ PlaneSurface() [10/13]

Trk::PlaneSurface::PlaneSurface ( const Amg::Transform3D & htrans)

Constructor for planar Surface without Bounds , reference.

Definition at line 107 of file PlaneSurface.cxx.

108 : Trk::Surface(htrans)
109 , m_bounds(nullptr)
110{}

◆ PlaneSurface() [11/13]

Trk::PlaneSurface::PlaneSurface ( const Amg::Transform3D & htrans,
double halephi,
double haleta )

Constructor for Rectangular Planes.

Definition at line 113 of file PlaneSurface.cxx.

114 : Trk::Surface(htrans)
115 , m_bounds(std::make_shared<const Trk::RectangleBounds>(halephi, haleta))
116{}

◆ PlaneSurface() [12/13]

Trk::PlaneSurface::PlaneSurface ( const Amg::Transform3D & htrans,
double minhalephi,
double maxhalephi,
double haleta )

Constructor for Trapezoidal Planes.

Definition at line 119 of file PlaneSurface.cxx.

120 : Trk::Surface(htrans)
121 , m_bounds(std::make_shared<const Trk::TrapezoidBounds>(minhalephi, maxhalephi, haleta))
122{}

◆ PlaneSurface() [13/13]

Trk::PlaneSurface::PlaneSurface ( const Amg::Transform3D & htrans,
std::shared_ptr< const Trk::SurfaceBounds > sbounds )

Constructor for Planes with shared object.

Definition at line 125 of file PlaneSurface.cxx.

128 : Trk::Surface(htrans), m_bounds(std::move(tbounds)) {}

◆ ~PlaneSurface()

Trk::PlaneSurface::~PlaneSurface ( )
virtualdefault

Destructor.

Member Function Documentation

◆ bounds()

virtual const SurfaceBounds & Trk::PlaneSurface::bounds ( ) const
finaloverridevirtual

This method returns the bounds by reference, static NoBounds in case of no boundaries.

Implements Trk::Surface.

◆ clone()

virtual PlaneSurface * Trk::PlaneSurface::clone ( ) const
overridevirtual

Virtual constructor.

Implements Trk::Surface.

◆ createUniqueNeutralParameters() [1/2]

Trk::Surface::NeutralTrackParametersUniquePtr Trk::PlaneSurface::createUniqueNeutralParameters ( const Amg::Vector3D & position,
const Amg::Vector3D & momentum,
double charge = 0.,
std::optional< AmgSymMatrix(5)> cov = std::nullopt ) const
finaloverridevirtual

Use the Surface as a ParametersBase constructor, from global parameters.

Use the Surface as a ParametersBase constructor, from global parameters - neutral.

  • neutral

Implements Trk::Surface.

Definition at line 188 of file PlaneSurface.cxx.

193{
194 return std::make_unique<ParametersT<5, Neutral, PlaneSurface>>(
195 position, momentum, charge, *this, std::move(cov));
196}
double charge(const T &p)
Definition AtlasPID.h:997

◆ createUniqueNeutralParameters() [2/2]

Trk::Surface::NeutralTrackParametersUniquePtr Trk::PlaneSurface::createUniqueNeutralParameters ( double l1,
double l2,
double phi,
double theta,
double oop,
std::optional< AmgSymMatrix(5)> cov = std::nullopt ) const
finaloverridevirtual

Use the Surface as a ParametersBase constructor, from local parameters - neutral.

Implements Trk::Surface.

Definition at line 173 of file PlaneSurface.cxx.

180{
181 return std::make_unique<ParametersT<5, Neutral, PlaneSurface>>(
182 l1, l2, phi, theta, oop, *this, std::move(cov));
183}
Scalar phi() const
phi method
Scalar theta() const
theta method

◆ createUniqueParameters() [1/2]

template<int DIM, class T>
std::unique_ptr< ParametersT< DIM, T, PlaneSurface > > Trk::PlaneSurface::createUniqueParameters ( const Amg::Vector3D & position,
const Amg::Vector3D & momentum,
double charge,
std::optional< AmgSymMatrix(DIM)> cov = std::nullopt ) const

Use the Surface as a ParametersBase constructor, from global parameters.

◆ createUniqueParameters() [2/2]

template<int DIM, class T>
std::unique_ptr< ParametersT< DIM, T, PlaneSurface > > Trk::PlaneSurface::createUniqueParameters ( double l1,
double l2,
double phi,
double theta,
double qop,
std::optional< AmgSymMatrix(DIM)> cov = std::nullopt ) const

Use the Surface as a ParametersBase constructor, from local parameters.

◆ createUniqueTrackParameters() [1/2]

Trk::Surface::ChargedTrackParametersUniquePtr Trk::PlaneSurface::createUniqueTrackParameters ( const Amg::Vector3D & position,
const Amg::Vector3D & momentum,
double charge,
std::optional< AmgSymMatrix(5)> cov = std::nullopt ) const
finaloverridevirtual

Use the Surface as a ParametersBase constructor, from global parameters - charged.

Implements Trk::Surface.

Definition at line 160 of file PlaneSurface.cxx.

165{
166 return std::make_unique<ParametersT<5, Charged, PlaneSurface>>(
167 position, momentum, charge, *this, std::move(cov));
168}

◆ createUniqueTrackParameters() [2/2]

Trk::Surface::ChargedTrackParametersUniquePtr Trk::PlaneSurface::createUniqueTrackParameters ( double l1,
double l2,
double phi,
double theta,
double qop,
std::optional< AmgSymMatrix(5)> cov = std::nullopt ) const
finaloverridevirtual

Use the Surface as a ParametersBase constructor, from local parameters - charged.

Implements Trk::Surface.

Definition at line 146 of file PlaneSurface.cxx.

153{
154 return std::make_unique<ParametersT<5, Charged, PlaneSurface>>(
155 l1, l2, phi, theta, qop, *this, std::move(cov));
156}

◆ globalToLocal()

ATH_FLATTEN bool Trk::PlaneSurface::globalToLocal ( const Amg::Vector3D & glob,
const Amg::Vector3D & mom,
Amg::Vector2D & loc ) const
finaloverridevirtual

Specified for PlaneSurface: GlobalToLocal method without dynamic memory allocation - boolean checks if on surface.

Implements Trk::Surface.

Definition at line 212 of file PlaneSurface.cxx.

214 {
215 Amg::Vector3D loc3Dframe = inverseTransformMultHelper(glopos);
216 locpos = Amg::Vector2D(loc3Dframe.x(), loc3Dframe.y());
217 return (loc3Dframe.z() * loc3Dframe.z() <=
219}
static constexpr double s_onSurfaceTolerance
Tolerance for being on Surface.
Amg::Vector3D inverseTransformMultHelper(const Amg::Vector3D &glopos) const
Eigen::Matrix< double, 2, 1 > Vector2D
Eigen::Matrix< double, 3, 1 > Vector3D

◆ globalToLocalDirection()

void Trk::PlaneSurface::globalToLocalDirection ( const Amg::Vector3D & glodir,
Trk::LocalDirection & locdir ) const

This method transforms the global direction to a local direction wrt the plane.

Definition at line 259 of file PlaneSurface.cxx.

260{
261 // bring the global direction into the surface frame
263 ldir = Trk::LocalDirection(std::atan2(d.z(), d.x()), std::atan2(d.z(), d.y()));
264}
Amg::Transform3D inverseTransformHelper() const
Helper method to factorize in one place common operations calculate inverse transofrm and multiply wi...

◆ insideBounds()

virtual bool Trk::PlaneSurface::insideBounds ( const Amg::Vector2D & locpos,
double tol1 = 0.,
double tol2 = 0. ) const
overridevirtual

This method calls the inside() method of the Bounds.

Implements Trk::Surface.

◆ insideBoundsCheck()

virtual bool Trk::PlaneSurface::insideBoundsCheck ( const Amg::Vector2D & locpos,
const BoundaryCheck & bchk ) const
finaloverridevirtual

◆ isOnSurface()

ATH_FLATTEN bool Trk::PlaneSurface::isOnSurface ( const Amg::Vector3D & glopo,
const BoundaryCheck & bchk = true,
double tol1 = 0.,
double tol2 = 0. ) const
finaloverridevirtual

This method returns true if the GlobalPosition is on the Surface for both, within or without check of whether the local position is inside boundaries or not.

Reimplemented from Trk::Surface.

Definition at line 268 of file PlaneSurface.cxx.

271{
272 Amg::Vector3D loc3Dframe = inverseTransformMultHelper(glopo);
273 if (std::abs(loc3Dframe(2)) > (s_onSurfaceTolerance + tol1)){
274 return false;
275 }
276 return (bchk ? bounds().inside(Amg::Vector2D(loc3Dframe(0), loc3Dframe(1)), tol1, tol2) : true);
277}
virtual const SurfaceBounds & bounds() const override final
This method returns the bounds by reference, static NoBounds in case of no boundaries.

◆ localToGlobal()

ATH_FLATTEN void Trk::PlaneSurface::localToGlobal ( const Amg::Vector2D & locp,
const Amg::Vector3D & mom,
Amg::Vector3D & glob ) const
finaloverridevirtual

Specified for PlaneSurface: LocalToGlobal method without dynamic memory allocation.

Implements Trk::Surface.

Definition at line 201 of file PlaneSurface.cxx.

204{
205 Amg::Vector3D loc3Dframe(locpos[Trk::locX], locpos[Trk::locY], 0.);
206 glopos = transform() * loc3Dframe;
207}
@ locY
local cartesian
Definition ParamDefs.h:38
@ locX
Definition ParamDefs.h:37

◆ localToGlobalDirection()

void Trk::PlaneSurface::localToGlobalDirection ( const Trk::LocalDirection & locdir,
Amg::Vector3D & globdir ) const

This method transforms a local direction wrt the plane to a global direction.

Definition at line 241 of file PlaneSurface.cxx.

242{
243
244 CxxUtils::sincos scXZ(locdir.angleXZ());
245 CxxUtils::sincos scYZ(locdir.angleYZ());
246
247 double norm = 1. / std::sqrt(scYZ.cs * scYZ.cs * scXZ.sn * scXZ.sn + scYZ.sn * scYZ.sn);
248
249 // decide on the sign
250 double sign = (scXZ.sn < 0.) ? -1. : 1.;
251
252 // now calculate the GlobalDirection in the global frame
253 globdir =
254 transform().linear() *
255 Amg::Vector3D(sign * scXZ.cs * scYZ.sn * norm, sign * scXZ.sn * scYZ.cs * norm, sign * scXZ.sn * scYZ.sn * norm);
256}
int sign(int a)
double angleXZ() const
access method for angle of local XZ projection
double angleYZ() const
access method for angle of local YZ projection

◆ name()

virtual std::string Trk::PlaneSurface::name ( ) const
overridevirtual

Return properly formatted class name for screen output.

Implements Trk::Surface.

◆ operator=() [1/2]

PlaneSurface & Trk::PlaneSurface::operator= ( const PlaneSurface & psf)
default

Assignment operator.

◆ operator=() [2/2]

PlaneSurface & Trk::PlaneSurface::operator= ( PlaneSurface && psf)
defaultnoexcept

Move assignment operator.

◆ operator==() [1/2]

bool Trk::PlaneSurface::operator== ( const PlaneSurface & cf) const

◆ operator==() [2/2]

bool Trk::PlaneSurface::operator== ( const Surface & sf) const
overridevirtual

Equality operator.

Implements Trk::Surface.

Definition at line 134 of file PlaneSurface.cxx.

135{
136 // first check the type not to compare apples with oranges
137 if (sf.type()!=Trk::SurfaceType::Plane){
138 return false;
139 }
140 return (*this) == static_cast<const Trk::PlaneSurface&>(sf);
141}

◆ straightLineDistanceEstimate() [1/2]

Trk::DistanceSolution Trk::PlaneSurface::straightLineDistanceEstimate ( const Amg::Vector3D & pos,
const Amg::Vector3D & dir ) const
finaloverridevirtual

fast straight line distance evaluation to Surface

distance to surface

Implements Trk::Surface.

Definition at line 281 of file PlaneSurface.cxx.

282{
283 static const double tol = 0.001;
284
285 const Amg::Vector3D& N = normal();
286
287 const double d = (pos - center()).dot(N);
288
289 const double A = dir.dot(N); // ignore sign
290 if (A == 0.) { // direction parallel to surface
291 if (std::abs(d) < tol) {
292 return {1, 0., true, 0.};
293 }
294 return {0, d, true, 0.};
295
296 }
297
298 return {1, d, true, -d / A};
299}
virtual const Amg::Vector3D & normal() const
Returns the normal vector of the Surface (i.e.
const Amg::Vector3D & center() const
Returns the center position of the Surface.
dot(G, fn, nodesToHighlight=[])
Definition dot.py:5

◆ straightLineDistanceEstimate() [2/2]

Trk::DistanceSolution Trk::PlaneSurface::straightLineDistanceEstimate ( const Amg::Vector3D & pos,
const Amg::Vector3D & dir,
bool Bound ) const
finaloverridevirtual

fast straight line distance evaluation to Surface - with bound option

Implements Trk::Surface.

Definition at line 302 of file PlaneSurface.cxx.

303{
304 const Amg::Transform3D& T = transform();
305 double Az[3] = { T(0, 2), T(1, 2), T(2, 2) };
306
307 // Transformation to plane system coordinates
308 //
309 double dx = pos[0] - T(0, 3);
310 double dy = pos[1] - T(1, 3);
311 double dz = pos[2] - T(2, 3);
312 double z = dx * Az[0] + dy * Az[1] + dz * Az[2];
313 double az = dir[0] * Az[0] + dir[1] * Az[1] + dir[2] * Az[2];
314
315 // Step to surface
316 //
317 int ns = 0;
318 double s = 0.;
319 if (az != 0.) {
320 s = -z / az;
321 ns = 1;
322 }
323 double dist = std::abs(z);
324 if (!bound)
325 return {ns, std::abs(z), true, s};
326
327 // Min distance to surface
328 //
329 double x = dx * T(0, 0) + dy * T(1, 0) + dz * T(2, 0);
330 double y = dx * T(0, 1) + dy * T(1, 1) + dz * T(2, 1);
331
332 Amg::Vector2D lp(x, y);
333
334 double d = bounds().minDistance(lp);
335 if (d > 0.)
336 dist = std::sqrt(dist * dist + d * d);
337
338 return {ns, dist, true, s};
339}
#define y
#define x
#define z
unsigned long long T

◆ straightLineIntersection()

Trk::Intersection Trk::PlaneSurface::straightLineIntersection ( const Amg::Vector3D & pos,
const Amg::Vector3D & dir,
bool forceDir,
Trk::BoundaryCheck bchk ) const
finaloverridevirtual

fast straight line intersection schema - standard: provides closest intersection and (signed) path length forceDir is to provide the closest forward solution

mathematical motivation:

the equation of the plane is given by:
\( \vec n \cdot \vec x = \vec n \cdot \vec p,\)
where \( \vec n = (n_{x}, n_{y}, n_{z})\) denotes the normal vector of the plane, \( \vec p = (p_{x}, p_{y}, p_{z})\) one specific point on the plane and \( \vec x = (x,y,z) \) all possible points on the plane.
Given a line with:
\( \vec l(u) = \vec l_{1} + u \cdot \vec v \),
the solution for \( u \) can be written: \( u = \frac{\vec n (\vec p - \vec l_{1})}{\vec n \vec v}\)
If the denominator is 0 then the line lies:

  • either in the plane
  • perpenticular to the normal of the plane

Implements Trk::Surface.

Definition at line 222 of file PlaneSurface.cxx.

225 {
226 double denom = dir.dot(normal());
227 if (denom) {
228 double u = (normal().dot((center() - pos))) / (denom);
229 Amg::Vector3D intersectPoint(pos + u * dir);
230 // evaluate the intersection in terms of direction
231 bool isValid = forceDir ? (u > 0.) : true;
232 // evaluate (if necessary in terms of boundaries)
233 isValid = bchk ? (isValid && isOnSurface(intersectPoint)) : isValid;
234 // return the result
235 return Trk::Intersection(intersectPoint, u, isValid);
236 }
237 return Trk::Intersection(pos, 0., false);
238}
bool isValid(const T &p)
Av: we implement here an ATLAS-sepcific convention: all particles which are 99xxxxx are fine.
Definition AtlasPID.h:878
virtual bool isOnSurface(const Amg::Vector3D &glopo, const BoundaryCheck &bchk=true, double tol1=0., double tol2=0.) const override final
This method returns true if the GlobalPosition is on the Surface for both, within or without check of...
@ u
Enums for curvilinear frames.
Definition ParamDefs.h:77

◆ type()

virtual constexpr SurfaceType Trk::PlaneSurface::type ( ) const
constexprfinaloverridevirtual

Return the surface type.

Implements Trk::Surface.

◆ ::BoundSurfaceCnv_p1

template<class SURFACE, class BOUNDS_CNV>
friend class ::BoundSurfaceCnv_p1
friend

< data members

Definition at line 280 of file PlaneSurface.h.

◆ ::BoundSurfaceCnv_p2

template<class SURFACE, class BOUNDS_CNV>
friend class ::BoundSurfaceCnv_p2
friend

Definition at line 282 of file PlaneSurface.h.

Member Data Documentation

◆ m_bounds

std::shared_ptr<const SurfaceBounds> Trk::PlaneSurface::m_bounds
protected

bounds (shared)

NoBounds as return object when no bounds are declared

Definition at line 284 of file PlaneSurface.h.

◆ s_boundless

const Trk::NoBounds Trk::PlaneSurface::s_boundless
staticprotected

Definition at line 286 of file PlaneSurface.h.

◆ staticType

SurfaceType Trk::PlaneSurface::staticType = SurfaceType::Plane
staticconstexpr

The surface type static constexpr.

Definition at line 67 of file PlaneSurface.h.


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