PandaRoot
PndCATrackParamVector Class Reference

#include <PndCATrackParamVector.h>

Classes

struct  PndCATrackFitParam
 

Public Member Functions

 PndCATrackParamVector ()
 
void Reset ()
 
void InitCovMatrix (float_v d2QMom=0.f)
 
void InitByTarget (const PndCATarget &target)
 
void InitByHit (const PndCAHitV &hit, const PndCAParam &param, const float_v &dQMom)
 
void InitDirection (float_v r0, float_v r1, float_v r2)
 
float_v X () const
 
float_v Y () const
 
float_v Z () const
 
float_v SinPhi () const
 
float_v DzDs () const
 
float_v QPt () const
 
float_v Angle () const
 
int_v ISec () const
 
float_v X0 () const
 
float_v X1 () const
 
float_v X2 () const
 
float_v Tx1 () const
 
float_v Tx2 () const
 
float_v QMomentum () const
 
float_v SignCosPhi () const
 
float_v Chi2 () const
 
int_v NDF () const
 
float_v Err2Y () const
 
float_v Err2Z () const
 
float_v Err2SinPhi () const
 
float_v Err2DzDs () const
 
float_v Err2QPt () const
 
float_v GetX () const
 
float_v GetY () const
 
float_v GetZ () const
 
float_v GetSinPhi () const
 
float_v GetDzDs () const
 
float_v GetQPt () const
 
float_v GetSignCosPhi () const
 
float_v GetChi2 () const
 
int_v GetNDF () const
 
float_v GetKappa (const float_v &Bz) const
 
float_v GetCosPhiPositive () const
 
float_v GetCosPhi () const
 
float_v Err2X1 () const
 
float_v Err2X2 () const
 
float_v Err2QMomentum () const
 
const float_v & Par (int i) const
 
const float_v & Cov (int i) const
 
void SetTrackParam (const PndCATrackParamVector &param, const float_m &m=float_m(true))
 
void SetTrackParam (const PndCATrackParamVector &p, int k, int m)
 
void SetPar (int i, const float_v &v)
 
void SetPar (int i, const float_v &v, const float_m &m)
 
void SetCov (int i, const float_v &v)
 
void SetCov (int i, const float_v &v, const float_m &m)
 
void SetX (const float_v &v)
 
void SetY (const float_v &v)
 
void SetZ (const float_v &v)
 
void SetX (const float_v &v, const float_m &m)
 
void SetY (const float_v &v, const float_m &m)
 
void SetZ (const float_v &v, const float_m &m)
 
void SetSinPhi (const float_v &v)
 
void SetSinPhi (const float_v &v, const float_m &m)
 
void SetDzDs (const float_v &v)
 
void SetDzDs (const float_v &v, const float_m &m)
 
void SetQPt (const float_v &v)
 
void SetQMomentum (const float_v &v)
 
void SetQPt (const float_v &v, const float_m &m)
 
void SetSignCosPhi (const float_v &v)
 
void SetSignCosPhi (const float_v &v, const float_m &m)
 
void SetChi2 (const float_v &v)
 
void SetChi2 (const float_v &v, const float_m &m)
 
void SetNDF (int v)
 
void SetNDF (const int_v &v)
 
void SetNDF (const int_v &v, const int_m &m)
 
void SetAngle (const float_v &v)
 
void SetAngle (const float_v &v, const float_m &m)
 
void SetISec (const int_v &v)
 
void SetISec (const int_v &v, const int_m &m)
 
void SetErr2Y (float_v v)
 
void SetErr2Z (float_v v)
 
void SetErr2QPt (float_v v)
 
float_v GetDist2 (const PndCATrackParamVector &t) const
 
float_v GetDistXZ2 (const PndCATrackParamVector &t) const
 
float_v GetS (const float_v &x, const float_v &y, const float_v &Bz) const
 
void GetDCAPoint (const float_v &x, const float_v &y, const float_v &z, float_v *px, float_v *py, float_v *pz, const float_v &Bz) const
 
float_m Transport0 (const int_v &ista, const PndCAParam &param, const float_m &mask=float_m(true))
 
float_m Transport0 (const PndCAHitV &hit, const PndCAParam &param, const float_m &mask=float_m(true))
 
float_m Transport0 (const PndCAHit &hit, const PndCAParam &p, const float_m &mask=float_m(true))
 
float_m TransportToXWithMaterial (const float_v &x, const float_v &XOverX0, const float_v &XThimesRho, const float_v &Bz, const float maxSinPhi=.999f)
 
float_m TransportToX (const float_v &x, const float_v &Bz, const float maxSinPhi=.999f, const float_m &mask=float_m(true))
 
float_m TransportToX (const float_v &x, PndCATrackLinearisationVector &t0, const float_v &Bz, const float maxSinPhi=.999f, float_v *DL=0, const float_m &mask=float_m(true))
 
float_m Transport0ToX (const float_v &x, const float_v &Bz, const float_m &mask)
 
float_m TransportToX (const float_v &x, const float_v &sinPhi0, const float_v &Bz, const float_v maxSinPhi=.999f, const float_m &mask=float_m(true))
 
float_m TransportToXWithMaterial (const float_v &x, PndCATrackLinearisationVector &t0, PndCATrackFitParam &par, const float_v &XOverX0, const float_v &XThimesRho, const float_v &Bz, const float maxSinPhi=.999f, const float_m &mask=float_m(true))
 
float_m TransportToXWithMaterial (const float_v &x, PndCATrackFitParam &par, const float_v &XOverX0, const float_v &XThimesRho, const float_v &Bz, const float maxSinPhi=.999f)
 
float_m Rotate (const float_v &alpha, PndCATrackLinearisationVector &t0, const float maxSinPhi=.999f, const float_m &mask=float_m(true))
 
float_m Rotate (const float_v &alpha, const float maxSinPhi=.999f, const float_m &mask=float_m(true))
 
void RotateXY (float_v alpha, float_v &x, float_v &y, float_v &sin, const float_m &mask=float_m(true)) const
 
float_m FilterWithMaterial (const float_v &y, const float_v &z, float_v err2Y, float_v errYZ, float_v err2Z, float maxSinPhi=0.999f, const float_m &mask=float_m(true), const int_v &hitNDF=int_v(2), const float_v &chi2Cut=10e10f)
 
float_m FilterWithMaterial (const float_v &y, const float_v &z, const PndCAStripInfo &info, float_v err2, float maxSinPhi=0.999f, const float_m &mask=float_m(true), const float_v &chi2Cut=10e10f)
 
void CalculateFitParameters (PndCATrackFitParam &par, const float_v &mass=0.13957f)
 
float_m CorrectForMeanMaterial (const float_v &xOverX0, const float_v &xTimesRho, const PndCATrackFitParam &par, const float_m &_mask)
 
float_m FilterVtx (const float_v &xV, const float_v &yV, const PndCAX1X2MeasurementInfo &info, float_v &extrDx, float_v &extrDy, float_v J[], const float_m &active=float_m(true))
 
float_m TransportJ0ToX0 (const float_v &x0, const float_v &cBz, float_v &extrDy, float_v &extrDz, float_v J[], const float_m &active=float_m(true))
 
float_m Transport (const int_v &ista, const PndCAParam &param, const float_m &mask=float_m(true))
 
float_m Transport (const PndCAHitV &hit, const PndCAParam &p, const float_m &mask=float_m(true))
 
float_m Filter (const PndCAHitV &hit, const PndCAParam &param, const float_m &mask=float_m(true), const float_v &chi2Cut=10e10f)
 
float_m Transport (const PndCAHit &hit, const PndCAParam &p, const float_m &mask=float_m(true))
 
