From ea40024607f87fe893e92f9ef05fe18360615945 Mon Sep 17 00:00:00 2001
From: "se.gorbunov" <se.gorbunov@gsi.de>
Date: Fri, 29 Dec 2023 17:46:57 +0000
Subject: [PATCH] Ca: extra checks for the input hits
---
algo/ca/core/pars/CaParameters.h | 6 +-
algo/ca/core/pars/CaStation.h | 14 +++
algo/ca/core/tracking/CaTrackFinder.cxx | 4 +-
core/detectors/trd/CbmTrdTrackingInterface.h | 32 +++--
reco/L1/CbmCaTimeSliceReader.cxx | 117 +++++++++++++++++++
reco/L1/CbmCaTimeSliceReader.h | 104 -----------------
reco/L1/CbmL1.cxx | 22 ++--
reco/L1/catools/CaToolsHitRecord.cxx | 13 +++
reco/L1/catools/CaToolsHitRecord.h | 4 +
reco/L1/qa/CbmCaInputQaSetup.cxx | 56 ++++++++-
reco/L1/qa/CbmCaInputQaSetup.h | 52 ---------
reco/L1/qa/CbmCaOutputQa.cxx | 17 ++-
12 files changed, 258 insertions(+), 183 deletions(-)
diff --git a/algo/ca/core/pars/CaParameters.h b/algo/ca/core/pars/CaParameters.h
index b1adb5beac..3ec7888091 100644
--- a/algo/ca/core/pars/CaParameters.h
+++ b/algo/ca/core/pars/CaParameters.h
@@ -112,6 +112,9 @@ namespace cbm::algo::ca
/// \param detectorID ID of the detector subsystem
[[gnu::always_inline]] int GetStationIndexGeometry(int localIndex, EDetectorID detectorID) const
{
+ if (localIndex >= GetNstationsGeometry(detectorID)) {
+ return -1;
+ }
return fvLocalToGeoIdMap[fvFirstGeoId[static_cast<int>(detectorID)] + localIndex];
}
@@ -120,7 +123,8 @@ namespace cbm::algo::ca
/// \param detectorID ID of the detector subsystem
[[gnu::always_inline]] int GetStationIndexActive(int localIndex, EDetectorID detectorID) const
{
- return fvGeoToActiveMap[GetStationIndexGeometry(localIndex, detectorID)];
+ int geoIndex = GetStationIndexGeometry(localIndex, detectorID);
+ return (geoIndex < 0) ? -1 : fvGeoToActiveMap[geoIndex];
}
/// \brief Provides access to L1CAIteration vector (const)
diff --git a/algo/ca/core/pars/CaStation.h b/algo/ca/core/pars/CaStation.h
index b7a072530c..a694f2cdd5 100644
--- a/algo/ca/core/pars/CaStation.h
+++ b/algo/ca/core/pars/CaStation.h
@@ -90,6 +90,13 @@ namespace cbm::algo::ca
return utils::simd::Cast<DataT, DataOut>(Xmax);
}
+ /// \brief Gets limit of the station size in x-axis direction
+ template<typename DataOut = DataT>
+ DataOut GetXmin() const
+ {
+ return utils::simd::Cast<DataT, DataOut>(-Xmax);
+ }
+
/// \brief Gets limit of the station size in y-axis direction
template<typename DataOut = DataT>
DataOut GetYmax() const
@@ -97,6 +104,13 @@ namespace cbm::algo::ca
return utils::simd::Cast<DataT, DataOut>(Ymax);
}
+ /// \brief Gets limit of the station size in y-axis direction
+ template<typename DataOut = DataT>
+ DataOut GetYmin() const
+ {
+ return utils::simd::Cast<DataT, DataOut>(-Ymax);
+ }
+
/// \brief String representation of class contents
/// \param verbosityLevel Verbosity level of the output
/// \param indentLevel Number of indent characters in the output
diff --git a/algo/ca/core/tracking/CaTrackFinder.cxx b/algo/ca/core/tracking/CaTrackFinder.cxx
index 8c231c174b..edc473faac 100644
--- a/algo/ca/core/tracking/CaTrackFinder.cxx
+++ b/algo/ca/core/tracking/CaTrackFinder.cxx
@@ -127,10 +127,10 @@ void TrackFinder::FindTracks()
info.fEventTimeMax = 1.e10;
}
- if (info.fEventTimeMin > 500 * 1.e6) { // cut hits with bogus start time > 500 ms
+ if (info.fEventTimeMin > 500.e6 || info.fEventTimeMax < -500.) { // cut hits with bogus start time > 500 ms
frAlgo.fvHitKeyFlags[h.FrontKey()] = 1;
frAlgo.fvHitKeyFlags[h.BackKey()] = 1;
- LOG(warning) << "CA: skip bogus hit " << h.ToString();
+ LOG(error) << "CA: skip bogus hit " << h.ToString();
continue;
}
diff --git a/core/detectors/trd/CbmTrdTrackingInterface.h b/core/detectors/trd/CbmTrdTrackingInterface.h
index cdb4d23938..88355618d4 100644
--- a/core/detectors/trd/CbmTrdTrackingInterface.h
+++ b/core/detectors/trd/CbmTrdTrackingInterface.h
@@ -17,9 +17,7 @@
#include "CbmTrdAddress.h"
#include "CbmTrdParModDigi.h"
#include "CbmTrdParSetDigi.h"
-
#include "FairTask.h"
-
#include "TMath.h"
#include <iostream>
@@ -32,7 +30,7 @@
/// @note For TRD one tracking station is a TRD module!
///
class CbmTrdTrackingInterface : public FairTask, public CbmTrackingDetectorInterfaceBase {
-public:
+ public:
/// @brief Default constructor
CbmTrdTrackingInterface();
@@ -84,22 +82,38 @@ public:
/// @brief Gets upper edge of a station along the X-axis
/// @param stationId Tracking station ID in the setup (NOTE: must be in range [0..GetNstations()-1])
/// @return Size of station along the X-axis [cm]
- double GetXmax(int stationId) const { return GetTrdModulePar(stationId)->GetSizeX(); }
+ double GetXmax(int stationId) const
+ {
+ const auto module = GetTrdModulePar(stationId);
+ return module->GetX() + module->GetSizeX();
+ }
/// @brief Gets lower edge coordinate of a station along the X-axis
/// @param stationId Tracking station ID in the setup (NOTE: must be in range [0..GetNstations()-1])
/// @return Size of station along the X-axis [cm]
- double GetXmin(int stationId) const { return -GetTrdModulePar(stationId)->GetSizeX(); }
+ double GetXmin(int stationId) const
+ {
+ const auto module = GetTrdModulePar(stationId);
+ return module->GetX() - module->GetSizeX();
+ }
/// @brief Gets upper edge of a station along the Y-axis
/// @param stationId Tracking station ID in the setup (NOTE: must be in range [0..GetNstations()-1])
/// @return Size of station along the Y-axis [cm]
- double GetYmax(int stationId) const { return GetTrdModulePar(stationId)->GetSizeY(); }
+ double GetYmax(int stationId) const
+ {
+ const auto module = GetTrdModulePar(stationId);
+ return module->GetY() + module->GetSizeY();
+ }
/// @brief Gets lower edge of a station along the Y-axis
/// @param stationId Tracking station ID in the setup (NOTE: must be in range [0..GetNstations()-1])
/// @return Size of station along the Y-axis [cm]
- double GetYmin(int stationId) const { return -GetTrdModulePar(stationId)->GetSizeY(); }
+ double GetYmin(int stationId) const
+ {
+ const auto module = GetTrdModulePar(stationId);
+ return module->GetY() - module->GetSizeY();
+ }
/// @brief Gets reference z of the station
/// @param stationId Tracking station ID in the setup (NOTE: must be in range [0..GetNstations()-1])
@@ -138,7 +152,7 @@ public:
/// @brief FairTask: sets parameter containers up
void SetParContainers();
-private:
+ private:
/// @brief Gets pointer to the TRD module
/// @param moduleId Id of the Trd module
/// @return Pointer to the particular CbmTrdParModDigi object
@@ -149,7 +163,7 @@ private:
static CbmTrdTrackingInterface* fpInstance; ///< Instance of the class
- CbmTrdParSetDigi* fTrdDigiPar {nullptr};
+ CbmTrdParSetDigi* fTrdDigiPar{nullptr};
//CbmTrdParModDigi* fTrdModuleInfo {nullptr};
ClassDef(CbmTrdTrackingInterface, 0);
diff --git a/reco/L1/CbmCaTimeSliceReader.cxx b/reco/L1/CbmCaTimeSliceReader.cxx
index c0ba1b46bd..d21ef8bd26 100644
--- a/reco/L1/CbmCaTimeSliceReader.cxx
+++ b/reco/L1/CbmCaTimeSliceReader.cxx
@@ -13,6 +13,7 @@
#include "CaInputData.h"
#include "CaParameters.h"
#include "CbmGlobalTrack.h"
+#include "CbmL1.h"
#include "CbmL1Util.h" // for CopyTrackParam2TC
#include "CbmMuchTrack.h"
#include "CbmStsTrack.h"
@@ -27,6 +28,7 @@
using cbm::ca::tools::HitRecord;
using namespace cbm::algo::ca::constants;
+using namespace cbm::algo;
using cbm::ca::TimeSliceReader;
@@ -350,6 +352,121 @@ void TimeSliceReader::SortQaHits()
fpvQaHits->SetName(name);
}
+
+// ---------------------------------------------------------------------------------------------------------------------
+//
+template<ca::EDetectorID DetID>
+int cbm::ca::TimeSliceReader::ReadHitsForDetector()
+{
+ using Hit_t = HitTypes_t::at<DetID>;
+
+ if (!fvbUseDet[DetID]) {
+ return 0;
+ } // Detector is entirelly not used
+
+ const auto* pDetInterface = fvpDetInterface[DetID];
+ int nHitsTot = fpEvent ? fpEvent->GetNofData(kCbmHitType[DetID]) : fvpBrHits[DetID]->GetEntriesFast();
+ int nHitsStored = 0; // number of hits used in tracking
+
+ fFirstHitKey = fNofHitKeys;
+
+ for (int iH = 0; iH < nHitsTot; ++iH) {
+ int iHext = fpEvent ? fpEvent->GetIndex(kCbmHitType[DetID], iH) : iH;
+ if (iHext < 0) {
+ LOG(warn) << "TimeSliceReader: hit index stored in the event is negative: " << iHext;
+ continue;
+ }
+ tools::HitRecord hitRecord;
+
+ auto* pHit = static_cast<Hit_t*>(fvpBrHits[DetID]->At(iHext));
+ int iStGeom = pDetInterface->GetTrackingStationIndex(pHit);
+
+ // Fill out detector specific data
+ if constexpr (ca::EDetectorID::kSts == DetID) {
+ hitRecord.fStripF = fFirstHitKey + pHit->GetFrontClusterId();
+ hitRecord.fStripB = fFirstHitKey + pHit->GetBackClusterId();
+ }
+ else if constexpr (ca::EDetectorID::kTof == DetID) {
+ // *** Additional cuts for TOF ***
+ // Skip Bmon hits
+ if (5 == CbmTofAddress::GetSmType(pHit->GetAddress())) {
+ continue;
+ }
+
+ // FIXME: Figure it out, if this cut is still needed (introduced a year ago for mCBM)
+ if (ECbmCaTrackingMode::kMCBM == fTrackingMode && pHit->GetZ() > 400) {
+ continue;
+ }
+ }
+
+ int iStActive = fpParameters->GetStationIndexActive(iStGeom, DetID);
+ if (iStActive < 0) {
+ continue;
+ } // Cut off inactive stations
+
+ // Fill out data common for all the detectors
+ hitRecord.fStaId = iStActive;
+ hitRecord.fX = pHit->GetX();
+ hitRecord.fY = pHit->GetY();
+ hitRecord.fZ = pHit->GetZ();
+ hitRecord.fDx2 = pHit->GetDx() * pHit->GetDx();
+ hitRecord.fDy2 = pHit->GetDy() * pHit->GetDy();
+ hitRecord.fDxy = pHit->GetDxy();
+ hitRecord.fT = pHit->GetTime();
+ hitRecord.fDt2 = pHit->GetTimeError() * pHit->GetTimeError();
+
+ // Apply hit error smearing according to misalignment
+ hitRecord.fDx2 += fpParameters->GetMisalignmentXsq(DetID);
+ hitRecord.fDy2 += fpParameters->GetMisalignmentYsq(DetID);
+ hitRecord.fDt2 += fpParameters->GetMisalignmentTsq(DetID);
+
+ std::tie(hitRecord.fRangeX, hitRecord.fRangeY, hitRecord.fRangeT) = pDetInterface->GetHitRanges(*pHit);
+
+ hitRecord.fDet = static_cast<int>(DetID);
+ hitRecord.fDataStream = (static_cast<int64_t>(hitRecord.fDet) << 60) | pHit->GetAddress();
+ hitRecord.fExtId = iHext;
+
+ // Check if hit values are reasonable
+ {
+ const ca::StationV& station = fpParameters->GetStation(iStActive);
+ if (pHit->GetX() < station.GetXmin<fscal>() * 1.2 || pHit->GetX() > station.GetXmax<fscal>() * 1.2
+ || pHit->GetY() < station.GetYmin<fscal>() * 1.2 || pHit->GetY() > station.GetYmax<fscal>() * 1.2
+ || fabs(pHit->GetZ() - station.GetZ<fscal>()) > 150) {
+ LOG(error) << "Ca:TimeSliceReader: " << CbmL1::GetDetectorName(DetID)
+ << " hit is outside of the station boundaries: " << hitRecord.ToString() << "; station X "
+ << station.GetXmin<fscal>() << ".." << station.GetXmax<fscal>() << " Y " << station.GetYmin<fscal>()
+ << ".." << station.GetYmax<fscal>() << " Z " << station.GetZ<fscal>();
+ continue;
+ }
+ }
+
+ // Update number of hit keys
+ if constexpr (ca::EDetectorID::kSts == DetID) {
+ if (fNofHitKeys <= hitRecord.fStripF) {
+ fNofHitKeys = hitRecord.fStripF + 1;
+ }
+ if (fNofHitKeys <= hitRecord.fStripB) {
+ fNofHitKeys = hitRecord.fStripB + 1;
+ }
+ }
+ else {
+ hitRecord.fStripF = fFirstHitKey + iHext;
+ hitRecord.fStripB = hitRecord.fStripF;
+ if (fNofHitKeys <= hitRecord.fStripF) {
+ fNofHitKeys = hitRecord.fStripF + 1;
+ }
+ }
+
+ // Save hit to data structures
+ if (hitRecord.Accept()) {
+ this->StoreHitRecord(hitRecord);
+ ++nHitsStored;
+ }
+ } // iH
+
+ return nHitsStored;
+}
+
// ---------------------------------------------------------------------------------------------------------------------
//
void TimeSliceReader::ReadHits()
diff --git a/reco/L1/CbmCaTimeSliceReader.h b/reco/L1/CbmCaTimeSliceReader.h
index e9b3b1c662..dc744a8b42 100644
--- a/reco/L1/CbmCaTimeSliceReader.h
+++ b/reco/L1/CbmCaTimeSliceReader.h
@@ -205,107 +205,3 @@ namespace cbm::ca
std::array<int, constants::size::MaxNdetectors + 1> fvHitFirstIndexDet = {{0}}; ///< First hit index in detector
};
} // namespace cbm::ca
-
-
-// **************************************************
-// ** Inline and template function implementations **
-// **************************************************
-
-// ---------------------------------------------------------------------------------------------------------------------
-//
-template<ca::EDetectorID DetID>
-int cbm::ca::TimeSliceReader::ReadHitsForDetector()
-{
- using Hit_t = HitTypes_t::at<DetID>;
-
- if (!fvbUseDet[DetID]) {
- return 0;
- } // Detector is entirelly not used
-
- const auto* pDetInterface = fvpDetInterface[DetID];
- int nHitsTot = fpEvent ? fpEvent->GetNofData(kCbmHitType[DetID]) : fvpBrHits[DetID]->GetEntriesFast();
- int nHitsStored = 0; // number of hits used in tracking
-
- fFirstHitKey = fNofHitKeys;
-
- for (int iH = 0; iH < nHitsTot; ++iH) {
- int iHext = fpEvent ? fpEvent->GetIndex(kCbmHitType[DetID], iH) : iH;
- if (iHext < 0) {
- continue;
- }
- tools::HitRecord hitRecord;
-
- auto* pHit = static_cast<Hit_t*>(fvpBrHits[DetID]->At(iHext));
- int iStGeom = pDetInterface->GetTrackingStationIndex(pHit);
-
- // Fill out detector specific data
- if constexpr (ca::EDetectorID::kSts == DetID) {
- hitRecord.fStripF = fFirstHitKey + pHit->GetFrontClusterId();
- hitRecord.fStripB = fFirstHitKey + pHit->GetBackClusterId();
- }
- else if constexpr (ca::EDetectorID::kTof == DetID) {
- // *** Additional cuts for TOF ***
- // Skip Bmon hits
- if (5 == CbmTofAddress::GetSmType(pHit->GetAddress())) {
- continue;
- }
-
- // FIXME: Figure it out, if this cut is still needed (introduced a year ago for mCBM)
- if (ECbmCaTrackingMode::kMCBM == fTrackingMode && pHit->GetZ() > 400) {
- continue;
- }
- }
-
- int iStActive = fpParameters->GetStationIndexActive(iStGeom, DetID);
- if (iStActive == -1) {
- continue;
- } // Cut off inactive stations
-
- // Fill out data common for all the detectors
- hitRecord.fStaId = iStActive;
- hitRecord.fX = pHit->GetX();
- hitRecord.fY = pHit->GetY();
- hitRecord.fZ = pHit->GetZ();
- hitRecord.fDx2 = pHit->GetDx() * pHit->GetDx();
- hitRecord.fDy2 = pHit->GetDy() * pHit->GetDy();
- hitRecord.fDxy = pHit->GetDxy();
- hitRecord.fT = pHit->GetTime();
- hitRecord.fDt2 = pHit->GetTimeError() * pHit->GetTimeError();
-
- // Apply hit error smearing according to misalignment
- hitRecord.fDx2 += fpParameters->GetMisalignmentXsq(DetID);
- hitRecord.fDy2 += fpParameters->GetMisalignmentYsq(DetID);
- hitRecord.fDt2 += fpParameters->GetMisalignmentTsq(DetID);
-
- std::tie(hitRecord.fRangeX, hitRecord.fRangeY, hitRecord.fRangeT) = pDetInterface->GetHitRanges(*pHit);
-
- hitRecord.fDet = static_cast<int>(DetID);
- hitRecord.fDataStream = (static_cast<int64_t>(hitRecord.fDet) << 60) | pHit->GetAddress();
- hitRecord.fExtId = iHext;
-
- // Update number of hit keys
- if constexpr (ca::EDetectorID::kSts == DetID) {
- if (fNofHitKeys <= hitRecord.fStripF) {
- fNofHitKeys = hitRecord.fStripF + 1;
- }
- if (fNofHitKeys <= hitRecord.fStripB) {
- fNofHitKeys = hitRecord.fStripB + 1;
- }
- }
- else {
- hitRecord.fStripF = fFirstHitKey + iHext;
- hitRecord.fStripB = hitRecord.fStripF;
- if (fNofHitKeys <= hitRecord.fStripF) {
- fNofHitKeys = hitRecord.fStripF + 1;
- }
- }
-
- // Save hit to data structures
- if (hitRecord.Accept()) {
- this->StoreHitRecord(hitRecord);
- ++nHitsStored;
- }
- } // iH
-
- return nHitsStored;
-}
diff --git a/reco/L1/CbmL1.cxx b/reco/L1/CbmL1.cxx
index b4afda8e6d..3c429ff45e 100644
--- a/reco/L1/CbmL1.cxx
+++ b/reco/L1/CbmL1.cxx
@@ -393,8 +393,8 @@ InitStatus CbmL1::Init()
stationInfo.SetZref(mvdInterface->GetZref(iSt));
stationInfo.SetZmin(mvdInterface->GetZmin(iSt));
stationInfo.SetZmax(mvdInterface->GetZmax(iSt));
- stationInfo.SetXmax(mvdInterface->GetXmax(iSt));
- stationInfo.SetYmax(mvdInterface->GetYmax(iSt));
+ stationInfo.SetXmax(std::max(fabs(mvdInterface->GetXmax(iSt)), fabs(mvdInterface->GetXmin(iSt))));
+ stationInfo.SetYmax(std::max(fabs(mvdInterface->GetYmax(iSt)), fabs(mvdInterface->GetYmin(iSt))));
stationInfo.SetTrackingStatus(true);
if (fvmDisabledStationIDs[ca::EDetectorID::kMvd].find(iSt)
!= fvmDisabledStationIDs[ca::EDetectorID::kMvd].cend()) {
@@ -416,8 +416,8 @@ InitStatus CbmL1::Init()
stationInfo.SetZref(stsInterface->GetZref(iSt));
stationInfo.SetZmin(stsInterface->GetZmin(iSt));
stationInfo.SetZmax(stsInterface->GetZmax(iSt));
- stationInfo.SetXmax(stsInterface->GetXmax(iSt));
- stationInfo.SetYmax(stsInterface->GetYmax(iSt));
+ stationInfo.SetXmax(std::max(fabs(stsInterface->GetXmax(iSt)), fabs(stsInterface->GetXmin(iSt))));
+ stationInfo.SetYmax(std::max(fabs(stsInterface->GetYmax(iSt)), fabs(stsInterface->GetYmin(iSt))));
stationInfo.SetTrackingStatus(true);
if (fvmDisabledStationIDs[ca::EDetectorID::kSts].find(iSt)
@@ -440,8 +440,8 @@ InitStatus CbmL1::Init()
stationInfo.SetZref(muchInterface->GetZref(iSt));
stationInfo.SetZmin(muchInterface->GetZmin(iSt));
stationInfo.SetZmax(muchInterface->GetZmax(iSt));
- stationInfo.SetXmax(muchInterface->GetXmax(iSt));
- stationInfo.SetYmax(muchInterface->GetYmax(iSt));
+ stationInfo.SetXmax(std::max(fabs(muchInterface->GetXmax(iSt)), fabs(muchInterface->GetXmin(iSt))));
+ stationInfo.SetYmax(std::max(fabs(muchInterface->GetYmax(iSt)), fabs(muchInterface->GetYmin(iSt))));
stationInfo.SetTrackingStatus(true);
if (fvmDisabledStationIDs[ca::EDetectorID::kMuch].find(iSt)
@@ -464,8 +464,8 @@ InitStatus CbmL1::Init()
stationInfo.SetZref(trdInterface->GetZref(iSt));
stationInfo.SetZmin(trdInterface->GetZmin(iSt));
stationInfo.SetZmax(trdInterface->GetZmax(iSt));
- stationInfo.SetXmax(trdInterface->GetXmax(iSt));
- stationInfo.SetYmax(trdInterface->GetYmax(iSt));
+ stationInfo.SetXmax(std::max(fabs(trdInterface->GetXmax(iSt)), fabs(trdInterface->GetXmin(iSt))));
+ stationInfo.SetYmax(std::max(fabs(trdInterface->GetYmax(iSt)), fabs(trdInterface->GetYmin(iSt))));
if (ca::Framework::TrackingMode::kGlobal == fTrackingMode) {
stationInfo.SetTimeInfo(false);
@@ -491,8 +491,8 @@ InitStatus CbmL1::Init()
stationInfo.SetZref(tofInterface->GetZref(iSt));
stationInfo.SetZmin(tofInterface->GetZmin(iSt));
stationInfo.SetZmax(tofInterface->GetZmax(iSt));
- stationInfo.SetXmax(tofInterface->GetXmax(iSt));
- stationInfo.SetYmax(tofInterface->GetYmax(iSt));
+ stationInfo.SetXmax(std::max(fabs(tofInterface->GetXmax(iSt)), fabs(tofInterface->GetXmin(iSt))));
+ stationInfo.SetYmax(std::max(fabs(tofInterface->GetYmax(iSt)), fabs(tofInterface->GetYmin(iSt))));
stationInfo.SetTrackingStatus(true);
if (fvmDisabledStationIDs[ca::EDetectorID::kTof].find(iSt)
@@ -818,7 +818,7 @@ void CbmL1::GenerateMaterialMaps()
for (unsigned int ist = 0; ist < vpActiveStations.size(); ++ist) {
auto* pStation = vpActiveStations[ist];
double z1 = pStation->GetZmax();
- double z2 = z1;
+ double z2 = z1;
if (ist < vpActiveStations.size() - 1) {
// split materials between the stations at the middle
auto* pStationNext = vpActiveStations[ist + 1];
diff --git a/reco/L1/catools/CaToolsHitRecord.cxx b/reco/L1/catools/CaToolsHitRecord.cxx
index 37be522ea2..e3781d4c92 100644
--- a/reco/L1/catools/CaToolsHitRecord.cxx
+++ b/reco/L1/catools/CaToolsHitRecord.cxx
@@ -67,3 +67,16 @@ bool HitRecord::Accept() const
}
return true;
}
+
+// ---------------------------------------------------------------------------------------------------------------------
+//
+std::string HitRecord::ToString() const
+{
+ std::stringstream msg;
+ msg << "HitRecord: det = " << fDet << ", addr = " << static_cast<int32_t>(fDataStream) << ", extId = " << fExtId
+ << ", staId = " << fStaId << ", pointId = " << fPointId << ", stripF = " << fStripF << ", stripB = " << fStripB
+ << ", x = " << fX << ", y = " << fY << ", z = " << fZ << ", t = " << fT << ", dx2 = " << fDx2
+ << ", dy2 = " << fDy2 << ", dt2 = " << fDt2 << ", dxy = " << fDxy << ", rangeX = " << fRangeX
+ << ", rangeY = " << fRangeY << ", rangeT = " << fRangeT;
+ return msg.str();
+}
diff --git a/reco/L1/catools/CaToolsHitRecord.h b/reco/L1/catools/CaToolsHitRecord.h
index d51d959404..f40c14ec03 100644
--- a/reco/L1/catools/CaToolsHitRecord.h
+++ b/reco/L1/catools/CaToolsHitRecord.h
@@ -11,6 +11,7 @@
#define CbmCaToolsHitRecord_h 1
#include <cstdint>
+#include <string>
namespace cbm::ca::tools
{
@@ -44,6 +45,9 @@ namespace cbm::ca::tools
/// @return false Hit is discarded
static constexpr bool kVerboseAccept = true;
bool Accept() const;
+
+ /// @brief Converts hit record to a string
+ std::string ToString() const;
};
} // namespace cbm::ca::tools
diff --git a/reco/L1/qa/CbmCaInputQaSetup.cxx b/reco/L1/qa/CbmCaInputQaSetup.cxx
index 4ff128997c..f264ff6f08 100644
--- a/reco/L1/qa/CbmCaInputQaSetup.cxx
+++ b/reco/L1/qa/CbmCaInputQaSetup.cxx
@@ -47,6 +47,60 @@ void InputQaSetup::CheckInit() const
}
}
+// -----------------------------------------------------------------------------------------------------------------
+//
+template<cbm::algo::ca::EDetectorID DetID>
+void InputQaSetup::FillHistogramsDet()
+{
+ using Hit_t = HitTypes_t::at<DetID>;
+ int nHits = fvpBrHits[DetID]->GetEntriesFast();
+
+ for (int iH = 0; iH < nHits; ++iH) {
+ const auto* pHit = dynamic_cast<const Hit_t*>(fvpBrHits[DetID]->At(iH));
+ if (!pHit) {
+ LOG(warn) << fName << ": hit with iH = " << iH << " for detector " << kDetName[DetID] << " is not found";
+ }
+ auto address = pHit->GetAddress();
+
+ // skip Bmon hits
+ if constexpr (ca::EDetectorID::kTof == DetID) {
+ if (5 == CbmTofAddress::GetSmType(address)) {
+ continue;
+ }
+ }
+
+ int iStLoc = fvpDetInterface[DetID]->GetTrackingStationIndex(address);
+ int iStGeo = fpParameters->GetStationIndexGeometry(iStLoc, DetID);
+ auto xHit = pHit->GetX();
+ auto yHit = pHit->GetY();
+ auto zHit = pHit->GetZ();
+ if (fvXmin[DetID][iStLoc] > xHit) {
+ fvXmin[DetID][iStLoc] = xHit;
+ }
+ if (fvXmax[DetID][iStLoc] < xHit) {
+ fvXmax[DetID][iStLoc] = xHit;
+ }
+ if (fvYmin[DetID][iStLoc] > yHit) {
+ fvYmin[DetID][iStLoc] = yHit;
+ }
+ if (fvYmax[DetID][iStLoc] < yHit) {
+ fvYmax[DetID][iStLoc] = yHit;
+ }
+ if (fvZmin[DetID][iStLoc] > zHit) {
+ fvZmin[DetID][iStLoc] = zHit;
+ }
+ if (fvZmax[DetID][iStLoc] < zHit) {
+ fvZmax[DetID][iStLoc] = zHit;
+ }
+ if (iStGeo >= 0) {
+ fph_hit_xz[iStGeo]->Fill(zHit, xHit);
+ fph_hit_yz[iStGeo]->Fill(zHit, yHit);
+ }
+ fph_hit_xz.back()->Fill(zHit, xHit);
+ fph_hit_yz.back()->Fill(zHit, yHit);
+ } // iH
+}
+
// ---------------------------------------------------------------------------------------------------------------------
//
void InputQaSetup::FillHistograms()
@@ -106,7 +160,7 @@ InitStatus InputQaSetup::InitCanvases()
int nSt = fvpDetInterface[iDet]->GetNtrackingStations();
for (int iSt = 0; iSt < nSt; ++iSt) {
int iStActive = fpParameters->GetStationIndexActive(iSt, static_cast<ca::EDetectorID>(iDet));
- Color_t boxColor = (iStActive == -1) ? kGray + 1 : kOrange + 2;
+ Color_t boxColor = (iStActive < 0) ? kGray + 1 : kOrange + 2;
const char* detName = Form("%s-%d", fvpDetInterface[iDet]->GetDetectorName().c_str(), iSt);
{
double zMin = LoBinEdge(fph_hit_xz.back()->GetXaxis(), fvZmin[iDet][iSt]);
diff --git a/reco/L1/qa/CbmCaInputQaSetup.h b/reco/L1/qa/CbmCaInputQaSetup.h
index c378e35641..e0b06e44df 100644
--- a/reco/L1/qa/CbmCaInputQaSetup.h
+++ b/reco/L1/qa/CbmCaInputQaSetup.h
@@ -110,56 +110,4 @@ namespace cbm::ca
std::vector<TH2F*> fph_point_yz = {}; ///< MC point occupancy in y-z plane vs. station
};
- // -----------------------------------------------------------------------------------------------------------------
- //
- template<ca::EDetectorID DetID>
- void InputQaSetup::FillHistogramsDet()
- {
- using Hit_t = HitTypes_t::at<DetID>;
- int nHits = fvpBrHits[DetID]->GetEntriesFast();
-
- for (int iH = 0; iH < nHits; ++iH) {
- const auto* pHit = dynamic_cast<const Hit_t*>(fvpBrHits[DetID]->At(iH));
- if (!pHit) {
- LOG(warn) << fName << ": hit with iH = " << iH << " for detector " << kDetName[DetID] << " is not found";
- }
- auto address = pHit->GetAddress();
-
- // skip Bmon hits
- if constexpr (ca::EDetectorID::kTof == DetID) {
- if (5 == CbmTofAddress::GetSmType(address)) {
- continue;
- }
- }
-
- int iStLoc = fvpDetInterface[DetID]->GetTrackingStationIndex(address);
- int iStGeo = fpParameters->GetStationIndexGeometry(iStLoc, DetID);
- auto xHit = pHit->GetX();
- auto yHit = pHit->GetY();
- auto zHit = pHit->GetZ();
- if (fvXmin[DetID][iStLoc] > xHit) {
- fvXmin[DetID][iStLoc] = xHit;
- }
- if (fvXmax[DetID][iStLoc] < xHit) {
- fvXmax[DetID][iStLoc] = xHit;
- }
- if (fvYmin[DetID][iStLoc] > yHit) {
- fvYmin[DetID][iStLoc] = yHit;
- }
- if (fvYmax[DetID][iStLoc] < yHit) {
- fvYmax[DetID][iStLoc] = yHit;
- }
- if (fvZmin[DetID][iStLoc] > zHit) {
- fvZmin[DetID][iStLoc] = zHit;
- }
- if (fvZmax[DetID][iStLoc] < zHit) {
- fvZmax[DetID][iStLoc] = zHit;
- }
-
- fph_hit_xz[iStGeo]->Fill(zHit, xHit);
- fph_hit_yz[iStGeo]->Fill(zHit, yHit);
- fph_hit_xz.back()->Fill(zHit, xHit);
- fph_hit_yz.back()->Fill(zHit, yHit);
- } // iH
- }
} // namespace cbm::ca
diff --git a/reco/L1/qa/CbmCaOutputQa.cxx b/reco/L1/qa/CbmCaOutputQa.cxx
index 88c0551657..363b38027d 100644
--- a/reco/L1/qa/CbmCaOutputQa.cxx
+++ b/reco/L1/qa/CbmCaOutputQa.cxx
@@ -432,7 +432,21 @@ InitStatus OutputQa::InitCanvases()
//
InitStatus OutputQa::InitDataBranches()
{
+
+ LOG_IF(fatal, !fpParameters.get())
+ << fName << ": CA parameters object is not defined. Please, provide initializer or read parameters from binary "
+ << "via OutputQa::ReadParameters(filename) from the qa macro";
+
+ // Turn off detectors that are not used in the reconstruction setup
+
+ fbUseMvd = fbUseMvd && (fpParameters->GetNstationsActive(ca::EDetectorID::kMvd) > 0);
+ fbUseSts = fbUseSts && (fpParameters->GetNstationsActive(ca::EDetectorID::kSts) > 0);
+ fbUseMuch = fbUseMuch && (fpParameters->GetNstationsActive(ca::EDetectorID::kMuch) > 0);
+ fbUseTrd = fbUseTrd && (fpParameters->GetNstationsActive(ca::EDetectorID::kTrd) > 0);
+ fbUseTof = fbUseTof && (fpParameters->GetNstationsActive(ca::EDetectorID::kTof) > 0);
+
// Turn off detectors, which hits are not presented in input tree
+
auto* fairManager = FairRootManager::Instance();
fbUseMvd = fbUseMvd && fairManager->GetObject("MvdHit");
fbUseSts = fbUseSts && fairManager->GetObject("StsHit");
@@ -449,9 +463,6 @@ InitStatus OutputQa::InitDataBranches()
LOG(info) << "\tTOF: " << (fbUseTof ? "ON" : "OFF");
LOG(info) << fName << ": Initializing data branches";
- LOG_IF(fatal, !fpParameters.get())
- << fName << ": CA parameters object is not defined. Please, provide initializer or read parameters from binary "
- << "via OutputQa::ReadParameters(filename) from the qa macro";
// Initialize IO data manager
if (!fpDataManager.get()) {
--
GitLab