diff --git a/algo/ca/core/pars/CaDefs.h b/algo/ca/core/pars/CaDefs.h
index 89b8d79a9fbbac7decdd8b899178968c13fbde2a..43507f8cfaf94e169419f52b6e809f93786090a0 100644
--- a/algo/ca/core/pars/CaDefs.h
+++ b/algo/ca/core/pars/CaDefs.h
@@ -17,13 +17,10 @@
namespace cbm::algo::ca
{
- using KfFramework_t = cbm::algo::kf::Framework<float>;
- using KfParameters_t = cbm::algo::kf::Parameters<float>;
- using KfSetup_t = cbm::algo::kf::Setup<float>;
+ using KfFramework_t = cbm::algo::kf::Framework<fvec>;
using cbm::algo::kf::TrackParamBase;
using cbm::algo::kf::TrackParamV;
-
} // namespace cbm::algo::ca
/// Namespace contains compile-time constants definition for the CA tracking algorithm
diff --git a/algo/ca/core/pars/CaMaterialMap.h b/algo/ca/core/pars/CaMaterialMap.h
index 24f8afedbb1775cb7f00fe462c9d57804cbc3e19..f72ebb8d31dacda8ffec447dd86851b69a37a938 100644
--- a/algo/ca/core/pars/CaMaterialMap.h
+++ b/algo/ca/core/pars/CaMaterialMap.h
@@ -105,6 +105,10 @@ namespace cbm::algo::ca
/// \brief String representation of the object
std::string ToString() const;
+ /// \brief Comparison operator
+ /// \note To be used for the material map ordering along the z-axis
+ bool operator<(const MaterialMap& r) const { return fZref < r.fZref; }
+
private:
int fNbins = kfdefs::Undef<int>; ///< Number of rows (== N columns) in the material budget table
float fXYmax = kfdefs::Undef<float>; ///< Size of the station in x and y dimensions [cm]
diff --git a/algo/kf/core/CMakeLists.txt b/algo/kf/core/CMakeLists.txt
index 55646bf8c5e58f8105264a5aa4678c65fccac269..5644088641ff9575ab1d53b5bad4b58d5b3ff49f 100644
--- a/algo/kf/core/CMakeLists.txt
+++ b/algo/kf/core/CMakeLists.txt
@@ -11,13 +11,13 @@ set(SRCS
${CMAKE_CURRENT_SOURCE_DIR}/KfFramework.cxx
${CMAKE_CURRENT_SOURCE_DIR}/data/KfTrackParam.cxx
${CMAKE_CURRENT_SOURCE_DIR}/geo/KfMaterialMap.cxx
- ${CMAKE_CURRENT_SOURCE_DIR}/geo/KfGeoLayer.cxx
${CMAKE_CURRENT_SOURCE_DIR}/geo/KfTarget.cxx
${CMAKE_CURRENT_SOURCE_DIR}/geo/KfField.cxx
${CMAKE_CURRENT_SOURCE_DIR}/geo/KfFieldValue.cxx
${CMAKE_CURRENT_SOURCE_DIR}/geo/KfFieldSlice.cxx
${CMAKE_CURRENT_SOURCE_DIR}/geo/KfFieldRegion.cxx
${CMAKE_CURRENT_SOURCE_DIR}/geo/KfSetup.cxx
+ ${CMAKE_CURRENT_SOURCE_DIR}/geo/KfSetupSerializer.cxx
${CMAKE_CURRENT_SOURCE_DIR}/pars/KfParameters.cxx
${CMAKE_CURRENT_SOURCE_DIR}/utils/KfUtils.cxx
)
@@ -68,8 +68,11 @@ install(
KfFramework.h
KfDefs.h
data/KfTrackParam.h
+ geo/KfField.h
+ geo/KfFieldRegion.h
+ geo/KfFieldSlice.h
+ geo/KfFieldValue.h
geo/KfMaterialMap.h
- geo/KfGeoLayer.h
geo/KfSetup.h
geo/KfTarget.h
pars/KfParameters.h
diff --git a/algo/kf/core/KfDefs.h b/algo/kf/core/KfDefs.h
index 96588834c608aef76760318b618fef4def0f610f..612bc9ebe0f99181f055eb551ac6ecd706b86b9c 100644
--- a/algo/kf/core/KfDefs.h
+++ b/algo/kf/core/KfDefs.h
@@ -74,8 +74,8 @@ namespace cbm::algo::kf
namespace cbm::algo::kf::defs
{
// ----- Array sizes -------------------------------------------------------------------------------------------------
- constexpr int MaxNofFieldSlices = 64; ///< Max number of field slices
-
+ constexpr int MaxNofFieldSlices = 64; ///< Max number of field slices
+ constexpr int MaxNofMaterialLayers = 32; ///< Max number of the material layers
// ----- Control -----------------------------------------------------------------------------------------------------
constexpr int DebugLvl = 0; ///< Level of debug output
diff --git a/algo/kf/core/KfFramework.cxx b/algo/kf/core/KfFramework.cxx
index a11fe973de397e57c9495756612d5439e33d0565..21334b755857c0bbb7e005c87469e15ea38f7cbc 100644
--- a/algo/kf/core/KfFramework.cxx
+++ b/algo/kf/core/KfFramework.cxx
@@ -24,12 +24,12 @@ try {
errMsg << "\tParameters initialization errors: " << err.what() << '\n';
}
- try {
- fSetup.Init();
- }
- catch (const std::runtime_error& err) {
- errMsg << "\tSetup initialization errors: " << err.what() << '\n';
- }
+ //try {
+ // fSetup.Init();
+ //}
+ //catch (const std::runtime_error& err) {
+ // errMsg << "\tSetup initialization errors: " << err.what() << '\n';
+ //}
if (!errMsg.str().empty()) {
throw std::runtime_error(errMsg.str());
@@ -43,6 +43,10 @@ catch (const std::exception& err) {
return false;
}
-template class cbm::algo::kf::Framework<float>;
-template class cbm::algo::kf::Framework<double>;
-// template class cbm::algo::kf::Framework<fvec>;
+
+namespace cbm::algo::kf
+{
+ template class Framework<float>;
+ template class Framework<double>;
+ template class Framework<fvec>;
+} // namespace cbm::algo::kf
diff --git a/algo/kf/core/KfFramework.h b/algo/kf/core/KfFramework.h
index 902126523d425470db2005ff7e1411f79cf06638..ea5309a051dfd21534d00fef968a80788151cea4 100644
--- a/algo/kf/core/KfFramework.h
+++ b/algo/kf/core/KfFramework.h
@@ -17,8 +17,8 @@ namespace cbm::algo::kf
{
/// \class Framework
/// \brief Main class of the KfCore library
- /// \tparam DataT Underlying data-type
- template<typename DataT>
+ /// \tparam T Underlying floating point data-type
+ template<typename T>
class Framework {
public:
/// \brief Default constructor
@@ -46,20 +46,20 @@ namespace cbm::algo::kf
bool IsInitialized() const { return fbInitialized; }
/// \brief Parameters access
- Parameters<DataT>& Pars() { return fPars; }
+ kf::Parameters<T>& Pars() { return fPars; }
/// \brief Parameters access
- const kf::Parameters<DataT>& Pars() const { return fPars; }
+ const kf::Parameters<T>& Pars() const { return fPars; }
/// \brief Setup access (mutable)
- kf::Setup<DataT>& Setup() { return fSetup; }
+ kf::Setup<T>& Setup() { return fSetup; }
/// \brief Setup access
- const kf::Setup<DataT>& Setup() const { return fSetup; }
+ const kf::Setup<T>& Setup() const { return fSetup; }
private:
- kf::Parameters<DataT> fPars; ///< KF parameters
- kf::Setup<DataT> fSetup; ///< KF setup
+ kf::Parameters<T> fPars; ///< KF parameters
+ kf::Setup<T> fSetup; ///< KF setup
bool fbInitialized = false; ///< Initialization status
};
} // namespace cbm::algo::kf
diff --git a/algo/kf/core/geo/KfField.cxx b/algo/kf/core/geo/KfField.cxx
index c9ab45a6b6b70f69f18d845a2b6fcdb2c917c1ec..fc5ae8cd980673c5d8bc84e6329ce85d290172dc 100644
--- a/algo/kf/core/geo/KfField.cxx
+++ b/algo/kf/core/geo/KfField.cxx
@@ -14,6 +14,7 @@
#include <tuple>
using cbm::algo::kf::EFieldMode;
+using cbm::algo::kf::Field;
using cbm::algo::kf::FieldFactory;
using cbm::algo::kf::FieldFnFair_t;
using cbm::algo::kf::detail::FieldBase;
@@ -34,31 +35,39 @@ FieldBaseIntrpl_t<T>::FieldBase(const FieldBaseIntrpl_t<I>& other)
// ---------------------------------------------------------------------------------------------------------------------
//
-template<typename T>
-template<typename I>
-FieldBaseIntrpl_t<T>& FieldBaseIntrpl_t<T>::operator=(const FieldBaseIntrpl_t<I>& other)
+template<typename T, EFieldMode FldMode>
+std::string Field<T, FldMode>::ToString(int indentLevel) const
{
- if (this != &other) {
- for (size_t i = 0; i < defs::MaxNofFieldSlices; ++i) {
- fvFieldSlices[i] = utils::simd::Cast<I, T>(other.fvFieldSlices[i]);
+ constexpr char IndentChar = '\t';
+ std::stringstream msg;
+ std::string indent(indentLevel, IndentChar);
+ msg << indent << "Field near primary vertex:\n" << this->fPrimVertexField.ToString(indentLevel + 1) << '\n';
+ msg << indent << "Field type: " << static_cast<int>(this->fFieldType) << '\n';
+ if constexpr (FldMode == EFieldMode::Orig) {
+ msg << indent << "Original field function";
+ }
+ else if constexpr (FldMode == EFieldMode::Intrpl) {
+ msg << indent << "Field slices:";
+ for (const auto& fldSlice : this->fvFieldSlices) {
+ msg << indent << "\n - " << fldSlice.ToString(indentLevel + 1);
}
}
- return *this;
+ return msg.str();
}
// ---------------------------------------------------------------------------------------------------------------------
//
void FieldFactory::AddSliceReference(double halfSizeX, double halfSizeY, double zRef)
{
- if (!fSliceReferences.emplace().second) {
+ if (!fSliceReferences.emplace(halfSizeX, halfSizeY, zRef).second) {
std::stringstream msg;
- msg << "FieldContainerFactory::AddReference: attempt of adding another slice reference with zRef = " << zRef
+ msg << "FieldFactory::AddReference: attempt of adding another slice reference with zRef = " << zRef
<< "(halfSizeX = " << halfSizeX << ", halfSizeY = " << halfSizeY << ").\nThe next slice references were "
<< "added:";
for (const auto& el : fSliceReferences) {
msg << "\n\t- halfSizeX = " << el.fHalfSizeX << ", halfSizeY = " << el.fHalfSizeY << ", zRef = " << el.fRefZ;
}
- throw std::runtime_error(msg.str());
+ throw std::logic_error(msg.str());
}
}
diff --git a/algo/kf/core/geo/KfField.h b/algo/kf/core/geo/KfField.h
index 5a18714ebaec9ce08cb955f87ad486db63289270..98882d7c2286c4415860c7e1d292c3034aaf7896 100644
--- a/algo/kf/core/geo/KfField.h
+++ b/algo/kf/core/geo/KfField.h
@@ -60,9 +60,12 @@ namespace cbm::algo::kf
template<typename I>
FieldBase& operator=(const FieldBase<I, EFieldMode::Orig>& other)
{
- if (this != &other) {
- fFieldFn = other.fFieldFn;
+ if constexpr (std::is_same_v<I, T>) {
+ if (this == &other) {
+ return *this;
+ }
}
+ fFieldFn = other.fFieldFn;
return *this;
}
@@ -130,7 +133,18 @@ namespace cbm::algo::kf
/// \brief Copy assignment operator
/// \tparam I Underlying floating type of the source
template<typename I>
- FieldBase& operator=(const FieldBase<I, EFieldMode::Intrpl>& other);
+ FieldBase& operator=(const FieldBase<I, EFieldMode::Intrpl>& other)
+ {
+ if constexpr (std::is_same_v<I, T>) {
+ if (this == &other) {
+ return *this;
+ }
+ }
+ for (size_t i = 0; i < defs::MaxNofFieldSlices; ++i) {
+ fvFieldSlices[i] = utils::simd::Cast<I, T>(other.fvFieldSlices[i]);
+ }
+ return *this;
+ }
SlicesContainer_t fvFieldSlices; ///< Array of field slices
EFieldType fFieldType{EFieldType::Null}; ///< Field type
@@ -159,6 +173,7 @@ namespace cbm::algo::kf
static constexpr EFieldMode FieldType = FldMode;
using FieldRegion_t = typename FieldBase_t::FieldRegion_t;
+ using Real_t = T;
/// \brief Default constructor
Field() = default;
@@ -178,8 +193,10 @@ namespace cbm::algo::kf
template<typename I>
Field& operator=(const Field<I, FldMode>& other)
{
- if (this == &other) {
- return *this;
+ if constexpr (std::is_same_v<I, T>) {
+ if (this == &other) {
+ return *this;
+ }
}
fPrimVertexField = other.fPrimVertexField;
return FieldBase_t::operator=(other);
@@ -202,6 +219,10 @@ namespace cbm::algo::kf
/// \brief Gets field region near primary vertex
const FieldRegion_t& GetPrimVertexField() const { return fPrimVertexField; }
+ /// \brief String representation of the class
+ /// \param indentLevel Indent level of the string output
+ std::string ToString(int indentLevel) const;
+
private:
/// \brief Serialization function
friend class boost::serialization::access;
@@ -226,12 +247,6 @@ namespace cbm::algo::kf
double fHalfSizeY{defs::Undef<double>}; ///< Half-size of the slice in y-direction [cm]
double fRefZ{defs::Undef<double>}; ///< Reference z-position of the slice [cm]
- /// \brief Default constructor
- SliceRef() = default;
-
- /// \brief Destructor
- //~SliceRef() = default;
-
/// \brief Constructor
/// \param halfX Half-size of the slice in x-direction [cm]
/// \param halfY Half-size of the slice in y-direction [cm]
@@ -239,7 +254,7 @@ namespace cbm::algo::kf
SliceRef(double halfX, double halfY, double refZ) : fHalfSizeX(halfX), fHalfSizeY(halfY), fRefZ(refZ) {}
/// \brief Comparision operator
- friend constexpr bool operator<(const SliceRef& l, const SliceRef& r) { return l.fRefZ < r.fRefZ; }
+ bool operator<(const SliceRef& r) const { return fRefZ < r.fRefZ; }
};
public:
@@ -249,18 +264,12 @@ namespace cbm::algo::kf
/// \brief Copy constructor
FieldFactory(const FieldFactory&) = default;
- /// \brief Move constructor
- FieldFactory(FieldFactory&&) = default;
-
/// \brief Destructor
~FieldFactory() = default;
/// \brief Copy assignment operator
FieldFactory& operator=(const FieldFactory&) = default;
- /// \brief Move assignment operator
- FieldFactory& operator=(FieldFactory&&) = default;
-
/// \brief Adds a slice reference
/// \param halfSizeX Half-size of the slice in x-direction [cm]
/// \param halfSizeY Half-size of the slice in y-direction [cm]
@@ -320,19 +329,19 @@ namespace cbm::algo::kf
// Check initialization
if (std::any_of(fTarget.begin(), fTarget.end(), [](double x) { return !utils::IsFinite(x); })) {
- throw std::runtime_error("FieldFactory::MakeField: target is undefined");
+ throw std::logic_error("FieldFactory::MakeField: target is undefined");
}
if (!fSliceReferences.size()) { // TODO: Remove requirement of slice references
- throw std::runtime_error("FieldFactory::MakeField: no slice references were provided");
+ throw std::logic_error("FieldFactory::MakeField: no slice references were provided");
}
else if (fSliceReferences.size() > defs::MaxNofFieldSlices) {
std::stringstream msg;
- msg << "FieldFactory::MakeField: to many slice references are provided (" << fSliceReferences.size() << "), "
- << "the maximum allowed number of field slices is " << defs::MaxNofFieldSlices;
- throw std::runtime_error(msg.str());
+ msg << "FieldFactory::MakeField: too many slice references are provided (" << fSliceReferences.size() << "), "
+ << "the maximum allowed number of field slices is kf::defs::MaxNofFieldSlices = " << defs::MaxNofFieldSlices;
+ throw std::logic_error(msg.str());
}
if (!fFieldFn) {
- throw std::runtime_error("FieldFactory::CreateField: no field function is provided");
+ throw std::logic_error("FieldFactory::CreateField: no field function is provided");
}
// Initialize the Field object
diff --git a/algo/kf/core/geo/KfFieldRegion.h b/algo/kf/core/geo/KfFieldRegion.h
index 6906a1f9bb832ad76200124ca43a0ff4f89fc540..d7f40d39b117084427671f1f039770572c7cde14 100644
--- a/algo/kf/core/geo/KfFieldRegion.h
+++ b/algo/kf/core/geo/KfFieldRegion.h
@@ -98,9 +98,12 @@ namespace cbm::algo::kf
template<typename I>
FieldRegionBase& operator=(const FieldRegionBase<I, EFieldMode::Orig>& other)
{
- if (this != &other) {
- fFieldFn = other.fFieldFn;
+ if constexpr (std::is_same_v<I, T>) {
+ if (this == &other) {
+ return *this;
+ }
}
+ fFieldFn = other.fFieldFn;
return *this;
}
@@ -159,10 +162,13 @@ namespace cbm::algo::kf
template<typename I>
FieldRegionBase& operator=(const FieldRegionBase<I, EFieldMode::Orig>& other)
{
- if (this != &other) {
- this->fCoeff = utils::simd::Cast<I, T>(other.fCoeff);
- this->fZfirst = utils::simd::Cast<I, T>(other.fZfirst);
+ if constexpr (std::is_same_v<I, T>) {
+ if (this == &other) {
+ return *this;
+ }
}
+ this->fCoeff = utils::simd::Cast<I, T>(other.fCoeff);
+ this->fZfirst = utils::simd::Cast<I, T>(other.fZfirst);
return *this;
}
@@ -213,6 +219,10 @@ namespace cbm::algo::kf
static constexpr EFieldMode FieldType = FldMode;
using FieldValue_t = FieldValue<T>;
+ using Real_t = T;
+
+ /// \brief Default constructor
+ FieldRegion() = default;
/// \brief Copy constructor
/// \tparam I Underlying floating point type of the source
@@ -238,8 +248,10 @@ namespace cbm::algo::kf
template<typename I>
FieldRegion& operator=(const FieldRegion<I, FldMode>& other)
{
- if (this == &other) {
- return *this;
+ if constexpr (std::is_same_v<I, T>) {
+ if (this == &other) {
+ return *this;
+ }
}
fFieldType = other.fFieldType;
return FieldRegionBase_t::operator=(other);
@@ -260,15 +272,12 @@ namespace cbm::algo::kf
std::string ToString(int indentLevel = 0) const;
private:
- /// \brief Default constructor
- FieldRegion() = default;
-
/// \brief Serialization method
friend class boost::serialization::access;
template<class Archive>
void serialize(Archive& ar, const unsigned int)
{
- ar& boost::serialization::base_object<FieldRegionBase_t>;
+ ar& boost::serialization::base_object<FieldRegionBase_t>(*this);
ar& fFieldType;
}
diff --git a/algo/kf/core/geo/KfGeoLayer.cxx b/algo/kf/core/geo/KfGeoLayer.cxx
deleted file mode 100644
index eaa8e9214b6ca46b71fdbaa166bd7044598fc86e..0000000000000000000000000000000000000000
--- a/algo/kf/core/geo/KfGeoLayer.cxx
+++ /dev/null
@@ -1,30 +0,0 @@
-/* Copyright (C) 2024 GSI Helmholtzzentrum fuer Schwerionenforschung, Darmstadt
- SPDX-License-Identifier: GPL-3.0-only
- Authors: Sergei Zharko [committer] */
-
-/// @file KfGeoLayer.cxx
-/// @brief A logical layer of the KF-setup (implementation)
-/// @since 18.07.2024
-/// @author Sergei Zharko <s.zharko@gsi.de>
-
-#include "KfGeoLayer.h"
-
-using cbm::algo::kf::GeoLayer;
-
-// ---------------------------------------------------------------------------------------------------------------------
-//
-GeoLayer::GeoLayer(kf::MaterialMap&& materialMap) : fMaterialMap(std::move(materialMap)) {}
-
-// ---------------------------------------------------------------------------------------------------------------------
-//
-std::string GeoLayer::ToString(int verbose)
-{
- if (verbose < 1) {
- return std::string{};
- }
-
- // TODO: SZh. 19.07.2024: Implement
- std::stringstream msg;
- msg << "";
- return msg.str();
-}
diff --git a/algo/kf/core/geo/KfGeoLayer.h b/algo/kf/core/geo/KfGeoLayer.h
deleted file mode 100644
index 60bd954410c0eb33585fc6450a5155e1ac4ccd12..0000000000000000000000000000000000000000
--- a/algo/kf/core/geo/KfGeoLayer.h
+++ /dev/null
@@ -1,69 +0,0 @@
-/* Copyright (C) 2024 GSI Helmholtzzentrum fuer Schwerionenforschung, Darmstadt
- SPDX-License-Identifier: GPL-3.0-only
- Authors: Sergei Zharko [committer] */
-
-/// @file KfGeoLayer.h
-/// @brief A logical layer of the KF-setup (header)
-/// @since 18.07.2024
-/// @author Sergei Zharko <s.zharko@gsi.de>
-
-#pragma once
-
-#include "KfMaterialMap.h"
-
-namespace cbm::algo::kf
-{
- /// \class GeoLayer
- /// \brief A logical layer of the KF-setup for a given [zMin, zMax] interval
- ///
- /// A layer includes a material budget map and a magnetic field, estimated for a reference z-position.
- class GeoLayer {
- public:
- /// \brief Default constructor
- GeoLayer() = default;
-
- /// \brief Constructor from the parameters
- /// \param materialMap Instance of the material budget map
- ///
- /// \note The constructor moves the material budget map from the source
- GeoLayer(kf::MaterialMap&& materialMap);
-
- /// \brief Copy constructor
- GeoLayer(const GeoLayer&) = default;
-
- /// \brief Move constructor
- GeoLayer(GeoLayer&&) = default;
-
- /// \brief Copy assignment operator
- GeoLayer& operator=(const GeoLayer&) = default;
-
- /// \brief Move assignment operator
- GeoLayer& operator=(GeoLayer&&) = default;
-
- /// \brief Destructor
- ~GeoLayer() = default;
-
- /// \brief Material map accessor (mutable)
- kf::MaterialMap& MaterialMap() { return fMaterialMap; }
-
- /// \brief Material map accessor (constant)
- const kf::MaterialMap& MaterialMap() const { return fMaterialMap; }
-
- /// \brief String representation of the class
- /// \param verbose Verbosity level
- std::string ToString(int verbose = 1);
-
- private:
- /// \brief Serialization method
- friend class boost::serialization::access;
- template<class Archive>
- void serialize(Archive& ar, const unsigned int /*version*/)
- {
- ar& fMaterialMap;
- }
-
- kf::MaterialMap fMaterialMap; ///< Material map instance for the geometry layer
- };
-
-
-} // namespace cbm::algo::kf
diff --git a/algo/kf/core/geo/KfMaterialMap.cxx b/algo/kf/core/geo/KfMaterialMap.cxx
index 9128811ff379a415d1772a2b84f0401a13acbf3d..e6f69ce164128064de0ad604d36060374053e48d 100644
--- a/algo/kf/core/geo/KfMaterialMap.cxx
+++ b/algo/kf/core/geo/KfMaterialMap.cxx
@@ -15,6 +15,20 @@
using cbm::algo::kf::fvec;
using cbm::algo::kf::MaterialMap;
+// ---------------------------------------------------------------------------------------------------------------------
+//
+MaterialMap::MaterialMap(int nBins, float xyMax, float zRef, float zMin, float zMax)
+ : fNbins(nBins)
+ , fXYmax(xyMax)
+ , fFactor(0.5 * fNbins / fXYmax)
+ , fZref(zRef)
+ , fZmin(zMin)
+ , fZmax(zMax)
+{
+ this->CheckConsistency();
+ fTable.resize(fNbins * fNbins);
+}
+
// ---------------------------------------------------------------------------------------------------------------------
//
MaterialMap::MaterialMap(MaterialMap&& other) noexcept { this->Swap(other); }
@@ -77,39 +91,6 @@ int MaterialMap::GetBin(float x, float y) const
return i + j * fNbins;
}
-// ---------------------------------------------------------------------------------------------------------------------
-//
-void MaterialMap::Initialize(int nBins, float xyMax, float zRef, float zMin, float zMax)
-{
- fNbins = nBins;
- fXYmax = xyMax;
- fZref = zRef;
- fZmin = zMin;
- fZmax = zMax;
-
- if (fNbins < 1) {
- std::stringstream aStream;
- aStream << "MaterialMap: object cannot be initialized with non-positive nBins = " << fNbins;
- throw std::logic_error(aStream.str());
- }
-
- if (fXYmax < 0.) {
- std::stringstream aStream;
- aStream << "MaterialMap: object cannot be initialized with non-positive XYmax = " << fXYmax << " [cm]";
- throw std::logic_error(aStream.str());
- }
-
- if (!((fZmin <= fZref) && (fZref <= zMax))) {
- std::stringstream aStream;
- aStream << "MaterialMap: object cannot be initialized with inconsistent Z: min " << fZmin << " ref " << fZref
- << " max " << fZmax << " [cm]";
- throw std::logic_error(aStream.str());
- }
-
- fFactor = 0.5 * fNbins / fXYmax;
- fTable.resize(fNbins * fNbins);
-}
-
// ---------------------------------------------------------------------------------------------------------------------
//
void MaterialMap::Swap(MaterialMap& other) noexcept
@@ -133,3 +114,27 @@ std::string MaterialMap::ToString() const
<< setw(12) << fZref << ' ' << setw(12) << fZmin << ' ' << setw(12) << fZmax;
return msg.str();
}
+
+// ---------------------------------------------------------------------------------------------------------------------
+//
+void MaterialMap::CheckConsistency() const
+{
+ if (fNbins < 1) {
+ std::stringstream aStream;
+ aStream << "MaterialMap: object cannot be initialized with non-positive nBins = " << fNbins;
+ throw std::logic_error(aStream.str());
+ }
+
+ if (fXYmax < 0.) {
+ std::stringstream aStream;
+ aStream << "MaterialMap: object cannot be initialized with non-positive XYmax = " << fXYmax << " [cm]";
+ throw std::logic_error(aStream.str());
+ }
+
+ if (!((fZmin <= fZref) && (fZref <= fZmax))) {
+ std::stringstream aStream;
+ aStream << "MaterialMap: object cannot be initialized with inconsistent Z: min " << fZmin << " ref " << fZref
+ << " max " << fZmax << " [cm]";
+ throw std::logic_error(aStream.str());
+ }
+}
diff --git a/algo/kf/core/geo/KfMaterialMap.h b/algo/kf/core/geo/KfMaterialMap.h
index d336d84eadbbe602a86b3f95f489d93858990bea..30289c882bdc9628c2bd2d945f86a3f2f550a26f 100644
--- a/algo/kf/core/geo/KfMaterialMap.h
+++ b/algo/kf/core/geo/KfMaterialMap.h
@@ -29,6 +29,14 @@ namespace cbm::algo::kf
/// \brief Default constructor
MaterialMap() = default;
+ /// \brief Constructor from parameters
+ /// \param nBins Number of rows or columns
+ /// \param xyMax Size of station in x and y dimensions [cm]
+ /// \param zRef Reference z-coordinate of the material layer [cm]
+ /// \param zMin Lower boundary z-coordinate for the material layer [cm]
+ /// \param zMax Upper boundary z-coordinate for the material layer [cm]
+ MaterialMap(int nBins, float xyMax, float zRef, float zMin, float zMax);
+
/// \brief Copy constructor
MaterialMap(const MaterialMap& other) = default;
@@ -98,8 +106,8 @@ namespace cbm::algo::kf
}
}
- /// \brief Checks, if the fields are NaN
- bool IsNaN() const
+ /// \brief Function to test the instance for NaN
+ bool IsUndefined() const
{
return utils::IsUndefined(fNbins) || utils::IsUndefined(fXYmax * fFactor * fZref * fZmin * fZmax);
}
@@ -113,11 +121,6 @@ namespace cbm::algo::kf
/// because iBinX = 0 and iBinY = 0 in the TH1::SetBinContent method of usually defines the underflow bin.
void SetRadThickBin(int iBinX, int iBinY, float thickness) { fTable[iBinX + fNbins * iBinY] = thickness; }
- /// \brief Sets properties of the material table -- number of rows or columnts and the size of station in XY plane
- /// \param nBins Number of rows or columns
- /// \param xyMax Size of station in x and y dimensions [cm]
- void Initialize(int nBins, float xyMax, float zRef, float zMin, float zMax);
-
/// \brief Swap method
void Swap(MaterialMap& other) noexcept;
@@ -127,7 +130,14 @@ namespace cbm::algo::kf
/// \brief String representation of the object
std::string ToString() const;
+ /// \brief Comparison operator (material map ordering by fZref)
+ friend bool operator<(const MaterialMap& lhs, const MaterialMap& rhs) { return lhs.fZref < rhs.fZref; }
+
private:
+ /// \brief Checks the object consistency
+ /// \throw std::logic_error If the object is in non-valid mode
+ void CheckConsistency() const;
+
int fNbins = defs::Undef<int>; ///< Number of rows (== N columns) in the material budget table
float fXYmax = defs::Undef<float>; ///< Size of the station in x and y dimensions [cm]
float fFactor = defs::Undef<float>; ///< Util. var. for the conversion of point coordinates to row/column id
@@ -150,5 +160,4 @@ namespace cbm::algo::kf
ar& fTable;
}
};
-
} // namespace cbm::algo::kf
diff --git a/algo/kf/core/geo/KfSetup.cxx b/algo/kf/core/geo/KfSetup.cxx
index bce84eece389e15cbf4ed00c23cd9ba090b48e9b..e2ac6a105fa3ce59d3c4a8a2a3f45e24286c65c6 100644
--- a/algo/kf/core/geo/KfSetup.cxx
+++ b/algo/kf/core/geo/KfSetup.cxx
@@ -12,50 +12,112 @@
#include <sstream>
using cbm::algo::kf::Setup;
+using cbm::algo::kf::SetupFactory;
+
+//
+// Setup methods definition
+//
// ---------------------------------------------------------------------------------------------------------------------
//
-template<typename DataT>
-void Setup<DataT>::CheckConsistency() const
+template<typename T>
+std::string Setup<T>::ToString(int verbosity, int indentLevel) const
{
- std::stringstream errMsg;
-
- // checks go here... the errMsg is to be filled on each checking stage
-
- if (!errMsg.str().empty()) {
- throw std::runtime_error(errMsg.str());
+ std::stringstream msg;
+ if (verbosity > 0) {
+ constexpr char indentCh = '\t';
+ std::string indent(indentLevel, indentCh);
+ msg << indent << "----- KF setup -----\n";
+ msg << indent << "Target:\n" << fTarget.ToString(indentLevel + 1);
+ msg << indent << "Material layers:\n";
+ for (const auto& layer : fvMaterialLayers) {
+ msg << indent << indentCh << layer.ToString() << '\n';
+ }
+ msg << indent << "Interpolated magnetic field:\n" << fIntrplField.ToString(indentLevel + 1);
+ if (foOrigField) {
+ msg << indent << "Original magnetic field:\n" << foOrigField->ToString(indentLevel + 1);
+ }
+ else {
+ msg << indent << "Original magnetic field: not provided";
+ }
}
+ return msg.str();
}
+namespace cbm::algo::kf
+{
+ template class Setup<float>;
+ template class Setup<double>;
+ template class Setup<fvec>;
+} // namespace cbm::algo::kf
+
+
+//
+// SetupFactory methods definition
+//
+
// ---------------------------------------------------------------------------------------------------------------------
//
-template<typename DataT>
-void Setup<DataT>::Init()
+void SetupFactory::AddMaterial(const MaterialMap& material)
{
- fbInitialized = true;
- CheckConsistency();
+ if (!fMaterialLayers.emplace(material).second) {
+ std::stringstream msg;
+ msg << "SetupFactory::AddMaterial: attempt of adding a duplicating material layer " << material.ToString()
+ << ".\nThe next material layers were already added:";
+ for (const auto& el : fMaterialLayers) {
+ msg << "\n\t- " << el.ToString();
+ }
+ throw std::logic_error(msg.str());
+ }
}
// ---------------------------------------------------------------------------------------------------------------------
//
-template<typename DataT>
-std::string Setup<DataT>::ToString(int verbosity, int indentLevel) const
+template<typename T>
+Setup<T> SetupFactory::MakeSetup(bool bProvideOrigField) const
{
- std::stringstream msg;
+ Setup<T> setup;
+ // Target initialization
+ setup.SetTarget(fTarget);
- if (verbosity > 0) {
- constexpr char indentCh = '\t';
- std::string indent(indentLevel, indentCh);
- msg << indent << " ----- KF setup -----\n";
- msg << indent << "MATERIAL LAYERS:\n";
- for (const auto& layer : fvMatLayers) {
- msg << indent << indentCh << layer.ToString() << '\n';
- }
+ // Field initialization
+ setup.SetField(fFieldFactory.MakeField<T, EFieldMode::Intrpl>());
+ if (bProvideOrigField) {
+ setup.SetField(fFieldFactory.MakeField<T, EFieldMode::Orig>());
}
- return msg.str();
+ // Material layers initialization
+ if (fMaterialLayers.size() > defs::MaxNofMaterialLayers) {
+ std::stringstream msg;
+ msg << "kf::SetupFactory::MakeSetup(): too many material layers are provided (" << fMaterialLayers.size()
+ << "), the maximum allowed number of the layers kf::defs::MaxNofMaterialLayers = "
+ << defs::MaxNofMaterialLayers;
+ throw std::logic_error(msg.str());
+ }
+ int layerId = 0;
+ for (const auto& material : fMaterialLayers) {
+ setup.SetMaterial(layerId, material);
+ }
+
+ return setup;
+}
+
+// ---------------------------------------------------------------------------------------------------------------------
+//
+void SetupFactory::Reset()
+{
+ fFieldFactory.Reset();
+ fMaterialLayers.clear();
}
-template class cbm::algo::kf::Setup<float>;
-template class cbm::algo::kf::Setup<double>;
-//template class cbm::algo::kf::Setup<fvec>;
+// ---------------------------------------------------------------------------------------------------------------------
+//
+void SetupFactory::SetTargetProperties(double x, double y, double z, double fieldInitStep, const MaterialMap& material)
+{
+ fTarget.SetX(x);
+ fTarget.SetY(y);
+ fTarget.SetZ(z);
+ fTarget.SetMaterial(material);
+ fFieldFactory.SetTarget(x, y, z);
+ fFieldFactory.SetStep(fieldInitStep);
+}
diff --git a/algo/kf/core/geo/KfSetup.h b/algo/kf/core/geo/KfSetup.h
index 755b02ab5fd97142466fdee938346d6c27a4dfb6..994c056a04b584656f0bfb1e3ca23c601fb761ec 100644
--- a/algo/kf/core/geo/KfSetup.h
+++ b/algo/kf/core/geo/KfSetup.h
@@ -10,24 +10,31 @@
#pragma once // include this header only once per compilation unit
#include "KfDefs.h"
-#include "KfGeoLayer.h"
+#include "KfField.h"
+#include "KfMaterialMap.h"
+#include "KfTarget.h"
#include "KfVector.h"
#include <boost/serialization/access.hpp>
+#include <boost/serialization/split_free.hpp>
+#include <optional>
#include <string>
#include <vector>
namespace cbm::algo::kf
{
- using MaterialContainer_t = Vector<MaterialMap>;
-
/// \class Setup
/// \brief KF-framework representation of the detector setup
- /// \tparam T Underlying floating point data-type (float, double or fvec)
- // TODO: Should the class have a template parameter?
+ /// \tparam T Underlying floating-point data type
template<typename T>
class alignas(VcMemAlign) Setup {
+ friend class boost::serialization::access; // Boost serializer methods
+ template<typename I>
+ friend class Setup;
+
+ using MaterialContainer_t = std::array<MaterialMap, defs::MaxNofMaterialLayers>;
+
public:
/// \brief Default constructor
Setup() = default;
@@ -35,65 +42,209 @@ namespace cbm::algo::kf
/// \brief Destructor
~Setup() = default;
- /// \brief Move constructor
- Setup(Setup&&) = default;
-
- /// \brief Move assignment operator
- Setup& operator=(Setup&&) = default;
-
- /// \brief Adds material layer
- void AddMaterialLayer(const MaterialMap& layer) { fvMatLayers.push_back(layer); }
-
- /// \brief Adds material layer, moving it from the source
- void AddMaterialLayer(MaterialMap&& layer) { fvMatLayers.push_back(std::move(layer)); }
+ /// \brief Copy constructor
+ /// \tparam I Underlying floating-point data type of the source
+ template<typename I>
+ Setup(const Setup<I>& other)
+ : fvMaterialLayers(other.fvMaterialLayers)
+ , foOrigField(other.foOrigField.has_value() ? std::make_optional(*other.foOrigField) : std::nullopt)
+ , fIntrplField(other.fIntrplField)
+ , fTarget(other.fTarget)
+ {
+ }
- /// \brief Check method
- /// \note Throws std::runtime_error, if the class inconsistent
- void CheckConsistency() const;
+ /// \brief Move constructor
+ //Setup(Setup&&) noexcept;
- /// \brief Gets material layer
- [[gnu::always_inline]] const MaterialMap& GetLayer(int iLayer) const
+ /// \brief Copy assignment operator
+ /// \tparam I Underlying floating-point data type of the source
+ template<typename I>
+ Setup& operator=(const Setup<I>& other)
{
- return fvMatLayers.at<defs::GetterCheck>(iLayer);
+ if constexpr (std::is_same_v<I, T>) {
+ if (this == &other) {
+ return *this;
+ }
+ }
+ fvMaterialLayers = other.fvMaterialLayers;
+ foOrigField = other.foOrigField;
+ fIntrplField = other.fIntrplField;
+ fTarget = other.fTarget;
+ return *this;
}
- /// \brief Initializes setup
- void Init();
+ /// \brief Move assignment operator
+ //Setup& operator=(Setup&&) noexcept;
+
+ /// \brief Makes an instance of the field depending on the template parameter
+ /// \tparam FldMode Field mode of the returned instance
+ /// \throw std::logic_error If the particular field member is undefined (nullopt)
+ template<EFieldMode FldMode>
+ const Field<T, FldMode>& GetField() const;
- /// \brief Returns initialization status
- bool IsInitialized() const { return fbInitialized; }
+ /// \brief Gets material layer
+ /// \param iLayer Index of layer
+ const MaterialMap& GetMaterial(int iLayer) const { return fvMaterialLayers[iLayer]; }
- /// \brief Resets initialization
- void Reset() { *this = std::move(Setup{}); }
+ /// \brief Gets target
+ const Target<T>& GetTarget() const { return fTarget; }
/// \brief String representation of the class contents
/// \param verbosity A verbose level for output
/// \param indentLevel Indent level of the string output
std::string ToString(int verbosity = 0, int indentLevel = 0) const;
- private:
- /// \brief Copy constructor
- /// \tparam I Underlying floating-point data type of the source
- //template<typename I>
- Setup(const Setup& other) = default;
+ /// \brief Sets field
+ template<EFieldMode FldMode>
+ void SetField(const Field<double, FldMode>& field);
- /// \brief Copy assignment operator
- /// \tparam I Underlying floating-point data type of the source
- //template<typename I>
- Setup& operator=(const Setup& other) = default;
+ /// \brief Sets material layer
+ /// \param iLayer Index of layer
+ /// \param map Material map
+ void SetMaterial(int iLayer, const MaterialMap& map) { fvMaterialLayers[iLayer] = map; }
+
+ /// \brief Sets target (moves from the source)
+ /// \param target Target instance
+ void SetTarget(const Target<double>& target) { fTarget = target; }
+
+ /// \brief BOOST serialization load method
+ template<class Archive>
+ void load(Archive& ar, const unsigned int /*version*/)
+ {
+ ar >> fvMaterialLayers;
+ ar >> fTarget;
+ ar >> fIntrplField;
+ foOrigField = std::nullopt;
+
+ // NOTE: The original field cannot be serialized, because it contains a pointer to the external magnetic field
+ // function.
+ }
+
+ template<class Archive>
+ void save(Archive& ar, const unsigned int /*version*/) const
+ {
+ ar << fvMaterialLayers;
+ ar << fTarget;
+ ar << fIntrplField;
+ }
- /// \brief Serialization method
- friend class boost::serialization::access;
template<class Archive>
- void serialize(Archive& ar, const unsigned int /*version*/)
+ void serialize(Archive& ar, const unsigned int version)
{
- ar& fvMatLayers;
- ar& fbInitialized;
+ boost::serialization::split_member(ar, *this, version);
}
- MaterialContainer_t fvMatLayers = {}; ///< Material layers map
+ MaterialContainer_t fvMaterialLayers = {}; ///< Container of the material maps
+ Field<T, EFieldMode::Intrpl> fIntrplField; ///< Interpolated field (NOTE: maybe make optional)
+ std::optional<Field<T, EFieldMode::Orig>> foOrigField{std::nullopt}; ///< Optional original field
+ Target<T> fTarget; ///< Target layer
+ };
- bool fbInitialized{false}; ///< Is instance initialized
- bool fbOrigFieldAvailable{false}; ///< Is original field available
+
+ /// \class SetupFactory
+ /// \brief Creates a valid initialized Setup instance
+ class SetupFactory {
+ public:
+ /// \brief Default constructor
+ SetupFactory() = default;
+
+ /// \brief Copy constructor
+ SetupFactory(const SetupFactory&) = delete;
+
+ /// \brief Move constructor
+ SetupFactory(SetupFactory&&) = delete;
+
+ /// \brief Destructor
+ ~SetupFactory() = default;
+
+ /// \brief Copy assignment operator
+ SetupFactory& operator=(const SetupFactory&) = delete;
+
+ /// \brief Move assignment operator
+ SetupFactory& operator=(SetupFactory&&) = delete;
+
+ /// \brief Adds magnetic field slice reference
+ /// \param halfSizeX Half-size of the slice in x-direction [cm]
+ /// \param halfSizeY Half-size of the slice in y-direction [cm]
+ /// \param refZ Reference z-position of the slice [cm]
+ void AddFieldSliceReference(double halfSizeX, double halfSizeY, double refZ)
+ {
+ fFieldFactory.AddSliceReference(halfSizeX, halfSizeY, refZ);
+ }
+
+ /// \brief Adds material layer
+ /// \param material Material layer
+ void AddMaterial(const MaterialMap& material);
+
+ /// \brief Creates a setup instance
+ /// \param bProvideOrigField Flag: true - defines optional variable of the original field
+ template<typename T>
+ Setup<T> MakeSetup(bool bProvideOrigField = true) const;
+
+ /// \brief Resets the instance
+ void Reset();
+
+ /// \brief Sets magnetic field function
+ /// \param fieldFn Magnetic field function (either KF, or FairRoot format)
+ /// \param fieldType Magnetic field type
+ template<typename FieldFn>
+ void SetFieldFunction(const FieldFn& fieldFn, EFieldType fieldType)
+ {
+ fFieldFactory.SetFieldFunction(fieldFn);
+ fFieldFactory.SetFieldType(fieldType);
+ }
+
+ /// \brief Sets target initialization properties
+ /// \param x Target x-coordinate [cm]
+ /// \param y Target y-coordinate [cm]
+ /// \param z Target z-coordinate [cm]
+ /// \param fieldInitStep Step between the nodal points in the target field region initialization [cm]
+ /// \param material Target material layer
+ void SetTargetProperties(double x, double y, double z, double fieldInitStep, const MaterialMap& material);
+
+ private:
+ std::set<MaterialMap> fMaterialLayers{}; ///< Set of material maps
+ FieldFactory fFieldFactory; ///< Instance of field factory
+ Target<double> fTarget; ///< Field to keep target properties
};
+
+ //
+ // Template method definition
+ //
+
+ // -------------------------------------------------------------------------------------------------------------------
+ //
+ template<typename T>
+ template<EFieldMode FldMode>
+ const Field<T, FldMode>& Setup<T>::GetField() const
+ {
+ if constexpr (FldMode == EFieldMode::Orig) {
+ if (foOrigField) {
+ return *foOrigField;
+ }
+ else {
+ throw std::logic_error("kf::Setup::MakeField(): attempt to create a Field instance in the original field mode "
+ "in the unsuitable scenario (for example, in the online reconstruction code). Please, "
+ "use the interpolated magnetic field calling Setup::MakeField<T, EFieldMode::Intrpl>() "
+ "instead.");
+ }
+ }
+ else if constexpr (FldMode == EFieldMode::Intrpl) {
+ return fIntrplField;
+ }
+ }
+
+ // -------------------------------------------------------------------------------------------------------------------
+ //
+ template<typename T>
+ template<EFieldMode FldMode>
+ void Setup<T>::SetField(const Field<double, FldMode>& field)
+ {
+ if constexpr (FldMode == EFieldMode::Intrpl) {
+ fIntrplField = field;
+ }
+ else if constexpr (FldMode == EFieldMode::Orig) {
+ foOrigField = std::make_optional(field);
+ }
+ }
} // namespace cbm::algo::kf
diff --git a/algo/kf/core/geo/KfSetupSerializer.cxx b/algo/kf/core/geo/KfSetupSerializer.cxx
new file mode 100644
index 0000000000000000000000000000000000000000..52cb6658334b1df83ff0a65df39a05e536583ff0
--- /dev/null
+++ b/algo/kf/core/geo/KfSetupSerializer.cxx
@@ -0,0 +1,27 @@
+/* Copyright (C) 2024 GSI Helmholtzzentrum fuer Schwerionenforschung, Darmstadt
+ SPDX-License-Identifier: GPL-3.0-only
+ Authors: Sergei Zharko [committer] */
+
+/// @file KfSetupSerializer.cxx
+/// @brief Helper class for a kf::Setup loading/storing (implementation)
+/// @since 26.08.2024
+/// @author Sergei Zharko <s.zharko@gsi.de>
+
+#include "KfSetupSerializer.h"
+
+using cbm::algo::kf::Setup;
+using cbm::algo::kf::SetupSerializer;
+
+// ---------------------------------------------------------------------------------------------------------------------
+//
+void SetupSerializer::Store(const Setup<double>& setup, const std::string& fileName)
+{
+ std::ofstream ofs(fileName, std::ios::binary);
+ if (!ofs) {
+ std::stringstream msg;
+ msg << "kf::SetupSerializer::Store: failed openning file \"" << fileName << "\" to store the setup";
+ throw std::runtime_error(msg.str());
+ }
+ boost::archive::binary_oarchive oa(ofs);
+ oa << setup;
+}
diff --git a/algo/kf/core/geo/KfSetupSerializer.h b/algo/kf/core/geo/KfSetupSerializer.h
new file mode 100644
index 0000000000000000000000000000000000000000..98af29fc96f7969705600ef6270bea2e45b117aa
--- /dev/null
+++ b/algo/kf/core/geo/KfSetupSerializer.h
@@ -0,0 +1,78 @@
+/* Copyright (C) 2024 GSI Helmholtzzentrum fuer Schwerionenforschung, Darmstadt
+ SPDX-License-Identifier: GPL-3.0-only
+ Authors: Sergei Zharko [committer] */
+
+/// @file KfSetupSerializer.h
+/// @brief Helper class for a kf::Setup loading/storing (header)
+/// @since 26.08.2024
+/// @author Sergei Zharko <s.zharko@gsi.de>
+
+#pragma once
+
+#include "KfSetup.h"
+
+#include <boost/archive/binary_iarchive.hpp>
+#include <boost/archive/binary_oarchive.hpp>
+
+#include <fstream>
+#include <sstream>
+#include <string>
+
+namespace cbm::algo::kf
+{
+ /// \class SetupSerializer
+ /// \brief A helper class for storing/loading a kf::Setup instance
+ ///
+ /// A particular setup instance is dependent from a template parameter, which represents an underlying floating
+ /// point data-type. Only a Setup with T = double can be serialized, so if a particular setup instance has
+ /// another floating-point type, it will be converted to the double-precision version.
+ /// Moreover, a Setup instance has two magnetic field representation, an interpolated magnetic field and an
+ /// optional original field. In the serialization only the interpolated magnetic field variant is stored, the
+ /// original field version will be set to nullopt during the loading from the file.
+ ///
+ class SetupSerializer {
+ public:
+ /// \brief Stores a serialized setup to a file
+ /// \param setup A reference to the Setup (NOTE: a double-Setup must be passed explicitly)
+ /// \param fileName Output file name
+ static void Store(const Setup<double>& setup, const std::string& fileName);
+
+ /// \brief Loads a serialized setup from a file
+ /// \tparam T Underlying floating-point type for the output setup object
+ /// \param fileName Input file name
+ /// \throw std::runtime_error If the fileName cannot be opened
+ template<typename T>
+ static Setup<T> Load(const std::string& fileName);
+ };
+
+
+ // -------------------------------------------------------------------------------------------------------------------
+ //
+ template<typename T>
+ Setup<T> SetupSerializer::Load(const std::string& fileName)
+ {
+ Setup<double> setup;
+
+ std::ifstream ifs(fileName, std::ios::binary);
+ if (!ifs) {
+ std::stringstream msg;
+ msg << "kf::SetupSerializer::Load: intput setup file \"" << fileName << "\" was not found";
+ throw std::runtime_error(msg.str());
+ }
+
+ try {
+ boost::archive::binary_iarchive ia(ifs);
+ ia >> setup;
+ }
+ catch (const std::exception& err) {
+ std::stringstream msg;
+ msg << "kf::SetupSerializer::Load: input setup file \"" << fileName
+ << "\" has inconsistent format or was "
+ "corrupted. The exception message: "
+ << err.what();
+ throw std::runtime_error(msg.str());
+ }
+
+ return Setup<T>(setup);
+ }
+} // namespace cbm::algo::kf
diff --git a/algo/kf/core/geo/KfTarget.cxx b/algo/kf/core/geo/KfTarget.cxx
index cb564fbd11f854b8bac9c1afab76f790ab676252..1df0856d6dee7072171f69d1cc9b84afad73e81e 100644
--- a/algo/kf/core/geo/KfTarget.cxx
+++ b/algo/kf/core/geo/KfTarget.cxx
@@ -9,16 +9,50 @@
#include "KfTarget.h"
-using cbm::algo::kf::GeoLayer;
+#include <sstream>
+
using cbm::algo::kf::MaterialMap;
using cbm::algo::kf::Target;
// ---------------------------------------------------------------------------------------------------------------------
//
-Target::Target(kf::MaterialMap&& materialMap, float x, float y, float z)
- : GeoLayer(std::move(materialMap))
- , fX(x)
- , fY(y)
- , fZ(z)
+template<typename T>
+Target<T>::Target(const MaterialMap& material, double x, double y, double z)
+ : fMaterial(material)
+ , fX(utils::simd::Cast<double, T>(x))
+ , fY(utils::simd::Cast<double, T>(y))
+ , fZ(utils::simd::Cast<double, T>(z))
+{
+}
+
+// ---------------------------------------------------------------------------------------------------------------------
+//
+template<typename T>
+void Target<T>::SetMaterial(const MaterialMap& material)
+{
+ if (material.IsUndefined()) {
+ throw std::logic_error("Target:ReceiveMaterial(): attempt to pass an undefined instance of the material map");
+ }
+ fMaterial = material;
+}
+
+// ---------------------------------------------------------------------------------------------------------------------
+//
+template<typename T>
+std::string Target<T>::ToString(int indentLevel) const
{
+ constexpr char IndentChar = '\t';
+ std::stringstream msg;
+ std::string indent(indentLevel, IndentChar);
+ msg << indent << "position: {" << fX << ", " << fY << ", " << fZ << "} [cm]\n";
+ msg << indent << "material: " << fMaterial.ToString();
+ return msg.str();
}
+
+
+namespace cbm::algo::kf
+{
+ template class Target<float>;
+ template class Target<double>;
+ template class Target<fvec>;
+} // namespace cbm::algo::kf
diff --git a/algo/kf/core/geo/KfTarget.h b/algo/kf/core/geo/KfTarget.h
index f37bca6eac969c53db01ba6fb5b2c4c3b722030d..78008cc4501d81e6503ddf2c45d82dcccfb2c97e 100644
--- a/algo/kf/core/geo/KfTarget.h
+++ b/algo/kf/core/geo/KfTarget.h
@@ -10,62 +10,95 @@
#pragma once
#include "KfDefs.h"
-#include "KfGeoLayer.h"
+#include "KfMaterialMap.h"
#include <boost/serialization/access.hpp>
-#include <boost/serialization/base_object.hpp>
+
+#include <string>
namespace cbm::algo::kf
{
- /// \class Target
- /// \brief A geometry layer in the target region
- class Target : public GeoLayer {
+ /// \class Target
+ /// \brief A geometry layer in the target region
+ /// \tparam T Underlying floating-point type
+ template<typename T>
+ class alignas(VcMemAlign) Target {
+ template<typename I>
+ friend class Target;
+
public:
/// \brief Default constructor
Target() = default;
/// \brief Constructor from parameters
- /// \param materialMap Instance of the material budget map
- /// \param x x-coordinate of the nominal target center [cm]
- /// \param y y-coordinate of the nominal target center [cm]
- /// \param z z-coordinate of the nominal target center [cm]
- Target(kf::MaterialMap&& materialMap, float x, float y, float z);
+ /// \param material Instance of the material budget map near the target region
+ /// \param x x-coordinate of the nominal target center [cm]
+ /// \param y y-coordinate of the nominal target center [cm]
+ /// \param z z-coordinate of the nominal target center [cm]
+ Target(const MaterialMap& material, double x, double y, double z);
/// \brief Copy constructor
- Target(const Target&) = default;
-
- /// \brief Move constructor
- Target(Target&&) = default;
-
- /// \brief Copy assignment operator
- Target& operator=(const Target&) = default;
-
- /// \brief Move assignment operator
- Target& operator=(Target&&) = default;
+ /// \tparam I Underlying floating-point type of the source
+ template<typename I>
+ Target(const Target<I>& other)
+ : fMaterial(other.fMaterial)
+ , fX(utils::simd::Cast<I, T>(other.fX))
+ , fY(utils::simd::Cast<I, T>(other.fY))
+ , fZ(utils::simd::Cast<I, T>(other.fZ))
+ {
+ }
/// \brief Destructor
~Target() = default;
+ /// \brief Copy assignment operator
+ /// \tparam I Underlying floating-point type of the source
+ template<typename I>
+ Target& operator=(const Target<I>& other)
+ {
+ if constexpr (std::is_same_v<I, T>) {
+ if (this == &other) {
+ return *this;
+ }
+ }
+ fMaterial = other.fMaterial;
+ fX = utils::simd::Cast<I, T>(other.fX);
+ fY = utils::simd::Cast<I, T>(other.fY);
+ fZ = utils::simd::Cast<I, T>(other.fZ);
+ return *this;
+ }
+
/// \brief Gets x-coordinate of the nominal target center
- float GetX() const { return fX; }
+ const T& GetX() const { return fX; }
/// \brief Gets x-coordinate of the nominal target center
- float GetY() const { return fY; }
+ const T& GetY() const { return fY; }
/// \brief Gets x-coordinate of the nominal target center
- float GetZ() const { return fZ; }
+ const T& GetZ() const { return fZ; }
+
+ /// \brief Gets material map
+ const MaterialMap& GetMaterial() const { return fMaterial; }
+
+ /// \brief Sets material map
+ /// \param material Material map
+ void SetMaterial(const MaterialMap& material);
/// \brief Sets x-coordinate of the nominal target center
/// \param x x-coordinate [cm]
- void SetX(float x) { fX = x; }
+ void SetX(const T& x) { fX = x; }
/// \brief Sets y-coordinate of the nominal target center
/// \param y y-coordinate [cm]
- void SetY(float y) { fY = y; }
+ void SetY(const T& y) { fY = y; }
/// \brief Sets x-coordinate of the nominal target center
/// \param z x-coordinate [cm]
- void SetZ(float z) { fZ = z; }
+ void SetZ(const T& z) { fZ = z; }
+
+ /// \brief String representation of the class
+ /// \param indentLevel Indent level of the string output
+ std::string ToString(int indentLevel = 0) const;
private:
/// \brief Serialization method
@@ -73,14 +106,15 @@ namespace cbm::algo::kf
template<class Archive>
void serialize(Archive& ar, const unsigned int /*version*/)
{
- ar& boost::serialization::base_object<GeoLayer>(*this);
+ ar& fMaterial;
ar& fX;
ar& fY;
ar& fZ;
}
- float fX = defs::Undef<float>; ///< x-coordinate of the nominal target center [cm]
- float fY = defs::Undef<float>; ///< y-coordinate of the nominal target center [cm]
- float fZ = defs::Undef<float>; ///< z-coordinate of the nominal target center [cm]
+ MaterialMap fMaterial; ///< Material map in the target region
+ T fX = defs::Undef<T>; ///< x-coordinate of the nominal target center [cm]
+ T fY = defs::Undef<T>; ///< y-coordinate of the nominal target center [cm]
+ T fZ = defs::Undef<T>; ///< z-coordinate of the nominal target center [cm]
};
} // namespace cbm::algo::kf
diff --git a/algo/kf/core/pars/KfParameters.cxx b/algo/kf/core/pars/KfParameters.cxx
index 999aec0b9a62980baa577b4ec3d76a9d75f6b928..ce76c88c748f8c900b3819ef64bdff40d1d3bb97 100644
--- a/algo/kf/core/pars/KfParameters.cxx
+++ b/algo/kf/core/pars/KfParameters.cxx
@@ -54,6 +54,9 @@ std::string Parameters<DataT>::ToString(int verbosity, int indentLevel) const
return msg.str();
}
-template class cbm::algo::kf::Parameters<float>;
-template class cbm::algo::kf::Parameters<double>;
-// template class cbm::algo::kf::Parameters<fvec>;
+namespace cbm::algo::kf
+{
+ template class Parameters<float>;
+ template class Parameters<double>;
+ template class Parameters<fvec>;
+} // namespace cbm::algo::kf