float_m Filter (const PndCAHit &hit, const PndCAParam &param, const float_m &mask=float_m(true), const float_v &chi2Cut=10e10f)
 
float_m Accept (const PndCAHit &hit, const PndCAParam &param, const float_m &mask=float_m(true), const float_v &chi2Cut=10e10f) const
 
float_m AcceptWithMaterial (const float_v &y, const float_v &z, float_v err2Y, float_v errYZ, float_v err2Z, float maxSinPhi=0.999f, const float_m &mask=float_m(true), const int_v &hitNDF=int_v(2), const float_v &chi2Cut=10e10f) const
 
float_m AcceptWithMaterial (const float_v &y, const float_v &z, const PndCAStripInfo &info, float_v err2, float maxSinPhi=0.999f, const float_m &mask=float_m(true), const float_v &chi2Cut=10e10f) const
 
float_m AddTarget (const PndCATarget &target, const float_m &mask=float_m(true))
 

Static Public Member Functions

static float_v ApproximateBetheBloch (const float_v &beta2)
 
static float_v BetheBlochGeant (const float_v &bg, const float_v &kp0=2.33f, const float_v &kp1=0.20f, const float_v &kp2=3.00f, const float_v &kp3=173e-9f, const float_v &kp4=0.49848f)
 
static float_v BetheBlochSolid (const float_v &bg)
 
static float_v BetheBlochGas (const float_v &bg)
 

Public Attributes

float_v fX
 
float_v fSignCosPhi
 
float_v fP [5]
 
float_v fC [15]
 
float_v fChi2
 
int_v fNDF
 
float_v fAlpha
 
int_v fISec
 

Friends

class PndCATrackParam
 

Detailed Description

PndCATrackParamVector class describes the track parametrisation which is used by the PndCATracker slice tracker.

Definition at line 48 of file PndCATrackParamVector.h.

Constructor & Destructor Documentation

◆ PndCATrackParamVector()

PndCATrackParamVector::PndCATrackParamVector ( )
inline

Definition at line 52 of file PndCATrackParamVector.h.

References fC, fP, and i.

52  : fX(Vc::Zero), fSignCosPhi(Vc::Zero), fChi2(Vc::Zero), fNDF(Vc::Zero), fISec(Vc::Zero)
53  {
54  for (int i = 0; i < 5; ++i)
55  fP[i].setZero();
56  for (int i = 0; i < 15; ++i)
57  fC[i].setZero();
58  }
unsigned int i
Definition: P4_F32vec4.h:33

Member Function Documentation

◆ Accept()

float_m PndCATrackParamVector::Accept ( const PndCAHit hit,
const PndCAParam param,
const float_m &  mask = float_m(true),
const float_v &  chi2Cut = 10e10f 
) const

Referenced by SetErr2QPt().

◆ AcceptWithMaterial() [1/2]

float_m PndCATrackParamVector::AcceptWithMaterial ( const float_v &  y,
const float_v &  z,
float_v  err2Y,
float_v  errYZ,
float_v  err2Z,
float  maxSinPhi = 0.999f,
const float_m &  mask = float_m(true),
const int_v &  hitNDF = int_v(2),
const float_v &  chi2Cut = 10e10f 
) const

Referenced by SetErr2QPt().

◆ AcceptWithMaterial() [2/2]

float_m PndCATrackParamVector::AcceptWithMaterial ( const float_v &  y,
const float_v &  z,
const PndCAStripInfo info,
float_v  err2,
float  maxSinPhi = 0.999f,
const float_m &  mask = float_m(true),
const float_v &  chi2Cut = 10e10f 
) const

◆ AddTarget()

float_m PndCATrackParamVector::AddTarget ( const PndCATarget target,
const float_m &  mask = float_m(true) 
)

Referenced by SetErr2QPt().

◆ Angle()

float_v PndCATrackParamVector::Angle ( ) const
inline

Definition at line 105 of file PndCATrackParamVector.h.

References fAlpha.

105 { return fAlpha; }

◆ ApproximateBetheBloch()

static float_v PndCATrackParamVector::ApproximateBetheBloch ( const float_v &  beta2)
static

Referenced by SetErr2QPt().

◆ BetheBlochGas()

static float_v PndCATrackParamVector::BetheBlochGas ( const float_v &  bg)
static

Referenced by SetErr2QPt().

◆ BetheBlochGeant()

static float_v PndCATrackParamVector::BetheBlochGeant ( const float_v &  bg,
const float_v &  kp0 = 2.33f,
const float_v &  kp1 = 0.20f,
const float_v &  kp2 = 3.00f,
const float_v &  kp3 = 173e-9f,
const float_v &  kp4 = 0.49848f 
)
static

Referenced by SetErr2QPt().

◆ BetheBlochSolid()

static float_v PndCATrackParamVector::BetheBlochSolid ( const float_v &  bg)
static

Referenced by SetErr2QPt().

◆ CalculateFitParameters()

void PndCATrackParamVector::CalculateFitParameters ( PndCATrackFitParam par,
const float_v &  mass = 0.13957f 
)

Referenced by SetErr2QPt().

◆ Chi2()

float_v PndCATrackParamVector::Chi2 ( ) const
inline

Definition at line 120 of file PndCATrackParamVector.h.

References fChi2.

120 { return fChi2; }

◆ CorrectForMeanMaterial()

float_m PndCATrackParamVector::CorrectForMeanMaterial ( const float_v &  xOverX0,
const float_v &  xTimesRho,
const PndCATrackFitParam par,
const float_m &  _mask 
)

Referenced by SetErr2QPt().

◆ Cov()

const float_v& PndCATrackParamVector::Cov ( int  i) const
inline

Definition at line 150 of file PndCATrackParamVector.h.

References fC, and i.

Referenced by PndCATrackParam::PndCATrackParam(), PndFTSCATrackParam::PndFTSCATrackParam(), and SetTrackParam().

150 { return fC[i]; }
unsigned int i
Definition: P4_F32vec4.h:33

◆ DzDs()

float_v PndCATrackParamVector::DzDs ( ) const
inline

Definition at line 102 of file PndCATrackParamVector.h.

References fP.

Referenced by TransportJ0ToX0(), TransportToX(), and Tx2().

102 { return fP[3]; }

◆ Err2DzDs()

float_v PndCATrackParamVector::Err2DzDs ( ) const
inline

Definition at line 126 of file PndCATrackParamVector.h.

References fC.

126 { return fC[9]; }

◆ Err2QMomentum()

float_v PndCATrackParamVector::Err2QMomentum ( ) const
inline

Definition at line 147 of file PndCATrackParamVector.h.

References fC.

147 { return fC[14]; }

◆ Err2QPt()

float_v PndCATrackParamVector::Err2QPt ( ) const
inline

Definition at line 127 of file PndCATrackParamVector.h.

References fC.

127 { return fC[14]; }

◆ Err2SinPhi()

float_v PndCATrackParamVector::Err2SinPhi ( ) const
inline

Definition at line 125 of file PndCATrackParamVector.h.

References fC.

125 { return fC[5]; }

◆ Err2X1()

float_v PndCATrackParamVector::Err2X1 ( ) const
inline

Definition at line 143 of file PndCATrackParamVector.h.

References fC.

143 { return fC[0]; }

◆ Err2X2()

float_v PndCATrackParamVector::Err2X2 ( ) const
inline

Definition at line 144 of file PndCATrackParamVector.h.

References fC.

144 { return fC[2]; }

◆ Err2Y()

float_v PndCATrackParamVector::Err2Y ( ) const
inline

Definition at line 123 of file PndCATrackParamVector.h.

References fC.

