CbmRoot
CbmStsKFTrackFitter.cxx
Go to the documentation of this file.
1 #include "CbmStsKFTrackFitter.h"
2 
3 #include "CbmKFMath.h"
4 #include "CbmKFStsHit.h"
5 #include "CbmKFTrack.h"
6 #include "CbmKFVertex.h"
7 
8 #include "CbmMvdHit.h"
9 #include "CbmStsTrack.h"
10 #include "CbmVertex.h"
11 #include "FairRootManager.h"
12 
13 #include "TClonesArray.h"
14 #include "TMath.h"
15 #include "TMatrixTSym.h"
16 
17 #include "math.h"
18 #include <iostream>
19 
20 using std::cout;
21 using std::endl;
22 
23 
24 ClassImp(CbmStsKFTrackFitter);
25 
26 CbmStsKFTrackFitter::CbmStsKFTrackFitter()
27  : fHits(), fMvdHitsArray(0), fStsHitsArray(0), fIsInitialised(0) {}
28 
29 void CbmStsKFTrackFitter::Init() {
30  // Initialisation
31  FairRootManager* rootMgr = FairRootManager::Instance();
32  if (NULL == rootMgr) {
33  cout << "-E- CbmStsKFTrackFitter::Init(): "
34  << "ROOT manager is not instantiated!" << endl;
35  return;
36  }
37  fStsHitsArray = reinterpret_cast<TClonesArray*>(rootMgr->GetObject("StsHit"));
38  if (!fStsHitsArray) {
39  cout << "-W- CbmStsKFTrackFitter::Init: "
40  << "no STS hits array" << endl;
41  //return;
42  }
43  fMvdHitsArray = reinterpret_cast<TClonesArray*>(rootMgr->GetObject("MvdHit"));
44  if (!fMvdHitsArray) {
45  cout << "-W- CbmStsKFTrackFitter::Init: "
46  << "no MVD hits array" << endl;
47  //return;
48  }
49  fIsInitialised = 1;
50 };
51 
52 void CbmStsKFTrackFitter::SetKFHits(CbmKFTrack& T, CbmStsTrack* track) {
53 
54  T.fHits.clear();
55 
56  if (!fIsInitialised) Init();
57 
58  Int_t NStsHits = (fStsHitsArray) ? track->GetNofHits() : 0;
59  Int_t NMvdHits = (fMvdHitsArray) ? track->GetNofMvdHits() : 0;
60 
61  fHits.resize(NMvdHits + NStsHits);
62  if (NMvdHits > 0) {
63  for (Int_t i = 0; i < NMvdHits; i++) {
64  Int_t j = track->GetMvdHitIndex(i);
65  fHits[i].Create(reinterpret_cast<CbmMvdHit*>(fMvdHitsArray->At(j)));
66  T.fHits.push_back(&(fHits[i]));
67  }
68  }
69  if (NStsHits > 0 && fStsHitsArray) {
70  for (Int_t i = 0; i < NStsHits; i++) {
71  Int_t j = track->GetHitIndex(i);
72  fHits[NMvdHits + i].Create(
73  reinterpret_cast<CbmStsHit*>(fStsHitsArray->At(j)));
74  T.fHits.push_back(&(fHits[NMvdHits + i]));
75  }
76  }
77 }
78 
79 Int_t CbmStsKFTrackFitter::DoFit(CbmStsTrack* track, Int_t pidHypo) {
80  track->SetPidHypo(pidHypo);
81 
82  CbmKFTrack T;
83  T.SetPID(pidHypo);
84  SetKFHits(T, track);
85  for (Int_t i = 0; i < 6; i++)
86  T.GetTrack()[i] = 0.; // no guess taken
87  T.Fit(1); // fit downstream
88  CheckTrack(T);
89  T.Fit(0); // fit upstream
90  CheckTrack(T);
91  Int_t err = T.Fit(1); // fit downstream
92  Bool_t ok = (!err) && CheckTrack(T);
93  if (ok) {
94  FairTrackParam par;
95  //T.GetTrackParam( *track->GetParamLast() ); // store fitted track & cov.matrix
96  T.GetTrackParam(par);
97  track->SetParamLast(&par);
98  err = T.Fit(0); // fit upstream
99  ok = ok && (!err) && CheckTrack(T);
100  if (ok)
101  T.GetStsTrack(*track, 1); // store fitted track & cov.matrix & chi2 & NDF
102  }
103  if (!ok) {
104  Double_t* t = T.GetTrack();
105  Double_t* c = T.GetCovMatrix();
106  for (int i = 0; i < 6; i++)
107  t[i] = 0;
108  for (int i = 0; i < 15; i++)
109  c[i] = 0;
110  c[0] = c[2] = c[5] = c[9] = c[14] = 100.;
111  T.GetRefChi2() = 100.;
112  T.GetRefNDF() = 0;
113  T.GetStsTrack(*track, 0);
114  T.GetStsTrack(*track, 1);
115  track->SetFlag(1);
116  } else {
117  track->SetFlag(0);
118  }
119  return !ok;
120 }
121 
122 
123 void CbmStsKFTrackFitter::Extrapolate(FairTrackParam* track,
124  Double_t z,
125  FairTrackParam* e_track) {
126  if (!track) return;
127  CbmKFTrack T;
128  T.SetTrackParam(*track);
129  T.Extrapolate(z);
130  if (e_track) T.GetTrackParam(*e_track);
131 }
132 
133 
134 void CbmStsKFTrackFitter::Extrapolate(CbmStsTrack* track,
135  Double_t z,
136  FairTrackParam* e_track) {
137  if (!track) return;
138  CbmKFTrack T;
139  T.SetPID(track->GetPidHypo());
140  const FairTrackParam *fpar = track->GetParamFirst(),
141  *lpar = track->GetParamLast();
142 
143  if (z <= fpar->GetZ()) { // extrapolate first parameters
144  T.SetTrackParam(*fpar);
145  T.Extrapolate(z);
146  } else if (z < fpar->GetZ() + 0.1) { // extrapolate first parameters
147  T.SetTrackParam(*fpar);
148  T.Propagate(z);
149  } else if (lpar->GetZ() <= z) { // extrapolate last parameters
150  T.SetTrackParam(*lpar);
151  T.Extrapolate(z);
152  } else if (lpar->GetZ() - 0.1 < z) { // extrapolate last parameters
153  T.SetTrackParam(*lpar);
154  T.Propagate(z);
155  } else { // refit with smoother
156  SetKFHits(T, track);
157  T.SetTrackParam(*fpar);
158  T.Smooth(z);
159  }
160  if (e_track) T.GetTrackParam(*e_track);
161 }
162 
163 
164 Double_t CbmStsKFTrackFitter::GetChiToVertex(CbmStsTrack* track,
165  CbmVertex* vtx) {
166  if (!vtx) {
167  FairRootManager* fManger = FairRootManager::Instance();
168  // Get pointer to PrimaryVertex object from IOManager if it exists
169  // The old name for the object is "PrimaryVertex" the new one
170  // "PrimaryVertex." Check first for the new name
171  // TODO: don't use reinterpret_cast
172  vtx = reinterpret_cast<CbmVertex*>(fManger->GetObject("PrimaryVertex."));
173  if (nullptr == vtx) {
174  vtx = reinterpret_cast<CbmVertex*>(fManger->GetObject("PrimaryVertex"));
175  }
176  if (!vtx) {
177  cout
178  << "-W- CbmStsKFTrackFitter::GetChiToVertex: No Primary Vertex found!"
179  << endl;
180  return 100.;
181  }
182  }
183  CbmKFTrack T;
184  T.SetStsTrack(*track, 1);
185  T.Extrapolate(vtx->GetZ());
186 
187  TMatrixFSym tmp(3);
188  vtx->CovMatrix(tmp);
189  Double_t Cv[3] = {tmp(0, 0), tmp(0, 1), tmp(1, 1)};
190 
191  return CbmKFMath::getDeviation(T.GetTrack()[0],
192  T.GetTrack()[1],
193  T.GetCovMatrix(),
194  vtx->GetX(),
195  vtx->GetY(),
196  Cv);
197 }
198 
199 
200 Double_t CbmStsKFTrackFitter::FitToVertex(CbmStsTrack* track,
201  CbmVertex* vtx,
202  FairTrackParam* v_track) {
203  Double_t ret = 100.;
204  if (!track || !vtx || !v_track) return ret;
205  CbmKFTrack T(*track);
206  CbmKFVertex V(*vtx);
207  T.Fit2Vertex(V);
208  T.GetTrackParam(*v_track);
209  if (T.GetRefNDF() > 0 && T.GetRefChi2() >= 0) {
210  ret = T.GetRefChi2() / T.GetRefNDF();
211  if (finite(ret)) ret = sqrt(ret);
212  }
213  return ret;
214 }
215 
216 Bool_t CbmStsKFTrackFitter::CheckTrack(CbmKFTrack& T) {
217  Bool_t ok = 1;
218  Double_t* t = T.GetTrack();
219  Double_t* c = T.GetCovMatrix();
220  for (int i = 0; i < 6; i++)
221  ok = ok && finite(t[i]) && TMath::Abs(t[i]) < 1.e5;
222  for (int i = 0; i < 15; i++)
223  ok = ok && finite(c[i]);
224  ok = ok && finite(T.GetMass()) && finite(T.GetRefChi2());
225  if (ok) {
226  ok = ok && (c[0] > 0);
227  ok = ok && (c[2] > 0);
228  ok = ok && (c[5] > 0);
229  ok = ok && (c[9] > 0);
230  ok = ok && (c[14] > 0);
231  }
232  if (ok) { // correct the cov matrix
233  Double_t c00 = TMath::Sqrt(c[0]);
234  Double_t c11 = TMath::Sqrt(c[2]);
235  Double_t c22 = TMath::Sqrt(c[5]);
236  Double_t c33 = TMath::Sqrt(c[9]);
237  Double_t c44 = TMath::Sqrt(c[14]);
238  Double_t a = c11 * c00;
239  if (c[1] > a) c[1] = a;
240  if (c[1] < -a) c[1] = -a;
241  a = c22 * c00;
242  if (c[3] > a) c[3] = a;
243  if (c[3] < -a) c[3] = -a;
244  a = c22 * c11;
245  if (c[4] > a) c[4] = a;
246  if (c[4] < -a) c[4] = -a;
247  a = c33 * c00;
248  if (c[6] > a) c[6] = a;
249  if (c[6] < -a) c[6] = -a;
250  a = c33 * c11;
251  if (c[7] > a) c[7] = a;
252  if (c[7] < -a) c[7] = -a;
253  a = c33 * c22;
254  if (c[8] > a) c[8] = a;
255  if (c[8] < -a) c[8] = -a;
256  a = c44 * c00;
257  if (c[10] > a) c[10] = a;
258  if (c[10] < -a) c[10] = -a;
259  a = c44 * c11;
260  if (c[11] > a) c[11] = a;
261  if (c[11] < -a) c[11] = -a;
262  a = c44 * c22;
263  if (c[12] > a) c[12] = a;
264  if (c[12] < -a) c[12] = -a;
265  a = c44 * c33;
266  if (c[13] > a) c[13] = a;
267  if (c[13] < -a) c[13] = -a;
268  }
269  if (!ok) {
270  for (int i = 0; i < 6; i++)
271  t[i] = 0;
272  for (int i = 0; i < 15; i++)
273  c[i] = 0;
274  c[0] = c[2] = c[5] = c[9] = c[14] = 100.;
275  T.GetRefChi2() = 100.;
276  T.GetRefNDF() = 0;
277  }
278  return ok;
279 }
CbmVertex.h
CbmTrack::GetParamLast
const FairTrackParam * GetParamLast() const
Definition: CbmTrack.h:62
CbmStsKFTrackFitter::CbmStsKFTrackFitter
CbmStsKFTrackFitter()
Definition: CbmStsKFTrackFitter.cxx:26
CbmStsKFTrackFitter::Extrapolate
void Extrapolate(CbmStsTrack *track, Double_t z, FairTrackParam *e_track)
Definition: CbmStsKFTrackFitter.cxx:134
CbmStsKFTrackFitter::SetKFHits
void SetKFHits(CbmKFTrack &T, CbmStsTrack *track)
Definition: CbmStsKFTrackFitter.cxx:52
CbmKFStsHit.h
sqrt
friend F32vec4 sqrt(const F32vec4 &a)
Definition: L1/vectors/P4_F32vec4.h:41
CbmKFTrack::SetPID
void SetPID(Int_t pidHypo)
Definition: CbmKFTrack.cxx:58
rootMgr
static FairRootManager * rootMgr
Definition: CbmDeviceHitBuilderTof.cxx:72
CbmKFTrack::GetTrack
Double_t * GetTrack()
Is it electron.
Definition: CbmKFTrack.h:58
CbmStsKFTrackFitter::GetChiToVertex
Double_t GetChiToVertex(CbmStsTrack *track, CbmVertex *vtx=0)
Definition: CbmStsKFTrackFitter.cxx:164
CbmKFTrack::GetTrackParam
void GetTrackParam(FairTrackParam &track)
Definition: CbmKFTrack.cxx:45
i
int i
Definition: L1/vectors/P4_F32vec4.h:25
CbmKFTrackInterface::Fit2Vertex
void Fit2Vertex(CbmKFVertexInterface &vtx)
Definition: CbmKFTrackInterface.cxx:350
CbmKFTrack::GetMass
Double_t GetMass()
Definition: CbmKFTrack.h:62
CbmTrack::SetParamLast
void SetParamLast(const FairTrackParam *par)
Definition: CbmTrack.h:76
CbmStsKFTrackFitter::DoFit
Int_t DoFit(CbmStsTrack *track, Int_t pidHypo=211)
Definition: CbmStsKFTrackFitter.cxx:79
CbmVertex::GetX
Double_t GetX() const
Definition: CbmVertex.h:68
CbmStsTrack::GetNofHits
virtual Int_t GetNofHits() const
Definition: CbmStsTrack.h:76
CbmMvdHit
Definition: CbmMvdHit.h:29
CbmTrack::SetFlag
void SetFlag(Int_t flag)
Definition: CbmTrack.h:69
CbmStsTrack::GetNofMvdHits
Int_t GetNofMvdHits() const
Definition: CbmStsTrack.h:84
CbmKFMath.h
CbmKFTrack::GetRefChi2
Double_t & GetRefChi2()
array[15] of covariance matrix
Definition: CbmKFTrack.h:60
CbmStsHit
data class for a reconstructed 3-d hit in the STS
Definition: CbmStsHit.h:31
CbmKFMath::getDeviation
static Double_t getDeviation(Double_t x, Double_t y, Double_t C[], Double_t vx, Double_t vy, Double_t Cv[]=0)
Definition: CbmKFMath.cxx:324
finite
T finite(T x)
Definition: CbmL1Def.h:21
CbmKFTrack::SetStsTrack
void SetStsTrack(CbmStsTrack &track, bool first=1)
Definition: CbmKFTrack.cxx:38
CbmStsTrack.h
Data class for STS tracks.
CbmKFTrack::GetRefNDF
Int_t & GetRefNDF()
Chi^2 after fit.
Definition: CbmKFTrack.h:61
ClassImp
ClassImp(CbmStsKFTrackFitter)
CbmStsKFTrackFitter
Definition: CbmStsKFTrackFitter.h:14
CbmTrack::GetHitIndex
Int_t GetHitIndex(Int_t iHit) const
Definition: CbmTrack.h:54
CbmKFTrack::GetCovMatrix
Double_t * GetCovMatrix()
array[6] of track parameters(x,y,tx,ty,qp,z)
Definition: CbmKFTrack.h:59
CbmVertex
Definition: CbmVertex.h:26
CbmTrack::SetPidHypo
void SetPidHypo(Int_t pid)
Definition: CbmTrack.h:68
CbmVertex::GetZ
Double_t GetZ() const
Definition: CbmVertex.h:70
CbmKFTrack.h
CbmStsKFTrackFitter::CheckTrack
Bool_t CheckTrack(CbmKFTrack &T)
Definition: CbmStsKFTrackFitter.cxx:216
CbmStsKFTrackFitter::fMvdHitsArray
TClonesArray * fMvdHitsArray
Definition: CbmStsKFTrackFitter.h:72
CbmTrack::GetParamFirst
const FairTrackParam * GetParamFirst() const
Definition: CbmTrack.h:61
CbmKFTrackInterface::Smooth
void Smooth(Double_t Z)
Definition: CbmKFTrackInterface.cxx:214
CbmStsKFTrackFitter::fHits
std::vector< CbmKFStsHit > fHits
Definition: CbmStsKFTrackFitter.h:70
CbmStsKFTrackFitter::Init
void Init()
Definition: CbmStsKFTrackFitter.cxx:29
CbmStsKFTrackFitter::fStsHitsArray
TClonesArray * fStsHitsArray
Definition: CbmStsKFTrackFitter.h:73
CbmVertex::GetY
Double_t GetY() const
Definition: CbmVertex.h:69
CbmMvdHit.h
CbmKFTrack::GetStsTrack
void GetStsTrack(CbmStsTrack &track, bool first=1)
Definition: CbmKFTrack.cxx:49
CbmStsKFTrackFitter::FitToVertex
Double_t FitToVertex(CbmStsTrack *track, CbmVertex *vtx, FairTrackParam *v_track)
Definition: CbmStsKFTrackFitter.cxx:200
CbmStsTrack::GetMvdHitIndex
Int_t GetMvdHitIndex(Int_t iHit) const
Definition: CbmStsTrack.h:70
CbmTrack::GetPidHypo
Int_t GetPidHypo() const
Definition: CbmTrack.h:56
CbmKFTrackInterface::Propagate
Int_t Propagate(Double_t z_out, Double_t QP0)
Definition: CbmKFTrackInterface.cxx:449
CbmStsTrack
Definition: CbmStsTrack.h:37
CbmVertex::CovMatrix
void CovMatrix(TMatrixFSym &covMat) const
Definition: CbmVertex.cxx:104
CbmKFTrack
Definition: CbmKFTrack.h:21
CbmKFTrackInterface::Fit
Int_t Fit(Bool_t downstream=1)
Definition: CbmKFTrackInterface.cxx:101
CbmKFVertex.h
CbmKFTrack::SetTrackParam
void SetTrackParam(const FairTrackParam &track)
Definition: CbmKFTrack.cxx:34
CbmStsKFTrackFitter::fIsInitialised
Bool_t fIsInitialised
Definition: CbmStsKFTrackFitter.h:75
CbmStsKFTrackFitter.h
CbmKFTrack::fHits
std::vector< CbmKFHit * > fHits
Definition: CbmKFTrack.h:29
CbmKFVertex
Definition: CbmKFVertex.h:6
CbmKFTrackInterface::Extrapolate
Int_t Extrapolate(Double_t z, Double_t *QP0=0)
Access to i-th hit.
Definition: CbmKFTrackInterface.cxx:39