LArSoft  v09_90_00
Liquid Argon Software toolkit - https://larsoft.org/
trkf::SurfYZPlane Class Reference

#include "SurfYZPlane.h"

Inheritance diagram for trkf::SurfYZPlane:
trkf::SurfPlane trkf::Surface trkf::SurfWireX

Public Types

enum  TrackDirection { FORWARD, BACKWARD, UNKNOWN }
 Track direction enum. More...
 

Public Member Functions

 SurfYZPlane ()
 Default constructor. More...
 
 SurfYZPlane (double x0, double y0, double z0, double phi)
 Initializing constructor. More...
 
virtual ~SurfYZPlane ()
 Destructor. More...
 
double x0 () const
 X origin. More...
 
double y0 () const
 Y origin. More...
 
double z0 () const
 Z origin. More...
 
double phi () const
 Rotation angle about x-axis. More...
 
virtual Surfaceclone () const
 Clone method. More...
 
virtual bool isTrackValid (const TrackVector &vec) const
 Surface-specific tests of validity of track parameters. More...
 
virtual void toLocal (const double xyz[3], double uvw[3]) const
 Transform global to local coordinates. More...
 
virtual void toGlobal (const double uvw[3], double xyz[3]) const
 Transform local to global coordinates. More...
 
virtual void getPosition (const TrackVector &vec, double xyz[3]) const
 Get position of track. More...
 
virtual void getMomentum (const TrackVector &vec, double mom[3], TrackDirection dir=UNKNOWN) const
 Get momentum vector of track. More...
 
virtual bool isParallel (const Surface &surf) const
 Test whether two surfaces are parallel, within tolerance. More...
 
virtual double distanceTo (const Surface &surf) const
 Find perpendicular forward distance to a parallel surface. More...
 
virtual bool isEqual (const Surface &surf) const
 Test two surfaces for equality, within tolerance. More...
 
virtual std::ostream & Print (std::ostream &out) const
 Printout. More...
 
double PointingError (const TrackVector &vec, const TrackError &err) const
 Get pointing error of track. More...
 
void getStartingError (TrackError &err) const
 Get starting error matrix for Kalman filter. More...
 
virtual TrackVector getDiff (const TrackVector &vec1, const TrackVector &vec2) const
 Calculate difference of two track parameter vectors. More...
 
virtual TrackDirection getDirection (const TrackVector &, TrackDirection dir=UNKNOWN) const
 Get direction of track (default UNKNOWN). More...
 

Private Attributes

double fX0
 X origin. More...
 
double fY0
 Y origin. More...
 
double fZ0
 Z origin. More...
 
double fPhi
 Rotation angle about x-axis. More...
 

Static Private Attributes

static double fPhiTolerance = 1.e-10
 Phi tolerance for parallel. More...
 
static double fSepTolerance = 1.e-6
 Separation tolerance for equal. More...
 

Detailed Description

Definition at line 46 of file SurfYZPlane.h.

Member Enumeration Documentation

Track direction enum.

Enumerator
FORWARD 
BACKWARD 
UNKNOWN 

Definition at line 54 of file Surface.h.

Constructor & Destructor Documentation

trkf::SurfYZPlane::SurfYZPlane ( )

Default constructor.

Definition at line 24 of file SurfYZPlane.cxx.

Referenced by clone(), and trkf::SurfWireX::SurfWireX().

24 : fX0(0.), fY0(0.), fZ0(0.), fPhi(0.) {}
double fPhi
Rotation angle about x-axis.
Definition: SurfYZPlane.h:106
double fZ0
Z origin.
Definition: SurfYZPlane.h:105
double fX0
X origin.
Definition: SurfYZPlane.h:103
double fY0
Y origin.
Definition: SurfYZPlane.h:104
trkf::SurfYZPlane::SurfYZPlane ( double  x0,
double  y0,
double  z0,
double  phi 
)

Initializing constructor.

Initializing constructor.