123 { return fC[0]; }

◆ Err2Z()

float_v PndCATrackParamVector::Err2Z ( ) const
inline

Definition at line 124 of file PndCATrackParamVector.h.

References fC.

124 { return fC[2]; }

◆ Filter() [1/2]

float_m PndCATrackParamVector::Filter ( const PndCAHitV hit,
const PndCAParam param,
const float_m &  mask = float_m(true),
const float_v &  chi2Cut = 10e10f 
)

Referenced by SetErr2QPt().

◆ Filter() [2/2]

float_m PndCATrackParamVector::Filter ( const PndCAHit hit,
const PndCAParam param,
const float_m &  mask = float_m(true),
const float_v &  chi2Cut = 10e10f 
)

◆ FilterVtx()

float_m PndCATrackParamVector::FilterVtx ( const float_v &  xV,
const float_v &  yV,
const PndCAX1X2MeasurementInfo info,
float_v &  extrDx,
float_v &  extrDy,
float_v  J[],
const float_m &  active = float_m(true) 
)
inline

Definition at line 446 of file PndCATrackParamVector.h.

References PndCAX1X2MeasurementInfo::C00(), PndCAX1X2MeasurementInfo::C10(), PndCAX1X2MeasurementInfo::C11(), fC, fChi2, fP, Y(), and Z().

Referenced by SetErr2QPt().

447 {
448  float_v zeta0, zeta1, S00, S10, S11, si;
449  float_v F00, F10, F20, F30, F40, F01, F11, F21, F31, F41;
450  float_v K00, K10, K20, K30, K40, K01, K11, K21, K31, K41;
451  float_v &c00 = fC[0];
452  float_v &c10 = fC[1];
453  float_v &c11 = fC[2];
454  float_v &c20 = fC[3];
455  float_v &c21 = fC[4];
456  float_v &c22 = fC[5];
457  float_v &c30 = fC[6];
458  float_v &c31 = fC[7];
459  float_v &c32 = fC[8];
460  float_v &c33 = fC[9];
461  float_v &c40 = fC[10];
462  float_v &c41 = fC[11];
463  float_v &c42 = fC[12];
464  float_v &c43 = fC[13];
465  float_v &c44 = fC[14];
466 
467  zeta0 = Y() + extrDy - yV;
468  zeta1 = Z() + extrDz - zV;
469 
470  // H = 1 0 J[0] J[1] J[2]
471  // 0 1 J[3] J[4] J[5]
472 
473  // F = CH'
474  F00 = c00;
475  F01 = c10;
476  F10 = c10;
477  F11 = c11;
478  F20 = J[0] * c22;
479  F21 = J[3] * c22;
480  F30 = J[1] * c33;
481  F31 = J[4] * c33;
482  F40 = J[2] * c44;
483  F41 = J[5] * c44;
484 
485  S00 = info.C00() + F00 + J[0] * F20 + J[1] * F30 + J[2] * F40;
486  S10 = info.C10() + F10 + J[3] * F20 + J[4] * F30 + J[5] * F40;
487  S11 = info.C11() + F11 + J[3] * F21 + J[4] * F31 + J[5] * F41;
488 
489  si = 1.f / (S00 * S11 - S10 * S10);
490  float_v S00tmp = S00;
491  S00 = si * S11;
492  S10 = -si * S10;
493  S11 = si * S00tmp;
494 
495  fChi2(active) += zeta0 * zeta0 * S00 + 2.f * zeta0 * zeta1 * S10 + zeta1 * zeta1 * S11;
496 
497  K00 = F00 * S00 + F01 * S10;
498  K01 = F00 * S10 + F01 * S11;
499  K10 = F10 * S00 + F11 * S10;
500  K11 = F10 * S10 + F11 * S11;
501  K20 = F20 * S00 + F21 * S10;
502  K21 = F20 * S10 + F21 * S11;
503  K30 = F30 * S00 + F31 * S10;
504  K31 = F30 * S10 + F31 * S11;
505  K40 = F40 * S00 + F41 * S10;
506  K41 = F40 * S10 + F41 * S11;
507 
508  fP[0](active) -= K00 * zeta0 + K01 * zeta1;
509  fP[1](active) -= K10 * zeta0 + K11 * zeta1;
510  fP[2](active) -= K20 * zeta0 + K21 * zeta1;
511  fP[3](active) -= K30 * zeta0 + K31 * zeta1;
512  fP[4](active) -= K40 * zeta0 + K41 * zeta1;
513 
514  c00(active) -= (K00 * F00 + K01 * F01);
515  c10(active) -= (K10 * F00 + K11 * F01);
516  c11(active) -= (K10 * F10 + K11 * F11);
517  c20(active) = -(K20 * F00 + K21 * F01);
518  c21(active) = -(K20 * F10 + K21 * F11);
519  c22(active) -= (K20 * F20 + K21 * F21);
520  c30(active) = -(K30 * F00 + K31 * F01);
521  c31(active) = -(K30 * F10 + K31 * F11);
522  c32(active) = -(K30 * F20 + K31 * F21);
523  c33(active) -= (K30 * F30 + K31 * F31);
524  c40(active) = -(K40 * F00 + K41 * F01);
525  c41(active) = -(K40 * F10 + K41 * F11);
526  c42(active) = -(K40 * F20 + K41 * F21);
527  c43(active) = -(K40 * F30 + K41 * F31);
528  c44(active) -= (K40 * F40 + K41 * F41);
529 
530  return active;
531 }
const float_v & C00() const
const float_v & C11() const
const float_v & C10() const

◆ FilterWithMaterial() [1/2]

float_m PndCATrackParamVector::FilterWithMaterial ( const float_v &  y,
const float_v &  z,
float_v  err2Y,
float_v  errYZ,
float_v  err2Z,
float  maxSinPhi = 0.999f,
const float_m &  mask = float_m(true),
const int_v &  hitNDF = int_v(2),
const float_v &  chi2Cut = 10e10f 
)

Referenced by SetErr2QPt().

◆ FilterWithMaterial() [2/2]

float_m PndCATrackParamVector::FilterWithMaterial ( const float_v &  y,
const float_v &  z,
const PndCAStripInfo info,
float_v  err2,
float  maxSinPhi = 0.999f,
const float_m &  mask = float_m(true),
const float_v &  chi2Cut = 10e10f 
)

◆ GetChi2()

float_v PndCATrackParamVector::GetChi2 ( ) const
inline

Definition at line 136 of file PndCATrackParamVector.h.

References fChi2.

Referenced by SetTrackParam().

136 { return fChi2; }

◆ GetCosPhi()

float_v PndCATrackParamVector::GetCosPhi ( ) const
inline

Definition at line 141 of file PndCATrackParamVector.h.

References fSignCosPhi, SinPhi(), and CAMath::Sqrt().

Referenced by Rotate(), and RotateXY().

141 { return fSignCosPhi * CAMath::Sqrt(float_v(Vc::One) - SinPhi() * SinPhi()); }
static T Sqrt(const T &x)
Definition: PndCAMath.h:57

◆ GetCosPhiPositive()

float_v PndCATrackParamVector::GetCosPhiPositive ( ) const
inline

Definition at line 140 of file PndCATrackParamVector.h.

References SinPhi(), and CAMath::Sqrt().

140 { return CAMath::Sqrt(float_v(Vc::One) - SinPhi() * SinPhi()); }
static T Sqrt(const T &x)
Definition: PndCAMath.h:57

◆ GetDCAPoint()

void PndCATrackParamVector::GetDCAPoint ( const float_v &  x,
const float_v &  y,
const float_v &  z,
float_v *  px,
float_v *  py,
float_v *  pz,
const float_v &  Bz 
) const

