From 35345d22417fdec1a51e70c8ff86049630b7b7bd Mon Sep 17 00:00:00 2001
From: "s.zharko@gsi.de" <s.zharko@gsi.de>
Date: Fri, 1 Nov 2024 13:34:40 +0100
Subject: [PATCH] CA: moving the material monitor to the KF-framework
---
algo/ca/core/CMakeLists.txt | 2 -
algo/ca/core/pars/CaMaterialMonitor.h | 112 ------------------
algo/kf/core/CMakeLists.txt | 2 +
algo/kf/core/geo/KfMaterialMap.h | 25 ++--
.../core/geo/KfMaterialMonitor.cxx} | 48 ++++----
algo/kf/core/geo/KfMaterialMonitor.h | 112 ++++++++++++++++++
reco/L1/CbmL1.cxx | 49 ++++----
reco/L1/CbmL1.h | 5 +-
8 files changed, 189 insertions(+), 166 deletions(-)
delete mode 100644 algo/ca/core/pars/CaMaterialMonitor.h
rename algo/{ca/core/pars/CaMaterialMonitor.cxx => kf/core/geo/KfMaterialMonitor.cxx} (70%)
create mode 100644 algo/kf/core/geo/KfMaterialMonitor.h
diff --git a/algo/ca/core/CMakeLists.txt b/algo/ca/core/CMakeLists.txt
index 3d512ff818..e54aeda076 100644
--- a/algo/ca/core/CMakeLists.txt
+++ b/algo/ca/core/CMakeLists.txt
@@ -20,7 +20,6 @@ set(SRCS
${CMAKE_CURRENT_SOURCE_DIR}/pars/CaInitManager.cxx
${CMAKE_CURRENT_SOURCE_DIR}/pars/CaIteration.cxx
${CMAKE_CURRENT_SOURCE_DIR}/pars/CaMaterialMap.cxx
- ${CMAKE_CURRENT_SOURCE_DIR}/pars/CaMaterialMonitor.cxx
${CMAKE_CURRENT_SOURCE_DIR}/pars/CaParameters.cxx
${CMAKE_CURRENT_SOURCE_DIR}/pars/CaSearchWindow.cxx
${CMAKE_CURRENT_SOURCE_DIR}/pars/CaStation.cxx
@@ -114,7 +113,6 @@ install(
pars/CaInitManager.h
pars/CaIteration.h
pars/CaMaterialMap.h
- pars/CaMaterialMonitor.h
pars/CaParameters.h
pars/CaSearchWindow.h
pars/CaStation.h
diff --git a/algo/ca/core/pars/CaMaterialMonitor.h b/algo/ca/core/pars/CaMaterialMonitor.h
deleted file mode 100644
index c220eef13d..0000000000
--- a/algo/ca/core/pars/CaMaterialMonitor.h
+++ /dev/null
@@ -1,112 +0,0 @@
-/* Copyright (C) 2023 GSI Helmholtzzentrum fuer Schwerionenforschung, Darmstadt
- SPDX-License-Identifier: GPL-3.0-only
- Authors: Sergey Gorbunov [committer] */
-
-/// \file MaterialMonitor.h
-/// \brief Class to collect statistics for ca::MaterialMap
-
-
-#pragma once // include this header only once per compilation unit
-
-
-#include "CaDefs.h"
-#include "CaMaterialMap.h"
-
-#include <string>
-#include <vector>
-
-namespace cbm::algo::ca
-{
-
- /// Class to collect statistics for ca::MaterialMap
- ///
- class MaterialMonitor {
- public:
- /// default constructor
- MaterialMonitor() : MaterialMonitor(nullptr) {}
-
- /// constructor
- /// \param materialMap Matareial map to be monitored
- /// \param name Name of the material map
- MaterialMonitor(const ca::MaterialMap* materialMap, std::string name = "") : fName(name)
- {
- SetMaterial(materialMap);
- }
-
- /// construction
- /// \param materialMap Matareial map to be monitored
- void SetMaterial(const ca::MaterialMap* materialMap);
-
- /// construction
- /// \param name Name of the material map
- void SetName(std::string name) { fName = name; }
-
- /// reset the map of active bins
- void ResetActiveBins()
- {
- for (auto& v : fActiveBinMap) {
- v = 0;
- }
- }
-
- /// mark a bin as active
- void MarkActiveBin(float x, float y);
-
- /// update values of statistical variables with respect to the active map
- void EvaluateStatistics();
-
- /// print statistics to a string
- std::string ToString();
-
- /// get number of active bins in the map
- int GetActiveNbins() const { return fActiveNbins; }
-
- /// get minimal radiation thickness among all active bins
- double GetActiveRadThickMin() const { return fPassiveRadThickMin; }
-
- /// get maximal radiation thickness among all active bins
- double GetActiveRadThickMax() const { return fPassiveRadThickMax; }
-
- /// get average radiation thickness among all active bins
- double GetActiveRadThickMean() const { return fPassiveRadThickMean; }
-
-
- /// get number of passive bins in the map
- int GetPassiveNbins() const { return fPassiveNbins; }
-
- /// get minimal radiation thickness among all passive bins
- double GetPassiveRadThickMin() const { return fPassiveRadThickMin; }
-
- /// get maximal radiation thickness among all passive bins
- double GetPassiveRadThickMax() const { return fPassiveRadThickMax; }
-
- /// get average radiation thickness among all passive bins
- double GetPassiveRadThickMean() const { return fPassiveRadThickMean; }
-
- /// get the ration of hits that show up outside the material map
- double GetRatioOfOutsideHits() const { return fNhitsOutside / (fNhitsTotal + 1.e-8); }
-
- /// get the number of processed hits
- double GetNhits() const { return fNhitsTotal; }
-
- private:
- const ca::MaterialMap* fMaterial{nullptr}; ///< Pointer to ca::MaterialMap
- std::string fName{}; ///< Name of the material map
-
- std::vector<char> fActiveBinMap{}; ///< Map of active bins in the material map (bins where hits appear)
-
- int fActiveNbins{constants::Undef<int>}; ///< Active material: number of bins
- double fActiveRadThickMin{constants::Undef<double>}; ///< Active material: minimal thickness
- double fActiveRadThickMax{constants::Undef<double>}; ///< Active material: maximal thickness
- double fActiveRadThickMean{constants::Undef<double>}; ///< Active material: average thickness
-
- int fPassiveNbins{constants::Undef<int>}; ///< Passive material: number of bins
- double fPassiveRadThickMin{constants::Undef<double>}; ///< Passive material: minimal thickness
- double fPassiveRadThickMax{constants::Undef<double>}; ///< Passive material: maximal thickness
- double fPassiveRadThickMean{constants::Undef<double>}; ///< Passive material: average thickness
-
- unsigned long fNhitsTotal{0}; ///< number of hits in statistics
- unsigned long fNhitsOutside{0}; ///< number of hits outside the material map
- };
-
-} // namespace cbm::algo::ca
diff --git a/algo/kf/core/CMakeLists.txt b/algo/kf/core/CMakeLists.txt
index a39e049f2f..0c2168e867 100644
--- a/algo/kf/core/CMakeLists.txt
+++ b/algo/kf/core/CMakeLists.txt
@@ -15,6 +15,7 @@ set(SRCS
${CMAKE_CURRENT_SOURCE_DIR}/data/KfMeasurementXy.cxx
${CMAKE_CURRENT_SOURCE_DIR}/data/KfMeasurementTime.cxx
${CMAKE_CURRENT_SOURCE_DIR}/geo/KfMaterialMap.cxx
+ ${CMAKE_CURRENT_SOURCE_DIR}/geo/KfMaterialMonitor.cxx
${CMAKE_CURRENT_SOURCE_DIR}/geo/KfTarget.cxx
${CMAKE_CURRENT_SOURCE_DIR}/geo/KfField.cxx
${CMAKE_CURRENT_SOURCE_DIR}/geo/KfFieldValue.cxx
@@ -109,6 +110,7 @@ install(
geo/KfFieldSlice.h
geo/KfFieldValue.h
geo/KfMaterialMap.h
+ geo/KfMaterialMonitor.h
geo/KfModuleIndexMap.h
geo/KfSetup.h
geo/KfSetupBuilder.h
diff --git a/algo/kf/core/geo/KfMaterialMap.h b/algo/kf/core/geo/KfMaterialMap.h
index fee65fe057..af1e742588 100644
--- a/algo/kf/core/geo/KfMaterialMap.h
+++ b/algo/kf/core/geo/KfMaterialMap.h
@@ -57,6 +57,9 @@ namespace cbm::algo::kf
/// \param zTarg z-coordinate of the target
void Add(const MaterialMap& other, float zTarg = defs::Undef<float>);
+ /// \brief Gets bin index for (x,y). Returns -1 when outside of the map
+ int GetBin(float x, float y) const;
+
/// \brief Gets number of bins (rows or columns) of the material table
int GetNbins() const { return fNbins; }
@@ -97,26 +100,30 @@ namespace cbm::algo::kf
}
}
-
/// \brief Gets material thickness in units of radiational length X0
- /// \tparam I Type of the x and y (floating point)
- /// \param iX Bin number along x axis
- /// \param iY Bin number along y axis
+ /// \tparam I Type of material thickness
+ /// \param iGlob Global bin number: iGlob = iX + iY * fNbins
template<typename I>
- I GetBinThicknessX0(int iX, int iY) const
- {
+ I GetBinThicknessX0(int iGlob) const {
if constexpr (std::is_same_v<I, fvec>) {
fvec res;
for (size_t i = 0; i < utils::simd::Size<I>(); ++i) {
- res[i] = GetBinThicknessX0<fscal>(iX, iY);
+ res[i] = GetBinThicknessX0<fscal>(iGlob);
}
return res;
}
else {
- return fTable[iX + iY * fNbins];
+ return fTable[iGlob];
}
}
+ /// \brief Gets material thickness in units of radiational length X0
+ /// \tparam I Type of material thickness
+ /// \param iX Bin number along x axis
+ /// \param iY Bin number along y axis
+ template<typename I>
+ I GetBinThicknessX0(int iX, int iY) const { return GetBinThicknessX0<I>(iX + iY * fNbins); }
+
/// \brief Function to test the instance for NaN
bool IsUndefined() const
{
@@ -152,8 +159,6 @@ namespace cbm::algo::kf
/// \throw std::logic_error If the object is in non-valid mode
void CheckConsistency() const;
- /// \brief Get bin index for (x,y). Returns -1 when outside of the map
- int GetBin(float x, float y) 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]
diff --git a/algo/ca/core/pars/CaMaterialMonitor.cxx b/algo/kf/core/geo/KfMaterialMonitor.cxx
similarity index 70%
rename from algo/ca/core/pars/CaMaterialMonitor.cxx
rename to algo/kf/core/geo/KfMaterialMonitor.cxx
index 1e9d02bd0e..57e41cdbec 100644
--- a/algo/ca/core/pars/CaMaterialMonitor.cxx
+++ b/algo/kf/core/geo/KfMaterialMonitor.cxx
@@ -1,30 +1,29 @@
-/* Copyright (C) 2023 GSI Helmholtzzentrum fuer Schwerionenforschung, Darmstadt
+/* Copyright (C) 2023-2024 GSI Helmholtzzentrum fuer Schwerionenforschung, Darmstadt
SPDX-License-Identifier: GPL-3.0-only
- Authors: Sergey Gorbunov [committer] */
+ Authors: Sergey Gorbunov [committer], Sergei Zharko */
-/// \file MaterialMonitor.cxx
-/// \brief Implementation of the MaterialMonitor class
+/// \file MaterialMonitor.cxx
+/// \brief Implementation of the MaterialMonitor class
+/// \author Sergey Gorbunov
-#include "CaMaterialMonitor.h"
+#include "KfMaterialMonitor.h"
#include "AlgoFairloggerCompat.h"
#include <iomanip>
#include <sstream>
-namespace cbm::algo::ca
+namespace cbm::algo::kf
{
-
namespace
{
using namespace cbm::algo;
}
- void MaterialMonitor::SetMaterial(const ca::MaterialMap* materialMap)
+ // -------------------------------------------------------------------------------------------------------------------
+ //
+ void MaterialMonitor::SetMaterial(const MaterialMap* materialMap)
{
- /// construction
- /// \param materialMap Matareial map to be monitored
-
fMaterial = materialMap;
fActiveBinMap.resize(0);
@@ -33,12 +32,15 @@ namespace cbm::algo::ca
fActiveBinMap.assign(fMaterial->GetNbins() * fMaterial->GetNbins(), 0);
}
+ // TODO: What do the hits mean here?
fNhitsTotal = 0;
fNhitsOutside = 0;
EvaluateStatistics();
}
+ // -------------------------------------------------------------------------------------------------------------------
+ //
void MaterialMonitor::MarkActiveBin(float x, float y)
{
/// mark a bin as active
@@ -56,6 +58,8 @@ namespace cbm::algo::ca
}
}
+ // -------------------------------------------------------------------------------------------------------------------
+ //
void MaterialMonitor::EvaluateStatistics()
{
/// update values of statistical variables with respect to the active map
@@ -63,13 +67,13 @@ namespace cbm::algo::ca
fActiveNbins = 0;
fPassiveNbins = 0;
- fActiveRadThickMin = constants::Undef<double>;
- fActiveRadThickMax = constants::Undef<double>;
- fActiveRadThickMean = constants::Undef<double>;
+ fActiveRadThickMin = defs::Undef<double>;
+ fActiveRadThickMax = defs::Undef<double>;
+ fActiveRadThickMean = defs::Undef<double>;
- fPassiveRadThickMin = constants::Undef<double>;
- fPassiveRadThickMax = constants::Undef<double>;
- fPassiveRadThickMean = constants::Undef<double>;
+ fPassiveRadThickMin = defs::Undef<double>;
+ fPassiveRadThickMax = defs::Undef<double>;
+ fPassiveRadThickMean = defs::Undef<double>;
if (!fMaterial) {
return;
@@ -78,13 +82,13 @@ namespace cbm::algo::ca
int nBins = fMaterial->GetNbins() * fMaterial->GetNbins();
if (nBins != (int) fActiveBinMap.size()) {
- LOG(fatal) << "MaterialMonitor: map of active bins is not consistent with the mterial map: nbins "
+ LOG(fatal) << "MaterialMonitor: map of active bins is not consistent with the material map: nbins "
<< fActiveBinMap.size() << " != " << nBins;
return;
}
for (int i = 0; i < nBins; i++) {
- double r = fMaterial->GetRadThickBin(i);
+ double r = fMaterial->GetBinThicknessX0<double>(i);
if (fActiveBinMap[i]) { // active material
if (fActiveNbins == 0) {
fActiveRadThickMin = r;
@@ -125,11 +129,15 @@ namespace cbm::algo::ca
}
+ // -------------------------------------------------------------------------------------------------------------------
+ //
std::string MaterialMonitor::ToString()
{
/// print statistics to a string
- EvaluateStatistics();
+ EvaluateStatistics();
+ // FIXME: a ToString method should not change data of its class, but only represent their current state. I would
+ // call EvaluateStatistics explicitly and provide a const ToString method inside it.
std::stringstream ss;
ss << std::setprecision(2) << std::fixed;
diff --git a/algo/kf/core/geo/KfMaterialMonitor.h b/algo/kf/core/geo/KfMaterialMonitor.h
new file mode 100644
index 0000000000..13d1f351d8
--- /dev/null
+++ b/algo/kf/core/geo/KfMaterialMonitor.h
@@ -0,0 +1,112 @@
+/* Copyright (C) 2023-2024 GSI Helmholtzzentrum fuer Schwerionenforschung, Darmstadt
+ SPDX-License-Identifier: GPL-3.0-only
+ Authors: Sergey Gorbunov [committer], Sergei Zharko */
+
+/// \file KfMaterialMonitor.h
+/// \brief A class to collect statistics for kf::MaterialMap
+/// \author Sergey Gorbunov
+
+#pragma once
+
+#include "KfMaterialMap.h"
+#include "KfDefs.h"
+
+#include <string>
+#include <vector>
+#include <algorithm>
+
+namespace cbm::algo::kf
+{
+ /// \class MaterialMonitor
+ /// \brief A class to collect statistics for a material budget map of the KF-framework
+ ///
+ class MaterialMonitor {
+ public:
+ /// \brief Default constructor
+ MaterialMonitor() : MaterialMonitor(nullptr) {}
+
+ /// \brief Constructor
+ /// \param materialMap Material map to be monitored
+ /// \param name Name of the material map
+ MaterialMonitor(const MaterialMap* materialMap, const std::string& name = "") : fName(name)
+ {
+ SetMaterial(materialMap);
+ }
+
+ /// \brief Set a material budget map
+ /// \param materialMap Material budget map to be monitored
+ void SetMaterial(const MaterialMap* materialMap);
+
+ /// \brief Set a name of the material budget map
+ /// \param name Name of the material map
+ void SetName(const std::string& name) { fName = name; }
+
+ /// \brief Reset the map of active bins
+ void ResetActiveBins()
+ {
+ std::fill(fActiveBinMap.begin(), fActiveBinMap.end(), 0);
+ }
+
+ /// \brief Mark a bin as active
+ /// \param x X-coordinate of the bin
+ /// \param y Y-coordinate of the bin
+ void MarkActiveBin(float x, float y);
+
+ /// \brief Update values of statistical variables with respect to the active map
+ void EvaluateStatistics();
+
+ /// \brief Print statistics to a string
+ std::string ToString();
+
+ /// \brief Get number of active bins in the map
+ int GetActiveNbins() const { return fActiveNbins; }
+
+ /// \brief Get minimal radiation thickness among all active bins
+ double GetActiveRadThickMin() const { return fPassiveRadThickMin; }
+
+ /// \brief Get maximal radiation thickness among all active bins
+ double GetActiveRadThickMax() const { return fPassiveRadThickMax; }
+
+ /// \brief Get average radiation thickness among all active bins
+ double GetActiveRadThickMean() const { return fPassiveRadThickMean; }
+
+ /// \brief Get number of passive bins in the map
+ int GetPassiveNbins() const { return fPassiveNbins; }
+
+ /// \brief Get minimal radiation thickness among all passive bins
+ double GetPassiveRadThickMin() const { return fPassiveRadThickMin; }
+
+ /// \brief Get maximal radiation thickness among all passive bins
+ double GetPassiveRadThickMax() const { return fPassiveRadThickMax; }
+
+ /// \brief Get average radiation thickness among all passive bins
+ double GetPassiveRadThickMean() const { return fPassiveRadThickMean; }
+
+ /// \brief Get the ration of hits that show up outside the material map
+ double GetRatioOfOutsideHits() const { return fNhitsOutside / (fNhitsTotal + 1.e-8); }
+
+ /// \brief Get the number of processed hits
+ double GetNhits() const { return fNhitsTotal; }
+
+ private:
+ std::string fName{}; ///< Name of the material map
+ std::vector<char> fActiveBinMap{}; ///< Map of active bins in the material map (bins where hits appear)
+
+ const MaterialMap* fMaterial{nullptr}; ///< Pointer to the material map
+
+ double fActiveRadThickMin{defs::Undef<double>}; ///< Active material: minimal thickness
+ double fActiveRadThickMax{defs::Undef<double>}; ///< Active material: maximal thickness
+ double fActiveRadThickMean{defs::Undef<double>}; ///< Active material: average thickness
+
+ double fPassiveRadThickMin{defs::Undef<double>}; ///< Passive material: minimal thickness
+ double fPassiveRadThickMax{defs::Undef<double>}; ///< Passive material: maximal thickness
+ double fPassiveRadThickMean{defs::Undef<double>}; ///< Passive material: average thickness
+
+ unsigned long fNhitsTotal{0}; ///< number of hits in statistics
+ unsigned long fNhitsOutside{0}; ///< number of hits outside the material map
+
+ int fActiveNbins{defs::Undef<int>}; ///< Active material: number of bins
+ int fPassiveNbins{defs::Undef<int>}; ///< Passive material: number of bins
+ };
+
+} // namespace cbm::algo::ca
diff --git a/reco/L1/CbmL1.cxx b/reco/L1/CbmL1.cxx
index 710bd4df91..e106cd4881 100644
--- a/reco/L1/CbmL1.cxx
+++ b/reco/L1/CbmL1.cxx
@@ -552,16 +552,17 @@ try {
fNStations = pParameters->GetNstationsActive();
- // monitor the material
+ // TODO: Probably move from CbmL1 to a more proper place (e.g. kf::Setup)
{
- LOG(info) << "\033[31;1m-------------------- L1 material -----------------------------\033[0m";
- fMaterialMonitor.clear();
- for (int i = 0; i < fNStations; i++) {
- ca::MaterialMonitor m(&(pParameters->GetThicknessMaps()[i]), Form("station %d", i));
- LOG(info) << m.ToString();
- fMaterialMonitor.push_back(m);
+ //std::stringstream msg;
+ //msg << " ----- Material budget map monitoring: active setup -----\n";
+ const auto& activeSetup = pParameters->GetActiveSetup();
+ for (int iSt = 0; iSt < activeSetup.GetNofLayers(); ++iSt) {
+ fMaterialMonitor.emplace_back(&(activeSetup.GetMaterial(iSt)), Form("Station %d", iSt));
+ //msg << monitor.ToString() << '\n';
}
- LOG(info) << "\033[31;1m-------------------------------------------------------------\033[0m";
+ //msg << " --------------------------------------------------------";
+ //LOG(info) << '\n' << msg.str();
}
DumpMaterialToFile("L1material.root");
@@ -683,20 +684,28 @@ void CbmL1::Finish()
EfficienciesPerformance(kTRUE);
}
- // monitor the material
- LOG(info) << "\033[31;1m ***************************\033[0m";
- LOG(info) << "\033[31;1m ** CA Tracking monitor **\033[0m";
- LOG(info) << "\033[31;1m ***************************\033[0m";
-
+ {
+ // monitor the material
+ std::stringstream msg;
+ msg << '\n';
+ msg << "\033[31;1m ***************************\033[0m\n";
+ msg << "\033[31;1m ** CA Tracking monitor **\033[0m\n";
+ msg << "\033[31;1m ***************************\033[0m\n";
+
+ // TODO: Probably move from CbmL1 to a more proper place (e.g. kf::Setup)
+ {
+ msg << '\n';
+ msg << " ----- Material budget map monitoring: active setup -----\n";
+ for (auto& monitor: fMaterialMonitor) {
+ msg << monitor.ToString() << '\n';
+ }
+ msg << " --------------------------------------------------------\n";
+ }
- LOG(info) << "\033[31;1m ----- Material budget -------------------------------------- \033[0m";
- for (int i = 0; i < fNStations; i++) {
- LOG(info) << fMaterialMonitor[i].ToString();
+ // monitor of the reconstructed tracks
+ msg << '\n' << fMonitor.ToString();
+ LOG(info) << msg.str();
}
- LOG(info) << "\033[31;1m -------------------------------------------------------------\033[0m";
-
- // monitor of the reconstructed tracks
- LOG(info) << fMonitor.ToString();
TDirectory* curr = gDirectory;
TFile* currentFile = gFile;
diff --git a/reco/L1/CbmL1.h b/reco/L1/CbmL1.h
index 2e78fc97df..fccff6814b 100644
--- a/reco/L1/CbmL1.h
+++ b/reco/L1/CbmL1.h
@@ -24,7 +24,7 @@
#include "CaDataManager.h"
#include "CaFramework.h"
#include "CaInitManager.h"
-#include "CaMaterialMonitor.h"
+#include "KfMaterialMonitor.h"
#include "CaMonitor.h"
#include "CaVector.h"
#include "CbmCaMCModule.h"
@@ -65,6 +65,7 @@ class TGeoNode;
namespace
{
namespace ca = cbm::algo::ca;
+ namespace kf = cbm::algo::kf;
}
/// Internal structure to handle link keys
@@ -533,7 +534,7 @@ class CbmL1 : public FairTask {
bool fExtrapolateToTheEndOfSTS{false};
- std::vector<ca::MaterialMonitor> fMaterialMonitor{}; ///< material monitoring
+ std::vector<kf::MaterialMonitor> fMaterialMonitor{}; ///< Material monitors for each material budget map
ClassDef(CbmL1, 0);
};
--
GitLab