Arguments:

x0, y0, z0 - Global coordinates of local origin. phi - Rotation angle about x-axis.

Definition at line 33 of file SurfYZPlane.cxx.

34  : fX0(x0), fY0(y0), fZ0(z0), fPhi(phi)
35  {}
double fPhi
Rotation angle about x-axis.
Definition: SurfYZPlane.h:106
double z0() const
Z origin.
Definition: SurfYZPlane.h:60
double x0() const
X origin.
Definition: SurfYZPlane.h:58
double y0() const
Y origin.
Definition: SurfYZPlane.h:59
double phi() const
Rotation angle about x-axis.
Definition: SurfYZPlane.h:61
double fZ0
Z origin.
Definition: SurfYZPlane.h:105
double fX0
X origin.
Definition: SurfYZPlane.h:103
double fY0
Y origin.
Definition: SurfYZPlane.h:104
trkf::SurfYZPlane::~SurfYZPlane ( )
virtual

Destructor.

Definition at line 38 of file SurfYZPlane.cxx.

38 {}

Member Function Documentation

Surface * trkf::SurfYZPlane::clone ( ) const
virtual

Clone method.

Implements trkf::Surface.

Definition at line 41 of file SurfYZPlane.cxx.

References SurfYZPlane().

Referenced by phi().

42  {
43  return new SurfYZPlane(*this);
44  }
SurfYZPlane()
Default constructor.
Definition: SurfYZPlane.cxx:24
double trkf::SurfYZPlane::distanceTo ( const Surface surf) const
virtual

Find perpendicular forward distance to a parallel surface.

Find perpendicular forward distance to a parallel surface.

Throw an exception if the other surface is not parallel.

Assuming the other surface is parallel, the distance is simply the w-coordinate of the other surface, and is signed.

Arguments:

surf - Other surface.

Returned value: Distance.

Implements trkf::Surface.

Definition at line 207 of file SurfYZPlane.cxx.

References isParallel(), trkf::Surface::toGlobal(), and toLocal().

Referenced by phi().

208  {
209  // Check if the other surface is parallel to this one.
210 
211  bool parallel = isParallel(surf);
212  if (!parallel)
213  throw cet::exception("SurfYZPlane") << "Attempt to find distance to non-parallel surface.\n";
214 
215  // Find the origin of the other surface in global coordinates,
216  // then convert to our local coordinates.
217 
218  double otheruvw[3] = {0., 0., 0.};
219  double xyz[3];
220  double myuvw[3];
221  surf.toGlobal(otheruvw, xyz);
222  toLocal(xyz, myuvw);
223 
224  // Distance is local w-coordinate of other surface origin.
225 
226  return myuvw[2];
227  }
virtual bool isParallel(const Surface &surf) const
Test whether two surfaces are parallel, within tolerance.
virtual void toLocal(const double xyz[3], double uvw[3]) const
Transform global to local coordinates.
Definition: SurfYZPlane.cxx:59
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
TrackVector trkf::Surface::getDiff ( const TrackVector vec1,
const TrackVector vec2 
) const
virtualinherited

Calculate difference of two track parameter vectors.

Calculate difference of two track parameter vectors. This method has a default implementation which is just the numeric difference. Surfaces that require a more sophisticated difference (e.g. phi-wrap difference) should override this method.

Reimplemented in trkf::SurfYZLine.

Definition at line 31 of file Surface.cxx.

32  {
33  return vec1 - vec2;
34  }
virtual TrackDirection trkf::Surface::getDirection ( const TrackVector ,
TrackDirection  dir = UNKNOWN 
) const
inlinevirtualinherited
void trkf::SurfYZPlane::getMomentum ( const TrackVector vec,
double  mom[3],
TrackDirection  dir = UNKNOWN 
) const
virtual

Get momentum vector of track.

Get momentum vector of track.

Arguments:

vec - Track state vector. mom - Momentum vector in global coordinate system. dir - Track direction.