Referenced by SetErr2QPt().

◆ GetDist2()

float_v PndCATrackParamVector::GetDist2 ( const PndCATrackParamVector t) const

Referenced by SetErr2QPt().

◆ GetDistXZ2()

float_v PndCATrackParamVector::GetDistXZ2 ( const PndCATrackParamVector t) const

Referenced by SetErr2QPt().

◆ GetDzDs()

float_v PndCATrackParamVector::GetDzDs ( ) const
inline

Definition at line 133 of file PndCATrackParamVector.h.

References fP.

133 { return fP[3]; }

◆ GetKappa()

float_v PndCATrackParamVector::GetKappa ( const float_v &  Bz) const
inline

Definition at line 139 of file PndCATrackParamVector.h.

References fP.

139 { return fP[4] * Bz; }

◆ GetNDF()

int_v PndCATrackParamVector::GetNDF ( ) const
inline

Definition at line 137 of file PndCATrackParamVector.h.

References fNDF.

Referenced by SetTrackParam().

137 { return fNDF; }

◆ GetQPt()

float_v PndCATrackParamVector::GetQPt ( ) const
inline

Definition at line 134 of file PndCATrackParamVector.h.

References fP.

134 { return fP[4]; }

◆ GetS()

float_v PndCATrackParamVector::GetS ( const float_v &  x,
const float_v &  y,
const float_v &  Bz 
) const

Referenced by SetErr2QPt().

◆ GetSignCosPhi()

float_v PndCATrackParamVector::GetSignCosPhi ( ) const
inline

Definition at line 135 of file PndCATrackParamVector.h.

References fSignCosPhi.

135 { return fSignCosPhi; }

◆ GetSinPhi()

float_v PndCATrackParamVector::GetSinPhi ( ) const
inline

Definition at line 132 of file PndCATrackParamVector.h.

References fP.

132 { return fP[2]; }

◆ GetX()

float_v PndCATrackParamVector::GetX ( ) const
inline

Definition at line 129 of file PndCATrackParamVector.h.

References fX.

129 { return fX; }

◆ GetY()

float_v PndCATrackParamVector::GetY ( ) const
inline

Definition at line 130 of file PndCATrackParamVector.h.

References fP.

130 { return fP[0]; }

◆ GetZ()

float_v PndCATrackParamVector::GetZ ( ) const
inline

Definition at line 131 of file PndCATrackParamVector.h.

References fP.

131 { return fP[1]; }

◆ InitByHit()

void PndCATrackParamVector::InitByHit ( const PndCAHitV hit,
const PndCAParam param,
const float_v &  dQMom 
)

Referenced by Reset().

◆ InitByTarget()

void PndCATrackParamVector::InitByTarget ( const PndCATarget target)

Referenced by Reset().

◆ InitCovMatrix()

void PndCATrackParamVector::InitCovMatrix ( float_v  d2QMom = 0.f)

Referenced by Reset().

◆ InitDirection()

void PndCATrackParamVector::InitDirection ( float_v  r0,
float_v  r1,
float_v  r2 
)
inline

Definition at line 78 of file PndCATrackParamVector.h.

References f, SetDzDs(), SetSignCosPhi(), SetSinPhi(), and sqrt().

79  {
80  const float_v r = sqrt(r0 * r0 + r1 * r1);
81  SetSinPhi(r1 / r);
82  SetSignCosPhi(1.f); // r0/abs(r0) );
83  SetDzDs(r2 / r);
84  }
void SetSignCosPhi(const float_v &v)
void SetSinPhi(const float_v &v)
friend F32vec4 sqrt(const F32vec4 &a)
Definition: P4_F32vec4.h:40
void SetDzDs(const float_v &v)
float f
Definition: P4_F32vec4.h:32

◆ ISec()

int_v PndCATrackParamVector::ISec ( ) const
inline

Definition at line 106 of file PndCATrackParamVector.h.

References fISec.

106 { return fISec; }

◆ NDF()

int_v PndCATrackParamVector::NDF ( ) const
inline

Definition at line 121 of file PndCATrackParamVector.h.

References fNDF.

121 { return fNDF; }

◆ Par()

const float_v& PndCATrackParamVector::Par ( int  i) const
inline

Definition at line 149 of file PndCATrackParamVector.h.

References fP, and i.

Referenced by PndCATrackParam::PndCATrackParam(), PndFTSCATrackParam::PndFTSCATrackParam(), and SetTrackParam().

149 { return fP[i]; }
unsigned int i
Definition: P4_F32vec4.h:33

◆ QMomentum()

float_v PndCATrackParamVector::QMomentum ( ) const
inline

Definition at line 113 of file PndCATrackParamVector.h.

References QPt().

113 { return QPt(); } // used for triplets comparison

◆ QPt()

float_v PndCATrackParamVector::QPt ( ) const
inline

Definition at line 103 of file PndCATrackParamVector.h.

References fP.

Referenced by QMomentum(), and TransportToX().

103 { return fP[4]; }

◆ Reset()

void PndCATrackParamVector::Reset ( )
inline

Definition at line 60 of file PndCATrackParamVector.h.

References f, fC, fChi2, fNDF, fP, fSignCosPhi, fX, i, InitByHit(), InitByTarget(), and InitCovMatrix().

61  {
62  fX.setZero();
63  fSignCosPhi = float_v(1.f);
64  fChi2.setZero();
65  fNDF.setZero();
66  // fISec.setZero();
67  {
68  for (int i = 0; i < 5; ++i)
69  fP[i].setZero();
70  for (int i = 0; i < 15; ++i)
71  fC[i].setZero();
72  }
73  }
unsigned int i
Definition: P4_F32vec4.h:33
float f
Definition: P4_F32vec4.h:32

◆ Rotate() [1/2]

float_m PndCATrackParamVector::Rotate ( const float_v &  alpha,
PndCATrackLinearisationVector t0,
const float  maxSinPhi = .999f,
const float_m &  mask = float_m(true) 
)

Referenced by SetErr2QPt(), and Transport0().

◆ Rotate() [2/2]

float_m PndCATrackParamVector::Rotate ( const float_v &  alpha,
const float  maxSinPhi = .999f,
const float_m &  mask = float_m(true) 
)
inline

Definition at line 381 of file PndCATrackParamVector.h.

References CAMath::Abs(), alpha, CAMath::Cos(), f, fAlpha, fC, GetCosPhi(), SetSinPhi(), SetX(), SetY(), CAMath::Sin(), SinPhi(), X(), and Y().

