diff --git a/algo/ca/core/pars/CaField.cxx b/algo/ca/core/pars/CaField.cxx
index 22401d4ea9d261d262a806f517d93128f845548c..de8f60531808f477d064ab056f01da358c1bc4e1 100644
--- a/algo/ca/core/pars/CaField.cxx
+++ b/algo/ca/core/pars/CaField.cxx
@@ -16,10 +16,11 @@ using namespace cbm::algo::ca;
//
// FieldValue methods
//
-
-bool FieldRegion::fgForceUseOfOriginalField = false;
+template<typename DataT>
+bool FieldRegion<DataT>::fgForceUseOfOriginalField = false;
+template<typename DataT>
std::function<std::tuple<double, double, double>(double x, double y, double z)>
- FieldRegion::fgOdiginalField([](double, double, double) {
+ FieldRegion<DataT>::fgOdiginalField([](double, double, double) {
assert(false && "ca::FieldRegion: The original field is not set, use cbm::ca::tools::SetCbmField()");
return std::tuple(constants::Undef<double>, constants::Undef<double>, constants::Undef<double>);
});
@@ -27,7 +28,8 @@ std::function<std::tuple<double, double, double>(double x, double y, double z)>
//----------------------------------------------------------------------------------------------------------------------
//
-void FieldValue::CheckConsistency() const
+template<typename DataT>
+void FieldValue<DataT>::CheckConsistency() const
{
/* Check SIMD data vectors for consistent initialization */
utils::CheckSimdVectorEquality(x, "FieldValue::x");
@@ -39,20 +41,24 @@ void FieldValue::CheckConsistency() const
//----------------------------------------------------------------------------------------------------------------------
// TODO:
-std::string FieldValue::ToString(int indentLevel) const
+template<typename DataT>
+std::string FieldValue<DataT>::ToString(int indentLevel) const
{
std::stringstream aStream{};
constexpr char indentChar = '\t';
std::string indent(indentLevel, indentChar);
- aStream << indent << "Bx [kG]: " << std::setw(12) << std::setfill(' ') << x[0] << '\n';
- aStream << indent << "By [kG]: " << std::setw(12) << std::setfill(' ') << y[0] << '\n';
- aStream << indent << "Bz [kG]: " << std::setw(12) << std::setfill(' ') << z[0];
+ aStream << indent << "Bx [kG]: " << std::setw(12) << std::setfill(' ') << utils::simd::Cast<DataT, float>(x) << '\n';
+ aStream << indent << "By [kG]: " << std::setw(12) << std::setfill(' ') << utils::simd::Cast<DataT, float>(y) << '\n';
+ aStream << indent << "Bz [kG]: " << std::setw(12) << std::setfill(' ') << utils::simd::Cast<DataT, float>(z);
return aStream.str();
}
-//----------------------------------------------------------------------------------------------------------------------
-//
-std::ostream& operator<<(std::ostream& out, const FieldValue& B) { return out << B.x << " | " << B.y << " | " << B.z; };
+namespace cbm::algo::ca
+{
+ template class FieldValue<fvec>;
+ template class FieldValue<float>;
+ template class FieldValue<double>;
+} // namespace cbm::algo::ca
//
// FieldSlice methods
@@ -60,7 +66,8 @@ std::ostream& operator<<(std::ostream& out, const FieldValue& B) { return out <<
//----------------------------------------------------------------------------------------------------------------------
//
-FieldSlice::FieldSlice()
+template<typename DataT>
+FieldSlice<DataT>::FieldSlice()
{
for (int i = 0; i < constants::size::MaxNFieldApproxCoefficients; ++i) {
cx[i] = constants::Undef<fscal>;
@@ -72,7 +79,8 @@ FieldSlice::FieldSlice()
//----------------------------------------------------------------------------------------------------------------------
//
-void FieldSlice::CheckConsistency() const
+template<typename DataT>
+void FieldSlice<DataT>::CheckConsistency() const
{
/* Check SIMD data vectors for consistent initialization */
for (int i = 0; i < constants::size::MaxNFieldApproxCoefficients; ++i) {
@@ -85,31 +93,32 @@ void FieldSlice::CheckConsistency() const
//----------------------------------------------------------------------------------------------------------------------
// TODO: Should it be inline? (S.Zharko)
-FieldValue FieldSlice::GetFieldValue(const fvec& x, const fvec& y) const
+template<typename DataT>
+FieldValue<DataT> FieldSlice<DataT>::GetFieldValue(const DataT& x, const DataT& y) const
{
- fvec x2 = x * x;
- fvec y2 = y * y;
- fvec xy = x * y;
-
- fvec x3 = x2 * x;
- fvec y3 = y2 * y;
- fvec xy2 = x * y2;
- fvec x2y = x2 * y;
-
- fvec x4 = x3 * x;
- fvec y4 = y3 * y;
- fvec xy3 = x * y3;
- fvec x2y2 = x2 * y2;
- fvec x3y = x3 * y;
-
- fvec x5 = x4 * x;
- fvec y5 = y4 * y;
- fvec xy4 = x * y4;
- fvec x2y3 = x2 * y3;
- fvec x3y2 = x3 * y2;
- fvec x4y = x4 * y;
-
- FieldValue B;
+ DataT x2 = x * x;
+ DataT y2 = y * y;
+ DataT xy = x * y;
+
+ DataT x3 = x2 * x;
+ DataT y3 = y2 * y;
+ DataT xy2 = x * y2;
+ DataT x2y = x2 * y;
+
+ DataT x4 = x3 * x;
+ DataT y4 = y3 * y;
+ DataT xy3 = x * y3;
+ DataT x2y2 = x2 * y2;
+ DataT x3y = x3 * y;
+
+ DataT x5 = x4 * x;
+ DataT y5 = y4 * y;
+ DataT xy4 = x * y4;
+ DataT x2y3 = x2 * y3;
+ DataT x3y2 = x3 * y2;
+ DataT x4y = x4 * y;
+
+ FieldValue<DataT> B;
B.x = cx[0] + cx[1] * x + cx[2] * y + cx[3] * x2 + cx[4] * xy + cx[5] * y2 + cx[6] * x3 + cx[7] * x2y + cx[8] * xy2
+ cx[9] * y3 + cx[10] * x4 + cx[11] * x3y + cx[12] * x2y2 + cx[13] * xy3 + cx[14] * y4 + cx[15] * x5
@@ -125,15 +134,17 @@ FieldValue FieldSlice::GetFieldValue(const fvec& x, const fvec& y) const
return B;
}
-FieldValue FieldSlice::GetFieldValueForLine(const TrackParamV& t) const
+template<typename DataT>
+FieldValue<DataT> FieldSlice<DataT>::GetFieldValueForLine(const TrackParamBase<DataT>& t) const
{
- fvec dz = z - t.GetZ();
+ DataT dz = z - t.GetZ();
return GetFieldValue(t.GetX() + t.GetTx() * dz, t.GetY() + t.GetTy() * dz);
}
//----------------------------------------------------------------------------------------------------------------------
//
-std::string FieldSlice::ToString(int indentLevel) const
+template<typename DataT>
+std::string FieldSlice<DataT>::ToString(int indentLevel) const
{
std::stringstream aStream{};
// TODO: possibly it is better to place the indentChar into L1Parameters (S.Zharko)
@@ -143,20 +154,28 @@ std::string FieldSlice::ToString(int indentLevel) const
for (int i = 0; i < constants::size::MaxNFieldApproxCoefficients; ++i) {
aStream << '\n' << indent;
aStream << std::setw(3) << std::setfill(' ') << i << ' ';
- aStream << std::setw(12) << std::setfill(' ') << cx[i][0] << ' ';
- aStream << std::setw(12) << std::setfill(' ') << cy[i][0] << ' ';
- aStream << std::setw(12) << std::setfill(' ') << cz[i][0];
+ aStream << std::setw(12) << std::setfill(' ') << utils::simd::Cast<DataT, float>(cx[i], 0) << ' ';
+ aStream << std::setw(12) << std::setfill(' ') << utils::simd::Cast<DataT, float>(cy[i], 0) << ' ';
+ aStream << std::setw(12) << std::setfill(' ') << utils::simd::Cast<DataT, float>(cz[i], 0);
}
return aStream.str();
}
+namespace cbm::algo::ca
+{
+ template class FieldSlice<fvec>;
+ template class FieldSlice<float>;
+ template class FieldSlice<double>;
+} // namespace cbm::algo::ca
+
//
// FieldRegion methdos
//
//----------------------------------------------------------------------------------------------------------------------
//
-FieldRegion::FieldRegion(float reg[10]) noexcept
+template<typename DataT>
+FieldRegion<DataT>::FieldRegion(float reg[10]) noexcept
: cx0(reg[0])
, cx1(reg[1])
, cx2(reg[2])
@@ -172,7 +191,8 @@ FieldRegion::FieldRegion(float reg[10]) noexcept
//----------------------------------------------------------------------------------------------------------------------
//
-void FieldRegion::CheckConsistency() const
+template<typename DataT>
+void FieldRegion<DataT>::CheckConsistency() const
{
/* Check SIMD data vectors for consistent initialization */
utils::CheckSimdVectorEquality(cx0, "FieldRegion::cx0");
@@ -191,23 +211,24 @@ void FieldRegion::CheckConsistency() const
//----------------------------------------------------------------------------------------------------------------------
// TODO: Should it be inline? (S.Zharko)
-FieldValue FieldRegion::Get(const fvec& x, const fvec& y, const fvec& z) const
+template<typename DataT>
+FieldValue<DataT> FieldRegion<DataT>::Get(const DataT& x, const DataT& y, const DataT& z) const
{
- FieldValue B;
+ FieldValue<DataT> B;
if (fgForceUseOfOriginalField || fUseOriginalField) {
- for (size_t i = 0; i < fvec::size(); i++) {
- auto [bx, by, bz] = fgOdiginalField(x[i], y[i], z[i]);
- B.x[i] = bx;
- B.y[i] = by;
- B.z[i] = bz;
+ for (size_t i = 0; i < utils::simd::Size<DataT>(); i++) {
+ auto [bx, by, bz] =
+ fgOdiginalField(utils::simd::Cast<DataT, double>(x, i), utils::simd::Cast<DataT, double>(y, i),
+ utils::simd::Cast<DataT, double>(z, i));
+ B.SetSimdEntry(bx, by, bz, i);
}
}
else {
- fvec dz = z - z0;
- fvec dz2 = dz * dz;
- B.x = cx0 + cx1 * dz + cx2 * dz2;
- B.y = cy0 + cy1 * dz + cy2 * dz2;
- B.z = cz0 + cz1 * dz + cz2 * dz2;
+ DataT dz = z - z0;
+ DataT dz2 = dz * dz;
+ B.x = cx0 + cx1 * dz + cx2 * dz2;
+ B.y = cy0 + cy1 * dz + cy2 * dz2;
+ B.z = cz0 + cz1 * dz + cz2 * dz2;
}
return B;
}
@@ -215,24 +236,25 @@ FieldValue FieldRegion::Get(const fvec& x, const fvec& y, const fvec& z) const
//----------------------------------------------------------------------------------------------------------------------
// TODO: Should it be inline? (S.Zharko)
-void FieldRegion::Set(const FieldValue& b0, const fvec b0z, const FieldValue& b1, const fvec b1z, const FieldValue& b2,
- const fvec b2z)
+template<typename DataT>
+void FieldRegion<DataT>::Set(const FieldValue<DataT>& b0, const DataT b0z, const FieldValue<DataT>& b1, const DataT b1z,
+ const FieldValue<DataT>& b2, const DataT b2z)
{
- z0 = b0z;
- fvec dz1 = b1z - b0z;
- fvec dz2 = b2z - b0z;
- fvec det = fvec::One() / (dz1 * dz2 * (b2z - b1z));
-
- fvec w21 = -dz2 * det;
- fvec w22 = dz1 * det;
- fvec w11 = -dz2 * w21;
- fvec w12 = -dz1 * w22;
-
- fvec db1 = b1.x - b0.x;
- fvec db2 = b2.x - b0.x;
- cx0 = b0.x;
- cx1 = db1 * w11 + db2 * w12;
- cx2 = db1 * w21 + db2 * w22;
+ z0 = b0z;
+ DataT dz1 = b1z - b0z;
+ DataT dz2 = b2z - b0z;
+ DataT det = utils::simd::One<DataT>() / (dz1 * dz2 * (b2z - b1z));
+
+ DataT w21 = -dz2 * det;
+ DataT w22 = dz1 * det;
+ DataT w11 = -dz2 * w21;
+ DataT w12 = -dz1 * w22;
+
+ DataT db1 = b1.x - b0.x;
+ DataT db2 = b2.x - b0.x;
+ cx0 = b0.x;
+ cx1 = db1 * w11 + db2 * w12;
+ cx2 = db1 * w21 + db2 * w22;
db1 = b1.y - b0.y;
db2 = b2.y - b0.y;
@@ -249,30 +271,32 @@ void FieldRegion::Set(const FieldValue& b0, const fvec b0z, const FieldValue& b1
//----------------------------------------------------------------------------------------------------------------------
// TODO: Should it be inline? (S.Zharko)
-void FieldRegion::Set(const FieldValue& b0, const fvec b0z, const FieldValue& b1, const fvec b1z)
+template<typename DataT>
+void FieldRegion<DataT>::Set(const FieldValue<DataT>& b0, const DataT b0z, const FieldValue<DataT>& b1, const DataT b1z)
{
- z0 = b0z;
- fvec dzi = fvec::One() / (b1z - b0z);
- cx0 = b0.x;
- cy0 = b0.y;
- cz0 = b0.z;
- cx1 = (b1.x - b0.x) * dzi;
- cy1 = (b1.y - b0.y) * dzi;
- cz1 = (b1.z - b0.z) * dzi;
- cx2 = 0.f;
- cy2 = 0.f;
- cz2 = 0.f;
+ z0 = b0z;
+ DataT dzi = utils::simd::One<DataT>() / (b1z - b0z);
+ cx0 = b0.x;
+ cy0 = b0.y;
+ cz0 = b0.z;
+ cx1 = (b1.x - b0.x) * dzi;
+ cy1 = (b1.y - b0.y) * dzi;
+ cz1 = (b1.z - b0.z) * dzi;
+ cx2 = 0.f;
+ cy2 = 0.f;
+ cz2 = 0.f;
}
//----------------------------------------------------------------------------------------------------------------------
// TODO: Should it be inline? (S.Zharko)
-void FieldRegion::Shift(fvec z)
+template<typename DataT>
+void FieldRegion<DataT>::Shift(DataT z)
{
- fvec dz = z - z0;
- fvec cx2dz = cx2 * dz;
- fvec cy2dz = cy2 * dz;
- fvec cz2dz = cz2 * dz;
- z0 = z;
+ DataT dz = z - z0;
+ DataT cx2dz = cx2 * dz;
+ DataT cy2dz = cy2 * dz;
+ DataT cz2dz = cz2 * dz;
+ z0 = z;
cx0 += (cx1 + cx2dz) * dz;
cy0 += (cy1 + cy2dz) * dz;
cz0 += (cz1 + cz2dz) * dz;
@@ -281,41 +305,10 @@ void FieldRegion::Shift(fvec z)
cz1 += cz2dz + cz2dz;
}
-//----------------------------------------------------------------------------------------------------------------------
-// TODO: Should it be inline? (S.Zharko)
-void FieldRegion::SetOneEntry(const int i0, const FieldRegion& f1, const int i1)
-{
- cx0[i0] = f1.cx0[i1];
- cx1[i0] = f1.cx1[i1];
- cx2[i0] = f1.cx2[i1];
- cy0[i0] = f1.cy0[i1];
- cy1[i0] = f1.cy1[i1];
- cy2[i0] = f1.cy2[i1];
- cz0[i0] = f1.cz0[i1];
- cz1[i0] = f1.cz1[i1];
- cz2[i0] = f1.cz2[i1];
- z0[i0] = f1.z0[i1];
-}
-
-//----------------------------------------------------------------------------------------------------------------------
-// TODO: Should it be inline? (S.Zharko)
-void FieldRegion::SetOneEntry(const FieldRegion& f1, const int i1)
-{
- cx0 = f1.cx0[i1];
- cx1 = f1.cx1[i1];
- cx2 = f1.cx2[i1];
- cy0 = f1.cy0[i1];
- cy1 = f1.cy1[i1];
- cy2 = f1.cy2[i1];
- cz0 = f1.cz0[i1];
- cz1 = f1.cz1[i1];
- cz2 = f1.cz2[i1];
- z0 = f1.z0[i1];
-}
-
//----------------------------------------------------------------------------------------------------------------------
//
-std::string FieldRegion::ToString(int indentLevel) const
+template<typename DataT>
+std::string FieldRegion<DataT>::ToString(int indentLevel) const
{
std::stringstream aStream{};
// TODO: possibly it is better to place the indentChar into L1Parameters (S.Zharko)
@@ -336,3 +329,10 @@ std::string FieldRegion::ToString(int indentLevel) const
aStream << indent << "z0 [cm]: " << z0;
return aStream.str();
}
+
+namespace cbm::algo::ca
+{
+ template class FieldRegion<fvec>;
+ template class FieldRegion<float>;
+ template class FieldRegion<double>;
+} // namespace cbm::algo::ca
diff --git a/algo/ca/core/pars/CaField.h b/algo/ca/core/pars/CaField.h
index e744448f5458978364f45c91dd1a88c55bb7435f..f7fbe1e0e58000cda0d494d9ce4d9a88ec795727 100644
--- a/algo/ca/core/pars/CaField.h
+++ b/algo/ca/core/pars/CaField.h
@@ -17,30 +17,60 @@
namespace cbm::algo::ca
{
+ template<typename DataT>
class FieldValue {
public:
- fvec x{0.f}; //< x-component of the field
- fvec y{0.f}; //< y-component of the field
- fvec z{0.f}; //< z-component of the field
+ DataT x{0.f}; //< x-component of the field
+ DataT y{0.f}; //< y-component of the field
+ DataT z{0.f}; //< z-component of the field
+
+ FieldValue() = default;
+
+ /// \brief Copy constructor with type conversion
+ template<typename DataI>
+ FieldValue(const FieldValue<DataI>& other)
+ : x(other.template GetX<DataT>())
+ , y(other.template GetY<DataT>())
+ , z(other.template GetZ<DataT>())
+ {
+ }
/// Combines the magnetic field with another field value object using weight
/// \param B other field value to combine with
/// \param w weight from 0 to 1
- void Combine(const FieldValue& B, const fvec& w);
+ void Combine(const FieldValue<DataT>& B, const DataT w);
- void Combine(const FieldValue& B, const fmask& w);
+ void Combine(const FieldValue<DataT>& B, const bool w);
+
+ void Combine(const FieldValue<fvec>& B, const fmask w);
/// Consistency checker
void CheckConsistency() const;
/// Operator << overloading
- friend std::ostream& operator<<(std::ostream& out, const FieldValue& B);
+ friend std::ostream& operator<<(std::ostream& out, const FieldValue<DataT>& B)
+ {
+ return out << B.x << " | " << B.y << " | " << B.z;
+ };
/// String representation of class contents
/// \param indentLevel number of indent characters in the output
std::string ToString(int indentLevel) const;
+ /// \brief Sets the entries of the field components into the i-th elements of the SIMD vectors
+ /// \details If DataT is not a SIMD type, the function simply overwrites values ignoring the index i
+ /// \param bx The value to set for the x-component of the field
+ /// \param by The value to set for the y-component of the field
+ /// \param bz The value to set for the z-component of the field
+ /// \param i The index at which the SIMD entries are set (ignored if DataT is not a SIMD type)
+ [[gnu::always_inline]] inline void SetSimdEntry(double bx, double by, double bz, size_t i)
+ {
+ utils::simd::SetEntry(x, bx, i); //B.x[i] = bx;
+ utils::simd::SetEntry(y, by, i); //B.y[i] = by;
+ utils::simd::SetEntry(z, bz, i); //B.z[i] = bz;
+ }
+
/// Serialization function
friend class boost::serialization::access;
template<class Archive>
@@ -52,23 +82,67 @@ namespace cbm::algo::ca
}
} _fvecalignment;
- [[gnu::always_inline]] inline void FieldValue::Combine(const FieldValue& B, const fvec& w)
+ template<typename DataT>
+ [[gnu::always_inline]] inline void FieldValue<DataT>::Combine(const FieldValue<DataT>& B, const DataT w)
{
x += w * (B.x - x);
y += w * (B.y - y);
z += w * (B.z - z);
}
- [[gnu::always_inline]] inline void FieldValue::Combine(const FieldValue& B, const fmask& w)
+ template<typename DataT>
+ [[gnu::always_inline]] inline void FieldValue<DataT>::Combine(const FieldValue<DataT>& B, const bool w)
+ {
+ if (w) {
+ x = B.x;
+ y = B.y;
+ z = B.z;
+ }
+ }
+
+ // TODO: Add specific masking implementation for vir (G.Kozlov)
+ template<>
+ [[gnu::always_inline]] inline void FieldValue<fvec>::Combine(const FieldValue<fvec>& B, const fmask w)
{
x(w) = B.x;
y(w) = B.y;
z(w) = B.z;
}
+ /// \brief Inherited classes for type names renaming without altering functionality
+ // TODO: This is a temporary solution and should be removed in the future (G.Kozlov)
+ class FieldValueV : public FieldValue<fvec> {
+ using FieldValue<fvec>::FieldValue; /// inherit all the constructors
+
+ public:
+ template<typename DataI>
+ FieldValueV(const FieldValue<DataI>& other) : FieldValue<fvec>(other)
+ {
+ }
+ };
+ class FieldValueF : public FieldValue<float> {
+ using FieldValue<float>::FieldValue; /// inherit all the constructors
+
+ public:
+ template<typename DataI>
+ FieldValueF(const FieldValue<DataI>& other) : FieldValue<float>(other)
+ {
+ }
+ };
+ class FieldValueD : public FieldValue<double> {
+ using FieldValue<double>::FieldValue; /// inherit all the constructors
+
+ public:
+ template<typename DataI>
+ FieldValueD(const FieldValue<DataI>& other) : FieldValue<double>(other)
+ {
+ }
+ };
+
/// Class represents a set of magnetic field approximation coefficients
///
// TODO: Crosscheck the default content (S.Zharko)
+ template<typename DataT>
class FieldSlice {
public:
/// Default constructor
@@ -81,12 +155,12 @@ namespace cbm::algo::ca
/// \param x x-coordinate of input
/// \param y y-coordinate of input
/// \return B the FieldValue output
- FieldValue GetFieldValue(const fvec& x, const fvec& y) const;
+ FieldValue<DataT> GetFieldValue(const DataT& x, const DataT& y) const;
/// Gets field value for the intersection with a straight track
/// \param t straight track
/// \return B the FieldValue output
- FieldValue GetFieldValueForLine(const TrackParamV& t) const;
+ FieldValue<DataT> GetFieldValueForLine(const TrackParamBase<DataT>& t) const;
/// String representation of class contents
/// \param indentLevel number of indent characters in the output
@@ -97,11 +171,11 @@ namespace cbm::algo::ca
// if the underlying type (fvec) has a default constructor, but
// we are sure, that it can be initialized with a float. (S.Zharko)
static constexpr auto kMaxNFieldApproxCoefficients = constants::size::MaxNFieldApproxCoefficients;
- fvec cx[kMaxNFieldApproxCoefficients]; ///< Polynomial coefficients for x-component of the field value
- fvec cy[kMaxNFieldApproxCoefficients]; ///< Polynomial coefficients for y-component of the field value
- fvec cz[kMaxNFieldApproxCoefficients]; ///< Polynomial coefficients for z-component of the field value
+ DataT cx[kMaxNFieldApproxCoefficients]; ///< Polynomial coefficients for x-component of the field value
+ DataT cy[kMaxNFieldApproxCoefficients]; ///< Polynomial coefficients for y-component of the field value
+ DataT cz[kMaxNFieldApproxCoefficients]; ///< Polynomial coefficients for z-component of the field value
- fvec z{constants::Undef<fscal>}; ///< z coordinate of the slice
+ DataT z{constants::Undef<float>}; ///< z coordinate of the slice
/// Serialization function
friend class boost::serialization::access;
@@ -115,7 +189,38 @@ namespace cbm::algo::ca
}
} _fvecalignment;
+ /// \brief Inherited classes for type names renaming without altering functionality
+ // TODO: This is a temporary solution and should be removed in the future (G.Kozlov)
+ class FieldSliceV : public FieldSlice<fvec> {
+ using FieldSlice<fvec>::FieldSlice; /// inherit all the constructors
+
+ public:
+ template<typename DataI>
+ FieldSliceV(const FieldSlice<DataI>& other) : FieldSlice<fvec>(other)
+ {
+ }
+ };
+ class FieldSliceF : public FieldSlice<float> {
+ using FieldSlice<float>::FieldSlice; /// inherit all the constructors
+
+ public:
+ template<typename DataI>
+ FieldSliceF(const FieldSlice<DataI>& other) : FieldSlice<float>(other)
+ {
+ }
+ };
+ class FieldSliceD : public FieldSlice<double> {
+ using FieldSlice<double>::FieldSlice; /// inherit all the constructors
+ public:
+ template<typename DataI>
+ FieldSliceD(const FieldSlice<DataI>& other) : FieldSlice<double>(other)
+ {
+ }
+ };
+
+
+ template<typename DataT>
class FieldRegion {
public:
// NOTE: When a custom constructor is defined, default constructor also should be provided (S.Zharko)
@@ -131,7 +236,7 @@ namespace cbm::algo::ca
/// \param y y-coordinate of the point to calculate the field
/// \param z z-coordinate of the point to calculate the field
/// \return the magnetic field value
- FieldValue Get(const fvec& x, const fvec& y, const fvec& z) const;
+ FieldValue<DataT> Get(const DataT& x, const DataT& y, const DataT& z) const;
/// Interpolates the magnetic field by three nodal points and sets the result to this FieldRegion object
/// The result is a quadratic interpolation of the field as a function of z
@@ -142,8 +247,8 @@ namespace cbm::algo::ca
/// \param b2 field value in the third nodal point
/// \param b2z z-coordinate of the third nodal point
/// TODO: does the sequence of b0z, b1z and b2z matter? If yes, probalby we need an assertion (S.Zharko)
- void Set(const FieldValue& b0, const fvec b0z, const FieldValue& b1, const fvec b1z, const FieldValue& b2,
- const fvec b2z);
+ void Set(const FieldValue<DataT>& b0, const DataT b0z, const FieldValue<DataT>& b1, const DataT b1z,
+ const FieldValue<DataT>& b2, const DataT b2z);
/// Interpolates the magnetic field by thwo nodal points and sets the result to this FieldRegion object.
/// The result is a linear interpolation of the field as a function of z
@@ -152,24 +257,31 @@ namespace cbm::algo::ca
/// \param b1 field value in the second nodal point
/// \param b1z z-coordinate of the second nodal point
/// TODO: does the sequence of b0z and b1z matter? If yes, probalby we need an assertion (S.Zharko)
- void Set(const FieldValue& b0, const fvec b0z, const FieldValue& b1, const fvec b1z);
+ void Set(const FieldValue<DataT>& b0, const DataT b0z, const FieldValue<DataT>& b1, const DataT b1z);
/// Shifts the coefficients to new central point
/// \param z z-coordinate of the new central point
- void Shift(fvec z);
-
- /// Replaces the selected layer of the coefficients with one from another
- /// FieldRegion object
- /// \param i0 the index of the fvect layer in this FieldRegion object
- /// \param fl the other FieldRegion object, which layer is going to be replaced
- /// \param i1 the index of the source fvect layer to copy
- void SetOneEntry(const int i0, const FieldRegion& f1, const int i1);
+ void Shift(DataT z);
+
+ /// Set one SIMD entry from one SIMD entry of the other class
+ /// It only works when DataT is fvec, TdataB is fvec
+ /// \param fb the other FieldRegion object, which layer is going to be replaced
+ /// \param ib the index of the source fvec layer to copy
+ template<typename TdataB>
+ void SetSimdEntry(const FieldRegion<TdataB>& fb, const int ib)
+ {
+ CopyBase<TdataB, true, false>(0, fb, ib);
+ }
- /// Replaces all the layers of the coefficients with one selected layer from another
- /// FieldRegion object
- /// \param fl the other FieldRegion object, which layer is going to be replaced
- /// \param i1 the index of the source fvect layer to copy
- void SetOneEntry(const FieldRegion& f1, const int i1);
+ /// Set one SIMD entry from one SIMD entry of the other class
+ /// It only works when DataT is fvec, TdataB is fvec
+ /// \param ia the index of the fvec layer in this FieldRegion object
+ /// \param fb the other FieldRegion object, which layer is going to be replaced
+ /// \param ib the index of the source fvec layer to copy
+ void SetSimdEntry(const int ia, const FieldRegion<fvec>& fb, const int ib)
+ {
+ CopyBase<fvec, false, false>(ia, fb, ib);
+ }
/// String representation of class contents
/// \param indentLevel number of indent characters in the output
@@ -187,27 +299,37 @@ namespace cbm::algo::ca
fgOdiginalField = f;
}
+ /// \brief Copies all/one entries from the other class
+ /// \tparam TdataB Type of the other class
+ /// \tparam TDoAllA If true, all entries of the current class must be set
+ /// \tparam TDoAllB If true, all entries of the other class must be used
+ /// \param ia Index of SIMD vector element of the current class
+ /// \param Tb Other class
+ /// \param ib Index of SIMD vector element of the other class
+ template<typename TdataB, bool TDoAllA, bool TDoAllB>
+ void CopyBase(const int ia, const FieldRegion<TdataB>& Tb, const int ib);
+
public:
static bool fgForceUseOfOriginalField;
static std::function<std::tuple<double, double, double>(double x, double y, double z)> fgOdiginalField;
// TODO: Probably it's better to have arrays instead of separate fvec values? (S.Zharko)
// Bx(z) = cx0 + cx1*(z-z0) + cx2*(z-z0)^2
- fvec cx0{0.f};
- fvec cx1{0.f};
- fvec cx2{0.f};
+ DataT cx0{0.f};
+ DataT cx1{0.f};
+ DataT cx2{0.f};
// By(z) = cy0 + cy1*(z-z0) + cy2*(z-z0)^2
- fvec cy0{0.f};
- fvec cy1{0.f};
- fvec cy2{0.f};
+ DataT cy0{0.f};
+ DataT cy1{0.f};
+ DataT cy2{0.f};
// Bz(z) = cz0 + cz1*(z-z0) + cz2*(z-z0)^2
- fvec cz0{0.f};
- fvec cz1{0.f};
- fvec cz2{0.f};
+ DataT cz0{0.f};
+ DataT cz1{0.f};
+ DataT cz2{0.f};
- fvec z0{0.f}; ///< z-coordinate of the field region central point
+ DataT z0{0.f}; ///< z-coordinate of the field region central point
bool fUseOriginalField{false};
@@ -228,4 +350,64 @@ namespace cbm::algo::ca
}
} _fvecalignment;
+ // ---------------------------------------------------------------------------------------------------------------------
+ //
+ template<typename TdataA>
+ template<typename TdataB, bool TDoAllA, bool TDoAllB>
+ inline void FieldRegion<TdataA>::CopyBase(const int ia, const FieldRegion<TdataB>& Tb, const int ib)
+ {
+ auto copy = [&](TdataA& a, const TdataB& b) {
+ utils::VecCopy<TdataA, TdataB, TDoAllA, TDoAllB>::CopyEntries(a, ia, b, ib);
+ };
+
+ copy(cx0, Tb.cx0);
+ copy(cx1, Tb.cx1);
+ copy(cx2, Tb.cx2);
+
+ copy(cy0, Tb.cy0);
+ copy(cy1, Tb.cy1);
+ copy(cy2, Tb.cy2);
+
+ copy(cz0, Tb.cz0);
+ copy(cz1, Tb.cz1);
+ copy(cz2, Tb.cz2);
+
+ copy(z0, Tb.z0);
+ } // CopyBase
+
+ /// \brief Inherited classes for type names renaming without altering functionality
+ // TODO: This is a temporary solution and should be removed in the future (G.Kozlov)
+ class FieldRegionV : public FieldRegion<fvec> {
+ using FieldRegion<fvec>::FieldRegion; /// inherit all the constructors
+
+ public:
+ template<typename DataI>
+ FieldRegionV(const FieldRegion<DataI>& other) : FieldRegion<fvec>(other)
+ {
+ }
+ };
+ class FieldRegionF : public FieldRegion<float> {
+ using FieldRegion<float>::FieldRegion; /// inherit all the constructors
+
+ public:
+ template<typename DataI>
+ FieldRegionF(const FieldRegion<DataI>& other) : FieldRegion<float>(other)
+ {
+ }
+ };
+ class FieldRegionD : public FieldRegion<double> {
+ using FieldRegion<double>::FieldRegion; /// inherit all the constructors
+
+ public:
+ template<typename DataI>
+ FieldRegionD(const FieldRegion<DataI>& other) : FieldRegion<double>(other)
+ {
+ }
+ };
+
+ /// \brief Explicit instantiation declarations for FieldRegion class with specific template types to ensure proper fgOdiginalField definition
+ extern template class FieldRegion<fvec>;
+ extern template class FieldRegion<float>;
+ extern template class FieldRegion<double>;
+
} // namespace cbm::algo::ca
diff --git a/algo/ca/core/pars/CaInitManager.cxx b/algo/ca/core/pars/CaInitManager.cxx
index 9e573d36cc5855ead203bf6bf27abb77a38294bd..39ea1c7ea9133b742964b610b4c2daea68eafcb7 100644
--- a/algo/ca/core/pars/CaInitManager.cxx
+++ b/algo/ca/core/pars/CaInitManager.cxx
@@ -66,8 +66,8 @@ void InitManager::ClearSetupInfo()
fInitController.SetFlag(EInitKey::kActiveDetectorIDs, false);
// Clear field info
- fParameters.fVertexFieldRegion = ca::FieldRegion();
- fParameters.fVertexFieldValue = ca::FieldValue();
+ fParameters.fVertexFieldRegion = ca::FieldRegionV();
+ fParameters.fVertexFieldValue = ca::FieldValueV();
fInitController.SetFlag(EInitKey::kPrimaryVertexField, false);
// Clear target position
@@ -282,7 +282,7 @@ void InitManager::InitTargetField(double zStep)
constexpr int nPointsNodal{3};
std::array<double, nPointsNodal> inputNodalZ{fTargetZ, fTargetZ + zStep, fTargetZ + 2. * zStep};
- std::array<ca::FieldValue, nPointsNodal> B{};
+ std::array<ca::FieldValueV, nPointsNodal> B{};
std::array<fvec, nPointsNodal> z{};
// loop over nodal points
for (int idx = 0; idx < nPointsNodal; ++idx) {
diff --git a/algo/ca/core/pars/CaParameters.h b/algo/ca/core/pars/CaParameters.h
index 988941163c6864640a441336c1598704123fff6c..b1adb5beac4bc5c915112838968bea64dc25ea8b 100644
--- a/algo/ca/core/pars/CaParameters.h
+++ b/algo/ca/core/pars/CaParameters.h
@@ -180,10 +180,10 @@ namespace cbm::algo::ca
fvec GetTargetPositionZ() const { return fTargetPos[2]; }
/// \brief Gets ca::FieldRegion object at primary vertex
- const FieldRegion& GetVertexFieldRegion() const { return fVertexFieldRegion; }
+ const FieldRegionV& GetVertexFieldRegion() const { return fVertexFieldRegion; }
/// \brief Gets ca::FieldValue object at primary vertex
- const FieldValue& GetVertexFieldValue() const { return fVertexFieldValue; }
+ const FieldValueV& GetVertexFieldValue() const { return fVertexFieldValue; }
/// \brief Gets random seed
///
@@ -256,10 +256,10 @@ namespace cbm::algo::ca
alignas(constants::misc::Alignment) std::array<fvec, 3> fTargetPos{Undef<float>, Undef<float>, Undef<float>};
/// Field value object at primary vertex (between target and the first station)
- alignas(constants::misc::Alignment) ca::FieldValue fVertexFieldValue{};
+ alignas(constants::misc::Alignment) ca::FieldValueV fVertexFieldValue{};
/// Field region object at primary vertex (between target and the first station)
- alignas(constants::misc::Alignment) ca::FieldRegion fVertexFieldRegion{};
+ alignas(constants::misc::Alignment) ca::FieldRegionV fVertexFieldRegion{};
/// Array of stations
alignas(constants::misc::Alignment) StationsContainer_t fStations{};
diff --git a/algo/ca/core/pars/CaStation.h b/algo/ca/core/pars/CaStation.h
index f9984340db88f20ead65637fec84cff7ec725c12..b7a072530c19228d7b494b8c51ccfb26fc55f8bf 100644
--- a/algo/ca/core/pars/CaStation.h
+++ b/algo/ca/core/pars/CaStation.h
@@ -29,13 +29,13 @@ namespace cbm::algo::ca
DataT Xmax = constants::Undef<DataT>; ///< min radius of the station [cm]
DataT Ymax = constants::Undef<DataT>; ///< max radius of the station [cm]
- FieldSlice fieldSlice{}; ///< Magnetic field near the station
+ FieldSlice<DataT> fieldSlice{}; ///< Magnetic field near the station
Station() = default;
/// \brief Copy constructor with type conversion
- template<typename DataI>
- Station(const Station<DataI>& other)
+ template<typename DataIn>
+ Station(const Station<DataIn>& other)
: type(other.GetType())
, timeInfo(other.GetTimeStatus())
, fieldStatus(other.GetFieldStatus())
@@ -77,24 +77,24 @@ namespace cbm::algo::ca
int GetFieldStatus() const { return fieldStatus; }
/// \brief Gets z-position of the station
- template<typename DataO = DataT>
- DataO GetZ() const
+ template<typename DataOut = DataT>
+ DataOut GetZ() const
{
- return utils::GetValSimpl<DataT, DataO>(fZ);
+ return utils::simd::Cast<DataT, DataOut>(fZ);
}
/// \brief Gets limit of the station size in x-axis direction
- template<typename DataO = DataT>
- DataO GetXmax() const
+ template<typename DataOut = DataT>
+ DataOut GetXmax() const
{
- return utils::GetValSimpl<DataT, DataO>(Xmax);
+ return utils::simd::Cast<DataT, DataOut>(Xmax);
}
/// \brief Gets limit of the station size in y-axis direction
- template<typename DataO = DataT>
- DataO GetYmax() const
+ template<typename DataOut = DataT>
+ DataOut GetYmax() const
{
- return utils::GetValSimpl<DataT, DataO>(Ymax);
+ return utils::simd::Cast<DataT, DataOut>(Ymax);
}
/// \brief String representation of class contents
@@ -106,12 +106,30 @@ namespace cbm::algo::ca
class StationV : public Station<fvec> {
using Station<fvec>::Station; /// inherit all the constructors
+
+ public:
+ template<typename DataIn>
+ StationV(const Station<DataIn>& other) : Station<fvec>(other)
+ {
+ }
};
class StationF : public Station<float> {
using Station<float>::Station; /// inherit all the constructors
+
+ public:
+ template<typename DataIn>
+ StationF(const Station<DataIn>& other) : Station<float>(other)
+ {
+ }
};
class StationD : public Station<double> {
using Station<double>::Station; /// inherit all the constructors
+
+ public:
+ template<typename DataIn>
+ StationD(const Station<DataIn>& other) : Station<double>(other)
+ {
+ }
};
} // namespace cbm::algo::ca
diff --git a/algo/ca/core/tracking/CaCloneMerger.cxx b/algo/ca/core/tracking/CaCloneMerger.cxx
index 274f6123fb80bc53a83873844bf3a6f757d3990b..28bae59f2ca0de08b02b551e5b3e5dc744383b7d 100644
--- a/algo/ca/core/tracking/CaCloneMerger.cxx
+++ b/algo/ca/core/tracking/CaCloneMerger.cxx
@@ -89,8 +89,8 @@ void CloneMerger::Exec(Vector<Track>& extTracks, Vector<ca::HitIndex_t>& extReco
TrackParamV& Tb = fitB.Tr();
TrackParamV& Tf = fitF.Tr();
- ca::FieldValue fBm, fBb, fBf _fvecalignment;
- ca::FieldRegion fld _fvecalignment;
+ ca::FieldValueV fBm, fBb, fBf _fvecalignment;
+ ca::FieldRegionV fld _fvecalignment;
// Max length for merging
unsigned char maxLengthForMerge = static_cast<unsigned char>(frAlgo.GetParameters().GetNstationsActive() - 3);
diff --git a/algo/ca/core/tracking/CaFramework.cxx b/algo/ca/core/tracking/CaFramework.cxx
index e5fc12f6ec797e7183d3b6bfdf670e61750e77b2..173093f53b4e99b74c6638cbe4184b603289f30a 100644
--- a/algo/ca/core/tracking/CaFramework.cxx
+++ b/algo/ca/core/tracking/CaFramework.cxx
@@ -61,7 +61,7 @@ void Framework::ReceiveParameters(Parameters&& parameters)
fGhostSuppression = fParameters.GetGhostSuppression();
- ca::FieldRegion::ForceUseOfOriginalField(fParameters.DevIsUseOfOriginalField());
+ ca::FieldRegionV::ForceUseOfOriginalField(fParameters.DevIsUseOfOriginalField());
}
int Framework::GetMcTrackIdForCaHit(int /*iHit*/) const
diff --git a/algo/ca/core/tracking/CaFramework.h b/algo/ca/core/tracking/CaFramework.h
index 4f5344c315f0893146b0cca39f49c2d8e2a2d896..ab3019b5ca3b3865ecc746733056c5ec55895c3b 100644
--- a/algo/ca/core/tracking/CaFramework.h
+++ b/algo/ca/core/tracking/CaFramework.h
@@ -318,7 +318,7 @@ namespace cbm::algo::ca
bool fIsTargetField{false}; ///< is the magnetic field present at the target
- ca::FieldValue fTargB _fvecalignment{}; // field in the target point (modifiable, do not touch!!)
+ ca::FieldValueV fTargB _fvecalignment{}; // field in the target point (modifiable, do not touch!!)
ca::MeasurementXy<fvec> TargetMeasurement _fvecalignment{}; // target constraint [cm]
int fGhostSuppression{1}; // NOTE: Should be equal to 0 in TRACKS_FROM_TRIPLETS mode!
diff --git a/algo/ca/core/tracking/CaTrackExtender.cxx b/algo/ca/core/tracking/CaTrackExtender.cxx
index 100714260d76a1b29f13a11cbe2a60139f7791b6..8b04c4a2548c40371185b53920fa742be955dc6b 100644
--- a/algo/ca/core/tracking/CaTrackExtender.cxx
+++ b/algo/ca/core/tracking/CaTrackExtender.cxx
@@ -98,7 +98,7 @@ void TrackExtender::FitBranchFast(const ca::Branch& t, TrackParamV& Tout, const
T.C10() = hit0.dXY();
T.C11() = hit0.dY2();
- ca::FieldRegion fld _fvecalignment;
+ ca::FieldRegionV fld _fvecalignment;
fvec fldZ0 = sta1.fZ; // suppose field is smoth
fvec fldZ1 = sta2.fZ;
fvec fldZ2 = sta0.fZ;
@@ -181,7 +181,7 @@ void TrackExtender::FindMoreHits(ca::Branch& t, TrackParamV& Tout, const bool up
fvec x2 = hit2.X();
fvec y2 = hit2.Y();
- ca::FieldRegion fld _fvecalignment;
+ ca::FieldRegionV fld _fvecalignment;
fvec fldZ0 = sta1.fZ;
fvec fldZ1 = sta2.fZ;
fvec fldZ2 = sta0.fZ;
diff --git a/algo/ca/core/tracking/CaTrackFit.cxx b/algo/ca/core/tracking/CaTrackFit.cxx
index 7707ea27dafa780b78ef2b3b7fa48a35447ba29b..136c47d422eaa226d22789f5d5e7460dece11d6a 100644
--- a/algo/ca/core/tracking/CaTrackFit.cxx
+++ b/algo/ca/core/tracking/CaTrackFit.cxx
@@ -599,7 +599,7 @@ namespace cbm::algo::ca
void TrackFit::Extrapolate // extrapolates track parameters and returns jacobian for extrapolation of CovMatrix
(fvec z_out, // extrapolate to this z position
- const ca::FieldRegion& F)
+ const ca::FieldRegionV& F)
{
// use Q/p linearisation at fQp0
@@ -613,7 +613,7 @@ namespace cbm::algo::ca
void TrackFit::ExtrapolateStep // extrapolates track parameters and returns jacobian for extrapolation of CovMatrix
(fvec zOut, // extrapolate to this z position
- const ca::FieldRegion& Field)
+ const ca::FieldRegionV& Field)
{
// use Q/p linearisation at fQp0
// implementation of the Runge-Kutta method without optimization
@@ -809,7 +809,7 @@ namespace cbm::algo::ca
}
- void TrackFit::ExtrapolateLine(fvec z_out, const ca::FieldRegion& F)
+ void TrackFit::ExtrapolateLine(fvec z_out, const ca::FieldRegionV& F)
{
// extrapolate the track assuming fQp0 == 0
// TODO: write special simplified procedure
@@ -1040,7 +1040,7 @@ namespace cbm::algo::ca
}
- void TrackFit::GetExtrapolatedXYline(fvec z, const ca::FieldRegion& F, fvec& extrX, fvec& extrY,
+ void TrackFit::GetExtrapolatedXYline(fvec z, const ca::FieldRegionV& F, fvec& extrX, fvec& extrY,
std::array<fvec, ca::TrackParamV::kNtrackParam>& Jx,
std::array<fvec, ca::TrackParamV::kNtrackParam>& Jy) const
{
@@ -1094,7 +1094,7 @@ namespace cbm::algo::ca
}
- void TrackFit::FilterWithTargetAtLine(fvec targZ, const ca::MeasurementXy<fvec>& targXY, const ca::FieldRegion& F)
+ void TrackFit::FilterWithTargetAtLine(fvec targZ, const ca::MeasurementXy<fvec>& targXY, const ca::FieldRegionV& F)
{
// Add the target constraint to a straight line track
diff --git a/algo/ca/core/tracking/CaTrackFit.h b/algo/ca/core/tracking/CaTrackFit.h
index 363073dd3f08c2830cc0835a0bc0943ee844ec66..61950cb662891957fe3f58695abca97fe4f7dd12 100644
--- a/algo/ca/core/tracking/CaTrackFit.h
+++ b/algo/ca/core/tracking/CaTrackFit.h
@@ -117,14 +117,14 @@ namespace cbm::algo::ca
/// it can do several extrapolation steps if the Z is far away
/// \param z - Z coordinate to extrapolate to
/// \param F - field region
- void Extrapolate(fvec z, const ca::FieldRegion& F);
+ void Extrapolate(fvec z, const ca::FieldRegionV& F);
/// extrapolate the track to the given Z using the field F
/// it does extrapolation in one step
- void ExtrapolateStep(fvec z_out, const ca::FieldRegion& F);
+ void ExtrapolateStep(fvec z_out, const ca::FieldRegionV& F);
/// extrapolate the track to the given Z using linearization at the straight line
- void ExtrapolateLine(fvec z_out, const ca::FieldRegion& F);
+ void ExtrapolateLine(fvec z_out, const ca::FieldRegionV& F);
/// extrapolate the track to the given Z assuming no magnetic field
void ExtrapolateLineNoField(fvec z_out);
@@ -159,7 +159,7 @@ namespace cbm::algo::ca
/// special utilities needed by the combinatorial track finder
/// extrapolate track as a line, return the extrapolated X, Y and the Jacobians
- void GetExtrapolatedXYline(fvec z, const ca::FieldRegion& F, fvec& extrX, fvec& extrY,
+ void GetExtrapolatedXYline(fvec z, const ca::FieldRegionV& F, fvec& extrX, fvec& extrY,
std::array<fvec, ca::TrackParamV::kNtrackParam>& Jx,
std::array<fvec, ca::TrackParamV::kNtrackParam>& Jy) const;
@@ -189,7 +189,7 @@ namespace cbm::algo::ca
fvec ExtrapolateLineDxy(fvec z_out) const;
/// add target measuremet to the track using linearisation at a straight line
- void FilterWithTargetAtLine(fvec targZ, const ca::MeasurementXy<fvec>& targXYInfo, const ca::FieldRegion& F);
+ void FilterWithTargetAtLine(fvec targZ, const ca::MeasurementXy<fvec>& targXYInfo, const ca::FieldRegionV& F);
/// \brief Approximate mean energy loss with Bethe-Bloch formula
/// \param bg2 (beta*gamma)^2
diff --git a/algo/ca/core/tracking/CaTrackFitter.cxx b/algo/ca/core/tracking/CaTrackFitter.cxx
index df42286c8e635a1ffb76132946ec9a5e0126f082..b85f1cfaf58547efd1fa5fc8408fb6de9d099dc4 100644
--- a/algo/ca/core/tracking/CaTrackFitter.cxx
+++ b/algo/ca/core/tracking/CaTrackFitter.cxx
@@ -37,12 +37,12 @@ void TrackFitter::FitCaTracks()
// static ca::FieldValue fldB0, fldB1, fldB2 _fvecalignment;
// static ca::FieldRegion fld _fvecalignment;
- ca::FieldValue fldB0, fldB1, fldB2 _fvecalignment;
- ca::FieldRegion fld _fvecalignment;
+ ca::FieldValueV fldB0, fldB1, fldB2 _fvecalignment;
+ ca::FieldRegionV fld _fvecalignment;
- ca::FieldValue fldB01, fldB11, fldB21 _fvecalignment;
- ca::FieldRegion fld1 _fvecalignment;
+ ca::FieldValueV fldB01, fldB11, fldB21 _fvecalignment;
+ ca::FieldRegionV fld1 _fvecalignment;
const int nStations = frAlgo.GetParameters().GetNstationsActive();
int nTracks_SIMD = fvec::size();
@@ -95,7 +95,7 @@ void TrackFitter::FitCaTracks()
fvec z_start;
fvec z_end;
- ca::FieldValue fB[constants::size::MaxNstations], fB_temp _fvecalignment;
+ ca::FieldValueV fB[constants::size::MaxNstations], fB_temp _fvecalignment;
fvec ZSta[constants::size::MaxNstations];
@@ -287,7 +287,7 @@ void TrackFitter::FitCaTracks()
vtxInfo.SetDy2(1.e-8);
if (ca::Framework::TrackingMode::kGlobal == frAlgo.fTrackingMode) {
- ca::FieldRegion fldFull _fvecalignment;
+ ca::FieldRegionV fldFull _fvecalignment;
fldFull.SetUseOriginalField(true);
fitpv.SetMaxExtrapolationStep(1.);
for (int vtxIter = 0; vtxIter < 2; vtxIter++) {
diff --git a/algo/ca/core/tracking/CaTripletConstructor.cxx b/algo/ca/core/tracking/CaTripletConstructor.cxx
index 3d1d0b1763e91e5eb4d87a11a2417b4b3b63b08c..458d5519b20de70712cc107b47faf4609daa9743 100644
--- a/algo/ca/core/tracking/CaTripletConstructor.cxx
+++ b/algo/ca/core/tracking/CaTripletConstructor.cxx
@@ -104,8 +104,8 @@ void TripletConstructor::CreateTripletsForHit(int istal, int istam, int istar, c
fFit.SetQp0(fAlgo->fMaxInvMom);
- ca::FieldValue lB, mB, cB, rB _fvecalignment;
- ca::FieldValue l_B_global, centB_global _fvecalignment;
+ ca::FieldValueV lB, mB, cB, rB _fvecalignment;
+ ca::FieldValueV l_B_global, centB_global _fvecalignment;
// get the magnetic field along the track.
// (suppose track is straight line with origin in the target)
@@ -148,18 +148,18 @@ void TripletConstructor::CreateTripletsForHit(int istal, int istam, int istar, c
// field made by the left hit, the target and the station istac in-between.
// is used for extrapolation to the target and to the middle hit
- ca::FieldRegion fld0;
+ ca::FieldRegionV fld0;
{
- ca::FieldValue B0 = fFld0Sta[0]->fieldSlice.GetFieldValueForLine(T);
- ca::FieldValue B1 = fFld0Sta[1]->fieldSlice.GetFieldValueForLine(T);
+ ca::FieldValueV B0 = fFld0Sta[0]->fieldSlice.GetFieldValueForLine(T);
+ ca::FieldValueV B1 = fFld0Sta[1]->fieldSlice.GetFieldValueForLine(T);
fld0.Set(fAlgo->fTargB, fAlgo->fTargZ, B0, fFld0Sta[0]->fZ, B1, fFld0Sta[1]->fZ);
}
{ // field, made by the left hit, the middle station and the right station
// Will be used for extrapolation to the right hit
- ca::FieldValue B0 = fFld1Sta[0]->fieldSlice.GetFieldValueForLine(T);
- ca::FieldValue B1 = fFld1Sta[1]->fieldSlice.GetFieldValueForLine(T);
- ca::FieldValue B2 = fFld1Sta[2]->fieldSlice.GetFieldValueForLine(T);
+ ca::FieldValueV B0 = fFld1Sta[0]->fieldSlice.GetFieldValueForLine(T);
+ ca::FieldValueV B1 = fFld1Sta[1]->fieldSlice.GetFieldValueForLine(T);
+ ca::FieldValueV B2 = fFld1Sta[2]->fieldSlice.GetFieldValueForLine(T);
fFldL.Set(B0, fFld1Sta[0]->fZ, B1, fFld1Sta[1]->fZ, B2, fFld1Sta[2]->fZ);
}
@@ -445,9 +445,9 @@ void TripletConstructor::FitTriplets()
// find the field along the track
- ca::FieldValue B[3] _fvecalignment;
- ca::FieldRegion fld _fvecalignment;
- ca::FieldRegion fldTarget _fvecalignment;
+ ca::FieldValueV B[3] _fvecalignment;
+ ca::FieldRegionV fld _fvecalignment;
+ ca::FieldRegionV fldTarget _fvecalignment;
fvec tx[3] = {(x[1] - x[0]) / (z[1] - z[0]), (x[2] - x[0]) / (z[2] - z[0]), (x[2] - x[1]) / (z[2] - z[1])};
fvec ty[3] = {(y[1] - y[0]) / (z[1] - z[0]), (y[2] - y[0]) / (z[2] - z[0]), (y[2] - y[1]) / (z[2] - z[1])};
diff --git a/algo/ca/core/tracking/CaTripletConstructor.h b/algo/ca/core/tracking/CaTripletConstructor.h
index 68234e0682a94b18daeaccf3ad378741bcc31be9..b58ac7a28ec2861cf27636f4d8daaaa59dc6c84b 100644
--- a/algo/ca/core/tracking/CaTripletConstructor.h
+++ b/algo/ca/core/tracking/CaTripletConstructor.h
@@ -104,7 +104,7 @@ namespace cbm::algo::ca
ca::HitIndex_t fIhitL; ///< index of the left hit in fAlgo->fWindowHits
TrackParamV fTrL;
- ca::FieldRegion fFldL;
+ ca::FieldRegionV fFldL;
Vector<ca::HitIndex_t> fHitsM_2{"TripletConstructor::fHitsM_2"};
Vector<TrackParamV> fTracks_2{"TripletConstructor::fTracks_2"};
diff --git a/algo/ca/core/utils/CaUtils.h b/algo/ca/core/utils/CaUtils.h
index 1f14d08f0f2a43ad3e0477cb1846c8a60b700b90..d6f03ff314606b6938d3f3096df7b056233cdf80 100644
--- a/algo/ca/core/utils/CaUtils.h
+++ b/algo/ca/core/utils/CaUtils.h
@@ -104,29 +104,6 @@ namespace cbm::algo::ca::utils
template<>
void CheckSimdVectorEquality(fvec v, const char* name);
- /// \brief template based common fvec/float/double variable converter
- template<typename DataT, typename DataO>
- inline DataO GetValSimpl(const DataT& val)
- {
- return static_cast<DataO>(val);
- }
-
- /// \brief template based common fvec/float/double variable converter
- /// fvec to float implementation
- template<>
- inline float GetValSimpl(const fvec& val)
- {
- return val[0];
- }
-
- /// \brief template based common fvec/float/double variable converter
- /// fvec to double implementation
- template<>
- inline double GetValSimpl(const fvec& val)
- {
- return static_cast<double>(val[0]);
- }
-
template<typename TdataA, typename TdataB, bool TDoAllA, bool TDoAllB>
class VecCopySpec {
public:
@@ -208,3 +185,138 @@ namespace cbm::algo::ca::utils
}
};
} // namespace cbm::algo::ca::utils
+
+/// Namespace contains compile-time constants definition for SIMD operations
+/// in the CA tracking algorithm
+///
+namespace cbm::algo::ca::utils::simd
+{
+ /// \brief Converts a value of type DataT to type DataOut
+ /// \details This function is a generic template that provides a flexible way for fvec/float/double conversions
+ /// \tparam DataT Input value data type
+ /// \tparam DataOut Converted value data type
+ /// \param val Input value of type DataT to be converted
+ /// \return Return val as DataOut
+ template<typename DataT, typename DataOut>
+ inline DataOut Cast(const DataT& val)
+ {
+ return static_cast<DataOut>(val);
+ }
+
+ /// \brief Converts a value of type DataT to type DataOut
+ /// \details This specialization extracts the first element of the input SIMD float vector (fvec) and returns it as a float
+ /// \param val Input value of type fvec to be converted
+ /// \return Return val[0] as float
+ template<>
+ inline float Cast(const fvec& val)
+ {
+ return val[0];
+ }
+
+ /// \brief Converts a value of type DataT to type DataOut
+ /// \details This specialization extracts the first element of the input SIMD float vector (fvec) and returns it as a double
+ /// \param val Input value of type fvec to be converted
+ /// \return Return val[0] as double
+ template<>
+ inline double Cast(const fvec& val)
+ {
+ return static_cast<double>(val[0]);
+ }
+
+ /// \brief Converts a value of type DataT at a specific index to type DataOut
+ /// \details This function is a generic template that provides a flexible way for fvec/float/double conversions
+ /// \tparam DataT Input value data type
+ /// \tparam DataOut Converted value data type
+ /// \param val Input value of type DataT to be converted
+ /// \param i The index indicating which element to convert, ignored for non-SIMD DataT
+ /// \return Return val as DataOut
+ template<typename DataT, typename DataOut>
+ inline DataOut Cast(const DataT& val, size_t /*i*/)
+ {
+ return static_cast<DataOut>(val);
+ }
+
+ /// \brief Converts a value of type DataT at a specific index to type DataOut
+ /// \details This specialization extracts the element at the specified index of the input SIMD float vector (fvec) and returns it as a float
+ /// \param val Input value of type fvec to be converted
+ /// \param i The index indicating which element to convert
+ /// \return Return val[i] as float
+ template<>
+ inline float Cast(const fvec& val, size_t i)
+ {
+ return val[i];
+ }
+
+ /// \brief Converts a value of type DataT at a specific index to type DataOut
+ /// \details This specialization extracts the element at the specified index of the input SIMD float vector (fvec) and returns it as a double
+ /// \param val Input value of type fvec to be converted
+ /// \param i The index indicating which element to convert.
+ /// \return Return val[i] as double
+ template<>
+ inline double Cast(const fvec& val, size_t i)
+ {
+ return static_cast<double>(val[i]);
+ }
+
+ /// \brief Returns the number of elements available for a given data type
+ /// \details This function is a template that provides a default implementation returning 1
+ /// \tparam DataT The data type for which the size is determined
+ /// \return Return 1
+ template<typename DataT>
+ inline size_t Size()
+ {
+ return 1;
+ }
+
+ /// \brief Returns the number of elements available for a given data type
+ /// \details This specialization returns the size of the SIMD float vector (fvec) using its static member function `size()`
+ /// \return Return number of elements in the SIMD
+ template<>
+ inline size_t Size<fvec>()
+ {
+ return fvec::size();
+ }
+
+ /// \brief Returns a value of the data type set to one
+ /// \details This function is a template that provides a default implementation returning 1
+ /// \tparam DataT The data type for which the value is generated
+ /// \return Return 1 as DataT
+ template<typename DataT>
+ inline DataT One()
+ {
+ return static_cast<DataT>(1);
+ }
+
+ /// \brief Returns a value of the data type set to one
+ /// \details This specialization returns a SIMD float vector (fvec) with all elements set to one using its static member function `One()`
+ /// \return Return SIMD vector of 1
+ template<>
+ inline fvec One<fvec>()
+ {
+ return fvec::One();
+ }
+
+ /// \brief Sets a value at a specific index in the output data
+ /// \details This function is a template that provides a default implementation for setting a value at a specific index in the output data
+ /// \tparam DataT The data type of the output data
+ /// \tparam DataIn The data type of the input data
+ /// \param out Reference to the output data where the value is to be set
+ /// \param in The input data from which the value is obtained
+ /// \param i The index at which the value is set, ignored for non-SIMD DataT or if both DataT and DataIn are fvec
+ template<typename DataT, typename DataIn>
+ inline void SetEntry(DataT& out, DataIn in, size_t /*i*/)
+ {
+ out = static_cast<DataT>(in);
+ }
+
+ /// \brief Sets a value at a specific index in the output data
+ /// \details This specialization sets a value at the specified index in the output SIMD float vector (fvec) using the input float value
+ /// \param out Reference to the output float vector (fvec) where the value is to be set
+ /// \param in The input non-SIMD value from which the value is obtained
+ /// \param i The index at which the value is set
+ template<typename DataIn>
+ inline void SetEntry(fvec& out, DataIn in, size_t i)
+ {
+ out[i] = in;
+ }
+} // namespace cbm::algo::ca::utils::simd
diff --git a/reco/KF/CbmKfTrackFitter.cxx b/reco/KF/CbmKfTrackFitter.cxx
index e14e9b3a7d1726ed86da3d0ccb7a327e86f9fba7..827f6cc87c7be8edfba9ec72fbe7ac6d086f08d3 100644
--- a/reco/KF/CbmKfTrackFitter.cxx
+++ b/reco/KF/CbmKfTrackFitter.cxx
@@ -480,7 +480,7 @@ void CbmKfTrackFitter::FitTrack(CbmKfTrackFitter::Track& t)
fFit.SetParticleMass(fMass);
t.fIsMsQp0Set = false;
- ca::FieldRegion field _fvecalignment;
+ ca::FieldRegionV field _fvecalignment;
field.SetUseOriginalField();
int nNodes = t.fNodes.size();
diff --git a/reco/KF/CbmKfTrackFitter.h b/reco/KF/CbmKfTrackFitter.h
index f9ba682f713a29aca8a82adb71bd5b59a3cd2ee7..58b618ad870a4f7066a01a9d800637e4219b7b01 100644
--- a/reco/KF/CbmKfTrackFitter.h
+++ b/reco/KF/CbmKfTrackFitter.h
@@ -29,7 +29,7 @@ class TClonesArray;
namespace cbm::algo::ca
{
- class FieldRegion;
+ class FieldRegionV;
} // namespace cbm::algo::ca
namespace
diff --git a/reco/KF/ParticleFitter/CbmL1PFFitter.cxx b/reco/KF/ParticleFitter/CbmL1PFFitter.cxx
index f5eba8ed2cbe727718503073ec55b8011711d697..1d2f7ee5f845994c2c90629f74bdc53b217ab8c7 100644
--- a/reco/KF/ParticleFitter/CbmL1PFFitter.cxx
+++ b/reco/KF/ParticleFitter/CbmL1PFFitter.cxx
@@ -56,7 +56,7 @@ namespace NS_L1TrackFitter
using namespace NS_L1TrackFitter;
-inline void CbmL1PFFitter::PFFieldRegion::setFromL1FieldRegion(const ca::FieldRegion& fld, int i)
+inline void CbmL1PFFitter::PFFieldRegion::setFromL1FieldRegion(const ca::FieldRegionV& fld, int i)
{
fP[0] = fld.cx0[i];
fP[1] = fld.cx1[i];
@@ -73,7 +73,7 @@ inline void CbmL1PFFitter::PFFieldRegion::setFromL1FieldRegion(const ca::FieldRe
fP[9] = fld.z0[i];
}
-inline void CbmL1PFFitter::PFFieldRegion::getL1FieldRegion(ca::FieldRegion& fld, int i)
+inline void CbmL1PFFitter::PFFieldRegion::getL1FieldRegion(ca::FieldRegionV& fld, int i)
{
fld.cx0[i] = fP[0];
fld.cx1[i] = fP[1];
@@ -90,7 +90,7 @@ inline void CbmL1PFFitter::PFFieldRegion::getL1FieldRegion(ca::FieldRegion& fld,
fld.z0[i] = fP[9];
}
-inline CbmL1PFFitter::PFFieldRegion::PFFieldRegion(const ca::FieldRegion& fld, int i) { setFromL1FieldRegion(fld, i); }
+inline CbmL1PFFitter::PFFieldRegion::PFFieldRegion(const ca::FieldRegionV& fld, int i) { setFromL1FieldRegion(fld, i); }
CbmL1PFFitter::CbmL1PFFitter() {}
@@ -165,8 +165,8 @@ void CbmL1PFFitter::Fit(std::vector<CbmStsTrack>& Tracks, const std::vector<CbmM
Initialize();
- ca::FieldValue b0, b1, b2 _fvecalignment;
- ca::FieldRegion fld _fvecalignment;
+ ca::FieldValueV b0, b1, b2 _fvecalignment;
+ ca::FieldRegionV fld _fvecalignment;
// fld.SetUseOriginalField();
static int nHits = CbmL1::Instance()->fpAlgo->GetParameters().GetNstationsActive();
@@ -196,8 +196,8 @@ void CbmL1PFFitter::Fit(std::vector<CbmStsTrack>& Tracks, const std::vector<CbmM
fvec dt2_first, dt2_last;
fvec z0, z1, z2, dz, z_start, z_end;
- ca::FieldValue fB[nHits];
- ca::FieldValue fB_temp _fvecalignment;
+ ca::FieldValueV fB[nHits];
+ ca::FieldValueV fB_temp _fvecalignment;
unsigned short N_vTracks = Tracks.size();
@@ -506,8 +506,8 @@ void CbmL1PFFitter::GetChiToVertex(vector<CbmStsTrack>& Tracks, vector<PFFieldRe
field.reserve(Tracks.size());
- ca::FieldRegion fld _fvecalignment;
- ca::FieldValue fB[3], fB_temp _fvecalignment;
+ ca::FieldRegionV fld _fvecalignment;
+ ca::FieldValueV fB[3], fB_temp _fvecalignment;
fvec zField[3];
unsigned short N_vTracks = Tracks.size();
@@ -644,7 +644,7 @@ void CbmL1PFFitter::CalculateFieldRegion(vector<CbmStsTrack>& Tracks, vector<PFF
field.reserve(Tracks.size());
- ca::FieldRegion fld _fvecalignment;
+ ca::FieldRegionV fld _fvecalignment;
int nTracks_SIMD = fvec::size();
TrackParamV T; // fitting parametr coresponding to current track
@@ -653,7 +653,7 @@ void CbmL1PFFitter::CalculateFieldRegion(vector<CbmStsTrack>& Tracks, vector<PFF
int ista;
const ca::StationV* sta = CbmL1::Instance()->fpAlgo->GetParameters().GetStations().begin();
- ca::FieldValue fB[3], fB_temp _fvecalignment;
+ ca::FieldValueV fB[3], fB_temp _fvecalignment;
fvec zField[3];
unsigned short N_vTracks = Tracks.size();
@@ -715,7 +715,7 @@ void CbmL1PFFitter::CalculateFieldRegionAtLastPoint(vector<CbmStsTrack>& Tracks,
field.reserve(Tracks.size());
- ca::FieldRegion fld _fvecalignment;
+ ca::FieldRegionV fld _fvecalignment;
int nTracks_SIMD = fvec::size();
TrackParamV T; // fitting parametr coresponding to current track
@@ -724,7 +724,7 @@ void CbmL1PFFitter::CalculateFieldRegionAtLastPoint(vector<CbmStsTrack>& Tracks,
int ista;
const ca::StationV* sta = CbmL1::Instance()->fpAlgo->GetParameters().GetStations().begin();
- ca::FieldValue fB[3], fB_temp _fvecalignment;
+ ca::FieldValueV fB[3], fB_temp _fvecalignment;
fvec zField[3];
unsigned short N_vTracks = Tracks.size();
diff --git a/reco/KF/ParticleFitter/CbmL1PFFitter.h b/reco/KF/ParticleFitter/CbmL1PFFitter.h
index 94617f8e29bff2f422733fb8c0c2fb625e8f6900..a7433653925ec4c337e314a7992181d82a28b30d 100644
--- a/reco/KF/ParticleFitter/CbmL1PFFitter.h
+++ b/reco/KF/ParticleFitter/CbmL1PFFitter.h
@@ -27,7 +27,7 @@ class CbmStsHit;
class CbmStsTrack;
namespace cbm::algo::ca
{
- class FieldRegion;
+ class FieldRegionV;
}
using namespace cbm::algo;
@@ -40,9 +40,9 @@ class CbmL1PFFitter {
// A container for parameters of ca::FieldRegion
struct PFFieldRegion {
PFFieldRegion() {}
- PFFieldRegion(const ca::FieldRegion&, int i);
- void setFromL1FieldRegion(const ca::FieldRegion&, int i);
- void getL1FieldRegion(ca::FieldRegion&, int i);
+ PFFieldRegion(const ca::FieldRegionV&, int i);
+ void setFromL1FieldRegion(const ca::FieldRegionV&, int i);
+ void getL1FieldRegion(ca::FieldRegionV&, int i);
float fP[10]{0.};
};
diff --git a/reco/L1/CbmL1Performance.cxx b/reco/L1/CbmL1Performance.cxx
index f55571a4465c1bb058ed590181cd66bbd8e437aa..1a15402da2426e41dfced6b68e2f16d2e3ed05cf 100644
--- a/reco/L1/CbmL1Performance.cxx
+++ b/reco/L1/CbmL1Performance.cxx
@@ -1299,7 +1299,7 @@ void CbmL1::TrackFitPerformance()
if (!mc.IsPrimary()) { // secondary
{ // extrapolate to vertex
- ca::FieldRegion fld _fvecalignment;
+ ca::FieldRegionV fld _fvecalignment;
fld.SetUseOriginalField();
fit.Extrapolate(mc.GetStartZ(), fld);
// add material
@@ -1367,7 +1367,7 @@ void CbmL1::TrackFitPerformance()
else { // primary
{ // extrapolate to vertex
- ca::FieldRegion fld _fvecalignment;
+ ca::FieldRegionV fld _fvecalignment;
fld.SetUseOriginalField();
// add material
@@ -1516,7 +1516,7 @@ void CbmL1::FillFitHistos(TrackParamV& track, const cbm::ca::tools::MCPoint& mcP
//fit.SetMaxExtrapolationStep(10.);
fit.SetDoFitVelocity(true);
fit.SetTrack(track);
- ca::FieldRegion fld _fvecalignment;
+ ca::FieldRegionV fld _fvecalignment;
fld.SetUseOriginalField();
fit.Extrapolate(mcP.GetZOut(), fld);
track = fit.Tr();
@@ -1597,8 +1597,8 @@ void CbmL1::FieldApproxCheck()
static_cast<int>(NbinsY + 1), -(Ymax + ddy / 2.), (Ymax + ddy / 2.));
Double_t r[3], B[3];
- ca::FieldSlice FSl;
- ca::FieldValue B_L1;
+ ca::FieldSliceV FSl;
+ ca::FieldValueV B_L1;
Double_t bbb, bbb_L1;
const int M = 5; // polinom order
diff --git a/reco/L1/catools/CaToolsField.h b/reco/L1/catools/CaToolsField.h
index c3b04a3e328bd93fbf2cf48817270af91865b6dc..765920772cfd17c250d769da8b8dd98b7468a569 100644
--- a/reco/L1/catools/CaToolsField.h
+++ b/reco/L1/catools/CaToolsField.h
@@ -26,7 +26,7 @@ namespace cbm::ca::tools
FairRunAna::Instance()->GetField()->GetFieldValue(pos, B);
return std::tuple(B[0], B[1], B[2]);
};
- ca::FieldRegion::SetOriginalField(fld);
+ ca::FieldRegionV::SetOriginalField(fld);
}
} // namespace cbm::ca::tools
diff --git a/reco/L1/qa/CbmCaTrackFitQa.cxx b/reco/L1/qa/CbmCaTrackFitQa.cxx
index e7b612b521b22aa6b384310d47844d1be6bdfe6d..5c4de0006d15fc6bf152a412d9bd690526145fd3 100644
--- a/reco/L1/qa/CbmCaTrackFitQa.cxx
+++ b/reco/L1/qa/CbmCaTrackFitQa.cxx
@@ -164,7 +164,7 @@ void TrackFitQa::Fill(const TrackParamV& trPar, const tools::MCPoint& mcPoint, b
fitter.SetMask(fmask::One());
fitter.SetDoFitVelocity(true);
fitter.SetTrack(trPar);
- cbm::algo::ca::FieldRegion fieldRegion;
+ cbm::algo::ca::FieldRegionV fieldRegion;
fieldRegion.SetUseOriginalField(); // Precised magnetic field is used
fitter.Extrapolate(mcPoint.GetZ(), fieldRegion);
diff --git a/reco/L1/qa/CbmCaTrackTypeQa.h b/reco/L1/qa/CbmCaTrackTypeQa.h
index 82c34998ffcdf3872cebcd984abd29790425dd28..1ebce85fcbdd108ac9731b6adf08b8f832732cab 100644
--- a/reco/L1/qa/CbmCaTrackTypeQa.h
+++ b/reco/L1/qa/CbmCaTrackTypeQa.h
@@ -277,7 +277,7 @@ namespace cbm::ca
TProfile* fph_stations_hit = nullptr; ///< Average number of stations with hit
ca::TrackFit fTrackFit; ///< Track fitter
- ca::FieldRegion fFieldRegion; ///< Magnetic field
+ ca::FieldRegionV fFieldRegion; ///< Magnetic field
int fCounterMC = 0; ///< Counter of MC tracks
int fCounterClones = 0; ///< Counter of clone tracks