Heavy ion Analysis Libriares
Loading...
Searching...
No Matches
HelixBase.cxx
1/*
2 * HalHelixBase.cxx
3 *
4 * Created on: 7 wrz 2018
5 * Author: Daniel Wielanek
6 * E-mail: daniel.wielanek@gmail.com
7 * Warsaw University of Technology, Faculty of Physics
8 */
9#include "HelixBase.h"
10
11#include <RtypesCore.h>
12#include <TMath.h>
13#include <TMathBase.h>
14#include <TString.h>
15#include <TVector2.h>
16#include <cmath>
17#include <iostream>
18
19
20namespace Hal {
21 Double_t HelixBase::fgHelixBz = 0.5;
22
23 HelixBase::~HelixBase() {}
24
25 HelixBase::HelixBase() :
26 fX0(0),
27 fY0(0),
28 fZ0(0),
29 fPhi0(0),
30 fPhiCos(0),
31 fPhiSin(0),
32 fLambda(0),
33 fLambdaCos(0),
34 fLambdaSin(0),
35 fCurv(0),
36 fXc(0),
37 fYc(0),
38 fH(0) {}
39
40 HelixBase::HelixBase(const HelixBase& t) :
41 TObject(t),
42 fX0(t.fX0),
43 fY0(t.fY0),
44 fZ0(t.fZ0),
45 fPhi0(t.fPhi0),
46 fPhiCos(t.fPhiCos),
47 fPhiSin(t.fPhiSin),
48 fLambda(t.fLambda),
49 fLambdaCos(t.fLambdaCos),
50 fLambdaSin(t.fLambdaSin),
51 fCurv(t.fCurv),
52 fXc(t.fXc),
53 fYc(t.fYc),
54 fH(t.fH) {}
55
56 HelixBase::HelixBase(const TVector3& x, const TVector3& mom, Double_t charge, Double_t conversion) : HelixBase() {
57 Double_t pt = mom.Pt();
58 if (TMath::Abs(conversion) < 0.00000001) conversion = -1.0 / 0.299792458 / GetBz();
59 fH = -TMath::Sign(1.0, GetBz() * charge);
60 fCurv = charge / (conversion * pt) * 0.01; // C [cm]
61 if (conversion * pt == 0) fCurv = 0; // singularity !
62 fCurv = TMath::Abs(fCurv);
63 fLambda = TMath::ATan2(mom.Pz(), pt);
64 fLambdaCos = TMath::Cos(fLambda);
65 fLambdaSin = TMath::Sin(fLambda);
66
67 Double_t Psi = mom.Phi();
68 fPhi0 = TVector2::Phi_mpi_pi(Psi - fH * TMath::PiOver2());
69 fPhiSin = TMath::Sin(fPhi0);
70 fPhiCos = TMath::Cos(fPhi0);
71 fXc = x.X() - fPhiCos / fCurv;
72 fYc = x.Y() - fPhiSin / fCurv;
73 //
74
75
76 fX0 = x.X(); // x0
77 fY0 = x.Y(); // y0
78 fZ0 = x.Z(); // z0
79 }
80
81 Int_t HelixBase::BaseIntersection(Double_t R, TVector3& x1, TVector3& x2) const {
82 Double_t s1, s2;
83 BasePathLength(R, s1, s2);
84 if (s1 == s2) {
85 x1 = BasePosition(s1);
86 return 1;
87 } else if (s1 == MaxPath()) {
88 return 0;
89 };
90 if (s1 < s2) {
91 x1 = BasePosition(s1);
92 x2 = BasePosition(s2);
93 } else {
94 x2 = BasePosition(s1);
95 x1 = BasePosition(s2);
96 }
97 return 2;
98 }
99
100 Double_t HelixBase::GetPeriod() const {
101 if (fCurv == 0) {
102 return 1E+9;
103 } else { //
104 return TMath::Abs(TMath::TwoPi() / (fH * GetCurv() * fLambdaCos));
105 }
106 }
107
108 void HelixBase::BasePrint() const {
109 std::cout << "*********** HALHELIX ***********" << std::endl;
110 std::cout << Form("X0 %4.4f %4.4f %4.4f", fX0, fY0, fZ0) << std::endl;
111 std::cout << Form("Xc = %4.4f Yc=%4.4f L=%4.4f ", fXc, fYc, fLambda) << std::endl;
112 std::cout << Form("Phi0=%4.2f Curv=%4.4f", fPhi0 / TMath::DegToRad(), fCurv) << std::endl;
113 }
114
115 void HelixBase::BaseSetParams(const TVector3& x, const TVector3& mom, Double_t charge, Double_t conversion) {
116 Double_t pt = mom.Pt();
117 if (TMath::Abs(conversion) < 0.00000001) conversion = -1 / 0.299792458 / GetBz();
118 fH = -TMath::Sign(1.0, GetBz() * charge);
119 fCurv = charge / (conversion * pt); // C
120 if (conversion * pt == 0) fCurv = 0; // singularity !
121 fCurv = TMath::Abs(fCurv * 0.01); // goto cm
122 fLambda = TMath::ATan2(mom.Pz(), pt);
123 fLambdaCos = TMath::Cos(fLambda);
124 fLambdaSin = TMath::Sin(fLambda);
125 Double_t Psi = mom.Phi();
126 fPhi0 = Psi - fH * TMath::PiOver2();
127 if (fPhi0 < -TMath::Pi()) fPhi0 += TMath::TwoPi();
128 if (fPhi0 > TMath::Pi()) fPhi0 -= TMath::TwoPi();
129 fPhiSin = TMath::Sin(fPhi0);
130 fPhiCos = TMath::Cos(fPhi0);
131 fXc = x.X() - fPhiCos / fCurv;
132 fYc = x.Y() - fPhiSin / fCurv;
133 //
134 fX0 = x.X(); // x0
135 fY0 = x.Y(); // y0
136 fZ0 = x.Z(); // z0
137 }
138
139 void HelixBase::BaseSetParams(Double_t x,
140 Double_t y,
141 Double_t z,
142 Double_t px,
143 Double_t py,
144 Double_t pz,
145 Double_t charge,
146 Double_t conversion) {
147 Double_t pt = TMath::Sqrt(px * px + py * py);
148 if (TMath::Abs(conversion) < 0.00000001) conversion = -1 / 0.299792458 / GetBz();
149 fH = -TMath::Sign(1.0, GetBz() * charge);
150 fCurv = charge / (conversion * pt); // C
151 if (conversion * pt == 0) fCurv = 0; // singularity !
152 fCurv = TMath::Abs(fCurv * 0.01); // goto cm
153 fLambda = TMath::ATan2(pz, pt);
154 fLambdaCos = TMath::Cos(fLambda);
155 fLambdaSin = TMath::Sin(fLambda);
156
157 Double_t Psi = 0;
158 if (pt != 0 && pz != 0) Psi = TMath::ATan2(py, px);
159 fPhi0 = Psi - fH * TMath::PiOver2();
160 if (fPhi0 < -TMath::Pi()) fPhi0 += TMath::TwoPi();
161 if (fPhi0 > TMath::Pi()) fPhi0 -= TMath::TwoPi();
162 fPhiSin = TMath::Sin(fPhi0);
163 fPhiCos = TMath::Cos(fPhi0);
164 fXc = x - fPhiCos / fCurv;
165 fYc = y - fPhiSin / fCurv;
166 //
167 fX0 = x; // x0
168 fY0 = y; // y0
169 fZ0 = z; // z0
170 }
171
172 void HelixBase::SetMagField(Double_t Bz) { fgHelixBz = Bz; }
173
174 Double_t HelixBase::GetBz() { return fgHelixBz; }
175
176 HelixBase& HelixBase::operator=(const HelixBase& helix) {
177 if (this != &helix) {
178 fY0 = helix.fY0;
179 fZ0 = helix.fZ0;
180 fX0 = helix.fX0;
181 fPhi0 = helix.fPhi0;
182 fLambda = helix.fLambda;
183 fCurv = helix.fCurv;
184 fXc = helix.fXc;
185 fYc = helix.fYc;
186 fH = helix.fH;
187 }
188 return *this;
189 }
190
191 TVector3 HelixBase::BasePosition(Double_t s) const {
192 if (fCurv == 0) { return TVector3(fX0 - s * fLambdaCos * fPhiSin, fY0 + s * fLambdaCos * fPhiCos, fZ0 + s * fLambdaSin); }
193 Double_t overC = 1.0 / fCurv;
194 return TVector3(fX0 + (TMath::Cos(fPhi0 + s * fH * fCurv * fLambdaCos) - fPhiCos) * overC,
195 fY0 + (TMath::Sin(fPhi0 + s * fH * fCurv * fLambdaCos) - fPhiSin) * overC,
196 fZ0 + s * fLambdaSin);
197 }
198
199 TVector3 HelixBase::BaseGetStartPoint() const { return TVector3(fX0, fY0, fZ0); }
200
201 Double_t HelixBase::FudgePathLenght(const TVector3& vec) const {
202 Double_t s;
203 Double_t dx = vec.X() - fX0;
204 Double_t dy = vec.Y() - fY0;
205
206 if (fCurv == 0) {
207 s = (dy * fPhiCos - dx * fPhiSin) / fLambdaCos;
208 } else {
209 s = TMath::ATan2(dy * fPhiCos - dx * fPhiSin, 1. / fCurv + dx * fPhiCos + dy * fPhiSin) / (fH * fCurv * fLambdaCos);
210 }
211 return s;
212 }
213
214 void HelixBase::BasePathLength(Double_t r, Double_t& s1, Double_t& s2) const {
215 //
216 // The math is taken from Maple with C(expr,optimized) and
217 // some hand-editing. It is not very nice but efficient.
218 // 'first' is the smallest of the two solutions (may be negative)
219 // 'second' is the other.
220 //
221 if (fCurv == 0) {
222 Double_t t1 = fLambdaCos * (fX0 * fPhiSin - fY0 * fPhiCos);
223 Double_t t12 = fY0 * fY0;
224 Double_t t13 = fPhiCos * fPhiCos;
225 Double_t t15 = r * r;
226 Double_t t16 = fX0 * fX0;
227 Double_t t20 = -fLambdaCos * fLambdaCos * (2.0 * fX0 * fPhiSin * fY0 * fPhiCos + t12 - t12 * t13 - t15 + t13 * t16);
228 if (t20 < 0.) {
229 s1 = s2 = 1E+9;
230 return;
231 }
232 t20 = TMath::Sqrt(t20);
233 s1 = (t1 - t20) / (fLambdaCos * fLambdaCos);
234 s2 = (t1 + t20) / (fLambdaCos * fLambdaCos);
235 } else {
236 Double_t t1 = fY0 * fCurv;
237 Double_t t2 = fPhiSin;
238 Double_t t3 = fCurv * fCurv;
239 Double_t t4 = fY0 * t2;
240 Double_t t5 = fPhiCos;
241 Double_t t6 = fX0 * t5;
242 Double_t t8 = fX0 * fX0;
243 Double_t t11 = fY0 * fY0;
244 Double_t t14 = r * r;
245 Double_t t15 = t14 * fCurv;
246 Double_t t17 = t8 * t8;
247 Double_t t19 = t11 * t11;
248 Double_t t21 = t11 * t3;
249 Double_t t23 = t5 * t5;
250 Double_t t32 = t14 * t14;
251 Double_t t35 = t14 * t3;
252 Double_t t38 = 8.0 * t4 * t6 - 4.0 * t1 * t2 * t8 - 4.0 * t11 * fCurv * t6 + 4.0 * t15 * t6 + t17 * t3 + t19 * t3
253 + 2.0 * t21 * t8 + 4.0 * t8 * t23 - 4.0 * t8 * fX0 * fCurv * t5 - 4.0 * t11 * t23
254 - 4.0 * t11 * fY0 * fCurv * t2 + 4.0 * t11 - 4.0 * t14 + t32 * t3 + 4.0 * t15 * t4 - 2.0 * t35 * t11
255 - 2.0 * t35 * t8;
256 Double_t t40 = (-t3 * t38);
257 if (t40 < 0.) {
258 s1 = s2 = 1E+9;
259 return;
260 }
261 t40 = TMath::Sqrt(t40);
262
263 Double_t t43 = fX0 * fCurv;
264 Double_t t45 = 2.0 * t5 - t35 + t21 + 2.0 - 2.0 * t1 * t2 - 2.0 * t43 - 2.0 * t43 * t5 + t8 * t3;
265 Double_t t46 = fH * fLambdaCos * fCurv;
266
267 s1 = (-fPhi0 + 2.0 * TMath::ATan((-2.0 * t1 + 2.0 * t2 + t40) / t45)) / t46;
268 s2 = -(fPhi0 + 2.0 * TMath::ATan((2.0 * t1 - 2.0 * t2 + t40) / t45)) / t46;
269
270 //
271 // Solution can be off by +/- one period, select smallest
272 //
273 Double_t p = GetPeriod();
274 if (!TMath::IsNaN(s1)) {
275 if (TMath::Abs(s1 - p) < TMath::Abs(s1))
276 s1 = s1 - p;
277 else if (TMath::Abs(s1 + p) < TMath::Abs(s1))
278 s1 = s1 + p;
279 }
280 if (!TMath::IsNaN(s2)) {
281 if (TMath::Abs(s2 - p) < TMath::Abs(s2))
282 s2 = s2 - p;
283 else if (TMath::Abs(s2 + p) < TMath::Abs(s2))
284 s2 = s2 + p;
285 }
286 }
287 if (s1 > s2) {
288 Double_t temp = s1;
289 s1 = s2;
290 s2 = temp;
291 }
292 }
293
294 void HelixBase::BasePathLength(Double_t r, Double_t x, Double_t y, Double_t& s1, Double_t& s2) const {
295 HelixBase helix(*this);
296 helix.fX0 = fX0 - x;
297 helix.fY0 = fY0 - y;
298 helix.BasePathLength(r, s1, s2);
299 helix.fX0 = fX0 + x;
300 helix.fY0 = fY0 + y;
301 }
302
303 Double_t HelixBase::BasePathLength(const TVector3& point, Bool_t scanPeriods) const {
304 Double_t dx = point.X() - fX0;
305 Double_t dy = point.Y() - fY0;
306 Double_t dz = point.Z() - fZ0;
307 Double_t s = 0;
308 Double_t Phase = TMath::ATan2(fYc - fY0, fXc - fX0) + TMath::Pi();
309 if (Phase < -TMath::Pi()) { Phase += TMath::TwoPi(); }
310 if (Phase > TMath::Pi()) { Phase -= TMath::TwoPi(); }
311 Double_t CosPhase = TMath::Cos(Phase);
312 Double_t SinPhase = TMath::Sin(Phase);
313 Double_t period = 0;
314 // Bool_t Singularity = kFALSE;
315
316 if (fCurv == 0) {
317 s = (dy * CosPhase - dx * SinPhase) / fLambdaCos;
318 period = 1E+9;
319 } else { //
320 const Double_t overC = 1.0 / fCurv;
321 s = atan2(dy * CosPhase - dx * SinPhase, overC + dx * CosPhase + dy * SinPhase) / (fH * fCurv * fLambdaCos);
322 period = TMath::Abs(TMath::TwoPi() / (fH * fCurv * fLambdaCos));
323
324 // using namespace units;
325 // #endif
326 // const Double_t MaxPrecisionNeeded = 0.000001; // micrometer;
327 const Double_t MaxPrecisionNeeded = 0.000001;
328 const int MaxIterations = 100;
329
330 //
331 // The math is taken from Maple with C(expr,optimized) and
332 // some hand-editing. It is not very nice but efficient.
333 //
334 Double_t t34 = fCurv * fLambdaCos * fLambdaCos;
335 Double_t t41 = fLambdaSin * fLambdaSin;
336 Double_t t6, t7, t11, t12, t19;
337
338 //
339 // Get a first guess by using the dca in 2D. Since
340 // in some extreme cases we might be off by n periods
341 // we add (subtract) periods in case we get any closer.
342 //
343
344 if (scanPeriods) {
345 Double_t ds = period;
346 int j, jmin = 0;
347
348 Double_t d, dmin = (BasePosition(s) - point).Mag2();
349 for (j = 1; j < MaxIterations; j++) {
350 d = (BasePosition(s + j * ds) - point).Mag2();
351 if (d < dmin) {
352 dmin = d;
353 jmin = j;
354 } else
355 break;
356 }
357 for (j = -1; - j < MaxIterations; j--) {
358 d = (BasePosition(s + j * ds) - point).Mag2();
359 if (d < dmin) {
360 dmin = d;
361 jmin = j;
362 } else
363 break;
364 }
365 if (jmin) s += jmin * ds;
366 }
367
368 //
369 // Newtons method:
370 // Stops after MaxIterations iterations or if the required
371 // precision is obtained. Whatever comes first.
372 //
373 Double_t sOld = s;
374
375 for (int i = 0; i < MaxIterations; i++) {
376 t6 = Phase + s * fH * fCurv * fLambdaCos;
377 t7 = cos(t6);
378 t11 = dx - overC * (t7 - CosPhase);
379 t12 = sin(t6);
380 t19 = dy - overC * (t12 - SinPhase);
381 s -= (t11 * t12 * fH * fLambdaCos - t19 * t7 * fH * fLambdaCos - (dz - s * fLambdaSin) * fLambdaSin)
382 / (t12 * t12 * fLambdaCos * fLambdaCos + t11 * t7 * t34 + t7 * t7 * fLambdaCos * fLambdaCos + t19 * t12 * t34 + t41);
383 if (TMath::Abs(sOld - s) < MaxPrecisionNeeded) break;
384 sOld = s;
385 }
386 }
387 return s;
388 }
389
390 Double_t HelixBase::BasePathLength(Double_t x, Double_t y) const { return FudgePathLenght(TVector3(x, y, 0)); }
391
392 void HelixBase::PathLengths(const HelixBase& h, Double_t& s1, Double_t& s2) const {
393 if ((fCurv == 0 && h.fCurv != 0) || ((fCurv != 0 && h.fCurv == 0))) {
394 s1 = s2 = 1E+9;
395 return;
396 }
397
398
399 if (fCurv == 0) {
400 //
401 // Analytic solution
402 //
403 TVector3 dv = h.BaseGetStartPoint() - BaseGetStartPoint();
404 TVector3 a(-fLambdaCos * fPhiSin, fLambdaCos * fPhiCos, fLambdaSin);
405 TVector3 b(-h.fLambdaCos * h.fPhiSin, h.fLambdaCos * h.fPhiCos, h.fLambdaSin);
406 Double_t ab = a * b;
407 Double_t g = dv * a;
408 Double_t k = dv * b;
409 s2 = (k - ab * g) / (ab * ab - 1.);
410 s1 = g + s2 * ab;
411 return;
412 } else {
413 //
414 // First step: get dca in the xy-plane as start value
415 //
416 Double_t dx = h.BaseGetXcenter() - BaseGetXcenter();
417 Double_t dy = h.BaseGetYcenter() - BaseGetYcenter();
418 Double_t dd = ::sqrt(dx * dx + dy * dy);
419 Double_t r1 = 1. / fCurv;
420 Double_t r2 = 1. / h.GetCurv();
421
422 Double_t cosAlpha = (r1 * r1 + dd * dd - r2 * r2) / (2 * r1 * dd);
423
424 Double_t s;
425 Double_t xx, yy;
426 if (TMath::Abs(cosAlpha) < 1) { // two solutions
427 const Double_t overDD = 1.0 / dd;
428 Double_t sinAlpha = TMath::Sin(TMath::ACos(cosAlpha));
429 xx = BaseGetXcenter() + r1 * (cosAlpha * dx - sinAlpha * dy) * overDD;
430 yy = BaseGetYcenter() + r1 * (sinAlpha * dx + cosAlpha * dy) * overDD;
431 s = BasePathLength(xx, yy);
432 xx = BaseGetXcenter() + r1 * (cosAlpha * dx + sinAlpha * dy) * overDD;
433 yy = BaseGetYcenter() + r1 * (cosAlpha * dy - sinAlpha * dx) * overDD;
434 Double_t a = BasePathLength(xx, yy);
435 if (h.Distance(BasePosition(a)) < h.Distance(BasePosition(s))) s = a;
436 } else { // no intersection (or exactly one)
437 int rsign = ((r2 - r1) > dd ? -1 : 1); // set -1 when *this* helix is
438 // completely contained in the other
439 xx = BaseGetXcenter() + rsign * r1 * dx / dd;
440 yy = BaseGetYcenter() + rsign * r1 * dy / dd;
441 s = BasePathLength(xx, yy);
442 }
443
444 //
445 // Second step: scan in decreasing intervals around seed 's'
446 //
447 // const Double_t MinStepSize = 0.001; // 10*micrometer;
448 // const Double_t MinRange = 10.; // 10*centimeter;
449 const Double_t MinStepSize = 10 * 0.000001;
450 const Double_t MinRange = 10 * 0.01;
451 Double_t dmin = h.Distance(BasePosition(s));
452 Double_t range = TMath::Max(2 * dmin, MinRange);
453 Double_t ds = range / 10;
454 Double_t slast = -999999, ss, d;
455 s1 = s - range / 2.;
456 s2 = s + range / 2.;
457
458 while (ds > MinStepSize) {
459 for (ss = s1; ss < s2 + ds; ss += ds) {
460 d = h.Distance(BasePosition(ss));
461 if (d < dmin) {
462 dmin = d;
463 s = ss;
464 }
465 slast = ss;
466 }
467 //
468 // In the rare cases where the minimum is at the
469 // the border of the current range we shift the range
470 // and start all over, i.e we do not decrease 'ds'.
471 // Else we decrease the search intervall around the
472 // current minimum and redo the scan in smaller steps.
473 //
474 if (s == s1) {
475 d = 0.8 * (s2 - s1);
476 s1 -= d;
477 s2 -= d;
478 } else if (s == slast) {
479 d = 0.8 * (s2 - s1);
480 s1 += d;
481 s2 += d;
482 } else {
483 s1 = s - ds;
484 s2 = s + ds;
485 ds /= 10;
486 }
487 }
488 s1 = s;
489 s2 = h.BasePathLength(BasePosition(s));
490 }
491 }
492
493 Double_t HelixBase::Distance(const TVector3& p, Bool_t scanPeriods) const {
494 return (BasePosition(BasePathLength(p, scanPeriods)) - p).Mag();
495 }
496
497 HelixBase::HelixBase(Double_t x,
498 Double_t y,
499 Double_t z,
500 Double_t px,
501 Double_t py,
502 Double_t pz,
503 Double_t charge,
504 Double_t conversion) :
505 HelixBase(TVector3(x, y, z), TVector3(px, py, pz), charge, conversion) {}
506
507 void HelixBase::BaseShift(Double_t x, Double_t y, Double_t z) {
508 fXc += x;
509 fYc += y;
510 fX0 += x;
511 fY0 += y;
512 fZ0 += z;
513 }
514
515 Double_t HelixBase::GetRotationDirection() const { return TMath::Sign(1, GetH() * TMath::Cos(GetDipAngle())); }
516
517 TVector3 HelixBase::BaseMomentum(Double_t s) const {
518 TVector3 origin = BasePosition(s);
519 if (fCurv != 0) {
520 Double_t newPhase = TMath::ATan2(origin.Y() - fYc, origin.X() - fXc);
521 const Double_t conversion = -1.0 / 0.299792458 / GetBz();
522 const Double_t charge = -TMath::Sign(1.0, GetBz() * fH);
523 Double_t pt = TMath::Abs(charge / (conversion * fCurv) * 0.01);
524 // set origin mOroigin
525 // set phase
526 if (newPhase < -TMath::Pi()) newPhase += TMath::TwoPi();
527 if (newPhase > TMath::Pi()) newPhase -= TMath::TwoPi();
528 Double_t shiftPhase = newPhase + fH * TMath::PiOver2();
529 return TVector3(pt * TMath::Cos(shiftPhase), pt * TMath::Sin(shiftPhase), pt * TMath::Tan(fLambda));
530 }
531 return TVector3(0, 0, 0);
532 }
533
534 Double_t HelixBase::GetPathAbsMin(Double_t s1, Double_t s2) const {
535 const Double_t period = GetPeriod();
536 if (period <= 0) return TMath::Min(s1, s2);
537 while (s1 < 0)
538 s1 += period;
539 while (s2 < 0)
540 s2 += period;
541 return TMath::Min(s1, s2);
542 }
543
544 void HelixBase::BaseFullEval(Double_t s, TVector3& mom, TVector3& pos) const {
545 pos = BasePosition(s);
546 if (fCurv != 0) {
547 Double_t newPhase = TMath::ATan2(pos.Y() - fYc, pos.X() - fXc);
548 const Double_t conversion = -1.0 / 0.299792458 / GetBz();
549 const Double_t charge = -TMath::Sign(1.0, GetBz() * fH);
550 Double_t pt = TMath::Abs(charge / (conversion * fCurv) * 0.01);
551 if (newPhase < -TMath::Pi()) newPhase += TMath::TwoPi();
552 if (newPhase > TMath::Pi()) newPhase -= TMath::TwoPi();
553 Double_t shiftPhase = newPhase + fH * TMath::PiOver2();
554 mom.SetXYZ(pt * TMath::Cos(shiftPhase), pt * TMath::Sin(shiftPhase), pt * TMath::Tan(fLambda));
555 return;
556 }
557 mom.SetXYZ(0, 0, 0);
558 }
559} // namespace Hal
Hal::HelixBase::MaxPath
static Double_t MaxPath()
Definition HelixBase.h:245
Hal::HelixBase::BaseGetYcenter
Double_t BaseGetYcenter() const
Definition HelixBase.h:97
Hal::HelixBase::GetCurv
Double_t GetCurv() const
Definition HelixBase.h:190
Hal::HelixBase::GetRotationDirection
Double_t GetRotationDirection() const
Definition HelixBase.cxx:515
Hal::HelixBase::GetPeriod
Double_t GetPeriod() const
Definition HelixBase.cxx:100
Hal::HelixBase::BaseShift
void BaseShift(Double_t x, Double_t y, Double_t z)
Definition HelixBase.cxx:507
Hal::HelixBase::GetH
Double_t GetH() const
Definition HelixBase.h:205
Hal::HelixBase::BaseIntersection
Int_t BaseIntersection(Double_t R, TVector3 &x1, TVector3 &x2) const
Definition HelixBase.cxx:81
Hal::HelixBase::GetPathAbsMin
Double_t GetPathAbsMin(Double_t s1, Double_t s2) const
Definition HelixBase.cxx:534
Hal::HelixBase::GetDipAngle
Double_t GetDipAngle() const
Definition HelixBase.h:200
Hal::HelixBase::BasePosition
TVector3 BasePosition(Double_t s) const
Definition HelixBase.cxx:191
Hal::HelixBase::operator=
HelixBase & operator=(const HelixBase &helix)
Definition HelixBase.cxx:176
Hal::HelixBase::BaseMomentum
TVector3 BaseMomentum(Double_t s) const
Definition HelixBase.cxx:517
Hal::HelixBase::BaseSetParams
void BaseSetParams(const TVector3 &x, const TVector3 &mom, Double_t charge, Double_t conversion=0.)
Definition HelixBase.cxx:115
Hal::HelixBase::BaseGetXcenter
Double_t BaseGetXcenter() const
Definition HelixBase.h:92
Hal::HelixBase::SetMagField
static void SetMagField(Double_t Bz)
Definition HelixBase.cxx:172
Hal::HelixBase::PathLengths
void PathLengths(const HelixBase &h, Double_t &s1, Double_t &s2) const
Definition HelixBase.cxx:392
Hal::HelixBase::BasePathLength
void BasePathLength(Double_t r, Double_t x, Double_t y, Double_t &s1, Double_t &s2) const
Definition HelixBase.cxx:294
Hal::HelixBase::BaseGetStartPoint
TVector3 BaseGetStartPoint() const
Definition HelixBase.cxx:199
Hal::HelixBase::BaseFullEval
void BaseFullEval(Double_t s, TVector3 &mom, TVector3 &pos) const
Definition HelixBase.cxx:544