Belle II Software  release-05-01-25
VertexFitKFit.cc
1 /**************************************************************************
2  * Copyright(C) 1997 - J. Tanaka *
3  * *
4  * Author: J. Tanaka *
5  * Contributor: J. Tanaka and *
6  * conversion to Belle II structure by T. Higuchi *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 
12 #include <algorithm>
13 #include <cstdio>
14 
15 #include <TMatrixFSym.h>
16 
17 #include <analysis/VertexFitting/KFit/MakeMotherKFit.h>
18 #include <analysis/VertexFitting/KFit/VertexFitKFit.h>
19 #include <analysis/utility/CLHEPToROOT.h>
20 
21 
22 using namespace std;
23 using namespace Belle2;
24 using namespace Belle2::analysis;
25 using namespace CLHEP;
26 
27 VertexFitKFit::VertexFitKFit():
28  m_BeforeVertex(HepPoint3D(0, 0, 0)),
29  m_AfterVertexError(HepSymMatrix(3, 0)),
30  m_BeamError(HepSymMatrix(3, 0))
31 {
32  m_CorrelationMode = false;
33  m_FlagFitted = false;
34  m_FlagKnownVertex = false;
35  m_FlagBeam = false;
36  m_FlagTube = false;
37  m_iTrackTube = -1;
38  m_CHIsqVertex = 0;
39  m_NecessaryTrackCount = 2;
40  m_V_E = HepMatrix(3, 3, 0);
41  m_v = HepMatrix(3, 1, 0);
42  m_v_a = HepMatrix(3, 1, 0);
43 
44  fill_n(m_EachCHIsq, KFitConst::kMaxTrackCount2, 0);
45 }
46 
47 VertexFitKFit::~VertexFitKFit() = default;
48 
49 
50 enum KFitError::ECode
51 VertexFitKFit::setInitialVertex(const HepPoint3D& v) {
52  m_BeforeVertex = v;
53 
54  return m_ErrorCode = KFitError::kNoError;
55 }
56 
57 enum KFitError::ECode VertexFitKFit::setInitialVertex(const TVector3& v)
58 {
59  m_BeforeVertex = HepPoint3D(v.X(), v.Y(), v.Z());
60  m_ErrorCode = KFitError::kNoError;
61  return m_ErrorCode;
62 }
63 
64 enum KFitError::ECode
65 VertexFitKFit::setIpProfile(const HepPoint3D& ip, const HepSymMatrix& ipe) {
66  if (m_FlagTube)
67  {
68  char buf[1024];
69  sprintf(buf, "%s:%s(): already constrained to IPtube", __FILE__, __func__);
70  B2FATAL(buf);
71  }
72 
73  m_FlagBeam = true;
74  m_BeforeVertex = ip;
75  m_BeamError = ipe;
76 
77  return m_ErrorCode = KFitError::kNoError;
78 }
79 
80 
81 enum KFitError::ECode
82 VertexFitKFit::setIpTubeProfile(const HepLorentzVector& p, const HepPoint3D& x, const HepSymMatrix& e, const double q) {
83  if (m_FlagBeam)
84  {
85  char buf[1024];
86  sprintf(buf, "%s:%s(): already constrained to IP", __FILE__, __func__);
87  B2FATAL(buf);
88  }
89 
90  m_FlagTube = true;
91  m_TubeTrack = KFitTrack(p, x, e, q);
92 
93  return m_ErrorCode = KFitError::kNoError;
94 }
95 
96 
97 enum KFitError::ECode
98 VertexFitKFit::setKnownVertex(const bool flag) {
99  m_FlagKnownVertex = flag;
100 
101  return m_ErrorCode = KFitError::kNoError;
102 }
103 
104 
105 enum KFitError::ECode
106 VertexFitKFit::setCorrelationMode(const bool m) {
107  m_CorrelationMode = m;
108 
109  return m_ErrorCode = KFitError::kNoError;
110 }
111 
112 
113 const HepPoint3D
114 VertexFitKFit::getVertex(const int flag) const
115 {
116  if (flag == KFitConst::kAfterFit && !isFitted()) return HepPoint3D();
117 
118  switch (flag) {
119  case KFitConst::kBeforeFit:
120  return m_BeforeVertex;
121 
122  case KFitConst::kAfterFit:
123  return m_AfterVertex;
124 
125  default:
126  KFitError::displayError(__FILE__, __LINE__, __func__, KFitError::kOutOfRange);
127  return HepPoint3D();
128  }
129 }
130 
131 
132 const HepSymMatrix
133 VertexFitKFit::getVertexError() const
134 {
135  return m_AfterVertexError;
136 }
137 
138 
139 double
140 VertexFitKFit::getCHIsq() const
141 {
142  return KFitBase::getCHIsq();
143 }
144 
145 
146 double
147 VertexFitKFit::getCHIsqVertex() const
148 {
149  // only m_FlagBeam = 1
150  return m_CHIsqVertex;
151 }
152 
153 
154 const HepMatrix
155 VertexFitKFit::getTrackVertexError(const int id) const
156 {
157  if (!isTrackIDInRange(id)) return HepMatrix(3, KFitConst::kNumber7, 0);
158  return m_AfterTrackVertexError[id];
159 }
160 
161 
162 double
163 VertexFitKFit::getTrackCHIsq(const int id) const
164 {
165  if (!isTrackIDInRange(id)) return -1;
166 
167  if (!m_FlagKnownVertex && !m_FlagBeam && m_CorrelationMode) {
168  return KFitBase::getTrackCHIsq(id);
169  }
170 
171  return m_EachCHIsq[id];
172 }
173 
174 
175 double
176 VertexFitKFit::getTrackPartCHIsq() const
177 {
178  if (m_ErrorCode != KFitError::kNoError) {
179  KFitError::displayError(__FILE__, __LINE__, __func__, m_ErrorCode);
180  return -1;
181  }
182 
183  if (m_TrackCount == 0) {
184  KFitError::displayError(__FILE__, __LINE__, __func__, KFitError::kBadTrackSize);
185  return -1;
186  }
187 
188  double chisq = 0.0;
189  for (int i = 0; i < m_TrackCount; i++) {
190  const double i_chisq = this->getTrackCHIsq(i);
191  chisq += i_chisq;
192  }
193 
194  return chisq;
195 }
196 
197 
198 int
199 VertexFitKFit::getTrackPartNDF() const
200 {
201  if (m_ErrorCode != KFitError::kNoError) {
202  KFitError::displayError(__FILE__, __LINE__, __func__, m_ErrorCode);
203  return 0;
204  }
205 
206  if (m_TrackCount == 0) {
207  KFitError::displayError(__FILE__, __LINE__, __func__, KFitError::kBadTrackSize);
208  return 0;
209  }
210 
211  return m_TrackCount * 2 - 2;
212 }
213 
214 
215 enum KFitError::ECode
216 VertexFitKFit::doFit() {
217  if (m_FlagTube) this->appendTube();
218 
219  if (m_FlagBeam) m_ErrorCode = this->doFit4();
220  else if (m_FlagKnownVertex) m_ErrorCode = this->doFit5();
221  else if (m_CorrelationMode) m_ErrorCode = this->doFit2();
222  else
223  m_ErrorCode = this->doFit3();
224 
225  const enum KFitError::ECode tmp_ErrorCode = m_ErrorCode;
226 
227  if (m_FlagTube) this->deleteTube();
228 
229  if (tmp_ErrorCode == KFitError::kNoError) m_FlagFitted = true;
230 
231  return m_ErrorCode = tmp_ErrorCode;
232 }
233 
234 
235 enum KFitError::ECode
236 VertexFitKFit::doFit3() {
237  // use small Matrix --> No Correlation
238  if (m_ErrorCode != KFitError::kNoError) return m_ErrorCode;
239 
240  if (m_TrackCount < m_NecessaryTrackCount)
241  {
242  m_ErrorCode = KFitError::kBadTrackSize;
243  KFitError::displayError(__FILE__, __LINE__, __func__, m_ErrorCode);
244  return m_ErrorCode;
245  }
246 
247  if (prepareInputMatrix() != KFitError::kNoError) return m_ErrorCode;
248  if (calculateNDF() != KFitError::kNoError) return m_ErrorCode;
249 
250 
251  double chisq = 0;
252  double tmp2_chisq = KFitConst::kInitialCHIsq;
253  int err_inverse = 0;
254 
255  m_al_a = m_al_0;
256  HepMatrix tmp_al_a(m_al_a);
257 
258  HepMatrix tmp_D(m_D), tmp_E(m_E);
259  HepMatrix tmp_V_D(m_V_D), tmp_V_E(m_V_E);
260  HepMatrix tmp_lam0(m_lam0), tmp_v_a(m_v_a);
261 
262  HepMatrix tmp2_D(m_D), tmp2_E(m_E);
263  HepMatrix tmp2_V_D(m_V_D), tmp2_V_E(m_V_E);
264  HepMatrix tmp2_lam0(m_lam0), tmp2_v_a(m_v_a);
265 
266  double tmp_each_chisq[KFitConst::kMaxTrackCount2];
267  double tmp2_each_chisq[KFitConst::kMaxTrackCount2];
268 
269  for (int j = 0; j < KFitConst::kMaxIterationCount; j++) // j'th loop start
270  {
271 
272  double tmp_chisq = KFitConst::kInitialCHIsq;
273 
274  for (int i = 0; i < KFitConst::kMaxIterationCount; i++) { // i'th loop start
275 
276  if (prepareInputSubMatrix() != KFitError::kNoError) return m_ErrorCode;
277  if (makeCoreMatrix() != KFitError::kNoError) return m_ErrorCode;
278 
279  HepMatrix tV_Ein(3, 3, 0);
280  chisq = 0;
281 
282  for (int k = 0; k < m_TrackCount; k++) { // k'th loop start
283 
284  HepMatrix tD = m_D.sub(2 * k + 1, 2 * (k + 1), KFitConst::kNumber6 * k + 1, KFitConst::kNumber6 * (k + 1)); // 2x6
285  HepMatrix tV_D = ((m_V_al_0.sub(KFitConst::kNumber6 * k + 1,
286  (int)(KFitConst::kNumber6 * (k + 1)))).similarity(tD)).inverse(err_inverse); // 2x2
287  if (err_inverse) {
288  m_ErrorCode = KFitError::kCannotGetMatrixInverse;
289  KFitError::displayError(__FILE__, __LINE__, __func__, m_ErrorCode);
290  return m_ErrorCode;
291  }
292 
293  m_V_D.sub(2 * k + 1, 2 * k + 1, tV_D);
294  HepMatrix tE = m_E.sub(2 * k + 1, 2 * (k + 1), 1, 3); // 2x3
295  tV_Ein += (tE.T()) * tV_D * tE; // 3x3
296  HepMatrix tDeltaAl = (m_al_0 - m_al_1).sub(KFitConst::kNumber6 * k + 1, KFitConst::kNumber6 * (k + 1), 1, 1); // 6x1
297  HepMatrix td = m_d.sub(2 * k + 1, 2 * (k + 1), 1, 1); // 2x1
298  HepMatrix tlam0 = tV_D * (tD * tDeltaAl + td); // 2x2x(2x6x6x1+2x1) = 2x1
299  m_lam0.sub(2 * k + 1, 1, tlam0);
300  m_EachCHIsq[k] = ((tlam0.T()) * (tD * tDeltaAl + tE * (m_v - m_v_a) + td))(1, 1); // 1x2x(2x6x6x1+2x3x3x1+2x1)
301  chisq += m_EachCHIsq[k];
302  } // k'th loop over
303 
304  m_V_E = tV_Ein.inverse(err_inverse);
305  if (err_inverse) {
306  m_ErrorCode = KFitError::kCannotGetMatrixInverse;
307  KFitError::displayError(__FILE__, __LINE__, __func__, m_ErrorCode);
308  return m_ErrorCode;
309  }
310 
311  m_v_a = m_v_a - m_V_E * (m_E.T()) * m_lam0;
312 
313  if (tmp_chisq <= chisq) {
314  if (i == 0) {
315  m_ErrorCode = KFitError::kBadInitialCHIsq;
316  } else {
317  for (int k = 0; k < m_TrackCount; k++) m_EachCHIsq[k] = tmp_each_chisq[k];
318  chisq = tmp_chisq;
319  m_v_a = tmp_v_a;
320  m_V_E = tmp_V_E;
321  m_V_D = tmp_V_D;
322  m_lam0 = tmp_lam0;
323  m_E = tmp_E;
324  m_D = tmp_D;
325  }
326  break;
327  } else {
328  for (int k = 0; k < m_TrackCount; k++) tmp_each_chisq[k] = m_EachCHIsq[k];
329  tmp_chisq = chisq;
330  tmp_v_a = m_v_a;
331  tmp_V_E = m_V_E;
332  tmp_V_D = m_V_D;
333  tmp_lam0 = m_lam0;
334  tmp_E = m_E;
335  tmp_D = m_D;
336  if (i == KFitConst::kMaxIterationCount - 1) {
337  m_FlagOverIteration = true;
338  }
339  }
340  } // i'th loop over
341 
342  m_al_a = m_al_1;
343  m_lam = m_lam0 - m_V_D * m_E * m_V_E * (m_E.T()) * m_lam0;
344  m_al_1 = m_al_0 - m_V_al_0 * (m_D.T()) * m_lam;
345 
346  if (j == 0) {
347 
348  for (int k = 0; k < m_TrackCount; k++) tmp2_each_chisq[k] = m_EachCHIsq[k];
349  tmp2_chisq = chisq;
350  tmp2_v_a = m_v_a;
351  tmp2_V_E = m_V_E;
352  tmp2_V_D = m_V_D;
353  tmp2_lam0 = m_lam0;
354  tmp2_E = m_E;
355  tmp2_D = m_D;
356  tmp_al_a = m_al_a;
357 
358  } else {
359 
360  if (tmp2_chisq <= chisq) {
361  for (int k = 0; k < m_TrackCount; k++) m_EachCHIsq[k] = tmp2_each_chisq[k];
362  chisq = tmp2_chisq;
363  m_v_a = tmp2_v_a;
364  m_V_E = tmp2_V_E;
365  m_V_D = tmp2_V_D;
366  m_lam0 = tmp2_lam0;
367  m_E = tmp2_E;
368  m_D = tmp2_D;
369  m_al_a = tmp_al_a;
370  break;
371  } else {
372  for (int k = 0; k < m_TrackCount; k++) tmp2_each_chisq[k] = m_EachCHIsq[k];
373  tmp2_chisq = chisq;
374  tmp2_v_a = m_v_a;
375  tmp2_V_E = m_V_E;
376  tmp2_V_D = m_V_D;
377  tmp2_lam0 = m_lam0;
378  tmp2_E = m_E;
379  tmp2_D = m_D;
380  tmp_al_a = m_al_a;
381  }
382  }
383  } // j'th loop over
384 
385 
386  if (m_ErrorCode != KFitError::kNoError) return m_ErrorCode;
387 
388  m_lam = m_lam0 - m_V_D * m_E * m_V_E * (m_E.T()) * m_lam0;
389  m_al_1 = m_al_0 - m_V_al_0 * (m_D.T()) * m_lam;
390  m_V_Dt = m_V_D - m_V_D * m_E * m_V_E * (m_E.T()) * m_V_D;
391  m_V_al_1 = m_V_al_0 - m_V_al_0 * (m_D.T()) * m_V_Dt * m_D * m_V_al_0;
392  m_Cov_v_al_1 = -m_V_E * (m_E.T()) * m_V_D * m_D * m_V_al_0;
393 
394  if (prepareOutputMatrix() != KFitError::kNoError) return m_ErrorCode;
395 
396  m_CHIsq = chisq;
397 
398  return m_ErrorCode = KFitError::kNoError;
399 }
400 
401 
402 enum KFitError::ECode
403 VertexFitKFit::doFit4() {
404  // included beam position constraint (only no correlation)
405  if (m_ErrorCode != KFitError::kNoError) return m_ErrorCode;
406 
407  if (m_TrackCount < m_NecessaryTrackCount - 1)
408  {
409  m_ErrorCode = KFitError::kBadTrackSize;
410  KFitError::displayError(__FILE__, __LINE__, __func__, m_ErrorCode);
411  return m_ErrorCode;
412  }
413 
414  if (prepareInputMatrix() != KFitError::kNoError) return m_ErrorCode;
415  if (calculateNDF() != KFitError::kNoError) return m_ErrorCode;
416 
417 
418  double chisq = 0;
419  double tmp_chisq = KFitConst::kInitialCHIsq;
420  int err_inverse = 0;
421 
422  m_al_a = m_al_0;
423  HepMatrix tmp_al_a(m_al_a);
424 
425  HepMatrix tmp_D(m_D), tmp_E(m_E);
426  HepMatrix tmp_lam(m_lam);
427 
428  // vertex
429  m_v[0][0] = m_BeforeVertex.x();
430  m_v[1][0] = m_BeforeVertex.y();
431  m_v[2][0] = m_BeforeVertex.z();
432 
433  double tmp_each_chisq[KFitConst::kMaxTrackCount2];
434  double tmp_vertex_chisq = 1.e+30; // An init-value is not needed but the C++ complier requires the init-value.
435 
436  // to avoid overestimation of vertex-z error.
437  bool it_flag = false;
438 
439  for (int i = 0; i < KFitConst::kMaxIterationCount ; i++)
440  {
441 
442  if (makeCoreMatrix() != KFitError::kNoError) return m_ErrorCode;
443 
444  chisq = 0;
445 
446  HepMatrix tV_Dtin = m_V_al_0.similarity(m_D) + m_BeamError.similarity(m_E);
447  HepMatrix tV_Dt = tV_Dtin.inverse(err_inverse);
448  if (err_inverse) {
449  m_ErrorCode = KFitError::kCannotGetMatrixInverse;
450  KFitError::displayError(__FILE__, __LINE__, __func__, m_ErrorCode);
451  return m_ErrorCode;
452  }
453  m_lam = tV_Dt * (m_D * (m_al_0 - m_al_1) + m_E * (m_v - m_v_a) + m_d); // (2*nTrk)x1
454  for (int k = 0; k < m_TrackCount; k++) {
455  HepMatrix tD = m_D.sub(2 * k + 1, 2 * (k + 1), KFitConst::kNumber6 * k + 1, KFitConst::kNumber6 * (k + 1)); // 2x6
456  HepMatrix tDeltaAl = (m_al_0 - m_al_1).sub(KFitConst::kNumber6 * k + 1, KFitConst::kNumber6 * (k + 1), 1, 1); // 6x1
457  HepMatrix td = m_d.sub(2 * k + 1, 2 * (k + 1), 1, 1); // 2x1
458  HepMatrix tE = m_E.sub(2 * k + 1, 2 * (k + 1), 1, 3); // 2x3
459  chisq += ((m_lam.sub(2 * k + 1, 2 * (k + 1), 1, 1).T()) * (tD * tDeltaAl + tE * (m_v - m_v_a) + td))(1,
460  1); // 1x2x(2x6x6x1+2x3x3x1+2x1)
461  m_EachCHIsq[k] = (m_lam.sub(2 * k + 1, 2 * (k + 1), 1, 1).T() * tD * m_V_al_0.sub(KFitConst::kNumber6 * k + 1,
462  KFitConst::kNumber6 * (k + 1)) * (tD.T()) * m_lam.sub(2 * k + 1, 2 * (k + 1), 1, 1))(1, 1);
463  }
464 
465  m_CHIsqVertex = (m_lam.T() * m_E * m_BeamError * (m_E.T()) * m_lam)(1, 1);
466  m_al_a = m_al_1;
467  m_v_a = m_v - m_BeamError * (m_E.T()) * m_lam;
468  m_al_1 = m_al_0 - m_V_al_0 * (m_D.T()) * m_lam;
469 
470  if (tmp_chisq <= chisq && it_flag) {
471  if (i == 0) {
472  m_ErrorCode = KFitError::kBadInitialCHIsq;
473  } else {
474  for (int k = 0; k < m_TrackCount; k++) m_EachCHIsq[k] = tmp_each_chisq[k];
475  m_CHIsqVertex = tmp_vertex_chisq;
476  chisq = tmp_chisq;
477  m_lam = tmp_lam;
478  m_E = tmp_E;
479  m_D = tmp_D;
480  m_al_a = tmp_al_a;
481  }
482  break;
483  } else {
484  if (tmp_chisq <= chisq) it_flag = true;
485  for (int k = 0; k < m_TrackCount; k++) tmp_each_chisq[k] = m_EachCHIsq[k];
486  tmp_vertex_chisq = m_CHIsqVertex;
487  tmp_chisq = chisq;
488  tmp_lam = m_lam;
489  tmp_E = m_E;
490  tmp_D = m_D;
491  tmp_al_a = m_al_a;
492  if (i == KFitConst::kMaxIterationCount - 1) {
493  m_FlagOverIteration = true;
494  }
495  }
496  }
497 
498 
499  if (m_ErrorCode != KFitError::kNoError) return m_ErrorCode;
500 
501  m_al_1 = m_al_0 - m_V_al_0 * (m_D.T()) * m_lam;
502  m_v_a = m_v - m_BeamError * (m_E.T()) * m_lam;
503  HepMatrix tV_Dtin = m_V_al_0.similarity(m_D) + m_BeamError.similarity(m_E);
504  m_V_Dt = tV_Dtin.inverse(err_inverse);
505  if (err_inverse)
506  {
507  m_ErrorCode = KFitError::kCannotGetMatrixInverse;
508  KFitError::displayError(__FILE__, __LINE__, __func__, m_ErrorCode);
509  return m_ErrorCode;
510  }
511 
512  m_V_al_1 = m_V_al_0 - m_V_al_0 * (m_D.T()) * m_V_Dt * m_D * m_V_al_0;
513  m_Cov_v_al_1 = -m_BeamError * (m_E.T()) * m_V_Dt * m_D * m_V_al_0;
514  // m_V_v is m_V_E
515  // --> need to replace m_V_E for my implimentaion.
516  m_V_E = m_BeamError - m_BeamError * (m_E.T()) * m_V_Dt * m_E * m_BeamError;
517 
518  if (prepareOutputMatrix() != KFitError::kNoError) return m_ErrorCode;
519 
520  m_CHIsq = chisq;
521 
522  return m_ErrorCode = KFitError::kNoError;
523 }
524 
525 
526 enum KFitError::ECode
527 VertexFitKFit::doFit5() {
528  // known vertex --> do not find vertex. (only no correlation)
529  if (m_ErrorCode != KFitError::kNoError) return m_ErrorCode;
530 
531  if (m_TrackCount < m_NecessaryTrackCount - 1)
532  {
533  m_ErrorCode = KFitError::kBadTrackSize;
534  KFitError::displayError(__FILE__, __LINE__, __func__, m_ErrorCode);
535  return m_ErrorCode;
536  }
537 
538  if (prepareInputMatrix() != KFitError::kNoError) return m_ErrorCode;
539  if (calculateNDF() != KFitError::kNoError) return m_ErrorCode;
540 
541 
542  double chisq = 0;
543  double tmp_chisq = KFitConst::kInitialCHIsq;
544  int err_inverse = 0;
545 
546  m_al_a = m_al_0;
547  HepMatrix tmp_al_a(m_al_a);
548 
549  HepMatrix tmp_al_0(m_al_1);
550  HepMatrix tmp_V_al_0(m_V_al_1);
551 
552  double tmp_each_chisq[KFitConst::kMaxTrackCount2];
553 
554  for (int i = 0; i < KFitConst::kMaxIterationCount; i++)
555  {
556 
557  if (makeCoreMatrix() != KFitError::kNoError) return m_ErrorCode;
558 
559  chisq = 0;
560 
561  for (int k = 0; k < m_TrackCount; k++) {
562  HepMatrix tD = m_D.sub(2 * k + 1, 2 * (k + 1), KFitConst::kNumber6 * k + 1, KFitConst::kNumber6 * (k + 1)); // 2x6
563  HepMatrix tV_D = ((m_V_al_0.sub(KFitConst::kNumber6 * k + 1,
564  (int)(KFitConst::kNumber6 * (k + 1)))).similarity(tD)).inverse(err_inverse); // 2x2
565  if (err_inverse) {
566  m_ErrorCode = KFitError::kCannotGetMatrixInverse;
567  KFitError::displayError(__FILE__, __LINE__, __func__, m_ErrorCode);
568  return m_ErrorCode;
569  }
570  m_V_D.sub(2 * k + 1, 2 * k + 1, tV_D);
571 
572  HepMatrix tDeltaAl = (m_al_0 - m_al_1).sub(KFitConst::kNumber6 * k + 1, KFitConst::kNumber6 * (k + 1), 1, 1); // 6x1
573  HepMatrix td = m_d.sub(2 * k + 1, 2 * (k + 1), 1, 1); // 2x1
574  HepMatrix tlam = tV_D * (tD * tDeltaAl + td); // 2x2x(2x6x6x1+2x1) = 2x1
575  m_lam.sub(2 * k + 1, 1, tlam);
576  m_EachCHIsq[k] = ((tlam.T()) * (tD * tDeltaAl + td))(1, 1); // 1x2x(2x6x6x1+2x1)
577  chisq += m_EachCHIsq[k];
578  }
579 
580  m_al_a = m_al_1;
581  m_al_1 = m_al_0 - m_V_al_0 * (m_D.T()) * m_lam;
582  m_V_al_1 = m_V_al_0 - m_V_al_0 * (m_D.T()) * m_V_D * m_D * m_V_al_0;
583 
584  if (tmp_chisq <= chisq) {
585  if (i == 0) {
586  m_ErrorCode = KFitError::kBadInitialCHIsq;
587  } else {
588  for (int k = 0; k < m_TrackCount; k++) m_EachCHIsq[k] = tmp_each_chisq[k];
589  chisq = tmp_chisq;
590  m_al_1 = tmp_al_0;
591  m_V_al_1 = tmp_V_al_0;
592  m_al_a = tmp_al_a;
593  }
594  break;
595  } else {
596  for (int k = 0; k < m_TrackCount; k++) tmp_each_chisq[k] = m_EachCHIsq[k];
597  tmp_chisq = chisq;
598  tmp_al_0 = m_al_1;
599  tmp_V_al_0 = m_V_al_1;
600  tmp_al_a = m_al_a;
601  if (i == KFitConst::kMaxIterationCount - 1) {
602  m_FlagOverIteration = true;
603  }
604  }
605  }
606 
607 
608  if (m_ErrorCode != KFitError::kNoError) return m_ErrorCode;
609 
610  if (prepareOutputMatrix() != KFitError::kNoError) return m_ErrorCode;
611 
612  m_CHIsq = chisq;
613 
614  return m_ErrorCode = KFitError::kNoError;
615 }
616 
617 
618 enum KFitError::ECode
619 VertexFitKFit::prepareInputMatrix() {
620  if (!m_CorrelationMode || m_FlagBeam || m_FlagKnownVertex)
621  {
622  if (m_TrackCount > KFitConst::kMaxTrackCount2) {
623  m_ErrorCode = KFitError::kBadTrackSize;
624  KFitError::displayError(__FILE__, __LINE__, __func__, m_ErrorCode);
625  return m_ErrorCode;
626  }
627  } else {
628  if (m_TrackCount > KFitConst::kMaxTrackCount)
629  {
630  m_ErrorCode = KFitError::kBadTrackSize;
631  KFitError::displayError(__FILE__, __LINE__, __func__, m_ErrorCode);
632  return m_ErrorCode;
633  }
634  }
635 
636 
637  int index = 0;
638  HepMatrix tmp_al_0(KFitConst::kNumber6 * m_TrackCount, 1, 0);
639  HepSymMatrix tmp_V_al_0(KFitConst::kNumber6 * m_TrackCount, 0);
640  HepMatrix tmp_property(m_TrackCount, 3, 0);
641 
642 
643  for (auto& track : m_Tracks)
644  {
645  // momentum x,y,z and position x,y,z
646  for (int j = 0; j < KFitConst::kNumber6; j++)
647  tmp_al_0[index * KFitConst::kNumber6 + j][0] = track.getFitParameter(j, KFitConst::kBeforeFit);
648  // these error
649  tmp_V_al_0.sub(index * KFitConst::kNumber6 + 1, track.getFitError(KFitConst::kBeforeFit));
650  // charge , mass , a
651  tmp_property[index][0] = track.getCharge();
652  tmp_property[index][1] = track.getMass();
653  const double c = KFitConst::kLightSpeed; // C++ bug?
654  // tmp_property[index][2] = -KFitConst::kLightSpeed * m_MagneticField * it->getCharge();
655  tmp_property[index][2] = -c * m_MagneticField * track.getCharge();
656  index++;
657  }
658 
659  // error between tarck and track
660  m_V_al_0 = tmp_V_al_0;
661  if (!m_FlagKnownVertex && !m_FlagBeam && m_CorrelationMode)
662  {
663  if (m_FlagCorrelation) {
664  this->prepareCorrelation();
665  if (m_ErrorCode != KFitError::kNoError) {
666  KFitError::displayError(__FILE__, __LINE__, __func__, m_ErrorCode);
667  return m_ErrorCode;
668  }
669  }
670  }
671 
672  // vertex
673  m_v_a[0][0] = m_BeforeVertex.x();
674  m_v_a[1][0] = m_BeforeVertex.y();
675  m_v_a[2][0] = m_BeforeVertex.z();
676 
677  // set member matrix
678  m_al_0 = tmp_al_0;
679  m_al_1 = m_al_0;
680  m_property = tmp_property;
681 
682  // define size of matrix
683  m_V_al_1 = HepMatrix(KFitConst::kNumber6 * m_TrackCount, KFitConst::kNumber6 * m_TrackCount, 0);
684  m_D = m_V_al_1.sub(1, m_TrackCount * 2, 1, KFitConst::kNumber6 * m_TrackCount);
685  m_E = m_V_al_1.sub(1, m_TrackCount * 2, 1, 3);
686  m_d = m_V_al_1.sub(1, m_TrackCount * 2, 1, 1);
687  m_V_D = m_V_al_1.sub(1, m_TrackCount * 2, 1, m_TrackCount * 2);
688  m_lam = m_V_al_1.sub(1, m_TrackCount * 2, 1, 1);
689  m_lam0 = m_V_al_1.sub(1, m_TrackCount * 2, 1, 1);
690  m_V_Dt = m_V_al_1.sub(1, m_TrackCount * 2, 1, m_TrackCount * 2);
691  m_Cov_v_al_1 = m_V_al_1.sub(1, 3, 1, KFitConst::kNumber6 * m_TrackCount);
692 
693  return m_ErrorCode = KFitError::kNoError;
694 }
695 
696 
697 enum KFitError::ECode
698 VertexFitKFit::prepareInputSubMatrix() {
699  // vertex
700  for (int i = 0; i < 3; i++)
701  {
702  m_v[i][0] = m_v_a[i][0];
703  }
704  return m_ErrorCode = KFitError::kNoError;
705 }
706 
707 
708 enum KFitError::ECode
709 VertexFitKFit::prepareOutputMatrix() {
710  Hep3Vector h3v;
711  unsigned index = 0;
712 
713  for (auto& pdata : m_Tracks)
714  {
715  // tracks
716  // momentum
717  h3v.setX(m_al_1[index * KFitConst::kNumber6 + 0][0]);
718  h3v.setY(m_al_1[index * KFitConst::kNumber6 + 1][0]);
719  h3v.setZ(m_al_1[index * KFitConst::kNumber6 + 2][0]);
720  pdata.setMomentum(HepLorentzVector(h3v, sqrt(h3v.mag2() + pdata.getMass()*pdata.getMass())), KFitConst::kAfterFit);
721  // position
722  pdata.setPosition(HepPoint3D(
723  m_al_1[index * KFitConst::kNumber6 + 3][0],
724  m_al_1[index * KFitConst::kNumber6 + 4][0],
725  m_al_1[index * KFitConst::kNumber6 + 5][0]), KFitConst::kAfterFit);
726  // error of the tracks
727  pdata.setError(makeError1(pdata.getMomentum(),
728  m_V_al_1.sub(
729  index * KFitConst::kNumber6 + 1,
730  (index + 1)*KFitConst::kNumber6,
731  index * KFitConst::kNumber6 + 1,
732  (index + 1)*KFitConst::kNumber6)),
733  KFitConst::kAfterFit);
734  if (m_ErrorCode != KFitError::kNoError) break;
735  index++;
736  }
737 
738  // vertex
739  m_AfterVertex.setX(m_v_a[0][0]);
740  m_AfterVertex.setY(m_v_a[1][0]);
741  m_AfterVertex.setZ(m_v_a[2][0]);
742 
743  // error of the vertex
744  for (int i = 0; i < 3; i++) for (int j = i; j < 3; j++)
745  {
746  m_AfterVertexError[i][j] = m_V_E[i][j];
747  }
748 
749  // error between vertex and tracks
750  for (int i = 0; i < m_TrackCount; i++)
751  {
752  HepMatrix hm(3, KFitConst::kNumber6, 0);
753  for (int j = 0; j < 3; j++) for (int k = 0; k < KFitConst::kNumber6; k++) {
754  hm[j][k] = m_Cov_v_al_1[j][KFitConst::kNumber6 * i + k];
755  }
756  m_AfterTrackVertexError.push_back(makeError2(m_Tracks[i].getMomentum(), hm));
757  }
758 
759  return m_ErrorCode = KFitError::kNoError;
760 }
761 
762 
763 enum KFitError::ECode
764 VertexFitKFit::makeCoreMatrix() {
765  // vertex fit
766  for (int i = 0; i < m_TrackCount; i++)
767  {
768  double B, S, U;
769  double sininv, sqrtag;
770 
771  double px = m_al_1[i * KFitConst::kNumber6 + 0][0];
772  double py = m_al_1[i * KFitConst::kNumber6 + 1][0];
773  double pz = m_al_1[i * KFitConst::kNumber6 + 2][0];
774  double x = m_al_1[i * KFitConst::kNumber6 + 3][0];
775  double y = m_al_1[i * KFitConst::kNumber6 + 4][0];
776  double z = m_al_1[i * KFitConst::kNumber6 + 5][0];
777  double a = m_property[i][2];
778 
779  double pt = sqrt(px * px + py * py);
780  if (pt == 0) {
781  m_ErrorCode = KFitError::kDivisionByZero;
782  KFitError::displayError(__FILE__, __LINE__, __func__, m_ErrorCode);
783  return m_ErrorCode;
784  }
785 
786  double invPt = 1. / pt;
787  double invPt2 = invPt * invPt;
788  double dlx = m_v_a[0][0] - x;
789  double dly = m_v_a[1][0] - y;
790  double dlz = m_v_a[2][0] - z;
791  double a1 = -dlx * py + dly * px;
792  double a2 = dlx * px + dly * py;
793  double r2d2 = dlx * dlx + dly * dly;
794  double Rx = dlx - 2.*px * a2 * invPt2;
795  double Ry = dly - 2.*py * a2 * invPt2;
796 
797  if (a != 0) { // charged
798 
799  B = a * a2 * invPt2;
800  if (fabs(B) > 1) {
801  m_ErrorCode = KFitError::kCannotGetARCSIN;
802  B2DEBUG(10, "KFitError: Cannot calculate arcsin");
803  //KFitError::displayError(__FILE__, __LINE__, __func__, m_ErrorCode);
804  return m_ErrorCode;
805  }
806  // sin^(-1)(B)
807  sininv = asin(B);
808  double tmp0 = 1.0 - B * B;
809  if (tmp0 == 0) {
810  m_ErrorCode = KFitError::kDivisionByZero;
811  KFitError::displayError(__FILE__, __LINE__, __func__, m_ErrorCode);
812  return m_ErrorCode;
813  }
814  // 1/sqrt(1-B^2)
815  sqrtag = 1.0 / sqrt(tmp0);
816  S = sqrtag * invPt2;
817  U = dlz - pz * sininv / a;
818 
819  } else { // neutral
820  B = 0.0;
821  sininv = 0.0;
822  sqrtag = 1.0;
823  S = invPt2;
824  U = dlz - pz * a2 * invPt2;
825  }
826 
827  // d
828  m_d[i * 2 + 0][0] = a1 - 0.5 * a * r2d2;
829  m_d[i * 2 + 1][0] = U * pt;
830 
831  // D
832  m_D[i * 2 + 0][i * KFitConst::kNumber6 + 0] = dly;
833  m_D[i * 2 + 0][i * KFitConst::kNumber6 + 1] = -dlx;
834  m_D[i * 2 + 0][i * KFitConst::kNumber6 + 2] = 0.0;
835  m_D[i * 2 + 0][i * KFitConst::kNumber6 + 3] = py + a * dlx;
836  m_D[i * 2 + 0][i * KFitConst::kNumber6 + 4] = -px + a * dly;
837  m_D[i * 2 + 0][i * KFitConst::kNumber6 + 5] = 0.0;
838  m_D[i * 2 + 1][i * KFitConst::kNumber6 + 0] = -pz * pt * S * Rx + U * px * invPt;
839  m_D[i * 2 + 1][i * KFitConst::kNumber6 + 1] = -pz * pt * S * Ry + U * py * invPt;
840  m_D[i * 2 + 1][i * KFitConst::kNumber6 + 2] = a != 0 ? -sininv * pt / a : -a2 * invPt;
841  m_D[i * 2 + 1][i * KFitConst::kNumber6 + 3] = px * pz * pt * S;
842  m_D[i * 2 + 1][i * KFitConst::kNumber6 + 4] = py * pz * pt * S;
843  m_D[i * 2 + 1][i * KFitConst::kNumber6 + 5] = -pt;
844 
845  // E
846  m_E[i * 2 + 0][0] = -py - a * dlx;
847  m_E[i * 2 + 0][1] = px - a * dly;
848  m_E[i * 2 + 0][2] = 0.0;
849  m_E[i * 2 + 1][0] = -px * pz * pt * S;
850  m_E[i * 2 + 1][1] = -py * pz * pt * S;
851  m_E[i * 2 + 1][2] = pt;
852  }
853 
854  return m_ErrorCode = KFitError::kNoError;
855 }
856 
857 
858 enum KFitError::ECode
859 VertexFitKFit::calculateNDF() {
860  if (m_FlagBeam) m_NDF = 2 * m_TrackCount;
861  else if (m_FlagTube) m_NDF = 2 * (m_TrackCount - 1) - 1;
862  else if (m_FlagKnownVertex) m_NDF = 2 * m_TrackCount;
863  else m_NDF = 2 * m_TrackCount - 3;
864 
865  return m_ErrorCode = KFitError::kNoError;
866 }
867 
868 
869 enum KFitError::ECode
870 VertexFitKFit::appendTube() {
871  if (!m_FlagTube) return m_ErrorCode = KFitError::kNoError;
872 
873  if (m_iTrackTube != -1)
874  {
875  char buf[1024];
876  sprintf(buf, "%s:%s(): internal error; duplicated appendTube() call?", __FILE__, __func__);
877  B2FATAL(buf);
878  }
879 
880  m_Tracks.push_back(m_TubeTrack);
881  m_TrackCount = m_Tracks.size();
882  m_iTrackTube = m_TrackCount;
883 
884  return m_ErrorCode = KFitError::kNoError;
885 }
886 
887 
888 enum KFitError::ECode
889 VertexFitKFit::deleteTube() {
890  if (!m_FlagTube) return m_ErrorCode = KFitError::kNoError;
891 
892  if (m_iTrackTube == -1)
893  {
894  char buf[1024];
895  sprintf(buf, "%s:%s(): internal error; duplicated deleteTube() call?", __FILE__, __func__);
896  B2FATAL(buf);
897  }
898 
899  m_Tracks.pop_back();
900  m_TrackCount = m_Tracks.size();
901  m_iTrackTube = -1;
902 
903  return m_ErrorCode = KFitError::kNoError;
904 }
905 
906 enum KFitError::ECode VertexFitKFit::updateMother(Particle* mother)
907 {
908  MakeMotherKFit kmm;
909  kmm.setMagneticField(m_MagneticField);
910  unsigned n = getTrackCount();
911  for (unsigned i = 0; i < n; ++i) {
912  kmm.addTrack(getTrackMomentum(i), getTrackPosition(i), getTrackError(i),
913  getTrack(i).getCharge());
914  kmm.setTrackVertexError(getTrackVertexError(i));
915  for (unsigned j = i + 1; j < n; ++j) {
916  kmm.setCorrelation(getCorrelation(i, j));
917  }
918  }
919  kmm.setVertex(getVertex());
920  kmm.setVertexError(getVertexError());
921  m_ErrorCode = kmm.doMake();
922  if (m_ErrorCode != KFitError::kNoError)
923  return m_ErrorCode;
924  double chi2 = getCHIsq();
925  int ndf = getNDF();
926  double prob = TMath::Prob(chi2, ndf);
927  //
928  bool haschi2 = mother->hasExtraInfo("chiSquared");
929  if (haschi2) {
930  mother->setExtraInfo("chiSquared", chi2);
931  mother->setExtraInfo("ndf", ndf);
932  } else {
933  mother->addExtraInfo("chiSquared", chi2);
934  mother->addExtraInfo("ndf", ndf);
935  }
936 
937  mother->updateMomentum(
938  CLHEPToROOT::getTLorentzVector(kmm.getMotherMomentum()),
939  CLHEPToROOT::getTVector3(kmm.getMotherPosition()),
940  CLHEPToROOT::getTMatrixFSym(kmm.getMotherError()),
941  prob);
942  m_ErrorCode = KFitError::kNoError;
943  return m_ErrorCode;
944 }
Belle2::analysis::KFitBase::m_MagneticField
double m_MagneticField
Magnetic field.
Definition: KFitBase.h:311
Belle2::analysis::KFitBase::m_v
CLHEP::HepMatrix m_v
See J.Tanaka Ph.D (2001) p137 for definition.
Definition: KFitBase.h:285
Belle2::analysis::KFitError::kCannotGetARCSIN
@ kCannotGetARCSIN
Cannot get arcsin (bad track property or internal error)
Definition: KFitError.h:67
Belle2::analysis::VertexFitKFit::getVertexError
const CLHEP::HepSymMatrix getVertexError(void) const
Get a fitted vertex error matrix.
Definition: VertexFitKFit.cc:133
Belle2::analysis::KFitConst::kMaxIterationCount
static const int kMaxIterationCount
Maximum iteration step (internal use)
Definition: KFitConst.h:53
Belle2::analysis::VertexFitKFit::setIpTubeProfile
enum KFitError::ECode setIpTubeProfile(const CLHEP::HepLorentzVector &p, const HepPoint3D &x, const CLHEP::HepSymMatrix &e, const double q)
Set a virtual IP-tube track for the vertex constraint fit.
Definition: VertexFitKFit.cc:82
Belle2::analysis::KFitBase::m_CHIsq
double m_CHIsq
chi-square of the fit.
Definition: KFitBase.h:297
Belle2::analysis::KFitBase::m_NDF
int m_NDF
NDF of the fit.
Definition: KFitBase.h:295
Belle2::analysis::VertexFitKFit::~VertexFitKFit
~VertexFitKFit(void)
Destruct the object.
Belle2::analysis::VertexFitKFit::m_AfterVertexError
CLHEP::HepSymMatrix m_AfterVertexError
Vertex error matrix after the fit.
Definition: VertexFitKFit.h:172
Belle2::analysis::VertexFitKFit::getCHIsq
double getCHIsq(void) const override
Get a chi-square of the fit.
Definition: VertexFitKFit.cc:140
Belle2::analysis::VertexFitKFit::m_EachCHIsq
double m_EachCHIsq[KFitConst::kMaxTrackCount2]
Container of chi-square's of the input tracks.
Definition: VertexFitKFit.h:162
Belle2::analysis::KFitBase::m_ErrorCode
enum KFitError::ECode m_ErrorCode
Error code.
Definition: KFitBase.h:243
Belle2::analysis::KFitBase::m_FlagCorrelation
bool m_FlagCorrelation
Flag whether a correlation among tracks exists.
Definition: KFitBase.h:306
Belle2::analysis::KFitError::kBadInitialCHIsq
@ kBadInitialCHIsq
Bad initial chi-square (internal error)
Definition: KFitError.h:60
Belle2::analysis::KFitBase::m_D
CLHEP::HepMatrix m_D
See J.Tanaka Ph.D (2001) p137 for definition.
Definition: KFitBase.h:266
Belle2::analysis::VertexFitKFit::getVertex
const HepPoint3D getVertex(const int flag=KFitConst::kAfterFit) const
Get a vertex position.
Definition: VertexFitKFit.cc:114
Belle2::analysis::KFitConst::kBeforeFit
static const int kBeforeFit
Input parameter to specify before-fit when setting/getting a track attribute.
Definition: KFitConst.h:43
Belle2::analysis::MakeMotherKFit::setTrackVertexError
enum KFitError::ECode setTrackVertexError(const CLHEP::HepMatrix &e)
Set a vertex error matrix of the child particle in the addTrack'ed order.
Definition: MakeMotherKFit.cc:95
Belle2::analysis::KFitBase::getCorrelation
virtual const CLHEP::HepMatrix getCorrelation(const int id1, const int id2, const int flag=KFitConst::kAfterFit) const
Get a correlation matrix between two tracks.
Definition: KFitBase.cc:184
Belle2::analysis::VertexFitKFit::doFit3
enum KFitError::ECode doFit3(void)
Perform a standard vertex-constraint fit including IP-tube constraint.
Definition: VertexFitKFit.cc:236
Belle2::analysis::MakeMotherKFit::setVertex
enum KFitError::ECode setVertex(const HepPoint3D &v)
Set a vertex position of the mother particle.
Definition: MakeMotherKFit.cc:71
Belle2::analysis::KFitBase::isFitted
bool isFitted(void) const
Return false if fit is not performed yet or performed fit is failed; otherwise true.
Definition: KFitBase.cc:728
Belle2::analysis::VertexFitKFit::m_TubeTrack
KFitTrack m_TubeTrack
Entity of the virtual IP-tube track.
Definition: VertexFitKFit.h:189
Belle2::analysis::KFitBase::doFit2
enum KFitError::ECode doFit2(void)
Perform a fit (used in VertexFitKFit::doFit() and MassVertexFitKFit::doFit()).
Definition: KFitBase.cc:578
Belle2::analysis::VertexFitKFit::setKnownVertex
enum KFitError::ECode setKnownVertex(const bool flag=true)
Tell the object to perform a fit with vertex position fixed.
Definition: VertexFitKFit.cc:98
Belle2::analysis::KFitBase::m_NecessaryTrackCount
int m_NecessaryTrackCount
Number needed tracks to perform fit.
Definition: KFitBase.h:303
Belle2::analysis::KFitBase::getTrackError
const CLHEP::HepSymMatrix getTrackError(const int id) const
Get an error matrix of the track.
Definition: KFitBase.cc:169
Belle2::analysis::KFitBase::m_lam0
CLHEP::HepMatrix m_lam0
See J.Tanaka Ph.D (2001) p138 for definition.
Definition: KFitBase.h:283
Belle2::analysis::KFitBase::m_Cov_v_al_1
CLHEP::HepMatrix m_Cov_v_al_1
See J.Tanaka Ph.D (2001) p137 for definition.
Definition: KFitBase.h:291
Belle2::analysis::KFitBase::m_Tracks
std::vector< KFitTrack > m_Tracks
Container of input tracks.
Definition: KFitBase.h:249
Belle2::analysis::KFitBase::m_property
CLHEP::HepMatrix m_property
Container of charges and masses.
Definition: KFitBase.h:263
Belle2::analysis::KFitBase::prepareCorrelation
virtual enum KFitError::ECode prepareCorrelation(void)
Build a grand correlation matrix from input-track properties.
Definition: KFitBase.cc:459
Belle2::analysis::KFitBase::makeError1
const CLHEP::HepSymMatrix makeError1(const CLHEP::HepLorentzVector &p, const CLHEP::HepMatrix &e) const
Rebuild an error matrix from a Lorentz vector and an error matrix.
Belle2::analysis::KFitBase::m_lam
CLHEP::HepMatrix m_lam
See J.Tanaka Ph.D (2001) p137 for definition.
Definition: KFitBase.h:276
Belle2::analysis::KFitBase::m_V_D
CLHEP::HepMatrix m_V_D
See J.Tanaka Ph.D (2001) p138 for definition.
Definition: KFitBase.h:271
Belle2::analysis::KFitBase::m_al_1
CLHEP::HepMatrix m_al_1
See J.Tanaka Ph.D (2001) p136 for definition.
Definition: KFitBase.h:259
Belle2::analysis::VertexFitKFit::deleteTube
enum KFitError::ECode deleteTube(void)
Delete the virtual tube track to m_Tracks just after the internal minimization call.
Definition: VertexFitKFit.cc:889
Belle2::analysis::VertexFitKFit::doFit5
enum KFitError::ECode doFit5(void)
Perform a fixed-vertex-position fit mainly for slow pion.
Definition: VertexFitKFit.cc:527
Belle2::analysis::KFitError::displayError
static void displayError(const char *file, const int line, const char *func, const enum ECode code)
Display a description of error and its location.
Definition: KFitError.h:79
Belle2::analysis::KFitBase::m_E
CLHEP::HepMatrix m_E
See J.Tanaka Ph.D (2001) p137 for definition.
Definition: KFitBase.h:279
Belle2::analysis::MakeMotherKFit::getMotherError
const CLHEP::HepSymMatrix getMotherError(void) const
Get an error matrix of the mother particle.
Definition: MakeMotherKFit.cc:185
Belle2::analysis::KFitBase::m_FlagFitted
bool m_FlagFitted
Flag to indicate if the fit is performed and succeeded.
Definition: KFitBase.h:245
Belle2::analysis::VertexFitKFit::setIpProfile
enum KFitError::ECode setIpProfile(const HepPoint3D &ip, const CLHEP::HepSymMatrix &ipe)
Set an IP-ellipsoid shape for the vertex constraint fit.
Definition: VertexFitKFit.cc:65
Belle2::analysis::MakeMotherKFit
MakeMotherKFit is a class to build mother particle from kinematically fitted daughters.
Definition: MakeMotherKFit.h:41
Belle2::analysis::VertexFitKFit::m_BeforeVertex
HepPoint3D m_BeforeVertex
Vertex position before the fit.
Definition: VertexFitKFit.h:167
Belle2::analysis::KFitConst::kInitialCHIsq
static constexpr double kInitialCHIsq
Initial chi-square value (internal use)
Definition: KFitConst.h:56
Belle2::analysis::MakeMotherKFit::getMotherMomentum
const CLHEP::HepLorentzVector getMotherMomentum(void) const
Get a Lorentz vector of the mother particle.
Definition: MakeMotherKFit.cc:171
Belle2::analysis::KFitError::kNoError
@ kNoError
No error.
Definition: KFitError.h:45
Belle2::analysis::KFitBase::getTrackCount
int getTrackCount(void) const
Get the number of added tracks.
Definition: KFitBase.cc:108
Belle2::analysis::KFitConst::kAfterFit
static const int kAfterFit
Input parameter to specify after-fit when setting/getting a track attribute.
Definition: KFitConst.h:45
Belle2::analysis::KFitBase::m_v_a
CLHEP::HepMatrix m_v_a
See J.Tanaka Ph.D (2001) p137 for definition.
Definition: KFitBase.h:287
Belle2::analysis::KFitConst::kMaxTrackCount
static const int kMaxTrackCount
Maximum track size.
Definition: KFitConst.h:48
Belle2::analysis::KFitConst::kMaxTrackCount2
static const int kMaxTrackCount2
Maximum track size (internal use)
Definition: KFitConst.h:50
Belle2::analysis::VertexFitKFit::m_iTrackTube
int m_iTrackTube
ID of the virtual tube track in the m_Tracks.
Definition: VertexFitKFit.h:187
Belle2::analysis::VertexFitKFit::doFit4
enum KFitError::ECode doFit4(void)
Perform a IP-ellipsoid and vertex-constraint fit.
Definition: VertexFitKFit.cc:403
Belle2::analysis::KFitBase::m_V_al_1
CLHEP::HepMatrix m_V_al_1
See J.Tanaka Ph.D (2001) p138 for definition.
Definition: KFitBase.h:274
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::analysis::VertexFitKFit::setCorrelationMode
enum KFitError::ECode setCorrelationMode(const bool m)
Tell the object to perform a fit with track correlations.
Definition: VertexFitKFit.cc:106
Belle2::analysis::VertexFitKFit::updateMother
enum KFitError::ECode updateMother(Particle *mother)
Update mother particle.
Definition: VertexFitKFit.cc:906
Belle2::analysis::KFitBase::getNDF
virtual int getNDF(void) const
Get an NDF of the fit.
Definition: KFitBase.cc:115
Belle2::analysis::KFitError::kDivisionByZero
@ kDivisionByZero
Division by zero (bad track property or internal error)
Definition: KFitError.h:63
Belle2::analysis::VertexFitKFit::m_AfterVertex
HepPoint3D m_AfterVertex
Vertex position after the fit.
Definition: VertexFitKFit.h:170
Belle2::analysis::VertexFitKFit::m_AfterTrackVertexError
std::vector< CLHEP::HepMatrix > m_AfterTrackVertexError
Array of vertex error matrices after the fit.
Definition: VertexFitKFit.h:174
Belle2::analysis::VertexFitKFit::getCHIsqVertex
double getCHIsqVertex(void) const
Get a chi-square of the fit excluding IP-constraint part.
Definition: VertexFitKFit.cc:147
Belle2::analysis::VertexFitKFit::getTrackCHIsq
double getTrackCHIsq(const int id) const override
Get a chi-square of the track.
Definition: VertexFitKFit.cc:163
Belle2::analysis::MakeMotherKFit::setMagneticField
enum KFitError::ECode setMagneticField(const double mf)
Change a magnetic field from the default value KFitConst::kDefaultMagneticField.
Definition: MakeMotherKFit.cc:63
Belle2::analysis::VertexFitKFit::m_CorrelationMode
bool m_CorrelationMode
Flag controlled by setCorrelationMode().
Definition: VertexFitKFit.h:160
Belle2::analysis::VertexFitKFit::getTrackPartNDF
int getTrackPartNDF(void) const
Get an NDF relevant to the getTrackPartCHIsq().
Definition: VertexFitKFit.cc:199
Belle2::analysis::KFitTrack
KFitTrack is a container of the track information (Lorentz vector, position, and error matrix),...
Definition: KFitTrack.h:39
Belle2::analysis::KFitBase::m_TrackCount
int m_TrackCount
Number of tracks.
Definition: KFitBase.h:301
Belle2::analysis::KFitError::kCannotGetMatrixInverse
@ kCannotGetMatrixInverse
Cannot calculate matrix inverse (bad track property or internal error)
Definition: KFitError.h:65
Belle2::analysis::KFitBase::m_al_0
CLHEP::HepMatrix m_al_0
See J.Tanaka Ph.D (2001) p136 for definition.
Definition: KFitBase.h:257
Belle2::analysis::VertexFitKFit::m_FlagBeam
bool m_FlagBeam
Flag if to perform IP-ellipsoid constraint fit.
Definition: VertexFitKFit.h:177
Belle2::analysis::KFitConst::kNumber6
static const int kNumber6
Constant 6 to check matrix size (internal use)
Definition: KFitConst.h:38
Belle2::analysis::VertexFitKFit::prepareInputSubMatrix
enum KFitError::ECode prepareInputSubMatrix(void) override
Build sub-matrices for minimum search from input-track properties.
Definition: VertexFitKFit.cc:698
Belle2::analysis::KFitConst::kLightSpeed
static constexpr double kLightSpeed
Speed of light.
Definition: KFitConst.h:62
Belle2::analysis::VertexFitKFit::getTrackVertexError
const CLHEP::HepMatrix getTrackVertexError(const int id) const
Get a vertex error matrix of the track.
Definition: VertexFitKFit.cc:155
Belle2::analysis::VertexFitKFit::doFit
enum KFitError::ECode doFit(void)
Perform a vertex-constraint fit.
Definition: VertexFitKFit.cc:216
Belle2::analysis::MakeMotherKFit::setCorrelation
enum KFitError::ECode setCorrelation(const CLHEP::HepMatrix &e)
Set a correlation matrix.
Definition: MakeMotherKFit.cc:119
Belle2::analysis::VertexFitKFit::prepareInputMatrix
enum KFitError::ECode prepareInputMatrix(void) override
Build grand matrices for minimum search from input-track properties.
Definition: VertexFitKFit.cc:619
Belle2::analysis::MakeMotherKFit::setVertexError
enum KFitError::ECode setVertexError(const CLHEP::HepSymMatrix &e)
Set a vertex error matrix of the mother particle.
Definition: MakeMotherKFit.cc:79
Belle2::analysis::KFitBase::getTrack
const KFitTrack getTrack(const int id) const
Get a specified track object.
Definition: KFitBase.cc:176
Belle2::analysis::VertexFitKFit::getTrackPartCHIsq
double getTrackPartCHIsq(void) const
Get a sum of the chi-square associated to the input tracks.
Definition: VertexFitKFit.cc:176
Belle2::analysis::KFitBase::m_d
CLHEP::HepMatrix m_d
See J.Tanaka Ph.D (2001) p137 for definition.
Definition: KFitBase.h:268
Belle2::analysis::MakeMotherKFit::getMotherPosition
const HepPoint3D getMotherPosition(void) const
Get a position of the mother particle.
Definition: MakeMotherKFit.cc:178
Belle2::analysis::KFitBase::getCHIsq
virtual double getCHIsq(void) const
Get a chi-square of the fit.
Definition: KFitBase.cc:122
Belle2::Particle
Class to store reconstructed particles.
Definition: Particle.h:77
Belle2::analysis::KFitBase::m_V_Dt
CLHEP::HepMatrix m_V_Dt
See J.Tanaka Ph.D (2001) p138 for definition.
Definition: KFitBase.h:289
Belle2::analysis::VertexFitKFit::appendTube
enum KFitError::ECode appendTube(void)
Add the virtual tube track to m_Tracks just before the internal minimization call.
Definition: VertexFitKFit.cc:870
Belle2::analysis::KFitBase::makeError2
const CLHEP::HepMatrix makeError2(const CLHEP::HepLorentzVector &p, const CLHEP::HepMatrix &e) const
Rebuild an error matrix from a Lorentz vector and an error matrix.
Definition: KFitBase.cc:297
Belle2::analysis::MakeMotherKFit::addTrack
enum KFitError::ECode addTrack(const KFitTrack &kp)
Add a track to the make-mother object.
Definition: MakeMotherKFit.cc:40
Belle2::analysis::KFitBase::isTrackIDInRange
bool isTrackIDInRange(const int id) const
Check if the id is in the range.
Definition: KFitBase.cc:739
Belle2::analysis::KFitBase::getTrackPosition
const HepPoint3D getTrackPosition(const int id) const
Get a position of the track.
Definition: KFitBase.cc:162
Belle2::analysis::VertexFitKFit::m_CHIsqVertex
double m_CHIsqVertex
chi-square of the fit excluding IP-constraint part.
Definition: VertexFitKFit.h:164
Belle2::analysis::VertexFitKFit::prepareOutputMatrix
enum KFitError::ECode prepareOutputMatrix(void) override
Build an output error matrix.
Definition: VertexFitKFit.cc:709
Belle2::analysis::KFitError::ECode
ECode
ECode is a error code enumerate.
Definition: KFitError.h:43
Belle2::analysis::MakeMotherKFit::doMake
enum KFitError::ECode doMake(void)
Perform a reconstruction of mother particle.
Definition: MakeMotherKFit.cc:192
Belle2::analysis::VertexFitKFit::calculateNDF
enum KFitError::ECode calculateNDF(void) override
Calculate an NDF of the fit.
Definition: VertexFitKFit.cc:859
Belle2::analysis::KFitError::kOutOfRange
@ kOutOfRange
Specified track-id out of range.
Definition: KFitError.h:51
HepGeom::Point3D< double >
Belle2::analysis::KFitConst::kNumber7
static const int kNumber7
Constant 7 to check matrix size (internal use)
Definition: KFitConst.h:40
Belle2::analysis::VertexFitKFit::m_FlagTube
bool m_FlagTube
Flag if to perform IP-tube constraint fit.
Definition: VertexFitKFit.h:185
Belle2::analysis::KFitBase::m_V_E
CLHEP::HepMatrix m_V_E
See J.Tanaka Ph.D (2001) p138 for definition.
Definition: KFitBase.h:281
Belle2::analysis::VertexFitKFit::m_BeamError
CLHEP::HepSymMatrix m_BeamError
Error matrix modeling the IP ellipsoid.
Definition: VertexFitKFit.h:179
Belle2::analysis::KFitError::kBadTrackSize
@ kBadTrackSize
Track count too small to perform fit.
Definition: KFitError.h:54
Belle2::analysis::KFitBase::m_V_al_0
CLHEP::HepSymMatrix m_V_al_0
See J.Tanaka Ph.D (2001) p137 for definition.
Definition: KFitBase.h:255
Belle2::analysis::KFitBase::m_FlagOverIteration
bool m_FlagOverIteration
Flag whether the iteration count exceeds the limit.
Definition: KFitBase.h:308
Belle2::analysis::VertexFitKFit::setInitialVertex
enum KFitError::ECode setInitialVertex(const HepPoint3D &v)
Set an initial vertex point for the vertex-vertex constraint fit.
Definition: VertexFitKFit.cc:51
Belle2::analysis::KFitBase::m_al_a
CLHEP::HepMatrix m_al_a
See J.Tanaka Ph.D (2001) p137 for definition.
Definition: KFitBase.h:261
Belle2::analysis::VertexFitKFit::makeCoreMatrix
enum KFitError::ECode makeCoreMatrix(void) override
Build matrices using the kinematical constraint.
Definition: VertexFitKFit.cc:764
Belle2::analysis::KFitBase::getTrackCHIsq
virtual double getTrackCHIsq(const int id) const
Get a chi-square of the track.
Definition: KFitBase.cc:136
Belle2::analysis::VertexFitKFit::m_FlagKnownVertex
bool m_FlagKnownVertex
Flag controlled by setKnownVertex().
Definition: VertexFitKFit.h:182
Belle2::analysis::KFitBase::getTrackMomentum
const CLHEP::HepLorentzVector getTrackMomentum(const int id) const
Get a Lorentz vector of the track.
Definition: KFitBase.cc:155