diff --git a/reco/L1/CMakeLists.txt b/reco/L1/CMakeLists.txt
index 6b07858575b42d36875dc8402d4c656c2b27307f..df1bb9c9ab1bc339ca8ccba004b24c36109bc3e5 100644
--- a/reco/L1/CMakeLists.txt
+++ b/reco/L1/CMakeLists.txt
@@ -46,7 +46,7 @@ set(SRCS
CbmL1ReadEvent.cxx
CbmCaPerformance.cxx
L1Algo/L1Station.cxx
- L1Algo/L1TrackParFit.cxx
+ L1Algo/L1Fit.cxx
L1Algo/L1MCEvent.cxx
CbmL1MCTrack.cxx
L1Algo/L1Material.cxx
diff --git a/reco/L1/CbmL1Performance.cxx b/reco/L1/CbmL1Performance.cxx
index 07d4b2f03fd3ce5dd36f2cdbd01f7c4bff5b6fd5..f18fb3dc8ae6e1f42be47e0bcdcc49e355d98873 100644
--- a/reco/L1/CbmL1Performance.cxx
+++ b/reco/L1/CbmL1Performance.cxx
@@ -56,7 +56,7 @@
#include "L1Algo/L1Algo.h"
#include "L1Algo/L1Def.h"
-#include "L1Algo/L1TrackParFit.h"
+#include "L1Algo/L1Fit.h"
using std::cout;
using std::endl;
@@ -1101,7 +1101,7 @@ void CbmL1::TrackFitPerformance()
static bool first_call = 1;
- L1TrackParFit fit;
+ L1Fit fit;
fit.SetParticleMass(fpAlgo->GetDefaultParticleMass());
if (first_call) {
diff --git a/reco/L1/L1Algo/L1Algo.h b/reco/L1/L1Algo/L1Algo.h
index 3aa9afa64d3bc4aa9ecfd58da5af24ccbe5e0c1e..7d6a8b96828c69401b71855523d5ab2d4f43c103 100644
--- a/reco/L1/L1Algo/L1Algo.h
+++ b/reco/L1/L1Algo/L1Algo.h
@@ -31,6 +31,7 @@ class L1AlgoDraw;
#include "L1Branch.h"
#include "L1CloneMerger.h"
#include "L1Field.h"
+#include "L1Fit.h"
#include "L1Grid.h"
#include "L1Hit.h"
#include "L1HitPoint.h"
@@ -41,7 +42,6 @@ class L1AlgoDraw;
#include "L1Station.h"
#include "L1Track.h"
#include "L1TrackPar.h"
-#include "L1TrackParFit.h"
#include "L1Triplet.h"
#include "L1Utils.h"
#include "L1Vector.h"
@@ -340,12 +340,12 @@ public:
/// ------ Subroutines used by L1Algo::KFTrackFitter() ------
void GuessVec(L1TrackPar& t, fvec* xV, fvec* yV, fvec* zV, fvec* Sy, fvec* wV, int NHits, fvec* zCur = 0);
- void GuessVec(L1TrackParFit& t, fvec* xV, fvec* yV, fvec* zV, fvec* Sy, fvec* wV, int NHits, fvec* zCur = 0,
- fvec* timeV = 0, fvec* w_time = 0);
- void GuessVecNoField(L1TrackParFit& t, fvec& x_last, fvec& x_2last, fvec& y_last, fvec& y_2last, fvec& z_end,
- fvec& z_2last, fvec& time_last, fvec* w_time, fvec& dt2_last);
+ void GuessVec(L1Fit& t, fvec* xV, fvec* yV, fvec* zV, fvec* Sy, fvec* wV, int NHits, fvec* zCur = 0, fvec* timeV = 0,
+ fvec* w_time = 0);
+ void GuessVecNoField(L1Fit& t, fvec& x_last, fvec& x_2last, fvec& y_last, fvec& y_2last, fvec& z_end, fvec& z_2last,
+ fvec& time_last, fvec* w_time, fvec& dt2_last);
- void FilterFirst(L1TrackParFit& track, fvec& x, fvec& y, fvec& t, fvec& dt2, fvec& dx2, fvec& dy2, fvec& dxy);
+ void FilterFirst(L1Fit& track, fvec& x, fvec& y, fvec& t, fvec& dt2, fvec& dx2, fvec& dy2, fvec& dxy);
#ifdef DRAW
diff --git a/reco/L1/L1Algo/L1BranchExtender.cxx b/reco/L1/L1Algo/L1BranchExtender.cxx
index 7434fdebe4df1008632e1f87ff107ca98162cb9f..28306b65357208a5ea60c0050f0dc0a24ff60706 100644
--- a/reco/L1/L1Algo/L1BranchExtender.cxx
+++ b/reco/L1/L1Algo/L1BranchExtender.cxx
@@ -27,7 +27,7 @@ void L1Algo::BranchFitterFast(const L1Branch& t, L1TrackPar& Tout, const bool di
{
L1_assert(t.NHits >= 3);
- L1TrackParFit fit;
+ L1Fit fit;
fit.SetParticleMass(GetDefaultParticleMass());
fit.fTr = Tout;
L1TrackPar& T = fit.fTr;
@@ -163,7 +163,7 @@ void L1Algo::FindMoreHits(L1Branch& t, L1TrackPar& Tout, const bool dir,
L1Vector<L1HitIndex_t> newHits {"L1TrackExtender::newHits"};
newHits.reserve(fParameters.GetNstationsActive());
- L1TrackParFit fit;
+ L1Fit fit;
fit.SetParticleMass(GetDefaultParticleMass());
fit.fQp0 = qp0;
fit.fTr = Tout;
diff --git a/reco/L1/L1Algo/L1CATrackFinder.cxx b/reco/L1/L1Algo/L1CATrackFinder.cxx
index a64086dce0dffdcfa747ee1c74dc94e195dc29aa..9dd1c86d8aa1a99746d4dc62107d5c0ba70958a0 100644
--- a/reco/L1/L1Algo/L1CATrackFinder.cxx
+++ b/reco/L1/L1Algo/L1CATrackFinder.cxx
@@ -25,13 +25,13 @@
#include "L1Algo.h"
#include "L1Assert.h"
#include "L1Branch.h"
+#include "L1Fit.h"
#include "L1Grid.h"
#include "L1HitArea.h"
#include "L1HitPoint.h"
#include "L1Portion.h"
#include "L1Track.h"
#include "L1TrackPar.h"
-#include "L1TrackParFit.h"
#ifdef _OPENMP
#include "omp.h"
@@ -211,7 +211,7 @@ inline void L1Algo::findSingletsStep1( /// input 1st stage of singlet search
const L1Station& fld1Sta1 = fParameters.GetStation(fld1Ista1);
const L1Station& fld1Sta2 = fParameters.GetStation(fld1Ista2);
- L1TrackParFit fit;
+ L1Fit fit;
fit.SetParticleMass(GetDefaultParticleMass());
fit.fQp0 = fMaxInvMom;
@@ -325,7 +325,7 @@ inline void L1Algo::findDoubletsStep0(
n2 = 0; // number of doublets
- L1TrackParFit fit;
+ L1Fit fit;
fit.SetParticleMass(GetDefaultParticleMass());
for (Tindex i1 = 0; i1 < n1; ++i1) // for each singlet
@@ -427,13 +427,13 @@ inline void L1Algo::findDoubletsStep0(
fvec chi2 = T1.chi2;
- L1TrackParFit::FilterChi2XYC00C10C11(stam.frontInfo, x, y, C00, C10, C11, chi2, hitm.U(), hitm.dU2());
+ L1Fit::FilterChi2XYC00C10C11(stam.frontInfo, x, y, C00, C10, C11, chi2, hitm.U(), hitm.dU2());
if (!fpCurrentIteration->GetTrackFromTripletsFlag()) {
if (chi2[i1_4] > fDoubletChi2Cut) continue;
}
- L1TrackParFit::FilterChi2(stam.backInfo, x, y, C00, C10, C11, chi2, hitm.V(), hitm.dV2());
+ L1Fit::FilterChi2(stam.backInfo, x, y, C00, C10, C11, chi2, hitm.V(), hitm.dV2());
// FilterTime(T1, hitm.T(), hitm.dT2());
@@ -540,7 +540,7 @@ inline void L1Algo::findTripletsStep0( // input
bool isMomentumFitted = (fIsTargetField || (stal.fieldStatus != 0) || (stam.fieldStatus != 0));
//bool isTimeFitted = ((stal.timeInfo != 0) || (stam.timeInfo != 0));
- L1TrackParFit fit;
+ L1Fit fit;
fit.SetParticleMass(fDefaultMass);
fit.fQp0 = fvec(0.);
@@ -681,7 +681,7 @@ inline void L1Algo::findTripletsStep0( // input
auto [xr, yr] = star.ConvUVtoXY<fscal>(hitr.U(), hitr.V());
- L1TrackParFit fit3;
+ L1Fit fit3;
fit3.SetParticleMass(fDefaultMass);
fit3.fQp0 = T2.qp;
@@ -725,9 +725,9 @@ inline void L1Algo::findTripletsStep0( // input
fvec chi2 = T2.chi2;
- L1TrackParFit::FilterChi2XYC00C10C11(star.frontInfo, x, y, C00, C10, C11, chi2, hitr.U(), hitr.dU2());
+ L1Fit::FilterChi2XYC00C10C11(star.frontInfo, x, y, C00, C10, C11, chi2, hitr.U(), hitr.dU2());
- L1TrackParFit::FilterChi2(star.backInfo, x, y, C00, C10, C11, chi2, hitr.V(), hitr.dV2());
+ L1Fit::FilterChi2(star.backInfo, x, y, C00, C10, C11, chi2, hitr.V(), hitr.dV2());
fit3.FilterTime(hitr.T(), hitr.dT2(), fvec::One(), star.timeInfo);
@@ -807,7 +807,7 @@ inline void L1Algo::findTripletsStep1( // input
// L1TrackPar *T_3
L1Vector<L1TrackPar>& T_3)
{
- L1TrackParFit fit;
+ L1Fit fit;
fit.SetParticleMass(fDefaultMass);
for (Tindex i3_V = 0; i3_V < n3_V; ++i3_V) {
@@ -838,7 +838,7 @@ inline void L1Algo::findTripletsStep2(Tindex n3, int istal, int istam, int istar
ca::tools::Debugger::Instance().AddNtuple("triplets", "ev:mc:sta:p:vx:vy:vz:chi2:ndf");
- L1TrackParFit fit;
+ L1Fit fit;
if (fpCurrentIteration->GetElectronFlag()) { fit.SetParticleMass(L1Constants::phys::kElectronMass); }
else {
diff --git a/reco/L1/L1Algo/L1CloneMerger.cxx b/reco/L1/L1Algo/L1CloneMerger.cxx
index 78cead9bb5f463415d67687a9bdcb0872ed1a72b..ed665648f901da66414aad087215b7a62960a061 100644
--- a/reco/L1/L1Algo/L1CloneMerger.cxx
+++ b/reco/L1/L1Algo/L1CloneMerger.cxx
@@ -12,9 +12,9 @@
#include <iostream>
#include "L1Algo.h"
+#include "L1Fit.h"
#include "L1Parameters.h"
#include "L1Track.h"
-#include "L1TrackParFit.h"
// ---------------------------------------------------------------------------------------------------------------------
//
@@ -76,11 +76,11 @@ void L1CloneMerger::Exec(L1Vector<L1Track>& extTracks, L1Vector<L1HitIndex_t>& e
isDownstreamNeighbour[iTr] = false;
}
- L1TrackParFit fitB;
+ L1Fit fitB;
fitB.SetParticleMass(frAlgo.GetDefaultParticleMass());
fitB.fQp0 = fvec(0.);
- L1TrackParFit fitF;
+ L1Fit fitF;
fitF.SetParticleMass(frAlgo.GetDefaultParticleMass());
fitF.fQp0 = fvec(0.);
diff --git a/reco/L1/L1Algo/L1Filtration.h b/reco/L1/L1Algo/L1Filtration.h
deleted file mode 100644
index f22e7c9e4e146a8b04b36127ea9d75691ac480df..0000000000000000000000000000000000000000
--- a/reco/L1/L1Algo/L1Filtration.h
+++ /dev/null
@@ -1,355 +0,0 @@
-/* Copyright (C) 2007-2018 GSI Helmholtzzentrum fuer Schwerionenforschung, Darmstadt
- SPDX-License-Identifier: GPL-3.0-only
- Authors: Sergey Gorbunov [committer], Maksym Zyzak, Igor Kulakov, Valentina Akishina */
-
-#ifndef L1Filtration_h
-#define L1Filtration_h
-
-#include "L1Def.h"
-#include "L1TrackPar.h"
-#include "L1UMeasurementInfo.h"
-#include "L1XYMeasurementInfo.h"
-
-//#define cnst static const fvec
-#define cnst const fvec
-
-
-inline void FilterTime(L1TrackPar& T, fvec t, fvec dt2, fvec timeInfo = fvec::One(), fvec w = fvec::One())
-{
- // filter track with a time measurement
-
- // F = CH'
- fvec F0 = T.C50;
- fvec F1 = T.C51;
- fvec F2 = T.C52;
- fvec F3 = T.C53;
- fvec F4 = T.C54;
- fvec F5 = T.C55;
-
- fvec HCH = T.C55;
-
- w.setZero(timeInfo <= fvec::Zero());
-
- // when dt0 is much smaller than current time error,
- // set track time exactly to the measurement value without filtering
- // it helps to keep the initial time errors reasonably small
- // the calculations in the covariance matrix are not affected
-
- const fmask maskDoFilter = (HCH < dt2 * 16.f);
- //const fvec maskDoFilter = _f32vec4_true;
-
- fvec wi = w / (dt2 + 1.0000001f * HCH);
- fvec zeta = T.t - t;
- fvec zetawi = w * zeta / (iif(maskDoFilter, dt2, fvec::Zero()) + HCH);
-
- //T.chi2 += maskDoFilter & (zeta * zetawi);
- T.chi2 += zeta * zeta * wi;
- T.NDF += w;
-
- fvec K1 = F1 * wi;
- fvec K2 = F2 * wi;
- fvec K3 = F3 * wi;
- fvec K4 = F4 * wi;
- fvec K5 = F5 * wi;
-
- T.x -= F0 * zetawi;
- T.y -= F1 * zetawi;
- T.tx -= F2 * zetawi;
- T.ty -= F3 * zetawi;
- T.qp -= F4 * zetawi;
- T.t -= F5 * zetawi;
-
- T.C00 -= F0 * F0 * wi;
- T.C10 -= K1 * F0;
- T.C11 -= K1 * F1;
- T.C20 -= K2 * F0;
- T.C21 -= K2 * F1;
- T.C22 -= K2 * F2;
- T.C30 -= K3 * F0;
- T.C31 -= K3 * F1;
- T.C32 -= K3 * F2;
- T.C33 -= K3 * F3;
- T.C40 -= K4 * F0;
- T.C41 -= K4 * F1;
- T.C42 -= K4 * F2;
- T.C43 -= K4 * F3;
- T.C44 -= K4 * F4;
- T.C50 -= K5 * F0;
- T.C51 -= K5 * F1;
- T.C52 -= K5 * F2;
- T.C53 -= K5 * F3;
- T.C54 -= K5 * F4;
- T.C55 -= K5 * F5;
-}
-
-
-inline void L1Filter(L1TrackPar& T, const L1UMeasurementInfo& info, fvec u, fvec du2, fvec w)
-{
- // filter the track T with an 1-D measurement u
-
- fvec zeta = info.cos_phi * T.x + info.sin_phi * T.y - u;
-
- // F = CH'
- fvec F0 = info.cos_phi * T.C00 + info.sin_phi * T.C10;
- fvec F1 = info.cos_phi * T.C10 + info.sin_phi * T.C11;
- fvec F2 = info.cos_phi * T.C20 + info.sin_phi * T.C21;
- fvec F3 = info.cos_phi * T.C30 + info.sin_phi * T.C31;
- fvec F4 = info.cos_phi * T.C40 + info.sin_phi * T.C41;
- fvec F5 = info.cos_phi * T.C50 + info.sin_phi * T.C51;
-
- fvec HCH = (F0 * info.cos_phi + F1 * info.sin_phi);
-
- // when the measurement error du2 is much smaller than the current track error HCH,
- // move track exactly to the measurement without filtering
- // it helps to keep the initial track errors reasonably small
- // the calculations in the covariance matrix are not affected
-
- const fmask maskDoFilter = (HCH < du2 * 16.f);
- //const fvec maskDoFilter = _f32vec4_true;
-
- // correction to HCH is needed for the case when du2 is so small
- // with respect to HCH that it disappears due to the roundoff error
- //
- fvec wi = w / (du2 + 1.0000001f * HCH);
- fvec zetawi = w * zeta / (iif(maskDoFilter, du2, fvec::Zero()) + HCH);
-
- // T.chi2 += iif( maskDoFilter, zeta * zetawi, fvec::Zero() );
- T.chi2 += zeta * zeta * wi;
- T.NDF += w;
-
- fvec K1 = F1 * wi;
- fvec K2 = F2 * wi;
- fvec K3 = F3 * wi;
- fvec K4 = F4 * wi;
- fvec K5 = F5 * wi;
-
- T.x -= F0 * zetawi;
- T.y -= F1 * zetawi;
- T.tx -= F2 * zetawi;
- T.ty -= F3 * zetawi;
- T.qp -= F4 * zetawi;
- T.t -= F5 * zetawi;
-
- T.C00 -= F0 * F0 * wi;
- T.C10 -= K1 * F0;
- T.C11 -= K1 * F1;
- T.C20 -= K2 * F0;
- T.C21 -= K2 * F1;
- T.C22 -= K2 * F2;
- T.C30 -= K3 * F0;
- T.C31 -= K3 * F1;
- T.C32 -= K3 * F2;
- T.C33 -= K3 * F3;
- T.C40 -= K4 * F0;
- T.C41 -= K4 * F1;
- T.C42 -= K4 * F2;
- T.C43 -= K4 * F3;
- T.C44 -= K4 * F4;
- T.C50 -= K5 * F0;
- T.C51 -= K5 * F1;
- T.C52 -= K5 * F2;
- T.C53 -= K5 * F3;
- T.C54 -= K5 * F4;
- T.C55 -= K5 * F5;
-}
-
-
-inline void L1FilterChi2(const L1UMeasurementInfo& info, const fvec& x, const fvec& y, const fvec& C00, const fvec& C10,
- const fvec& C11, fvec& chi2, const fvec& u, const fvec& du2)
-{
- fvec zeta, HCH;
- fvec F0, F1;
-
- zeta = info.cos_phi * x + info.sin_phi * y - u;
-
- // F = CH'
- F0 = info.cos_phi * C00 + info.sin_phi * C10;
- F1 = info.cos_phi * C10 + info.sin_phi * C11;
-
- HCH = (F0 * info.cos_phi + F1 * info.sin_phi);
-
- chi2 += zeta * zeta / (du2 + HCH);
-}
-
-inline void L1FilterChi2XYC00C10C11(const L1UMeasurementInfo& info, fvec& x, fvec& y, fvec& C00, fvec& C10, fvec& C11,
- fvec& chi2, const fvec& u, const fvec& du2)
-{
- fvec wi, zeta, zetawi, HCH;
- fvec F0, F1;
- fvec K1;
-
- zeta = info.cos_phi * x + info.sin_phi * y - u;
-
- // F = CH'
- F0 = info.cos_phi * C00 + info.sin_phi * C10;
- F1 = info.cos_phi * C10 + info.sin_phi * C11;
-
- HCH = (F0 * info.cos_phi + F1 * info.sin_phi);
-
- wi = fvec(1.) / (du2 + HCH);
- zetawi = zeta * wi;
- chi2 += zeta * zetawi;
-
- K1 = F1 * wi;
-
- x -= F0 * zetawi;
- y -= F1 * zetawi;
-
- C00 -= F0 * F0 * wi;
- C10 -= K1 * F0;
- C11 -= K1 * F1;
-}
-
-inline void L1FilterVtx(L1TrackPar& T, fvec x, fvec y, const L1XYMeasurementInfo& info, fvec extrDx, fvec extrDy,
- fvec J[])
-{
- cnst TWO(2.);
-
- fvec zeta0, zeta1, S00, S10, S11, si;
- fvec F00, F10, F20, F30, F40, F01, F11, F21, F31, F41;
- fvec K00, K10, K20, K30, K40, K01, K11, K21, K31, K41;
-
- //zeta0 = T.x + J[0]*T.tx + J[1]*T.ty + J[2]*T.qp - x;
- //zeta1 = T.y + J[3]*T.tx + J[4]*T.ty + J[5]*T.qp - y;
-
- zeta0 = T.x + extrDx - x;
- zeta1 = T.y + extrDy - y;
-
- // H = 1 0 J[0] J[1] J[2]
- // 0 1 J[3] J[4] J[5]
-
- // F = CH'
- F00 = T.C00;
- F01 = T.C10;
- F10 = T.C10;
- F11 = T.C11;
- F20 = J[0] * T.C22;
- F21 = J[3] * T.C22;
- F30 = J[1] * T.C33;
- F31 = J[4] * T.C33;
- F40 = J[2] * T.C44;
- F41 = J[5] * T.C44;
-
- S00 = info.C00 + F00 + J[0] * F20 + J[1] * F30 + J[2] * F40;
- S10 = info.C10 + F10 + J[3] * F20 + J[4] * F30 + J[5] * F40;
- S11 = info.C11 + F11 + J[3] * F21 + J[4] * F31 + J[5] * F41;
-
- si = fvec(1.) / (S00 * S11 - S10 * S10);
- //si = fvec(rcp(fvec((S00*S11 - S10*S10)[0])));
- fvec S00tmp = S00;
- S00 = si * S11;
- S10 = -si * S10;
- S11 = si * S00tmp;
-
- T.chi2 += zeta0 * zeta0 * S00 + fvec(2.) * zeta0 * zeta1 * S10 + zeta1 * zeta1 * S11;
- T.NDF += TWO;
-
- K00 = F00 * S00 + F01 * S10;
- K01 = F00 * S10 + F01 * S11;
- K10 = F10 * S00 + F11 * S10;
- K11 = F10 * S10 + F11 * S11;
- K20 = F20 * S00 + F21 * S10;
- K21 = F20 * S10 + F21 * S11;
- K30 = F30 * S00 + F31 * S10;
- K31 = F30 * S10 + F31 * S11;
- K40 = F40 * S00 + F41 * S10;
- K41 = F40 * S10 + F41 * S11;
-
- T.x -= K00 * zeta0 + K01 * zeta1;
- T.y -= K10 * zeta0 + K11 * zeta1;
- T.tx -= K20 * zeta0 + K21 * zeta1;
- T.ty -= K30 * zeta0 + K31 * zeta1;
- T.qp -= K40 * zeta0 + K41 * zeta1;
-
- T.C00 -= (K00 * F00 + K01 * F01);
- T.C10 -= (K10 * F00 + K11 * F01);
- T.C11 -= (K10 * F10 + K11 * F11);
- T.C20 = -(K20 * F00 + K21 * F01);
- T.C21 = -(K20 * F10 + K21 * F11);
- T.C22 -= (K20 * F20 + K21 * F21);
- T.C30 = -(K30 * F00 + K31 * F01);
- T.C31 = -(K30 * F10 + K31 * F11);
- T.C32 = -(K30 * F20 + K31 * F21);
- T.C33 -= (K30 * F30 + K31 * F31);
- T.C40 = -(K40 * F00 + K41 * F01);
- T.C41 = -(K40 * F10 + K41 * F11);
- T.C42 = -(K40 * F20 + K41 * F21);
- T.C43 = -(K40 * F30 + K41 * F31);
- T.C44 -= (K40 * F40 + K41 * F41);
-}
-
-
-inline void L1FilterXY(L1TrackPar& T, fvec x, fvec y, const L1XYMeasurementInfo& info)
-{
- cnst TWO = 2.f;
-
- fvec zeta0, zeta1, S00, S10, S11, si;
- fvec F00, F10, F20, F30, F40, F01, F11, F21, F31, F41;
- fvec K00, K10, K20, K30, K40, K01, K11, K21, K31, K41;
-
- zeta0 = T.x - x;
- zeta1 = T.y - y;
-
- // F = CH'
- F00 = T.C00;
- F10 = T.C10;
- F20 = T.C20;
- F30 = T.C30;
- F40 = T.C40;
- F01 = T.C10;
- F11 = T.C11;
- F21 = T.C21;
- F31 = T.C31;
- F41 = T.C41;
-
- S00 = F00 + info.C00;
- S10 = F10 + info.C10;
- S11 = F11 + info.C11;
-
- si = 1.f / (S00 * S11 - S10 * S10);
- fvec S00tmp = S00;
- S00 = si * S11;
- S10 = -si * S10;
- S11 = si * S00tmp;
-
- T.chi2 += zeta0 * zeta0 * S00 + 2.f * zeta0 * zeta1 * S10 + zeta1 * zeta1 * S11;
- T.NDF += TWO;
-
- K00 = F00 * S00 + F01 * S10;
- K01 = F00 * S10 + F01 * S11;
- K10 = F10 * S00 + F11 * S10;
- K11 = F10 * S10 + F11 * S11;
- K20 = F20 * S00 + F21 * S10;
- K21 = F20 * S10 + F21 * S11;
- K30 = F30 * S00 + F31 * S10;
- K31 = F30 * S10 + F31 * S11;
- K40 = F40 * S00 + F41 * S10;
- K41 = F40 * S10 + F41 * S11;
-
- T.x -= K00 * zeta0 + K01 * zeta1;
- T.y -= K10 * zeta0 + K11 * zeta1;
- T.tx -= K20 * zeta0 + K21 * zeta1;
- T.ty -= K30 * zeta0 + K31 * zeta1;
- T.qp -= K40 * zeta0 + K41 * zeta1;
-
- T.C00 -= K00 * F00 + K01 * F01;
- T.C10 -= K10 * F00 + K11 * F01;
- T.C11 -= K10 * F10 + K11 * F11;
- T.C20 -= K20 * F00 + K21 * F01;
- T.C21 -= K20 * F10 + K21 * F11;
- T.C22 -= K20 * F20 + K21 * F21;
- T.C30 -= K30 * F00 + K31 * F01;
- T.C31 -= K30 * F10 + K31 * F11;
- T.C32 -= K30 * F20 + K31 * F21;
- T.C33 -= K30 * F30 + K31 * F31;
- T.C40 -= K40 * F00 + K41 * F01;
- T.C41 -= K40 * F10 + K41 * F11;
- T.C42 -= K40 * F20 + K41 * F21;
- T.C43 -= K40 * F30 + K41 * F31;
- T.C44 -= K40 * F40 + K41 * F41;
-}
-
-
-#undef cnst
-
-#endif
diff --git a/reco/L1/L1Algo/L1TrackParFit.cxx b/reco/L1/L1Algo/L1Fit.cxx
similarity index 93%
rename from reco/L1/L1Algo/L1TrackParFit.cxx
rename to reco/L1/L1Algo/L1Fit.cxx
index d595db3e3e35e3ac8eb0a6f41edec2a7bebdada2..cd79f57cb417a7df48d1d82c04426b17760fb550 100644
--- a/reco/L1/L1Algo/L1TrackParFit.cxx
+++ b/reco/L1/L1Algo/L1Fit.cxx
@@ -2,11 +2,11 @@
SPDX-License-Identifier: GPL-3.0-only
Authors: Maksym Zyzak [committer], Valentina Akishina */
-#include "L1TrackParFit.h"
+#include "L1Fit.h"
#define cnst const fvec
-void L1TrackParFit::Filter(const L1UMeasurementInfo& info, const fvec& u, const fvec& sigma2, const fvec& w)
+void L1Fit::Filter(const L1UMeasurementInfo& info, const fvec& u, const fvec& sigma2, const fvec& w)
{
fvec zeta, HCH;
fvec F0, F1, F2, F3, F4, F5;
@@ -75,7 +75,7 @@ void L1TrackParFit::Filter(const L1UMeasurementInfo& info, const fvec& u, const
}
-void L1TrackParFit::FilterTime(fvec t, fvec dt2, fvec w, fvec timeInfo)
+void L1Fit::FilterTime(fvec t, fvec dt2, fvec w, fvec timeInfo)
{
// filter track with a time measurement
@@ -142,7 +142,7 @@ void L1TrackParFit::FilterTime(fvec t, fvec dt2, fvec w, fvec timeInfo)
}
-void L1TrackParFit::FilterXY(const L1XYMeasurementInfo& info, fvec x, fvec y)
+void L1Fit::FilterXY(const L1XYMeasurementInfo& info, fvec x, fvec y)
{
cnst TWO(2.);
@@ -222,8 +222,8 @@ void L1TrackParFit::FilterXY(const L1XYMeasurementInfo& info, fvec x, fvec y)
fTr.C55 -= K50 * F50 + K51 * F51;
}
-void L1TrackParFit::FilterExtrapolatedXY(const fvec& x, const fvec& y, const L1XYMeasurementInfo& info,
- const fvec& extrX, const fvec& extrY, const fvec Jx[6], const fvec Jy[6])
+void L1Fit::FilterExtrapolatedXY(const fvec& x, const fvec& y, const L1XYMeasurementInfo& info, const fvec& extrX,
+ const fvec& extrY, const fvec Jx[6], const fvec Jy[6])
{
// add a 2-D measurenent (x,y) at some z, that differs from fTr.z
// extrX, extrY are extrapolated track parameters at z, Jx, Jy are derivatives of the extrapolation
@@ -303,8 +303,8 @@ void L1TrackParFit::FilterExtrapolatedXY(const fvec& x, const fvec& y, const L1X
}
-void L1TrackParFit::Extrapolate // extrapolates track parameters and returns jacobian for extrapolation of CovMatrix
- (fvec z_out, // extrapolate to this z position
+void L1Fit::Extrapolate // extrapolates track parameters and returns jacobian for extrapolation of CovMatrix
+ (fvec z_out, // extrapolate to this z position
const L1FieldRegion& F, const fvec& w)
{
// use Q/p linearisation at fQp0
@@ -317,9 +317,8 @@ void L1TrackParFit::Extrapolate // extrapolates track parameters and returns ja
}
}
-void
- L1TrackParFit::ExtrapolateStep // extrapolates track parameters and returns jacobian for extrapolation of CovMatrix
- (fvec z_out, // extrapolate to this z position
+void L1Fit::ExtrapolateStep // extrapolates track parameters and returns jacobian for extrapolation of CovMatrix
+ (fvec z_out, // extrapolate to this z position
const L1FieldRegion& F, const fvec& w)
{
// use Q/p linearisation at fQp0
@@ -657,9 +656,9 @@ void
}
}
-void L1TrackParFit::
- ExtrapolateStepAnalytic // extrapolates track parameters and returns jacobian for extrapolation of CovMatrix
- (fvec z_out, // extrapolate to this z position
+void
+ L1Fit::ExtrapolateStepAnalytic // extrapolates track parameters and returns jacobian for extrapolation of CovMatrix
+ (fvec z_out, // extrapolate to this z position
const L1FieldRegion& F, const fvec& w)
{
//
@@ -924,7 +923,7 @@ void L1TrackParFit::
//cout<<"Extrapolation ok"<<endl;
}
-void L1TrackParFit::ExtrapolateLine(fvec z_out, const L1FieldRegion& F, const fvec& w)
+void L1Fit::ExtrapolateLine(fvec z_out, const L1FieldRegion& F, const fvec& w)
{
// extrapolate the track assuming fQp0 == 0
// TODO: write special simplified procedure
@@ -935,7 +934,7 @@ void L1TrackParFit::ExtrapolateLine(fvec z_out, const L1FieldRegion& F, const fv
fQp0 = qp0;
}
-void L1TrackParFit::ExtrapolateLine(fvec z_out, const fvec& w)
+void L1Fit::ExtrapolateLine(fvec z_out, const fvec& w)
{
// extrapolate the track assuming fQp0 == 0
//
@@ -993,7 +992,7 @@ void L1TrackParFit::ExtrapolateLine(fvec z_out, const fvec& w)
}
-void L1TrackParFit::AddMsInMaterial(const fvec& radThick, fvec w)
+void L1Fit::AddMsInMaterial(const fvec& radThick, fvec w)
{
cnst ONE = 1.;
@@ -1018,7 +1017,7 @@ void L1TrackParFit::AddMsInMaterial(const fvec& radThick, fvec w)
fTr.C33 += w * (ONE + tyty) * a;
}
-void L1TrackParFit::AddMsInThickMaterial(fvec radThick, fvec w, fvec thickness, bool fDownstream)
+void L1Fit::AddMsInThickMaterial(fvec radThick, fvec w, fvec thickness, bool fDownstream)
{
cnst ONE = 1.;
@@ -1058,7 +1057,7 @@ void L1TrackParFit::AddMsInThickMaterial(fvec radThick, fvec w, fvec thickness,
}
-void L1TrackParFit::EnergyLossCorrection(const fvec& radThick, fvec upstreamDirection, fvec w)
+void L1Fit::EnergyLossCorrection(const fvec& radThick, fvec upstreamDirection, fvec w)
{
const fvec qp2cut(1. / (10. * 10.)); // 10 GeV cut
const fvec qp02 = max(fQp0 * fQp0, qp2cut);
@@ -1089,8 +1088,7 @@ void L1TrackParFit::EnergyLossCorrection(const fvec& radThick, fvec upstreamDire
}
template<int atomicZ>
-void L1TrackParFit::EnergyLossCorrection(float atomicA, float rho, float radLen, const fvec& radThick, fvec direction,
- fvec w)
+void L1Fit::EnergyLossCorrection(float atomicA, float rho, float radLen, const fvec& radThick, fvec direction, fvec w)
{
const fvec qp2cut(1. / (10. * 10.)); // 10 GeV cut
const fvec qp02 = max(fQp0 * fQp0, qp2cut);
@@ -1153,7 +1151,7 @@ void L1TrackParFit::EnergyLossCorrection(float atomicA, float rho, float radLen,
}
-void L1TrackParFit::EnergyLossCorrectionIron(const fvec& radThick, fvec direction, fvec w)
+void L1Fit::EnergyLossCorrectionIron(const fvec& radThick, fvec direction, fvec w)
{
constexpr int atomicZ = 26;
constexpr float atomicA = 55.845f;
@@ -1163,7 +1161,7 @@ void L1TrackParFit::EnergyLossCorrectionIron(const fvec& radThick, fvec directio
}
-void L1TrackParFit::EnergyLossCorrectionCarbon(const fvec& radThick, fvec direction, fvec w)
+void L1Fit::EnergyLossCorrectionCarbon(const fvec& radThick, fvec direction, fvec w)
{
constexpr int atomicZ = 6;
constexpr float atomicA = 12.011f;
@@ -1172,7 +1170,7 @@ void L1TrackParFit::EnergyLossCorrectionCarbon(const fvec& radThick, fvec direct
EnergyLossCorrection<atomicZ>(atomicA, rho, radLen, radThick, direction, w);
}
-void L1TrackParFit::EnergyLossCorrectionAl(const fvec& radThick, fvec direction, fvec w)
+void L1Fit::EnergyLossCorrectionAl(const fvec& radThick, fvec direction, fvec w)
{
constexpr int atomicZ = 13;
constexpr float atomicA = 26.981f;
@@ -1181,8 +1179,8 @@ void L1TrackParFit::EnergyLossCorrectionAl(const fvec& radThick, fvec direction,
EnergyLossCorrection<atomicZ>(atomicA, rho, radLen, radThick, direction, w);
}
-void L1TrackParFit::GetExtrapolatedXYline(const fvec& z, const L1FieldRegion& F, fvec& extrX, fvec& extrY, fvec Jx[6],
- fvec Jy[6]) const
+void L1Fit::GetExtrapolatedXYline(const fvec& z, const L1FieldRegion& F, fvec& extrX, fvec& extrY, fvec Jx[6],
+ fvec Jy[6]) const
{
// extrapolate track assuming it is straight (qp==0)
// return the extrapolated X, Y and the derivatives of the extrapolated X and Y
@@ -1229,8 +1227,8 @@ void L1TrackParFit::GetExtrapolatedXYline(const fvec& z, const L1FieldRegion& F,
}
-void L1TrackParFit::AddTargetToLine(const fvec& targX, const fvec& targY, const fvec& targZ,
- const L1XYMeasurementInfo& targXYInfo, const L1FieldRegion& F)
+void L1Fit::AddTargetToLine(const fvec& targX, const fvec& targY, const fvec& targZ,
+ const L1XYMeasurementInfo& targXYInfo, const L1FieldRegion& F)
{
// Add the target constraint to a straight line track
@@ -1239,7 +1237,7 @@ void L1TrackParFit::AddTargetToLine(const fvec& targX, const fvec& targY, const
FilterExtrapolatedXY(targX, targY, targXYInfo, eX, eY, Jx, Jy);
}
-fvec L1TrackParFit::ApproximateBetheBloch(const fvec& bg2)
+fvec L1Fit::ApproximateBetheBloch(const fvec& bg2)
{
//
// This is the parameterization of the Bethe-Bloch formula inspired by Geant.
@@ -1285,8 +1283,8 @@ fvec L1TrackParFit::ApproximateBetheBloch(const fvec& bg2)
* (0.5f * log(_2me * bg2 * maxT / (mI * mI)) - bg2 / (fvec(1.f) + bg2) - d2);
}
-fvec L1TrackParFit::ApproximateBetheBloch(const fvec& bg2, const fvec& kp0, const fvec& kp1, const fvec& kp2,
- const fvec& kp3, const fvec& kp4)
+fvec L1Fit::ApproximateBetheBloch(const fvec& bg2, const fvec& kp0, const fvec& kp1, const fvec& kp2, const fvec& kp3,
+ const fvec& kp4)
{
//
// This is the parameterization of the Bethe-Bloch formula inspired by Geant.
@@ -1333,8 +1331,8 @@ fvec L1TrackParFit::ApproximateBetheBloch(const fvec& bg2, const fvec& kp0, cons
}
-void L1TrackParFit::FilterChi2XYC00C10C11(const L1UMeasurementInfo& info, fvec& x, fvec& y, fvec& C00, fvec& C10,
- fvec& C11, fvec& chi2, const fvec& u, const fvec& du2)
+void L1Fit::FilterChi2XYC00C10C11(const L1UMeasurementInfo& info, fvec& x, fvec& y, fvec& C00, fvec& C10, fvec& C11,
+ fvec& chi2, const fvec& u, const fvec& du2)
{
fvec wi, zeta, zetawi, HCH;
fvec F0, F1;
@@ -1363,8 +1361,8 @@ void L1TrackParFit::FilterChi2XYC00C10C11(const L1UMeasurementInfo& info, fvec&
}
-void L1TrackParFit::FilterChi2(const L1UMeasurementInfo& info, const fvec& x, const fvec& y, const fvec& C00,
- const fvec& C10, const fvec& C11, fvec& chi2, const fvec& u, const fvec& du2)
+void L1Fit::FilterChi2(const L1UMeasurementInfo& info, const fvec& x, const fvec& y, const fvec& C00, const fvec& C10,
+ const fvec& C11, fvec& chi2, const fvec& u, const fvec& du2)
{
fvec zeta, HCH;
fvec F0, F1;
diff --git a/reco/L1/L1Algo/L1TrackParFit.h b/reco/L1/L1Algo/L1Fit.h
similarity index 86%
rename from reco/L1/L1Algo/L1TrackParFit.h
rename to reco/L1/L1Algo/L1Fit.h
index 49bcc41c4da063f6916c2309e68cdaea6330c9c0..c057cf94e1e748116962ea4846fcdd3b4a277ed3 100644
--- a/reco/L1/L1Algo/L1TrackParFit.h
+++ b/reco/L1/L1Algo/L1Fit.h
@@ -2,8 +2,8 @@
SPDX-License-Identifier: GPL-3.0-only
Authors: Maksym Zyzak [committer], Sergey Gorbunov */
-#ifndef L1TrackParFit_h
-#define L1TrackParFit_h
+#ifndef L1Fit_h
+#define L1Fit_h
#include "L1Def.h"
#include "L1Field.h"
@@ -11,15 +11,15 @@
#include "L1UMeasurementInfo.h"
#include "L1XYMeasurementInfo.h"
-class L1TrackParFit {
+class L1Fit {
public:
- L1TrackParFit() = default;
+ L1Fit() = default;
- L1TrackParFit(const L1TrackPar& t) { SetTrack(t); }
+ L1Fit(const L1TrackPar& t) { SetTrack(t); }
template<typename T>
- L1TrackParFit(const T* tr, const T* C)
+ L1Fit(const T* tr, const T* C)
{
SetTrack(tr, C);
}
@@ -37,7 +37,7 @@ public:
fQp0 = fTr.qp;
}
- void SetOneEntry(const int i0, const L1TrackParFit& T1, const int i1);
+ void SetOneEntry(const int i0, const L1Fit& T1, const int i1);
void Print(int i = -1);
@@ -121,30 +121,30 @@ public:
// =============================================================================================
-inline void L1TrackParFit::Print(int i) { fTr.Print(i); }
+inline void L1Fit::Print(int i) { fTr.Print(i); }
-inline void L1TrackParFit::SetOneEntry(const int i0, const L1TrackParFit& T1, const int i1)
+inline void L1Fit::SetOneEntry(const int i0, const L1Fit& T1, const int i1)
{
fTr.SetOneEntry(i0, T1.fTr, i1);
fQp0[i0] = T1.fQp0[i1];
}
-inline void L1TrackParFit::ExtrapolateXC00Line(fvec z_out, fvec& x, fvec& C00) const
+inline void L1Fit::ExtrapolateXC00Line(fvec z_out, fvec& x, fvec& C00) const
{
const fvec dz = (z_out - fTr.z);
x = fTr.x + fTr.tx * dz;
C00 = fTr.C00 + dz * (2 * fTr.C20 + dz * fTr.C22);
}
-inline void L1TrackParFit::ExtrapolateYC11Line(fvec z_out, fvec& y, fvec& C11) const
+inline void L1Fit::ExtrapolateYC11Line(fvec z_out, fvec& y, fvec& C11) const
{
const fvec dz = (z_out - fTr.z);
y = fTr.y + fTr.ty * dz;
C11 = fTr.C11 + dz * (2 * fTr.C31 + dz * fTr.C33);
}
-inline void L1TrackParFit::ExtrapolateC10Line(fvec z_out, fvec& C10) const
+inline void L1Fit::ExtrapolateC10Line(fvec z_out, fvec& C10) const
{
const fvec dz = (z_out - fTr.z);
C10 = fTr.C10 + dz * (fTr.C21 + fTr.C30 + dz * fTr.C32);
diff --git a/reco/L1/L1Algo/L1TrackFitter.cxx b/reco/L1/L1Algo/L1TrackFitter.cxx
index d9bed492b92e930d9ca1c406a2f74c3333d70ead..c3bf044a6cc22a255b39bac22e2f63335d8bf18d 100644
--- a/reco/L1/L1Algo/L1TrackFitter.cxx
+++ b/reco/L1/L1Algo/L1TrackFitter.cxx
@@ -6,8 +6,8 @@
#include <vector>
#include "L1Algo.h"
+#include "L1Fit.h"
#include "L1TrackPar.h"
-#include "L1TrackParFit.h"
using std::cout;
using std::endl;
@@ -41,7 +41,7 @@ void L1Algo::L1KFTrackFitter()
const int nStations = fParameters.GetNstationsActive();
int nTracks_SIMD = fvec::size();
- L1TrackParFit fit; // fitting parametr coresponding to current track
+ L1Fit fit; // fitting parametr coresponding to current track
L1TrackPar& tr = fit.fTr;
fit.SetParticleMass(GetDefaultParticleMass());
@@ -273,7 +273,7 @@ void L1Algo::L1KFTrackFitter()
// extrapolate to the PV region
- L1TrackParFit fitpv = fit;
+ L1Fit fitpv = fit;
{
L1UMeasurementInfo vtxInfoX;
vtxInfoX.cos_phi = 1.;
@@ -452,7 +452,7 @@ void L1Algo::L1KFTrackFitter()
} // iter
}
}
-void L1Algo::GuessVecNoField(L1TrackParFit& t, fvec& x_last, fvec& x_2last, fvec& y_last, fvec& y_2last, fvec& z_end,
+void L1Algo::GuessVecNoField(L1Fit& t, fvec& x_last, fvec& x_2last, fvec& y_last, fvec& y_2last, fvec& z_end,
fvec& z_2last, fvec& time_last, fvec* /*w_time*/, fvec& dt2_last)
{
fvec dzi = fvec(1.) / (z_end - z_2last);
@@ -543,8 +543,8 @@ void L1Algo::GuessVec(L1TrackPar& t, fvec* xV, fvec* yV, fvec* zV, fvec* Sy, fve
t.z = z0;
}
-void L1Algo::GuessVec(L1TrackParFit& t, fvec* xV, fvec* yV, fvec* zV, fvec* Sy, fvec* wV, int NHits, fvec* zCur,
- fvec* timeV, fvec* w_time)
+void L1Algo::GuessVec(L1Fit& t, fvec* xV, fvec* yV, fvec* zV, fvec* Sy, fvec* wV, int NHits, fvec* zCur, fvec* timeV,
+ fvec* w_time)
// gives nice initial approximation for x,y,tx,ty - almost same as KF fit. qp - is shifted by 4%, resid_ual - ~3.5% (KF fit resid_ual - 1%).
{
fvec A0, A1 = ZERO, A2 = ZERO, A3 = ZERO, A4 = ZERO, A5 = ZERO, a0, a1 = ZERO, a2 = ZERO, b0, b1 = ZERO, b2 = ZERO;
@@ -628,7 +628,7 @@ void L1Algo::GuessVec(L1TrackParFit& t, fvec* xV, fvec* yV, fvec* zV, fvec* Sy,
}
-void L1Algo::FilterFirst(L1TrackParFit& track, fvec& x, fvec& y, fvec& t, fvec& dt2, fvec& d_xx, fvec& d_yy, fvec& d_xy)
+void L1Algo::FilterFirst(L1Fit& track, fvec& x, fvec& y, fvec& t, fvec& dt2, fvec& d_xx, fvec& d_yy, fvec& d_xy)
{
track.fTr.C00 = d_xx;
track.fTr.C10 = d_xy;
diff --git a/reco/L1/ParticleFinder/CbmL1PFFitter.cxx b/reco/L1/ParticleFinder/CbmL1PFFitter.cxx
index 1ef9732af5d5d214684a9d2d17094ceb5db51c9a..ff5f029fe6253aa3c07fbd601028bd2293f66fb0 100644
--- a/reco/L1/ParticleFinder/CbmL1PFFitter.cxx
+++ b/reco/L1/ParticleFinder/CbmL1PFFitter.cxx
@@ -39,9 +39,9 @@
#include "L1Algo.h" // Also defines L1Parameters
#include "L1Def.h"
#include "L1Field.h"
+#include "L1Fit.h"
#include "L1Station.h"
#include "L1TrackPar.h"
-#include "L1TrackParFit.h"
//typedef L1Fit1 L1Fit;
@@ -97,7 +97,7 @@ inline void CbmL1PFFitter::PFFieldRegion::getL1FieldRegion(L1FieldRegion& fld, i
inline CbmL1PFFitter::PFFieldRegion::PFFieldRegion(const L1FieldRegion& fld, int i) { setFromL1FieldRegion(fld, i); }
-void FilterFirst(L1TrackParFit& fit, fvec& x, fvec& y, fvec& dxx, fvec& dxy, fvec& dyy)
+void FilterFirst(L1Fit& fit, fvec& x, fvec& y, fvec& dxx, fvec& dxy, fvec& dyy)
{
L1TrackPar& tr = fit.fTr;
tr.C00 = dxx;
@@ -139,7 +139,7 @@ void CbmL1PFFitter::Fit(vector<CbmStsTrack>& Tracks, vector<int>& pidHypo)
int iVec = 0, i = 0;
int nTracks_SIMD = fvec::size();
- L1TrackParFit fit;
+ L1Fit fit;
fit.SetParticleMass(CbmL1::Instance()->fpAlgo->GetDefaultParticleMass());
L1TrackPar& T = fit.fTr; // fitting parametr coresponding to current track
@@ -453,7 +453,7 @@ void CbmL1PFFitter::GetChiToVertex(vector<CbmStsTrack>& Tracks, vector<PFFieldRe
int nTracks_SIMD = fvec::size();
- L1TrackParFit fit;
+ L1Fit fit;
L1TrackPar& T = fit.fTr; // fitting parametr coresponding to current track
CbmStsTrack* tr[fvec::size()] {nullptr};