Implements trkf::Surface.

Definition at line 126 of file SurfYZPlane.cxx.

References util::abs(), trkf::Surface::BACKWARD, e, trkf::Surface::FORWARD, fPhi, and trkf::Surface::getDirection().

Referenced by phi().

127  {
128 
129  // Get momentum.
130 
131  double invp = std::abs(vec(4));
132  double p = 1. / std::max(invp, 1.e-3); // Capped at 1000. GeV/c.
133 
134  // Get track slope parameters.
135 
136  double dudw = vec(2);
137  double dvdw = vec(3);
138 
139  // Calculate dw/ds.
140 
141  double dwds = 1. / std::sqrt(1. + dudw * dudw + dvdw * dvdw);
142  TrackDirection realdir = getDirection(vec, dir); // Should be same as original direction.
143  if (realdir == BACKWARD)
144  dwds = -dwds;
145  else if (realdir != FORWARD)
146  throw cet::exception("SurfYZPlane") << "Track direction not specified.\n";
147 
148  // Calculate momentum vector in local coordinate system.
149 
150  double pu = p * dudw * dwds;
151  double pv = p * dvdw * dwds;
152  double pw = p * dwds;
153 
154  // Rotate momentum to global coordinte system.
155 
156  double sinphi = std::sin(fPhi);
157  double cosphi = std::cos(fPhi);
158 
159  mom[0] = pu;
160  mom[1] = pv * cosphi - pw * sinphi;
161  mom[2] = pv * sinphi + pw * cosphi;
162 
163  return;
164  }
TrackDirection
Track direction enum.
Definition: Surface.h:54
double fPhi
Rotation angle about x-axis.
Definition: SurfYZPlane.h:106
constexpr auto abs(T v)
Returns the absolute value of the argument.
TDirectory * dir
Definition: macro.C:5
virtual TrackDirection getDirection(const TrackVector &, TrackDirection dir=UNKNOWN) const
Get direction of track (default UNKNOWN).
Definition: Surface.h:83
Float_t e
Definition: plot.C:35
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
void trkf::SurfYZPlane::getPosition ( const TrackVector vec,
double  xyz[3] 
) const
virtual

Get position of track.

Get position of track.

Arguments:

vec - Track state vector. xyz - Position in global coordinate system.

Implements trkf::Surface.

Definition at line 103 of file SurfYZPlane.cxx.

References toGlobal().

Referenced by phi().

104  {
105  // Get position in local coordinate system.
106 
107  double uvw[3];
108  uvw[0] = vec(0);
109  uvw[1] = vec(1);
110  uvw[2] = 0.;
111 
112  // Transform to global coordinate system.
113 
114  toGlobal(uvw, xyz);
115  return;
116  }
virtual void toGlobal(const double uvw[3], double xyz[3]) const
Transform local to global coordinates.
Definition: SurfYZPlane.cxx:81
void trkf::SurfPlane::getStartingError ( TrackError err) const
virtualinherited

Get starting error matrix for Kalman filter.

Get starting error matrix for Kalman filter.

Arguments:

err - Error matrix.

Implements trkf::Surface.

Definition at line 83 of file SurfPlane.cxx.

84  {
85  err.resize(5, false);
86  err.clear();
87  err(0, 0) = 1000.;
88  err(1, 1) = 1000.;
89  err(2, 2) = 0.25;
90  err(3, 3) = 0.25;
91  err(4, 4) = 10.;
92  }
bool trkf::SurfYZPlane::isEqual ( const Surface surf) const
virtual

Test two surfaces for equality, within tolerance.

Test two surfaces for equality, within tolerance. Here equal is defined as having all surface parameters the same, not just having the surfaces coincide spatially, so that the local coordinate systems are the same between the two surfaces.

Arguments:

surf - Other surface.

Returned values: true if equal.

Implements trkf::Surface.

Definition at line 240 of file SurfYZPlane.cxx.