382 {
383  // Rotate the coordinate system in XY on the angle alpha
384  if ((CAMath::Abs(alpha) < 1e-6f || !mask).isFull())
385  return mask;
386 
387  const float_v cA = CAMath::Cos(alpha);
388  const float_v sA = CAMath::Sin(alpha);
389  const float_v x0 = X(), y0 = Y(), sP = SinPhi(), cP = GetCosPhi();
390  const float_v cosPhi = cP * cA + sP * sA;
391  const float_v sinPhi = -cP * sA + sP * cA;
392 
393  float_m ReturnMask(mask);
394  ReturnMask &= (!(CAMath::Abs(sinPhi) > maxSinPhi || CAMath::Abs(cosPhi) < 1.e-4f || CAMath::Abs(cP) < 1.e-4f));
395 
396  float_v tmp = alpha * 0.15915494f; // 1/(2.f*3.1415f);
397  ReturnMask &= abs(tmp - round(tmp)) < 0.167f; // allow turn by 60 degree only TODO scalar
398 
399  // float J[5][5] = { { j0, 0, 0, 0, 0 }, // Y
400  // { 0, 1, 0, 0, 0 }, // Z
401  // { 0, 0, j2, 0, 0 }, // SinPhi
402  // { 0, 0, 0, 1, 0 }, // DzDs
403  // { 0, 0, 0, 0, 1 } }; // Kappa
404 
405  const float_v j0 = cP / cosPhi;
406  const float_v j2 = cosPhi / cP;
407  const float_v d = SinPhi() - sP;
408 
409  SetX(x0 * cA + y0 * sA, ReturnMask);
410  SetY(-x0 * sA + y0 * cA, ReturnMask);
411 
412  SetSinPhi(sinPhi + j2 * d, ReturnMask);
413 
414  fC[0](ReturnMask) *= j0 * j0;
415  fC[1](ReturnMask) *= j0;
416  fC[3](ReturnMask) *= j0;
417  fC[6](ReturnMask) *= j0;
418  fC[10](ReturnMask) *= j0;
419 
420  fC[3](ReturnMask) *= j2;
421  fC[4](ReturnMask) *= j2;
422  fC[5](ReturnMask) *= j2 * j2;
423  fC[8](ReturnMask) *= j2;
424  fC[12](ReturnMask) *= j2;
425 
426  fAlpha(ReturnMask) += alpha;
427 
428  return ReturnMask;
429 }
void SetSinPhi(const float_v &v)
static T Sin(const T &x)
Definition: PndCAMath.h:83
static T Cos(const T &x)
Definition: PndCAMath.h:88
static T Abs(const T &x)
Definition: PndCAMath.h:68
void SetY(const float_v &v)
void SetX(const float_v &v)
float f
Definition: P4_F32vec4.h:32
double alpha
Definition: f_Init.h:31

◆ RotateXY()

void PndCATrackParamVector::RotateXY ( float_v  alpha,
float_v &  x,
float_v &  y,
float_v &  sin,
const float_m &  mask = float_m(true) 
) const
inline

Definition at line 431 of file PndCATrackParamVector.h.

References CAMath::Cos(), f, GetCosPhi(), CAMath::Sin(), sin(), SinPhi(), X(), and Y().

Referenced by SetErr2QPt().

432 {
433  //* Rotate the coordinate system in XY on the angle alpha
434  if ((abs(alpha) < 1e-6f || !mask).isFull())
435  return;
436 
437  const float_v cA = CAMath::Cos(alpha);
438  const float_v sA = CAMath::Sin(alpha);
439 
440  x(mask) = (X() * cA + Y() * sA);
441  y(mask) = (-X() * sA + Y() * cA);
442  sin(mask) = -GetCosPhi() * sA + SinPhi() * cA;
443 }
friend F32vec4 sin(const F32vec4 &a)
Definition: P4_F32vec4.h:130
static T Sin(const T &x)
Definition: PndCAMath.h:83
static T Cos(const T &x)
Definition: PndCAMath.h:88
float f
Definition: P4_F32vec4.h:32
double alpha
Definition: f_Init.h:31

◆ SetAngle() [1/2]

void PndCATrackParamVector::SetAngle ( const float_v &  v)
inline

Definition at line 213 of file PndCATrackParamVector.h.

References fAlpha, and v.

213 { fAlpha = v; }
__m128 v
Definition: P4_F32vec4.h:15

◆ SetAngle() [2/2]

void PndCATrackParamVector::SetAngle ( const float_v &  v,
const float_m &  m 
)
inline

Definition at line 214 of file PndCATrackParamVector.h.

References fAlpha, and v.

214 { fAlpha(m) = v; }
__m128 m
Definition: P4_F32vec4.h:38
__m128 v
Definition: P4_F32vec4.h:15

◆ SetChi2() [1/2]

void PndCATrackParamVector::SetChi2 ( const float_v &  v)
inline

Definition at line 207 of file PndCATrackParamVector.h.

References fChi2, and v.

207 { fChi2 = v; }
__m128 v
Definition: P4_F32vec4.h:15

◆ SetChi2() [2/2]

void PndCATrackParamVector::SetChi2 ( const float_v &  v,
const float_m &  m 
)
inline

Definition at line 208 of file PndCATrackParamVector.h.

References fChi2, and v.

208 { fChi2(m) = v; }
__m128 m
Definition: P4_F32vec4.h:38
__m128 v
Definition: P4_F32vec4.h:15

◆ SetCov() [1/2]

void PndCATrackParamVector::SetCov ( int  i,
const float_v &  v 
)
inline

Definition at line 189 of file PndCATrackParamVector.h.

References fC, i, and v.

189 { fC[i] = v; }
__m128 v
Definition: P4_F32vec4.h:15
unsigned int i
Definition: P4_F32vec4.h:33

◆ SetCov() [2/2]

void PndCATrackParamVector::SetCov ( int  i,
const float_v &  v,
const float_m &  m 
)
inline

Definition at line 190 of file PndCATrackParamVector.h.

References fC, i, and m.

190 { fC[i](m) = v; }
__m128 m
Definition: P4_F32vec4.h:38
__m128 v
Definition: P4_F32vec4.h:15
unsigned int i
Definition: P4_F32vec4.h:33

◆ SetDzDs() [1/2]

void PndCATrackParamVector::SetDzDs ( const float_v &  v)
inline

Definition at line 200 of file PndCATrackParamVector.h.

References fP, and v.

Referenced by InitDirection().

200 { fP[3] = v; }
__m128 v
Definition: P4_F32vec4.h:15

◆ SetDzDs() [2/2]

void PndCATrackParamVector::SetDzDs ( const float_v &  v,
const float_m &  m 
)
inline

Definition at line 201 of file PndCATrackParamVector.h.

References fP, and m.

201 { fP[3](m) = v; }
__m128 m
Definition: P4_F32vec4.h:38
__m128 v
Definition: P4_F32vec4.h:15

◆ SetErr2QPt()

◆ SetErr2Y()

void PndCATrackParamVector::SetErr2Y ( float_v  v)
inline

Definition at line 219 of file PndCATrackParamVector.h.

References fC, and v.

219 { fC[0] = v; }
__m128 v
Definition: P4_F32vec4.h:15

◆ SetErr2Z()

void PndCATrackParamVector::SetErr2Z ( float_v  v)
inline

Definition at line 220 of file PndCATrackParamVector.h.

References fC, and v.

220 { fC[2] = v; }
__m128 v
Definition: P4_F32vec4.h:15

◆ SetISec() [1/2]

void PndCATrackParamVector::SetISec ( const int_v &  v)
inline

Definition at line 216 of file PndCATrackParamVector.h.

References fISec, and v.

216 { fISec = v; }
__m128 v
Definition: P4_F32vec4.h:15

◆ SetISec() [2/2]

void PndCATrackParamVector::SetISec ( const int_v &  v,
const int_m &  m 
)
inline

Definition at line 217 of file PndCATrackParamVector.h.

References fISec, and v.

217 { fISec(m) = v; }
__m128 m
Definition: P4_F32vec4.h:38
__m128 v
Definition: P4_F32vec4.h:15

◆ SetNDF() [1/3]

void PndCATrackParamVector::SetNDF ( int  v)
inline

Definition at line 209 of file PndCATrackParamVector.h.

References fNDF, and v.

209 { fNDF = v; }
__m128 v
Definition: P4_F32vec4.h:15

◆ SetNDF() [2/3]

void PndCATrackParamVector::SetNDF ( const int_v &  v)
inline

Definition at line 210 of file PndCATrackParamVector.h.

References fNDF, and v.

210 { fNDF = v; }
__m128 v
Definition: P4_F32vec4.h:15

◆ SetNDF() [3/3]

