PandaRoot
PndCATrackParam.h
Go to the documentation of this file.
1 //****************************************************************************
2 //* This file is part of PandaRoot. *
3 //* *
4 //* PandaRoot is distributed under the terms of the *
5 //* GNU General Public License (GPL) version 3, *
6 //* copied verbatim in the file "LICENSE". *
7 //* *
8 //* Copyright (C) 2006 - 2024 FAIR GmbH and copyright holders of PandaRoot *
9 //* The copyright holders are listed in the file "COPYRIGHTHOLDERS". *
10 //* The authors are listed in the file "AUTHORS". *
11 //****************************************************************************
12 
13 //-*- Mode: C++ -*-
14 // $Id: PndCATrackParam.h,v 1.5 2011/05/20 16:28:05 fisyak Exp $
15 // ************************************************************************
16 // This file is property of and copyright by the ALICE HLT Project *
17 // ALICE Experiment at CERN, All rights reserved. *
18 // See cxx source for full Copyright notice *
19 // *
20 //*************************************************************************
21 
22 #ifndef PNDCATRACKPARAM_H
23 #define PNDCATRACKPARAM_H
24 #include <string.h>
25 #include "PndCADef.h"
26 #include "PndCAMath.h"
27 #include "PndCATrackParamVector.h"
28 #include "PndCAHits.h"
29 
31 
40  // friend std::istream &operator>>( std::istream &, PndCATrackParam & );
41  // friend std::ostream &operator<<( std::ostream &, const PndCATrackParam & );
42  public:
44  PndCATrackParam(const TrackParamVector &v, int i) : fX(v.X()[i]), fSignCosPhi(v.SignCosPhi()[i]), fChi2(v.Chi2()[i]), fNDF(v.NDF()[i]), fAlpha(v.Angle()[i]), fISec(v.ISec()[i])
45  {
46  for (int j = 0; j < 5; ++j)
47  fP[j] = v.Par()[j][i];
48  for (int j = 0; j < 15; ++j)
49  fC[j] = v.Cov()[j][i];
50  }
51 
53  float fBethe;
54  float fE;
55  float fTheta2;
56  float fEP2;
57  float fSigmadE2;
58  float fK22;
59  float fK33;
60  float fK43;
61  float fK44;
62  };
63 
64  float X() const { return fX; }
65  float Y() const { return fP[0]; }
66  float Z() const { return fP[1]; }
67  float SinPhi() const { return fP[2]; }
68  float DzDs() const { return fP[3]; }
69  float QPt() const { return fP[4]; }
70  float QMomentum() const { return fP[4]; }
71 
76  float SignCosPhi() const { return fSignCosPhi; }
77  float Chi2() const { return fChi2; }
78  int NDF() const { return fNDF; }
79 
80  float Err2Y() const { return fC[0]; }
81  float Err2Z() const { return fC[2]; }
82  float Err2SinPhi() const { return fC[5]; }
83  float Err2DzDs() const { return fC[9]; }
84  float Err2QPt() const { return fC[14]; }
85 
86  float Angle() const { return fAlpha; }
87  int ISec() const { return fISec; }
88 
89  float GetX() const { return fX; }
90  float GetY() const { return fP[0]; }
91  float GetZ() const { return fP[1]; }
92  float GetSinPhi() const { return fP[2]; }
93  float GetDzDs() const { return fP[3]; }
94  float GetQPt() const { return fP[4]; }
95  float GetSignCosPhi() const { return fSignCosPhi; }
96  float GetChi2() const { return fChi2; }
97  int GetNDF() const { return fNDF; }
98 
99  float GetKappa(float Bz) const { return fP[4] * Bz; }
100  float GetCosPhiPositive() const { return CAMath::Sqrt(1 - SinPhi() * SinPhi()); }
101  float GetCosPhi() const { return fSignCosPhi * CAMath::Sqrt(1 - SinPhi() * SinPhi()); }
102 
103  float GetErr2Y() const { return fC[0]; }
104  float GetErr2Z() const { return fC[2]; }
105  float GetErr2SinPhi() const { return fC[5]; }
106  float GetErr2DzDs() const { return fC[9]; }
107  float GetErr2QPt() const { return fC[14]; }
108 
109  float Err2X1() const { return fC[0]; }
110  float Err2X2() const { return fC[2]; }
111  float Err2QMomentum() const { return fC[14]; }
112 
113  const float *Par() const { return fP; }
114  const float *Cov() const { return fC; }
115 
116  const float *GetPar() const { return fP; }
117  const float *GetCov() const { return fC; }
118 
119  void SetPar(int i, float v) { fP[i] = v; }
120  void SetCov(int i, float v) { fC[i] = v; }
121 
122  void SetX(float v) { fX = v; }
123  void SetY(float v) { fP[0] = v; }
124  void SetZ(float v) { fP[1] = v; }
125  void SetSinPhi(float v) { fP[2] = v; }
126  void SetDzDs(float v) { fP[3] = v; }
127  void SetQPt(float v) { fP[4] = v; }
128  void SetSignCosPhi(float v) { fSignCosPhi = v; }
129  void SetChi2(float v) { fChi2 = v; }
130  void SetNDF(int v) { fNDF = v; }
131 
132  void SetAngle(float v) { fAlpha = v; }
133  void SetISec(int v) { fISec = v; }
134 
135  void SetErr2QPt(float v) { fC[14] = v; }
136 
137  void InitDirection(float r0, float r1, float r2) // initialize direction parameters according to a given tangent vector
138  {
139  const float r = sqrt(r0 * r0 + r1 * r1);
140  SetSinPhi(r1 / r);
141  SetSignCosPhi(sign(r0));
142  SetDzDs(r2 / r);
143  }
144 
145  float GetDist2(const PndCATrackParam &t) const;
146  float GetDistXZ2(const PndCATrackParam &t) const;
147 
148  float GetS(float x, float y, float Bz) const;
149 
150  void GetDCAPoint(float x, float y, float z, float &px, float &py, float &pz, float Bz) const;
151 
152  bool TransportToX(float x, float Bz, float maxSinPhi = .999);
153  bool TransportToXWithMaterial(float x, float Bz, float maxSinPhi = .999);
154 
155  bool TransportToX(float x, PndCATrackLinearisation &t0, float Bz, float maxSinPhi = .999, float *DL = nullptr);
156 
157  bool TransportToX(float x, float sinPhi0, float cosPhi0, float Bz, float maxSinPhi = .999);
158 
159  bool TransportToXWithMaterial(float x, PndCATrackLinearisation &t0, PndCATrackFitParam &par, float Bz, float maxSinPhi = .999);
160 
161  bool TransportToXWithMaterial(float x, PndCATrackFitParam &par, float Bz, float maxSinPhi = .999);
162 
163  static float ApproximateBetheBloch(float beta2);
164  static float BetheBlochGeant(float bg, float kp0 = 2.33, float kp1 = 0.20, float kp2 = 3.00, float kp3 = 173e-9, float kp4 = 0.49848);
165  static float BetheBlochSolid(float bg);
166  static float BetheBlochGas(float bg);
167 
168  void CalculateFitParameters(PndCATrackFitParam &par, float mass = 0.13957);
169  bool CorrectForMeanMaterial(float xOverX0, float xTimesRho, const PndCATrackFitParam &par);
170 
171  bool Rotate(float alpha, float maxSinPhi = .999);
172  bool Rotate(float alpha, PndCATrackLinearisation &t0, float maxSinPhi = .999);
173 
174  void RotateXY(float alpha, float &x, float &y, float &sin) const;
175  bool Filter(float y, float z, float err2Y, float errYZ, float err2Z, float maxSinPhi = .999);
176 
177  void Print() const;
178 
180  {
181  fC[0] = 10.f;
182  fC[1] = 0.f;
183  fC[2] = 10.f;
184  fC[3] = 0.f;
185  fC[4] = 0.f;
186  fC[5] = 1.f;
187  fC[6] = 0.f;
188  fC[7] = 0.f;
189  fC[8] = 0.f;
190  fC[9] = 1.f;
191  fC[10] = 0.f;
192  fC[11] = 0.f;
193  fC[12] = 0.f;
194  fC[13] = 0.f;
195  fC[14] = 10.f;
196  }
197 
198  PndCATrackParam GetGlobalParam(float alpha) const; // alpha - angle of the current slice
199 
200  void Reset()
201  {
202  fX = 0;
203  fSignCosPhi = 0;
204  for (int i = 0; i < 5; i++)
205  fP[i] = 0;
206  for (int i = 0; i < 15; i++)
207  fC[i] = 0;
208  fChi2 = 0;
209  fNDF = 0;
210  fISec = 0;
211  }
212 
213  bool Transport(const PndCAHit &hit, float Bz);
214  bool Filter(const PndCAHit &hit);
215 
216  bool IsValid() const { return fChi2 != -1; }
217  void SetAsInvalid() { fChi2 = -1; }
218 
219  private:
220  float fX; // x position
221  float fSignCosPhi; // sign of cosPhi // phi = arctg (Dy/Dx)
222  float fP[5]; // 'active' track parameters: Y, Z, SinPhi, Dz/Ds (ds = sqrt( dx^2 + dy^2 )), q/Pt
223  float fC[15]; // the covariance matrix for Y,Z,SinPhi,..
224  float fChi2; // the chi^2 value
225  int fNDF; // the Number of Degrees of Freedom
226 
227  float fAlpha; // coor system
228  int fISec;
229 };
230 
231 inline void PndCATrackParam::RotateXY(float alpha, float &x, float &y, float &sin) const
232 {
233  //* Rotate the coordinate system in XY on the angle alpha
234 
235  const float cA = CAMath::Cos(alpha);
236  const float sA = CAMath::Sin(alpha);
237 
238  x = (X() * cA + Y() * sA);
239  y = (-X() * sA + Y() * cA);
240  sin = -GetCosPhi() * sA + SinPhi() * cA;
241 }
242 
243 inline bool PndCATrackParam::Filter(float y, float z, float err2Y, float errYZ, float err2Z, float maxSinPhi)
244 {
245  assert(maxSinPhi > 0.f);
246  //* Add the y,z measurement with the Kalman filter
247 
248  const float c00 = fC[0];
249  const float c10 = fC[1];
250  const float c11 = fC[2];
251  const float c20 = fC[3];
252  const float c21 = fC[4];
253  // float c22 = fC[5];
254  const float c30 = fC[6];
255  const float c31 = fC[7];
256  // float c32 = fC[8];
257  // float c33 = fC[9];
258  const float c40 = fC[10];
259  const float c41 = fC[11];
260  // float c42 = fC[12];
261  // float c43 = fC[13];
262  // float c44 = fC[14];
263 
264  float d = 1.f / (err2Y * err2Z + err2Y * c11 + err2Z * c00 + c00 * c11 - c10 * c10 - 2 * errYZ * c10 - errYZ * errYZ);
265  err2Y += c00;
266  err2Z += c11;
267  errYZ += c10;
268 
269  const float z0 = y - fP[0], z1 = z - fP[1];
270 
271  if (ISUNLIKELY(err2Y < 1.e-8f) || ISUNLIKELY(err2Z < 1.e-8f))
272  return 0;
273 
274  const float mS0 = err2Z * d;
275  const float mS1 = -errYZ * d;
276  const float mS2 = err2Y * d;
277 
278  // K = CHtS
279 
280  const float k00 = c00 * mS0 + c10 * mS1, k01 = c00 * mS1 + c10 * mS2, k10 = c10 * mS0 + c11 * mS1, k11 = c10 * mS1 + c11 * mS2, k20 = c20 * mS0 + c21 * mS1,
281  k21 = c20 * mS1 + c21 * mS2, k30 = c30 * mS0 + c31 * mS1, k31 = c30 * mS1 + c31 * mS2, k40 = c40 * mS0 + c41 * mS1, k41 = c40 * mS1 + c41 * mS2;
282 
283  const float sinPhi = fP[2] + k20 * z0 + k21 * z1;
284 
285  if (ISUNLIKELY(CAMath::Abs(sinPhi) >= maxSinPhi))
286  return 0;
287 
288  fNDF += 2;
289  fChi2 += mS0 * z0 * z0 + mS2 * z1 * z1 + 2 * z0 * z1 * mS1;
290 
291  fP[0] += k00 * z0 + k01 * z1;
292  fP[1] += k10 * z0 + k11 * z1;
293  fP[2] = sinPhi;
294  fP[3] += k30 * z0 + k31 * z1;
295  fP[4] += k40 * z0 + k41 * z1;
296 
297  fC[0] -= (k00 * c00 + k01 * c10); // c00
298 
299  fC[1] -= (k10 * c00 + k11 * c10); // c10
300  fC[2] -= (k10 * c10 + k11 * c11); // c11
301 
302  fC[3] -= (k20 * c00 + k21 * c10); // c20
303  fC[4] -= (k20 * c10 + k21 * c11); // c21
304  fC[5] -= (k20 * c20 + k21 * c21); // c22
305 
306  fC[6] -= (k30 * c00 + k31 * c10); // c30
307  fC[7] -= (k30 * c10 + k31 * c11); // c31
308  fC[8] -= (k30 * c20 + k31 * c21); // c32
309  fC[9] -= (k30 * c30 + k31 * c31); // c33
310 
311  fC[10] -= (k40 * c00 + k41 * c10); // c40
312  fC[11] -= (k40 * c10 + k41 * c11); // c41
313  fC[12] -= (k40 * c20 + k41 * c21); // c42
314  fC[13] -= (k40 * c30 + k41 * c31); // c43
315  fC[14] -= (k40 * c40 + k41 * c41); // c44
316 
317  return 1;
318 }
319 
320 inline float PndCATrackParam::ApproximateBetheBloch(float beta2)
321 {
322  //------------------------------------------------------------------
323  // This is an approximation of the Bethe-Bloch formula with
324  // the density effect taken into account at beta*gamma > 3.5
325  // (the approximation is reasonable only for solid materials)
326  //------------------------------------------------------------------
327  if (beta2 >= 1)
328  return 0;
329  else {
330  const float beta2_beta21i = beta2 / (1 - beta2);
331  if (beta2_beta21i > 12.25) // 3.5^2 = 12.25
332  return 0.153e-3 / beta2 * (9.94223 + 0.5 * log(beta2_beta21i) - beta2); // log( 3.5*5940 ) = 9.94223
333  else
334  return 0.153e-3 / beta2 * (8.6895 + log(beta2_beta21i) - beta2); // log( 5940 ) = 8.6895
335  }
336 }
337 
339 {
340  const float p2 = (1. + fP[3] * fP[3]);
341  const float k2 = fP[4] * fP[4];
342  const float mass2 = mass * mass;
343 
344  const float beta2 = p2 / (p2 + mass2 * k2);
345 
346  const float pp2 = (k2 > 1.e-8) ? p2 / k2 : 10000; // impuls 2
347 
348  // par.fBethe = BetheBlochGas( pp2/mass2);
349  par.fBethe = ApproximateBetheBloch(pp2 / mass2);
350  par.fE = CAMath::Sqrt(pp2 + mass2);
351  par.fTheta2 = 198.81e-6 / (beta2 * pp2); // 14.1^2 * 1e-6
352  par.fEP2 = par.fE / pp2; // have tried reduce number of "/", but it was slower. (may be bacause of additional of constants = memory)
353 
354  // Approximate energy loss fluctuation (M.Ivanov)
355 
356  const float knst = 0.07; // To be tuned.
357  par.fSigmadE2 = knst * par.fEP2 * fP[4];
358  par.fSigmadE2 = par.fSigmadE2 * par.fSigmadE2;
359 
360  par.fK22 = p2;
361  par.fK33 = par.fK22 * par.fK22;
362  par.fK43 = fP[3] * fP[4] * par.fK22;
363  par.fK44 = (p2 - 1.f) * k2;
364 }
365 
366 #include "PndCAMath.h"
367 #include "PndCATrackLinearisation.h"
368 
369 inline bool PndCATrackParam::TransportToX(float x, PndCATrackLinearisation &t0, float Bz, float maxSinPhi, float *DL)
370 {
371  //* Transport the track parameters to X=x, using linearization at t0, and the field value Bz
372  //* maxSinPhi is the max. allowed value for |t0.SinPhi()|
373  //* linearisation of trajectory t0 is also transported to X=x,
374  //* returns 1 if OK
375  //*
376 
377  const float ex = t0.CosPhi();
378  const float ey = t0.SinPhi();
379  const float k = t0.QPt() * Bz;
380  const float dx = x - X();
381 
382  const float ey1 = k * dx + ey;
383 
384  // check for intersection with X=x
385 
386  if (CAMath::Abs(ey1) > maxSinPhi)
387  return 0;
388 
389  float ex1 = CAMath::Sqrt(1.f - ey1 * ey1);
390  if (ex < 0)
391  ex1 = -ex1;
392 
393  const float dx2 = dx * dx;
394  const float ss = ey + ey1;
395  const float cc = ex + ex1;
396 
397  if ((CAMath::Abs(cc) < 1.e-4 || CAMath::Abs(ex) < 1.e-4 || CAMath::Abs(ex1) < 1.e-4))
398  return 0;
399 
400  const float cci = 1.f / cc;
401  const float exi = 1.f / ex;
402  const float ex1i = 1.f / ex1;
403 
404  const float tg = ss * cci; // tan((phi1+phi)/2)
405 
406  const float dy = dx * tg;
407  float dl = dx * CAMath::Sqrt(1.f + tg * tg);
408 
409  if (cc < 0)
410  dl = -dl;
411  float dSin = dl * k * 0.5;
412  if (dSin > 1.f)
413  dSin = 1.f;
414  if (dSin < -1.f)
415  dSin = -1.f;
416  const float dS = (CAMath::Abs(k) > 1.e-4) ? (2 * CAMath::ASin(dSin) / k) : dl;
417  const float dz = dS * t0.DzDs();
418 
419  if (DL)
420  *DL = -dS * CAMath::Sqrt(1.f + t0.DzDs() * t0.DzDs());
421 
422  const float d[3] = {fP[2] - t0.SinPhi(), fP[3] - t0.DzDs(), fP[4] - t0.QPt()};
423 
424  // float H0[5] = { 1,0, h2, 0, h4 };
425  // float H1[5] = { 0, 1, 0, dS, 0 };
426  // float H2[5] = { 0, 0, 1, 0, dxBz };
427  // float H3[5] = { 0, 0, 0, 1, 0 };
428  // float H4[5] = { 0, 0, 0, 0, 1 };
429 
430  const float h2 = dx * (1.f + ey * ey1 + ex * ex1) * exi * ex1i * cci;
431  const float h4 = dx2 * (cc + ss * ey1 * ex1i) * cci * cci * Bz;
432  const float dxBz = dx * Bz;
433 
434  t0.SetCosPhi(ex1);
435  t0.SetSinPhi(ey1);
436 
437  fX = X() + dx;
438  fP[0] = Y() + dy + h2 * d[0] + h4 * d[2];
439  fP[1] = Z() + dz + dS * d[1];
440  fP[2] = t0.SinPhi() + d[0] + dxBz * d[2];
441  if (CAMath::Abs(fP[2]) > maxSinPhi)
442  fP[2] = t0.SinPhi();
443 
444 #if 1
445  const float c00 = fC[0];
446  const float c10 = fC[1];
447  const float c11 = fC[2];
448  const float c20 = fC[3];
449  const float c21 = fC[4];
450  const float c22 = fC[5];
451  const float c30 = fC[6];
452  const float c31 = fC[7];
453  const float c32 = fC[8];
454  const float c33 = fC[9];
455  const float c40 = fC[10];
456  const float c41 = fC[11];
457  const float c42 = fC[12];
458  const float c43 = fC[13];
459  const float c44 = fC[14];
460 
461  fC[0] = (c00 + h2 * h2 * c22 + h4 * h4 * c44 + 2.f * (h2 * c20 + h4 * c40 + h2 * h4 * c42));
462 
463  fC[1] = c10 + h2 * c21 + h4 * c41 + dS * (c30 + h2 * c32 + h4 * c43);
464  fC[2] = c11 + 2.f * dS * c31 + dS * dS * c33;
465 
466  fC[3] = c20 + h2 * c22 + h4 * c42 + dxBz * (c40 + h2 * c42 + h4 * c44);
467  fC[4] = c21 + dS * c32 + dxBz * (c41 + dS * c43);
468  fC[5] = c22 + 2.f * dxBz * c42 + dxBz * dxBz * c44;
469 
470  fC[6] = c30 + h2 * c32 + h4 * c43;
471  fC[7] = c31 + dS * c33;
472  fC[8] = c32 + dxBz * c43;
473  fC[9] = c33;
474 
475  fC[10] = c40 + h2 * c42 + h4 * c44;
476  fC[11] = c41 + dS * c43;
477  fC[12] = c42 + dxBz * c44;
478  fC[13] = c43;
479  fC[14] = c44;
480 #else
481  fC[0] = (fC[0] + h2 * h2 * fC[5] + h4 * h4 * fC[14] + 2 * (h2 * fC[3] + h4 * fC[10] + h2 * h4 * fC[12]));
482 
483  fC[1] = fC[1] + h2 * fC[4] + h4 * fC[11] + dS * (fC[6] + h2 * fC[8] + h4 * fC[13]);
484  fC[2] = fC[2] + 2 * dS * fC[7] + dS * dS * fC[9];
485 
486  fC[3] = fC[3] + h2 * fC[5] + h4 * fC[12] + dxBz * (fC[10] + h2 * fC[12] + h4 * fC[14]);
487  fC[4] = fC[4] + dS * fC[8] + dxBz * (fC[11] + dS * fC[13]);
488  fC[5] = fC[5] + 2 * dxBz * fC[12] + dxBz * dxBz * fC[14];
489 
490  fC[6] = fC[6] + h2 * fC[8] + h4 * fC[13];
491  fC[7] = fC[7] + dS * fC[9];
492  fC[8] = fC[8] + dxBz * fC[13];
493  fC[9] = fC[9];
494 
495  fC[10] = fC[10] + h2 * fC[12] + h4 * fC[14];
496  fC[11] = fC[11] + dS * fC[13];
497  fC[12] = fC[12] + dxBz * fC[14];
498  fC[13] = fC[13];
499  fC[14] = fC[14];
500 #endif
501 
502  // std::cout << fC[0] << " "<<fC1[0]<<" "<<fC[2] << " "<<fC1[2]<<" "<<fC[5] << " "<<fC1[5]<<" "<<fC[9] << " "<<fC1[9]<<" "<<fC[14] << "
503  // "<<fC1[14]<<std::endl;
504  return 1;
505 }
506 
507 inline bool PndCATrackParam::CorrectForMeanMaterial(float xOverX0, float xTimesRho, const PndCATrackFitParam &par)
508 {
509  //------------------------------------------------------------------
510  // This function corrects the track parameters for the crossed material.
511  // "xOverX0" - X/X0, the thickness in units of the radiation length.
512  // "xTimesRho" - is the product length*density (g/cm^2).
513  //------------------------------------------------------------------
514  // float &fC22 = fC[5];
515  // float &fC33 = fC[9];
516  // float &fC40 = fC[10];
517  // float &fC41 = fC[11];
518  // float &fC42 = fC[12];
519  // float &fC43 = fC[13];
520  // float &fC44 = fC[14];
521 
522  // Energy losses************************
523 
524  const float dE = par.fBethe * xTimesRho;
525  if (CAMath::Abs(dE) > 0.3 * par.fE)
526  return 0; // 30% energy loss is too much!
527  const float corr = (1. - par.fEP2 * dE);
528  if (corr < 0.3 || corr > 1.3)
529  return 0;
530 
531  fP[4] *= corr;
532  fC[10] *= corr;
533  fC[11] *= corr;
534  fC[12] *= corr;
535  fC[13] *= corr;
536  fC[14] *= corr * corr;
537  fC[14] += par.fSigmadE2 * CAMath::Abs(dE);
538  // std::cout << "dE "<<dE<<" corr "<<corr<<" fBethe " <<par.fBethe<<" XxRo "<<xTimesRho<<std::endl;
539 
540  // Multiple scattering******************
541 
542  const float theta2 = par.fTheta2 * CAMath::Abs(xOverX0);
543  fC[5] += theta2 * par.fK22 * (1. - fP[2] * fP[2]);
544  fC[9] += theta2 * par.fK33;
545  fC[13] += theta2 * par.fK43;
546  fC[14] += theta2 * par.fK44;
547 
548  return 1;
549 }
550 
551 inline bool PndCATrackParam::TransportToXWithMaterial(float x, PndCATrackLinearisation &t0, PndCATrackFitParam &par, float Bz, float maxSinPhi)
552 {
553  //* Transport the track parameters to X=x taking into account material budget
554 
555  const float kRho = 1.54e-3; // 1.025e-3 ;//0.9e-3;
556  // const float kRadLen = 29.532;//28.94;
557  // const float kRhoOverRadLen = kRho / kRadLen;
558  const float kRhoOverRadLen = 7.68e-5;
559  float dl;
560 
561  if (!TransportToX(x, t0, Bz, maxSinPhi, &dl))
562  return 0;
563 
564  UNUSED_PARAM3(kRho, kRhoOverRadLen, par); // TODO
565  // CorrectForMeanMaterial( dl*kRhoOverRadLen, dl*kRho, par );
566  return 1;
567 }
568 
569 inline bool PndCATrackParam::TransportToXWithMaterial(float x, PndCATrackFitParam &par, float Bz, float maxSinPhi)
570 {
571  //* Transport the track parameters to X=x taking into account material budget
572 
573  PndCATrackLinearisation t0(*this);
574  return TransportToXWithMaterial(x, t0, par, Bz, maxSinPhi);
575 }
576 
577 inline bool PndCATrackParam::Transport(const PndCAHit &hit, float Bz)
578 {
579  // TODO material. See Vector part
580  PndCATrackFitParam fitPar;
581  CalculateFitParameters(fitPar);
582  PndCATrackLinearisation tR(*this);
583  const bool rotated = Rotate(-fAlpha + hit.Angle(), tR, .999f);
584  PndCATrackLinearisation tE(*this);
585  const bool transported = TransportToXWithMaterial(hit.X0(), tE, fitPar, Bz, 0.999f);
586  return rotated & transported;
587 }
588 
589 inline bool PndCATrackParam::Filter(const PndCAHit &hit)
590 {
591  return Filter(hit.X1(), hit.X2(), hit.Err2X1(), hit.ErrX12(), hit.Err2X2(), 0.999f);
592 }
593 
595 
596 // std::istream &operator>>( std::istream &in, PndCATrackParam &ot );
597 // std::ostream &operator<<( std::ostream &out, const PndCATrackParam &ot );
598 
599 #endif
static T ASin(const T &x)
float DzDs() const
bool CorrectForMeanMaterial(float xOverX0, float xTimesRho, const PndCATrackFitParam &par)
float GetErr2QPt() const
float GetDist2(const PndCATrackParam &t) const
void SetZ(float v)
float GetY() const
bool Filter(float y, float z, float err2Y, float errYZ, float err2Z, float maxSinPhi=.999)
float GetS(float x, float y, float Bz) const
const float * GetCov() const
const float * Cov() const
void SetCov(int i, float v)
bool Transport(const PndCAHit &hit, float Bz)
bool Rotate(float alpha, float maxSinPhi=.999)
friend F32vec4 sqrt(const F32vec4 &a)
Definition: P4_F32vec4.h:40
static float BetheBlochGeant(float bg, float kp0=2.33, float kp1=0.20, float kp2=3.00, float kp3=173e-9, float kp4=0.49848)
float Chi2() const
static T Sqrt(const T &x)
Definition: PndCAMath.h:57
bool IsValid() const
float Angle() const
Definition: PndCAHits.h:64
friend F32vec4 sin(const F32vec4 &a)
Definition: P4_F32vec4.h:130
static T Sin(const T &x)
Definition: PndCAMath.h:83
PndCATrackParam TrackParam
void SetSignCosPhi(float v)
float GetQPt() const
friend F32vec4 log(const F32vec4 &a)
Definition: P4_F32vec4.h:129
void SetX(float v)
float Err2QPt() const
void GetDCAPoint(float x, float y, float z, float &px, float &py, float &pz, float Bz) const
float Err2X1() const
Definition: PndCAHits.h:51
const float * Par() const
static T Cos(const T &x)
Definition: PndCAMath.h:88
bool TransportToXWithMaterial(float x, float Bz, float maxSinPhi=.999)
float SinPhi() const
__m128 v
Definition: P4_F32vec4.h:15
const float * GetPar() const
float Z() const
unsigned int i
Definition: P4_F32vec4.h:33
float Err2DzDs() const
static T Abs(const T &x)
Definition: PndCAMath.h:68
float Y() const
void InitDirection(float r0, float r1, float r2)
void SetSinPhi(float v)
void SetDzDs(float v)
float GetErr2SinPhi() const
float Err2Y() const
static float BetheBlochSolid(float bg)
float GetX() const
void SetISec(int v)
int NDF() const
void SetPar(int i, float v)
float GetSignCosPhi() const
float GetDistXZ2(const PndCATrackParam &t) const
float GetErr2Y() const
void SetQPt(float v)
float X() const
int ISec() const
float Err2X2() const
Definition: PndCAHits.h:53
float Err2QMomentum() const
float Err2X1() const
float GetDzDs() const
float SignCosPhi() const
float Angle() const
float GetChi2() const
float GetKappa(float Bz) const
bool TransportToX(float x, float Bz, float maxSinPhi=.999)
float X1() const
Definition: PndCAHits.h:48
float GetErr2Z() const
float f
Definition: P4_F32vec4.h:32
void Print() const
#define ISUNLIKELY(x)
Definition: PndCADef.h:169
static float BetheBlochGas(float bg)
PndCATrackParam GetGlobalParam(float alpha) const
int sign(T val)
Definition: PndCADef.h:61
float QPt() const
double alpha
Definition: f_Init.h:31
float GetCosPhi() const
static float ApproximateBetheBloch(float beta2)
float GetErr2DzDs() const
void SetErr2QPt(float v)
float X2() const
Definition: PndCAHits.h:49
float Err2X2() const
const float_v & Cov(int i) const
float Err2SinPhi() const
int GetNDF() const
float ErrX12() const
Definition: PndCAHits.h:52
float GetSinPhi() const
void SetAngle(float v)
float QMomentum() const
void SetY(float v)
void CalculateFitParameters(PndCATrackFitParam &par, float mass=0.13957)
float GetCosPhiPositive() const
float GetZ() const
PndCATrackParam(const TrackParamVector &v, int i)
double pz[39]
Definition: pipisigmas.h:25
float Err2Z() const
void SetChi2(float v)
float X0() const
Definition: PndCAHits.h:47
const float_v & Par(int i) const
void RotateXY(float alpha, float &x, float &y, float &sin) const
void SetNDF(int v)