References util::abs(), fPhi, fPhiTolerance, fSepTolerance, fX0, fY0, fZ0, phi(), x0(), y0(), and z0().

Referenced by trkf::KHitWireX::KHitWireX(), and phi().

241  {
242  bool result = false;
243 
244  // Test if the other surface is a SurfYZPlane.
245 
246  const SurfYZPlane* psurf = dynamic_cast<const SurfYZPlane*>(&surf);
247  if (psurf != 0) {
248 
249  // Test whether surface parameters are the same within tolerance.
250 
251  double delta_phi = TVector2::Phi_mpi_pi(fPhi - psurf->phi());
252  double dx = fX0 - psurf->x0();
253  double dy = fY0 - psurf->y0();
254  double dz = fZ0 - psurf->z0();
255  if (std::abs(delta_phi) <= fPhiTolerance && std::abs(dx) <= fSepTolerance &&
256  std::abs(dy) <= fSepTolerance && std::abs(dz) <= fSepTolerance)
257  result = true;
258  }
259  return result;
260  }
static double fSepTolerance
Separation tolerance for equal.
Definition: SurfYZPlane.h:99
double fPhi
Rotation angle about x-axis.
Definition: SurfYZPlane.h:106
static double fPhiTolerance
Phi tolerance for parallel.
Definition: SurfYZPlane.h:98
constexpr auto abs(T v)
Returns the absolute value of the argument.
SurfYZPlane()
Default constructor.
Definition: SurfYZPlane.cxx:24
double fZ0
Z origin.
Definition: SurfYZPlane.h:105
double fX0
X origin.
Definition: SurfYZPlane.h:103
double fY0
Y origin.
Definition: SurfYZPlane.h:104
bool trkf::SurfYZPlane::isParallel ( const Surface surf) const
virtual

Test whether two surfaces are parallel, within tolerance.

Test whether two surfaces are parallel, within tolerance. This method will only return true if the other surface is a SurfYZPlane.

Arguments:

surf - Other surface.

Returned value: true if parallel.

Implements trkf::Surface.

Definition at line 176 of file SurfYZPlane.cxx.

References util::abs(), fPhi, fPhiTolerance, and phi().

Referenced by distanceTo(), and phi().

177  {
178  bool result = false;
179 
180  // Test if the other surface is a SurfYZPlane.
181 
182  const SurfYZPlane* psurf = dynamic_cast<const SurfYZPlane*>(&surf);
183  if (psurf != 0) {
184 
185  // Test whether surface angle parameters are the same
186  // within tolerance.
187 
188  double delta_phi = TVector2::Phi_mpi_pi(fPhi - psurf->phi());
189  if (std::abs(delta_phi) <= fPhiTolerance) result = true;
190  }
191  return result;
192  }
double fPhi
Rotation angle about x-axis.
Definition: SurfYZPlane.h:106
static double fPhiTolerance
Phi tolerance for parallel.
Definition: SurfYZPlane.h:98
constexpr auto abs(T v)
Returns the absolute value of the argument.
SurfYZPlane()
Default constructor.
Definition: SurfYZPlane.cxx:24
bool trkf::SurfYZPlane::isTrackValid ( const TrackVector vec) const
virtual

Surface-specific tests of validity of track parameters.

Implements trkf::Surface.

Definition at line 47 of file SurfYZPlane.cxx.

Referenced by phi().

48  {
49  return true;
50  }
double trkf::SurfPlane::PointingError ( const TrackVector vec,
const TrackError err 
) const
virtualinherited

Get pointing error of track.

Get pointing error of track.

Arguments:

vec - Track parameters. err - Track error matrix.

Returns: Pointing error.

This method calculates the track pointing error based on the slope track paramers and errors (parameters 2 and 3).

Implements trkf::Surface.

Definition at line 34 of file SurfPlane.cxx.

References den.

