diff --git a/algo/kf/core/KfTrackKalmanFilter.cxx b/algo/kf/core/KfTrackKalmanFilter.cxx
index 1f0966a3303182ab19db4d4107717d08e35c4440..ad4da0068a14424984fd631fc470653417fb9567 100644
--- a/algo/kf/core/KfTrackKalmanFilter.cxx
+++ b/algo/kf/core/KfTrackKalmanFilter.cxx
@@ -997,15 +997,15 @@ namespace cbm::algo::kf
DataT s0 = lg * qp * t;
DataT a = (DataT(1.) + fMass2 * qp * qp) * s0 * s0 * t * radThick;
-
+
// Approximate formula
// DataT h = txtx + tyty;
// DataT h2 = h * h;
- //cnst c1 = 0.0136f, c2 = c1 * 0.038f, c3 = c2 * 0.5f, c4 = -c3 / 2.0f, c5 = c3 / 3.0f, c6 = -c3 / 4.0f;
- //DataT s0 = (c1 + c2 * log(radThick) + c3 * h + h2 * (c4 + c5 * h + c6 * h2)) * qp0t;
- //DataT a = ( (kONE+mass2*qp0*qp0t)*radThick*s0*s0 );
- //DataT a = ((t + fMass2 * qp0 * qp0t) * radThick * s0 * s0);
+ // cnst c1 = 0.0136f, c2 = c1 * 0.038f, c3 = c2 * 0.5f, c4 = -c3 / 2.0f, c5 = c3 / 3.0f, c6 = -c3 / 4.0f;
+ // DataT s0 = (c1 + c2 * log(radThick) + c3 * h + h2 * (c4 + c5 * h + c6 * h2)) * qp0t;
+ // DataT a = ( (kONE+mass2*qp0*qp0t)*radThick*s0*s0 );
+ // DataT a = ((t + fMass2 * qp0 * qp0t) * radThick * s0 * s0);
fTr.C22() += kf::utils::iif(fMask, txtx1 * a, DataT(0.));
fTr.C32() += kf::utils::iif(fMask, tx * ty * a, DataT(0.));
diff --git a/reco/KF/KF.cmake b/reco/KF/KF.cmake
index 1143c637c3e6d2ddbf0fcc62a3dfb2bc9796619f..7eccc20caa32e128c866f299d1897e55527d27f1 100644
--- a/reco/KF/KF.cmake
+++ b/reco/KF/KF.cmake
@@ -18,8 +18,6 @@ set(SRCS
CbmKFPrimaryVertexFinder.cxx
CbmKFTrackInterface.cxx
CbmKFVertexInterface.cxx
- CbmKfTrackFitter.cxx
- CbmKfFitTracksTask.cxx
Interface/CbmKFStsHit.cxx
Interface/CbmKFTrack.cxx
diff --git a/reco/KF/KFLinkDef.h b/reco/KF/KFLinkDef.h
index 3e05ec2fe7c2cd031b2ea0d0768ac98727ff0083..b8dd970b728a08ea0f09375ea8ae2090102e6b8d 100644
--- a/reco/KF/KFLinkDef.h
+++ b/reco/KF/KFLinkDef.h
@@ -27,8 +27,8 @@
#pragma link C++ class CbmStsKFTrackFitter + ;
#pragma link C++ class CbmPVFinderKF + ;
#pragma link C++ class CbmPVFinderKFGlobal + ;
-#pragma link C++ class CbmKfTrackFitter + ;
-#pragma link C++ class CbmKfFitTracksTask + ;
+//#pragma link C++ class CbmKfTrackFitter + ;
+//#pragma link C++ class CbmKfFitTracksTask + ;
//#pragma link C++ class CbmTrackingDetectorInterfaceInit + ;
//KFQA
diff --git a/reco/L1/CMakeLists.txt b/reco/L1/CMakeLists.txt
index 4a0da6c73a52643eb24dff5ac87629790ee84c94..8c61bf88f11d7fdec12c422ab765c00b9b659dff 100644
--- a/reco/L1/CMakeLists.txt
+++ b/reco/L1/CMakeLists.txt
@@ -34,7 +34,6 @@ set(SRCS
OffLineInterface/CbmL1GlobalFindTracksEvents.cxx
OffLineInterface/CbmGenerateMaterialMaps.cxx
- CbmL1Util.cxx
CbmL1Performance.cxx
CbmL1MCTrack.cxx
CbmL1Track.cxx
@@ -71,7 +70,6 @@ set(NO_DICT_SRCS
set(HEADERS
CbmL1Constants.h
- CbmL1Util.h
CbmL1DetectorID.h
CbmL1MCPoint.h
CbmL1Hit.h
diff --git a/reco/L1/CbmCaTimeSliceReader.cxx b/reco/L1/CbmCaTimeSliceReader.cxx
index bd367294c28a613bd13a0729bf6461accc2f2c3a..156e0387f890175a58c341cf56fd2985edc16d0c 100644
--- a/reco/L1/CbmCaTimeSliceReader.cxx
+++ b/reco/L1/CbmCaTimeSliceReader.cxx
@@ -13,8 +13,8 @@
#include "CaInputData.h"
#include "CaParameters.h"
#include "CbmGlobalTrack.h"
+#include "CbmKfUtil.h" // for CopyTrackParam2TC
#include "CbmL1.h"
-#include "CbmL1Util.h" // for CopyTrackParam2TC
#include "CbmMuchTrack.h"
#include "CbmStsTrack.h"
#include "CbmTofTrack.h"
@@ -217,8 +217,8 @@ void TimeSliceReader::ReadRecoTracks()
auto* pInputTrack = static_cast<CbmStsTrack*>(fpBrRecoTracks->At(iText));
auto& track = (*fpvTracks)[iT];
- track.Set(cbm::L1Util::ConvertTrackParam(*pInputTrack->GetParamFirst()));
- track.TLast = cbm::L1Util::ConvertTrackParam(*pInputTrack->GetParamLast());
+ track.Set(cbm::kf::ConvertTrackParam(*pInputTrack->GetParamFirst()));
+ track.TLast = cbm::kf::ConvertTrackParam(*pInputTrack->GetParamLast());
track.ChiSq() = pInputTrack->GetChiSq();
track.Ndf() = pInputTrack->GetNDF();
track.Tpv.Time() = pInputTrack->GetStartTime();
@@ -250,8 +250,8 @@ void TimeSliceReader::ReadRecoTracks()
int iText = fpEvent ? fpEvent->GetIndex(ECbmDataType::kGlobalTrack, iT) : iT;
auto* pInputTrack = static_cast<CbmGlobalTrack*>(fpBrRecoTracks->At(iText));
auto& track = (*fpvTracks)[iT];
- track.Set(cbm::L1Util::ConvertTrackParam(*pInputTrack->GetParamFirst()));
- track.TLast = cbm::L1Util::ConvertTrackParam(*pInputTrack->GetParamLast());
+ track.Set(cbm::kf::ConvertTrackParam(*pInputTrack->GetParamFirst()));
+ track.TLast = cbm::kf::ConvertTrackParam(*pInputTrack->GetParamLast());
track.ChiSq() = pInputTrack->GetChi2();
track.Ndf() = pInputTrack->GetNDF();
diff --git a/reco/L1/OffLineInterface/CbmL1GlobalTrackFinder.cxx b/reco/L1/OffLineInterface/CbmL1GlobalTrackFinder.cxx
index fa8b55181b54ccc42badb203b8347bcb22d92655..c1f15963b2f5b3ce19368dfb8e5b742776878cf1 100644
--- a/reco/L1/OffLineInterface/CbmL1GlobalTrackFinder.cxx
+++ b/reco/L1/OffLineInterface/CbmL1GlobalTrackFinder.cxx
@@ -21,7 +21,7 @@
#include "CbmEvent.h"
#include "CbmGlobalTrack.h"
-#include "CbmL1Util.h"
+#include "CbmKfUtil.h"
#include "CbmMuchTrack.h"
#include "CbmStsHit.h"
#include "CbmStsTrack.h"
@@ -84,8 +84,8 @@ Int_t CbmL1GlobalTrackFinder::CopyL1Tracks(CbmEvent* event)
if (event) event->AddData(ECbmDataType::kGlobalTrack, globalTrackIndex);
CbmGlobalTrack* t = dynamic_cast<CbmGlobalTrack*>(fGlobalTracks->At(globalTrackIndex++));
t->SetFlag(0);
- t->SetParamFirst(cbm::L1Util::ConvertTrackParam(T));
- t->SetParamLast(cbm::L1Util::ConvertTrackParam(T.TLast));
+ t->SetParamFirst(cbm::kf::ConvertTrackParam(T));
+ t->SetParamLast(cbm::kf::ConvertTrackParam(T.TLast));
t->SetChi2(T.GetChiSq());
// t->SetLength(T.length);
t->SetNDF(T.GetNdf());
@@ -156,8 +156,8 @@ void CbmL1GlobalTrackFinder::CbmL1TrackToCbmTrack(CbmL1Track l1track, CbmTrack*
//track->SetPreviousTrackId(l1track.GetPreviousTrackId());//???
//track->SetFlag(l1track.GetQuality());//???
- track->SetParamFirst(cbm::L1Util::ConvertTrackParam(l1track));
- track->SetParamLast(cbm::L1Util::ConvertTrackParam(l1track.TLast));
+ track->SetParamFirst(cbm::kf::ConvertTrackParam(l1track));
+ track->SetParamLast(cbm::kf::ConvertTrackParam(l1track.TLast));
}
// -------------------------------------------------------------------------
@@ -188,8 +188,8 @@ void CbmL1GlobalTrackFinder::CbmL1TrackToCbmStsTrack(CbmL1Track l1track, CbmStsT
track->SetLastHitTime(T.TLast.GetTime());
track->SetLastHitTimeError(T.TLast.GetTimeError());
- track->SetParamFirst(cbm::L1Util::ConvertTrackParam(T));
- track->SetParamLast(cbm::L1Util::ConvertTrackParam(T.TLast));
+ track->SetParamFirst(cbm::kf::ConvertTrackParam(T));
+ track->SetParamLast(cbm::kf::ConvertTrackParam(T.TLast));
}
// -------------------------------------------------------------------------
@@ -211,8 +211,8 @@ void CbmL1GlobalTrackFinder::CbmL1TrackToCbmMuchTrack(CbmL1Track l1track, CbmMuc
//track->SetPreviousTrackId(T->GetPreviousTrackId());//???
//track->SetFlag(T->GetQuality());//???
- track->SetParamFirst(cbm::L1Util::ConvertTrackParam(T));
- track->SetParamLast(cbm::L1Util::ConvertTrackParam(T.TLast));
+ track->SetParamFirst(cbm::kf::ConvertTrackParam(T));
+ track->SetParamLast(cbm::kf::ConvertTrackParam(T.TLast));
}
// -------------------------------------------------------------------------
@@ -234,8 +234,8 @@ void CbmL1GlobalTrackFinder::CbmL1TrackToCbmTrdTrack(CbmL1Track l1track, CbmTrdT
//track->SetPreviousTrackId(T->GetPreviousTrackId());//???
//track->SetFlag(T->GetQuality());//???
- track->SetParamFirst(cbm::L1Util::ConvertTrackParam(T.Tpv));
- track->SetParamLast(cbm::L1Util::ConvertTrackParam(T.TLast));
+ track->SetParamFirst(cbm::kf::ConvertTrackParam(T.Tpv));
+ track->SetParamLast(cbm::kf::ConvertTrackParam(T.TLast));
}
// -------------------------------------------------------------------------
@@ -256,8 +256,8 @@ void CbmL1GlobalTrackFinder::CbmL1TrackToCbmTofTrack(CbmL1Track l1track, CbmTofT
//track->SetPreviousTrackId(T->GetPreviousTrackId());//???
//track->SetFlag(T->GetQuality());//???
- track->SetParamFirst(cbm::L1Util::ConvertTrackParam(T));
- track->SetParamLast(cbm::L1Util::ConvertTrackParam(T.TLast));
+ track->SetParamFirst(cbm::kf::ConvertTrackParam(T));
+ track->SetParamLast(cbm::kf::ConvertTrackParam(T.TLast));
}
// -------------------------------------------------------------------------
diff --git a/reco/L1/OffLineInterface/CbmL1StsTrackFinder.cxx b/reco/L1/OffLineInterface/CbmL1StsTrackFinder.cxx
index f0b3f5ea29eda9e9dda664581cf7edb381ade7af..3986a33d3cb18be575f6e5c51a3e8929f94065a6 100644
--- a/reco/L1/OffLineInterface/CbmL1StsTrackFinder.cxx
+++ b/reco/L1/OffLineInterface/CbmL1StsTrackFinder.cxx
@@ -20,7 +20,7 @@
#include "CbmL1StsTrackFinder.h"
#include "CbmEvent.h"
-#include "CbmL1Util.h"
+#include "CbmKfUtil.h"
#include "CbmStsHit.h"
#include "CbmStsTrack.h"
#include "FairHit.h"
@@ -70,8 +70,8 @@ Int_t CbmL1StsTrackFinder::CopyL1Tracks(CbmEvent* event)
if (event) event->AddData(ECbmDataType::kStsTrack, trackIndex);
CbmStsTrack* t = dynamic_cast<CbmStsTrack*>(fTracks->At(trackIndex++));
t->SetFlag(0);
- t->SetParamFirst(cbm::L1Util::ConvertTrackParam(T));
- t->SetParamLast(cbm::L1Util::ConvertTrackParam(T.TLast));
+ t->SetParamFirst(cbm::kf::ConvertTrackParam(T));
+ t->SetParamLast(cbm::kf::ConvertTrackParam(T.TLast));
t->SetChiSq(T.GetChiSq());
t->SetNDF(T.GetNdf());
t->SetChiSqTime(T.GetChiSqTime());
diff --git a/reco/alignment/CbmBbaAlignmentTask.cxx b/reco/alignment/CbmBbaAlignmentTask.cxx
index 3e854451a176d63b5e4f531154bb347126935d74..17268f7b28d6aa3ddf12ea53127b610ae79a31b1 100644
--- a/reco/alignment/CbmBbaAlignmentTask.cxx
+++ b/reco/alignment/CbmBbaAlignmentTask.cxx
@@ -318,12 +318,12 @@ void CbmBbaAlignmentTask::Exec(Option_t* /*opt*/)
n.fIsTimeSet = false;
}
- if (!fFitter.FitTrack(t.fUnalignedTrack)) {
+ if (!fFitter.FitTrajectory(t.fUnalignedTrack)) {
LOG(warning) << "failed to fit the track! ";
continue;
}
- const auto& par = t.fUnalignedTrack.fNodes.front().fTrack;
+ const auto& par = t.fUnalignedTrack.fNodes.front().fParamUp;
if (t.fNstsHits < l1Par.GetNstationsActive(ca::EDetectorID::kSts)) continue;
if (!std::isfinite(par.GetQp())) continue;
@@ -364,13 +364,12 @@ void CbmBbaAlignmentTask::Exec(Option_t* /*opt*/)
n.fMxy.SetNdfX(0);
}
}
- n.fTrack.SetQp(0.);
}
if (t.fNstsHits < 2) continue;
if (t.fNtofHits < 2) continue;
if (t.fNtrd1dHits + t.fNtrd2dHits < 2) continue;
- if (!fFitter.FitTrack(t.fUnalignedTrack)) {
+ if (!fFitter.FitTrajectory(t.fUnalignedTrack)) {
LOG(warning) << "failed to fit the track! ";
continue;
}
@@ -385,7 +384,7 @@ void CbmBbaAlignmentTask::Exec(Option_t* /*opt*/)
for (TrackContainer& t : fTracks) {
if (!t.fIsActive) continue;
for (unsigned int in = 0; in < t.fUnalignedTrack.fNodes.size(); in++) {
- CbmKfTrackFitter::FitNode& node = t.fUnalignedTrack.fNodes[in];
+ CbmKfTrackFitter::TrajectoryNode& node = t.fUnalignedTrack.fNodes[in];
if (!node.fIsFitted) {
t.fIsActive = false;
break;
@@ -409,20 +408,20 @@ void CbmBbaAlignmentTask::Exec(Option_t* /*opt*/)
const double cutX = 8.;
const double cutY = 8.;
for (unsigned int in = 0; in < t.fUnalignedTrack.fNodes.size(); in++) {
- CbmKfTrackFitter::Track tr = t.fUnalignedTrack;
- CbmKfTrackFitter::FitNode& node = tr.fNodes[in];
+ CbmKfTrackFitter::Trajectory tr = t.fUnalignedTrack;
+ CbmKfTrackFitter::TrajectoryNode& node = tr.fNodes[in];
if (!node.fIsXySet) {
continue;
}
node.fIsXySet = false; // exclude the node from the fit
- if (!fFitter.FitTrack(tr)) {
+ if (!fFitter.FitTrajectory(tr)) {
t.fIsActive = false;
break;
}
- double x_residual = node.fMxy.X() - node.fTrack.X(); // this should be the difference vector
- double y_residual = node.fMxy.Y() - node.fTrack.Y();
- double x_pull = x_residual / sqrt(node.fMxy.Dx2() + node.fTrack.GetCovariance(0, 0));
- double y_pull = y_residual / sqrt(node.fMxy.Dy2() + node.fTrack.GetCovariance(1, 1));
+ double x_residual = node.fMxy.X() - node.fParamDn.X(); // this should be the difference vector
+ double y_residual = node.fMxy.Y() - node.fParamDn.Y();
+ double x_pull = x_residual / sqrt(node.fMxy.Dx2() + node.fParamDn.GetCovariance(0, 0));
+ double y_pull = y_residual / sqrt(node.fMxy.Dy2() + node.fParamDn.GetCovariance(1, 1));
if (!node.fIsFitted || (node.fMxy.NdfX() > 0. && fabs(x_pull) > cutX)
|| (node.fMxy.NdfY() > 0. && fabs(y_pull) > cutY)) {
t.fIsActive = false;
@@ -488,12 +487,12 @@ void CbmBbaAlignmentTask::ApplyAlignment(const std::vector<double>& par)
for (TrackContainer& t : fTracks) {
for (unsigned int in = 0; in < t.fUnalignedTrack.fNodes.size(); in++) {
- CbmKfTrackFitter::FitNode& node = t.fUnalignedTrack.fNodes[in];
+ CbmKfTrackFitter::TrajectoryNode& node = t.fUnalignedTrack.fNodes[in];
if (!node.fIsXySet) {
continue;
}
- CbmKfTrackFitter::FitNode& nodeAligned = t.fAlignedTrack.fNodes[in];
- int iSensor = node.fReference1;
+ CbmKfTrackFitter::TrajectoryNode& nodeAligned = t.fAlignedTrack.fNodes[in];
+ int iSensor = node.fReference1;
assert(iSensor >= 0 && iSensor < (int) fSensors.size());
assert(nodeAligned.fReference1 == node.fReference1);
@@ -510,9 +509,18 @@ void CbmBbaAlignmentTask::ApplyAlignment(const std::vector<double>& par)
nodeAligned.fZ = node.fZ + dz;
// shift the fitted track to the Z position of the hit
- nodeAligned.fTrack.SetX(node.fTrack.X() + node.fTrack.Tx() * dz);
- nodeAligned.fTrack.SetY(node.fTrack.Y() + node.fTrack.Ty() * dz);
- nodeAligned.fTrack.SetZ(node.fTrack.Z() + dz);
+ {
+ auto& p = nodeAligned.fParamDn;
+ p.SetX(p.X() + p.Tx() * dz);
+ p.SetY(p.Y() + p.Ty() * dz);
+ p.SetZ(nodeAligned.fZ);
+ }
+ {
+ auto& p = nodeAligned.fParamUp;
+ p.SetX(p.X() + p.Tx() * dz);
+ p.SetY(p.Y() + p.Ty() * dz);
+ p.SetZ(nodeAligned.fZ);
+ }
}
}
@@ -555,14 +563,14 @@ double CbmBbaAlignmentTask::CostFunction(const std::vector<double>& par)
// fit.SetDebugInfo(" track " + std::to_string(itr));
for (unsigned int in = 0; in < t.fAlignedTrack.fNodes.size(); in++) {
- CbmKfTrackFitter::FitNode& node = t.fAlignedTrack.fNodes[in];
+ CbmKfTrackFitter::TrajectoryNode& node = t.fAlignedTrack.fNodes[in];
if (!node.fIsFitted) {
LOG(fatal) << "BBA: Cost function: aligned node is not fitted! ";
}
}
- bool ok = fit.FitTrack(t.fAlignedTrack);
- double chi2 = t.fAlignedTrack.fNodes.back().fTrack.GetChiSq();
- double ndf = t.fAlignedTrack.fNodes.back().fTrack.GetNdf();
+ bool ok = fit.FitTrajectory(t.fAlignedTrack);
+ double chi2 = t.fAlignedTrack.fNodes.back().fParamDn.GetChiSq();
+ double ndf = t.fAlignedTrack.fNodes.back().fParamDn.GetNdf();
if (!ok || !std::isfinite(chi2) || chi2 <= 0. || ndf <= 0.) {
// TODO: deactivate the track and continue
LOG(fatal) << "BBA: fit failed! ";
@@ -836,8 +844,8 @@ void CbmBbaAlignmentTask::Finish()
fCostInitial = CostFunction(parMisaligned);
for (auto& t : fTracks) {
- double ndf = t.fAlignedTrack.fNodes.back().fTrack.GetNdf();
- double chi2 = t.fAlignedTrack.fNodes.back().fTrack.GetChiSq();
+ double ndf = t.fAlignedTrack.fNodes.back().fParamDn.GetNdf();
+ double chi2 = t.fAlignedTrack.fNodes.back().fParamDn.GetChiSq();
if (ndf < 0. || chi2 < 0. || !std::isfinite(ndf) || !std::isfinite(chi2)) {
t.fIsActive = false;
}
@@ -850,17 +858,17 @@ void CbmBbaAlignmentTask::Finish()
for (const auto& t : fTracks) {
if (!t.fIsActive) continue;
- // calculate the residuals for all tracks at each fitNode before alignment
+ // calculate the residuals for all tracks at each TrajectoryNode before alignment
for (unsigned int in = 0; in < t.fAlignedTrack.fNodes.size(); in++) {
- CbmKfTrackFitter::Track tr = t.fAlignedTrack;
- CbmKfTrackFitter::FitNode& node = tr.fNodes[in];
+ CbmKfTrackFitter::Trajectory tr = t.fAlignedTrack;
+ CbmKfTrackFitter::TrajectoryNode& node = tr.fNodes[in];
if (!node.fIsXySet) {
continue;
}
node.fIsXySet = false;
- fFitter.FitTrack(tr);
+ fFitter.FitTrajectory(tr);
if (!node.fIsFitted) {
continue;
@@ -868,10 +876,10 @@ void CbmBbaAlignmentTask::Finish()
Sensor& s = fSensors[node.fReference1];
- double x_residual = node.fMxy.X() - node.fTrack.X(); // this should be the difference vector
- double y_residual = node.fMxy.Y() - node.fTrack.Y();
- double x_pull = x_residual / sqrt(node.fMxy.Dx2() + node.fTrack.GetCovariance(0, 0));
- double y_pull = y_residual / sqrt(node.fMxy.Dy2() + node.fTrack.GetCovariance(1, 1));
+ double x_residual = node.fMxy.X() - node.fParamDn.X(); // this should be the difference vector
+ double y_residual = node.fMxy.Y() - node.fParamDn.Y();
+ double x_pull = x_residual / sqrt(node.fMxy.Dx2() + node.fParamDn.GetCovariance(0, 0));
+ double y_pull = y_residual / sqrt(node.fMxy.Dy2() + node.fParamDn.GetCovariance(1, 1));
if (node.fMxy.NdfX() > 0) {
hResidualsBeforeAlignmentX[0]->Fill(x_residual);
@@ -1010,22 +1018,22 @@ void CbmBbaAlignmentTask::Finish()
for (auto& t : fTracks) {
if (!t.fIsActive) continue;
- // calculate the residuals for all tracks at each fitNode before alignment
+ // calculate the residuals for all tracks at each TrajectoryNode before alignment
for (unsigned int in = 0; in < t.fAlignedTrack.fNodes.size(); in++) {
- CbmKfTrackFitter::Track tr = t.fAlignedTrack;
- CbmKfTrackFitter::FitNode& node = tr.fNodes[in];
+ CbmKfTrackFitter::Trajectory tr = t.fAlignedTrack;
+ CbmKfTrackFitter::TrajectoryNode& node = tr.fNodes[in];
if (!node.fIsXySet) {
continue;
}
node.fIsXySet = false;
- fFitter.FitTrack(tr);
+ fFitter.FitTrajectory(tr);
Sensor& s = fSensors[node.fReference1];
- double x_residual = node.fMxy.X() - node.fTrack.X();
- double y_residual = node.fMxy.Y() - node.fTrack.Y();
- double x_pull = x_residual / sqrt(node.fMxy.Dx2() + node.fTrack.GetCovariance(0, 0));
- double y_pull = y_residual / sqrt(node.fMxy.Dy2() + node.fTrack.GetCovariance(1, 1));
+ double x_residual = node.fMxy.X() - node.fParamDn.X();
+ double y_residual = node.fMxy.Y() - node.fParamDn.Y();
+ double x_pull = x_residual / sqrt(node.fMxy.Dx2() + node.fParamDn.GetCovariance(0, 0));
+ double y_pull = y_residual / sqrt(node.fMxy.Dy2() + node.fParamDn.GetCovariance(1, 1));
if (node.fMxy.NdfX() > 0) {
hResidualsAfterAlignmentX[0]->Fill(x_residual);
diff --git a/reco/alignment/CbmBbaAlignmentTask.h b/reco/alignment/CbmBbaAlignmentTask.h
index e89dd022dac1bae2bc5bb87daf6298394bb6f590..f7bb9ffd35ad892df0b3c3b12c0d08bddd8186d4 100644
--- a/reco/alignment/CbmBbaAlignmentTask.h
+++ b/reco/alignment/CbmBbaAlignmentTask.h
@@ -65,8 +65,8 @@ class CbmBbaAlignmentTask : public FairTask {
/// aligned and unaligned hits are stored in the corresponding CbmKfTrackFitter::Track objects
struct TrackContainer {
- CbmKfTrackFitter::Track fUnalignedTrack; // track before alignment
- CbmKfTrackFitter::Track fAlignedTrack; // track after alignment
+ CbmKfTrackFitter::Trajectory fUnalignedTrack; // track before alignment
+ CbmKfTrackFitter::Trajectory fAlignedTrack; // track after alignment
int fNmvdHits{0}; // number of MVD hits
int fNstsHits{0}; // number of STS hits
diff --git a/reco/kfnew/CMakeLists.txt b/reco/kfnew/CMakeLists.txt
index 4b2b4a5c7ba0c4b287bdb84692b82b03048a6be1..1bb0e057ea82fb3ab1cea07a2ec4d51bb5bf29aa 100644
--- a/reco/kfnew/CMakeLists.txt
+++ b/reco/kfnew/CMakeLists.txt
@@ -6,7 +6,10 @@ set(INCLUDE_DIRECTORIES
set(SRCS
${CMAKE_CURRENT_SOURCE_DIR}/CbmKfTarget.cxx
+ ${CMAKE_CURRENT_SOURCE_DIR}/CbmKfUtil.cxx
${CMAKE_CURRENT_SOURCE_DIR}/CbmKfTrackingSetupBuilder.cxx
+ ${CMAKE_CURRENT_SOURCE_DIR}/CbmKfTrackFitter.cxx
+ ${CMAKE_CURRENT_SOURCE_DIR}/CbmKfFitTracksTask.cxx
)
SET_SOURCE_FILES_PROPERTIES(${SRCS} PROPERTIES COMPILE_FLAGS "-O3")
@@ -73,6 +76,7 @@ install(
FILES
CbmKfOriginalField.h
CbmKfTrackingSetupBuilder.h
+ CbmKfUtil.h
CbmKfTarget.h
DESTINATION
include/
diff --git a/reco/KF/CbmKfFitTracksTask.cxx b/reco/kfnew/CbmKfFitTracksTask.cxx
similarity index 84%
rename from reco/KF/CbmKfFitTracksTask.cxx
rename to reco/kfnew/CbmKfFitTracksTask.cxx
index 07467a4a6c5b395de671f95bb12d06062e8c9bd9..be74bf398711351f47cc28fb41487a083cf4be3b 100644
--- a/reco/KF/CbmKfFitTracksTask.cxx
+++ b/reco/kfnew/CbmKfFitTracksTask.cxx
@@ -10,10 +10,8 @@
#include "CbmKfFitTracksTask.h"
-#include "CbmCaTimeSliceReader.h"
#include "CbmGlobalTrack.h"
-#include "CbmL1.h"
-#include "CbmL1Util.h"
+#include "CbmKfUtil.h"
#include "CbmMuchTrack.h"
#include "CbmMuchTrackingInterface.h"
#include "CbmMvdHit.h"
@@ -33,7 +31,7 @@
#include <iostream>
-ClassImp(CbmKfFitTracksTask);
+// ClassImp(CbmKfFitTracksTask);
namespace
{
@@ -100,24 +98,24 @@ void CbmKfFitTracksTask::Exec(Option_t* /*opt*/)
LOG(fatal) << "CbmKfFitTracksTask: null pointer to the sts track!";
return;
}
- CbmKfTrackFitter::Track t;
+ CbmKfTrackFitter::Trajectory t;
if (!fFitter.CreateMvdStsTrack(t, iTr)) {
LOG(fatal) << "CbmKfFitTracksTask: can not create the sts track for the fit! ";
return;
}
- fFitter.FitTrack(t);
+ fFitter.FitTrajectory(t);
{
- const auto& parV = t.fNodes.front().fTrack;
+ const auto& parV = t.fNodes.front().fParamUp;
cbm::algo::kf::TrackParamD parD;
parD.Set(parV, 0);
- FairTrackParam trackFirst = cbm::L1Util::ConvertTrackParam(parD);
+ FairTrackParam trackFirst = cbm::kf::ConvertTrackParam(parD);
stsTrack->SetParamFirst(&trackFirst);
}
{
- const auto& parV = t.fNodes.back().fTrack;
+ const auto& parV = t.fNodes.back().fParamDn;
cbm::algo::kf::TrackParamD parD;
parD.Set(parV, 0);
- FairTrackParam trackLast = cbm::L1Util::ConvertTrackParam(parD);
+ FairTrackParam trackLast = cbm::kf::ConvertTrackParam(parD);
stsTrack->SetParamLast(&trackLast);
}
}
@@ -134,24 +132,24 @@ void CbmKfFitTracksTask::Exec(Option_t* /*opt*/)
LOG(fatal) << "CbmKfFitTracksTask: null pointer to the global track!";
return;
}
- CbmKfTrackFitter::Track t;
+ CbmKfTrackFitter::Trajectory t;
if (!fFitter.CreateGlobalTrack(t, *globalTrack)) {
LOG(fatal) << "CbmKfFitTracksTask: can not create the global track for the fit! ";
return;
}
- fFitter.FitTrack(t);
+ fFitter.FitTrajectory(t);
{
- const auto& parV = t.fNodes.front().fTrack;
+ const auto& parV = t.fNodes.front().fParamUp;
cbm::algo::kf::TrackParamD parD;
parD.Set(parV, 0);
- FairTrackParam trackFirst = cbm::L1Util::ConvertTrackParam(parD);
+ FairTrackParam trackFirst = cbm::kf::ConvertTrackParam(parD);
globalTrack->SetParamFirst(&trackFirst);
}
{
- const auto& parV = t.fNodes.back().fTrack;
+ const auto& parV = t.fNodes.back().fParamDn;
cbm::algo::kf::TrackParamD parD;
parD.Set(parV, 0);
- FairTrackParam trackLast = cbm::L1Util::ConvertTrackParam(parD);
+ FairTrackParam trackLast = cbm::kf::ConvertTrackParam(parD);
globalTrack->SetParamLast(&trackLast);
}
}
diff --git a/reco/KF/CbmKfFitTracksTask.h b/reco/kfnew/CbmKfFitTracksTask.h
similarity index 97%
rename from reco/KF/CbmKfFitTracksTask.h
rename to reco/kfnew/CbmKfFitTracksTask.h
index 4b43d73529dca12cb8a854611770367591f8c375..dfbf2f8b3f4107dcc8edd430b680f8f181f645be 100644
--- a/reco/KF/CbmKfFitTracksTask.h
+++ b/reco/kfnew/CbmKfFitTracksTask.h
@@ -58,5 +58,5 @@ class CbmKfFitTracksTask : public FairTask {
Int_t fNeventsProcessed{0}; ///< number of processed events
- ClassDefOverride(CbmKfFitTracksTask, 0);
+ // ClassDefOverride(CbmKfFitTracksTask, 0);
};
diff --git a/reco/KF/CbmKfTrackFitter.cxx b/reco/kfnew/CbmKfTrackFitter.cxx
similarity index 61%
rename from reco/KF/CbmKfTrackFitter.cxx
rename to reco/kfnew/CbmKfTrackFitter.cxx
index 805a1cd82c0882b7a90c9c59abb21b7ed85b8bb6..14a67f0dfba7962319c4194edfeea6e1d80f8317 100644
--- a/reco/KF/CbmKfTrackFitter.cxx
+++ b/reco/kfnew/CbmKfTrackFitter.cxx
@@ -4,13 +4,8 @@
#include "CbmKfTrackFitter.h"
-#include "CaDefs.h"
-#include "CaFramework.h"
-#include "CaSimd.h"
-#include "CaStation.h"
#include "CbmGlobalTrack.h"
-#include "CbmL1.h"
-#include "CbmL1Util.h"
+#include "CbmKfTrackingSetupBuilder.h"
#include "CbmMuchPixelHit.h"
#include "CbmMuchTrack.h"
#include "CbmMuchTrackingInterface.h"
@@ -29,6 +24,7 @@
#include "CbmTrdTrackingInterface.h"
#include "FairRootManager.h"
#include "KFParticleDatabase.h"
+#include "KfDefs.h"
#include "KfTrackKalmanFilter.h"
#include "KfTrackParam.h"
#include "TClonesArray.h"
@@ -36,15 +32,12 @@
using std::vector;
using namespace std;
-using ca::fmask;
-using ca::fvec;
-
using namespace cbm::algo;
-void CbmKfTrackFitter::Track::OrderNodesInZ()
+void CbmKfTrackFitter::Trajectory::OrderNodesInZ()
{
// sort the nodes in z
- std::sort(fNodes.begin(), fNodes.end(), [](const FitNode& a, const FitNode& b) { return a.fZ < b.fZ; });
+ std::sort(fNodes.begin(), fNodes.end(), [](const TrajectoryNode& a, const TrajectoryNode& b) { return a.fZ < b.fZ; });
}
@@ -58,16 +51,20 @@ void CbmKfTrackFitter::Init()
return;
}
- if (!CbmL1::Instance() || !CbmL1::Instance()->fpAlgo) {
- LOG(fatal) << "CbmKfTrackFitter: no CbmL1 task initialized ";
- }
-
FairRootManager* ioman = FairRootManager::Instance();
if (!ioman) {
LOG(fatal) << "CbmKfTrackFitter: no FairRootManager";
}
+ // TODO: better check if all the interfaces are properly initialized
+ if (!CbmStsTrackingInterface::Instance() || !CbmTrdTrackingInterface::Instance()
+ || !CbmTofTrackingInterface::Instance()) {
+ LOG(fatal) << "CbmTrackingDetectorInterface instance was not found. Please, add it as a task to your "
+ "reco macro right before the KF and L1 tasks:\n"
+ << "\033[1;30mrun->AddTask(new CbmTrackingDetectorInterfaceInit());\033[0m";
+ }
+
// Get hits
fInputMvdHits = dynamic_cast<TClonesArray*>(ioman->GetObject("MvdHit"));
@@ -85,6 +82,8 @@ void CbmKfTrackFitter::Init()
fInputTrdTracks = dynamic_cast<TClonesArray*>(ioman->GetObject("TrdTrack"));
fInputTofTracks = dynamic_cast<TClonesArray*>(ioman->GetObject("TofTrack"));
+ fKfSetup = cbm::kf::TrackingSetupBuilder::Instance()->GetSharedGeoSetup();
+
fIsInitialized = true;
}
@@ -109,7 +108,7 @@ void CbmKfTrackFitter::SetMassHypothesis(double mass)
void CbmKfTrackFitter::SetElectronFlag(bool isElectron) { fIsElectron = isElectron; }
-bool CbmKfTrackFitter::CreateMvdStsTrack(CbmKfTrackFitter::Track& kfTrack, int stsTrackIndex)
+bool CbmKfTrackFitter::CreateMvdStsTrack(CbmKfTrackFitter::Trajectory& kfTrack, int stsTrackIndex)
{
Init();
if (!fIsInitialized) {
@@ -138,7 +137,7 @@ bool CbmKfTrackFitter::CreateMvdStsTrack(CbmKfTrackFitter::Track& kfTrack, int s
}
-bool CbmKfTrackFitter::CreateGlobalTrack(CbmKfTrackFitter::Track& kfTrack, int globalTrackIndex)
+bool CbmKfTrackFitter::CreateGlobalTrack(CbmKfTrackFitter::Trajectory& kfTrack, int globalTrackIndex)
{
Init();
if (!fIsInitialized) {
@@ -165,7 +164,7 @@ bool CbmKfTrackFitter::CreateGlobalTrack(CbmKfTrackFitter::Track& kfTrack, int g
}
-bool CbmKfTrackFitter::CreateGlobalTrack(CbmKfTrackFitter::Track& kfTrack, const CbmGlobalTrack& globalTrack)
+bool CbmKfTrackFitter::CreateGlobalTrack(CbmKfTrackFitter::Trajectory& kfTrack, const CbmGlobalTrack& globalTrack)
{
kfTrack = {};
Init();
@@ -173,39 +172,36 @@ bool CbmKfTrackFitter::CreateGlobalTrack(CbmKfTrackFitter::Track& kfTrack, const
return false;
}
- const ca::Parameters<ca::fvec>& caPar = CbmL1::Instance()->fpAlgo->GetParameters();
-
- int nStations = caPar.GetNstationsActive();
-
- CbmKfTrackFitter::Track t;
+ CbmKfTrackFitter::Trajectory t;
{
- t.fNodes.resize(nStations);
- for (int i = 0; i < nStations; i++) {
- CbmKfTrackFitter::FitNode& n = t.fNodes[i];
- n = {};
- n.fMaterialLayer = i;
- n.fZ = caPar.GetStation(i).GetZ<float>();
- n.fRadThick = 0.;
- n.fIsRadThickFixed = false;
- n.fIsFitted = false;
- n.fIsXySet = false;
- n.fIsTimeSet = false;
- n.fHitSystemId = ECbmModuleId::kNotExist;
- n.fHitAddress = -1;
- n.fHitIndex = -1;
+ t.fNodes.resize(fKfSetup->GetNofLayers());
+ for (int i = 0; i < fKfSetup->GetNofLayers(); i++) {
+ CbmKfTrackFitter::TrajectoryNode& n = t.fNodes[i];
+ n = {};
+ n.fMaterialLayer = i;
+ n.fZ = fKfSetup->GetMaterial(i).GetZref();
+ n.fRadThick = 0.;
+ n.fIsRadThickFixed = false;
+ n.fIsFitted = false;
+ n.fIsXySet = false;
+ n.fIsTimeSet = false;
+ n.fHitSystemId = ECbmModuleId::kNotExist;
+ n.fHitAddress = -1;
+ n.fHitIndex = -1;
}
}
// lambda to set the node from the pixel hit
- auto setNode = [&](const CbmPixelHit& h, int stIdx, int hitIdx, ca::EDetectorID detId) -> CbmKfTrackFitter::FitNode* {
- int ista = caPar.GetStationIndexActive(stIdx, detId);
- if (ista < 0) {
+ auto setNode = [&](const CbmPixelHit& h, int stIdx, int hitIdx,
+ ca::EDetectorID detId) -> CbmKfTrackFitter::TrajectoryNode* {
+ int iLayer = fKfSetup->GetIndexMap().LocalToGlobal(detId, stIdx);
+ if (iLayer < 0) {
return nullptr;
}
- assert(ista < nStations);
- CbmKfTrackFitter::FitNode& n = t.fNodes[ista];
- n.fZ = h.GetZ();
+ assert(iLayer >= 0 && iLayer < fKfSetup->GetNofLayers());
+ CbmKfTrackFitter::TrajectoryNode& n = t.fNodes[iLayer];
+ n.fZ = h.GetZ();
n.fMxy.SetX(h.GetX());
n.fMxy.SetY(h.GetY());
n.fMxy.SetDx2(h.GetDx() * h.GetDx());
@@ -220,7 +216,7 @@ bool CbmKfTrackFitter::CreateGlobalTrack(CbmKfTrackFitter::Track& kfTrack, const
n.fIsTimeSet = (detId != ca::EDetectorID::kMvd);
n.fRadThick = 0.;
n.fIsRadThickFixed = false;
- n.fHitSystemId = cbm::L1Util::ConvertDetSystemId(detId);
+ n.fHitSystemId = ca::ToCbmModuleId(detId);
n.fHitAddress = h.GetAddress();
n.fHitIndex = hitIdx;
return &n;
@@ -254,6 +250,10 @@ bool CbmKfTrackFitter::CreateGlobalTrack(CbmKfTrackFitter::Track& kfTrack, const
LOG(error) << "CbmKfTrackFitter: Mvd hit array not found!";
return false;
}
+ if (!CbmMvdTrackingInterface::Instance()) {
+ LOG(error) << "CbmKfTrackFitter: Mvd tracking interface not found!";
+ return false;
+ }
for (int ih = 0; ih < nMvdHits; ih++) {
int hitIndex = stsTrack->GetMvdHitIndex(ih);
if (hitIndex >= fInputMvdHits->GetEntriesFast()) {
@@ -274,6 +274,10 @@ bool CbmKfTrackFitter::CreateGlobalTrack(CbmKfTrackFitter::Track& kfTrack, const
LOG(error) << "CbmKfTrackFitter: Sts hit array not found!";
return false;
}
+ if (!CbmStsTrackingInterface::Instance()) {
+ LOG(error) << "CbmKfTrackFitter: Sts tracking interface not found!";
+ return false;
+ }
for (int ih = 0; ih < nStsHits; ih++) {
int hitIndex = stsTrack->GetStsHitIndex(ih);
if (hitIndex >= fInputStsHits->GetEntriesFast()) {
@@ -311,6 +315,10 @@ bool CbmKfTrackFitter::CreateGlobalTrack(CbmKfTrackFitter::Track& kfTrack, const
LOG(error) << "CbmKfTrackFitter: Much hit array not found!";
return false;
}
+ if (!CbmMuchTrackingInterface::Instance()) {
+ LOG(error) << "CbmKfTrackFitter: Much tracking interface not found!";
+ return false;
+ }
for (int ih = 0; ih < nHits; ih++) {
int hitIndex = track->GetHitIndex(ih);
if (hitIndex >= fInputMuchHits->GetEntriesFast()) {
@@ -347,6 +355,10 @@ bool CbmKfTrackFitter::CreateGlobalTrack(CbmKfTrackFitter::Track& kfTrack, const
LOG(error) << "CbmKfTrackFitter: Trd hit array not found!";
return false;
}
+ if (!CbmTrdTrackingInterface::Instance()) {
+ LOG(error) << "CbmKfTrackFitter: Trd tracking interface not found!";
+ return false;
+ }
for (int ih = 0; ih < nHits; ih++) {
int hitIndex = track->GetHitIndex(ih);
if (hitIndex >= fInputTrdHits->GetEntriesFast()) {
@@ -389,6 +401,10 @@ bool CbmKfTrackFitter::CreateGlobalTrack(CbmKfTrackFitter::Track& kfTrack, const
LOG(error) << "CbmKfTrackFitter: Tof hit array not found!";
return false;
}
+ if (!CbmTofTrackingInterface::Instance()) {
+ LOG(error) << "CbmKfTrackFitter: Tof tracking interface not found!";
+ return false;
+ }
for (int ih = 0; ih < nHits; ih++) {
int hitIndex = track->GetHitIndex(ih);
if (hitIndex >= fInputTofHits->GetEntriesFast()) {
@@ -404,11 +420,6 @@ bool CbmKfTrackFitter::CreateGlobalTrack(CbmKfTrackFitter::Track& kfTrack, const
t.OrderNodesInZ();
- // set initial track parameters
- // kf::TrackParamD tmp = cbm::L1Util::ConvertTrackParam(*globalTrack.GetParamFirst());
- // t.fNodes[0].fTrack.Set(tmp);
- // t.fNodes[0].fIsFitted = 1;
-
kfTrack = t;
return true;
@@ -416,7 +427,7 @@ bool CbmKfTrackFitter::CreateGlobalTrack(CbmKfTrackFitter::Track& kfTrack, const
}
-void CbmKfTrackFitter::FilterFirstMeasurement(const FitNode& n)
+void CbmKfTrackFitter::FilterFirstMeasurement(const TrajectoryNode& n)
{
// a special routine to filter the first measurement.
// the measurement errors are simply copied to the track covariance matrix
@@ -428,7 +439,7 @@ void CbmKfTrackFitter::FilterFirstMeasurement(const FitNode& n)
auto& tr = fFit.Tr();
- if (n.fIsFitted && (fabs(tr.GetZ() - n.fZ) > 1.e-10)) {
+ if (fabs(tr.GetZ() - n.fZ) > 1.e-10) {
LOG(fatal) << "CbmKfTrackFitter: Z mismatch: fitted track " << tr.GetZ() << " != node " << n.fZ;
}
@@ -438,13 +449,13 @@ void CbmKfTrackFitter::FilterFirstMeasurement(const FitNode& n)
tr.SetY(mxy.Y());
tr.SetZ(n.fZ);
- if (fSkipUnmeasuredCoordinates && n.fIsFitted) {
+ if (fSkipUnmeasuredCoordinates) { // TODO: take C10 correlation into account
if (mxy.NdfX() == 0) {
- tr.SetX(n.fTrack.GetX());
+ tr.SetX(n.fParamDn.GetX());
tr.SetC00(1.e4);
}
if (mxy.NdfY() == 0) {
- tr.SetY(n.fTrack.GetY());
+ tr.SetY(n.fParamDn.GetY());
tr.SetC11(1.e4);
tr.SetC10(0.);
}
@@ -452,7 +463,7 @@ void CbmKfTrackFitter::FilterFirstMeasurement(const FitNode& n)
tr.SetChiSq(0.);
tr.SetChiSqTime(0.);
- tr.SetNdf(-5 + mxy.NdfX() + mxy.NdfY());
+ tr.SetNdf(-5. + mxy.NdfX() + mxy.NdfY());
if (n.fIsTimeSet) {
tr.SetTime(mt.T());
tr.SetC55(mt.Dt2());
@@ -461,43 +472,52 @@ void CbmKfTrackFitter::FilterFirstMeasurement(const FitNode& n)
else {
tr.SetNdfTime(-2);
}
- tr.Vi() = constants::phys::SpeedOfLightInv;
+ tr.Vi() = cbm::algo::kf::defs::SpeedOfLightInv<>;
tr.InitVelocityRange(0.5);
}
-void CbmKfTrackFitter::AddMaterialEffects(CbmKfTrackFitter::FitNode& n, kf::FitDirection direction)
+void CbmKfTrackFitter::AddMaterialEffects(CbmKfTrackFitter::TrajectoryNode& n, const LinearizationAtNode& l,
+ kf::FitDirection direction)
{
// add material effects
if (n.fMaterialLayer < 0) {
return;
}
- auto& tr = n.fIsFitted ? n.fTrack : fFit.Tr();
+ double tx = 0.5 * (l.fParamDn.GetTx() + l.fParamUp.GetTx());
+ double ty = 0.5 * (l.fParamDn.GetTy() + l.fParamUp.GetTy());
+ double msQp = 0.5 * (l.fParamDn.GetQp() + l.fParamUp.GetQp());
- double msQp0 = n.fIsFitted ? n.fTrack.GetQp() : fFit.Qp0();
if (fIsQpForMsFixed) {
- msQp0 = fDefaultQpForMs;
+ msQp = fDefaultQpForMs;
}
// calculate the radiation thickness from the current track
if (!n.fIsRadThickFixed) {
- n.fRadThick =
- CbmL1::Instance()->fpAlgo->GetParameters().GetMaterialThicknessScal(n.fMaterialLayer, tr.GetX(), tr.GetY());
+ const kf::MaterialMap& map = fKfSetup->GetMaterial(n.fMaterialLayer);
+ n.fRadThick = map.GetThicknessX0(l.fParamDn.GetX(), l.fParamDn.GetY());
}
- fFit.MultipleScattering(n.fRadThick, tr.GetTx(), tr.GetTy(), msQp0);
+ fFit.MultipleScattering(n.fRadThick, tx, ty, msQp);
+
if (!fIsQpForMsFixed) {
+ if (direction == kf::FitDirection::kDownstream) {
+ fFit.SetQp0(l.fParamUp.GetQp());
+ }
+ else {
+ fFit.SetQp0(l.fParamDn.GetQp());
+ }
fFit.EnergyLossCorrection(n.fRadThick, direction);
}
}
-bool CbmKfTrackFitter::FitTrack(CbmKfTrackFitter::Track& t)
+bool CbmKfTrackFitter::FitTrajectory(CbmKfTrackFitter::Trajectory& t)
{
// (re)fit the track
if (fVerbosityLevel > 0) {
- std::cout << "FitTrack ... " << std::endl;
+ std::cout << "FitTrajectory ... " << std::endl;
}
bool ok = true;
@@ -518,16 +538,17 @@ bool CbmKfTrackFitter::FitTrack(CbmKfTrackFitter::Track& t)
{ // find the first and the last hit nodes. Check if the nodes are ordered in Z
double zOld = -1.e10;
for (int i = 0; i < nNodes; i++) {
- if (t.fNodes[i].fIsXySet) {
+ const auto& n = t.fNodes[i];
+ if (n.fIsXySet) {
if (firstHitNode < 0) {
firstHitNode = i;
}
lastHitNode = i;
}
- if (t.fNodes[i].fZ < zOld) {
+ if (n.fZ < zOld) {
isOrdered = false;
}
- zOld = t.fNodes[i].fZ;
+ zOld = n.fZ;
}
}
@@ -546,164 +567,214 @@ bool CbmKfTrackFitter::FitTrack(CbmKfTrackFitter::Track& t)
kf::FieldRegion<double> field(kf::GlobalField::fgOriginalFieldType, kf::GlobalField::fgOriginalField);
- { // fit downstream
-
- { // initialisation at the first measurement
- FitNode& n = t.fNodes[firstHitNode];
- assert(n.fIsXySet);
-
- if (n.fIsFitted) { // The initial parameters are taken from the previous fit
- // Linearisation at the previously fitted momentum
- //if (n.fZ != n.fTrack.GetZ()) {
- //LOG(ERROR) << "Z " << n.fZ << " track Z " << n.fTrack.GetZ() << " dz " << n.fZ - n.fTrack.GetZ();
- //}
- //assert(n.fZ == n.fTrack.GetZ());
- fFit.SetTrack(n.fTrack);
- }
- else { // The initial parameters are calculated from the first and the last measurements
- if (!fDoSmooth) { //SG!!!
- LOG(fatal) << "CbmKfTrackFitter: Debug: downstream: the first measurement is not fitted";
+ std::vector<LinearizationAtNode> linearization(nNodes);
+
+ if (t.fIsFitted) { // take the linearization from the previously fitted trajectory
+ for (int i = firstHitNode; i <= lastHitNode; i++) {
+ if (!t.fNodes[i].fIsFitted) {
+ LOG(fatal) << "CbmKfTrackFitter: node " << i << " in the measured region is not fitted";
+ }
+ linearization[i].fParamDn = t.fNodes[i].fParamDn;
+ linearization[i].fParamUp = t.fNodes[i].fParamUp;
+ }
+ }
+ else { // first approximation with straight line segments connecting the nodes
+ for (int i1 = firstHitNode, i2 = firstHitNode + 1; i2 <= lastHitNode; i2++) {
+ auto& n1 = t.fNodes[i1];
+ auto& n2 = t.fNodes[i2];
+ if (!n2.fIsXySet) {
+ continue;
+ }
+ double dz = n2.fZ - n1.fZ;
+ double dzi = (fabs(dz) > 1.e-4) ? 1. / dz : 0.;
+ double tx = (n2.fMxy.X() - n1.fMxy.X()) * dzi;
+ double ty = (n2.fMxy.Y() - n1.fMxy.Y()) * dzi;
+ for (int i = i1; i <= i2; i++) { // fill the nodes inbetween
+ auto& n = t.fNodes[i];
+ auto& l = linearization[i];
+ double x = n1.fMxy.X() + tx * (n.fZ - n1.fZ);
+ double y = n1.fMxy.Y() + ty * (n.fZ - n1.fZ);
+ if (i < i2
+ || i == lastHitNode) { // downstream parameters for i2 will be set later, except for the last hit node
+ l.fParamDn.SetX(x);
+ l.fParamDn.SetY(y);
+ l.fParamDn.SetZ(n.fZ);
+ l.fParamDn.SetTx(tx);
+ l.fParamDn.SetTy(ty);
+ l.fParamDn.SetQp(0.);
+ l.fParamDn.SetTime(0.);
+ l.fParamDn.SetVi(cbm::algo::kf::defs::SpeedOfLightInv<>);
+ }
+ if (i > i1 || i == firstHitNode) { // upstream parameters for i1 are already set, except for the first hit node
+ l.fParamUp.SetX(x);
+ l.fParamUp.SetY(y);
+ l.fParamUp.SetZ(n.fZ);
+ l.fParamUp.SetTx(tx);
+ l.fParamUp.SetTy(ty);
+ l.fParamUp.SetQp(0.);
+ l.fParamUp.SetTime(0.);
+ l.fParamUp.SetVi(cbm::algo::kf::defs::SpeedOfLightInv<>);
}
- auto& tr = fFit.Tr();
- auto& n1 = t.fNodes[lastHitNode];
- tr.X() = n.fMxy.X();
- tr.Y() = n.fMxy.Y();
- tr.Z() = n.fZ;
- double dz = n1.fZ - n.fZ;
- tr.Tx() = (n1.fMxy.X() - n.fMxy.X()) / dz;
- tr.Ty() = (n1.fMxy.Y() - n.fMxy.Y()) / dz;
- tr.Qp() = 0.;
- tr.Time() = 0.;
- tr.Vi() = constants::phys::SpeedOfLightInv;
- fFit.SetQp0(0.); // Linearisation at the infinite momentum
}
+ i1 = i2;
+ }
+ }
- FilterFirstMeasurement(n);
+ t.fIsFitted = false;
+ for (auto& n : t.fNodes) {
+ n.fIsFitted = false;
+ }
- if (fVerbosityLevel > 0) {
- std::cout << " fit downstream: node " << firstHitNode << " chi2 " << fFit.Tr().GetChiSq() << " x "
- << fFit.Tr().GetX() << " y " << fFit.Tr().GetY() << " z " << fFit.Tr().GetZ() << " tx "
- << fFit.Tr().GetTx() << " ty " << fFit.Tr().GetTy() << std::endl;
- }
+ auto printNode = [&](std::string str, int node) {
+ if (fVerbosityLevel > 0) {
+ LOG(info) << str << ": node " << node << " chi2 " << fFit.Tr().GetChiSq() << " x " << fFit.Tr().GetX() << " y "
+ << fFit.Tr().GetY() << " z " << fFit.Tr().GetZ() << " tx " << fFit.Tr().GetTx() << " ty "
+ << fFit.Tr().GetTy();
}
+ };
- for (int iNode = firstHitNode + 1; iNode < nNodes; iNode++) {
-
- FitNode& n = t.fNodes[iNode];
- if (!fDoSmooth && !n.fIsFitted) { //SG!!!
- LOG(fatal) << "CbmKfTrackFitter: Debug: downstream: node " << iNode << " is not fitted";
- }
+ { // fit downstream up to the last measurement
- fFit.Extrapolate(n.fZ, field);
+ { // initialisation at the first measurement
+ TrajectoryNode& n = t.fNodes[firstHitNode];
+ assert(n.fIsXySet);
+ fFit.SetTrack(linearization[firstHitNode].fParamDn);
+ FilterFirstMeasurement(n);
+ printNode("fit downstream", firstHitNode);
+ }
- if (n.fIsFitted && (iNode <= lastHitNode)) { // linearisation is taken from the previous fit
- fFit.SetQp0(n.fTrack.GetQp());
- }
+ for (int iNode = firstHitNode + 1; iNode <= lastHitNode; iNode++) {
- if (fDoSmooth) { // store the partially fitted track
- n.fTrack = fFit.Tr();
- n.fIsFitted = false; // the stored track is not fully fitted
- }
+ TrajectoryNode& n = t.fNodes[iNode];
- AddMaterialEffects(n, kf::FitDirection::kDownstream);
+ // transport the partially fitted track to the current node
+ fFit.SetQp0(linearization[iNode - 1].fParamDn.GetQp());
+ fFit.Extrapolate(n.fZ, field);
+ // measure the track at the current node
if (n.fIsXySet) {
- assert(iNode <= lastHitNode);
fFit.FilterXY(n.fMxy, fSkipUnmeasuredCoordinates);
}
if (n.fIsTimeSet) {
fFit.FilterTime(n.fMt);
}
+ printNode("fit downstream", iNode);
- if (fVerbosityLevel > 0) {
- std::cout << " fit downstream: node " << iNode << " chi2 " << fFit.Tr().GetChiSq() << " x " << fFit.Tr().GetX()
- << " y " << fFit.Tr().GetY() << " z " << fFit.Tr().GetZ() << " tx " << fFit.Tr().GetTx() << " ty "
- << fFit.Tr().GetTy() << std::endl;
- }
- if (fDoSmooth && iNode >= lastHitNode) { // the track is fully fitted, store it
- n.fTrack = fFit.Tr();
- n.fIsFitted = true;
- fFit.SetQp0(fFit.Tr().GetQp()); // set the linearisation of q/p to the fitted value
- }
+ // store the partially fitted track before the scattering
+ n.fParamUp = fFit.Tr();
+
+ // add material effects
+ AddMaterialEffects(n, linearization[iNode], kf::FitDirection::kDownstream);
+
+ // store the partially fitted track after the scattering
+ n.fParamDn = fFit.Tr();
}
}
+ // store the chi2 for debugging purposes
double dnChi2 = fFit.Tr().GetChiSq();
- double dnNdf = fFit.Tr().GetNdf();
-
- { // fit upstream
- {
- FitNode& n = t.fNodes[lastHitNode];
- if (!fDoSmooth && !n.fIsFitted) {
- LOG(fatal) << "CbmKfTrackFitter: Debug: upstream: the last measurement is not fitted";
- }
+ { // fit upstream with the Kalman Filter smoothing
- assert(n.fIsFitted);
- fFit.SetTrack(n.fTrack);
+ { // initialisation at the last measurement
+ TrajectoryNode& n = t.fNodes[lastHitNode];
+ assert(n.fIsXySet);
+ n.fIsFitted = true;
+ fFit.SetTrack(linearization[lastHitNode].fParamUp);
FilterFirstMeasurement(n);
- if (fVerbosityLevel > 0) {
- std::cout << "\n\n fit upstream: node " << lastHitNode << " chi2 " << fFit.Tr().GetChiSq() << " x "
- << fFit.Tr().GetX() << " y " << fFit.Tr().GetY() << " z " << fFit.Tr().GetZ() << " tx "
- << fFit.Tr().GetTx() << " ty " << fFit.Tr().GetTy() << std::endl;
- }
+ printNode("fit upstream", lastHitNode);
}
- for (int iNode = lastHitNode - 1; iNode >= 0; iNode--) {
- FitNode& n = t.fNodes[iNode];
+ for (int iNode = lastHitNode - 1; ok && (iNode > firstHitNode); iNode--) {
- if (!fDoSmooth && !n.fIsFitted) {
- LOG(fatal) << "CbmKfTrackFitter: Debug: upstream: the node " << iNode << "is not fitted";
- }
+ TrajectoryNode& n = t.fNodes[iNode];
+ // transport the partially fitted track to the current node
+ fFit.SetQp0(linearization[iNode + 1].fParamUp.GetQp());
fFit.Extrapolate(n.fZ, field);
+ // combine partially fitted downstream and upstream tracks
+ ok = ok && Smooth(n.fParamDn, fFit.Tr());
+
+ AddMaterialEffects(n, linearization[iNode], kf::FitDirection::kUpstream);
+
+ // combine partially fitted downstream and upstream tracks
+ ok = ok && Smooth(n.fParamUp, fFit.Tr());
+
+ n.fIsFitted = true;
+
+ // measure the track at the current node
if (n.fIsXySet) {
fFit.FilterXY(n.fMxy, fSkipUnmeasuredCoordinates);
}
if (n.fIsTimeSet) {
fFit.FilterTime(n.fMt);
}
+ printNode("fit upstream", iNode);
+ }
- AddMaterialEffects(n, kf::FitDirection::kUpstream);
+ if (ok) {
+ int iNode = firstHitNode;
- if (fVerbosityLevel > 0) {
- std::cout << " fit upstream: node " << iNode << " chi2 " << fFit.Tr().GetChiSq() << " x " << fFit.Tr().GetX()
- << " y " << fFit.Tr().GetY() << " z " << fFit.Tr().GetZ() << " tx " << fFit.Tr().GetTx() << " ty "
- << fFit.Tr().GetTy() << std::endl;
- }
+ TrajectoryNode& n = t.fNodes[iNode];
- // combine partially fitted downstream and upstream tracks
+ // transport the partially fitted track to the current node
+ fFit.SetQp0(linearization[iNode + 1].fParamUp.GetQp());
+ fFit.Extrapolate(n.fZ, field);
- if (iNode > firstHitNode) {
- if (fDoSmooth) {
- n.fIsFitted = Smooth(n.fTrack, fFit.Tr());
- }
+ // measure the track at the current node
+ if (n.fIsXySet) {
+ fFit.FilterXY(n.fMxy, fSkipUnmeasuredCoordinates);
}
- else {
- if (fDoSmooth) {
- n.fIsFitted = true;
- n.fTrack = fFit.Tr();
- }
+ if (n.fIsTimeSet) {
+ fFit.FilterTime(n.fMt);
}
+ printNode("fit upstream", iNode);
+ n.fParamDn = fFit.Tr();
- if (fDoSmooth && !n.fIsFitted) {
- ok = false;
- }
- if (n.fIsFitted) {
- fFit.SetQp0(n.fTrack.GetQp());
- }
+ AddMaterialEffects(n, linearization[iNode], kf::FitDirection::kUpstream);
+
+ n.fParamUp = fFit.Tr();
+
+ n.fIsFitted = true;
}
}
- double upChi2 = fFit.Tr().GetChiSq();
- double upNdf = fFit.Tr().GetNdf();
+ if (ok) { // propagate downstream
+ fFit.SetTrack(t.fNodes[lastHitNode].fParamDn);
+ for (int i = lastHitNode + 1; i < nNodes; i++) {
+ auto& n = t.fNodes[i];
+ fFit.Extrapolate(n.fZ, field);
+ n.fParamDn = fFit.Tr();
+ LinearizationAtNode l;
+ l.fParamDn = fFit.Tr();
+ l.fParamUp = fFit.Tr();
+ AddMaterialEffects(n, l, kf::FitDirection::kDownstream);
+ n.fParamUp = fFit.Tr();
+ n.fIsFitted = true;
+ }
+ }
- if (!fDoSmooth) {
- if (upNdf != dnNdf) {
- LOG(fatal) << "CbmKfTrackFitter: ndf mismatch: dn " << dnNdf << " != up " << upNdf;
+ if (ok) { // propagate upstream
+ fFit.SetTrack(t.fNodes[firstHitNode].fParamUp);
+ for (int i = firstHitNode - 1; i >= 0; i--) {
+ auto& n = t.fNodes[i];
+ fFit.Extrapolate(n.fZ, field);
+ n.fParamUp = fFit.Tr();
+ LinearizationAtNode l;
+ l.fParamDn = fFit.Tr();
+ l.fParamUp = fFit.Tr();
+ AddMaterialEffects(n, l, kf::FitDirection::kUpstream);
+ n.fParamDn = fFit.Tr();
+ n.fIsFitted = true;
}
+ }
+
+ assert(ok);
+
+ if (!fDoSmooth) {
+ double upChi2 = fFit.Tr().GetChiSq();
if (fabs(upChi2 - dnChi2) > 1.e-1) {
//if (upChi2 > dnChi2 + 1.e-2) {
LOG(error) << "CbmKfTrackFitter: " << fDebugInfo << ": chi2 mismatch: dn " << dnChi2 << " != up " << upChi2
@@ -717,17 +788,15 @@ bool CbmKfTrackFitter::FitTrack(CbmKfTrackFitter::Track& t)
const auto& tt = fFit.Tr();
for (int iNode = 0; iNode < nNodes; iNode++) {
- FitNode& n = t.fNodes[iNode];
- if (iNode == lastHitNode) {
- //continue;
- }
- if (!fDoSmooth && iNode != lastHitNode) {
- //n.fIsFitted = false;
- }
- n.fTrack.SetNdf(tt.GetNdf());
- n.fTrack.SetNdfTime(tt.GetNdfTime());
- n.fTrack.SetChiSq(tt.GetChiSq());
- n.fTrack.SetChiSqTime(tt.GetChiSqTime());
+ TrajectoryNode& n = t.fNodes[iNode];
+ n.fParamDn.SetNdf(tt.GetNdf());
+ n.fParamDn.SetNdfTime(tt.GetNdfTime());
+ n.fParamDn.SetChiSq(tt.GetChiSq());
+ n.fParamDn.SetChiSqTime(tt.GetChiSqTime());
+ n.fParamUp.SetNdf(tt.GetNdf());
+ n.fParamUp.SetNdfTime(tt.GetNdfTime());
+ n.fParamUp.SetChiSq(tt.GetChiSq());
+ n.fParamUp.SetChiSqTime(tt.GetChiSqTime());
}
return ok;
diff --git a/reco/KF/CbmKfTrackFitter.h b/reco/kfnew/CbmKfTrackFitter.h
similarity index 67%
rename from reco/KF/CbmKfTrackFitter.h
rename to reco/kfnew/CbmKfTrackFitter.h
index 7449b49423fd1f856e8cadf5c77fe1e7c931d8da..7a484d9346635917fc7ccbab7e974994a29bbd4b 100644
--- a/reco/KF/CbmKfTrackFitter.h
+++ b/reco/kfnew/CbmKfTrackFitter.h
@@ -5,11 +5,11 @@
#pragma once // include this header only once per compilation unit
-#include "CaDefs.h"
-#include "CaSimd.h"
#include "CbmDefs.h"
+#include "KfDefs.h"
#include "KfMeasurementTime.h"
#include "KfMeasurementXy.h"
+#include "KfSetup.h"
#include "KfTrackKalmanFilter.h"
#include "KfTrackParam.h"
@@ -27,26 +27,20 @@ class CbmTofHit;
class FairTrackParam;
class TClonesArray;
-namespace cbm::algo::ca
-{
- class FieldRegionV;
-} // namespace cbm::algo::ca
-
-namespace
-{
- using namespace cbm::algo;
-}
/// A fitter for the Cbm tracks
///
class CbmKfTrackFitter {
+
public:
+ /// The class to fit the tracks with the Kalman Filter with explicit linearisation of the equations
+ ///
/// A node on the trajectory where the track parameters are:
/// a) measured and / or
/// b) scattered and / or
/// c) need to be estimated
/// The nodes are expected to be ordered by increasing Z.
- /// It can be done via CbmKfTrackFitter::Track::MakeConsistent() method.
+ /// It can be done via CbmKfTrackFitter::Trajectory::MakeConsistent() method.
///
/// When fIsFitted flag is set (e.g. by the fitter), the node track parameters are used
/// for the trajectory linearisation during the fit.
@@ -55,11 +49,12 @@ class CbmKfTrackFitter {
///
/// TODO: proper interface and description
///
- struct FitNode {
+ struct TrajectoryNode {
- double fZ{0.}; ///< Z coordinate
+ double fZ{0.}; ///< Z coordinate of the node
- kf::TrackParamD fTrack{}; ///< fitted track
+ cbm::algo::kf::TrackParamD fParamUp{}; ///< fitted track parameters upstream the node
+ cbm::algo::kf::TrackParamD fParamDn{}; ///< fitted track parameters downstream the node
/// == Material information (if present)
@@ -72,9 +67,9 @@ class CbmKfTrackFitter {
/// == Hit information ( if present )
- kf::MeasurementXy<double> fMxy{}; ///< XY-measurement at fZ
+ cbm::algo::kf::MeasurementXy<double> fMxy{}; ///< XY-measurement at fZ
- kf::MeasurementTime<double> fMt{}; ///< time measurement at fZ
+ cbm::algo::kf::MeasurementTime<double> fMt{}; ///< time measurement at fZ
/// == Flags etc
bool fIsXySet{false}; ///< true if the XY measurement is set
@@ -91,11 +86,12 @@ class CbmKfTrackFitter {
int fReference2{-1}; ///< some reference can be set by the user
};
- /// A track to be fitted
+ /// A trajectory to be fitted
// TODO: proper interface
- struct Track {
- std::vector<FitNode> fNodes; ///< nodes on the track
- void OrderNodesInZ(); // sort the nodes in Z and set fFirstHitNode and fLastHitNode
+ struct Trajectory {
+ std::vector<TrajectoryNode> fNodes; ///< nodes on the trajectory
+ bool fIsFitted{false}; ///< true if the trajectory is fitted
+ void OrderNodesInZ(); // sort the nodes in Z and set fFirstHitNode and fLastHitNode
};
CbmKfTrackFitter();
@@ -130,12 +126,12 @@ class CbmKfTrackFitter {
void FixMomentumForMs(bool fix = true) { fIsQpForMsFixed = fix; }
/// set the input data arrays
- bool CreateMvdStsTrack(Track& kfTrack, int stsTrackIndex);
- bool CreateGlobalTrack(Track& kfTrack, int globalTrackIndex);
- bool CreateGlobalTrack(Track& kfTrack, const CbmGlobalTrack& globalTrack);
+ bool CreateMvdStsTrack(Trajectory& kfTrack, int stsTrackIndex);
+ bool CreateGlobalTrack(Trajectory& kfTrack, int globalTrackIndex);
+ bool CreateGlobalTrack(Trajectory& kfTrack, const CbmGlobalTrack& globalTrack);
/// fit the track
- bool FitTrack(CbmKfTrackFitter::Track& t);
+ bool FitTrajectory(CbmKfTrackFitter::Trajectory& t);
/// do the KF-smoothing to define track pars at all the nodes
void SetDoSmooth(bool doSmooth) { fDoSmooth = doSmooth; }
@@ -144,17 +140,24 @@ class CbmKfTrackFitter {
void SetVerbosityLevel(int level) { fVerbosityLevel = level; }
/// set information about the track for debug output
- void SetDebugInfo(const std::string &info) { fDebugInfo = info; }
+ void SetDebugInfo(const std::string& info) { fDebugInfo = info; }
private:
- void FilterFirstMeasurement(const FitNode& n);
+ struct LinearizationAtNode {
+ cbm::algo::kf::TrackParamD fParamUp{}; ///< fitted track parameters upstream the node
+ cbm::algo::kf::TrackParamD fParamDn{}; ///< fitted track parameters downstream the node
+ };
+
+ void FilterFirstMeasurement(const TrajectoryNode& n);
- void AddMaterialEffects(FitNode& n, kf::FitDirection direction);
+ void AddMaterialEffects(TrajectoryNode& n, const LinearizationAtNode& l, cbm::algo::kf::FitDirection direction);
// combine two tracks
- bool Smooth(kf::TrackParamD& t1, const kf::TrackParamD& t2);
+ bool Smooth(cbm::algo::kf::TrackParamD& t1, const cbm::algo::kf::TrackParamD& t2);
private:
+ std::shared_ptr<const cbm::algo::kf::Setup<double>> fKfSetup;
+
// input data arrays
TClonesArray* fInputMvdHits{nullptr};
TClonesArray* fInputStsHits{nullptr};
@@ -172,7 +175,7 @@ class CbmKfTrackFitter {
bool fIsInitialized = {false}; // is the fitter initialized
bool fSkipUnmeasuredCoordinates{false};
- kf::TrackKalmanFilter<double> fFit; // track fit object
+ cbm::algo::kf::TrackKalmanFilter<double> fFit; // track fit object
/// externally defined inverse momentum for the Multiple Scattering calculation.
/// It is used for the tracks in field-free regions.
@@ -181,9 +184,9 @@ class CbmKfTrackFitter {
double fDefaultQpForMs{1. / 0.1};
bool fIsQpForMsFixed{false};
- double fMass{kf::defs::PionMass<double>}; // mass hypothesis for the fit
- bool fIsElectron{false}; // fit track as an electron (with the bermsstrallung effect)
- bool fDoSmooth{true}; // do the KF-smoothing to define track pars at all the nodes
- int fVerbosityLevel{0}; // verbosity level
- std::string fDebugInfo{}; // debug info
+ double fMass{cbm::algo::kf::defs::PionMass<double>}; // mass hypothesis for the fit
+ bool fIsElectron{false}; // fit track as an electron (with the bermsstrallung effect)
+ bool fDoSmooth{true}; // do the KF-smoothing to define track pars at all the nodes
+ int fVerbosityLevel{0}; // verbosity level
+ std::string fDebugInfo{}; // debug info
};
diff --git a/reco/kfnew/CbmKfTrackingSetupBuilder.cxx b/reco/kfnew/CbmKfTrackingSetupBuilder.cxx
index 356242cbcfea78e67d2780e936c83550ec486347..1b4a0d62d675531ae2e0faf6d31dd1871780a49a 100644
--- a/reco/kfnew/CbmKfTrackingSetupBuilder.cxx
+++ b/reco/kfnew/CbmKfTrackingSetupBuilder.cxx
@@ -60,7 +60,7 @@ void TrackingSetupBuilder::CheckDetectorPresence()
// ---------------------------------------------------------------------------------------------------------------------
//
-const std::shared_ptr<cbm::algo::kf::Setup<double>> TrackingSetupBuilder::GetSharedGeoSetup()
+std::shared_ptr<const cbm::algo::kf::Setup<double>> TrackingSetupBuilder::GetSharedGeoSetup()
{
using cbm::algo::kf::EFieldMode;
using cbm::algo::kf::Setup;
diff --git a/reco/kfnew/CbmKfTrackingSetupBuilder.h b/reco/kfnew/CbmKfTrackingSetupBuilder.h
index 72f240b75f10764df6aba8e85eb4a7ddce082808..00ab71c5a5f5097796591ce43c6617cade085ef4 100644
--- a/reco/kfnew/CbmKfTrackingSetupBuilder.h
+++ b/reco/kfnew/CbmKfTrackingSetupBuilder.h
@@ -28,7 +28,7 @@ namespace cbm::kf
/// \brief Gets a shared pointer to the geometry setup
/// \note The original magnetic field is defined.
/// \note Use-cases: precise fit in physical analyses and QA.
- const std::shared_ptr<cbm::algo::kf::Setup<double>> GetSharedGeoSetup();
+ std::shared_ptr<const cbm::algo::kf::Setup<double>> GetSharedGeoSetup();
/// \brief Makes setup object
/// \param fldMode Field mode (kf::EFiledMode)
diff --git a/reco/L1/CbmL1Util.cxx b/reco/kfnew/CbmKfUtil.cxx
similarity index 88%
rename from reco/L1/CbmL1Util.cxx
rename to reco/kfnew/CbmKfUtil.cxx
index 745f6db483592ea59087bedf0390d5ee2106e439..8b627b8e53b4919b2839c75d7e29c18cef8e1b0c 100644
--- a/reco/L1/CbmL1Util.cxx
+++ b/reco/kfnew/CbmKfUtil.cxx
@@ -2,12 +2,12 @@
SPDX-License-Identifier: GPL-3.0-only
Authors: Sergey Gorbunov [committer] */
-#include "CbmL1Util.h"
+#include "CbmKfUtil.h"
-#include "CaDefs.h"
#include "FairTrackParam.h"
+#include "KfDefs.h"
-namespace cbm::L1Util
+namespace cbm::kf
{
cbm::algo::kf::TrackParamD ConvertTrackParam(const FairTrackParam& par)
{
@@ -19,7 +19,7 @@ namespace cbm::L1Util
t.Qp() = par.GetQp();
t.Z() = par.GetZ();
t.Time() = 0.;
- t.Vi() = cbm::algo::ca::constants::phys::SpeedOfLightInv;
+ t.Vi() = cbm::algo::kf::defs::SpeedOfLightInv<>;
t.CovMatrix().fill(0.);
@@ -58,4 +58,4 @@ namespace cbm::L1Util
}
-} // namespace cbm::L1Util
+} // namespace cbm::kf
diff --git a/reco/L1/CbmL1Util.h b/reco/kfnew/CbmKfUtil.h
similarity index 51%
rename from reco/L1/CbmL1Util.h
rename to reco/kfnew/CbmKfUtil.h
index aefd392a0ed2c2fe0362a772def9cab7ae905569..b81a8612d116431a94744c976994748a3a3228b7 100644
--- a/reco/L1/CbmL1Util.h
+++ b/reco/kfnew/CbmKfUtil.h
@@ -2,11 +2,10 @@
SPDX-License-Identifier: GPL-3.0-only
Authors: Sergey Gorbunov [committer] */
-#ifndef CbmL1Util_H
-#define CbmL1Util_H 1
+#ifndef CbmKfUtil_H
+#define CbmKfUtil_H 1
#include "CbmDefs.h"
-#include "CbmL1DetectorID.h"
#include "CbmMuchTrackingInterface.h"
#include "CbmMvdTrackingInterface.h"
#include "CbmStsTrackingInterface.h"
@@ -18,9 +17,9 @@
class FairTrackParam;
///
-/// Collection of useful utilites for CbmL1
+/// Collection of useful utilites for CbmKf
///
-namespace cbm::L1Util
+namespace cbm::kf
{
/// copy fair track param to Ca track param
@@ -29,32 +28,6 @@ namespace cbm::L1Util
/// copy Ca track param to fair track param
FairTrackParam ConvertTrackParam(const cbm::algo::kf::TrackParamD& t);
- /// convert Ca detector Id to Cbm detector Id
- inline ECbmModuleId ConvertDetSystemId(const cbm::algo::ca::EDetectorID caDetId)
- {
- switch (caDetId) {
- case cbm::algo::ca::EDetectorID::kMvd: return ECbmModuleId::kMvd;
- case cbm::algo::ca::EDetectorID::kSts: return ECbmModuleId::kSts;
- case cbm::algo::ca::EDetectorID::kMuch: return ECbmModuleId::kMuch;
- case cbm::algo::ca::EDetectorID::kTrd: return ECbmModuleId::kTrd;
- case cbm::algo::ca::EDetectorID::kTof: return ECbmModuleId::kTof;
- default: return ECbmModuleId::kNotExist;
- }
- }
-
- /// convert Cbm detector Id to Ca detector Id
- inline cbm::algo::ca::EDetectorID ConvertDetSystemId(const ECbmModuleId cbmDetId)
- {
- switch (cbmDetId) {
- case ECbmModuleId::kMvd: return cbm::algo::ca::EDetectorID::kMvd;
- case ECbmModuleId::kSts: return cbm::algo::ca::EDetectorID::kSts;
- case ECbmModuleId::kMuch: return cbm::algo::ca::EDetectorID::kMuch;
- case ECbmModuleId::kTrd: return cbm::algo::ca::EDetectorID::kTrd;
- case ECbmModuleId::kTof: return cbm::algo::ca::EDetectorID::kTof;
- default: return cbm::algo::ca::EDetectorID::END;
- }
- }
-
inline const CbmTrackingDetectorInterfaceBase* GetTrackingInterface(const cbm::algo::ca::EDetectorID caDetId)
{
switch (caDetId) {
@@ -69,10 +42,10 @@ namespace cbm::L1Util
inline const CbmTrackingDetectorInterfaceBase* GetTrackingInterface(const ECbmModuleId cbmDetId)
{
- cbm::algo::ca::EDetectorID caDetId = ConvertDetSystemId(cbmDetId);
+ cbm::algo::ca::EDetectorID caDetId = cbm::algo::ca::ToCaDetectorID(cbmDetId);
return GetTrackingInterface(caDetId);
}
-} // namespace cbm::L1Util
+} // namespace cbm::kf
#endif
diff --git a/reco/qa/CbmRecoQaTask.cxx b/reco/qa/CbmRecoQaTask.cxx
index bbc881af9db2537f22baca29e8003a6390e5b33c..e5617fe623cca27f4f483a913d0ff366b28474f5 100644
--- a/reco/qa/CbmRecoQaTask.cxx
+++ b/reco/qa/CbmRecoQaTask.cxx
@@ -39,6 +39,8 @@
#include <regex>
using namespace std;
+using namespace cbm::algo;
+
std::bitset<kRecoQaNConfigs> CbmRecoQaTask::fuRecoConfig = {};
//_____________________________________________________________________
@@ -373,11 +375,12 @@ bool CbmRecoQaTask::Detector::View::Load(const CbmHit* h, const FairMCPoint* poi
return true;
}
-bool CbmRecoQaTask::Detector::View::Load(const CbmKfTrackFitter::FitNode* n, const FairMCPoint* point)
+bool CbmRecoQaTask::Detector::View::Load(const CbmKfTrackFitter::TrajectoryNode* n, const FairMCPoint* point)
{
- double dx = n->fMxy.X() - n->fTrack.X(), dy = n->fMxy.Y() - n->fTrack.Y(), dt = n->fMt.T() - n->fTrack.Time(),
- pullx = dx / sqrt(n->fMxy.Dx2() + n->fTrack.GetCovariance(0, 0)),
- pully = dy / sqrt(n->fMxy.Dy2() + n->fTrack.GetCovariance(1, 1));
+ const kf::TrackParamD& t = n->fParamUp;
+ double dx = n->fMxy.X() - t.X(), dy = n->fMxy.Y() - t.Y(), dt = n->fMt.T() - t.Time(),
+ pullx = dx / sqrt(n->fMxy.Dx2() + t.GetCovariance(0, 0)),
+ pully = dy / sqrt(n->fMxy.Dy2() + t.GetCovariance(1, 1));
// pullt = dt / sqrt(n->fMt.Dt2() + n->fTrack.GetCovariance(5, 5));
for (auto& projection : fProjection) {
@@ -387,7 +390,7 @@ bool CbmRecoQaTask::Detector::View::Load(const CbmKfTrackFitter::FitNode* n, con
switch (projection.first) {
case eViewProjection::kXYa: hh->Fill(n->fMxy.X(), n->fMxy.Y()); break;
- case eViewProjection::kXYp: hh->Fill(n->fTrack.X(), n->fTrack.Y()); break;
+ case eViewProjection::kXYp: hh->Fill(t.X(), t.Y()); break;
case eViewProjection::kXdX: hh->Fill(n->fMxy.X(), scale * dx); break;
case eViewProjection::kYdY: hh->Fill(n->fMxy.Y(), scale * dy); break;
case eViewProjection::kWdT: hh->Fill(n->fMxy.X(), scale * dt); break;
@@ -397,7 +400,7 @@ bool CbmRecoQaTask::Detector::View::Load(const CbmKfTrackFitter::FitNode* n, con
}
if (point) {
- double dxMC = point->GetX() - n->fTrack.X(), dyMC = point->GetY() - n->fTrack.Y();
+ double dxMC = point->GetX() - t.X(), dyMC = point->GetY() - t.Y();
switch (projection.first) {
case eViewProjection::kXdXMC: hh->Fill(point->GetX(), scale * dxMC); break;
@@ -499,7 +502,7 @@ void CbmRecoQaTask::Exec(Option_t*)
if (!FilterTrack(trackObj)) continue;
auto track = dynamic_cast<const CbmGlobalTrack*>(trackObj);
- CbmKfTrackFitter::Track trkKf;
+ CbmKfTrackFitter::Trajectory trkKf;
if (!fFitter.CreateGlobalTrack(trkKf, *track)) {
LOG(fatal) << GetName() << "::Exec: can not create the track for the fit!";
break;
@@ -521,7 +524,7 @@ void CbmRecoQaTask::Exec(Option_t*)
}
n.fIsTimeSet = n.fIsXySet = false;
- fFitter.FitTrack(trkKf);
+ fFitter.FitTrajectory(trkKf);
// match to MC
const CbmHit* hit = dynamic_cast<CbmHit*>(fHits[n.fHitSystemId]->At(n.fHitIndex));
@@ -552,14 +555,14 @@ void CbmRecoQaTask::Exec(Option_t*)
// Fit track with all hits ON for vx projection
int iprj(0);
for (auto& zprj : fProjs) {
- CbmKfTrackFitter::FitNode n = CbmKfTrackFitter::FitNode();
+ CbmKfTrackFitter::TrajectoryNode n = CbmKfTrackFitter::TrajectoryNode();
n.fZ = zprj.first;
n.fIsXySet = true;
n.fReference1 = iprj++;
trkKf.fNodes.push_back(n);
}
trkKf.OrderNodesInZ();
- fFitter.FitTrack(trkKf);
+ fFitter.FitTrajectory(trkKf);
for (auto& n : trkKf.fNodes) {
if (n.fReference1 < 0) continue;
if (fProjs.find(n.fReference1 + nnodes) == fProjs.end()) {
diff --git a/reco/qa/CbmRecoQaTask.h b/reco/qa/CbmRecoQaTask.h
index 197264f9d4b87cac65928d5df4ab15aefcab621f..f36c458da682e6125e4269081831c6ae655c1929 100644
--- a/reco/qa/CbmRecoQaTask.h
+++ b/reco/qa/CbmRecoQaTask.h
@@ -88,7 +88,7 @@ class CbmRecoQaTask : public FairTask {
bool HasAddress(const CbmHit* h, double& x, double& y, double& dx, double& dy) const;
template<class Hit>
bool Load(const CbmHit* h, const FairMCPoint* point, const CbmEvent* ev);
- bool Load(const CbmKfTrackFitter::FitNode* n, const FairMCPoint* point);
+ bool Load(const CbmKfTrackFitter::TrajectoryNode* n, const FairMCPoint* point);
std::string ToString() const;
static std::string ToString(eViewProjection prj);