void PndCATrackParamVector::SetNDF ( const int_v &  v,
const int_m &  m 
)
inline

Definition at line 211 of file PndCATrackParamVector.h.

References fNDF, and v.

211 { fNDF(m) = v; }
__m128 m
Definition: P4_F32vec4.h:38
__m128 v
Definition: P4_F32vec4.h:15

◆ SetPar() [1/2]

void PndCATrackParamVector::SetPar ( int  i,
const float_v &  v 
)
inline

Definition at line 187 of file PndCATrackParamVector.h.

References fP, i, and v.

187 { fP[i] = v; }
__m128 v
Definition: P4_F32vec4.h:15
unsigned int i
Definition: P4_F32vec4.h:33

◆ SetPar() [2/2]

void PndCATrackParamVector::SetPar ( int  i,
const float_v &  v,
const float_m &  m 
)
inline

Definition at line 188 of file PndCATrackParamVector.h.

References fP, i, and m.

188 { fP[i](m) = v; }
__m128 m
Definition: P4_F32vec4.h:38
__m128 v
Definition: P4_F32vec4.h:15
unsigned int i
Definition: P4_F32vec4.h:33

◆ SetQMomentum()

void PndCATrackParamVector::SetQMomentum ( const float_v &  v)
inline

Definition at line 203 of file PndCATrackParamVector.h.

References SetQPt().

203 { SetQPt(v); }
__m128 v
Definition: P4_F32vec4.h:15
void SetQPt(const float_v &v)

◆ SetQPt() [1/2]

void PndCATrackParamVector::SetQPt ( const float_v &  v)
inline

Definition at line 202 of file PndCATrackParamVector.h.

References fP, and v.

Referenced by SetQMomentum().

202 { fP[4] = v; }
__m128 v
Definition: P4_F32vec4.h:15

◆ SetQPt() [2/2]

void PndCATrackParamVector::SetQPt ( const float_v &  v,
const float_m &  m 
)
inline

Definition at line 204 of file PndCATrackParamVector.h.

References fP, and m.

204 { fP[4](m) = v; }
__m128 m
Definition: P4_F32vec4.h:38
__m128 v
Definition: P4_F32vec4.h:15

◆ SetSignCosPhi() [1/2]

void PndCATrackParamVector::SetSignCosPhi ( const float_v &  v)
inline

Definition at line 205 of file PndCATrackParamVector.h.

References fSignCosPhi, and v.

Referenced by InitDirection().

205 { fSignCosPhi = v; }
__m128 v
Definition: P4_F32vec4.h:15

◆ SetSignCosPhi() [2/2]

void PndCATrackParamVector::SetSignCosPhi ( const float_v &  v,
const float_m &  m 
)
inline

Definition at line 206 of file PndCATrackParamVector.h.

References fSignCosPhi, and v.

206 { fSignCosPhi(m) = v; }
__m128 m
Definition: P4_F32vec4.h:38
__m128 v
Definition: P4_F32vec4.h:15

◆ SetSinPhi() [1/2]

void PndCATrackParamVector::SetSinPhi ( const float_v &  v)
inline

Definition at line 198 of file PndCATrackParamVector.h.

References fP, and v.

Referenced by InitDirection(), and Rotate().

198 { fP[2] = v; }
__m128 v
Definition: P4_F32vec4.h:15

◆ SetSinPhi() [2/2]

void PndCATrackParamVector::SetSinPhi ( const float_v &  v,
const float_m &  m 
)
inline

Definition at line 199 of file PndCATrackParamVector.h.

References fP, and m.

199 { fP[2](m) = v; }
__m128 m
Definition: P4_F32vec4.h:38
__m128 v
Definition: P4_F32vec4.h:15

◆ SetTrackParam() [1/2]

void PndCATrackParamVector::SetTrackParam ( const PndCATrackParamVector param,
const float_m &  m = float_m(true) 
)
inline

Definition at line 160 of file PndCATrackParamVector.h.

References Cov(), fC, fChi2, fISec, fNDF, fP, fSignCosPhi, fX, GetChi2(), GetNDF(), i, m, Par(), SignCosPhi(), and X().

161  {
162  for (int i = 0; i < 5; i++)
163  fP[i](m) = param.Par()[i];
164  for (int i = 0; i < 15; i++)
165  fC[i](m) = param.Cov()[i];
166  fX(m) = param.X();
167  fSignCosPhi(m) = param.SignCosPhi();
168  fChi2(m) = param.GetChi2();
169  fNDF(static_cast<int_m>(m)) = param.GetNDF();
170  fISec(m) = param.fISec;
171  }
__m128 m
Definition: P4_F32vec4.h:38
unsigned int i
Definition: P4_F32vec4.h:33
const float_v & Cov(int i) const
const float_v & Par(int i) const

◆ SetTrackParam() [2/2]

void PndCATrackParamVector::SetTrackParam ( const PndCATrackParamVector p,
int  k,
int  m 
)
inline

Definition at line 173 of file PndCATrackParamVector.h.

References Cov(), f, fAlpha, fC, fChi2, fISec, fNDF, fP, fSignCosPhi, fX, i, m, and Par().

174  {
175  fX[k] = p.fX[m];
176  fSignCosPhi = float_v(1.f);
177  for (int i = 0; i < 5; i++)
178  fP[i][k] = p.Par()[i][m];
179  for (int i = 0; i < 15; i++)
180  fC[i][k] = p.Cov()[i][m];
181  fChi2[k] = p.fChi2[m];
182  fNDF[k] = p.fNDF[m];
183  fAlpha[k] = p.fAlpha[m];
184  fISec[k] = p.fISec[m];
185  }
__m128 m
Definition: P4_F32vec4.h:38
unsigned int i
Definition: P4_F32vec4.h:33
float f
Definition: P4_F32vec4.h:32
const float_v & Cov(int i) const
const float_v & Par(int i) const

◆ SetX() [1/2]

void PndCATrackParamVector::SetX ( const float_v &  v)
inline

Definition at line 192 of file PndCATrackParamVector.h.

References fX, and v.

Referenced by Rotate().

192 { fX = v; }
__m128 v
Definition: P4_F32vec4.h:15

◆ SetX() [2/2]

void PndCATrackParamVector::SetX ( const float_v &  v,
const float_m &  m 
)
inline

Definition at line 195 of file PndCATrackParamVector.h.

References fX, and v.

195 { fX(m) = v; }
__m128 m
Definition: P4_F32vec4.h:38
__m128 v
Definition: P4_F32vec4.h:15

◆ SetY() [1/2]

void PndCATrackParamVector::SetY ( const float_v &  v)
inline

Definition at line 193 of file PndCATrackParamVector.h.

References fP, and v.

Referenced by Rotate().

193 { fP[0] = v; }
__m128 v
Definition: P4_F32vec4.h:15

◆ SetY() [2/2]

void PndCATrackParamVector::SetY ( const float_v &  v,
const float_m &  m 
)
inline

Definition at line 196 of file PndCATrackParamVector.h.

References fP, and m.

196 { fP[0](m) = v; }
__m128 m
Definition: P4_F32vec4.h:38
__m128 v
Definition: P4_F32vec4.h:15

◆ SetZ() [1/2]

void PndCATrackParamVector::SetZ ( const float_v &  v)
inline

Definition at line 194 of file PndCATrackParamVector.h.

References fP, and v.

194 { fP[1] = v; }
__m128 v
Definition: P4_F32vec4.h:15

◆ SetZ() [2/2]

void PndCATrackParamVector::SetZ ( const float_v &  v,
const float_m &  m 
)
inline

Definition at line 197 of file PndCATrackParamVector.h.

References fP, and m.

