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