Belle II Software  release-05-01-25
Helix.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2013 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Martin Heck *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 #include <framework/dataobjects/Helix.h>
11 
12 #include <framework/gearbox/Const.h>
13 
14 #include <boost/math/special_functions/sign.hpp>
15 #include <boost/math/special_functions/sinc.hpp>
16 #include <boost/math/tools/precision.hpp>
17 
18 #include <TMatrixD.h>
19 
20 #include <cassert>
21 
22 using namespace Belle2;
23 using namespace HelixParameterIndex;
24 
26  m_d0(0.0),
27  m_phi0(0.0),
28  m_omega(0.0),
29  m_z0(0.0),
30  m_tanLambda(0.0)
31 {
32 }
33 
34 Helix::Helix(const TVector3& position,
35  const TVector3& momentum,
36  const short int charge,
37  const double bZ)
38 {
39  setCartesian(position, momentum, charge, bZ);
40 }
41 
42 Helix::Helix(const double& d0,
43  const double& phi0,
44  const double& omega,
45  const double& z0,
46  const double& tanLambda)
47  : m_d0(d0),
48  m_phi0(phi0),
49  m_omega(omega),
50  m_z0(z0),
51  m_tanLambda(tanLambda)
52 {
53 }
54 
55 double Helix::getPerigeeX() const
56 {
57  return getD0() * getSinPhi0();
58 }
59 
60 double Helix::getPerigeeY() const
61 {
62  return -getD0() * getCosPhi0();
63 }
64 
65 double Helix::getPerigeeZ() const
66 {
67  return getZ0();
68 }
69 
70 TVector3 Helix::getPerigee() const
71 {
72  return TVector3(getPerigeeX(), getPerigeeY(), getPerigeeZ());
73 }
74 
75 double Helix::getMomentumX(const double bZ) const
76 {
77  return getCosPhi0() * getTransverseMomentum(bZ);
78 }
79 
80 double Helix::getMomentumY(const double bZ) const
81 {
82  return getSinPhi0() * getTransverseMomentum(bZ);
83 }
84 
85 double Helix::getMomentumZ(const double bZ) const
86 {
87  return getTanLambda() * getTransverseMomentum(bZ);
88 }
89 
90 TVector3 Helix:: getMomentum(const double bZ) const
91 {
92  return TVector3(getMomentumX(bZ), getMomentumY(bZ), getMomentumZ(bZ));
93 }
94 
95 TVector3 Helix::getDirection() const
96 {
97  return TVector3(getCosPhi0(), getSinPhi0(), getTanLambda());
98 }
99 
100 double Helix::getTransverseMomentum(const double bZ) const
101 {
102  return 1 / std::fabs(getAlpha(bZ) * getOmega());
103 }
104 
105 double Helix::getKappa(const double bZ) const
106 {
107  return getOmega() * getAlpha(bZ);
108 }
109 
110 double Helix::getAlpha(const double bZ)
111 {
112  return 1.0 / (bZ * Const::speedOfLight) * 1E4;
113 }
114 
115 short Helix::getChargeSign() const
116 {
117  return boost::math::sign(getOmega());
118 }
119 
120 void Helix::reverse()
121 {
122  // All except z0 have to be taken to their opposites
123  m_d0 = -m_d0; //d0
124  m_phi0 = reversePhi(m_phi0);
125  m_omega = -m_omega;
126  m_tanLambda = -m_tanLambda;
127 }
128 
129 double Helix::getArcLength2DAtCylindricalR(const double& cylindricalR) const
130 {
131  // Slight trick here
132  // Since the sought point is on the helix we treat it as the perigee
133  // and the origin as the point to extrapolate to.
134  // We know the distance of the origin to the circle, which is just d0
135  // The direct distance from the origin to the imaginary perigee is just the given cylindricalR.
136  const double dr = getD0();
137  const double deltaCylindricalR = cylindricalR;
138  const double absArcLength2D = calcArcLength2DAtDeltaCylindricalRAndDr(deltaCylindricalR, dr);
139  return absArcLength2D;
140 }
141 
142 double Helix::getArcLength2DAtXY(const double& x, const double& y) const
143 {
144  double dr = 0;
145  double arcLength2D = 0;
146  calcArcLength2DAndDrAtXY(x, y, arcLength2D, dr);
147  return arcLength2D;
148 }
149 
150 double Helix::getArcLength2DAtNormalPlane(const double& byX, const double& byY,
151  const double& nX, const double& nY) const
152 {
153  // Construct the tangential vector to the plan in xy space
154  const double tX = nY;
155  const double tY = -nX;
156 
157  // Fetch the helix parameters
158  const double omega = getOmega();
159  const double cosPhi0 = getCosPhi0();
160  const double sinPhi0 = getSinPhi0();
161  const double d0 = getD0();
162 
163  // Prepare a delta vector, which is the vector from the perigee point to the support point of the plane
164  // Split it in component parallel and a component orthogonal to tangent at the perigee.
165  const double deltaParallel = byX * cosPhi0 + byY * sinPhi0;
166  const double deltaOrthogonal = byY * cosPhi0 - byX * sinPhi0 + d0;
167  const double deltaCylindricalR = hypot(deltaOrthogonal, deltaParallel);
168  const double deltaCylindricalRSquared = deltaCylindricalR * deltaCylindricalR;
169 
170  const double UOrthogonal = 1 + omega * deltaOrthogonal; // called nu in the Karimaki paper.
171  const double UParallel = omega * deltaParallel;
172 
173  // Some commonly used terms - compare Karimaki 1990
174  const double A = 2 * deltaOrthogonal + omega * deltaCylindricalRSquared;
175 
176  const double tParallel = tX * cosPhi0 + tY * sinPhi0;
177  const double tOrthogonal = tY * cosPhi0 - tX * sinPhi0;
178  const double tCylindricalR = (tX * tX + tY * tY);
179 
180  const double c = A / 2;
181  const double b = UParallel * tParallel + UOrthogonal * tOrthogonal;
182  const double a = omega / 2 * tCylindricalR;
183 
184  const double discriminant = ((double)b) * b - 4 * a * c;
185  const double root = sqrt(discriminant);
186  const double bigSum = (b > 0) ? -b - root : -b + root;
187 
188  const double bigOffset = bigSum / 2 / a;
189  const double smallOffset = 2 * c / bigSum;
190 
191  const double distance1 = hypot(byX + bigOffset * tX, byY + bigOffset * tY);
192  const double distance2 = hypot(byX + smallOffset * tX, byY + smallOffset * tY);
193 
194  if (distance1 < distance2) {
195  return getArcLength2DAtXY(byX + bigOffset * tX, byY + bigOffset * tY);
196  } else {
197  return getArcLength2DAtXY(byX + smallOffset * tX, byY + smallOffset * tY);
198  }
199 }
200 
201 
202 TVector3 Helix::getPositionAtArcLength2D(const double& arcLength2D) const
203 {
204  /*
205  / \ / \ / \
206  | x | | cos phi0 -sin phi0 | | - sin(chi) / omega |
207  | | = | | * | |
208  | y | | sin phi0 cos phi0 | | -(1 - cos(chi)) / omega - d0 |
209  \ / \ / \ /
210 
211  and
212 
213  z = tanLambda * arclength + z0;
214 
215  where chi = -arcLength2D * omega
216 
217  Here arcLength2D means the arc length of the circle in the xy projection
218  traversed in the forward direction.
219  */
220 
221  // First calculating the position assuming the circle center lies on the y axes (phi0 = 0)
222  // Rotate to the right phi position afterwards
223  // Using the sinus cardinalis yields expressions that are smooth in the limit of omega -> 0
224 
225  // Do calculations in double
226  const double chi = -arcLength2D * getOmega();
227  const double chiHalf = chi / 2.0;
228 
229  using boost::math::sinc_pi;
230 
231  const double x = arcLength2D * sinc_pi(chi);
232  const double y = arcLength2D * sinc_pi(chiHalf) * sin(chiHalf) - getD0();
233 
234  // const double z = s * tan lambda + z0
235  const double z = fma((double)arcLength2D, getTanLambda(), getZ0());
236 
237  // Unrotated position
238  TVector3 position(x, y, z);
239 
240  // Rotate to the right phi0 position
241  position.RotateZ(getPhi0());
242 
243  return position;
244 }
245 
246 TVector3 Helix::getTangentialAtArcLength2D(const double& arcLength2D) const
247 {
248  const double omega = getOmega();
249  const double phi0 = getPhi0();
250  const double tanLambda = getTanLambda();
251 
252  const double chi = - omega * arcLength2D;
253 
254  const double tx = cos(chi + phi0);
255  const double ty = sin(chi + phi0);
256  const double tz = tanLambda;
257 
258  TVector3 tangential(tx, ty, tz);
259  return tangential;
260 }
261 
262 TVector3 Helix::getUnitTangentialAtArcLength2D(const double& arcLength2D) const
263 {
264  TVector3 unitTangential = getTangentialAtArcLength2D(arcLength2D);
265  const double norm = hypot(1, getTanLambda());
266  const double invNorm = 1 / norm;
267  unitTangential *= invNorm;
268  return unitTangential;
269 }
270 
271 TVector3 Helix::getMomentumAtArcLength2D(const double& arcLength2D, const double& bz) const
272 {
273  TVector3 momentum = getTangentialAtArcLength2D(arcLength2D);
274  const double pr = getTransverseMomentum(bz);
275  momentum *= pr;
276 
277  return momentum;
278 }
279 
280 double Helix::passiveMoveBy(const double& byX,
281  const double& byY,
282  const double& byZ)
283 {
284  // First calculate the distance of the new origin to the helix in the xy projection
285  double new_d0 = 0;
286  double arcLength2D = 0;
287  calcArcLength2DAndDrAtXY(byX, byY, arcLength2D, new_d0);
288 
289  // Third the new phi0 and z0 can be calculated from the arc length
290  double chi = -arcLength2D * getOmega();
291  double new_phi0 = m_phi0 + chi;
292  double new_z0 = getZ0() - byZ + getTanLambda() * arcLength2D;
293 
295  m_d0 = new_d0;
296  m_phi0 = new_phi0;
297  m_z0 = new_z0;
298 
299  return arcLength2D;
300 }
301 
302 TMatrixD Helix::calcPassiveMoveByJacobian(const TVector3& by, const double expandBelowChi) const
303 {
304  TMatrixD jacobian(5, 5);
305  calcPassiveMoveByJacobian(by.X(), by.Y(), jacobian, expandBelowChi);
306  return jacobian;
307 }
308 
309 
310 void Helix::calcPassiveMoveByJacobian(const double& byX,
311  const double& byY,
312  TMatrixD& jacobian,
313  const double expandBelowChi) const
314 {
315  // 0. Preparations
316  // Initialise the return value to a unit matrix
317  jacobian.UnitMatrix();
318  assert(jacobian.GetNrows() == jacobian.GetNcols());
319  assert(jacobian.GetNrows() == 5 or jacobian.GetNrows() == 6);
320 
321  // Fetch the helix parameters
322  const double omega = getOmega();
323  const double cosPhi0 = getCosPhi0();
324  const double sinPhi0 = getSinPhi0();
325  const double d0 = getD0();
326  const double tanLambda = getTanLambda();
327 
328  // Prepare a delta vector, which is the vector from the perigee point to the new origin
329  // Split it in component parallel and a component orthogonal to tangent at the perigee.
330  const double deltaParallel = byX * cosPhi0 + byY * sinPhi0;
331  const double deltaOrthogonal = byY * cosPhi0 - byX * sinPhi0 + d0;
332  const double deltaCylindricalR = hypot(deltaOrthogonal, deltaParallel);
333  const double deltaCylindricalRSquared = deltaCylindricalR * deltaCylindricalR;
334 
335  // Some commonly used terms - compare Karimaki 1990
336  const double A = 2 * deltaOrthogonal + omega * deltaCylindricalRSquared;
337  const double USquared = 1 + omega * A;
338  const double U = sqrt(USquared);
339  const double UOrthogonal = 1 + omega * deltaOrthogonal; // called nu in the Karimaki paper.
340  const double UParallel = omega * deltaParallel;
341  // Note U is a vector pointing from the middle of the projected circle scaled by a factor omega.
342  const double u = 1 + omega * d0; // just called u in the Karimaki paper.
343 
344  // ---------------------------------------------------------------------------------------
345  // 1. Set the parts related to the xy coordinates
346  // Fills the upper left 3x3 corner of the jacobain
347 
348  // a. Calculate the row related to omega
349  // The row related to omega is not different from the unit jacobian initialised above.
350  // jacobian(iOmega, iOmega) = 1.0; // From UnitMatrix above.
351 
352  // b. Calculate the row related to d0
353  const double new_d0 = A / (1 + U);
354 
355  jacobian(iD0, iOmega) = (deltaCylindricalR + new_d0) * (deltaCylindricalR - new_d0) / 2 / U;
356  jacobian(iD0, iPhi0) = - u * deltaParallel / U;
357  jacobian(iD0, iD0) = UOrthogonal / U;
358 
359  // c. Calculate the row related to phi0
360  // Also calculate the derivatives of the arc length
361  // which are need for the row related to z.
362  const double dArcLength2D_dD0 = - UParallel / USquared;
363  const double dArcLength2D_dPhi0 = (omega * deltaCylindricalRSquared + deltaOrthogonal - d0 * UOrthogonal) / USquared;
364 
365  jacobian(iPhi0, iD0) = - dArcLength2D_dD0 * omega;
366  jacobian(iPhi0, iPhi0) = u * UOrthogonal / USquared;
367  jacobian(iPhi0, iOmega) = -deltaParallel / USquared;
368 
369  // For jacobian(iPhi0, iOmega) we have to dig deeper
370  // since the normal equations have a divergence for omega -> 0.
371  // This hinders not only the straight line case but extrapolations between points which are close together,
372  // like it happens when the current perigee is very close the new perigee.
373  // To have a smooth transition in this limit we have to carefully
374  // factor out the divergent quotients and approximate them with their taylor series.
375 
376  const double chi = -atan2(UParallel, UOrthogonal);
377  double arcLength2D = 0;
378  double dArcLength2D_dOmega = 0;
379 
380  if (fabs(chi) < std::min(expandBelowChi, M_PI / 2.0)) {
381  // Never expand for the far side of the circle by limiting the expandBelow to maximally half a pi.
382  // Close side of the circle
383  double principleArcLength2D = deltaParallel / UOrthogonal;
384  const double dPrincipleArcLength2D_dOmega = - principleArcLength2D * deltaOrthogonal / UOrthogonal;
385 
386  const double x = principleArcLength2D * omega;
387  const double f = calcATanXDividedByX(x);
388  const double df_dx = calcDerivativeOfATanXDividedByX(x);
389 
390  arcLength2D = principleArcLength2D * f;
391  dArcLength2D_dOmega = (f + x * df_dx) * dPrincipleArcLength2D_dOmega + principleArcLength2D * df_dx * principleArcLength2D;
392  } else {
393  // Far side of the circle
394  // If the far side of the circle is a well definied concept, omega is high enough that
395  // we can divide by it.
396  // Otherwise nothing can rescue us since the far side of the circle is so far away that no reasonable extrapolation can be made.
397  arcLength2D = - chi / omega;
398  dArcLength2D_dOmega = (-arcLength2D - jacobian(iPhi0, iOmega)) / omega;
399  }
400 
401  // ---------------------------------------------------------------------------------------
402  // 2. Set the parts related to the z coordinate
403  // Since tanLambda stays the same there are no entries for tanLambda in the jacobian matrix
404  // For the new z0' = z0 + arcLength2D * tanLambda
405  jacobian(iZ0, iD0) = tanLambda * dArcLength2D_dD0;
406  jacobian(iZ0, iPhi0) = tanLambda * dArcLength2D_dPhi0;
407  jacobian(iZ0, iOmega) = tanLambda * dArcLength2D_dOmega;
408  jacobian(iZ0, iZ0) = 1.0; // From UnitMatrix above.
409  jacobian(iZ0, iTanLambda) = arcLength2D;
410 
411  if (jacobian.GetNrows() == 6) {
412  // Also write the derivates of arcLength2D to the jacobian if the matrix size has space for them.
413  jacobian(iArcLength2D, iD0) = dArcLength2D_dD0;
414  jacobian(iArcLength2D, iPhi0) = dArcLength2D_dPhi0;
415  jacobian(iArcLength2D, iOmega) = dArcLength2D_dOmega;
416  jacobian(iArcLength2D, iZ0) = 0.0;
417  jacobian(iArcLength2D, iTanLambda) = 0;
418  jacobian(iArcLength2D, iArcLength2D) = 1.0;
419  }
420 }
421 
422 double Helix::reversePhi(const double& phi)
423 {
424  return std::remainder(phi + M_PI, 2 * M_PI);
425 }
426 
427 double Helix::calcArcLength2DFromSecantLength(const double& secantLength) const
428 {
429  return secantLength * calcSecantLengthToArcLength2DFactor(secantLength);
430 }
431 
432 
433 double Helix::calcSecantLengthToArcLength2DFactor(const double& secantLength) const
434 {
435  double x = secantLength * getOmega() / 2.0;
436  return calcASinXDividedByX(x);
437 }
438 
439 double Helix::calcASinXDividedByX(const double& x)
440 {
441  // Approximation of asin(x) / x
442  // Inspired by BOOST's sinc
443 
444  BOOST_MATH_STD_USING;
445 
446  auto const taylor_n_bound = boost::math::tools::forth_root_epsilon<double>();
447 
448  if (abs(x) >= taylor_n_bound) {
449  if (fabs(x) == 1) {
450  return M_PI / 2.0;
451 
452  } else {
453  return asin(x) / x;
454 
455  }
456 
457  } else {
458  // approximation by taylor series in x at 0 up to order 0
459  double result = 1.0;
460 
461  auto const taylor_0_bound = boost::math::tools::epsilon<double>();
462  if (abs(x) >= taylor_0_bound) {
463  double x2 = x * x;
464  // approximation by taylor series in x at 0 up to order 2
465  result += x2 / 6.0;
466 
467  auto const taylor_2_bound = boost::math::tools::root_epsilon<double>();
468  if (abs(x) >= taylor_2_bound) {
469  // approximation by taylor series in x at 0 up to order 4
470  result += x2 * x2 * (3.0 / 40.0);
471  }
472  }
473  return result;
474  }
475 
476 }
477 
478 double Helix::calcATanXDividedByX(const double& x)
479 {
480  // Approximation of atan(x) / x
481  // Inspired by BOOST's sinc
482 
483  BOOST_MATH_STD_USING;
484 
485  auto const taylor_n_bound = boost::math::tools::forth_root_epsilon<double>();
486 
487  if (abs(x) >= taylor_n_bound) {
488  return atan(x) / x;
489 
490  } else {
491  // approximation by taylor series in x at 0 up to order 0
492  double result = 1.0;
493 
494  auto const taylor_0_bound = boost::math::tools::epsilon<double>();
495  if (abs(x) >= taylor_0_bound) {
496  double x2 = x * x;
497  // approximation by taylor series in x at 0 up to order 2
498  result -= x2 / 3.0;
499 
500  auto const taylor_2_bound = boost::math::tools::root_epsilon<double>();
501  if (abs(x) >= taylor_2_bound) {
502  // approximation by taylor series in x at 0 up to order 4
503  result += x2 * x2 * (1.0 / 5.0);
504  }
505  }
506  return result;
507  }
508 }
509 
510 double Helix::calcDerivativeOfATanXDividedByX(const double& x)
511 {
512  // Approximation of atan(x) / x
513  // Inspired by BOOST's sinc
514 
515  BOOST_MATH_STD_USING;
516 
517  auto const taylor_n_bound = boost::math::tools::forth_root_epsilon<double>();
518 
519  const double x2 = x * x;
520  if (abs(x) >= taylor_n_bound) {
521  const double chi = atan(x);
522  return ((1 - chi / x) / x - chi) / (1 + x2);
523 
524  } else {
525  // approximation by taylor series in x at 0 up to order 0
526  double result = 1.0;
527 
528  auto const taylor_0_bound = boost::math::tools::epsilon<double>();
529  if (abs(x) >= taylor_0_bound) {
530  // approximation by taylor series in x at 0 up to order 2
531  result -= 2.0 * x / 3.0;
532 
533  auto const taylor_2_bound = boost::math::tools::root_epsilon<double>();
534  if (abs(x) >= taylor_2_bound) {
535  // approximation by taylor series in x at 0 up to order 4
536  result += x2 * x * (4.0 / 5.0);
537  }
538  }
539  return result;
540  }
541 }
542 
543 
544 
545 
546 void Helix::calcArcLength2DAndDrAtXY(const double& x, const double& y, double& arcLength2D, double& dr) const
547 {
548  // Prepare common variables
549  const double omega = getOmega();
550  const double cosPhi0 = getCosPhi0();
551  const double sinPhi0 = getSinPhi0();
552  const double d0 = getD0();
553 
554  const double deltaParallel = x * cosPhi0 + y * sinPhi0;
555  const double deltaOrthogonal = y * cosPhi0 - x * sinPhi0 + d0;
556  const double deltaCylindricalR = hypot(deltaOrthogonal, deltaParallel);
557  const double deltaCylindricalRSquared = deltaCylindricalR * deltaCylindricalR;
558 
559  const double A = 2 * deltaOrthogonal + omega * deltaCylindricalRSquared;
560  const double U = sqrt(1 + omega * A);
561  const double UOrthogonal = 1 + omega * deltaOrthogonal; // called nu in the Karimaki paper.
562  const double UParallel = omega * deltaParallel;
563  // Note U is a vector pointing from the middle of the projected circle scaled by a factor omega.
564 
565  // Calculate dr
566  dr = A / (1 + U);
567 
568  // Calculate the absolute value of the arc length
569  const double chi = -atan2(UParallel, UOrthogonal);
570 
571  if (fabs(chi) < M_PI / 8) { // Rough guess where the critical zone for approximations begins
572  // Close side of the circle
573  double principleArcLength2D = deltaParallel / UOrthogonal;
574  arcLength2D = principleArcLength2D * calcATanXDividedByX(principleArcLength2D * omega);
575  } else {
576  // Far side of the circle
577  // If the far side of the circle is a well definied concept meaning that we have big enough omega.
578  arcLength2D = -chi / omega;
579  }
580 }
581 
582 double Helix::calcArcLength2DAtDeltaCylindricalRAndDr(const double& deltaCylindricalR, const double& dr) const
583 {
584  const double omega = getOmega();
585  double secantLength = sqrt((deltaCylindricalR + dr) * (deltaCylindricalR - dr) / (1 + dr * omega));
586  return calcArcLength2DFromSecantLength(secantLength);
587 }
588 
589 void Helix::setCartesian(const TVector3& position,
590  const TVector3& momentum,
591  const short int charge,
592  const double bZ)
593 {
594  assert(abs(charge) <= 1); // Not prepared for doubly-charged particles.
595  const double alpha = getAlpha(bZ);
596 
597  // We allow for the case that position, momentum are not given
598  // exactly in the perigee. Therefore we have to transform the momentum
599  // with the position as the reference point and then move the coordinate system
600  // to the origin.
601 
602  const double x = position.X();
603  const double y = position.Y();
604  const double z = position.Z();
605 
606  const double px = momentum.X();
607  const double py = momentum.Y();
608  const double pz = momentum.Z();
609 
610  const double ptinv = 1 / hypot(px, py);
611  const double omega = charge * ptinv / alpha;
612  const double tanLambda = ptinv * pz;
613  const double phi0 = atan2(py, px);
614  const double z0 = z;
615  const double d0 = 0;
616 
617  m_omega = omega;
618  m_phi0 = phi0;
619  m_d0 = d0;
620  m_tanLambda = tanLambda;
621  m_z0 = z0;
622 
623  passiveMoveBy(-x, -y, 0);
624 }
625 
626 namespace Belle2 {
632  std::ostream& operator<<(std::ostream& output, const Helix& helix)
633  {
634  return output
635  << "Helix("
636  << "d0=" << helix.getD0() << ", "
637  << "phi0=" << helix.getPhi0() << ", "
638  << "omega=" << helix.getOmega() << ", "
639  << "z0=" << helix.getZ0() << ", "
640  << "tanLambda=" << helix.getTanLambda() << ")";
641  }
643 }
Belle2::operator<<
std::ostream & operator<<(std::ostream &output, const IntervalOfValidity &iov)
Definition: IntervalOfValidity.cc:196
Belle2::CDC::Helix::phi0
double phi0(void) const
Return helix parameter phi0.
Definition: Helix.h:391
Belle2::CDC::Helix::momentum
Hep3Vector momentum(double dPhi=0.) const
returns momentum vector after rotating angle dPhi in phi direction.
Definition: Helix.cc:261
Belle2::CDC::Helix::x
HepPoint3D x(double dPhi=0.) const
returns position after rotating angle dPhi in phi direction.
Definition: Helix.cc:199
Belle2::Const::speedOfLight
static const double speedOfLight
[cm/ns]
Definition: Const.h:568
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::CDC::Helix
Helix parameter class.
Definition: Helix.h:51
Belle2::CDC::Helix::a
const HepVector & a(void) const
Returns helix parameters.
Definition: Helix.h:431
Belle2::CDC::Helix::sinPhi0
double sinPhi0(void) const
Return sin phi0.
Definition: Helix.h:495
Belle2::CDC::Helix::Helix
Helix(const HepPoint3D &pivot, const HepVector &a, const HepSymMatrix &Ea)
Constructor with pivot, helix parameter a, and its error matrix.
Definition: Helix.cc:134
Belle2::CDC::Helix::cosPhi0
double cosPhi0(void) const
Return cos phi0.
Definition: Helix.h:503
Belle2::CDC::Helix::dr
double dr(void) const
Return helix parameter dr.
Definition: Helix.h:383