35  {
36  // Get slope parameters and error matrix.
37 
38  double xp = vec(2);
39  double yp = vec(3);
40  double exx = err(2, 2);
41  double eyy = err(3, 3);
42  double exy = err(3, 2);
43 
44  // Calculate error matrix of pointing unit vector in some coordinate system.
45 
46  double den = 1. + xp * xp + yp * yp;
47  double den3 = den * den * den;
48 
49  double vxx = ((1. + yp * yp) * (1. + yp * yp) * exx + xp * xp * yp * yp * eyy -
50  2. * xp * yp * (1. + yp * yp) * exy) /
51  den3;
52  double vyy = (xp * xp * yp * yp * exx + (1. + xp * xp) * (1. + xp * xp) * eyy -
53  2. * xp * yp * (1. + xp * xp) * exy) /
54  den3;
55  double vzz = (xp * xp * exx + yp * yp * eyy + 2. * xp * yp * exy) / den3;
56 
57  double vxy = (-xp * yp * (1. + yp * yp) * exx - xp * yp * (1. + xp * xp) * eyy +
58  (1. + xp * xp + yp * yp + 2. * xp * xp * yp * yp) * exy) /
59  den3;
60  double vyz =
61  (xp * xp * yp * exx - yp * (1. + xp * xp) * eyy - xp * (1. + xp * xp - yp * yp) * exy) / den3;
62  double vxz =
63  (-xp * (1. + yp * yp) * exx + xp * yp * yp * eyy - yp * (1. - xp * xp + yp * yp) * exy) /
64  den3;
65 
66  // Calculate square root of the largest eigenvalue of error matrix.
67 
68  double ddd2 = vxx * vxx + vyy * vyy + vzz * vzz - 2. * vxx * vyy - 2. * vxx * vzz -
69  2. * vyy * vzz + 4. * vxy * vxy + 4. * vyz * vyz + 4. * vxz * vxz;
70  double ddd = sqrt(ddd2 > 0. ? ddd2 : 0.);
71  double lambda2 = 0.5 * (vxx + vyy + vzz + ddd);
72  double lambda = sqrt(lambda2 > 0. ? lambda2 : 0.);
73 
74  return lambda;
75  }
Float_t den
Definition: plot.C:35
std::ostream & trkf::SurfYZPlane::Print ( std::ostream &  out) const
virtual

Printout.

Implements trkf::Surface.

Definition at line 263 of file SurfYZPlane.cxx.

References fPhi, fX0, fY0, and fZ0.

Referenced by phi().

264  {
265  out << "SurfYZPlane{ x0=" << fX0 << ", y0=" << fY0 << ", z0=" << fZ0 << ", phi=" << fPhi << "}";
266  return out;
267  }
double fPhi
Rotation angle about x-axis.
Definition: SurfYZPlane.h:106
double fZ0
Z origin.
Definition: SurfYZPlane.h:105
double fX0
X origin.
Definition: SurfYZPlane.h:103
double fY0
Y origin.
Definition: SurfYZPlane.h:104
void trkf::SurfYZPlane::toGlobal ( const double  uvw[3],
double  xyz[3] 
) const
virtual

Transform local to global coordinates.

Transform local to global coordinates.

Arguments:

uvw - Cartesian coordinates in local coordinate system. xyz - Cartesian coordinates in global coordinate system.

Implements trkf::Surface.

Definition at line 81 of file SurfYZPlane.cxx.

References fPhi, fX0, fY0, and fZ0.

Referenced by getPosition(), and phi().