197 { fP[1](m) = v; }
__m128 m
Definition: P4_F32vec4.h:38
__m128 v
Definition: P4_F32vec4.h:15

◆ SignCosPhi()

float_v PndCATrackParamVector::SignCosPhi ( ) const
inline

The sign of cos phi is always positive in the slice tracker. Only after coordinate transformation can the sign change to negative.

Definition at line 119 of file PndCATrackParamVector.h.

References fSignCosPhi.

Referenced by PndCATrackLinearisationVector::PndCATrackLinearisationVector(), SetTrackParam(), Tx1(), and Tx2().

119 { return fSignCosPhi; }

◆ SinPhi()

float_v PndCATrackParamVector::SinPhi ( ) const
inline

◆ Transport() [1/3]

float_m PndCATrackParamVector::Transport ( const int_v &  ista,
const PndCAParam param,
const float_m &  mask = float_m(true) 
)

Referenced by SetErr2QPt().

◆ Transport() [2/3]

float_m PndCATrackParamVector::Transport ( const PndCAHitV hit,
const PndCAParam p,
const float_m &  mask = float_m(true) 
)

◆ Transport() [3/3]

float_m PndCATrackParamVector::Transport ( const PndCAHit hit,
const PndCAParam p,
const float_m &  mask = float_m(true) 
)

◆ Transport0() [1/3]

float_m PndCATrackParamVector::Transport0 ( const int_v &  ista,
const PndCAParam param,
const float_m &  mask = float_m(true) 
)
inline

Definition at line 562 of file PndCATrackParamVector.h.

References PndCAParam::cBz(), PndCAParam::GetR(), and Transport0ToX().

Referenced by SetErr2QPt().

563 {
564  float_m active = mask;
565  // PndCATrackLinearisationVector tE( *this );
566  // active &= TransportToX( param.GetR( ista, active ), tE, param.cBz(), 0.999f, nullptr, active );
567  active &= Transport0ToX(param.GetR(ista, active), param.cBz(), active);
568  return active;
569 }
float_m Transport0ToX(const float_v &x, const float_v &Bz, const float_m &mask)
float GetR(short iSt) const
Definition: PndCAParam.h:78
float cBz() const
Definition: PndCAParam.h:60

◆ Transport0() [2/3]

float_m PndCATrackParamVector::Transport0 ( const PndCAHitV hit,
const PndCAParam param,
const float_m &  mask = float_m(true) 
)

◆ Transport0() [3/3]

float_m PndCATrackParamVector::Transport0 ( const PndCAHit hit,
const PndCAParam p,
const float_m &  mask = float_m(true) 
)
inline

Definition at line 571 of file PndCATrackParamVector.h.

References CAMath::Abs(), PndCAHit::Angle(), PndCAParam::cBz(), f, fAlpha, fX, Rotate(), Transport0ToX(), and PndCAHit::X0().

572 {
573  if (((CAMath::Abs(-fAlpha + hit.Angle()) < 1e-6f && CAMath::Abs(-fX + hit.X0()) < 2e-4f) || !mask).isFull())
574  return mask;
575  float_m active = mask;
576  // PndCATrackLinearisationVector tR( *this );
577  active &= Rotate(-fAlpha + hit.Angle(), .999f, active);
578  // active &= Rotate( -fAlpha + hit.Angle(), tR, .999f, active );
579  // active &= TransportToX( hit.X0(), tR, param.cBz(), 0.999f, nullptr, active );
580  active &= Transport0ToX(hit.X0(), param.cBz(), active);
581  return active;
582 }
float_m Rotate(const float_v &alpha, PndCATrackLinearisationVector &t0, const float maxSinPhi=.999f, const float_m &mask=float_m(true))
float_m Transport0ToX(const float_v &x, const float_v &Bz, const float_m &mask)
float Angle() const
Definition: PndCAHits.h:64
static T Abs(const T &x)
Definition: PndCAMath.h:68
float f
Definition: P4_F32vec4.h:32
float X0() const
Definition: PndCAHits.h:47

◆ Transport0ToX()

float_m PndCATrackParamVector::Transport0ToX ( const float_v &  x,
const float_v &  Bz,
const float_m &  mask 
)

Referenced by SetErr2QPt(), and Transport0().

◆ TransportJ0ToX0()

float_m PndCATrackParamVector::TransportJ0ToX0 ( const float_v &  x0,
const float_v &  cBz,
float_v &  extrDy,
float_v &  extrDz,
float_v  J[],
const float_m &  active = float_m(true) 
)
inline

Definition at line 533 of file PndCATrackParamVector.h.

References CAMath::Abs(), DzDs(), CAMath::RSqrt(), SinPhi(), and X().

Referenced by SetErr2QPt().

534 {
535  const float_v &ey = SinPhi();
536  const float_v &dx = x - X();
537  const float_v &exi = CAMath::RSqrt(float_v(Vc::One) - ey * ey); // RSqrt
538 
539  const float_v &dxBz = dx * cBz;
540  const float_v &dS = dx * exi;
541  const float_v &h2 = dS * exi * exi;
542  const float_v &h4 = .5f * h2 * dxBz;
543 
544  float_m mask = active && CAMath::Abs(exi) <= 1.e4f;
545 
546  extrDy(mask) = dS * ey;
547  extrDz(mask) = dS * DzDs();
548  J[0](mask) = dS;
549  J[1](mask) = 0;
550  J[2](mask) = h4;
551  J[3](mask) = dS;
552  J[4](mask) = dS;
553  J[5](mask) = 0;
554  return mask;
555 }
static T Abs(const T &x)
Definition: PndCAMath.h:68
static T RSqrt(const T &x)
Definition: PndCAMath.h:62

◆ TransportToX() [1/3]

float_m PndCATrackParamVector::TransportToX ( const float_v &  x,
const float_v &  Bz,
const float  maxSinPhi = .999f,
const float_m &  mask = float_m(true) 
)

Referenced by SetErr2QPt().

◆ TransportToX() [2/3]

float_m PndCATrackParamVector::TransportToX ( const float_v &  x,
PndCATrackLinearisationVector t0,
const float_v &  Bz,
const float  maxSinPhi = .999f,
float_v *  DL = 0,
const float_m &  mask = float_m(true) 
)

◆ TransportToX() [3/3]

float_m PndCATrackParamVector::TransportToX ( const float_v &  x,
const float_v &  sinPhi0,
const float_v &  Bz,
const float_v  maxSinPhi = .999f,
const float_m &  mask = float_m(true) 
)
inline

mvz start 23.01.2010

mvz end 23.01.2010

Definition at line 307 of file PndCATrackParamVector.h.

References CAMath::Abs(), DzDs(), f, fC, fP, fX, QPt(), CAMath::RSqrt(), SinPhi(), and X().

