diff --git a/core/base/CbmHistManager.h b/core/base/CbmHistManager.h
index 5294f3325adb7ef18df457729a592ecf681043cf..390a5a86b7ed56067981d29591f109e0e14f37a5 100644
--- a/core/base/CbmHistManager.h
+++ b/core/base/CbmHistManager.h
@@ -88,6 +88,21 @@ public:
Add(name, h);
}
+ /**
+ * \brief Helper function for creation of 1-dimensional histograms and profiles.
+ * Template argument is a real object type that has to be created, for example,
+ * Create1<TH1F>("name", "title", binsX);
+ * \param[in] name Object name.
+ * \param[in] title Object title.
+ * \param[in] binsX array of low-edges for each bin in X
+ */
+ template<class T>
+ void Create1(const std::string& name, const std::string& title, const std::vector<Double_t>& binsX)
+ {
+ T* h = new T(name.c_str(), title.c_str(), binsX.size() - 1, binsX.data());
+ Add(name, h);
+ }
+
/**
* \brief Helper function for creation of 2-dimensional histograms and profiles.
* Template argument is a real object type that has to be created, for example,
@@ -109,6 +124,23 @@ public:
Add(name, h);
}
+ /**
+ * \brief Helper function for creation of 2-dimensional histograms and profiles.
+ * Template argument is a real object type that has to be created, for example,
+ * Create2<TH2F>("name", "title", binsX, binsY);
+ * \param[in] name Object name.
+ * \param[in] title Object title.
+ * \param[in] binsX array of low-edges for each bin in X
+ * \param[in] binsY array of low-edges for each bin in Y
+ */
+ template<class T>
+ void Create2(const std::string& name, const std::string& title, const std::vector<Double_t>& binsX,
+ const std::vector<Double_t>& binsY)
+ {
+ T* h = new T(name.c_str(), title.c_str(), binsX.size() - 1, binsX.data(), binsY.size() - 1, binsY.data());
+ Add(name, h);
+ }
+
/**
* \brief Helper function for creation of 3-dimensional histograms and profiles.
* Template argument is a real object type that has to be created, for example,
diff --git a/core/base/utils/CbmUtils.cxx b/core/base/utils/CbmUtils.cxx
index c4448298030302f83ba26ee7ae2645e217939500..7fea819a9f491b82250728caa4d555c2c812e050 100644
--- a/core/base/utils/CbmUtils.cxx
+++ b/core/base/utils/CbmUtils.cxx
@@ -25,6 +25,14 @@ using std::vector;
namespace Cbm
{
+ CbmMCDataArray* InitOrFatalMc(const std::string& objName, const std::string& description)
+ {
+ CbmMCDataManager* mcManager = GetOrFatal<CbmMCDataManager>("MCDataManager", description);
+ CbmMCDataArray* array = mcManager->InitBranch(objName.c_str());
+ if (array == nullptr) { LOG(fatal) << description << " No MCTrack!"; }
+ return array;
+ }
+
void SaveCanvasAsImage(TCanvas* c, const string& dir, const string& option)
{
if (dir == "") return;
diff --git a/core/base/utils/CbmUtils.h b/core/base/utils/CbmUtils.h
index cb486b2ecc3948f389b438ddfc77e60fee575b40..63826c516ef81ab555edddfef4a025c1705bced0 100644
--- a/core/base/utils/CbmUtils.h
+++ b/core/base/utils/CbmUtils.h
@@ -5,6 +5,11 @@
#ifndef CBMUTILS_H_
#define CBMUTILS_H_
+#include "../CbmMCDataManager.h"
+
+#include "FairLogger.h"
+#include "FairRootManager.h"
+
#include <sstream> // for string, stringstream
#include <vector> // for vector
@@ -46,6 +51,28 @@ namespace Cbm
return (x > ZERO) ? 1 : ((x < ZERO) ? -1 : 0);
}
+ /*
+ \brief Tries to get object from FairRootManager. If object is not found create Fatal error.
+ \param[in] objName Name of the object.
+ \param[in] description Optional description for LOG, usually class name + method name.
+ */
+ template<typename T>
+ T* GetOrFatal(const std::string& objName, const std::string& description = "")
+ {
+ FairRootManager* ioman = FairRootManager::Instance();
+ if (ioman == nullptr) { LOG(fatal) << description << " No FairRootManager!"; }
+ T* obj = static_cast<T*>(ioman->GetObject(objName.c_str()));
+ if (obj == nullptr) { LOG(fatal) << description << " No " << objName << " object!"; }
+ return obj;
+ }
+
+ /*
+ \brief Tries to get CbmMCDataArray from CbmMCDataManager. If object is not found create Fatal error.
+ \param[in] objName Name of the object.
+ \param[in] description Optional description for LOG, usually class name + method name.
+ */
+ CbmMCDataArray* InitOrFatalMc(const std::string& objName, const std::string& description = "");
+
void SaveCanvasAsImage(TCanvas* c, const std::string& dir, const std::string& option = "eps;png;gif");
void SaveCanvasAsImageImpl(const std::string& imageType, TCanvas* c, const std::string& dir,
const std::string& option);
diff --git a/core/detectors/rich/CMakeLists.txt b/core/detectors/rich/CMakeLists.txt
index a1c1339fe075dde7e43b93d2f4fbb19119cd5ed1..bebe61a875dda0c41d9f565fdfcce83b4321eaca 100644
--- a/core/detectors/rich/CMakeLists.txt
+++ b/core/detectors/rich/CMakeLists.txt
@@ -4,7 +4,9 @@ set(INCLUDE_DIRECTORIES
${CBMDATA_DIR}
${CBMDATA_DIR}/rich
+ ${CBMDATA_DIR}/tof
${CBMDATA_DIR}/global
+ ${CBMDATA_DIR}/base
${CBMBASE_DIR}
${CBMBASE_DIR}/draw
@@ -39,6 +41,7 @@ set(NO_DICT_SRCS
CbmRichGeoManager.cxx
CbmRichDigiMapManager.cxx
CbmRichMCbmDigiMapManager.cxx
+ utils/CbmRichUtil.cxx
)
IF (SSE_FOUND)
diff --git a/core/detectors/rich/CbmRichDraw.h b/core/detectors/rich/CbmRichDraw.h
index 897be27f742cc6259a6ef21851a35084c5ec8b21..5303ceab3c7570ca2fd4eb4586e594046ef4a980 100644
--- a/core/detectors/rich/CbmRichDraw.h
+++ b/core/detectors/rich/CbmRichDraw.h
@@ -5,9 +5,8 @@
#ifndef RICH_CbmRichDraw
#define RICH_CbmRichDraw
-#include "CbmRichDetectorData.h" // for CbmRichPmtData, CbmRichPixelData
-#include "CbmRichDigiMapManager.h" // for CbmRichDigiMapManager
-#include "CbmRichGeoManager.h" // for CbmRichGeoManager
+
+#include "CbmRichUtil.h"
#include <RtypesCore.h> // for ROOT data types
#include <TCanvas.h> // for TCanvas
@@ -29,7 +28,7 @@ public:
TH2D* hUp = (TH2D*) h->Clone();
DrawH2(hUp);
if (usePmtBins) {
- std::vector<Double_t> yPmtBins = CbmRichDraw::GetPmtHistYbins();
+ std::vector<Double_t> yPmtBins = CbmRichUtil::GetPmtHistYbins();
hUp->GetYaxis()->SetRange(yPmtBins.size() / 2 + 1, yPmtBins.size());
}
else {
@@ -42,7 +41,7 @@ public:
c->cd(2);
TH2D* hDown = (TH2D*) h->Clone();
if (usePmtBins) {
- std::vector<Double_t> yPmtBins = CbmRichDraw::GetPmtHistYbins();
+ std::vector<Double_t> yPmtBins = CbmRichUtil::GetPmtHistYbins();
hDown->GetYaxis()->SetRange(0, yPmtBins.size() / 2 - 1);
}
else {
@@ -72,48 +71,6 @@ public:
gPad->SetLeftMargin(0.1);
gPad->SetRightMargin(0.15);
}
-
- static std::vector<Double_t> GetPmtHistXbins() { return CbmRichDraw::GetPmtHistBins(true); }
-
- static std::vector<Double_t> GetPmtHistYbins() { return CbmRichDraw::GetPmtHistBins(false); }
-
-private:
- static std::vector<Double_t> GetPmtHistBins(Bool_t isX)
- {
- std::vector<Double_t> initVec;
- std::vector<Int_t> pmts = CbmRichDigiMapManager::GetInstance().GetPmtIds();
- for (Int_t pmtId : pmts) {
- CbmRichPmtData* pmtData = CbmRichDigiMapManager::GetInstance().GetPmtDataById(pmtId);
- TVector3 inPos(pmtData->fX, pmtData->fY, pmtData->fZ);
- TVector3 outPos;
- CbmRichGeoManager::GetInstance().RotatePoint(&inPos, &outPos);
- if (isX) {
- initVec.push_back(outPos.X() - 0.5 * pmtData->fWidth);
- initVec.push_back(outPos.X() + 0.5 * pmtData->fWidth);
- }
- else {
- initVec.push_back(outPos.Y() - 0.5 * pmtData->fHeight);
- initVec.push_back(outPos.Y() + 0.5 * pmtData->fHeight);
- }
- }
- sort(initVec.begin(), initVec.end());
-
- std::vector<Double_t> uniVec;
- for (unsigned int i = 0; i < initVec.size(); i++) {
- if (i == 0) uniVec.push_back(initVec[i]);
- if (initVec[i] - uniVec[uniVec.size() - 1] > 0.000001) uniVec.push_back(initVec[i]);
- }
-
- // cout << "uniVec.size():" << uniVec.size() << endl;
- // for (int i = 0; i < uniVec.size(); i++) {
- // cout << std::setprecision(9);
- // cout << fixed;
- // cout << uniVec[i] << " " ;
- // }
- // cout << endl;
-
- return uniVec;
- }
};
#endif
diff --git a/core/detectors/rich/utils/CbmRichUtil.cxx b/core/detectors/rich/utils/CbmRichUtil.cxx
new file mode 100644
index 0000000000000000000000000000000000000000..af7622e14b94aa41a7af3f739794637344cb73d9
--- /dev/null
+++ b/core/detectors/rich/utils/CbmRichUtil.cxx
@@ -0,0 +1,136 @@
+/* Copyright (C) 2021 Justus-Liebig-Universitaet Giessen, Giessen
+ SPDX-License-Identifier: GPL-3.0-only
+ Authors: Semen Lebedev [committer]*/
+
+#include "CbmRichUtil.h"
+
+#include "CbmDigiManager.h"
+#include "CbmGlobalTrack.h"
+#include "CbmMCDataArray.h"
+#include "CbmMatchRecoToMC.h"
+#include "CbmRichDetectorData.h"
+#include "CbmRichDigi.h"
+#include "CbmRichDigiMapManager.h"
+#include "CbmRichGeoManager.h"
+#include "CbmRichHit.h"
+#include "CbmRichRing.h"
+
+#include "TClonesArray.h"
+
+using namespace std;
+
+map<pair<int, int>, int> CbmRichUtil::CreateNofHitsInRingMap(TClonesArray* richHits, CbmMCDataArray* richPoints,
+ CbmMCDataArray* mcTracks, CbmDigiManager* digiMan)
+{
+ map<pair<int, int>, int> nofHitsInRingMap;
+ int nofRichHits = richHits->GetEntriesFast();
+ for (int iHit = 0; iHit < nofRichHits; iHit++) {
+ const CbmRichHit* hit = static_cast<CbmRichHit*>(richHits->At(iHit));
+ if (hit == nullptr || hit->GetRefId() < 0) continue;
+ const CbmRichDigi* digi = digiMan->Get<CbmRichDigi>(hit->GetRefId());
+ if (digi == nullptr) continue;
+ const CbmMatch* digiMatch = digiMan->GetMatch(ECbmModuleId::kRich, hit->GetRefId());
+ if (digiMatch == nullptr) continue;
+ int eventId = digiMatch->GetMatchedLink().GetEntry();
+
+ vector<pair<int, int>> motherIds =
+ CbmMatchRecoToMC::GetMcTrackMotherIdsForRichHit(digiMan, hit, richPoints, mcTracks, eventId);
+ for (size_t i = 0; i < motherIds.size(); i++) {
+ nofHitsInRingMap[motherIds[i]]++;
+ }
+ }
+ return nofHitsInRingMap;
+}
+
+double CbmRichUtil::GetRingTrackDistance(Int_t globalTrackId) { return GetRingTrackDistanceImpl(globalTrackId)[0]; }
+
+double CbmRichUtil::GetRingTrackDistanceX(Int_t globalTrackId) { return GetRingTrackDistanceImpl(globalTrackId)[1]; }
+
+double CbmRichUtil::GetRingTrackDistanceY(Int_t globalTrackId) { return GetRingTrackDistanceImpl(globalTrackId)[2]; }
+
+vector<Double_t> CbmRichUtil::GetRingTrackDistanceImpl(Int_t globalTrackId)
+{
+ vector<Double_t> errorVec = {999., 999., 999.};
+ FairRootManager* ioman = FairRootManager::Instance();
+ if (ioman == NULL) return errorVec;
+ // Do we really need static here, depends on ioman->GetObject() method
+ static TClonesArray* globalTracks = NULL;
+ static TClonesArray* richRings = NULL;
+ static TClonesArray* richProjections = NULL;
+
+ if (globalTracks == NULL || richRings == NULL || richProjections == NULL) {
+ //cout << "globalTracks == NULL || richRings == NULL || richProjections == NULL" << endl;
+ globalTracks = (TClonesArray*) ioman->GetObject("GlobalTrack");
+ richRings = (TClonesArray*) ioman->GetObject("RichRing");
+ richProjections = (TClonesArray*) ioman->GetObject("RichProjection");
+ }
+ else {
+ //cout << "globalTracks, richRings, richProjections NOT NULL" << endl;
+ }
+
+ if (globalTracks == NULL || richRings == NULL || richProjections == NULL) {
+ LOG(error) << "CbmRichUtil::GetRingTrackDistance globalTracks, "
+ "richRings, richProjections NOT INITIALIZED"
+ << endl;
+ return errorVec;
+ }
+
+ const CbmGlobalTrack* globalTrack = static_cast<const CbmGlobalTrack*>(globalTracks->At(globalTrackId));
+ if (globalTrack == NULL) return errorVec;
+
+ Int_t stsId = globalTrack->GetStsTrackIndex();
+ if (stsId < 0) return errorVec;
+
+ const FairTrackParam* pTrack = static_cast<const FairTrackParam*>(richProjections->At(stsId));
+ if (pTrack == NULL) return errorVec;
+
+ if (pTrack->GetX() == 0 && pTrack->GetY() == 0) return errorVec;
+
+ Int_t richId = globalTrack->GetRichRingIndex();
+ if (richId < 0) return errorVec;
+
+ const CbmRichRing* richRing = static_cast<const CbmRichRing*>(richRings->At(richId));
+ if (richRing == NULL) return errorVec;
+
+ // Double_t xRing = richRing->GetCenterX();
+ // Double_t yRing = richRing->GetCenterY();
+ Double_t dx = richRing->GetCenterX() - pTrack->GetX();
+ Double_t dy = richRing->GetCenterY() - pTrack->GetY();
+
+ Double_t dist = TMath::Sqrt(dx * dx + dy * dy);
+
+ vector<Double_t> v = {dist, dx, dy};
+ return v;
+}
+
+vector<double> CbmRichUtil::GetPmtHistXbins() { return GetPmtHistBins(true); }
+
+vector<double> CbmRichUtil::GetPmtHistYbins() { return GetPmtHistBins(false); }
+
+vector<double> CbmRichUtil::GetPmtHistBins(bool isX)
+{
+ vector<double> initVec;
+ vector<int> pmts = CbmRichDigiMapManager::GetInstance().GetPmtIds();
+ for (int pmtId : pmts) {
+ CbmRichPmtData* pmtData = CbmRichDigiMapManager::GetInstance().GetPmtDataById(pmtId);
+ TVector3 inPos(pmtData->fX, pmtData->fY, pmtData->fZ);
+ TVector3 outPos;
+ CbmRichGeoManager::GetInstance().RotatePoint(&inPos, &outPos);
+ if (isX) {
+ initVec.push_back(outPos.X() - 0.5 * pmtData->fWidth);
+ initVec.push_back(outPos.X() + 0.5 * pmtData->fWidth);
+ }
+ else {
+ initVec.push_back(outPos.Y() - 0.5 * pmtData->fHeight);
+ initVec.push_back(outPos.Y() + 0.5 * pmtData->fHeight);
+ }
+ }
+ sort(initVec.begin(), initVec.end());
+
+ vector<double> uniVec;
+ for (size_t i = 0; i < initVec.size(); i++) {
+ if (i == 0) uniVec.push_back(initVec[i]);
+ if (initVec[i] - uniVec[uniVec.size() - 1] > 0.000001) uniVec.push_back(initVec[i]);
+ }
+ return uniVec;
+}
diff --git a/core/detectors/rich/utils/CbmRichUtil.h b/core/detectors/rich/utils/CbmRichUtil.h
index 47f564e7e5c58c6e4022e354c8ec0f318c48eadf..ab6de905f75e64e269d14a9a3bdcd8aad4963ea2 100644
--- a/core/detectors/rich/utils/CbmRichUtil.h
+++ b/core/detectors/rich/utils/CbmRichUtil.h
@@ -1,108 +1,42 @@
-/* Copyright (C) 2017-2019 GSI/JINR-LIT, Darmstadt/Dubna
+/* Copyright (C) 2017-2021 GSI/JINR-LIT, Darmstadt/Dubna
SPDX-License-Identifier: GPL-3.0-only
Authors: Andrey Lebedev [committer], Semen Lebedev */
#ifndef RICH_CbmRichUtil
#define RICH_CbmRichUtil
-#include "CbmGlobalTrack.h"
+#include <map>
+#include <vector>
-#include "FairRootManager.h"
-#include "FairTrackParam.h"
-#include <Logger.h>
-
-#include "TCanvas.h"
-#include "TClonesArray.h"
-#include "TH2.h"
-#include "TObject.h"
-
-#include <string>
-
-using namespace std;
+class CbmDigiManager;
+class TClonesArray;
+class CbmMCDataArray;
+class TH2D;
class CbmRichUtil {
public:
- static Double_t GetRingTrackDistance(Int_t globalTrackId)
- {
- vector<Double_t> v = GetRingTrackDistanceImpl(globalTrackId);
- return v[0];
- }
+ static double GetRingTrackDistance(int globalTrackId);
+ static double GetRingTrackDistanceX(int globalTrackId);
+ static double GetRingTrackDistanceY(int globalTrackId);
- static Double_t GetRingTrackDistanceX(Int_t globalTrackId)
- {
- vector<Double_t> v = GetRingTrackDistanceImpl(globalTrackId);
- return v[1];
- }
+ // Create PMT XY histograms
+ static std::vector<double> GetPmtHistXbins();
+ static std::vector<double> GetPmtHistYbins();
+ static std::vector<double> GetPmtHistBins(bool isX);
- static Double_t GetRingTrackDistanceY(Int_t globalTrackId)
- {
- vector<Double_t> v = GetRingTrackDistanceImpl(globalTrackId);
- return v[2];
- }
+ static std::map<std::pair<int, int>, int> CreateNofHitsInRingMap(TClonesArray* richHits, CbmMCDataArray* richPoints,
+ CbmMCDataArray* mcTracks, CbmDigiManager* digiMan);
- static uint16_t GetDirichId(Int_t Address) { return ((Address >> 16) & 0xFFFF); }
+ static uint16_t GetDirichId(int Address) { return ((Address >> 16) & 0xFFFF); }
- static uint16_t GetDirichChannel(Int_t Address) { return (Address & 0xFFFF); }
+ static uint16_t GetDirichChannel(int Address) { return (Address & 0xFFFF); }
private:
/**
* \brief Return a vector with total distance and x, y components. [0] - total distance, [1] - x component, [2] - y component
*/
- static vector<Double_t> GetRingTrackDistanceImpl(Int_t globalTrackId)
- {
- vector<Double_t> errorVec = {999., 999., 999.};
- FairRootManager* ioman = FairRootManager::Instance();
- if (ioman == NULL) return errorVec;
- // Do we really need static here, depends on ioman->GetObject() method
- static TClonesArray* globalTracks = NULL;
- static TClonesArray* richRings = NULL;
- static TClonesArray* richProjections = NULL;
-
- if (globalTracks == NULL || richRings == NULL || richProjections == NULL) {
- //cout << "globalTracks == NULL || richRings == NULL || richProjections == NULL" << endl;
- globalTracks = (TClonesArray*) ioman->GetObject("GlobalTrack");
- richRings = (TClonesArray*) ioman->GetObject("RichRing");
- richProjections = (TClonesArray*) ioman->GetObject("RichProjection");
- }
- else {
- //cout << "globalTracks, richRings, richProjections NOT NULL" << endl;
- }
-
- if (globalTracks == NULL || richRings == NULL || richProjections == NULL) {
- LOG(error) << "CbmRichUtil::GetRingTrackDistance globalTracks, "
- "richRings, richProjections NOT INITIALIZED"
- << endl;
- return errorVec;
- }
-
- const CbmGlobalTrack* globalTrack = static_cast<const CbmGlobalTrack*>(globalTracks->At(globalTrackId));
- if (globalTrack == NULL) return errorVec;
-
- Int_t stsId = globalTrack->GetStsTrackIndex();
- if (stsId < 0) return errorVec;
-
- const FairTrackParam* pTrack = static_cast<const FairTrackParam*>(richProjections->At(stsId));
- if (pTrack == NULL) return errorVec;
-
- if (pTrack->GetX() == 0 && pTrack->GetY() == 0) return errorVec;
-
- Int_t richId = globalTrack->GetRichRingIndex();
- if (richId < 0) return errorVec;
-
- const CbmRichRing* richRing = static_cast<const CbmRichRing*>(richRings->At(richId));
- if (richRing == NULL) return errorVec;
-
- // Double_t xRing = richRing->GetCenterX();
- // Double_t yRing = richRing->GetCenterY();
- Double_t dx = richRing->GetCenterX() - pTrack->GetX();
- Double_t dy = richRing->GetCenterY() - pTrack->GetY();
-
- Double_t dist = TMath::Sqrt(dx * dx + dy * dy);
-
- vector<Double_t> v = {dist, dx, dy};
- return v;
- }
+ static std::vector<double> GetRingTrackDistanceImpl(int globalTrackId);
};
#endif
diff --git a/macro/rich/geotest/run_many.py b/macro/rich/geotest/run_many.py
index 3253ed41dc5d92711bdd96f76a12efa262d0cfb1..29c0e0b740c8ba9c158ff6aace5e239df8384315 100755
--- a/macro/rich/geotest/run_many.py
+++ b/macro/rich/geotest/run_many.py
@@ -11,4 +11,4 @@ commandXterm3 = ('xterm -hold -e python3 run_one.py {} {}&').format(3, "10gev")
os.system(commandXterm3)
commandXterm4 = ('xterm -hold -e python3 run_one.py {} {}&').format(4, "12gev")
-os.system(commandXterm4)
\ No newline at end of file
+os.system(commandXterm4)
diff --git a/macro/rich/geotest/run_one.py b/macro/rich/geotest/run_one.py
index 5d93adab7daec29bf9ea465357b38a21aacb2930..143ee7bb80ac6c1df0c5e6475580a7da5d798400 100755
--- a/macro/rich/geotest/run_one.py
+++ b/macro/rich/geotest/run_one.py
@@ -32,9 +32,9 @@ def main():
recoCmd = ('root -l -b -q run_reco.C\(\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",{}\)').format(testType, traFile, parFile, digiFile, recoFile, geoSetup, nEvents)
qaCmd = ('root -l -b -q run_qa.C\(\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",{}\)').format(testType, traFile, parFile, digiFile, recoFile, qaFile, geoSetup, resultDir, nEvents)
- os.system((". /{} -a; {}").format(cbmrootConfigPath, traCmd))
- os.system((". /{} -a; {}").format(cbmrootConfigPath, digiCmd))
- os.system((". /{} -a; {}").format(cbmrootConfigPath, recoCmd))
+ # os.system((". /{} -a; {}").format(cbmrootConfigPath, traCmd))
+ # os.system((". /{} -a; {}").format(cbmrootConfigPath, digiCmd))
+ # os.system((". /{} -a; {}").format(cbmrootConfigPath, recoCmd))
os.system((". /{} -a; {}").format(cbmrootConfigPath, qaCmd))
# def make_args():
diff --git a/reco/detectors/rich/qa/CbmRichGeoTest.cxx b/reco/detectors/rich/qa/CbmRichGeoTest.cxx
index 28a08083ead04f23678fa552865b040cc7dbe3d2..b508221a207cb951378c611a9b2c420290032120 100644
--- a/reco/detectors/rich/qa/CbmRichGeoTest.cxx
+++ b/reco/detectors/rich/qa/CbmRichGeoTest.cxx
@@ -20,6 +20,7 @@
#include "CbmMCTrack.h"
#include "CbmReport.h"
#include "CbmRichConverter.h"
+#include "CbmRichDetectorData.h"
#include "CbmRichDigi.h"
#include "CbmRichDigiMapManager.h"
#include "CbmRichDraw.h"
diff --git a/reco/detectors/rich/qa/CbmRichUrqmdTest.cxx b/reco/detectors/rich/qa/CbmRichUrqmdTest.cxx
index 3e0995d89138fb1f61fc7277a3137ba07372ee9d..0c61bb252bb25ec2f64b008307664604838ba0e4 100644
--- a/reco/detectors/rich/qa/CbmRichUrqmdTest.cxx
+++ b/reco/detectors/rich/qa/CbmRichUrqmdTest.cxx
@@ -2,13 +2,6 @@
SPDX-License-Identifier: GPL-3.0-only
Authors: Semen Lebedev [committer], Andrey Lebedev */
-/**
- * \file CbmRichUrqmdTest.cxx
- *
- * \author Semen Lebedev <s.lebedev@gsi.de>
- * \date 2012
- **/
-
#include "CbmRichUrqmdTest.h"
#include "CbmDigiManager.h"
@@ -18,13 +11,15 @@
#include "CbmMCDataManager.h"
#include "CbmMCEventList.h"
#include "CbmMCTrack.h"
-#include "CbmMatchRecoToMC.h"
+#include "CbmRichDetectorData.h"
#include "CbmRichDigi.h"
+#include "CbmRichDigiMapManager.h"
#include "CbmRichDraw.h"
#include "CbmRichGeoManager.h"
#include "CbmRichHit.h"
#include "CbmRichPoint.h"
#include "CbmRichRing.h"
+#include "CbmRichUtil.h"
#include "CbmTrackMatchNew.h"
#include "CbmUtils.h"
@@ -38,70 +33,38 @@
#include "TStyle.h"
#include <TFile.h>
-#include <boost/assign/list_of.hpp>
-
#include <iostream>
#include <sstream>
#include <string>
using namespace std;
-using boost::assign::list_of;
-
-CbmRichUrqmdTest::CbmRichUrqmdTest()
- : FairTask("CbmRichUrqmdTest")
- , fHM(NULL)
- , fOutputDir("")
- , fRichHits(NULL)
- , fRichRings(NULL)
- , fRichPoints(NULL)
- , fMcTracks(NULL)
- , fRichRingMatches(NULL)
- , fRichProjections(NULL)
- , fDigiMan(nullptr)
- , fEventList(NULL)
- , fEventNum(0)
- , fMinNofHits(7)
- , fNofHitsInRingMap()
-{
-}
+using namespace Cbm;
+
+CbmRichUrqmdTest::CbmRichUrqmdTest() : FairTask("CbmRichUrqmdTest") {}
CbmRichUrqmdTest::~CbmRichUrqmdTest() {}
InitStatus CbmRichUrqmdTest::Init()
{
cout << "CbmRichUrqmdTest::Init" << endl;
- FairRootManager* ioman = FairRootManager::Instance();
- if (NULL == ioman) { Fatal("CbmRichUrqmdTest::Init", "RootManager not instantised!"); }
-
- CbmMCDataManager* mcManager = (CbmMCDataManager*) ioman->GetObject("MCDataManager");
- if (mcManager == nullptr) LOG(fatal) << "CbmRichUrqmdTest::Init() NULL MCDataManager.";
-
- fMcTracks = mcManager->InitBranch("MCTrack");
- if (NULL == fMcTracks) { LOG(fatal) << "CbmRichUrqmdTest::Init No MCTrack!"; }
-
- fRichPoints = mcManager->InitBranch("RichPoint");
- if (NULL == fRichPoints) { LOG(fatal) << "CbmRichUrqmdTest::Init No RichPoint!"; }
-
- fRichHits = (TClonesArray*) ioman->GetObject("RichHit");
- if (NULL == fRichHits) { LOG(fatal) << "CbmRichUrqmdTest::Init No RichHit!"; }
-
- fRichRings = (TClonesArray*) ioman->GetObject("RichRing");
- if (NULL == fRichRings) { LOG(fatal) << "CbmRichUrqmdTest::Init No RichRing!"; }
+ fMcTracks = InitOrFatalMc("MCTrack", "CbmRichUrqmdTest::Init");
+ fRichPoints = InitOrFatalMc("RichPoint", "CbmRichUrqmdTest::Init");
+ fRichHits = GetOrFatal<TClonesArray>("RichHit", "CbmRichUrqmdTest::Init");
+ fRichRings = GetOrFatal<TClonesArray>("RichRing", "CbmRichUrqmdTest::Init");
+ fRichRingMatches = GetOrFatal<TClonesArray>("RichRingMatch", "CbmRichUrqmdTest::Init");
+ fRichProjections = GetOrFatal<TClonesArray>("RichProjection", "CbmRichUrqmdTest::Init");
+ fEventList = GetOrFatal<CbmMCEventList>("MCEventList.", "CbmRichUrqmdTest::Init");
fDigiMan = CbmDigiManager::Instance();
fDigiMan->Init();
- fRichRingMatches = (TClonesArray*) ioman->GetObject("RichRingMatch");
- if (NULL == fRichRingMatches) { LOG(fatal) << "CbmRichUrqmdTest::Init No RichRingMatch!"; }
-
- fRichProjections = (TClonesArray*) ioman->GetObject("RichProjection");
- if (NULL == fRichProjections) { LOG(fatal) << "CbmRichUrqmdTest::Init No fRichProjection !"; }
-
- fEventList = (CbmMCEventList*) ioman->GetObject("MCEventList.");
- if (NULL == fEventList) { LOG(fatal) << "CbmRichUrqmdTest::Init No MCEventList!"; }
+ fVertexZStsSlices = {make_pair(0., 5.), make_pair(5., 15.), make_pair(15., 25.), make_pair(25., 35.),
+ make_pair(35., 45.), make_pair(45., 55.), make_pair(55., 65.), make_pair(65., 75.),
+ make_pair(75., 85.), make_pair(85., 95.), make_pair(95., 105.)};
InitHistograms();
+
return kSUCCESS;
}
@@ -111,7 +74,7 @@ void CbmRichUrqmdTest::Exec(Option_t* /*option*/)
cout << "CbmRichUrqmdTest, event No. " << fEventNum << endl;
- FillRichRingNofHits();
+ fNofHitsInRingMap = CbmRichUtil::CreateNofHitsInRingMap(fRichHits, fRichPoints, fMcTracks, fDigiMan);
NofRings();
NofHitsAndPoints();
NofProjections();
@@ -123,53 +86,33 @@ void CbmRichUrqmdTest::InitHistograms()
{
fHM = new CbmHistManager();
- fHM->Create1<TH1D>("fh_vertex_z", "fh_vertex_z;z [cm];# vertices per event", 350, -1., 350);
- fHM->Create1<TH1D>("fh_vertex_z_sts", "fh_vertex_z_sts;z [cm];# vertices per event", 222, -1., 110.);
- fHM->Create2<TH2D>("fh_vertex_xy", "fh_vertex_xy;x [cm];y [cm];# vertices per event", 100, -200., 200., 100, -200.,
- 200.);
- fHM->Create2<TH2D>("fh_vertex_zy", "fh_vertex_zy;z [cm];y [cm];# vertices per event", 350, -1., 350, 100, -200.,
- 200.);
- fHM->Create2<TH2D>("fh_vertex_zx", "fh_vertex_zx;z [cm];x [cm];# vertices per event", 350, -1., 350, 100, -200.,
- 200.);
- fHM->Create2<TH2D>("fh_vertex_xy_z100_180", "fh_vertex_xy_z100_180;x [cm];y [cm];# vertices per event", 100, -200.,
- 200., 100, -200., 200.);
- fHM->Create2<TH2D>("fh_vertex_xy_z180_370", "fh_vertex_xy_z180_370;x [cm];y [cm];# vertices per event", 100, -200.,
- 200., 100, -200., 200.);
- fHM->Create2<TH2D>("fh_vertex_xy_z180_230", "fh_vertex_xy_z180_230;x [cm];y [cm];# vertices per event", 100, -200.,
- 200., 100, -200., 200.);
-
- fHM->Create2<TH2D>("fh_vertex_xy_z5", "fh_vertex_xy_z5;x [cm];y [cm];# vertices per event", 100, -100., 100., 100,
- -100., 100.);
- fHM->Create2<TH2D>("fh_vertex_xy_z5_15", "fh_vertex_xy_z5_15;x [cm];y [cm];# vertices per event", 100, -100., 100.,
- 100, -100., 100.);
- fHM->Create2<TH2D>("fh_vertex_xy_z15_25", "fh_vertex_xy_z15_25;x [cm];y [cm];# vertices per event", 100, -100., 100.,
- 100, -100., 100.);
- fHM->Create2<TH2D>("fh_vertex_xy_z25_35", "fh_vertex_xy_z25_35;x [cm];y [cm];# vertices per event", 100, -100., 100.,
- 100, -100., 100.);
- fHM->Create2<TH2D>("fh_vertex_xy_z35_45", "fh_vertex_xy_z35_45;x [cm];y [cm];# vertices per event", 100, -100., 100.,
- 100, -100., 100.);
- fHM->Create2<TH2D>("fh_vertex_xy_z45_55", "fh_vertex_xy_z45_55;x [cm];y [cm];# vertices per event", 100, -100., 100.,
- 100, -100., 100.);
- fHM->Create2<TH2D>("fh_vertex_xy_z55_65", "fh_vertex_xy_z55_65;x [cm];y [cm];# vertices per event", 100, -100., 100.,
- 100, -100., 100.);
- fHM->Create2<TH2D>("fh_vertex_xy_z65_75", "fh_vertex_xy_z65_75;x [cm];y [cm];# vertices per event", 100, -100., 100.,
- 100, -100., 100.);
- fHM->Create2<TH2D>("fh_vertex_xy_z75_85", "fh_vertex_xy_z75_85;x [cm];y [cm];# vertices per event", 100, -100., 100.,
- 100, -100., 100.);
- fHM->Create2<TH2D>("fh_vertex_xy_z85_95", "fh_vertex_xy_z85_95;x [cm];y [cm];# vertices per event", 100, -100., 100.,
- 100, -100., 100.);
- fHM->Create2<TH2D>("fh_vertex_xy_z95_105", "fh_vertex_xy_z95_105;x [cm];y [cm];# vertices per event", 100, -100.,
- 100., 100, -100., 100.);
-
- fHM->Create1<TH1D>("fh_nof_rings_1hit", "fh_nof_rings_1hit;# detected particles/event;Yield", 250, -.5, 249.5);
- fHM->Create1<TH1D>("fh_nof_rings_7hits", "fh_nof_rings_7hits;# detected particles/event;Yield", 250, -.5, 249.5);
- fHM->Create1<TH1D>("fh_nof_rings_prim_1hit", "fh_nof_rings_prim_1hit;# detected particles/event;Yield", 50, -.5,
- 69.5);
- fHM->Create1<TH1D>("fh_nof_rings_prim_7hits", "fh_nof_rings_prim_7hits;# detected particles/event;Yield", 50, -.5,
+ fHM->Create1<TH1D>("fh_vertex_z", "fh_vertex_z;z [cm];# vertices/ev.", 350, -1., 350);
+ fHM->Create1<TH1D>("fh_vertex_z_sts", "fh_vertex_z_sts;z [cm];# vertices/ev.", 222, -1., 110.);
+ fHM->Create2<TH2D>("fh_vertex_xy", "fh_vertex_xy;x [cm];y [cm];# vertices/ev.", 100, -200., 200., 100, -200., 200.);
+ fHM->Create2<TH2D>("fh_vertex_zy", "fh_vertex_zy;z [cm];y [cm];# vertices/ev.", 350, -1., 350, 100, -200., 200.);
+ fHM->Create2<TH2D>("fh_vertex_zx", "fh_vertex_zx;z [cm];x [cm];# vertices/ev.", 350, -1., 350, 100, -200., 200.);
+
+
+ fHM->Create2<TH2D>("fh_vertex_xy_z100_180", "fh_vertex_xy_z100_180;x [cm];y [cm];# vertices/ev.", 100, -200., 200.,
+ 100, -200., 200.);
+ fHM->Create2<TH2D>("fh_vertex_xy_z180_370", "fh_vertex_xy_z180_370;x [cm];y [cm];# vertices/ev.", 100, -200., 200.,
+ 100, -200., 200.);
+ fHM->Create2<TH2D>("fh_vertex_xy_z180_230", "fh_vertex_xy_z180_230;x [cm];y [cm];# vertices/ev.", 100, -200., 200.,
+ 100, -200., 200.);
+
+ for (auto pair : fVertexZStsSlices) {
+ string name = "fh_vertex_xy_z" + to_string(pair.first) + "_" + to_string(pair.second);
+ fHM->Create2<TH2D>(name, name + ";x [cm];y [cm];# vertices/ev.", 100, -100., 100., 100, -100., 100.);
+ }
+
+ fHM->Create1<TH1D>("fh_nof_rings_1hit", "fh_nof_rings_1hit;# detected particles/ev.;Yield", 250, -.5, 249.5);
+ fHM->Create1<TH1D>("fh_nof_rings_7hits", "fh_nof_rings_7hits;# detected particles/ev.;Yield", 250, -.5, 249.5);
+ fHM->Create1<TH1D>("fh_nof_rings_prim_1hit", "fh_nof_rings_prim_1hit;# detected particles/ev.;Yield", 50, -.5, 69.5);
+ fHM->Create1<TH1D>("fh_nof_rings_prim_7hits", "fh_nof_rings_prim_7hits;# detected particles/ev.;Yield", 50, -.5,
69.5);
- fHM->Create1<TH1D>("fh_nof_rings_target_1hit", "fh_nof_rings_target_1hit;# detected particles/event;Yield", 60, -.5,
+ fHM->Create1<TH1D>("fh_nof_rings_target_1hit", "fh_nof_rings_target_1hit;# detected particles/ev.;Yield", 60, -.5,
79.5);
- fHM->Create1<TH1D>("fh_nof_rings_target_7hits", "fh_nof_rings_target_7hits;# detected particles/event;Yield", 60, -.5,
+ fHM->Create1<TH1D>("fh_nof_rings_target_7hits", "fh_nof_rings_target_7hits;# detected particles/ev.;Yield", 60, -.5,
79.5);
fHM->Create1<TH1D>("fh_secel_mom", "fh_secel_mom;p [GeV/c];Number per event", 100, 0., 20);
@@ -179,80 +122,31 @@ void CbmRichUrqmdTest::InitHistograms()
fHM->Create1<TH1D>("fh_kaon_mom", "fh_kaon_mom;p [GeV/c];Number per event", 100, 0., 20);
fHM->Create1<TH1D>("fh_mu_mom", "fh_mu_mom;p [GeV/c];Number per event", 100, 0., 20);
- fHM->Create1<TH1D>("fh_nof_points_per_event", "fh_nof_points_per_event;Particle;# MC points per event", 7, .5, 7.5);
- fHM->Create1<TH1D>("fh_nof_hits_per_event", "fh_nof_hits_per_event;# hits per event;Yield", 100, 0, 2000);
+ fHM->Create1<TH1D>("fh_nof_points_per_event", "fh_nof_points_per_event;Particle;# MC points/ev.", 7, .5, 7.5);
+ fHM->Create1<TH1D>("fh_nof_hits_per_event", "fh_nof_hits_per_event;# hits per event;Yield", 200, 0, 2000);
fHM->Create1<TH1D>("fh_nof_hits_per_pmt", "fh_nof_hits_per_pmt;# hits per PMT;% of total", 65, -0.5, 64.5);
- vector<Double_t> xPmtBins = CbmRichDraw::GetPmtHistXbins();
- vector<Double_t> yPmtBins = CbmRichDraw::GetPmtHistYbins();
+ vector<double> xPmtBins = CbmRichUtil::GetPmtHistXbins();
+ vector<double> yPmtBins = CbmRichUtil::GetPmtHistYbins();
// before drawing must be normalized by 1/64
- TH2D* fh_hitrate_xy = new TH2D("fh_hitrate_xy", "fh_hitrate_xy;X [cm];Y [cm];# hits/pixel/s", xPmtBins.size() - 1,
- &xPmtBins[0], yPmtBins.size() - 1, &yPmtBins[0]);
- fHM->Add("fh_hitrate_xy", fh_hitrate_xy);
-
- TH2D* fh_hits_xy = new TH2D("fh_hits_xy", "fh_hits_xy;X [cm];Y [cm];# hits/PMT/event", xPmtBins.size() - 1,
- &xPmtBins[0], yPmtBins.size() - 1, &yPmtBins[0]);
- fHM->Add("fh_hits_xy", fh_hits_xy);
-
- TH2D* fh_points_xy = new TH2D("fh_points_xy", "fh_points_xy;X [cm];Y [cm];# MC points/PMT/event", xPmtBins.size() - 1,
- &xPmtBins[0], yPmtBins.size() - 1, &yPmtBins[0]);
- fHM->Add("fh_points_xy", fh_points_xy);
-
- TH2D* fh_points_xy_pions = new TH2D("fh_points_xy_pions", "fh_points_xy_pions;X [cm];Y [cm];# MC points/PMT/event",
- xPmtBins.size() - 1, &xPmtBins[0], yPmtBins.size() - 1, &yPmtBins[0]);
- fHM->Add("fh_points_xy_pions", fh_points_xy_pions);
-
- TH2D* fh_points_xy_gamma_target =
- new TH2D("fh_points_xy_gamma_target", "fh_points_xy_gamma_target;X [cm];Y [cm];# MC points/PMT/event",
- xPmtBins.size() - 1, &xPmtBins[0], yPmtBins.size() - 1, &yPmtBins[0]);
- fHM->Add("fh_points_xy_gamma_target", fh_points_xy_gamma_target);
-
- TH2D* fh_points_xy_gamma_nontarget =
- new TH2D("fh_points_xy_gamma_nontarget", "fh_points_xy_gamma_nontarget;X [cm];Y [cm];# MC points/PMT/event",
- xPmtBins.size() - 1, &xPmtBins[0], yPmtBins.size() - 1, &yPmtBins[0]);
- fHM->Add("fh_points_xy_gamma_nontarget", fh_points_xy_gamma_nontarget);
-
- TH2D* fh_skipped_pmt_10_xy =
- new TH2D("fh_skipped_pmt_10_xy", "fh_skipped_pmt_10_xy;X [cm];Y [cm];# skipped PMTs (>10 hits) [%]",
- xPmtBins.size() - 1, &xPmtBins[0], yPmtBins.size() - 1, &yPmtBins[0]);
- fHM->Add("fh_skipped_pmt_10_xy", fh_skipped_pmt_10_xy);
-
- TH2D* fh_skipped_pmt_20_xy =
- new TH2D("fh_skipped_pmt_20_xy", "fh_skipped_pmt_20_xy;X [cm];Y [cm];# skipped PMTs (>20 hits) [%]",
- xPmtBins.size() - 1, &xPmtBins[0], yPmtBins.size() - 1, &yPmtBins[0]);
- fHM->Add("fh_skipped_pmt_20_xy", fh_skipped_pmt_20_xy);
-
- TH2D* fh_skipped_pmt_30_xy =
- new TH2D("fh_skipped_pmt_30_xy", "fh_skipped_pmt_30_xy;X [cm];Y [cm];# skipped PMTs (>30 hits) [%]",
- xPmtBins.size() - 1, &xPmtBins[0], yPmtBins.size() - 1, &yPmtBins[0]);
- fHM->Add("fh_skipped_pmt_30_xy", fh_skipped_pmt_30_xy);
-
- fHM->Create1<TH1D>("fh_nof_proj_per_event", "fh_nof_proj_per_event;# tracks per event;Yield", 50, 0, 1000);
- fHM->Create2<TH2D>("fh_proj_xy", "fh_proj_xy;X [cm];Y [cm];# tracks/cm^{2}/event", 240, -120, 120, 420, -210, 210);
-}
-
-void CbmRichUrqmdTest::FillRichRingNofHits()
-{
- fNofHitsInRingMap.clear();
- Int_t nofRichHits = fRichHits->GetEntriesFast();
- for (Int_t iHit = 0; iHit < nofRichHits; iHit++) {
- CbmRichHit* hit = static_cast<CbmRichHit*>(fRichHits->At(iHit));
- if (NULL == hit) continue;
- Int_t digiIndex = hit->GetRefId();
- if (digiIndex < 0) continue;
- const CbmRichDigi* digi = fDigiMan->Get<CbmRichDigi>(digiIndex);
- if (NULL == digi) continue;
- const CbmMatch* digiMatch = fDigiMan->GetMatch(ECbmModuleId::kRich, digiIndex);
- if (NULL == digiMatch) continue;
- Int_t eventId = digiMatch->GetMatchedLink().GetEntry();
-
- vector<pair<Int_t, Int_t>> motherIds =
- CbmMatchRecoToMC::GetMcTrackMotherIdsForRichHit(fDigiMan, hit, fRichPoints, fMcTracks, eventId);
- for (UInt_t i = 0; i < motherIds.size(); i++) {
- fNofHitsInRingMap[motherIds[i]]++;
- }
- }
+ fHM->Create2<TH2D>("fh_hitrate_xy", "fh_hitrate_xy;X [cm];Y [cm];# hits/pixel/s", xPmtBins, yPmtBins);
+ fHM->Create2<TH2D>("fh_hits_xy", "fh_hits_xy;X [cm];Y [cm];# hits/PMT/ev.", xPmtBins, yPmtBins);
+ fHM->Create2<TH2D>("fh_points_xy", "fh_points_xy;X [cm];Y [cm];# MC points/PMT/ev.", xPmtBins, yPmtBins);
+ fHM->Create2<TH2D>("fh_points_xy_pions", "fh_points_xy_pions;X [cm];Y [cm];# MC points/PMT/ev.", xPmtBins, yPmtBins);
+ fHM->Create2<TH2D>("fh_points_xy_gamma_target", "fh_points_xy_gamma_target;X [cm];Y [cm];# MC points/PMT/ev.",
+ xPmtBins, yPmtBins);
+ fHM->Create2<TH2D>("fh_points_xy_gamma_nontarget", "fh_points_xy_gamma_nontarget;X [cm];Y [cm];# MC points/PMT/ev.",
+ xPmtBins, yPmtBins);
+ fHM->Create2<TH2D>("fh_skipped_pmt_10_xy", "fh_skipped_pmt_10_xy;X [cm];Y [cm];# skipped PMTs (>10 hits) [%]",
+ xPmtBins, yPmtBins);
+ fHM->Create2<TH2D>("fh_skipped_pmt_20_xy", "fh_skipped_pmt_20_xy;X [cm];Y [cm];# skipped PMTs (>20 hits) [%]",
+ xPmtBins, yPmtBins);
+ fHM->Create2<TH2D>("fh_skipped_pmt_30_xy", "fh_skipped_pmt_30_xy;X [cm];Y [cm];# skipped PMTs (>30 hits) [%]",
+ xPmtBins, yPmtBins);
+
+ fHM->Create1<TH1D>("fh_nof_proj_per_event", "fh_nof_proj_per_event;# tracks/ev.;Yield", 50, 0, 1000);
+ fHM->Create2<TH2D>("fh_proj_xy", "fh_proj_xy;X [cm];Y [cm];# tracks/cm^{2}/ev.", 240, -120, 120, 420, -210, 210);
}
void CbmRichUrqmdTest::NofRings()
@@ -487,17 +381,10 @@ void CbmRichUrqmdTest::Vertex()
if (v.Z() >= 180 && v.Z() <= 370) fHM->H2("fh_vertex_xy_z180_370")->Fill(v.X(), v.Y());
if (v.Z() >= 180 && v.Z() <= 230) fHM->H2("fh_vertex_xy_z180_230")->Fill(v.X(), v.Y());
- if (v.Z() <= 5) fHM->H2("fh_vertex_xy_z5")->Fill(v.X(), v.Y());
- if (v.Z() > 5 && v.Z() <= 15) fHM->H2("fh_vertex_xy_z5_15")->Fill(v.X(), v.Y());
- if (v.Z() > 15 && v.Z() <= 25) fHM->H2("fh_vertex_xy_z15_25")->Fill(v.X(), v.Y());
- if (v.Z() > 25 && v.Z() <= 35) fHM->H2("fh_vertex_xy_z25_35")->Fill(v.X(), v.Y());
- if (v.Z() > 35 && v.Z() <= 45) fHM->H2("fh_vertex_xy_z35_45")->Fill(v.X(), v.Y());
- if (v.Z() > 45 && v.Z() <= 55) fHM->H2("fh_vertex_xy_z45_55")->Fill(v.X(), v.Y());
- if (v.Z() > 55 && v.Z() <= 65) fHM->H2("fh_vertex_xy_z55_65")->Fill(v.X(), v.Y());
- if (v.Z() > 65 && v.Z() <= 75) fHM->H2("fh_vertex_xy_z65_75")->Fill(v.X(), v.Y());
- if (v.Z() > 75 && v.Z() <= 85) fHM->H2("fh_vertex_xy_z75_85")->Fill(v.X(), v.Y());
- if (v.Z() > 85 && v.Z() <= 95) fHM->H2("fh_vertex_xy_z85_95")->Fill(v.X(), v.Y());
- if (v.Z() > 95 && v.Z() <= 105) fHM->H2("fh_vertex_xy_z95_105")->Fill(v.X(), v.Y());
+ for (auto pair : fVertexZStsSlices) {
+ string name = "fh_vertex_xy_z" + to_string(pair.first) + "_" + to_string(pair.second);
+ if (v.Z() > pair.first && v.Z() <= pair.second) { fHM->H2(name)->Fill(v.X(), v.Y()); }
+ }
}
}
}
@@ -538,54 +425,16 @@ void CbmRichUrqmdTest::DrawHist()
{
gStyle->SetOptTitle(1);
- fHM->H2("fh_vertex_xy_z5")->Scale(1. / fEventNum);
- fHM->H2("fh_vertex_xy_z5_15")->Scale(1. / fEventNum);
- fHM->H2("fh_vertex_xy_z15_25")->Scale(1. / fEventNum);
- fHM->H2("fh_vertex_xy_z25_35")->Scale(1. / fEventNum);
- fHM->H2("fh_vertex_xy_z35_45")->Scale(1. / fEventNum);
- fHM->H2("fh_vertex_xy_z45_55")->Scale(1. / fEventNum);
- fHM->H2("fh_vertex_xy_z55_65")->Scale(1. / fEventNum);
- fHM->H2("fh_vertex_xy_z65_75")->Scale(1. / fEventNum);
- fHM->H2("fh_vertex_xy_z75_85")->Scale(1. / fEventNum);
- fHM->H2("fh_vertex_xy_z85_95")->Scale(1. / fEventNum);
- fHM->H2("fh_vertex_xy_z95_105")->Scale(1. / fEventNum);
-
-
TCanvas* c = fHM->CreateCanvas("rich_urqmd_vertex_sts_xyz", "rich_urqmd_vertex_sts_xyz", 1600, 1200);
c->Divide(4, 3);
- c->cd(1);
- DrawH2(fHM->H2("fh_vertex_xy_z5"), kLinear, kLinear, kLog);
- DrawTextOnPad("Z < 5 cm", 0.3, 0.9, 0.7, 0.98);
- c->cd(2);
- DrawH2(fHM->H2("fh_vertex_xy_z5_15"), kLinear, kLinear, kLog);
- DrawTextOnPad("5 < Z < 15 cm", 0.3, 0.9, 0.7, 0.98);
- c->cd(3);
- DrawH2(fHM->H2("fh_vertex_xy_z15_25"), kLinear, kLinear, kLog);
- DrawTextOnPad("15 < Z < 25 cm", 0.3, 0.9, 0.7, 0.98);
- c->cd(4);
- DrawH2(fHM->H2("fh_vertex_xy_z25_35"), kLinear, kLinear, kLog);
- DrawTextOnPad("25 < Z < 35 cm", 0.3, 0.9, 0.7, 0.98);
- c->cd(5);
- DrawH2(fHM->H2("fh_vertex_xy_z35_45"), kLinear, kLinear, kLog);
- DrawTextOnPad("35 < Z < 45 cm", 0.3, 0.9, 0.7, 0.98);
- c->cd(6);
- DrawH2(fHM->H2("fh_vertex_xy_z45_55"), kLinear, kLinear, kLog);
- DrawTextOnPad("45 < Z < 55 cm", 0.3, 0.9, 0.7, 0.98);
- c->cd(7);
- DrawH2(fHM->H2("fh_vertex_xy_z55_65"), kLinear, kLinear, kLog);
- DrawTextOnPad("55 < Z < 65 cm", 0.3, 0.9, 0.7, 0.98);
- c->cd(8);
- DrawH2(fHM->H2("fh_vertex_xy_z65_75"), kLinear, kLinear, kLog);
- DrawTextOnPad("65 < Z < 75 cm", 0.3, 0.9, 0.7, 0.98);
- c->cd(9);
- DrawH2(fHM->H2("fh_vertex_xy_z75_85"), kLinear, kLinear, kLog);
- DrawTextOnPad("75 < Z < 85 cm", 0.3, 0.9, 0.7, 0.98);
- c->cd(10);
- DrawH2(fHM->H2("fh_vertex_xy_z85_95"), kLinear, kLinear, kLog);
- DrawTextOnPad("85 < Z < 95 cm", 0.3, 0.9, 0.7, 0.98);
- c->cd(11);
- DrawH2(fHM->H2("fh_vertex_xy_z95_105"), kLinear, kLinear, kLog);
- DrawTextOnPad("95 < Z < 105 cm", 0.3, 0.9, 0.7, 0.98);
+ int i = 1;
+ for (auto pair : fVertexZStsSlices) {
+ string name = "fh_vertex_xy_z" + to_string(pair.first) + "_" + to_string(pair.second);
+ fHM->H2(name)->Scale(1. / fEventNum);
+ c->cd(i++);
+ DrawH2(fHM->H2(name), kLinear, kLinear, kLog);
+ DrawTextOnPad(to_string(pair.first) + " cm < Z < " + to_string(pair.second) + " cm", 0.3, 0.9, 0.7, 0.98);
+ }
gStyle->SetOptTitle(0);
}
@@ -677,8 +526,8 @@ void CbmRichUrqmdTest::DrawHist()
ss6 << "#mu^{#pm} (" << fHM->H1("fh_mu_mom")->GetEntries() / fEventNum << ")";
DrawH1({fHM->H1("fh_gamma_target_mom"), fHM->H1("fh_gamma_nontarget_mom"), fHM->H1("fh_secel_mom"),
fHM->H1("fh_pi_mom"), fHM->H1("fh_kaon_mom"), fHM->H1("fh_mu_mom")},
- list_of(ss1.str())(ss2.str())(ss3.str()) (ss4.str()) (ss5.str()) (ss6.str()), kLinear, kLog, true, 0.5, 0.7,
- 0.99, 0.99, "hist");
+ {ss1.str(), ss2.str(), ss3.str(), ss4.str(), ss5.str(), ss6.str()}, kLinear, kLog, true, 0.5, 0.7, 0.99,
+ 0.99, "hist");
}
{
diff --git a/reco/detectors/rich/qa/CbmRichUrqmdTest.h b/reco/detectors/rich/qa/CbmRichUrqmdTest.h
index 848e57dca6364df54e3874fd4cdefac76a7a943f..26f5b142ff0d6637f0ead922639b85e88c89d8df 100644
--- a/reco/detectors/rich/qa/CbmRichUrqmdTest.h
+++ b/reco/detectors/rich/qa/CbmRichUrqmdTest.h
@@ -1,16 +1,7 @@
-/* Copyright (C) 2012-2020 UGiessen/JINR-LIT, Giessen/Dubna
+/* Copyright (C) 2012-2021 UGiessen/JINR-LIT, Giessen/Dubna
SPDX-License-Identifier: GPL-3.0-only
Authors: Andrey Lebedev, Semen Lebedev [committer] */
-/**
- * \file CbmRichUrqmdTest.h
- *
- * \brief RICH geometry testing in Urqmd collisions.
- *
- * \author Semen Lebedev <s.lebedev@gsi.de>
- * \date 2012
- **/
-
#ifndef CBM_RICH_URQMD_TEST
#define CBM_RICH_URQMD_TEST
@@ -26,16 +17,7 @@ class CbmDigiManager;
#include <map>
#include <vector>
-using namespace std;
-/**
- * \class CbmRichUrqmdTest
- *
- * \brief RICH geometry testing in Urqmd collisions.
- *
- * \author Semen Lebedev <s.lebedev@gsi.de>
- * \date 2012
- **/
class CbmRichUrqmdTest : public FairTask {
public:
@@ -69,20 +51,35 @@ public:
* \brief Set output directory where you want to write results (figures and json).
* \param[in] dir Path to the output directory.
*/
- void SetOutputDir(const string& dir) { fOutputDir = dir; }
+ void SetOutputDir(const std::string& dir) { fOutputDir = dir; }
private:
+ CbmHistManager* fHM = nullptr;
+
+ std::string fOutputDir = ""; // output dir for results
+
+ TClonesArray* fRichHits = nullptr;
+ TClonesArray* fRichRings = nullptr;
+ CbmMCDataArray* fRichPoints = nullptr;
+ CbmMCDataArray* fMcTracks = nullptr;
+ TClonesArray* fRichRingMatches = nullptr;
+ TClonesArray* fRichProjections = nullptr;
+ CbmDigiManager* fDigiMan = nullptr;
+ CbmMCEventList* fEventList = nullptr;
+
+ Int_t fEventNum = 0;
+ Int_t fMinNofHits = 7; // Min number of hits in ring for detector acceptance calculation.
+
+ std::map<std::pair<Int_t, Int_t>, Int_t> fNofHitsInRingMap; // Number of hits in the MC RICH ring
+
+ std::vector<std::pair<int, int>> fVertexZStsSlices; // Slices (Zmin - Zmax) inside STS
+
/**
* \brief Initialize histograms.
*/
void InitHistograms();
- /**
- * \brief
- */
- void FillRichRingNofHits();
-
/**
* \brief
*/
@@ -121,25 +118,6 @@ private:
*/
CbmRichUrqmdTest& operator=(const CbmRichUrqmdTest&);
- CbmHistManager* fHM;
-
- string fOutputDir; // output dir for results
-
- TClonesArray* fRichHits;
- TClonesArray* fRichRings;
- CbmMCDataArray* fRichPoints;
- CbmMCDataArray* fMcTracks;
- TClonesArray* fRichRingMatches;
- TClonesArray* fRichProjections;
- CbmDigiManager* fDigiMan;
- CbmMCEventList* fEventList;
-
- Int_t fEventNum;
- Int_t fMinNofHits; // Min number of hits in ring for detector acceptance calculation.
-
- // Number of hits in the MC RICH ring
- std::map<pair<Int_t, Int_t>, Int_t> fNofHitsInRingMap;
-
ClassDef(CbmRichUrqmdTest, 1)
};