82  {
83  double sinphi = std::sin(fPhi);
84  double cosphi = std::cos(fPhi);
85 
86  // x = x0 + u
87  xyz[0] = fX0 + uvw[0];
88 
89  // y = y0 + v*cos(phi) - w*sin(phi)
90  xyz[1] = fY0 + uvw[1] * cosphi - uvw[2] * sinphi;
91 
92  // z = z0 + v*sin(phi) + w*cos(phi)
93  xyz[2] = fZ0 + uvw[1] * sinphi + uvw[2] * cosphi;
94  }
double fPhi
Rotation angle about x-axis.
Definition: SurfYZPlane.h:106
double fZ0
Z origin.
Definition: SurfYZPlane.h:105
double fX0
X origin.
Definition: SurfYZPlane.h:103
double fY0
Y origin.
Definition: SurfYZPlane.h:104
void trkf::SurfYZPlane::toLocal ( const double  xyz[3],
double  uvw[3] 
) const
virtual

Transform global to local coordinates.

Transform global to local coordinates.

Arguments:

xyz - Cartesian coordinates in global coordinate system. uvw - Cartesian coordinates in local coordinate system.

Implements trkf::Surface.

Definition at line 59 of file SurfYZPlane.cxx.

References fPhi, fX0, fY0, and fZ0.

Referenced by distanceTo(), and phi().

60  {
61  double sinphi = std::sin(fPhi);
62  double cosphi = std::cos(fPhi);
63 
64  // u = x-x0
65  uvw[0] = xyz[0] - fX0;
66 
67  // v = (y-y0)*cos(phi) + (z-z0)*sin(phi)
68  uvw[1] = (xyz[1] - fY0) * cosphi + (xyz[2] - fZ0) * sinphi;
69 
70  // w = -(y-y0)*sin(phi) + (z-z0)*cos(phi)
71  uvw[2] = -(xyz[1] - fY0) * sinphi + (xyz[2] - fZ0) * cosphi;
72  }
double fPhi
Rotation angle about x-axis.
Definition: SurfYZPlane.h:106
double fZ0
Z origin.
Definition: SurfYZPlane.h:105
double fX0
X origin.
Definition: SurfYZPlane.h:103
double fY0
Y origin.
Definition: SurfYZPlane.h:104
double trkf::SurfYZPlane::x0 ( ) const
inline

X origin.

Definition at line 58 of file SurfYZPlane.h.

References fX0.

Referenced by isEqual(), recob::tracking::makePlane(), and trkf::PropYZPlane::short_vec_prop().

double trkf::SurfYZPlane::y0 ( ) const
inline

Y origin.

Definition at line 59 of file SurfYZPlane.h.

References fY0.

Referenced by isEqual(), recob::tracking::makePlane(), and trkf::PropYZPlane::short_vec_prop().

double trkf::SurfYZPlane::z0 ( ) const
inline

Z origin.

Definition at line 60 of file SurfYZPlane.h.

References fZ0.

Referenced by isEqual(), recob::tracking::makePlane(), and trkf::PropYZPlane::short_vec_prop().

Member Data Documentation

double trkf::SurfYZPlane::fPhi
private

Rotation angle about x-axis.

Definition at line 106 of file SurfYZPlane.h.

Referenced by getMomentum(), isEqual(), isParallel(), phi(), Print(), toGlobal(), and toLocal().

double trkf::SurfYZPlane::fPhiTolerance = 1.e-10
staticprivate

Phi tolerance for parallel.

Definition at line 98 of file SurfYZPlane.h.

Referenced by isEqual(), and isParallel().

double trkf::SurfYZPlane::fSepTolerance = 1.e-6
staticprivate

Separation tolerance for equal.

Definition at line 99 of file SurfYZPlane.h.

Referenced by isEqual().

double trkf::SurfYZPlane::fX0
private

X origin.

Definition at line 103 of file SurfYZPlane.h.

Referenced by isEqual(), Print(), toGlobal(), toLocal(), and x0().

double trkf::SurfYZPlane::fY0
private

Y origin.

Definition at line 104 of file SurfYZPlane.h.

Referenced by isEqual(), Print(), toGlobal(), toLocal(), and y0().

double trkf::SurfYZPlane::fZ0
private

Z origin.

Definition at line 105 of file SurfYZPlane.h.

Referenced by isEqual(), Print(), toGlobal(), toLocal(), and z0().


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