diff --git a/core/detectors/much/CbmMuchTrackingInterface.h b/core/detectors/much/CbmMuchTrackingInterface.h
index 362288aeef77153c448ff96b5f216c2edec80e5e..467bd8d4c40cdd30b9dcb259cfc7dbf4dfb87df1 100644
--- a/core/detectors/much/CbmMuchTrackingInterface.h
+++ b/core/detectors/much/CbmMuchTrackingInterface.h
@@ -102,18 +102,14 @@ public:
/// Gets a tracking station of a CbmPixelHit
/// \param hit A pointer to CbmPixelHit
/// \return Local index of the tracking station
- int GetTrackingStationIndex(const CbmPixelHit* hit) const
- {
- return GetTrackingStationIndex(hit->GetAddress());
- }
+ int GetTrackingStationIndex(const CbmPixelHit* hit) const { return GetTrackingStationIndex(hit->GetAddress()); }
/// Gets a tracking station by the address of element
/// \param address Unique element address
/// \return Local index of the tracking station
int GetTrackingStationIndex(int address) const
{
- return CbmMuchGeoScheme::GetStationIndex(address) * 3
- + CbmMuchGeoScheme::GetLayerIndex(address);
+ return CbmMuchGeoScheme::GetStationIndex(address) * 3 + CbmMuchGeoScheme::GetLayerIndex(address);
}
/// Gets max size of a station along the X-axis
diff --git a/core/detectors/sts/CbmStsTrackingInterface.h b/core/detectors/sts/CbmStsTrackingInterface.h
index d279ba6d73a98fa38c4cf7f1c9cad570637e4912..830efa4ee2a1410798aef3cbe1e557d4fc863e01 100644
--- a/core/detectors/sts/CbmStsTrackingInterface.h
+++ b/core/detectors/sts/CbmStsTrackingInterface.h
@@ -107,18 +107,12 @@ public:
/// Gets a tracking station of a CbmPixelHit
/// \param hit A pointer to CbmPixelHit
/// \return Local index of the tracking station
- int GetTrackingStationIndex(const CbmPixelHit* hit) const
- {
- return GetTrackingStationIndex(hit->GetAddress());
- }
+ int GetTrackingStationIndex(const CbmPixelHit* hit) const { return GetTrackingStationIndex(hit->GetAddress()); }
/// Gets a tracking station by the address of element
/// \param address Unique element address
/// \return Local index of the tracking station
- int GetTrackingStationIndex(int address) const
- {
- return CbmStsSetup::Instance()->GetStationNumber(address);
- }
+ int GetTrackingStationIndex(int address) const { return CbmStsSetup::Instance()->GetStationNumber(address); }
/// Gets max size of a station along the X-axis
/// \param stationId Tracking station ID in the setup (NOTE: must be in range [0..GetNstations()-1])
diff --git a/core/detectors/tof/CbmTofTrackingInterface.h b/core/detectors/tof/CbmTofTrackingInterface.h
index 1a0fa5fdba3a00558246b3e5e93fabb88f3ded52..5e53e389d559633e2d8a2cd75dbb18c6ce26d9c1 100644
--- a/core/detectors/tof/CbmTofTrackingInterface.h
+++ b/core/detectors/tof/CbmTofTrackingInterface.h
@@ -121,7 +121,7 @@ public:
// CbmTofAddress::GetRpcId(address));
// NOTE: Implement, when the "mcbm_beam_2021_07_surveyed" parameters will be fixed
- return -1; // iSt;
+ return -1; // iSt;
}
diff --git a/core/detectors/trd/CbmTrdTrackingInterface.h b/core/detectors/trd/CbmTrdTrackingInterface.h
index 21108d9d0bedc4c3879bbcf452df0402fb3cdcda..411a42bb0b6575af24f55abb80ef2cc73cf67357 100644
--- a/core/detectors/trd/CbmTrdTrackingInterface.h
+++ b/core/detectors/trd/CbmTrdTrackingInterface.h
@@ -103,10 +103,7 @@ public:
/// Gets a tracking station by the address
/// \param address Unique element address
/// \return Local index of the tracking station
- int GetTrackingStationIndex(int address) const
- {
- return CbmTrdAddress::GetLayerId(address);
- }
+ int GetTrackingStationIndex(int address) const { return CbmTrdAddress::GetLayerId(address); }
/// Gets max size of a station along the X-axis
/// \param stationId Tracking station ID in the setup (NOTE: must be in range [0..GetNstations()-1])
diff --git a/core/qa/CbmQaHist.h b/core/qa/CbmQaHist.h
index edb6f195ce47cf47444339e24e92df2d26d69f5b..5f838c0a55301de5d4408e83641ba6a78b7e1084 100644
--- a/core/qa/CbmQaHist.h
+++ b/core/qa/CbmQaHist.h
@@ -84,7 +84,9 @@ public:
f->SetLineColor(kRed);
f->SetLineWidth(3);
TPaveStats* st = (TPaveStats*) this->FindObject("stats");
- if (!st) { LOG(fatal) << "CbmQaHist: can not access statistics of histogram with name \"" << this->GetName() << '\"'; }
+ if (!st) {
+ LOG(fatal) << "CbmQaHist: can not access statistics of histogram with name \"" << this->GetName() << '\"';
+ }
else {
st->SetX1NDC(0.6);
st->SetX2NDC(0.940);
diff --git a/macro/run/run_reco.C b/macro/run/run_reco.C
index 450e632bb3725ee004775bfbd4ddc2964303db03..8b9add96503cd1c8aed8c4c140003eaebafd6580 100644
--- a/macro/run/run_reco.C
+++ b/macro/run/run_reco.C
@@ -380,6 +380,7 @@ void run_reco(TString input = "", Int_t nTimeSlices = -1, Int_t firstTimeSlice =
// L1 tracking
auto l1 = (debugWithMC) ? new CbmL1("L1", 2, 3) : new CbmL1("L1", 0);
+ l1->SetInputConfigName(TString(gSystem->Getenv("VMCWORKDIR")) + "/reco/L1/L1Algo/L1ConfigExample.yaml");
// --- Material budget file names
TString mvdGeoTag;
diff --git a/mvd/CbmMvdTrackingInterface.h b/mvd/CbmMvdTrackingInterface.h
index 01b845e3ed29312ad822b01dd68096735a5bc7c9..f7604cbb4a7986b0fff88d209e0f1f0651031dca 100644
--- a/mvd/CbmMvdTrackingInterface.h
+++ b/mvd/CbmMvdTrackingInterface.h
@@ -12,8 +12,8 @@
#ifndef CbmMvdTrackingInterface_h
#define CbmMvdTrackingInterface_h 1
-#include "CbmMvdDetectorId.h"
#include "CbmMvdDetector.h"
+#include "CbmMvdDetectorId.h"
#include "CbmMvdHit.h"
#include "CbmMvdStationPar.h"
#include "CbmPixelHit.h"
@@ -116,14 +116,11 @@ public:
}();
return hitMvd->GetStationNr();
}
-
+
/// Gets a tracking station by the address of element (detectorID in terms of MVD)
/// \param detectorId Unique element address (detectorID in terms of MVD)
/// \return Local index of the tracking station
- int GetTrackingStationIndex(int detectorId) const
- {
- return StationNr(detectorId);
- }
+ int GetTrackingStationIndex(int detectorId) const { return StationNr(detectorId); }
/// Gets max size of a tracking station along the X-axis
/// \param stationId Tracking station ID in the setup (NOTE: must be in range [0..GetNstations()-1])
diff --git a/reco/L1/CMakeLists.txt b/reco/L1/CMakeLists.txt
index 4f43b55df77fa1c511b345239113b224c61be170..25d1c057ec53b09d3b6aae851fa3d0d8d97321bd 100644
--- a/reco/L1/CMakeLists.txt
+++ b/reco/L1/CMakeLists.txt
@@ -98,7 +98,7 @@ endif()
link_directories( ${LINK_DIRECTORIES})
-
+# ----- Compilation sources --------------------------------------------
set(SRCS
#L1Algo / L1CATrackFinder.cxx
#CbmL1Performance.cxx
@@ -149,17 +149,21 @@ L1Algo/L1CAIteration.cxx
L1Algo/L1BaseStationInfo.cxx
L1Algo/L1InitManager.cxx
L1Algo/L1Parameters.cxx
+L1Algo/L1ConfigRW.cxx
L1Algo/utils/L1AlgoDraw.cxx
L1Algo/utils/L1AlgoEfficiencyPerformance.cxx
L1Algo/utils/L1AlgoPulls.cxx
+
ParticleFinder/CbmL1PFFitter.cxx
ParticleFinder/CbmL1PFMCParticle.cxx
qa/CbmTrackerInputQaTrd.cxx
qa/CbmTrackingInputQaSts.cxx
)
+# ----------------------------------------------------------------------
+# ----- Headers --------------------------------------------------------
set(HEADERS
CbmL1CATrdTrackFinderSA.h
CbmL1.h
@@ -192,6 +196,7 @@ L1Algo/L1Vector.h
qa/CbmTrackerInputQaTrd.h
qa/CbmTrackingInputQaSts.h
)
+# ----------------------------------------------------------------------
@@ -247,6 +252,7 @@ Set(DEPENDENCIES
CbmRecoSts
CbmQaBase
boost_regex
+ external::yaml-cpp
)
if (OPENMP_FOUND AND APPLE)
@@ -280,6 +286,7 @@ Install(FILES CbmL1Counters.h
L1Algo/L1InitManager.h
L1Algo/L1CAIteration.h
L1Algo/L1Parameters.h
+ L1Algo/L1ConfigRW.h
L1Algo/L1Constants.h
L1Algo/L1Utils.h
L1Algo/L1NaN.h
diff --git a/reco/L1/CbmL1.cxx b/reco/L1/CbmL1.cxx
index d9bbb550d6b2d18f36e876166fe9d3b9ab56730b..fb8ca7941c32dab66e972d9bc963a534706e375d 100644
--- a/reco/L1/CbmL1.cxx
+++ b/reco/L1/CbmL1.cxx
@@ -617,8 +617,8 @@ InitStatus CbmL1::Init()
** **
****************************************/
- // TODO: Need to provide a selection: default iterations input (these hardcoded ones), input from file or input
- // from runing macro (S.Zharko)
+ // TODO: Need to provide a selection: default iterations input (these hard-coded ones), input from file or input
+ // from running macro (S.Zharko)
auto trackingIterFastPrim = L1CAIteration("FastPrimIter");
trackingIterFastPrim.SetTrackChi2Cut(10.f);
trackingIterFastPrim.SetTripletChi2Cut(23.4450f); // = 7.815 * 3; // prob = 0.05
@@ -1024,7 +1024,7 @@ void CbmL1::Reconstruct(CbmEvent* event)
if (fVerbose > 1) { cout << "L1 Track fitter ok" << endl; }
- // save recontstructed tracks
+ // save reconstructed tracks
if (fLegacyEventMode) vRTracksCur.clear();
int start_hit = 0;
diff --git a/reco/L1/CbmL1.h b/reco/L1/CbmL1.h
index 0d9741de880892ad110492b41559939fcaff240a..aee4c87c52a290214bdd85d3ff22fced3b94d3eb 100644
--- a/reco/L1/CbmL1.h
+++ b/reco/L1/CbmL1.h
@@ -175,65 +175,26 @@ public:
if (fileName != "") fMatBudgetFileName[L1DetectorID::kTof] = fileName;
}
- /// Correction function for the material budget map
- /// It fills bins with no statistics
- template<L1DetectorID detID>
- void ApplyCorrectionToMaterialMap(L1Material& material, const L1MaterialInfo& homogenious)
+ /// \brief Sets a name for the input configuration file
+ /// If the file is undefined, default tracking parameters will be used. Otherwise, the default parameters will be
+ /// overridden with ones from the configuration file
+ /// \param filename Name of the input tracking configuration file
+ void SetInputConfigName(const char* filename)
{
- // TODO: unify the correction function for all detectors
- [[maybe_unused]] float minVal = 0.;
- if constexpr (detID == L1DetectorID::kMuch) { minVal = 0.15f; }
- else if constexpr (detID == L1DetectorID::kTof || detID == L1DetectorID::kTrd) {
- minVal = 0.0015f;
- }
-
- // A bit ugly solution, but so can we avoid dependency on input maps file
- std::vector<float> keepRow {};
- if constexpr (detID != L1DetectorID::kSts) { keepRow.resize(material.GetNbins()); }
-
- for (int iBinX = 0; iBinX < material.GetNbins(); ++iBinX) {
- if constexpr (detID == L1DetectorID::kTof) { minVal = 0.0015f; }
- if constexpr (detID != L1DetectorID::kSts) {
- for (int iBinY = 0; iBinY < material.GetNbins(); ++iBinY) {
- keepRow[iBinY] = material.GetRadThick(iBinX, iBinY);
- }
- }
- for (int iBinY = 0; iBinY < material.GetNbins(); ++iBinY) {
- if constexpr (detID == L1DetectorID::kMvd) {
- // Correction for holes in the material map
- if (material.GetRadThick(iBinX, iBinY) < homogenious.RadThick[0]) {
- if (iBinY > 0 && iBinY < material.GetNbins() - 1) {
- material.SetRadThick(iBinX, iBinY, TMath::Min(keepRow[iBinY - 1], keepRow[iBinY + 1]));
- }
- }
+ if (fpInitManager) { fpInitManager->SetInputConfigName(std::string(filename)); }
+ }
- // Correction for the hardcodded value of RadThick of MVD stations
- if (material.GetRadThick(iBinX, iBinY) < 0.0015) { material.SetRadThick(iBinX, iBinY, 0.0015); }
- }
- else if constexpr (detID == L1DetectorID::kSts) {
- if (material.GetRadThick(iBinX, iBinY) < homogenious.RadThick[0]) {
- material.SetRadThick(iBinX, iBinY, homogenious.RadThick[0]);
- }
- }
- else if constexpr (detID == L1DetectorID::kMuch || detID == L1DetectorID::kTrd || detID == L1DetectorID::kTof) {
- // Correction for holes in the material map
- if (L1Algo::TrackingMode::kGlobal != fTrackingMode) {
- if ((iBinY > 0) && (iBinY < material.GetNbins() - 1)) {
- material.SetRadThick(iBinX, iBinY, TMath::Min(keepRow[iBinY - 1], keepRow[iBinY + 1]));
- }
- }
- float val = material.GetRadThick(iBinX, iBinY);
- if (val > 0.0015) { // remember last non-zero value
- minVal = val;
- }
- else { // empty bin with no statistics, fill it with the neoighbours value
- material.SetRadThick(iBinX, iBinY, minVal);
- }
- }
- }
- }
+ /// \brief Sets a name for the output configuration file
+ /// \param filename Name of the input tracking configuration file
+ void SetOutputConfigName(const char* filename)
+ {
+ if (fpInitManager) { fpInitManager->SetOutputConfigName(std::string(filename)); }
}
+ /// Correction function for the material budget map
+ /// It fills bins with no statistics
+ template<L1DetectorID detID>
+ void ApplyCorrectionToMaterialMap(L1Material& material, const L1MaterialInfo& homogenious);
/// Utility to map the L1DetectorID items into detector names
static constexpr const char* GetDetectorName(L1DetectorID detectorID)
@@ -294,7 +255,7 @@ private:
/// Read information about hits, mcPoints and mcTracks into L1 classes
/// Repacks data from the external TClonesArray objects to the internal L1 arrays
- /// \param fData_ Pointer to the target object containig L1Algo internal arrays of hits
+ /// \param fData_ Pointer to the target object containing L1Algo internal arrays of hits
/// \param TsStart Reference to the timeslice start time
/// \param TsLength Reference to the timeslice length
/// \param TsOverlap Reference to the timeslice overlap length (does not used at the moment)
@@ -308,7 +269,7 @@ private:
/// \param MC Pointer to a target CbmL1MCPoint object
/// \param iPoint Index of the point into the input MC points CbmMCDataArray object for the particular detector
/// \param MVD Index of the detector subsystem
- /// \return flag: false - success, true - some errors occured
+ /// \return flag: false - success, true - some errors occurred
bool ReadMCPoint(CbmL1MCPoint* MC, int iPoint, int MVD); // help procedure
/// Converts data from generic FairMCPoint based class to the CbmL1MCPoint
@@ -317,7 +278,7 @@ private:
/// \param file Index of the input file
/// \param event Index of the input event
/// \param MVD Index of the detector subsystem
- /// \return flag: false - success, true - some errors occured
+ /// \return flag: false - success, true - some errors occurred
// TODO: Probably, we should replace input parameter MVD with the template (S.Zharko)
bool ReadMCPoint(CbmL1MCPoint* MC, int iPoint, int file, int event, int MVD);
@@ -332,13 +293,13 @@ private:
*/
/// Procedure for match hits and MCPoints.
- /// Reads information about correspondence between hits and mcpoints and fill CbmL1MCPoint::hitIds and CbmL1Hit::mcPointIds arrays
+ /// Reads information about correspondence between hits and MC-points and fill CbmL1MCPoint::hitIds and CbmL1Hit::mcPointIds arrays
/// should be called after fill of algo
void HitMatch();
- void FieldApproxCheck(); // Build histos with difference between Field map and approximated field
- void FieldIntegralCheck(); // Build 2D histo: dependence of the field integral on phi and theta
- void InputPerformance(); // Build histos about input data, like hit pulls, etc.
+ void FieldApproxCheck(); // Build histograms with difference between Field map and approximated field
+ void FieldIntegralCheck(); // Build 2D histogram: dependence of the field integral on phi and theta
+ void InputPerformance(); // Build histograms about input data, like hit pulls, etc.
void TimeHist();
/// Reconstruction Performance
@@ -516,4 +477,65 @@ private:
ClassDef(CbmL1, 0);
};
+
+// ---------------------------------------------------------------------------------------------------------------------
+//
+template<L1DetectorID detID>
+void CbmL1::ApplyCorrectionToMaterialMap(L1Material& material, const L1MaterialInfo& homogenious)
+{
+ // TODO: unify the correction function for all detectors
+ [[maybe_unused]] float minVal = 0.;
+ if constexpr (detID == L1DetectorID::kMuch) { minVal = 0.15f; }
+ else if constexpr (detID == L1DetectorID::kTof || detID == L1DetectorID::kTrd) {
+ minVal = 0.0015f;
+ }
+
+ // A bit ugly solution, but so can we avoid dependency on input maps file
+ std::vector<float> keepRow {};
+ if constexpr (detID != L1DetectorID::kSts) { keepRow.resize(material.GetNbins()); }
+
+ for (int iBinX = 0; iBinX < material.GetNbins(); ++iBinX) {
+ if constexpr (detID == L1DetectorID::kTof) { minVal = 0.0015f; }
+ if constexpr (detID != L1DetectorID::kSts) {
+ for (int iBinY = 0; iBinY < material.GetNbins(); ++iBinY) {
+ keepRow[iBinY] = material.GetRadThick(iBinX, iBinY);
+ }
+ }
+ for (int iBinY = 0; iBinY < material.GetNbins(); ++iBinY) {
+ if constexpr (detID == L1DetectorID::kMvd) {
+ // Correction for holes in the material map
+ if (material.GetRadThick(iBinX, iBinY) < homogenious.RadThick[0]) {
+ if (iBinY > 0 && iBinY < material.GetNbins() - 1) {
+ material.SetRadThick(iBinX, iBinY, TMath::Min(keepRow[iBinY - 1], keepRow[iBinY + 1]));
+ }
+ }
+
+ // Correction for the hard-coded value of RadThick of MVD stations
+ if (material.GetRadThick(iBinX, iBinY) < 0.0015) { material.SetRadThick(iBinX, iBinY, 0.0015); }
+ }
+ else if constexpr (detID == L1DetectorID::kSts) {
+ if (material.GetRadThick(iBinX, iBinY) < homogenious.RadThick[0]) {
+ material.SetRadThick(iBinX, iBinY, homogenious.RadThick[0]);
+ }
+ }
+ else if constexpr (detID == L1DetectorID::kMuch || detID == L1DetectorID::kTrd || detID == L1DetectorID::kTof) {
+ // Correction for holes in the material map
+ if (L1Algo::TrackingMode::kGlobal != fTrackingMode) {
+ if ((iBinY > 0) && (iBinY < material.GetNbins() - 1)) {
+ material.SetRadThick(iBinX, iBinY, TMath::Min(keepRow[iBinY - 1], keepRow[iBinY + 1]));
+ }
+ }
+ float val = material.GetRadThick(iBinX, iBinY);
+ if (val > 0.0015) { // remember last non-zero value
+ minVal = val;
+ }
+ else { // empty bin with no statistics, fill it with the neighbours value
+ material.SetRadThick(iBinX, iBinY, minVal);
+ }
+ }
+ }
+ }
+}
+
+
#endif //_CbmL1_h_
diff --git a/reco/L1/CbmL1ReadEvent.cxx b/reco/L1/CbmL1ReadEvent.cxx
index 45ea39e4cf26b124a84301b2a25e1b27ed0c4991..2d85d7a05112d6bb8143db4dfd6e661610857ab5 100644
--- a/reco/L1/CbmL1ReadEvent.cxx
+++ b/reco/L1/CbmL1ReadEvent.cxx
@@ -1482,8 +1482,8 @@ void CbmL1::HitMatch()
static_cast<const CbmMatch*>(listStsClusterMatch->At(sh->GetBackClusterId()));
CbmMatch stsHitMatch;
- Float_t Sum_of_front = 0; // total weight of all front links
- Float_t Sum_of_back = 0; // total weight of all back links
+ Float_t Sum_of_front = 0; // total weight of all front links
+ Float_t Sum_of_back = 0; // total weight of all back links
for (Int_t iFrontLink = 0; iFrontLink < frontClusterMatch->GetNofLinks(); iFrontLink++) {
diff --git a/reco/L1/L1Algo/L1BaseStationInfo.cxx b/reco/L1/L1Algo/L1BaseStationInfo.cxx
index a978aaaa8c47bd0e8dba87b35df150780714dda1..0592d84b0d067df8a576403c0b7fca4f16b54eef 100644
--- a/reco/L1/L1Algo/L1BaseStationInfo.cxx
+++ b/reco/L1/L1Algo/L1BaseStationInfo.cxx
@@ -4,7 +4,7 @@
/************************************************************************************************************
* @file L1BaseStationInfo.cxx
- * @bried Realization of L1BaseStationInfo class large methods
+ * @brief Realization of L1BaseStationInfo class large methods
* @since 18.01.2022
*
* The class is implemented to connect concrete experimental setup (CBM or BMN in CbmL1
diff --git a/reco/L1/L1Algo/L1BaseStationInfo.h b/reco/L1/L1Algo/L1BaseStationInfo.h
index caf3e36130d49112a516371b82697aeabfeb9189..45910b0ec3166544fe0a62ffcfd921e92dfdfeaf 100644
--- a/reco/L1/L1Algo/L1BaseStationInfo.h
+++ b/reco/L1/L1Algo/L1BaseStationInfo.h
@@ -4,7 +4,7 @@
/************************************************************************************************************
* @file L1BaseStationInfo.h
- * @bried A base class for a L1 station interface
+ * @brief A base class for a L1 station interface
* @since 18.12.2021
*
* The class is implemented to connect concrete experimental setup (CBM or BMN in CbmL1
diff --git a/reco/L1/L1Algo/L1CAIteration.cxx b/reco/L1/L1Algo/L1CAIteration.cxx
index 022507e46d32ed84d9eabd4d9028d0ed38877e2c..359ebb1ecc9dd3e03828bab56d4a2bb4792f0c0c 100644
--- a/reco/L1/L1Algo/L1CAIteration.cxx
+++ b/reco/L1/L1Algo/L1CAIteration.cxx
@@ -132,7 +132,7 @@ std::string L1CAIteration::ToString(int indentLevel) const
aStream << indent << indentChar << "Doublet chi2 cut: " << fDoubletChi2Cut << '\n';
aStream << indent << indentChar << "Pick gather: " << fPickGather << '\n';
aStream << indent << indentChar << "Pick neighbour: " << fPickNeighbour << '\n';
- aStream << indent << indentChar << "Max invariant momentum: " << fMaxInvMom << '\n';
+ aStream << indent << indentChar << "Max inverse momentum: " << fMaxInvMom << '\n';
aStream << indent << indentChar << "Max slope at primary vertex: " << fMaxSlopePV << '\n';
aStream << indent << indentChar << "Max slope: " << fMaxSlope << '\n';
aStream << indent << indentChar << "Max DZ: " << fMaxDZ << '\n';
diff --git a/reco/L1/L1Algo/L1CAIteration.h b/reco/L1/L1Algo/L1CAIteration.h
index 42ffe9563de6a35443e39ad5dd633466389ab985..c4dcc3e82576be32f9b66fb2919df0c90be8b1a9 100644
--- a/reco/L1/L1Algo/L1CAIteration.h
+++ b/reco/L1/L1Algo/L1CAIteration.h
@@ -109,6 +109,7 @@ public:
void SetMaxDZ(float input) { fMaxDZ = input; }
/// Sets max considered q/p for tracks
+ /// TODO: Replace with minimum momentum setter (S.Zharko)
void SetMaxInvMom(float input) { fMaxInvMom = input; }
/// Sets max slope (tx\ty) in 3D hit position of a triplet
diff --git a/reco/L1/L1Algo/L1ConfigExample.yaml b/reco/L1/L1Algo/L1ConfigExample.yaml
new file mode 100644
index 0000000000000000000000000000000000000000..cdc1b2c7c8bb26b9b3b8e8f548309063910deb80
--- /dev/null
+++ b/reco/L1/L1Algo/L1ConfigExample.yaml
@@ -0,0 +1,74 @@
+################################################
+# #
+# L1 Tracking Configuration File (prototype) #
+# #
+################################################
+
+---
+l1:
+ setup: Null
+ algorithm:
+ # Parameters dependent on the CA track finder iteration
+ # Iterations will be used in the same sequence as they defined below
+ iterations:
+ - name: "TESTFastPrimIter"
+ track_chi2_cut: 10.
+ triplet_chi2_cut: 23.4450
+ doublet_chi2_cut: 7.56327
+ pick_gather: 3.0
+ pick_neighbour: 5.0
+ min_momentum: 0.5 # GeV/c
+ max_slope_pv: 1.1
+ max_slope: 2.748
+ max_dz: 0. # cm
+ min_start_triplet_lvl: 0
+ target_pos_sigma_x: 1. # cm
+ target_pos_sigma_y: 1. # cm
+ is_primary: true
+ is_electron: false
+ - name: "AllPrimIter"
+ track_chi2_cut: 10.
+ triplet_chi2_cut: 23.4450
+ doublet_chi2_cut: 7.56327
+ pick_gather: 4.0
+ pick_neighbour: 5.0
+ min_momentum: 0.05 # GeV/c
+ max_slope_pv: 1.1
+ max_slope: 2.748
+ max_dz: 0.1 # cm
+ min_start_triplet_lvl: 0
+ target_pos_sigma_x: 1. # cm
+ target_pos_sigma_y: 1. # cm
+ is_primary: true
+ is_electron: false
+ - name: "AllPrimJumpIter"
+ track_chi2_cut: 10.
+ triplet_chi2_cut: 18.7560
+ doublet_chi2_cut: 7.56327
+ pick_gather: 4.0
+ pick_neighbour: 5.0
+ min_momentum: 0.1 # GeV/c
+ max_slope_pv: 1.1
+ max_slope: 2.748
+ max_dz: 0.1 # cm
+ min_start_triplet_lvl: 0
+ target_pos_sigma_x: 5. # cm
+ target_pos_sigma_y: 5. # cm
+ is_primary: true
+ is_electron: false
+ - name: "AllSecIter"
+ track_chi2_cut: 10.
+ triplet_chi2_cut: 18.7560
+ doublet_chi2_cut: 7.56327
+ pick_gather: 4.0
+ pick_neighbour: 5.0
+ min_momentum: 0.1 # GeV/c
+ max_slope_pv: 1.5
+ max_slope: 2.748
+ max_dz: 0.1 # cm
+ min_start_triplet_lvl: 1
+ target_pos_sigma_x: 10. # cm
+ target_pos_sigma_y: 10. # cm
+ is_primary: false
+ is_electron: false
+...
diff --git a/reco/L1/L1Algo/L1ConfigRW.cxx b/reco/L1/L1Algo/L1ConfigRW.cxx
new file mode 100644
index 0000000000000000000000000000000000000000..7b89b925eb941c4a586b86a5f12fd7749d129a4d
--- /dev/null
+++ b/reco/L1/L1Algo/L1ConfigRW.cxx
@@ -0,0 +1,116 @@
+/* Copyright (C) 2022 GSI Helmholtzzentrum fuer Schwerionenforschung, Darmstadt
+ SPDX-License-Identifier: GPL-3.0-only
+ Authors: Sergei Zharko [committer] */
+
+/// @file L1ConfigRW.h
+/// @brief Configuration parameter file reader/writer for L1 tracking algorithm (implementation)
+/// @author S.Zharko <s.zharko@gsi.de>
+/// @since 18.07.2022
+
+#include "L1ConfigRW.h"
+
+#include <iostream>
+#include <numeric>
+#include <sstream>
+
+#include <yaml-cpp/yaml.h>
+
+#include "L1CAIteration.h"
+#include "L1InitManager.h"
+
+using namespace std::string_literals;
+
+// ---------------------------------------------------------------------------------------------------------------------
+//
+L1ConfigRW::L1ConfigRW(L1InitManager* pInitManager, int verbose) : fpInitManager(pInitManager), fVerbose(verbose) {}
+
+// ---------------------------------------------------------------------------------------------------------------------
+//
+std::vector<std::string> L1ConfigRW::GetNodeKeys(const YAML::Node& node) const
+{
+ std::vector<std::string> res;
+ res.reserve(node.size());
+ try {
+ for (const auto& item : node) {
+ res.push_back(item.first.as<std::string>());
+ }
+ }
+ catch (const YAML::InvalidNode& exc) { // It is impossible to
+ LOG(warn) << "L1 config: attempt to call L1ConfigRW::GetNodeKeys for node, which keys could not be represented "
+ << "with strings. An empty vector will be returned";
+ std::vector<std::string>().swap(res);
+ }
+
+ return res;
+}
+
+// ---------------------------------------------------------------------------------------------------------------------
+//
+void L1ConfigRW::ReadCAIterations(const YAML::Node& node)
+{
+ if (node) {
+ if (fVerbose > -1) {
+ LOG(info) << "L1 config: Reading CA tracking iterations sequence. Default iterations will be ignored";
+ }
+ fpInitManager->ClearCAIterations();
+ fpInitManager->SetCAIterationsNumberCrosscheck(node.size());
+ if (fVerbose > 2) { LOG(info) << "L1 config: " << fVerbose << " CA iterations were recorded"; }
+ if (fVerbose > 3) { LOG(info) << "L1 config: Recorded iterations:"; }
+
+ // Loop over input sub-nodes. Each sub-node corresponds to one L1CAIteration object
+ for (const auto& input : node) {
+ try {
+ auto caIter = L1CAIteration(input["name"].as<std::string>());
+ caIter.SetTrackChi2Cut(input["track_chi2_cut"].as<float>());
+ caIter.SetTripletChi2Cut(input["triplet_chi2_cut"].as<float>());
+ caIter.SetDoubletChi2Cut(input["doublet_chi2_cut"].as<float>());
+ caIter.SetPickGather(input["pick_gather"].as<float>());
+ caIter.SetPickNeighbour(input["pick_neighbour"].as<float>());
+ caIter.SetMaxInvMom(1. / input["min_momentum"].as<float>());
+ caIter.SetMaxSlopePV(input["max_slope_pv"].as<float>());
+ caIter.SetMaxSlope(input["max_slope"].as<float>());
+ caIter.SetMaxDZ(input["max_dz"].as<float>());
+ caIter.SetTargetPosSigmaXY(input["target_pos_sigma_x"].as<float>(), input["target_pos_sigma_y"].as<float>());
+ caIter.SetMinLevelTripletStart(input["min_start_triplet_lvl"].as<int>());
+ caIter.SetPrimaryFlag(input["is_primary"].as<bool>());
+ caIter.SetElectronFlag(input["is_electron"].as<bool>());
+ if (fVerbose > 3) { LOG(info) << "L1 config:\n" << caIter.ToString(1); }
+ fpInitManager->PushBackCAIteration(caIter);
+ }
+ catch (const YAML::InvalidNode& exc) {
+ const auto nodeKeys = this->GetNodeKeys(input);
+ const auto nodeKeysStr =
+ std::accumulate(nodeKeys.cbegin(), nodeKeys.cend(), std::string(""),
+ [](std::string lhs, std::string rhs) { return std::move(lhs) + "\n\t" + std::move(rhs); });
+ LOG(fatal) << "L1 config: attempt to access key which does not exist in the configuration file (message from "
+ << "YAML exception: " << exc.what() << "). Defined keys: " << nodeKeysStr;
+ }
+ } // Loop over input sub-nodes: end
+ }
+ else {
+ LOG(warn) << "L1 config: CA tracking iterations sequence was not found in the parameters file, default will be "
+ << "used. To define iterations please use the following path in the YAML file: l1/algorithm/iterations";
+ }
+}
+
+// ---------------------------------------------------------------------------------------------------------------------
+//
+void L1ConfigRW::ReadYaml(const std::string& filename)
+{
+ // Load YAML node from the file
+ YAML::Node config;
+ try {
+ config = YAML::LoadFile(filename);
+ }
+ catch (const YAML::BadFile& exc) {
+ LOG(error) << "L1 config: tracking parameters file \"" << filename << "\" does not exist, default will be used";
+ return;
+ }
+
+ // Tracking iterations
+ this->ReadCAIterations(config["l1"]["algorithm"]["iterations"]);
+}
+
+// ---------------------------------------------------------------------------------------------------------------------
+//
+void L1ConfigRW::WriteYaml(const std::string& /*filename*/) {}
diff --git a/reco/L1/L1Algo/L1ConfigRW.h b/reco/L1/L1Algo/L1ConfigRW.h
new file mode 100644
index 0000000000000000000000000000000000000000..ce83641362a1ad01b78a461c0c60065c4b5323bf
--- /dev/null
+++ b/reco/L1/L1Algo/L1ConfigRW.h
@@ -0,0 +1,56 @@
+/* Copyright (C) 2022 GSI Helmholtzzentrum fuer Schwerionenforschung, Darmstadt
+ SPDX-License-Identifier: GPL-3.0-only
+ Authors: Sergei Zharko [committer] */
+
+/// @file L1ConfigRW.h
+/// @brief Configuration parameter file reader/writer for L1 tracking algorithm (declaration)
+/// @author S.Zharko <s.zharko@gsi.de>
+/// @since 18.07.2022
+
+#ifndef L1ConfigRW_h
+#define L1ConfigRW_h 1
+
+#include <string>
+#include <vector>
+
+namespace YAML
+{
+ class Node;
+}
+class L1InitManager;
+
+/// Class L1ConfigRW provides internal read and write methods for L1 tracking algorithm parameters
+class L1ConfigRW {
+public:
+ /// Constructor
+ /// \param pInitManager Pointer to the L1InitManager instance
+ L1ConfigRW(L1InitManager* pInitManager, int verbose = 0);
+
+ /// Destructor
+ ~L1ConfigRW() {}
+
+ /// Reads parameters from the file
+ /// \param filename Name of the configuration file
+ void ReadYaml(const std::string& filename);
+
+ /// Writes parameters to the file
+ /// \param filename Name of the configuration file
+ void WriteYaml(const std::string& filename);
+
+private:
+ /// Reads CA track finder iterations from YAML node
+ /// \param node YAML node
+ void ReadCAIterations(const YAML::Node& node);
+
+ /// Gets parameters content of the node
+ /// FIXME: We assume, that all of the keys of node can be represented with strings. This assumption is not always true, so
+ /// one should be sure about the assumption when using this function.
+ /// \param node YAML node
+ /// \return Vector of key names
+ std::vector<std::string> GetNodeKeys(const YAML::Node& node) const;
+
+ L1InitManager* fpInitManager = nullptr; ///< Pointer to the L1InitManager instance
+ int fVerbose = 0; ///< Verbosity level
+};
+
+#endif // L1ConfigRW_h
diff --git a/reco/L1/L1Algo/L1Constants.h b/reco/L1/L1Algo/L1Constants.h
index 2f3f86f0a8c73cadeee81a90140fd0855dd8fb15..0e6bca4b4b33b0c4187b88ef17eb11183efb725d 100644
--- a/reco/L1/L1Algo/L1Constants.h
+++ b/reco/L1/L1Algo/L1Constants.h
@@ -49,10 +49,10 @@ namespace L1Constants
/// Flag for calling the CAMergeClones procedure ... TODO
constexpr bool kIfMergeClones {true};
- /// Flag: debug mode for analyzing the doublets pergormance efficiencies
+ /// Flag: debug mode for analyzing the doublets performance efficiencies
//constexpr bool kIfDebugDoubletsPerformance {false};
- /// Flag: debug mode for analyzing the tiplets pergormance efficiencies
+ /// Flag: debug mode for analyzing the tiplets performance efficiencies
//constexpr bool kIfDebugTripletsPerformance {false};
/// Flag: debug mode for creating pools for triplets.
@@ -75,7 +75,7 @@ namespace L1Constants
constexpr int kAlignment {16};
} // end namespace misc
- /// NoInit constants (aliasses)
+ /// NoInit constants (aliases)
namespace noin
{
constexpr float kF {L1NaN::SetNaN<float>()};
diff --git a/reco/L1/L1Algo/L1InitManager.cxx b/reco/L1/L1Algo/L1InitManager.cxx
index 84dc651a2e43f5cd794a3d2860c7123228b7bd40..78f230a6737d1eba1e5e095adb9a680dc0f2d608 100644
--- a/reco/L1/L1Algo/L1InitManager.cxx
+++ b/reco/L1/L1Algo/L1InitManager.cxx
@@ -4,7 +4,7 @@
/************************************************************************************************************
* @file L1InitManager.cxx
- * @bried Input data management class for L1Algo
+ * @brief Input data management class for L1Algo
* @since 19.01.2022
***********************************************************************************************************/
@@ -15,14 +15,14 @@
#include "L1Assert.h"
-//-----------------------------------------------------------------------------------------------------------------------
+// ----------------------------------------------------------------------------------------------------------------------
//
void L1InitManager::AddStation(const L1BaseStationInfo& inStation)
{
// Check if other fields were defined already
// Active detector IDs
L1MASSERT(0, fInitController.GetFlag(EInitKey::kActiveDetectorIDs),
- "Attempt to add a station info before the active detetors set had been initialized");
+ "Attempt to add a station info before the active detectors set had been initialized");
// Number of stations check
L1MASSERT(0, fInitController.GetFlag(EInitKey::kStationsNumberCrosscheck),
@@ -30,7 +30,7 @@ void L1InitManager::AddStation(const L1BaseStationInfo& inStation)
// Field function
L1MASSERT(0, fInitController.GetFlag(EInitKey::kFieldFunction),
- "Attempt to add a station info before the magnetic field function had been intialized");
+ "Attempt to add a station info before the magnetic field function had been initialized");
// Check activeness of this station type
bool isStationActive =
@@ -44,7 +44,7 @@ void L1InitManager::AddStation(const L1BaseStationInfo& inStation)
if (!inStationCopy.GetInitController().GetFlag(L1BaseStationInfo::EInitKey::kThicknessMap)) {
LOG(warn) << "Station material map was not set for detectorID = "
<< static_cast<int>(inStationCopy.GetDetectorID()) << ", stationID = " << inStationCopy.GetStationID()
- << ". Homogenious material budget will be used: " << inStationCopy.GetRadThick()[0];
+ << ". Homogeneous material budget will be used: " << inStationCopy.GetRadThick()[0];
L1Material material;
material.SetBins(1, 100);
material.SetRadThick(0, 0, inStationCopy.GetRadThick()[0]);
@@ -62,7 +62,7 @@ void L1InitManager::AddStation(const L1BaseStationInfo& inStation)
// insert the station in a set
auto insertionResult = fStationsInfo.insert(std::move(inStationCopy));
if (!insertionResult.second) {
- LOG(fatal) << "Attempt to add a dublicating station info object (detectorID = "
+ LOG(fatal) << "Attempt to add a duplicated station info object (detectorID = "
<< static_cast<int>(inStationCopy.GetDetectorID()) << ", stationID = " << inStationCopy.GetStationID()
<< ")";
}
@@ -83,7 +83,7 @@ void L1InitManager::AddStation(const L1BaseStationInfo& inStation)
<< ". Is active: " << isStationActive;
}
-//-----------------------------------------------------------------------------------------------------------------------
+// ----------------------------------------------------------------------------------------------------------------------
//
void L1InitManager::CheckInit()
{
@@ -91,10 +91,24 @@ void L1InitManager::CheckInit()
this->CheckStationsInfoInit();
}
-//-----------------------------------------------------------------------------------------------------------------------
+// ----------------------------------------------------------------------------------------------------------------------
+//
+void L1InitManager::ClearCAIterations()
+{
+ fParameters.fCAIterations.clear();
+ fCAIterationsNumberCrosscheck = -1;
+ fInitController.SetFlag(EInitKey::kCAIterations, false);
+}
+
+// ----------------------------------------------------------------------------------------------------------------------
// NOTE: this function should be called once in the TransferParametersContainer
void L1InitManager::FormParametersContainer()
{
+ // Read configuration file
+ // NOTE: We consider parameters from the configuration file as ones with a higher priority, so all the defined
+ // variables there would be rewritten by the configuration
+ if (fConfigInputName != "") { fConfigRW.ReadYaml(fConfigInputName); }
+
// Check initialization
this->CheckInit();
@@ -127,7 +141,7 @@ void L1InitManager::FormParametersContainer()
fParameters.CheckConsistency();
}
-//-----------------------------------------------------------------------------------------------------------------------
+// ----------------------------------------------------------------------------------------------------------------------
//
void L1InitManager::InitTargetField(double zStep)
{
@@ -139,7 +153,7 @@ void L1InitManager::InitTargetField(double zStep)
// Check for field function
L1MASSERT(0, fInitController.GetFlag(EInitKey::kFieldFunction),
- "Attempt to initialze the field value and field region near target before initializing field function");
+ "Attempt to initialize the field value and field region near target before initializing field function");
// Check for target defined
L1MASSERT(
@@ -168,25 +182,7 @@ void L1InitManager::InitTargetField(double zStep)
fInitController.SetFlag(EInitKey::kPrimaryVertexField);
}
-//-----------------------------------------------------------------------------------------------------------------------
-//
-void L1InitManager::PrintCAIterations(int verbosityLevel) const
-{
- for (const auto& iteration : fParameters.fCAIterations) {
- iteration.Print(verbosityLevel);
- }
-}
-
-//-----------------------------------------------------------------------------------------------------------------------
-//
-void L1InitManager::PrintStations(int verbosityLevel) const
-{
- for (const auto& station : fStationsInfo) {
- station.Print(verbosityLevel);
- }
-}
-
-//-----------------------------------------------------------------------------------------------------------------------
+// ----------------------------------------------------------------------------------------------------------------------
//
void L1InitManager::PushBackCAIteration(const L1CAIteration& iteration)
{
@@ -199,28 +195,28 @@ void L1InitManager::PushBackCAIteration(const L1CAIteration& iteration)
fParameters.fCAIterations.push_back(iteration);
}
-//-----------------------------------------------------------------------------------------------------------------------
+// ----------------------------------------------------------------------------------------------------------------------
//
void L1InitManager::SetActiveDetectorIDs(const L1DetectorIDSet_t& detectorIDs)
{
- // TODO: To think about redifinition possibilities: should it be allowed or not? (S.Zh.)
+ // TODO: To think about redefinition possibilities: should it be allowed or not? (S.Zh.)
fActiveDetectorIDs = detectorIDs;
fInitController.SetFlag(EInitKey::kActiveDetectorIDs);
}
-//-----------------------------------------------------------------------------------------------------------------------
+// ----------------------------------------------------------------------------------------------------------------------
//
void L1InitManager::SetCAIterationsNumberCrosscheck(int nIterations)
{
fCAIterationsNumberCrosscheck = nIterations;
L1Vector<L1CAIteration>& iterationsContainer = fParameters.fCAIterations;
- // NOTE: should be called to prevent multiple copyings of objects between the memory realocations
+ // NOTE: should be called to prevent multiple copies of objects between the memory reallocations
iterationsContainer.reserve(nIterations);
fInitController.SetFlag(EInitKey::kCAIterationsNumberCrosscheck);
}
-//-----------------------------------------------------------------------------------------------------------------------
+// ----------------------------------------------------------------------------------------------------------------------
//
void L1InitManager::SetFieldFunction(const L1FieldFunction_t& fieldFunction)
{
@@ -233,7 +229,7 @@ void L1InitManager::SetFieldFunction(const L1FieldFunction_t& fieldFunction)
}
}
-//-----------------------------------------------------------------------------------------------------------------------
+// ----------------------------------------------------------------------------------------------------------------------
//
void L1InitManager::SetGhostSuppression(int ghostSuppression)
{
@@ -245,7 +241,7 @@ void L1InitManager::SetGhostSuppression(int ghostSuppression)
fInitController.SetFlag(EInitKey::kGhostSuppression);
}
-//-----------------------------------------------------------------------------------------------------------------------
+// ----------------------------------------------------------------------------------------------------------------------
//
void L1InitManager::SetMomentumCutOff(float momentumCutOff)
{
@@ -257,7 +253,7 @@ void L1InitManager::SetMomentumCutOff(float momentumCutOff)
fInitController.SetFlag(EInitKey::kMomentumCutOff);
}
-//-----------------------------------------------------------------------------------------------------------------------
+// ----------------------------------------------------------------------------------------------------------------------
//
void L1InitManager::SetNstations(L1DetectorID detectorID, int nStations)
{
@@ -298,7 +294,7 @@ void L1InitManager::SetNstations(L1DetectorID detectorID, int nStations)
}
}
-//-----------------------------------------------------------------------------------------------------------------------
+// ----------------------------------------------------------------------------------------------------------------------
//
void L1InitManager::SetTargetPosition(double x, double y, double z)
{
@@ -316,7 +312,7 @@ void L1InitManager::SetTargetPosition(double x, double y, double z)
fInitController.SetFlag(EInitKey::kTargetPos);
}
-//-----------------------------------------------------------------------------------------------------------------------
+// ----------------------------------------------------------------------------------------------------------------------
//
void L1InitManager::SetTrackingLevel(int trackingLevel)
{
@@ -328,7 +324,7 @@ void L1InitManager::SetTrackingLevel(int trackingLevel)
fInitController.SetFlag(EInitKey::kTrackingLevel);
}
-//-----------------------------------------------------------------------------------------------------------------------
+// ----------------------------------------------------------------------------------------------------------------------
//
void L1InitManager::TransferParametersContainer(L1Parameters& destination)
{
@@ -341,7 +337,7 @@ void L1InitManager::TransferParametersContainer(L1Parameters& destination)
// INIT CHECKERS
//
-//-----------------------------------------------------------------------------------------------------------------------
+// ----------------------------------------------------------------------------------------------------------------------
//
void L1InitManager::CheckCAIterationsInit()
{
@@ -361,7 +357,7 @@ void L1InitManager::CheckCAIterationsInit()
fInitController.SetFlag(EInitKey::kCAIterations, ifInitPassed);
}
-//-----------------------------------------------------------------------------------------------------------------------
+// ----------------------------------------------------------------------------------------------------------------------
// TODO: REWRITE! and add const qualifier (S.Zharko)
void L1InitManager::CheckStationsInfoInit()
{
diff --git a/reco/L1/L1Algo/L1InitManager.h b/reco/L1/L1Algo/L1InitManager.h
index 8732fa5c4e48e291a87c73510eff1e120ef4e47b..86dad6330c4dbc7592635058cb88599e496a9f53 100644
--- a/reco/L1/L1Algo/L1InitManager.h
+++ b/reco/L1/L1Algo/L1InitManager.h
@@ -4,7 +4,7 @@
/************************************************************************************************************
* @file L1InitManager.h
- * @bried Input data management class for L1Algo
+ * @brief Input data management class for L1Algo
* @since 24.12.2021
***********************************************************************************************************/
#ifndef L1InitManager_h
@@ -12,6 +12,7 @@
#include "L1BaseStationInfo.h"
#include "L1CAIteration.h"
+#include "L1ConfigRW.h"
#include "L1Constants.h"
#include "L1Field.h"
#include "L1ObjectInitController.h"
@@ -28,6 +29,8 @@
#include <numeric>
#include <set>
+class L1ConfigRW;
+
/// Forward declaration of the tracking detectors scoped enumeration. Concrete realization of this enumeration must be
/// determined in the concrete setup class (i.e. CbmL1/BmnL1)
enum class L1DetectorID;
@@ -75,7 +78,7 @@ private:
// this->GetObjectInitController().ToString() method call (S.Zharko)
kActiveDetectorIDs, ///< 0) If the detector sequence is set
kStationsNumberCrosscheck, ///< 1) If the crosscheck station numbers were setup
- kFieldFunction, ///< 2) If magnetic field getter funciton is set
+ kFieldFunction, ///< 2) If magnetic field getter function is set
kTargetPos, ///< 3) If target position was defined
kPrimaryVertexField, ///< 4) If magnetic field value and region defined at primary vertex
kStationsInfo, ///< 5) If all the planned stations were added to the manager
@@ -121,22 +124,28 @@ public:
/// Adds another station of a given type using std::unique_ptr-wraped pointer to L1BaseStationInfo
void AddStation(const std::unique_ptr<L1BaseStationInfo>& puStation) { AddStation(*puStation); }
- /// Provides final checks of large fields initialization calling Check"Object"Init() privat methods,
- /// must be called in the begining of L1Algo::Init()
+ /// Provides final checks of large fields initialization calling Check"Object"Init() private methods,
+ /// must be called in the beginning of L1Algo::Init()
void CheckInit();
// NOTE: This method calls checkers of large fields initializations like a station or an iteration. The method must be
// called in the L1Algo class. (S.Zharko)
+ /// Clears vector of CA track finder iterations
+ void ClearCAIterations();
+
/// Gets ghost suppression flag
int GetGhostSuppression() const { return fGhostSuppression; }
+ /// Gets a name of the input configuration file
+ const std::string& GetInputConfigName() const { return fConfigInputName; }
+
/// Gets momentum cutoff
float GetMomentumCutOff() const { return fMomentumCutOff; }
/// Gets a const reference to L1ObjectInitController
const InitController_t& GetInitController() const { return fInitController; }
- /// Gets a pointer to L1Parameters instance with a posibility of its fields modification
+ /// Gets a pointer to L1Parameters instance with a possibility of its fields modification
//const L1Parameters* GetParameters() const { return fpParameters; }
/// Gets total number of active stations
@@ -157,20 +166,17 @@ public:
return fNstationsGeometry[static_cast<L1DetectorID_t>(detectorID)];
}
+ /// Gets a name of the output configuration file
+ const std::string& GetOutputConfigName() const { return fConfigOutputName; }
+
/// Gets tracking level
int GetTrackingLevel() const { return fTrackingLevel; }
/// Calculates L1FieldValue and L1FieldReference values for a selected step in z-axis from the target position
/// \param zStep step between nodal points
- // TODO: Consider posibility for linear approximation (S.Zh.)
+ // TODO: Consider possibility for linear approximation (S.Zh.)
void InitTargetField(double zStep);
- /// Prints a list of CA track finder iterations
- void PrintCAIterations(int verbosityLevel = 0) const;
-
- /// Prints a list of stations
- void PrintStations(int verbosityLevel = 0) const;
-
/// Pushes an CA track finder iteration into a sequence of iteration using reference
void PushBackCAIteration(const L1CAIteration& iteration);
@@ -189,35 +195,49 @@ public:
/// Sets a magnetic field function, which will be applied for all the stations
void SetFieldFunction(const L1FieldFunction_t& fieldFcn);
- ///
+ /// FIXME: ...
void SetGhostSuppression(int ghostSuppression);
- ///
- void SetMomentumCutOff(float momentumCutOff);
+ /// Sets a name of the input configuration file. If the file is undefined, default settings will be used. If the file is
+ /// defined, the default settings will be overwritten with
+ /// \param filename Name of the input L1 parameters configuration
+ void SetInputConfigName(const std::string& filename) { fConfigInputName = filename; }
///
- void SetTrackingLevel(int trackingLevel);
+ void SetMomentumCutOff(float momentumCutOff);
/// Sets a number of actual stations for a particular tracking detector ID to provide initialization cross-check
void SetNstations(L1DetectorID detectorID, int nStations);
- /// Sets target poisition
+ /// Sets a name of the output configuration file. The output file is created from the fields, saved in the resulted
+ /// L1Parameters object
+ /// \param filename Name of the output L1 parameters configuration
+ void SetOutputConfigName(const std::string& filename) { fConfigOutputName = filename; }
+
+ /// Sets target position
+ /// \param x Position X component [cm]
+ /// \param y Position Y component [cm]
+ /// \param z Position Z component [cm]
void SetTargetPosition(double x, double y, double z);
+ ///
+ void SetTrackingLevel(int trackingLevel);
+
/// Transfers L1Parameters object to the destination
+ /// \param destination Reference to the destination of the L1 object
void TransferParametersContainer(L1Parameters& destination);
/// ***************************
/// ** Flags for development **
/// ***************************
- /// Ignopre hit search areas
+ /// Ignore hit search areas
void DevSetIgnoreHitSearchAreas(bool value = true) { fParameters.fDevIsIgnoreHitSearchAreas = value; }
/// Start singlets fit at the target
void DevSetFitSingletsFromTarget(bool value = true) { fParameters.fDevIsFitSingletsFromTarget = value; }
- /// Flag to match doublets using Mc information
+ /// Flag to match doublets using MC information
void DevSetIsMatchDoubletsViaMc(bool value = true) { fParameters.fDevIsMatchDoubletsViaMc = value; }
private:
@@ -246,16 +266,16 @@ private:
std::set<L1BaseStationInfo> fStationsInfo {}; ///< Set of L1BaseStationInfo objects
/// Numbers of stations, which are active in tracking. Index of an array element (except the last one) corresponds to a given
- /// L1DetectorID of the detector subystem. The last array element corresponds to the total number of stations.
+ /// L1DetectorID of the detector subsystem. The last array element corresponds to the total number of stations.
std::array<int, L1Constants::size::kMaxNdetectors + 1> fNstationsActive {};
/// Actual numbers of stations, provided by geometry. Index of an array element (except the last one) corresponds to a given
- /// L1DetectorID of the detector subystem. The last array element corresponds to the total number of stations.
+ /// L1DetectorID of the detector subsystem. The last array element corresponds to the total number of stations.
std::array<int, L1Constants::size::kMaxNdetectors + 1> fNstationsGeometry {};
/// A function which returns magnetic field vector B in a radius-vector xyz
L1FieldFunction_t fFieldFunction {[](const double (&)[3], double (&)[3]) {}};
- // NOTE: Stations of daetectors which will not be assigned as active, will not be included in the tracking!
+ // NOTE: Stations of detectors which will not be assigned as active, will not be included in the tracking!
// * CA track finder iterations related *
@@ -265,7 +285,12 @@ private:
int fTrackingLevel {0}; ///< tracking level
int fGhostSuppression {0}; ///< flag: if true, ghost tracks are suppressed
- float fMomentumCutOff {0}; ///< minimum momentum of tracks TODO: ?
+ float fMomentumCutOff {0}; ///< minimum momentum of tracks
+
+ // * Configuration related *
+ std::string fConfigInputName {""}; ///< name for the input configuration file
+ std::string fConfigOutputName {""}; ///< name for the output configuration file
+ L1ConfigRW fConfigRW {this, /*verbose = */ 4}; ///< configuration file reader and writer
};
#endif
diff --git a/reco/L1/L1Algo/L1Parameters.h b/reco/L1/L1Algo/L1Parameters.h
index 8dbf3582b3670c16257a8a71661fee42ae059f01..d8476c0b16bacd0a0530ba2de73ac1129eb892c9 100644
--- a/reco/L1/L1Algo/L1Parameters.h
+++ b/reco/L1/L1Algo/L1Parameters.h
@@ -174,15 +174,14 @@ public:
/// Is singlets fit starts at the target
bool DevIsFitSingletsFromTarget() const { return fDevIsFitSingletsFromTarget; }
- /// Flag to match doublets using Mc information
+ /// Flag to match doublets using MC information
bool DevIsMatchDoubletsViaMc() const { return fDevIsMatchDoubletsViaMc; }
private:
unsigned int fMaxDoubletsPerSinglet {150}; ///< Upper-bound cut on max number of doublets per one singlet
unsigned int fMaxTripletPerDoublets {15}; ///< Upper-bound cut on max number of triplets per one doublet
- alignas(L1Constants::misc::kAlignment)
- L1IterationsContainer_t fCAIterations {}; ///< L1 Track finder iterations vector
+ alignas(L1Constants::misc::kAlignment) L1IterationsContainer_t fCAIterations {}; ///< L1 tracking iterations vector
/*************************
** Geometry parameters **
@@ -204,14 +203,14 @@ private:
alignas(L1Constants::misc::kAlignment) L1MaterialContainer_t fThickMap {};
/// Numbers of stations, which are active in tracking. Index of an array element (except the last one) corresponds to a given
- /// L1DetectorID of the detector subystem. The last array element corresponds to the total number of stations.
+ /// L1DetectorID of the detector subsystem. The last array element corresponds to the total number of stations.
alignas(L1Constants::misc::kAlignment) std::array<int, (L1Constants::size::kMaxNdetectors + 1)> fNstationsActive {};
/// Actual numbers of stations, provided by geometry. Index of an array element (except the last one) corresponds to a given
- /// L1DetectorID of the detector subystem. The last array element corresponds to the total number of stations.
+ /// L1DetectorID of the detector subsystem. The last array element corresponds to the total number of stations.
alignas(L1Constants::misc::kAlignment) std::array<int, (L1Constants::size::kMaxNdetectors + 1)> fNstationsGeometry {};
- /// Map of the actual detector indeces to the active detector indeces
+ /// Map of the actual detector indexes to the active detector indexes
/// The vector maps actual station index (which is defined by ) to the index of station in tracking. If the station is inactive,
/// its index is equal to -1. Example: let stations 1 and 4 be inactive. Then:
/// actual index: 0 1 2 3 4 5 6 7 8 9 0 0 0 0
@@ -226,13 +225,12 @@ private:
bool fDevIsFitSingletsFromTarget {false}; ///< Fit singlet starting from the target with the KF
- bool fDevIsMatchDoubletsViaMc {false}; ///< Flag to match doublets using Mc information
+ bool fDevIsMatchDoubletsViaMc {false}; ///< Flag to match doublets using MC information
/********************************
** Friend classes declaration **
********************************/
- /// Note provide friend access to L1InitManager (TODO: L1Algo auch?)
friend class L1InitManager;
};
diff --git a/reco/L1/qa/CbmTrackingInputQaSts.cxx b/reco/L1/qa/CbmTrackingInputQaSts.cxx
index 29a546ef9d45823b40ce74d01a269cbe83815718..2ffee72734fd26dc6dd02c9ac502c118d3653d7a 100644
--- a/reco/L1/qa/CbmTrackingInputQaSts.cxx
+++ b/reco/L1/qa/CbmTrackingInputQaSts.cxx
@@ -11,23 +11,23 @@
#include "CbmTrackingInputQaSts.h"
-#include "CbmMCDataArray.h"
-#include "CbmStsCluster.h"
#include "CbmDigiManager.h"
+#include "CbmMCDataArray.h"
#include "CbmMCEventList.h"
+#include "CbmMCTrack.h"
#include "CbmMatch.h"
-#include "CbmStsHit.h"
#include "CbmStsAddress.h" // TODO: TMP for tests, must be removed!!!! (S.Zharko)
-#include "CbmTimeSlice.h"
-#include "CbmMCTrack.h"
+#include "CbmStsCluster.h"
+#include "CbmStsHit.h"
#include "CbmStsPoint.h"
#include "CbmStsTrackingInterface.h" // Communicate via tracking detector interface
+#include "CbmTimeSlice.h"
#include "FairLogger.h"
#include "TClonesArray.h"
-#include "TParticlePDG.h"
#include "TDatabasePDG.h"
+#include "TParticlePDG.h"
#include <iostream>
#include <stdexcept>
@@ -36,49 +36,43 @@ ClassImp(CbmTrackingInputQaSts);
// ---------------------------------------------------------------------------------------------------------------------
//
-CbmTrackingInputQaSts::CbmTrackingInputQaSts(int verbosity) : FairTask("CbmTrackingInputQaSts", verbosity)
-{
-}
+CbmTrackingInputQaSts::CbmTrackingInputQaSts(int verbosity) : FairTask("CbmTrackingInputQaSts", verbosity) {}
// ---------------------------------------------------------------------------------------------------------------------
//
-CbmTrackingInputQaSts::~CbmTrackingInputQaSts()
-{
- DeInit();
-}
+CbmTrackingInputQaSts::~CbmTrackingInputQaSts() { DeInit(); }
// ---------------------------------------------------------------------------------------------------------------------
//
bool CbmTrackingInputQaSts::CheckDistributions()
{
+ std::cout << "CALL CbmTrackingInputQaSts::CheckDistributions()\n";
+
bool res = true;
-
+
const int nStations = fpDetectorInterface->GetNtrackingStations();
// ** Check pulls **
//
-
-
// Function to compare pulls distributions RMS with 1
// TODO: Choose proper selection criteria (S.Zharko)
auto CheckPullsDistribution = [&](const TH1& pHist) {
// Checks sigma of distribution fit with unity
auto* pFitFunc = pHist.GetFunction("gaus");
if (!pFitFunc) {
- throw std::runtime_error("STS tracking input QA: attempt to check sigma of unfitted pulls histogram");
+ throw std::runtime_error(TString("STS tracking input QA: attempt to check sigma of histogram \"")
+ + pHist.GetName() + "\" with undefined fit");
}
auto vSigma = pFitFunc->GetParameter(2);
auto eSigma = pFitFunc->GetParError(2);
// Select 3 sigma interval
- LOG(info) << "Checking histogram \"" << pHist.GetName() << "\": fit result = " << vSigma << " +/- " << eSigma;
- if (std::fabs(vSigma - 1.) < 3. * eSigma) {
- return true;
- }
+ LOG(info) << "Checking histogram \"" << pHist.GetName() << "\": fit result = " << vSigma << " +/- " << eSigma;
+ if (std::fabs(vSigma - 1.) < 3. * eSigma) { return true; }
else {
return false;
}
-
+
//auto vRms = pHist.GetRMS(); // RMS value
//auto eRms = pHist.GetRMSError(); // RMS error
//std::cout << std::fabs(vRms - 1.) << '\n';
@@ -120,14 +114,14 @@ void CbmTrackingInputQaSts::DeInit()
fpDigiManager = nullptr;
fpMcManager = nullptr;
- fpMcTracks = nullptr;
- fpMcPoints = nullptr;
+ fpMcTracks = nullptr;
+ fpMcPoints = nullptr;
- fpMcEventList = nullptr;
- fpClusters = nullptr;
- fpHits = nullptr;
+ fpMcEventList = nullptr;
+ fpClusters = nullptr;
+ fpHits = nullptr;
fpClusterMatches = nullptr;
- fpHitMatches = nullptr;
+ fpHitMatches = nullptr;
fOutFolder.Clear();
fpOutFolderHists = nullptr;
@@ -163,25 +157,44 @@ void CbmTrackingInputQaSts::Exec(Option_t*)
// --------------------------------------------------------------------------------------------------------------------
//
void CbmTrackingInputQaSts::Finish()
-{
+{
+ std::cout << "CALL CbmTrackingInputQaSts::Finish()\n";
+
+ // Fit histograms
+ this->FitHistograms();
+
// Add output to a sink
auto* pSink = FairRootManager::Instance()->GetSink();
- if (pSink) {
- pSink->WriteObject(&GetQa(), nullptr);
- }
+ if (pSink) { pSink->WriteObject(&GetQa(), nullptr); }
- // Check accumulated distributions
+ // Check accumulated distributions
bool areResolutionsOk = CheckDistributions();
// TODO: Collect all the flags in one place and make a decission here (S.Zharko)
- if (areResolutionsOk) {
- LOG(info) << this->GetName() << ": \033[1;32mtask succeeded\033[0m";
- }
+ if (areResolutionsOk) { LOG(info) << this->GetName() << ": \033[1;32mtask succeeded\033[0m"; }
else {
- LOG(info) << this->GetName() << ": \033[1;31mtask failed\033[0m";
+ LOG(info) << this->GetName() << ": \033[1;31mtask failed\033[0m";
}
+}
+
+// ---------------------------------------------------------------------------------------------------------------------
+//
+void CbmTrackingInputQaSts::FitHistograms()
+{
+ std::cout << "CALL CbmTrackingInputQaSts::FitHistograms()\n";
+
+ const int nStations = fpDetectorInterface->GetNtrackingStations();
+ for (int iSt = 0; iSt < nStations; ++iSt) {
+ // Fit histograms
+ fHistResidualX[iSt].Fit("gaus", "Q");
+ fHistResidualY[iSt].Fit("gaus", "Q");
+ fHistResidualT[iSt].Fit("gaus", "Q");
+ fHistPullX[iSt].Fit("gaus", "Q");
+ fHistPullY[iSt].Fit("gaus", "Q");
+ fHistPullT[iSt].Fit("gaus", "Q");
+ }
}
// ---------------------------------------------------------------------------------------------------------------------
@@ -198,19 +211,10 @@ InitStatus CbmTrackingInputQaSts::GeometryQa()
TFolder& CbmTrackingInputQaSts::GetQa()
{
//gStyle->SetPaperSize(20, 20);
-
+
// Loop over tracking stations
const int nStations = fpDetectorInterface->GetNtrackingStations();
for (int iSt = 0; iSt < nStations; ++iSt) {
- // Fit histograms
- fHistResidualX[iSt].Fit("gaus", "Q");
- fHistResidualY[iSt].Fit("gaus", "Q");
- fHistResidualT[iSt].Fit("gaus", "Q");
-
- fHistPullX[iSt].Fit("gaus", "Q");
- fHistPullY[iSt].Fit("gaus", "Q");
- fHistPullT[iSt].Fit("gaus", "Q");
-
// Draw histograms
fCanvResidualX.cd(iSt + 1);
fHistResidualX[iSt].DrawCopy("", "");
@@ -237,7 +241,7 @@ TFolder& CbmTrackingInputQaSts::GetQa()
fCanvEfficiencyXY.cd(iSt + 1);
fHistEfficiencyXY[iSt].DrawCopy("colz", "");
- }// Loop over tracking stations: end
+ } // Loop over tracking stations: end
return fOutFolder;
}
@@ -267,14 +271,14 @@ InitStatus CbmTrackingInputQaSts::Init()
fpOutFolderHists = fOutFolder.AddFolder("rawHist", "Raw histograms");
gStyle->SetOptStat(0);
- fNevents.SetVal(0); // redundant
+ fNevents.SetVal(0); // redundant
fpOutFolderHists->Add(&fNevents);
// Register histograms in the output folder
for (auto* histo : fHistList) {
fpOutFolderHists->Add(histo);
}
-
+
// Do QA of geometry in the end of initialization
return GeometryQa();
}
@@ -313,7 +317,7 @@ InitStatus CbmTrackingInputQaSts::InitCanvases()
fCanvEfficiencyR.Divide2D(fpDetectorInterface->GetNtrackingStations());
fCanvEfficiencyXY.Divide2D(fpDetectorInterface->GetNtrackingStations());
-
+
// Add the canvases to the output folder
fOutFolder.Add(&fCanvResidualX);
fOutFolder.Add(&fCanvResidualY);
@@ -325,10 +329,10 @@ InitStatus CbmTrackingInputQaSts::InitCanvases()
fOutFolder.Add(&fCanvPointsPerHit);
fOutFolder.Add(&fCanvHitsPerPoint);
-
+
fOutFolder.Add(&fCanvEfficiencyR);
fOutFolder.Add(&fCanvEfficiencyXY);
-
+
return kSUCCESS;
}
@@ -337,7 +341,7 @@ InitStatus CbmTrackingInputQaSts::InitCanvases()
InitStatus CbmTrackingInputQaSts::InitHistograms()
{
gStyle->SetOptStat();
-
+
fHistResidualX.clear();
fHistResidualY.clear();
fHistResidualT.clear();
@@ -348,7 +352,7 @@ InitStatus CbmTrackingInputQaSts::InitHistograms()
fHistPointsPerHit.clear();
fHistHitsPerPoint.clear();
-
+
fHistEfficiencyR.clear();
fHistEfficiencyXY.clear();
@@ -418,24 +422,24 @@ InitStatus CbmTrackingInputQaSts::InitHistograms()
histTitle = Form("STS: Station %d: Hit per MC points; N_{hits} per MC point", iSt);
std::tie(nBins, rMin, rMax) = fRangeHitsPerPoint;
fHistHitsPerPoint.emplace_back(histName, histTitle, nBins, rMin, rMax);
-
+
// Efficiencies
- double xMax = fpDetectorInterface->GetXmax(iSt); /// TODO: test ranges
+ double xMax = fpDetectorInterface->GetXmax(iSt); /// TODO: test ranges
double yMax = fpDetectorInterface->GetYmax(iSt); /// TODO: test ranges
-
+
// Efficiency vs. distance of point from center
- histName = Form("pr1_sts_EfficiencyR_st%d", iSt);
- histTitle = Form("STS: Station %d: Efficiency R; R [cm]", iSt);
- rMin = 0.;
- rMax = sqrt(xMax * xMax + yMax * yMax);
- nBins = 100;
+ histName = Form("pr1_sts_EfficiencyR_st%d", iSt);
+ histTitle = Form("STS: Station %d: Efficiency R; R [cm]", iSt);
+ rMin = 0.;
+ rMax = sqrt(xMax * xMax + yMax * yMax);
+ nBins = 100;
fHistEfficiencyR.emplace_back(histName, histTitle, nBins, rMin, rMax);
fHistEfficiencyR[iSt].SetOptStat(1110);
- // Efficiency vs. x and y
- histName = Form("pr1_sts_EfficiencyXY_st%d", iSt);
- histTitle = Form("STS: Station %d: Efficiency XY; X [cm]; Y [cm]", iSt);
- nBins = 50;
+ // Efficiency vs. x and y
+ histName = Form("pr1_sts_EfficiencyXY_st%d", iSt);
+ histTitle = Form("STS: Station %d: Efficiency XY; X [cm]; Y [cm]", iSt);
+ nBins = 50;
fHistEfficiencyXY.emplace_back(histName, histTitle, nBins, -xMax, xMax, nBins, -yMax, yMax);
fHistEfficiencyXY[iSt].SetOptStat(10);
}
@@ -472,7 +476,7 @@ CbmMatch CbmTrackingInputQaSts::MatchHits(const CbmStsHit* pHit, int iHit)
LOG(error) << "STD: hit (id = " << iHit << ") has incorrect front cluster index: " << iClusterF;
throw std::runtime_error("STS tracking input QA: wrong front cluster index");
}
-
+
const int iClusterB = pHit->GetBackClusterId();
if (iClusterB < 0) {
LOG(error) << "STD: hit (id = " << iHit << ") has incorrect back cluster index: " << iClusterF;
@@ -491,7 +495,7 @@ CbmMatch CbmTrackingInputQaSts::MatchHits(const CbmStsHit* pHit, int iHit)
LOG(error) << "STD: back cluster does not exist for hit (id = " << iHit << ')';
throw std::runtime_error("STS tracking input QA: back cluster not found");
}
-
+
const CbmMatch* pClusterMatchF = dynamic_cast<const CbmMatch*>(fpClusterMatches->At(iClusterF));
const CbmMatch* pClusterMatchB = dynamic_cast<const CbmMatch*>(fpClusterMatches->At(iClusterB));
@@ -519,7 +523,7 @@ CbmMatch CbmTrackingInputQaSts::MatchHits(const CbmStsHit* pHit, int iHit)
}
}
- return res; // Rely on NRVO
+ return res; // Rely on NRVO
}
// ---------------------------------------------------------------------------------------------------------------------
@@ -529,25 +533,25 @@ double CbmTrackingInputQaSts::ParticleMass(int pdg)
if (fabs(pdg) < 9999999 && ((TParticlePDG*) TDatabasePDG::Instance()->GetParticle(pdg))) {
return TDatabasePDG::Instance()->GetParticle(pdg)->Mass();
}
-
+
constexpr double kAmuToGev {0.93149410242}; // 1 amu in GeV/c2
-
+
// Define masses of several light nuclei by hands
- if (fabs(pdg) == 1000010020) {
- return kAmuToGev * 2.01410177812; // deutron
- }
- else if (fabs(pdg) == 1000010030) {
- return kAmuToGev * 3.01604927790; // triton
- }
- else if (fabs(pdg) == 1000020030) {
- return kAmuToGev * 3.01602932007; // He-3
- }
- else if (fabs(pdg) == 1000020040) {
- return kAmuToGev * 4.00260325413; // He-4
+ if (fabs(pdg) == 1000010020) {
+ return kAmuToGev * 2.01410177812; // deutron
+ }
+ else if (fabs(pdg) == 1000010030) {
+ return kAmuToGev * 3.01604927790; // triton
+ }
+ else if (fabs(pdg) == 1000020030) {
+ return kAmuToGev * 3.01602932007; // He-3
+ }
+ else if (fabs(pdg) == 1000020040) {
+ return kAmuToGev * 4.00260325413; // He-4
}
-
+
LOG(warn) << "\033[1;31m Found mass for pdg = " << pdg << " is undefined. "
- << "Please, provide data for this particle\033[0m\n";
+ << "Please, provide data for this particle\033[0m\n";
return 0.;
}
@@ -557,9 +561,7 @@ double CbmTrackingInputQaSts::ParticleMass(int pdg)
InitStatus CbmTrackingInputQaSts::ReadAndCreateDataBranches()
{
auto* pManager = FairRootManager::Instance();
- if (!pManager) {
- LOG(fatal) << "FairRootManager was not found";
- }
+ if (!pManager) { LOG(fatal) << "FairRootManager was not found"; }
fpDigiManager = CbmDigiManager::Instance();
fpDigiManager->Init(); // NOTE: is initialized only once
@@ -574,7 +576,8 @@ InitStatus CbmTrackingInputQaSts::ReadAndCreateDataBranches()
return kERROR;
}
- fpHits = dynamic_cast<TClonesArray*>(FairRootManager::Instance()->GetObject("StsHit"));;
+ fpHits = dynamic_cast<TClonesArray*>(FairRootManager::Instance()->GetObject("StsHit"));
+ ;
if (!fpHits) {
LOG(error) << "Hits input array was not found for STS";
return kERROR;
@@ -627,21 +630,21 @@ void CbmTrackingInputQaSts::ResolutionQa()
{
// Resolution QA is impossible without MC information
if (!fIsMcPresent) { return; }
-
+
const int nHits = fpHits->GetEntriesFast(); // Number of hits in event/ts
const int nStations = fpDetectorInterface->GetNtrackingStations(); // Number of stations in event/ts
const int nMcEvents = fpMcEventList->GetNofEvents(); // Number of MC events in event/ts
-
- std::vector<std::vector<int>> nHitsPerMcPoint; // Number of hits per one MC point in different MC events
+
+ std::vector<std::vector<int>> nHitsPerMcPoint; // Number of hits per one MC point in different MC events
nHitsPerMcPoint.resize(nMcEvents);
-
+
// ** Loop over MC events within event/ts **
for (int iE = 0; iE < nMcEvents; ++iE) {
const int fileId = fpMcEventList->GetFileIdByIndex(iE);
const int eventId = fpMcEventList->GetEventIdByIndex(iE);
const int nPoints = fpMcPoints->Size(fileId, eventId);
nHitsPerMcPoint[iE].resize(nPoints, 0);
- } // Loop over MC events within event/ts: end
+ } // Loop over MC events within event/ts: end
// ** Loop over hits within event/ts **
@@ -651,23 +654,23 @@ void CbmTrackingInputQaSts::ResolutionQa()
LOG(error) << "STS hit with index " << iHit << " does not exist in \"StsHit\" array";
throw std::runtime_error("STS tracking input QA: hit does not exist");
}
-
+
// Check station index of a hit
const int iStation = fpDetectorInterface->GetTrackingStationIndex(pHit);
if (iStation < 0 || iStation >= nStations) {
- LOG(error) << "STS hit with index " << iHit << ": tracking station index = "
- << iStation << " is out of range [0, " << nStations << ']';
+ LOG(error) << "STS hit with index " << iHit << ": tracking station index = " << iStation
+ << " is out of range [0, " << nStations << ']';
throw std::runtime_error("STS tracking input QA: hit has inconsistent station index");
}
// Get custom match for the hit
- const CbmMatch hitMatch = MatchHits(pHit, iHit); // throws std::runtime_error
+ const CbmMatch hitMatch = MatchHits(pHit, iHit); // throws std::runtime_error
// Fill number of hits per one MC point vector and update the number of hits per a given MC point
int nMcPoints = 0; // Number of non-noisy MC points per event
for (int iLink = 0; iLink < hitMatch.GetNofLinks(); ++iLink) {
const auto& link = hitMatch.GetLink(iLink);
- if (link.GetIndex() >= 0) { // Select only non-noisy links (non-noisy digis)
+ if (link.GetIndex() >= 0) { // Select only non-noisy links (non-noisy digis)
++nMcPoints;
const int iE = fpMcEventList->GetEventIndex(link);
if (iE < 0 || iE >= nMcEvents) {
@@ -686,23 +689,23 @@ void CbmTrackingInputQaSts::ResolutionQa()
// Select only hits with one MC point
if (nMcPoints != 1) { continue; }
-
+
// Take the best link in the match
const auto& bestLink = hitMatch.GetMatchedLink();
- // Skip noise
+ // Skip noise
if (bestLink.GetIndex() < 0) { continue; }
// Get MC point
const auto* pPoint = dynamic_cast<const CbmStsPoint*>(fpMcPoints->Get(bestLink));
if (!pPoint) {
LOG(error) << "STS: MC point does not exist (iHit = " << iHit << ')';
- throw std::runtime_error("STS tracking input QA: MC point does not exist"); // NOTE:
+ throw std::runtime_error("STS tracking input QA: MC point does not exist"); // NOTE:
}
// TODO: Study the distributions of |pointZIn - stationZ|, |pointZout - stationZ| and |hitZ - stationZ|
//const double stationZ = fpDetectorInterface->GetZ(iStation);
-
+
// ***********************************
// ** Calculate residuals and pools **
// ***********************************
@@ -713,21 +716,21 @@ void CbmTrackingInputQaSts::ResolutionQa()
LOG(error) << "STS: MC event time is undefined (iHit = " << iHit << ')';
throw std::runtime_error("STS tracking input QA: the MC event time is undefined");
}
-
+
// Get time, space position and momenta components for the MC point (values are taken for the "In" point)
- double mcX = pPoint->GetXIn(); // [cm]
- double mcY = pPoint->GetYIn(); // [cm]
- double mcZ = pPoint->GetZIn(); // [cm]
- double mcT = pPoint->GetTime() + t0; // [ns]
- double mcPx = pPoint->GetPx(); // [GeV/c]
- double mcPy = pPoint->GetPy(); // [GeV/c]
- double mcPz = pPoint->GetPz(); // [GeV/c]
+ double mcX = pPoint->GetXIn(); // [cm]
+ double mcY = pPoint->GetYIn(); // [cm]
+ double mcZ = pPoint->GetZIn(); // [cm]
+ double mcT = pPoint->GetTime() + t0; // [ns]
+ double mcPx = pPoint->GetPx(); // [GeV/c]
+ double mcPy = pPoint->GetPy(); // [GeV/c]
+ double mcPz = pPoint->GetPz(); // [GeV/c]
// Skip slow particles (TODO: why pZ, not p?
if (fabs(pPoint->GetPzOut()) < fMinMomentum) { continue; }
// Difference between z components of MC in point and the hit
- double dz = pHit->GetZ() - mcZ; // [cm]
+ double dz = pHit->GetZ() - mcZ; // [cm]
// Propagate MC-point x and y "In" coordinates to z-plane of the hit
mcX += dz * mcPx / mcPz;
@@ -736,14 +739,14 @@ void CbmTrackingInputQaSts::ResolutionQa()
// Propagete MC-point "In" time to z-plane of the hit
int pdgCode = pPoint->GetPid();
double mass = ParticleMass(pdgCode);
- constexpr double speedOfLight {29.9792458}; // cm/ns
+ constexpr double speedOfLight {29.9792458}; // cm/ns
TVector3 mom;
-
+
mcT += dz / (mcPz * speedOfLight) * sqrt(mass * mass + mom.Mag2());
- double dx = pHit->GetX() - mcX;
- double dy = pHit->GetY() - mcY;
- double dt = pHit->GetTime() - mcT;
+ double dx = pHit->GetX() - mcX;
+ double dy = pHit->GetY() - mcY;
+ double dt = pHit->GetTime() - mcT;
double rmsX = pHit->GetDx();
double rmsY = pHit->GetDy();
double rmsT = pHit->GetTimeError();
@@ -754,9 +757,9 @@ void CbmTrackingInputQaSts::ResolutionQa()
fHistPullX[iStation].Fill(dx / rmsX);
fHistPullY[iStation].Fill(dy / rmsY);
- fHistPullT[iStation].Fill(dt / rmsT);
- }// Loop over hits within event/ts: end
-
+ fHistPullT[iStation].Fill(dt / rmsT);
+ } // Loop over hits within event/ts: end
+
// *********************
// ** Fill efficiency **
// *********************
@@ -766,7 +769,7 @@ void CbmTrackingInputQaSts::ResolutionQa()
const int fileId = fpMcEventList->GetFileIdByIndex(iE);
const int eventId = fpMcEventList->GetEventIdByIndex(iE);
const int nPoints = fpMcPoints->Size(fileId, eventId);
-
+
// ** Loop over MC points **
for (int iP = 0; iP < nPoints; ++iP) {
const auto* pPoint = dynamic_cast<CbmStsPoint*>(fpMcPoints->Get(fileId, eventId, iP));
@@ -774,22 +777,22 @@ void CbmTrackingInputQaSts::ResolutionQa()
LOG(error) << "MC point does not exist for iE = " << iE << ", iP = " << iP;
throw std::runtime_error("STS tracking input QA: MC point does not exist");
}
-
+
int address = pPoint->GetDetectorID();
int iStation = fpDetectorInterface->GetTrackingStationIndex(address);
if (iStation < 0 || iStation >= nStations) {
- LOG(error) << "STS MC point with index " << iP << ": tracking station index = "
- << iStation << " is out of range [0, " << nStations << ']';
+ LOG(error) << "STS MC point with index " << iP << ": tracking station index = " << iStation
+ << " is out of range [0, " << nStations << ']';
throw std::runtime_error("STS tracking input QA: MC point has inconsistent station index");
}
fHistHitsPerPoint[iStation].Fill(nHitsPerMcPoint[iE][iP]);
-
- double pointDistance = sqrt(pPoint->GetXIn() * pPoint->GetXIn() + pPoint->GetYIn() * pPoint->GetYIn()); // [cm]
+
+ double pointDistance = sqrt(pPoint->GetXIn() * pPoint->GetXIn() + pPoint->GetYIn() * pPoint->GetYIn()); // [cm]
fHistEfficiencyR[iStation].Fill(pointDistance, (nHitsPerMcPoint[iE][iP] > 0));
fHistEfficiencyXY[iStation].Fill(pPoint->GetXIn(), pPoint->GetYIn(), (nHitsPerMcPoint[iE][iP] > 0));
}
- } // Loop over MC events within event/ts: end
+ } // Loop over MC events within event/ts: end
}
@@ -799,4 +802,3 @@ void CbmTrackingInputQaSts::SetParContainers()
{
std::cout << "\033[1;32mCALL:\033[0mCbmTrackingInputQaSts::SetParContainers()\n";
}
-
diff --git a/reco/L1/qa/CbmTrackingInputQaSts.h b/reco/L1/qa/CbmTrackingInputQaSts.h
index 7a0b324baa62fd71f1e81440076c36d43bc91585..c022f97b3290954a249037a9732c34cd5f15512f 100644
--- a/reco/L1/qa/CbmTrackingInputQaSts.h
+++ b/reco/L1/qa/CbmTrackingInputQaSts.h
@@ -27,9 +27,9 @@
#include "TParameter.h"
#include "TProfile.h"
#include "TString.h"
-#include <tuple>
#include <memory>
+#include <tuple>
#include <vector>
class CbmDigiManager;
@@ -58,7 +58,7 @@ public:
CbmTrackingInputQaSts(CbmTrackingInputQaSts&&) = delete;
CbmTrackingInputQaSts& operator=(const CbmTrackingInputQaSts&) = delete;
CbmTrackingInputQaSts& operator=(CbmTrackingInputQaSts&&) = delete;
-
+
/// TTask: Executes the task (processes a timeslice)
/// \param option Specify action for the TTask::Exec (not defined yet)
void Exec(Option_t* /*option*/);
@@ -66,19 +66,22 @@ public:
/// FairTask: Action at hte end of the run
void Finish();
+ /// Fits histograms
+ void FitHistograms();
+
/// Gets maximum allowed distance between z-components of hit/MC-point position and the station center [cm]
double GetMaxDistanceZ() const { return fMaxDistanceZ; }
-
+
/// Gets maximum allowed difference between pulls distributions RMS and the expected RMS value of RMS = 1
double GetMaxPullsRmsDiff() const { return fMaxPullsRmsDiff; }
-
+
/// Gets minimum particle momentum
double GetMinMomentum() const { return fMinMomentum; }
-
- /// Fits and draws histograms
+
+ /// Draws histograms
TFolder& GetQa();
-
- /// FairTask: Task initialization an the beginning of the run
+
+ /// FairTask: Task initialization in the beginning of the run
InitStatus Init();
/// FairTask: Task reinitialization
@@ -102,7 +105,7 @@ private:
// ***************
// ** Functions **
// ***************
-
+
/// Checks accumulated distributions
bool CheckDistributions();
@@ -129,12 +132,12 @@ private:
/// \param iHit An index of the hit in the hits array
/// \return Match object
CbmMatch MatchHits(const CbmStsHit* pHit, int iHit);
-
+
/// Initializes input data branches
InitStatus ReadAndCreateDataBranches();
/// Registers histogram
- /// \param pHisto
+ /// \param pHisto
void RegisterHist(TH1* pHisto)
{
fHistList.push_back(pHisto);
@@ -151,17 +154,18 @@ private:
// *****************
// ** Flags and cuts **
- bool fIsMcPresent {false}; ///< Flag: true - all procedures involving MC will be processed, otherwice data will be processed only
+ bool fIsMcPresent {
+ false}; ///< Flag: true - all procedures involving MC will be processed, otherwise data will be processed only
- bool fIsQaPassed {true}; ///< Flag: true - QA is successfull, false - Qa is failed
+ bool fIsQaPassed {true}; ///< Flag: true - QA is successful, false - QA is failed
/// Maximum allowed distance between z-components of hit/MC-point position and the station center [cm]
- double fMaxDistanceZ {1.0};
+ double fMaxDistanceZ {1.0};
+
+ double fMinMomentum {0.01}; ///< Minimum value of particle momentum z component [GeV/c]
- double fMinMomentum {0.01}; ///< Minimum value of particle momentum z component [GeV/c]
-
double fMaxPullsRmsDiff {0.05}; ///< Maximum difference of the pulls RMS from unity
-
+
// ** Input arrays **
CbmTrackingDetectorInterfaceBase* fpDetectorInterface {nullptr}; ///< Pointer to current tracking detector I/F
@@ -171,10 +175,10 @@ private:
CbmMCEventList* fpMcEventList {nullptr};
CbmMCDataManager* fpMcManager {nullptr};
- TClonesArray* fpClusters {nullptr}; ///< Clusters
- TClonesArray* fpHits {nullptr}; ///< Hits
- TClonesArray* fpClusterMatches {nullptr}; ///< Matches for clusters (used only in STS)
- TClonesArray* fpHitMatches {nullptr}; ///< Matches for hits, provided with CbmMatchRecoToMC task
+ TClonesArray* fpClusters {nullptr}; ///< Clusters
+ TClonesArray* fpHits {nullptr}; ///< Hits
+ TClonesArray* fpClusterMatches {nullptr}; ///< Matches for clusters (used only in STS)
+ TClonesArray* fpHitMatches {nullptr}; ///< Matches for hits, provided with CbmMatchRecoToMC task
CbmMCDataArray* fpMcTracks {nullptr};
CbmMCDataArray* fpMcPoints {nullptr};
@@ -182,15 +186,15 @@ private:
// ** Output folders **
- TFolder fOutFolder {"TrackingInputQaSts", "TrackingInputQaSts"}; ///< Output folder for this QA task
- TFolder* fpOutFolderHists {nullptr}; ///< Subfolder for raw histograms
+ TFolder fOutFolder {"TrackingInputQaSts", "TrackingInputQaSts"}; ///< Output folder for this QA task
+ TFolder* fpOutFolderHists {nullptr}; ///< Subfolder for raw histograms
std::vector<TH1*> fHistList; ///< List of the pointers to all histograms contained in the class
// ** Histograms **
TParameter<long> fNevents {"nEvents", 0}; ///< Number of processed events
-
+
std::vector<CbmQaHist<TH1D>> fHistResidualX; ///< Residual distributions for X vs index of a station
std::vector<CbmQaHist<TH1D>> fHistResidualY; ///< Residual distributions for X vs index of a station
std::vector<CbmQaHist<TH1D>> fHistResidualT; ///< Residual distributions for T vs index of a station
@@ -203,10 +207,10 @@ private:
std::vector<CbmQaHist<TH1D>> fHistHitsPerPoint; ///< Distribution of hits per one MC point vs. index of a station
// TODO: There is a class TEfficiency in the ROOT, probably it's a good idea to use it here (S.Zharko)
- std::vector<CbmQaHist<TProfile2D>> fHistEfficiencyXY; ///<
- std::vector<CbmQaHist<TProfile>> fHistEfficiencyR; ///<
+ std::vector<CbmQaHist<TProfile2D>> fHistEfficiencyXY; ///<
+ std::vector<CbmQaHist<TProfile>> fHistEfficiencyR; ///<
- static constexpr std::tuple<int, double, double> fRangeResidualX {100, -0.01, 0.01};
+ static constexpr std::tuple<int, double, double> fRangeResidualX {100, -0.01, 0.01};
static constexpr std::tuple<int, double, double> fRangeResidualY {100, -0.03, 0.03};
static constexpr std::tuple<int, double, double> fRangeResidualT {100, -25., 25.};
@@ -220,9 +224,9 @@ private:
// ** Canvases **
- CbmQaCanvas fCanvResidualX {"c_sts_residualX", "STS X Residual distributions", 1500, 1000 };
- CbmQaCanvas fCanvResidualY {"c_sts_residualY", "STS Y Residual distributions", 1500, 1000 };
- CbmQaCanvas fCanvResidualT {"c_sts_residualT", "STS T Residual distributions", 1500, 1000 };
+ CbmQaCanvas fCanvResidualX {"c_sts_residualX", "STS X Residual distributions", 1500, 1000};
+ CbmQaCanvas fCanvResidualY {"c_sts_residualY", "STS Y Residual distributions", 1500, 1000};
+ CbmQaCanvas fCanvResidualT {"c_sts_residualT", "STS T Residual distributions", 1500, 1000};
CbmQaCanvas fCanvPullX {"c_sts_pullX", "STS X Pull distributions", 1500, 1000};
CbmQaCanvas fCanvPullY {"c_sts_pullY", "STS Y Pull distributions", 1500, 1000};
@@ -230,7 +234,7 @@ private:
CbmQaCanvas fCanvPointsPerHit {"c_sts_PointsPerHit", "STS MC Points per Hit", 1500, 1000};
CbmQaCanvas fCanvHitsPerPoint {"c_sts_HitsPerPoint", "STS Hits per MC Point", 1500, 1000};
-
+
CbmQaCanvas fCanvEfficiencyR {"c_sts_EfficiencyR", "STS Efficiency vs. R", 1500, 1000};
CbmQaCanvas fCanvEfficiencyXY {"c_sts_EfficiencyXY", "STS Efficiency vs. X and Y", 1500, 1000};
@@ -238,8 +242,6 @@ private:
// ************
// ** Output **
// ************
-
-
};
#endif