308 {
309  //* Transport the track parameters to X=x, using linearization at phi0 with 0 curvature,
310  //* and the field value Bz
311  //* maxSinPhi is the max. allowed value for |t0.SinPhi()|
312  //* linearisation of trajectory t0 is also transported to X=x,
313  //* returns 1 if OK
314  //*
315 
316  // debugKF() << "Start TransportToX(" << x << ", " << _mask << ")\n" << *this << std::endl;
317 
318  const float_v &ey = sinPhi0;
319  const float_v &dx = x - X();
320  const float_v &exi = float_v(Vc::One) * CAMath::RSqrt(float_v(Vc::One) - ey * ey); // RSqrt
321 
322  const float_v &dxBz = dx * Bz;
323  const float_v &dS = dx * exi;
324  const float_v &h2 = dS * exi * exi;
325  const float_v &h4 = .5f * h2 * dxBz;
327  // const float_v &sinPhi = SinPhi() * (float_v( Vc::One ) - 0.5f * dxBz * QPt() *dxBz * QPt()/ ( float_v( Vc::One ) - SinPhi()*SinPhi() )) + dxBz * QPt();
328  const float_v &sinPhi = SinPhi() + dxBz * QPt();
330 
331  float_m mask = _mask && CAMath::Abs(exi) <= 1.e4f;
332  mask &= ((CAMath::Abs(sinPhi) <= maxSinPhi) || (maxSinPhi <= 0.f));
333 
334  fX(mask) += dx;
335  fP[0](mask) += dS * ey + h2 * (SinPhi() - ey) + h4 * QPt();
336  fP[1](mask) += dS * DzDs();
337  fP[2](mask) = sinPhi;
338 
339  // const float_v c00 = fC[0];
340  // const float_v c11 = fC[2];
341  const float_v c20 = fC[3];
342  // const float_v c21 = fC[4];
343  const float_v c22 = fC[5];
344  // const float_v c30 = fC[6];
345  const float_v c31 = fC[7];
346  // const float_v c32 = fC[8];
347  const float_v c33 = fC[9];
348  const float_v c40 = fC[10];
349  // const float_v c41 = fC[11];
350  const float_v c42 = fC[12];
351  // const float_v c43 = fC[13];
352  const float_v c44 = fC[14];
353 
354  const float_v two(2.f);
355 
356  fC[0](mask) += h2 * h2 * c22 + h4 * h4 * c44 + two * (h2 * c20 + h4 * (c40 + h2 * c42));
357 
358  // fC[1] ( mask ) += h2 * c21 + h4 * c41 + dS * ( c30 + h2 * c32 + h4 * c43 );
359  fC[2](mask) += dS * (two * c31 + dS * c33);
360 
361  fC[3](mask) += h2 * c22 + h4 * c42 + dxBz * (c40 + h2 * c42 + h4 * c44);
362  // fC[4] ( mask ) += dS * c32 + dxBz * ( c41 + dS * c43 );
363  const float_v &dxBz_c44 = dxBz * c44;
364  fC[12](mask) += dxBz_c44;
365  fC[5](mask) += dxBz * (two * c42 + dxBz_c44);
366 
367  // fC[6] ( mask ) += h2 * c32 + h4 * c43;
368  fC[7](mask) += dS * c33;
369  // fC[8] ( mask ) += dxBz * c43;
370  // fC[9] ( mask ) = c33;
371 
372  fC[10](mask) += h2 * c42 + h4 * c44;
373  // fC[11]( mask ) += dS * c43;
374  // fC[13]( mask ) = c43;
375  // fC[14]( mask ) = c44;
376 
377  // debugKF() << mask << "\n" << *this << std::endl;
378  return mask;
379 }
static T Abs(const T &x)
Definition: PndCAMath.h:68
static T RSqrt(const T &x)
Definition: PndCAMath.h:62
float f
Definition: P4_F32vec4.h:32

◆ TransportToXWithMaterial() [1/3]

float_m PndCATrackParamVector::TransportToXWithMaterial ( const float_v &  x,
const float_v &  XOverX0,
const float_v &  XThimesRho,
const float_v &  Bz,
const float  maxSinPhi = .999f 
)

Referenced by SetErr2QPt().

◆ TransportToXWithMaterial() [2/3]

float_m PndCATrackParamVector::TransportToXWithMaterial ( const float_v &  x,
PndCATrackLinearisationVector t0,
PndCATrackFitParam par,
const float_v &  XOverX0,
const float_v &  XThimesRho,
const float_v &  Bz,
const float  maxSinPhi = .999f,
const float_m &  mask = float_m(true) 
)

◆ TransportToXWithMaterial() [3/3]

float_m PndCATrackParamVector::TransportToXWithMaterial ( const float_v &  x,
PndCATrackFitParam par,
const float_v &  XOverX0,
const float_v &  XThimesRho,
const float_v &  Bz,
const float  maxSinPhi = .999f 
)

◆ Tx1()

float_v PndCATrackParamVector::Tx1 ( ) const
inline

Definition at line 111 of file PndCATrackParamVector.h.

References SignCosPhi(), SinPhi(), and sqrt().

111 { return SinPhi() / (SignCosPhi() * sqrt(1 - SinPhi() * SinPhi())); } // CHECKME
friend F32vec4 sqrt(const F32vec4 &a)
Definition: P4_F32vec4.h:40

◆ Tx2()

float_v PndCATrackParamVector::Tx2 ( ) const
inline

Definition at line 112 of file PndCATrackParamVector.h.

References DzDs(), SignCosPhi(), SinPhi(), and sqrt().

112 { return DzDs() / (SignCosPhi() * sqrt(1 - SinPhi() * SinPhi())); } // dx2/dx0 = dz/dx
friend F32vec4 sqrt(const F32vec4 &a)
Definition: P4_F32vec4.h:40

◆ X()

float_v PndCATrackParamVector::X ( ) const
inline

Definition at line 98 of file PndCATrackParamVector.h.

References fX.

Referenced by Rotate(), RotateXY(), SetTrackParam(), TransportJ0ToX0(), TransportToX(), and X0().

98 { return fX; }

◆ X0()

float_v PndCATrackParamVector::X0 ( ) const
inline

Definition at line 108 of file PndCATrackParamVector.h.

References X().

108 { return X(); }

◆ X1()

float_v PndCATrackParamVector::X1 ( ) const
inline

Definition at line 109 of file PndCATrackParamVector.h.

References Y().

109 { return Y(); }

◆ X2()

float_v PndCATrackParamVector::X2 ( ) const
inline

Definition at line 110 of file PndCATrackParamVector.h.

References Z().

110 { return Z(); }

◆ Y()

float_v PndCATrackParamVector::Y ( ) const
inline

Definition at line 99 of file PndCATrackParamVector.h.

References fP.

Referenced by FilterVtx(), Rotate(), RotateXY(), and X1().

99 { return fP[0]; }

◆ Z()

float_v PndCATrackParamVector::Z ( ) const
inline

Definition at line 100 of file PndCATrackParamVector.h.

References fP.

Referenced by FilterVtx(), and X2().

100 { return fP[1]; }

Friends And Related Function Documentation

◆ PndCATrackParam

friend class PndCATrackParam
friend

Definition at line 153 of file PndCATrackParamVector.h.

Member Data Documentation

◆ fAlpha

float_v PndCATrackParamVector::fAlpha

Definition at line 301 of file PndCATrackParamVector.h.

Referenced by Angle(), Rotate(), SetAngle(), SetTrackParam(), and Transport0().

◆ fC

◆ fChi2

float_v PndCATrackParamVector::fChi2

Definition at line 298 of file PndCATrackParamVector.h.

Referenced by Chi2(), FilterVtx(), GetChi2(), Reset(), SetChi2(), and SetTrackParam().

◆ fISec

int_v PndCATrackParamVector::fISec

Definition at line 302 of file PndCATrackParamVector.h.

Referenced by ISec(), SetISec(), and SetTrackParam().

◆ fNDF

int_v PndCATrackParamVector::fNDF

Definition at line 299 of file PndCATrackParamVector.h.

Referenced by GetNDF(), NDF(), Reset(), SetNDF(), and SetTrackParam().

◆ fP

◆ fSignCosPhi

float_v PndCATrackParamVector::fSignCosPhi

◆ fX

float_v PndCATrackParamVector::fX

Definition at line 294 of file PndCATrackParamVector.h.

Referenced by GetX(), Reset(), SetTrackParam(), SetX(), Transport0(), TransportToX(), and X().


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