From 018b2d30293977878dbc45fc64cb626ebf5adf20 Mon Sep 17 00:00:00 2001
From: "s.zharko@gsi.de" <s.zharko@gsi.de>
Date: Thu, 9 Jun 2022 17:30:17 +0200
Subject: [PATCH] L1: add consistency checkers for L1Parameters, L1Station,
L1UMeasurementInfo, L1XYMeasurementInfo, L1MaterialInfo, L1Material,
L1FieldSlice, L1FieldRegion, L1FieldValue classes
---
reco/L1/CbmL1.cxx | 23 +---
reco/L1/CbmMvdTrackerIF.cxx | 4 +-
reco/L1/L1Algo/L1Algo.h | 4 +-
reco/L1/L1Algo/L1Field.cxx | 69 ++++++++++--
reco/L1/L1Algo/L1Field.h | 22 +++-
reco/L1/L1Algo/L1InitManager.cxx | 7 +-
reco/L1/L1Algo/L1InitManager.h | 5 +-
reco/L1/L1Algo/L1MaterialInfo.cxx | 80 +++++++++++++-
reco/L1/L1Algo/L1MaterialInfo.h | 25 +++--
reco/L1/L1Algo/L1Parameters.cxx | 144 ++++++++++++++++++++-----
reco/L1/L1Algo/L1Parameters.h | 4 +-
reco/L1/L1Algo/L1Station.cxx | 83 +++++++++++++-
reco/L1/L1Algo/L1Station.h | 25 +++--
reco/L1/L1Algo/L1UMeasurementInfo.cxx | 21 +++-
reco/L1/L1Algo/L1UMeasurementInfo.h | 5 +
reco/L1/L1Algo/L1Utils.h | 27 ++++-
reco/L1/L1Algo/L1XYMeasurementInfo.cxx | 10 +-
reco/L1/L1Algo/L1XYMeasurementInfo.h | 5 +-
reco/L1/vectors/P4_F32vec4.h | 20 ++++
19 files changed, 484 insertions(+), 99 deletions(-)
diff --git a/reco/L1/CbmL1.cxx b/reco/L1/CbmL1.cxx
index fd5bb314f2..1c1a77072d 100644
--- a/reco/L1/CbmL1.cxx
+++ b/reco/L1/CbmL1.cxx
@@ -28,8 +28,6 @@
#include "CbmMuchModuleGem.h"
#include "CbmMuchPad.h"
#include "CbmMuchStation.h"
-#include "CbmMvdDetector.h"
-#include "CbmMvdStationPar.h"
#include "CbmSetup.h" // TODO: To be replaced to the CbmStsTrackerIF !! (S.Zharko)
#include "CbmMvdTrackerIF.h"
#include "CbmStsTrackerIF.h"
@@ -491,28 +489,11 @@ InitStatus CbmL1::Init()
/*** MVD and STS ***/
- CbmStsSetup* stsSetup = CbmStsSetup::Instance();
- if (!stsSetup->IsInit()) { stsSetup->Init(nullptr); }
- if (!stsSetup->IsModuleParsInit()) { stsSetup->SetModuleParameters(fStsParSetModule); }
- if (!stsSetup->IsSensorParsInit()) { stsSetup->SetSensorParameters(fStsParSetSensor); }
- if (!stsSetup->IsSensorCondInit()) { stsSetup->SetSensorConditions(fStsParSetSensorCond); }
-
- NMvdStationsGeom = 0;
- if (fUseMVD) {
- CbmMvdDetector* mvdDetector = CbmMvdDetector::Instance();
- if (mvdDetector) {
- CbmMvdStationPar* mvdStationPar = mvdDetector->GetParameterFile();
- assert(mvdStationPar);
- NMvdStationsGeom = mvdStationPar->GetStationCount();
- }
- }
- //NStsStationsGeom = (fUseSTS) ? CbmStsSetup::Instance()->GetNofStations() : 0;
-
auto mvdInterface = CbmMvdTrackerIF::Instance();
auto stsInterface = CbmStsTrackerIF::Instance();
- //NMvdStationsGeom = mvdInterface->GetNstations();
- NStsStationsGeom = stsInterface->GetNstations();
+ NMvdStationsGeom = (fUseMVD) ? mvdInterface->GetNstations() : 0;
+ NStsStationsGeom = (fUseSTS) ? stsInterface->GetNstations() : 0;
NStationGeom = NMvdStationsGeom + NStsStationsGeom + NMuchStationsGeom + NTrdStationsGeom + NTOFStationGeom;
// Provide crosscheck number of stations for the fpInitManagera
diff --git a/reco/L1/CbmMvdTrackerIF.cxx b/reco/L1/CbmMvdTrackerIF.cxx
index d465dd58d0..9a855edbb0 100644
--- a/reco/L1/CbmMvdTrackerIF.cxx
+++ b/reco/L1/CbmMvdTrackerIF.cxx
@@ -81,7 +81,7 @@ double CbmMvdTrackerIF::GetRmax(int stationId) const
//
int CbmMvdTrackerIF::GetNstations() const
{
- return fMvdMaterial->size();
+ return fMvdStationPar->GetStationCount();
}
//-------------------------------------------------------------------------------------------------------------------------------------
@@ -140,7 +140,7 @@ InitStatus CbmMvdTrackerIF::Init()
fMvdMaterial = &(CbmKF::Instance()->vMvdMaterial);
fMvdStationPar = CbmMvdDetector::Instance()->GetParameterFile();
-
+
if (!fMvdMaterial) { return kFATAL; }
if (!fMvdStationPar) { return kFATAL; }
diff --git a/reco/L1/L1Algo/L1Algo.h b/reco/L1/L1Algo/L1Algo.h
index 5a789389a9..ca2faa6fd2 100644
--- a/reco/L1/L1Algo/L1Algo.h
+++ b/reco/L1/L1Algo/L1Algo.h
@@ -395,8 +395,8 @@ public:
L1InitManager* GetInitManager() { return &fInitManager; }
private:
- L1Parameters fParameters {}; ///< Object of L1Algo parameters class
- L1InitManager fInitManager {&fParameters}; ///< Object of L1Algo initialization manager class
+ L1Parameters fParameters {}; ///< Object of L1Algo parameters class
+ L1InitManager fInitManager {}; ///< Object of L1Algo initialization manager class
/*********************************************************************************************/ /**
diff --git a/reco/L1/L1Algo/L1Field.cxx b/reco/L1/L1Algo/L1Field.cxx
index 8816d409ac..4f71d06cff 100644
--- a/reco/L1/L1Algo/L1Field.cxx
+++ b/reco/L1/L1Algo/L1Field.cxx
@@ -3,6 +3,7 @@
Authors: Sergey Gorbunov, Sergei Zharko [committer] */
#include "L1Field.h"
+#include "L1Utils.h"
#include <iomanip>
#include <iostream>
@@ -12,14 +13,16 @@
// L1FieldValue methods
//
-
//----------------------------------------------------------------------------------------------------------------------
-// TODO: Should it be inline? (S.Zharko)
-void L1FieldValue::Combine(L1FieldValue& B, fvec w)
+//
+void L1FieldValue::CheckConsistency() const
{
- x += w * (B.x - x);
- y += w * (B.y - y);
- z += w * (B.z - z);
+ /* Check SIMD data vectors for consistent initialization */
+ L1Utils::CheckSimdVectorEquality(x, "L1FieldValue::x");
+ L1Utils::CheckSimdVectorEquality(y, "L1FieldValue::y");
+ L1Utils::CheckSimdVectorEquality(z, "L1FieldValue::z");
+
+ // TODO: Any other checks? (S.Zharko)
}
//----------------------------------------------------------------------------------------------------------------------
@@ -35,6 +38,8 @@ std::string L1FieldValue::ToString(int indentLevel) const
return aStream.str();
}
+//----------------------------------------------------------------------------------------------------------------------
+//
std::ostream& operator<<(std::ostream& out, const L1FieldValue& B)
{
return out << B.x[0] << " | " << B.y[0] << " | " << B.z[0];
@@ -46,12 +51,38 @@ std::ostream& operator<<(std::ostream& out, const L1FieldValue& B)
//----------------------------------------------------------------------------------------------------------------------
//
-L1FieldSlice::L1FieldSlice() noexcept
+L1FieldSlice::L1FieldSlice()
+{
+ for (int i = 0; i < L1Constants::size::kMaxNFieldApproxCoefficients; ++i) {
+ cx[i] = L1Utils::kNaN;
+ cy[i] = L1Utils::kNaN;
+ cz[i] = L1Utils::kNaN;
+ }
+}
+
+//----------------------------------------------------------------------------------------------------------------------
+//
+void L1FieldSlice::CheckConsistency() const
{
+ /* Check SIMD data vectors for consistent initialization */
for (int i = 0; i < L1Constants::size::kMaxNFieldApproxCoefficients; ++i) {
- cx[i] = 0.f;
- cy[i] = 0.f;
- cz[i] = 0.f;
+ if (!cx[i].IsHorizontallyEqual()) {
+ std::stringstream msg;
+ msg << "L1FieldSlice: \"cx[" << i << "]\" SIMD vector is inconsistent not all of the words are equal each other: " << cx[i];
+ throw std::logic_error(msg.str());
+ }
+
+ if (!cy[i].IsHorizontallyEqual()) {
+ std::stringstream msg;
+ msg << "L1FieldSlice: \"cy[" << i << "]\" SIMD vector is inconsistent not all of the words are equal each other: " << cy[i];
+ throw std::logic_error(msg.str());
+ }
+
+ if (!cz[i].IsHorizontallyEqual()) {
+ std::stringstream msg;
+ msg << "L1FieldSlice: \"cz[" << i << "]\" SIMD vector is inconsistent not all of the words are equal each other: " << cz[i];
+ throw std::logic_error(msg.str());
+ }
}
}
@@ -133,6 +164,24 @@ L1FieldRegion::L1FieldRegion(float reg[10]) noexcept
{
}
+//----------------------------------------------------------------------------------------------------------------------
+//
+void L1FieldRegion::CheckConsistency() const
+{
+ /* Check SIMD data vectors for consistent initialization */
+ L1Utils::CheckSimdVectorEquality(cx0, "L1FieldRegion::cx0");
+ L1Utils::CheckSimdVectorEquality(cx1, "L1FieldRegion::cx1");
+ L1Utils::CheckSimdVectorEquality(cx2, "L1FieldRegion::cx2");
+ L1Utils::CheckSimdVectorEquality(cy0, "L1FieldRegion::cy0");
+ L1Utils::CheckSimdVectorEquality(cy1, "L1FieldRegion::cy1");
+ L1Utils::CheckSimdVectorEquality(cy2, "L1FieldRegion::cy2");
+ L1Utils::CheckSimdVectorEquality(cz0, "L1FieldRegion::cz0");
+ L1Utils::CheckSimdVectorEquality(cz1, "L1FieldRegion::cz1");
+ L1Utils::CheckSimdVectorEquality(cz2, "L1FieldRegion::cz2");
+ L1Utils::CheckSimdVectorEquality(z0, "L1FieldRegion::z0");
+ // TODO: Any other checks? (S.Zharko)
+}
+
//----------------------------------------------------------------------------------------------------------------------
// TODO: Should it be inline? (S.Zharko)
L1FieldValue L1FieldRegion::Get(const fvec z)
diff --git a/reco/L1/L1Algo/L1Field.h b/reco/L1/L1Algo/L1Field.h
index b80cd4b7a7..b69ce7b8ff 100644
--- a/reco/L1/L1Algo/L1Field.h
+++ b/reco/L1/L1Algo/L1Field.h
@@ -21,6 +21,9 @@ public:
/// \param w weight from 0 to 1 // TODO: Do we need any checks here? (S.Zharko)
void Combine(L1FieldValue& B, fvec w); // TODO: Shouldn't the B parameter be const? (S.Zharko)
+ /// Consistency checker
+ void CheckConsistency() const;
+
/// Operator << overloading
friend std::ostream& operator<<(std::ostream& out, const L1FieldValue& B);
@@ -29,18 +32,30 @@ public:
std::string ToString(int indentLevel) const;
} _fvecalignment;
+inline __attribute__((always_inline)) void L1FieldValue::Combine(L1FieldValue& B, fvec w)
+{
+ x += w * (B.x - x);
+ y += w * (B.y - y);
+ z += w * (B.z - z);
+}
/// Class represents a set of magnetic field approximation coefficients
///
// TODO: Crosscheck the default content (S.Zharko)
class L1FieldSlice {
public:
- L1FieldSlice() noexcept;
+ /// Default constructor
+ L1FieldSlice();
+
+ /// Consistency checker
+ void CheckConsistency() const;
+
/// Gets field value from (x, y) fvec point
/// \param x x-coordinate of input
/// \param y y-coordinate of input
/// \param B the L1FieldValue output
void GetFieldValue(const fvec& x, const fvec& y, L1FieldValue& B) const;
+
/// String representation of class contents
/// \param indentLevel number of indent characters in the output
std::string ToString(int indentLevel = 0) const;
@@ -61,7 +76,10 @@ public:
L1FieldRegion() = default;
L1FieldRegion(float reg[10]) noexcept;
-
+
+ /// Consistency checker
+ void CheckConsistency() const;
+
/// Gets the field vector at z
// TODO: Probably we need a const specifier here, because the method does not change the fields
L1FieldValue Get(const fvec z);
diff --git a/reco/L1/L1Algo/L1InitManager.cxx b/reco/L1/L1Algo/L1InitManager.cxx
index b4ac3de902..b621b83004 100644
--- a/reco/L1/L1Algo/L1InitManager.cxx
+++ b/reco/L1/L1Algo/L1InitManager.cxx
@@ -15,9 +15,6 @@
#include <sstream>
#include "L1Assert.h"
-//-----------------------------------------------------------------------------------------------------------------------
-//
-L1InitManager::L1InitManager(L1Parameters* pParameters) {}
//-----------------------------------------------------------------------------------------------------------------------
//
@@ -128,7 +125,9 @@ void L1InitManager::FormParametersContainer()
fStationsInfo.insert(std::move(node));
++thickMapIt;
}
- //LOG(info) << "L1InitManager: L1Material vector was successfully transfered to L1Algo core";
+
+ // Check the consistency of the parameters object. If object inconsistent, it throws std::logic_error
+ fParameters.CheckConsistency();
}
//-----------------------------------------------------------------------------------------------------------------------
diff --git a/reco/L1/L1Algo/L1InitManager.h b/reco/L1/L1Algo/L1InitManager.h
index f98a03d640..2925d7cceb 100644
--- a/reco/L1/L1Algo/L1InitManager.h
+++ b/reco/L1/L1Algo/L1InitManager.h
@@ -94,10 +94,7 @@ private:
public:
/// Default constructor
- L1InitManager() = delete;
-
- /// Constructor from ptr to L1Paramters object
- L1InitManager(L1Parameters* pParameters);
+ L1InitManager() = default;
/// Destructor
~L1InitManager() = default;
diff --git a/reco/L1/L1Algo/L1MaterialInfo.cxx b/reco/L1/L1Algo/L1MaterialInfo.cxx
index 7a367a66fb..915c367311 100644
--- a/reco/L1/L1Algo/L1MaterialInfo.cxx
+++ b/reco/L1/L1Algo/L1MaterialInfo.cxx
@@ -7,6 +7,7 @@
#include <iomanip>
#include <sstream>
#include <vector>
+#include <FairLogger.h>
/************************
* L1MaterialInfo class *
@@ -18,13 +19,49 @@ std::string L1MaterialInfo::ToString(int indentLevel) const
// TODO: possibly it is better to place indentChar into L1Parameters (S.Zharko)
constexpr char indentChar = '\t';
std::string indent(indentLevel, indentChar);
- aStream << indent << "Station thickness (X) [??]: " << std::setw(12) << std::setfill(' ') << thick[0] << '\n';
- aStream << indent << "Radiation length (X0) [??]: " << std::setw(12) << std::setfill(' ') << RL[0] << '\n';
+ aStream << indent << "Station thickness (X) [cm]: " << std::setw(12) << std::setfill(' ') << thick[0] << '\n';
+ aStream << indent << "Radiation length (X0) [cm]: " << std::setw(12) << std::setfill(' ') << RL[0] << '\n';
aStream << indent << "X/X0: " << std::setw(12) << std::setfill(' ') << RadThick[0] << '\n';
aStream << indent << "log(X/X0): " << std::setw(12) << std::setfill(' ') << logRadThick[0];
return aStream.str();
}
+
+
+//------------------------------------------------------------------------------------------------------------------------------------
+//
+void L1MaterialInfo::CheckConsistency() const
+{
+ /* (i) Checks for the horizontal equality of elements */
+ L1Utils::CheckSimdVectorEquality(thick, "L1MaterialInfo::thick");
+ L1Utils::CheckSimdVectorEquality(RL, "L1MaterialInfo::RL");
+ L1Utils::CheckSimdVectorEquality(RadThick, "L1MaterialInfo::RadThick");
+ L1Utils::CheckSimdVectorEquality(logRadThick, "L1MaterialInfo::logRadThick");
+
+ /* (ii) Checks for physical sence: thick and RL must be larger then 0. */
+ if (thick[0] < fscal(0.) || L1Utils::CmpFloats(thick[0], fscal(0.))) {
+ std::stringstream aStream;
+ aStream <<"L1MaterialInfo: illegal value for station thickness: (" << thick[0] << ", positive value expected) [cm]";
+ throw std::logic_error(aStream.str());
+ }
+
+ if (RL[0] < fscal(0.) || L1Utils::CmpFloats(RL[0], fscal(0.))) {
+ std::stringstream aStream;
+ aStream <<"L1MaterialInfo: illegal value for station radiation length: (" << RL[0] << ", positive value expected) [cm]";
+ throw std::logic_error(aStream.str());
+ }
+
+ /* (iii) Checks for RadThick and logRadThick */
+ if (!L1Utils::CmpFloats(RadThick[0] * RL[0], thick[0])) {
+ throw std::logic_error("L1MaterialInfo: illegal relation between thickness, radiation length and their ratio (RadThick)");
+ }
+
+ if (!L1Utils::CmpFloats(std::exp(logRadThick[0]), RadThick[0])) {
+ throw std::logic_error("L1MaterialInfo: illegal relation between RadThick and logRadThick data fields");
+ }
+}
+
+
/********************
* L1Material class *
********************/
@@ -33,7 +70,8 @@ L1Material::L1Material() {}
//------------------------------------------------------------------------------------------------------------------------------------
//
-L1Material::~L1Material() noexcept {}
+L1Material::~L1Material()
+{}
//------------------------------------------------------------------------------------------------------------------------------------
//
@@ -101,14 +139,48 @@ fvec L1Material::GetRadThick(fvec x, fvec y) const
return r;
}
+
+//------------------------------------------------------------------------------------------------------------------------------------
+//
+void L1Material::CheckConsistency() const
+{
+ /* (i) Check, if the object was initialized */
+ if (fNbins < 1 || std::isnan(fRmax)) {
+ throw std::logic_error("L1Material: class was not initialized");
+ }
+
+ /* (ii) Check if the thickness values correct (non-negative) */
+ bool isThicknessOk = true;
+ for (int i = 0; i < int(fTable.size()); ++i) {
+ if (fTable[i] < 0) {
+ isThicknessOk = false;
+ LOG(error) << "L1Material: found illegal negative station thickness value " << fTable[i] << " at iBinX = " << (i % fNbins)
+ << ", iBin = " << (i / fNbins);
+ }
+ }
+ if (!isThicknessOk) { throw std::logic_error("L1Material: incorrect station thickness values found in the thickness map"); }
+}
+
//------------------------------------------------------------------------------------------------------------------------------------
//
void L1Material::SetBins(int nBins, float stationSize)
{
fNbins = nBins;
fRmax = stationSize;
- fFactor = 0.5 * fNbins / fRmax;
+
+ if (fNbins < 1) {
+ std::stringstream aStream;
+ aStream <<"L1Material: object cannot be initialized with non-positive nBins = " << nBins;
+ throw std::logic_error(aStream.str());
+ }
+ if (stationSize < 0 || L1Utils::CmpFloats(stationSize, 0.f)) {
+ std::stringstream aStream;
+ aStream << "L1Material: object cannot be initialized with non-positive stationStation = " << stationSize << " [cm]";
+ throw std::logic_error(aStream.str());
+ }
+
+ fFactor = 0.5 * fNbins / fRmax;
fTable.resize(fNbins * fNbins);
}
diff --git a/reco/L1/L1Algo/L1MaterialInfo.h b/reco/L1/L1Algo/L1MaterialInfo.h
index 762d048281..85f17a6f70 100644
--- a/reco/L1/L1Algo/L1MaterialInfo.h
+++ b/reco/L1/L1Algo/L1MaterialInfo.h
@@ -11,15 +11,19 @@
#include <vector>
#include "L1Def.h"
+#include "L1Utils.h"
/// Class L1MaterialInfo contains SIMDized vector fields of the
/// The fields of the structure should ONLY be initialized within L1BaseStationInfo::SetMaterial(double, double) method, when the
/// stations sequence is initialized
struct L1MaterialInfo {
- fvec thick {0}; ///< Average thickness of the station in arbitary length units
- fvec RL {0}; ///< Average radiation length (X0) of the station material in THE SAME UNITS as the thickness
- fvec RadThick {0}; ///< Average thickness in units of radiation length (X/X0)
- fvec logRadThick {0};
+ fvec thick {L1Utils::kNaN}; ///< Average thickness of the station in arbitary length units
+ fvec RL {L1Utils::kNaN}; ///< Average radiation length (X0) of the station material in THE SAME UNITS as the thickness
+ fvec RadThick {L1Utils::kNaN}; ///< Average thickness in units of radiation length (X/X0)
+ fvec logRadThick {L1Utils::kNaN};
+
+ /// Verifies class invariant consistency
+ void CheckConsistency() const;
/// String representation of class contents
/// \param indentLevel number of indent characters in the output
@@ -29,7 +33,9 @@ struct L1MaterialInfo {
/// Class L1Material describes a map of station thickness in units of radiation length (X0) to the specific point in XY plane
class L1Material {
public:
- /// Default constructor
+ /// Constructor
+ /// \param nBins Number of rows or columns
+ /// \param stationSize Size of station in x and y dimensions [cm]
L1Material();
/// Copy constructor
@@ -64,6 +70,9 @@ public:
/// \param x X coordinate of the point [cm] (SIMDized vector)
/// \param y Y coordinate of the point [cm] (SIMDized veotor)
fvec GetRadThick(fvec x, fvec y) const;
+
+ /// Verifies class invariant consistency
+ void CheckConsistency() const;
/// Sets value of material thickness in units of X0 for a given cell of the material table
/// WARNING: Indeces of rows and columns in the table runs from 0 to nBins-1 inclusively, where nBins is the number both of rows
@@ -83,9 +92,9 @@ public:
void Swap(L1Material& other) noexcept;
private:
- int fNbins {0}; ///< Number of rows (columns) in the material budget table
- float fRmax {0.f}; ///< Size of the station in x and y dimensions [cm]
- float fFactor {0.f}; ///< Factor used in the recalculation of point coordinates to row/column id
+ int fNbins {-1}; ///< Number of rows (columns) in the material budget table
+ float fRmax {L1Utils::kNaN}; ///< Size of the station in x and y dimensions [cm]
+ float fFactor {L1Utils::kNaN}; ///< Factor used in the recalculation of point coordinates to row/column id
std::vector<float> fTable {}; ///< Material budget table
} _fvecalignment;
diff --git a/reco/L1/L1Algo/L1Parameters.cxx b/reco/L1/L1Algo/L1Parameters.cxx
index d89191e183..8f395aedc6 100644
--- a/reco/L1/L1Algo/L1Parameters.cxx
+++ b/reco/L1/L1Algo/L1Parameters.cxx
@@ -14,7 +14,10 @@
//------------------------------------------------------------------------------------------------------------------------------------
//
-L1Parameters::L1Parameters() {}
+L1Parameters::L1Parameters()
+{
+ fActiveStationGlobalIDs.fill(-1); // by default, all stations are inactive, thus all the IDs must be -1
+}
//------------------------------------------------------------------------------------------------------------------------------------
//
@@ -77,26 +80,118 @@ void L1Parameters::Swap(L1Parameters& other) noexcept
std::swap(fNstationsActive, other.fNstationsActive);
std::swap(fNstationsGeometry, other.fNstationsGeometry);
std::swap(fActiveStationGlobalIDs, other.fActiveStationGlobalIDs);
+}
+
+//------------------------------------------------------------------------------------------------------------------------------------
+//
+void L1Parameters::CheckConsistency() const {
+ LOG(info) << "Consistency test for L1 parameters object... ";
+ /*
+ * Arrays containing numbers of stations
+ *
+ * In the fNstationsActive and fNstationsGeometry array first L1Constants::size::kMaxNdetectors elements are the numbers of stations
+ * in particular detector subsystem. The last element in both arrays corresponds to the total number of stations (geometry or
+ * active). This fact applies restriction on the arrays: the last element must be equal to the sum of all previous elements.
+ *
+ */
+
+ if (std::accumulate(fNstationsGeometry.cbegin(), fNstationsGeometry.cend() - 1, 0) != *(fNstationsGeometry.cend() - 1)) {
+ throw std::logic_error("L1Parameters: invalid object condition: total number of stations provided by geometry "
+ "differs from the sum of the station numbers for individual detector subsystems");
+ };
+
+ if (std::accumulate(fNstationsActive.cbegin(), fNstationsActive.cend() - 1, 0) != *(fNstationsActive.cend() - 1)) {
+ throw std::logic_error("L1Parameters: invalid object condition: total number of stations active in tracking "
+ "differs from the sum of the station numbers for individual detector subsystems");
+ };
+
+ /*
+ * Array of active station IDs
+ *
+ * In the array of active station IDs,
+ * (i) sum of elements, which are not equal -1, must be equal the number of stations,
+ * (ii) subsequence of all the elements, which are not equal -1, must be a gapless subset of integer numbers starting with 0
+ */
+
+ auto filterInactiveStationIDs = [](int x) {return x != -1;};
+ int nStationsCheck = std::count_if(fActiveStationGlobalIDs.cbegin(), fActiveStationGlobalIDs.cend(), filterInactiveStationIDs);
+ if (nStationsCheck != *(fNstationsActive.cend() - 1)) {
+ std::stringstream msg;
+ msg << "L1Parameters: invalid object condition: array of active station IDs is not consistent "
+ << "with the total number of stations (" << nStationsCheck << " vs. " << *(fNstationsActive.cend() - 1) << ' '
+ << "expected)";
+ throw std::logic_error(msg.str());
+ }
+
+ // Check, if the subsequence of all the elements, which are not equal -1, must be a gapless subset of integer numbers
+ // starting from 0. If it is, the testValue in the end must be equal the number of non -1 elements
+
+ std::vector<int> idsCheck(nStationsCheck); // temporary vector containing a sequence of active station IDs without -1 elements
+ std::copy_if(fActiveStationGlobalIDs.cbegin(), fActiveStationGlobalIDs.cend(), idsCheck.begin(), filterInactiveStationIDs);
+ bool isStationIDsOk = true;
+ for (int id = 0; id < nStationsCheck; ++id) {
+ isStationIDsOk = isStationIDsOk && idsCheck[id] == id;
+ }
+ if (!isStationIDsOk) {
+ std::stringstream msg;
+ msg << "L1Parameters: invalid object condition: array of active station IDs is not a gapless subset "
+ << "of integer numbers starting from 0:\n\t";
+ for (auto id: fActiveStationGlobalIDs) {
+ msg << std::setw(3) << std::setfill(' ') << id << ' ';
+ }
+ throw std::logic_error(msg.str());
+ }
+
+ /*
+ * Check target position SIMD vector
+ */
+
+ L1Utils::CheckSimdVectorEquality(fTargetPos[0], "L1Parameters: target position x");
+ L1Utils::CheckSimdVectorEquality(fTargetPos[1], "L1Parameters: target position y");
+ L1Utils::CheckSimdVectorEquality(fTargetPos[2], "L1Parameters: target position z");
+
+ /*
+ * Check vertex field region and value objects at primary vertex
+ */
+
+ fVertexFieldValue.CheckConsistency();
+ fVertexFieldRegion.CheckConsistency();
+
+
+ /*
+ * Check if each station object is consistent itself, and if all of them are placed after the target
+ * NOTE: If a station was not set up, it is accounted inconsistent (uninitialized). In the array of stations there are uninitialized
+ * stations possible (with id > NstationsActiveTotal), thus one should NOT run the loop above over all the stations in array
+ * but only until *(fNstationsActive.cend() - 1) (== NstationsActiveTotal).
+ * TODO: Probably, we should introduce methods, which check the consistency of fully initialized objects such as L1Station,
+ * L1MaterialInfo, etc. (S.Zharko)
+ */
- //for (int i = 0; i < int(L1Constants::size::kMaxNstations); ++i) {
- // int tmp = fActiveStationGlobalIDs[i];
- // fActiveStationGlobalIDs[i] = other.fActiveStationGlobalIDs[i];
- // other.fActiveStationGlobalIDs[i] = tmp;
- //}
-
- //for (int i = 0; i < int(L1Constants::size::kMaxNdetectors + 1); ++i) {
- // int tmp = fNstationsGeometry[i];
- // fNstationsGeometry[i] = other.fNstationsGeometry[i];
- // other.fNstationsGeometry[i] = tmp;
- //}
-
- //for (int i = 0; i < int(L1Constants::size::kMaxNdetectors + 1); ++i) {
- // int tmp = fNstationsActive[i];
- // fNstationsActive[i] = other.fNstationsActive[i];
- // other.fNstationsActive[i] = tmp;
- //}
+ for (int iSt = 0; iSt < *(fNstationsActive.cend() - 1); ++iSt) {
+ fStations[iSt].CheckConsistency();
+ if (fStations[iSt].z[0] < fTargetPos[2][0]) {
+ std::stringstream msg;
+ msg << "L1Parameters: station with global ID = " << iSt << " is placed before target "
+ << "(z_st = " << fStations[iSt].z[0] << " [cm] < z_targ = " << fTargetPos[2][0] << " [cm])";
+ throw std::logic_error(msg.str());
+ }
+ }
+
+ /*
+ * Check thickness maps
+ * NOTE: If a L1Material map was not set up, it is accounted inconsistent (uninitialized). In the array of thickness maps for each
+ * there are uninitialized elements possible (with id > NstationsActiveTotal), thus one should NOT run the loop above over
+ * all the stations in array but only until *(fNstationsActive.cend() - 1) (== NstationsActiveTotal).
+ */
+
+ for (int iSt = 0; iSt < *(fNstationsActive.cend() - 1); ++iSt) {
+ fThickMap[iSt].CheckConsistency();
+ }
+
+ LOG(info) << "Consistency test for L1 parameters object... \033[1;32mpassed\033[0m";
}
+
//------------------------------------------------------------------------------------------------------------------------------------
//
void L1Parameters::Print(int /*verbosityLevel*/) const
@@ -130,26 +225,25 @@ std::string L1Parameters::ToString(int verbosity, int indentLevel) const
aStream << indent << "GEOMETRY:\n";
aStream << indent << indentChar << "TARGET:\n";
- aStream << indent << indentChar << indentChar << "Position (x [cm], y [cm], z[cm]): ";
+ aStream << indent << indentChar << indentChar << "Position:\n";
for (int dim = 0; dim < 3 /*nDimensions*/; ++dim ) {
- aStream << std::setw(12) << std::setfill(' ') << fTargetPos[dim][0] << ' ';
+ aStream << indent << indentChar << indentChar << indentChar << char(120 + dim) << " = " << fTargetPos[dim][0] << " [cm]\n";
}
- aStream << '\n';
aStream << indent << indentChar << "NUMBER OF STATIONS:\n";
aStream << indent << indentChar << indentChar << "Number of stations (Geometry): ";
for (int idx = 0; idx < int(fNstationsGeometry.size()); ++idx ) {
- if (int(fNstationsGeometry.size() - 2) == idx) { aStream << " | total = "; }
+ if (int(fNstationsGeometry.size() - 1) == idx) { aStream << " | total = "; }
aStream << std::setw(2) << std::setfill(' ') << fNstationsGeometry[idx] << ' ';
}
aStream << '\n';
- aStream << indent << indentChar << indentChar << "Number of stations (Active): ";
+ aStream << indent << indentChar << indentChar << "Number of stations (Active): ";
for (int idx = 0; idx < int(fNstationsActive.size()); ++idx ) {
- if (int(fNstationsActive.size() - 2) == idx) { aStream << " | total = "; }
+ if (int(fNstationsActive.size() - 1) == idx) { aStream << " | total = "; }
aStream << std::setw(2) << std::setfill(' ') << fNstationsActive[idx] << ' ';
}
aStream << '\n';
- aStream << indent << indentChar << indentChar << "Active station index: ";
+ aStream << indent << indentChar << indentChar << "Active station indeces: ";
for (int idx = 0; idx < *(fNstationsActive.end() - 1); ++idx ) {
aStream << std::setw(3) << std::setfill(' ') << fActiveStationGlobalIDs[idx] << ' ';
}
diff --git a/reco/L1/L1Algo/L1Parameters.h b/reco/L1/L1Algo/L1Parameters.h
index 5fc1873b66..59d869d81c 100644
--- a/reco/L1/L1Algo/L1Parameters.h
+++ b/reco/L1/L1Algo/L1Parameters.h
@@ -155,6 +155,8 @@ public:
/// Gets L1FieldValue object at primary vertex
const L1FieldValue& GetVertexFieldValue() const { return fVertexFieldValue; }
+ /// Class invariant checker
+ void CheckConsistency() const;
private:
@@ -172,7 +174,7 @@ private:
/// Field value object at primary vertex (between target and the first station)
alignas(L1Constants::misc::kAlignment) L1FieldValue fVertexFieldValue {};
- /// Field value object at primary vertex (between target and the first station)
+ /// Field region object at primary vertex (between target and the first station)
alignas(L1Constants::misc::kAlignment) L1FieldRegion fVertexFieldRegion {};
/// Array of stations
diff --git a/reco/L1/L1Algo/L1Station.cxx b/reco/L1/L1Algo/L1Station.cxx
index 30941e3a5d..3d242286fe 100644
--- a/reco/L1/L1Algo/L1Station.cxx
+++ b/reco/L1/L1Algo/L1Station.cxx
@@ -9,9 +9,87 @@
#include <iomanip>
#include <sstream>
+//------------------------------------------------------------------------------------------------------------------------------------
+//
+void L1Station::CheckConsistency() const
+{
+ /*
+ * Integer fields initialization checks
+ */
+
+ if (type < 0) {
+ std::stringstream msg;
+ msg << "L1Station: station type was not initialized (type = " << type << ", type > 0)";
+ throw std::logic_error(msg.str());
+ }
+
+ if (timeInfo != 0 && timeInfo != 1) {
+ std::stringstream msg;
+ msg << "L1Station: illegal time information flag (timeInfo = " << timeInfo << ", "
+ << "0 [time information is not used] or 1 [time information is used] expected)";
+ throw std::logic_error(msg.str());
+ }
+
+ if (fieldStatus != 0 && fieldStatus != 1) {
+ std::stringstream msg;
+ msg << "L1Station: illegal field status flag (fieldStatus = " << fieldStatus << ", "
+ << "0 [station is outside the field] or 1 [station is inside the field] expected";
+ throw std::logic_error(msg.str());
+ }
+
+ /*
+ * SIMD vector checks: all the words in a SIMD vector must be equal
+ */
+
+ L1Utils::CheckSimdVectorEquality(z, "L1Station::z");
+ L1Utils::CheckSimdVectorEquality(Rmin, "L1Station::Rmin");
+ L1Utils::CheckSimdVectorEquality(Rmax, "L1Station::Rmax");
+ L1Utils::CheckSimdVectorEquality(dt, "L1Station::dt");
+
+ /*
+ * Inner and outer radia checks:
+ * (i) both Rmin and Rmax must be >= 0
+ * (ii) Rmax cannot be smaller or equal to Rmin
+ */
+
+ if (Rmin[0] < 0 || L1Utils::CmpFloats(Rmin[0], Rmax[0]) || Rmax[0] < Rmin[0]) {
+ std::stringstream msg;
+ msg << "L1Station: " << this->ToString() << " has incorrect radia values: "
+ << "Rmin = " << Rmin[0] << " [cm], Rmax = " << Rmax[0] << " [cm] (0 <= Rmin < Rmax expected)";
+ throw std::logic_error(msg.str());
+ }
+
+ /*
+ * Time resolution cannot be smaller then 0
+ */
+
+ if (dt[0] < 0) {
+ std::stringstream msg;
+ msg << "L1Station: " << this->ToString() << " has incorrect time resolution value: "
+ << "dt = " << Rmin[0] << " [ns], Rmax = " << Rmax[0] << " [cm] (0 <= Rmin < Rmax expected)";
+ throw std::logic_error(msg.str());
+ }
+
+ /*
+ * Check consistency of other members
+ */
+
+ materialInfo.CheckConsistency();
+ fieldSlice.CheckConsistency();
+ frontInfo.CheckConsistency();
+ backInfo.CheckConsistency();
+ xInfo.CheckConsistency();
+ yInfo.CheckConsistency();
+ XYInfo.CheckConsistency();
+}
+
// TODO: Improve log style (S.Zh.)
+//------------------------------------------------------------------------------------------------------------------------------------
+//
void L1Station::Print(int verbosity) const { LOG(info) << ToString(verbosity); }
+//------------------------------------------------------------------------------------------------------------------------------------
+//
std::string L1Station::ToString(int verbosityLevel, int indentLevel) const
{
std::stringstream aStream {};
@@ -19,13 +97,12 @@ std::string L1Station::ToString(int verbosityLevel, int indentLevel) const
std::string indent(indentLevel, indentChar);
if (verbosityLevel > 1) {
// TODO: probably we can have verbosity level and address for each underlying object (S.Zharko)
- aStream << '\n' << indent << "Address: " << this << '\n';
}
if (verbosityLevel == 0) {
- aStream << indent << "- <z [cm], typeID> = " << std::setw(12) << std::setfill(' ') << z[0] << ", " << std::setw(4)
- << std::setfill(' ') << type;
+ aStream << "L1Station at z = " << z[0] << " [cm]" ;
}
else {
+ aStream << '\n' << indent << "Address: " << this << '\n';
aStream << indent << "Station type ID: " << std::setw(12) << std::setfill(' ') << type << '\n';
aStream << indent << "Is time info used: " << std::setw(12) << std::setfill(' ') << timeInfo << '\n';
aStream << indent << "Is station placed in field: " << std::setw(12) << std::setfill(' ') << fieldStatus << '\n';
diff --git a/reco/L1/L1Algo/L1Station.h b/reco/L1/L1Algo/L1Station.h
index 7b7d426288..8a8a51f706 100644
--- a/reco/L1/L1Algo/L1Station.h
+++ b/reco/L1/L1Algo/L1Station.h
@@ -11,20 +11,21 @@
#include "L1MaterialInfo.h"
#include "L1UMeasurementInfo.h"
#include "L1XYMeasurementInfo.h"
+#include "L1Utils.h"
/// Structure L1Station
/// Contains a set of geometry parameters for a particular station
///
struct L1Station {
- int type {0};
- int timeInfo {0}; ///< flag: if time information can be used
- int fieldStatus {0}; ///< flag: 1 - station is INSIDE the field, 0 - station is OUTSIDE the field
- fvec z {0}; ///< z position of station
- fvec Rmin {0}; ///< min radius of the station
- fvec Rmax {0}; ///< max radius of the station
- fvec dt {0};
- L1MaterialInfo materialInfo {};
- L1FieldSlice fieldSlice {};
+ int type {-1};
+ int timeInfo {-1}; ///< flag: if time information can be used
+ int fieldStatus {-1}; ///< flag: 1 - station is INSIDE the field, 0 - station is OUTSIDE the field
+ fvec z {L1Utils::kNaN}; ///< z position of station [cm]
+ fvec Rmin {L1Utils::kNaN}; ///< min radius of the station [cm]
+ fvec Rmax {L1Utils::kNaN}; ///< max radius of the station [cm]
+ fvec dt {L1Utils::kNaN}; ///< time resolution [ns]
+ L1MaterialInfo materialInfo {}; ///< structure containing station thickness(X), rad. length (X0), X/X0 and log(X/X0) values
+ L1FieldSlice fieldSlice {};
L1UMeasurementInfo frontInfo {};
L1UMeasurementInfo backInfo {};
L1UMeasurementInfo xInfo {};
@@ -34,11 +35,17 @@ struct L1Station {
/// Prints object fields
/// \param verbosity Verbosity level of the output
void Print(int verbosity = 0) const;
+
/// String representation of class contents
/// \param verbosityLevel Verbosity level of the output
/// \param indentLevel Number of indent characters in the output
std::string ToString(int verbosityLevel = 0, int indentLevel = 0) const;
+ /// Verifies class invariant consistency
+ /// \note Object is considered undefined in the creation time, so this function should be called after the object
+ /// initialization
+ void CheckConsistency() const;
+
} _fvecalignment;
#endif // L1Station_h
diff --git a/reco/L1/L1Algo/L1UMeasurementInfo.cxx b/reco/L1/L1Algo/L1UMeasurementInfo.cxx
index a53f7a5b61..2f8992f2be 100644
--- a/reco/L1/L1Algo/L1UMeasurementInfo.cxx
+++ b/reco/L1/L1Algo/L1UMeasurementInfo.cxx
@@ -7,6 +7,8 @@
#include <iomanip>
#include <sstream>
+//----------------------------------------------------------------------------------------------------------------------
+//
std::string L1UMeasurementInfo::ToString(int indentLevel) const
{
std::stringstream aStream {};
@@ -15,6 +17,23 @@ std::string L1UMeasurementInfo::ToString(int indentLevel) const
std::string indent(indentLevel, indentChar);
aStream << indent << "cos(phi): " << std::setw(12) << std::setfill(' ') << cos_phi[0] << '\n';
aStream << indent << "sin(phi): " << std::setw(12) << std::setfill(' ') << sin_phi[0] << '\n';
- aStream << indent << "sigma2 [??]: " << std::setw(12) << std::setfill(' ') << sigma2[0];
+ aStream << indent << "sigma2 [cm2]: " << std::setw(12) << std::setfill(' ') << sigma2[0];
return aStream.str();
}
+
+//----------------------------------------------------------------------------------------------------------------------
+//
+void L1UMeasurementInfo::CheckConsistency() const
+{
+ /* (i) Checks for the horizontal equality of SIMD vector elements */
+ L1Utils::CheckSimdVectorEquality(cos_phi, "L1UMeasurementsInfo::cos_phi");
+ L1Utils::CheckSimdVectorEquality(sin_phi, "L1UMeasurementsInfo::sin_phi");
+ L1Utils::CheckSimdVectorEquality(sigma2, "L1UMeasurementsInfo::sigma2");
+
+ /*(ii) Sigma2 possible values*/
+ if (sigma2[0] < 0) {
+ std::stringstream msg;
+ msg << "L1UMeasurementsInfo: illegal value of sigma2: " << sigma2[0] << " [cm2] (sigma2 >= 0 excepted)";
+ throw std::logic_error(msg.str());
+ }
+}
diff --git a/reco/L1/L1Algo/L1UMeasurementInfo.h b/reco/L1/L1Algo/L1UMeasurementInfo.h
index 3e3423eaf2..3a33617d2a 100644
--- a/reco/L1/L1Algo/L1UMeasurementInfo.h
+++ b/reco/L1/L1Algo/L1UMeasurementInfo.h
@@ -7,6 +7,7 @@
#include <string>
+#include "L1Utils.h"
#include "L1Def.h"
class L1UMeasurementInfo {
@@ -15,10 +16,14 @@ public:
fvec cos_phi {0.f};
fvec sin_phi {0.f};
fvec sigma2 {0.f};
+
/// String representation of class contents
/// \param indentLevel number of indent characters in the output
std::string ToString(int indentLevel = 0) const;
+ /// Check consistency
+ void CheckConsistency() const;
+
} _fvecalignment;
diff --git a/reco/L1/L1Algo/L1Utils.h b/reco/L1/L1Algo/L1Utils.h
index f47260906e..fee4f61e52 100644
--- a/reco/L1/L1Algo/L1Utils.h
+++ b/reco/L1/L1Algo/L1Utils.h
@@ -13,16 +13,41 @@
#include <iomanip>
#include <limits>
#include <map>
+#include <cmath>
#include <sstream>
#include <string>
#include <unordered_map>
+#include "vectors/P4_F32vec4.h"
+
/// Class contains some utility functions for L1Algo
struct L1Utils {
- /// Some constants
+ /// NaN value for float
static constexpr float kNaN {std::numeric_limits<float>::signaling_NaN()};
+ /// Comparison method for floats
+ /// \param lhs Left floating point to compare
+ /// \param rhs Right floating point to compare
+ /// \return Comparison result: true - equals within epsilon
+ template <typename T>
+ static bool CmpFloats(T lhs, T rhs)
+ {
+ static_assert(!std::numeric_limits<T>::is_integer, "L1Utils::CmpFloatingPoint does not work with integers");
+ return fabs(lhs - rhs) < 2. * std::numeric_limits<T>::epsilon() * fabs(lhs + rhs)
+ || fabs(lhs - rhs) < std::numeric_limits<T>::min();
+ }
+
+
+ static __attribute__((always_inline)) void CheckSimdVectorEquality(fvec v, const char* name)
+ {
+ if (!v.IsHorizontallyEqual()) {
+ std::stringstream msg;
+ msg << name << " SIMD vector is inconsistent, not all of the words are equal each other: " << v;
+ throw std::logic_error(msg.str());
+ }
+ }
+
/// Hash for unordered_map with enum class keys
struct EnumClassHash {
template<typename T>
diff --git a/reco/L1/L1Algo/L1XYMeasurementInfo.cxx b/reco/L1/L1Algo/L1XYMeasurementInfo.cxx
index 07e3cdc92e..b8eba31513 100644
--- a/reco/L1/L1Algo/L1XYMeasurementInfo.cxx
+++ b/reco/L1/L1Algo/L1XYMeasurementInfo.cxx
@@ -3,10 +3,18 @@
Authors: Sergey Gorbunov, Sergei Zharko [committer] */
#include "L1XYMeasurementInfo.h"
-
+#include "L1Utils.h"
#include <iomanip>
#include <sstream> // for stringstream
+void L1XYMeasurementInfo::CheckConsistency() const
+{
+ /* (i) Checks for the horizontal equality of SIMD vector elements */
+ L1Utils::CheckSimdVectorEquality(C00, "L1XYMeasurementsInfo::C00");
+ L1Utils::CheckSimdVectorEquality(C10, "L1XYMeasurementsInfo::C10");
+ L1Utils::CheckSimdVectorEquality(C11, "L1XYMeasurementsInfo::C11");
+}
+
std::string L1XYMeasurementInfo::ToString(int indentLevel) const
{
std::stringstream aStream {};
diff --git a/reco/L1/L1Algo/L1XYMeasurementInfo.h b/reco/L1/L1Algo/L1XYMeasurementInfo.h
index 569a80d6a2..d8b7ac0c6b 100644
--- a/reco/L1/L1Algo/L1XYMeasurementInfo.h
+++ b/reco/L1/L1Algo/L1XYMeasurementInfo.h
@@ -15,10 +15,13 @@ public:
fvec C00 {0};
fvec C10 {0};
fvec C11 {0};
+
+ /// Consistency checker
+ void CheckConsistency() const;
+
/// String representation of class contents
/// \param indentLevel number of indent characters in the output
std::string ToString(int indentLevel = 0) const;
-
} _fvecalignment;
diff --git a/reco/L1/vectors/P4_F32vec4.h b/reco/L1/vectors/P4_F32vec4.h
index cefaa3dc7c..1d22e5e241 100644
--- a/reco/L1/vectors/P4_F32vec4.h
+++ b/reco/L1/vectors/P4_F32vec4.h
@@ -228,6 +228,26 @@ public:
return strm;
}
+ /// Checks, if all bands are equal
+ /// NOTE: two values defined as signaling_NaN() are not equal, thus if there are all or one
+ /// of the words are kNaN, the function returns false
+ bool IsHorizontallyEqual() const
+ {
+ return v[0] == v[1] && v[1] == v[2] && v[2] == v[3];
+ }
+
+ /// Checks, if any of the bands is NaN
+ bool IsNanAny() const
+ {
+ return std::isnan(v[0]) || std::isnan(v[1]) || std::isnan(v[2]) || std::isnan(v[3]);
+ }
+
+ /// Checks, if all of the bands is NaN
+ bool IsNanAll() const
+ {
+ return std::isnan(v[0]) && std::isnan(v[1]) && std::isnan(v[2]) && std::isnan(v[3]);
+ }
+
} __attribute__((aligned(16)));
--
GitLab