From a8601a15398943a1f934fdc6ac3347c8784389ea Mon Sep 17 00:00:00 2001
From: "se.gorbunov" <se.gorbunov@gsi.de>
Date: Thu, 25 Nov 2021 12:34:12 +0000
Subject: [PATCH] tracker input QA for TRD
---
core/data/CbmLink.h | 2 +
core/data/CbmMCEventList.cxx | 11 +
core/data/CbmMCEventList.h | 18 +
core/qa/CMakeLists.txt | 6 +
core/qa/CbmQaBaseLinkDef.h | 12 +-
core/qa/CbmQaHist.cxx | 5 +-
core/qa/CbmQaHist.h | 94 ++--
macro/run/run_qa.C | 5 +-
reco/L1/CMakeLists.txt | 8 +-
reco/L1/L1LinkDef.h | 2 +
reco/L1/qa/CbmTrackerInputQaTrd.cxx | 739 ++++++++++++++++++++++++++++
reco/L1/qa/CbmTrackerInputQaTrd.h | 166 +++++++
12 files changed, 1024 insertions(+), 44 deletions(-)
create mode 100644 reco/L1/qa/CbmTrackerInputQaTrd.cxx
create mode 100644 reco/L1/qa/CbmTrackerInputQaTrd.h
diff --git a/core/data/CbmLink.h b/core/data/CbmLink.h
index 2823928c4c..1ca2563482 100644
--- a/core/data/CbmLink.h
+++ b/core/data/CbmLink.h
@@ -61,6 +61,8 @@ public:
return (lhs.GetFile() == rhs.GetFile() && lhs.GetEntry() == rhs.GetEntry() && lhs.GetIndex() == rhs.GetIndex());
}
+ friend Bool_t operator!=(const CbmLink& lhs, const CbmLink& rhs) { return !(lhs == rhs); }
+
/** Comparison operators by //Dr.Sys **/
friend Bool_t operator<(const CbmLink& l, const CbmLink& r)
{
diff --git a/core/data/CbmMCEventList.cxx b/core/data/CbmMCEventList.cxx
index 054fa8dc7a..8dd8e851ae 100644
--- a/core/data/CbmMCEventList.cxx
+++ b/core/data/CbmMCEventList.cxx
@@ -97,6 +97,17 @@ Double_t CbmMCEventList::GetEventTime(UInt_t eventId, UInt_t fileId)
// ----------------------------------------------------------------------------
+// ----- Get event index in the list -------------------------------------
+Int_t CbmMCEventList::GetEventIndex(UInt_t eventId, UInt_t fileId)
+{
+ if (!fIsSorted) Sort();
+ auto it = Find(fileId, eventId);
+ if (it == fEvents.end()) return -1.;
+ return (it - fEvents.begin());
+}
+// ----------------------------------------------------------------------------
+
+
// ----- Get event time for event at index in list ------------------------
Double_t CbmMCEventList::GetEventTimeByIndex(UInt_t index)
{
diff --git a/core/data/CbmMCEventList.h b/core/data/CbmMCEventList.h
index 82e4d99267..ae5ef757aa 100644
--- a/core/data/CbmMCEventList.h
+++ b/core/data/CbmMCEventList.h
@@ -10,6 +10,7 @@
#ifndef CBMMCEVENTLIST_H
#define CBMMCEVENTLIST_H 1
+#include "CbmLink.h"
#include "CbmMCEventInfo.h" // for CbmMCEventInfo
#include <Rtypes.h> // for THashConsistencyHolder, ClassDef
@@ -72,6 +73,10 @@ public:
**/
Double_t GetEventTime(UInt_t event, UInt_t file);
+ /** @brief Event start time
+ * * @param link link to the MC event
+ **/
+ Double_t GetEventTime(const CbmLink& link) { return GetEventTime(link.GetEntry(), link.GetFile()); }
/** @brief Event time by index
** @value Event time for event at given index in list
@@ -88,6 +93,19 @@ public:
**/
Int_t GetFileIdByIndex(UInt_t index);
+ /** @brief Event index
+ ** @param event MC event number
+ ** @param file MC input file number
+ ** @value eventindex in list
+ **
+ ** Returns -1. if the event is not present in the list.
+ **/
+ Int_t GetEventIndex(UInt_t event, UInt_t file);
+
+ /** @brief Event index
+ **/
+ Int_t GetEventIndex(const CbmLink& link) { return GetEventIndex(link.GetEntry(), link.GetFile()); }
+
/** @brief Number of events in the list
** @value Number of events
diff --git a/core/qa/CMakeLists.txt b/core/qa/CMakeLists.txt
index 80589aa5a5..aaa8e0a0f2 100644
--- a/core/qa/CMakeLists.txt
+++ b/core/qa/CMakeLists.txt
@@ -24,6 +24,12 @@ CbmQaPie.cxx
CbmQaHist.cxx
)
+set(HEADERS
+CbmQaCanvas.h
+CbmQaPie.h
+CbmQaHist.h
+)
+
set(LINKDEF CbmQaBaseLinkDef.h)
Set(LIBRARY_NAME CbmQaBase)
Set(DEPENDENCIES Hist Gpad)
diff --git a/core/qa/CbmQaBaseLinkDef.h b/core/qa/CbmQaBaseLinkDef.h
index 06d82d7b0e..c17d1514d6 100644
--- a/core/qa/CbmQaBaseLinkDef.h
+++ b/core/qa/CbmQaBaseLinkDef.h
@@ -11,10 +11,16 @@
#pragma link C++ nestedclasses;
#pragma link C++ nestedtypedef;
+// custom streamers are already implemented in the class body
#pragma link C++ class CbmQaCanvas - ;
-#pragma link C++ class CbmQaPieSlice + ;
#pragma link C++ class CbmQaPie - ;
-#pragma link C++ class CbmQaH1F + ;
-#pragma link C++ class CbmQaH1D + ;
+
+// create streamers automatically
+#pragma link C++ class CbmQaPieSlice + ;
+#pragma link C++ class CbmQaHist < TH1F> + ;
+#pragma link C++ class CbmQaHist < TH1D> + ;
+#pragma link C++ class CbmQaHist < TH1I> + ;
+#pragma link C++ class CbmQaHist < TProfile> + ;
+#pragma link C++ class CbmQaHist < TProfile2D> + ;
#endif
diff --git a/core/qa/CbmQaHist.cxx b/core/qa/CbmQaHist.cxx
index 0651387e7f..5e46f88cf0 100644
--- a/core/qa/CbmQaHist.cxx
+++ b/core/qa/CbmQaHist.cxx
@@ -13,6 +13,9 @@ templateClassImp(CbmQaHist);
// create definitions of specific CbmQaHist classes,
// otherwise they will not be linked by the ROOT linker
-//
+
template class CbmQaHist<TH1F>;
template class CbmQaHist<TH1D>;
+template class CbmQaHist<TH1I>;
+template class CbmQaHist<TProfile>;
+template class CbmQaHist<TProfile2D>;
diff --git a/core/qa/CbmQaHist.h b/core/qa/CbmQaHist.h
index e0e609f053..860e9d1e8b 100644
--- a/core/qa/CbmQaHist.h
+++ b/core/qa/CbmQaHist.h
@@ -12,13 +12,24 @@
#include "CbmQaCanvas.h"
-#include "TFitResultPtr.h"
-#include "TH1D.h"
-#include "TH1F.h"
-#include "TStyle.h"
-#include "TVirtualPad.h"
-/// A modification of TH* classes that keeps statistics & fit drawing options
-/// and resizes the stat window accordingly without actual drawing.
+#include <FairLogger.h>
+
+#include <TFitResultPtr.h>
+#include <TPaveStats.h>
+#include <TStyle.h>
+#include <TVirtualPad.h>
+
+// The TH* headers are needed here for the ROOT linker
+
+#include <TF1.h>
+#include <TH1D.h>
+#include <TH1F.h>
+#include <TH1I.h>
+#include <TProfile.h>
+#include <TProfile2D.h>
+
+/// A modification of TH* and TProfile* classes that keeps statistics & fit drawing options
+/// and resizes the stat window accordingly without the actual drawing.
/// In the original classes, hist->Draw() & canv->Update() must be called
/// for resetting Stat/Fit window.
///
@@ -34,6 +45,7 @@ public:
fOptStat = gStyle->GetOptStat();
fOptFit = gStyle->GetOptFit();
}
+ this->SetLineWidth(2);
}
/// Copy constructor
@@ -45,6 +57,7 @@ public:
CbmQaHist(Types... args) : HistTypeT(args...)
{
if (gStyle) { SetOptStatFit(gStyle->GetOptStat(), gStyle->GetOptFit()); }
+ this->SetLineWidth(2);
}
/// Destructor
@@ -57,54 +70,63 @@ public:
TFitResultPtr Fit(Types... args)
{
TVirtualPad* padsav = gPad;
- GetCanvas().cd();
+ GetDummyCanvas().cd();
+ this->Sumw2();
auto ret = HistTypeT::Fit(args...);
- GetCanvas().Clear();
+ GetDummyCanvas().Clear();
+
+ // make the output look nice
+
if (padsav) padsav->cd();
auto* f = this->GetFunction("gaus");
if (f) {
- f->SetParName(0, "Peak");
- f->SetParName(1, "#mu");
- f->SetParName(2, "#sigma");
+ f->SetParNames("Peak", "#mu", "#sigma");
f->SetLineColor(kRed);
f->SetLineWidth(3);
+ TPaveStats* st = (TPaveStats*) this->FindObject("stats");
+ if (!st) { LOG(fatal) << "CbmQaHist: can not access histogram statistics"; }
+ else {
+ st->SetX1NDC(0.6);
+ st->SetX2NDC(0.940);
+ st->SetY1NDC(0.5);
+ st->SetY2NDC(0.930);
+ st->SetOptStat(111110);
+ st->SetOptFit(10001);
+ }
}
return ret;
}
- /// Set stat drawing options and resize the stat window
+ /// Set stat drawing options and autoresize the stat window
void SetOptStat(Int_t stat = 1) { SetOptStatFit(stat, fOptFit); }
- /// Set fit drawing options and resize the stat window
+ /// Set fit drawing options and autoresize the stat window
void SetOptFit(Int_t fit = 1) { SetOptStatFit(fOptStat, fit); }
- /// Set stat & fit drawing options and resize the stat window
+ /// Set stat & fit drawing options and autoresize the stat window
void SetOptStatFit(int stat, int fit)
{
+ // the only way to create and auto-size the stat window is to draw the histogram
+
fOptStat = stat;
fOptFit = fit;
- if (!gStyle) { return; }
+ if (!gStyle) { return; } // should not happen
+
TVirtualPad* savePad = gPad;
int saveStat = gStyle->GetOptStat();
int saveFit = gStyle->GetOptFit();
- this->SetStats(0); // remove the old stat window
- this->SetStats(1); // set the flag to create thes stat window during Draw()
- CbmQaHist* tmp = (CbmQaHist*) this->Clone("myClone");
- GetCanvas().cd();
+ GetDummyCanvas().cd();
gStyle->SetOptStat(fOptStat);
gStyle->SetOptFit(fOptFit);
- // the only way to create and auto-size the stat window
- tmp->Draw();
- GetCanvas().Update();
- // move the stat window to *this
- TObject* obj = tmp->FindObject("stats");
- // the list of functions seems always to exist, but let's check it
- if (obj && this->GetListOfFunctions()) {
- tmp->GetListOfFunctions()->Remove(obj);
- this->GetListOfFunctions()->Add(obj);
- }
- GetCanvas().Clear();
+
+ this->SetStats(0); // remove the old stat window
+ this->SetStats(1); // set the flag to create the stat window during Draw()
+ this->Draw();
+ GetDummyCanvas().Update();
+ GetDummyCanvas().Clear();
+
+ // restore the environment
gStyle->SetOptStat(saveStat);
gStyle->SetOptFit(saveFit);
if (savePad) savePad->cd();
@@ -112,12 +134,12 @@ public:
private:
/// a static canvas for temporary drawing
- static CbmQaCanvas& GetCanvas()
+ static CbmQaCanvas& GetDummyCanvas()
{
/// the static variable will be initialised at the first call;
/// deleted at the application end (c++11)
- static CbmQaCanvas tmp("CbmQaTempCanvas", "CbmQaTempCanvas", 1, 1);
- return tmp;
+ static CbmQaCanvas dummy("CbmQaTempCanvas", "CbmQaTempCanvas", 1, 1);
+ return dummy;
}
int fOptStat = 1;
@@ -126,8 +148,4 @@ private:
ClassDefNV(CbmQaHist, 1);
};
-// shortcuts
-typedef CbmQaHist<TH1F> CbmQaH1F;
-typedef CbmQaHist<TH1D> CbmQaH1D;
-
#endif
diff --git a/macro/run/run_qa.C b/macro/run/run_qa.C
index b56a8a2007..1754454a3c 100644
--- a/macro/run/run_qa.C
+++ b/macro/run/run_qa.C
@@ -48,7 +48,9 @@ void run_qa(TString dataTra = "data/sis100_muon_jpsi_test", TString dataRaw = "d
// ----- Logger settings ----------------------------------------------
FairLogger::GetLogger()->SetLogScreenLevel("INFO");
- FairLogger::GetLogger()->SetLogVerbosityLevel("LOW");
+ fair::Logger::DefineVerbosity(
+ "user1", fair::VerbositySpec::Make(fair::VerbositySpec::Info::severity, fair::VerbositySpec::Info::file_line));
+ FairLogger::GetLogger()->SetLogVerbosityLevel("user1");
// ------------------------------------------------------------------------
// ----- Environment --------------------------------------------------
@@ -173,6 +175,7 @@ void run_qa(TString dataTra = "data/sis100_muon_jpsi_test", TString dataRaw = "d
//run->AddTask(new CbmTrdDigitizerPRFQa()); //works put currently doesn't do anything
//run->AddTask(new CbmTrdHitRateFastQa()); //opens lots of windows
run->AddTask(new CbmTrdHitProducerQa());
+ run->AddTask(new CbmTrackerInputQaTrd()); // Tracker requirements to TRD
}
// ------------------------------------------------------------------------
diff --git a/reco/L1/CMakeLists.txt b/reco/L1/CMakeLists.txt
index 6c4fb45a51..5c76eaf2a1 100644
--- a/reco/L1/CMakeLists.txt
+++ b/reco/L1/CMakeLists.txt
@@ -41,6 +41,7 @@ ${CBMROOT_SOURCE_DIR}/reco/L1
${CBMROOT_SOURCE_DIR}/reco/L1/L1Algo
${CBMROOT_SOURCE_DIR}/reco/L1/OffLineInterface
${CBMROOT_SOURCE_DIR}/reco/L1/ParticleFinder
+${CBMROOT_SOURCE_DIR}/reco/L1/qa
${CBMROOT_SOURCE_DIR}/reco/KF
${CBMROOT_SOURCE_DIR}/reco/KF/KFQA
@@ -52,7 +53,7 @@ ${CBMROOT_SOURCE_DIR}/sim/transport/geosetup
${CBMDATA_DIR}
${CBMDATA_DIR}/base
- ${CBMDETECTORBASE_DIR}/sts
+${CBMDETECTORBASE_DIR}/sts
${CBMROOT_SOURCE_DIR}/mvd
${CBMDETECTORBASE_DIR}/trd
@@ -144,6 +145,8 @@ L1Algo/utils/L1AlgoPulls.cxx
ParticleFinder/CbmL1PFFitter.cxx
ParticleFinder/CbmL1PFMCParticle.cxx
+
+qa/CbmTrackerInputQaTrd.cxx
)
set(HEADERS
@@ -175,6 +178,8 @@ OffLineInterface/CbmL1GlobalFindTracksEvents.h
#OffLineInterface / CbmL1SttTrack.h
L1Algo/L1Def.h
L1Algo/L1Vector.h
+
+qa/CbmTrackerInputQaTrd.h
)
@@ -231,6 +236,7 @@ Set(DEPENDENCIES
CbmStsBase
CbmRecoBase
CbmRecoSts
+ CbmQaBase
boost_regex
)
diff --git a/reco/L1/L1LinkDef.h b/reco/L1/L1LinkDef.h
index 4ddd908135..d7b233e93c 100644
--- a/reco/L1/L1LinkDef.h
+++ b/reco/L1/L1LinkDef.h
@@ -29,5 +29,7 @@
//#pragma link C++ class CbmL1SttHit+;
//#pragma link C++ class CbmL1SttTrackFinder+;
//#pragma link C++ class CbmL1SttTrack+;
+#pragma link C++ class CbmTrackerInputQaTrd + ;
+
#endif
diff --git a/reco/L1/qa/CbmTrackerInputQaTrd.cxx b/reco/L1/qa/CbmTrackerInputQaTrd.cxx
new file mode 100644
index 0000000000..5ebe1fae51
--- /dev/null
+++ b/reco/L1/qa/CbmTrackerInputQaTrd.cxx
@@ -0,0 +1,739 @@
+/* Copyright (C) 2021 Facility for Antiproton and Ion Research in Europe, Darmstadt
+ SPDX-License-Identifier: GPL-3.0-only
+ Authors: Sergey Gorbunov[committer]*/
+
+/// @file CbmTrackerInputQaTrd.cxx
+/// @author Sergey Gorbunov
+/// @date 02.11.2021
+
+#include "CbmTrackerInputQaTrd.h"
+
+#include "CbmDefs.h"
+#include "CbmDigiManager.h"
+#include "CbmLink.h"
+#include "CbmMCDataArray.h"
+#include "CbmMCDataManager.h"
+#include "CbmMCEventList.h"
+#include "CbmMCTrack.h"
+#include "CbmMatch.h"
+#include "CbmQaCanvas.h"
+//#include "CbmSetup.h"
+#include "CbmTimeSlice.h"
+#include "CbmTrdCluster.h"
+#include "CbmTrdDigi.h"
+#include "CbmTrdHit.h"
+#include "CbmTrdParModDigi.h" // for CbmTrdModule
+#include "CbmTrdParModGeo.h"
+#include "CbmTrdParSetDigi.h" // for CbmTrdParSetDigi
+#include "CbmTrdParSetGeo.h"
+#include "CbmTrdPoint.h"
+
+#include <FairRootManager.h>
+#include <FairRunAna.h>
+#include <FairRuntimeDb.h>
+#include <FairSink.h>
+#include <FairTask.h>
+#include <Logger.h>
+
+#include <TClonesArray.h>
+#include <TDatabasePDG.h>
+#include <TGenericClassInfo.h>
+#include <TGeoManager.h>
+#include <TGeoNode.h>
+#include <TMathBase.h>
+#include <TObjArray.h>
+#include <TParameter.h>
+#include <TParticlePDG.h>
+#include <TString.h>
+#include <TStyle.h>
+
+#include <algorithm>
+#include <iostream>
+#include <map>
+#include <vector>
+
+#include <stdio.h>
+
+ClassImp(CbmTrackerInputQaTrd);
+
+// -------------------------------------------------------------------------
+
+CbmTrackerInputQaTrd::CbmTrackerInputQaTrd(const char* name, Int_t verbose) : FairTask(name, verbose)
+{
+ // Constructor
+
+ // Create a list of histogramms
+
+ fHistList.clear();
+ fHistList.reserve(20);
+ fHistList.push_back(&fh1DresidualU);
+ fHistList.push_back(&fh1DresidualV);
+ fHistList.push_back(&fh1DresidualT);
+ fHistList.push_back(&fh2DresidualX);
+ fHistList.push_back(&fh2DresidualY);
+ fHistList.push_back(&fh2DresidualT);
+ fHistList.push_back(&fh1DpullU);
+ fHistList.push_back(&fh1DpullV);
+ fHistList.push_back(&fh1DpullT);
+ fHistList.push_back(&fh2DpullX);
+ fHistList.push_back(&fh2DpullY);
+ fHistList.push_back(&fh2DpullT);
+
+ // Keep the ownership on the histograms to avoid double deletion
+ for (unsigned int i = 0; i < fHistList.size(); i++) {
+ fHistList[i]->SetDirectory(0);
+ }
+}
+
+// -------------------------------------------------------------------------
+CbmTrackerInputQaTrd::~CbmTrackerInputQaTrd()
+{ /// Destructor
+ DeInit();
+}
+
+
+// -------------------------------------------------------------------------
+void CbmTrackerInputQaTrd::DeInit()
+{
+ fIsTrdInSetup = 0;
+ fIsMcPresent = false;
+ fNtrackingStations = 0;
+
+ fTimeSlice = nullptr;
+ fDigiManager = nullptr;
+
+ fMcManager = nullptr;
+ fMcTracks = nullptr;
+ fMcPoints = nullptr;
+
+ fClusters = nullptr;
+ fHits = nullptr;
+ fHitMatches = nullptr;
+
+ fOutFolder.Clear();
+
+ fHistFolder = nullptr;
+
+ fNevents.SetVal(0);
+
+ fhPointsPerHit.clear();
+ fhHitsPerPoint.clear();
+}
+// -------------------------------------------------------------------------
+
+void CbmTrackerInputQaTrd::SetParContainers()
+{
+ fTrdDigiPar = nullptr;
+ fTrdGeoPar = nullptr;
+
+ // Get run and runtime database
+ FairRunAna* ana = FairRunAna::Instance();
+ if (!ana) { LOG(fatal) << "No FairRunAna instance"; }
+
+ FairRuntimeDb* rtdb = ana->GetRuntimeDb();
+ if (!rtdb) { LOG(fatal) << "No FairRuntimeDb instance"; }
+
+ fTrdDigiPar = dynamic_cast<CbmTrdParSetDigi*>(rtdb->getContainer("CbmTrdParSetDigi"));
+ fTrdGeoPar = dynamic_cast<CbmTrdParSetGeo*>(rtdb->getContainer("CbmTrdParSetGeo"));
+}
+
+// -------------------------------------------------------------------------
+InitStatus CbmTrackerInputQaTrd::ReInit()
+{
+ // Initialize and check the setup
+
+ DeInit();
+
+ // Check if TRD is present in the CBM setup
+ // A straightforward way to do so is CbmSetup::Instance()->IsActive(ECbmModuleId::kTrd),
+ // but unfortunatelly CbmSetup class is unaccesible.
+ // ( CbmSimSteer library requires libfreetype that has linking problems on MacOsX )
+ // Therefore let's simply check if any TRD material is present in the geometry.
+ //
+ fIsTrdInSetup = 0;
+ {
+ TObjArray* topNodes = gGeoManager->GetTopNode()->GetNodes();
+ for (Int_t iTopNode = 0; iTopNode < topNodes->GetEntriesFast(); iTopNode++) {
+ TGeoNode* topNode = static_cast<TGeoNode*>(topNodes->At(iTopNode));
+ if (TString(topNode->GetName()).Contains("trd")) {
+ fIsTrdInSetup = 1;
+ break;
+ }
+ }
+ }
+
+ if (!fIsTrdInSetup) {
+ LOG(info) << "TRD is not in the setup, do nothing";
+ return kSUCCESS;
+ }
+
+ // check parameter containers
+
+ if (!fTrdDigiPar) {
+ LOG(error) << "No CbmTrdParSetDigi container in FairRuntimeDb";
+ return kERROR;
+ }
+
+ if (!fTrdGeoPar) {
+ LOG(error) << "No CbmTrdParSetGeo container in FairRuntimeDb";
+ return kERROR;
+ }
+
+ FairRootManager* manager = FairRootManager::Instance();
+ fDigiManager = CbmDigiManager::Instance();
+ fDigiManager->Init();
+
+ fTimeSlice = dynamic_cast<CbmTimeSlice*>(manager->GetObject("TimeSlice."));
+ if (!fTimeSlice) {
+ LOG(error) << "No time slice found";
+ return kERROR;
+ }
+
+ fHits = dynamic_cast<TClonesArray*>(manager->GetObject("TrdHit"));
+
+ if (!fHits) {
+ LOG(error) << "No TRD hit array found";
+ return kERROR;
+ }
+
+ fClusters = dynamic_cast<TClonesArray*>(manager->GetObject("TrdCluster"));
+
+ if (!fClusters) {
+ LOG(error) << "No TRD cluster array found";
+ return kERROR;
+ }
+
+ fMcManager = dynamic_cast<CbmMCDataManager*>(manager->GetObject("MCDataManager"));
+
+ fIsMcPresent = (fMcManager != nullptr);
+
+ if (fIsMcPresent) {
+
+ fMcEventList = dynamic_cast<CbmMCEventList*>(manager->GetObject("MCEventList."));
+ fMcTracks = fMcManager->InitBranch("MCTrack");
+ fMcPoints = fMcManager->InitBranch("TrdPoint");
+ fHitMatches = dynamic_cast<TClonesArray*>(manager->GetObject("TrdHitMatch"));
+
+ // we assume that when TRD is in the setup and the MC manager is present,
+ // the TRD MC data must be present too
+
+ if (!fMcEventList) {
+ LOG(error) << ": No MCEventList data!";
+ return kERROR;
+ }
+
+ if (!fMcTracks) {
+ LOG(error) << "No MC tracks found";
+ return kERROR;
+ }
+
+ if (!fMcPoints) {
+ LOG(error) << "No MC points found";
+ return kERROR;
+ }
+
+ if (!fHitMatches) {
+ LOG(error) << "No TRD hit matches found";
+ return kERROR;
+ }
+
+ if (!fDigiManager->IsMatchPresent(ECbmModuleId::kTrd)) {
+ LOG(error) << "No TRD digi matches found";
+ return kERROR;
+ }
+ }
+
+ // count TRD stations
+ // TODO: put the code below to some TRD geometry interface
+
+ fNtrackingStations = 0;
+ {
+ Int_t layerCounter = 0;
+ TObjArray* topNodes = gGeoManager->GetTopNode()->GetNodes();
+ for (Int_t iTopNode = 0; iTopNode < topNodes->GetEntriesFast(); iTopNode++) {
+ TGeoNode* topNode = static_cast<TGeoNode*>(topNodes->At(iTopNode));
+ if (TString(topNode->GetName()).Contains("trd")) {
+ TObjArray* layers = topNode->GetNodes();
+ for (Int_t iLayer = 0; iLayer < layers->GetEntriesFast(); iLayer++) {
+ TGeoNode* layer = static_cast<TGeoNode*>(layers->At(iLayer));
+ if (TString(layer->GetName()).Contains("layer")) { layerCounter++; }
+ }
+ }
+ }
+ fNtrackingStations = layerCounter;
+ }
+
+ LOG(debug) << "Init(): fNtrackingStations = " << fNtrackingStations;
+
+ if (fNtrackingStations <= 0) {
+ LOG(error) << "can not count TRD tracking stations";
+ return kERROR;
+ }
+
+ // initialise histograms
+ fOutFolder.SetOwner(false);
+ fHistFolder = fOutFolder.AddFolder("rawHist", "Raw Histograms");
+ gStyle->SetOptStat(0);
+
+ fNevents.SetVal(0);
+ fHistFolder->Add(&fNevents);
+
+ for (unsigned int i = 0; i < fHistList.size(); i++) {
+ fHistList[i]->Reset();
+ fHistFolder->Add(fHistList[i]);
+ }
+
+ fhPointsPerHit.clear();
+ fhHitsPerPoint.clear();
+ fhEfficiencyR.clear();
+ fhEfficiencyXY.clear();
+
+ fhPointsPerHit.reserve(fNtrackingStations);
+ fhHitsPerPoint.reserve(fNtrackingStations);
+ fhEfficiencyR.reserve(fNtrackingStations);
+ fhEfficiencyXY.reserve(fNtrackingStations);
+
+ for (Int_t i = 0; i < fNtrackingStations; i++) {
+ fhPointsPerHit.emplace_back(Form("hMcPointsPerHit%i", i),
+ Form("MC Points per Hit: Station %i;N mc points / hit", i), 10, -0.5, 9.5);
+ fhPointsPerHit[i].SetDirectory(0);
+ fhPointsPerHit[i].SetLineWidth(2);
+ fhPointsPerHit[i].SetOptStat(110);
+ fHistFolder->Add(&fhPointsPerHit[i]);
+ }
+
+ for (Int_t i = 0; i < fNtrackingStations; i++) {
+ fhHitsPerPoint.emplace_back(Form("hHitsPerMcPoint%i", i),
+ Form("Hits per MC Point: Station %i; N hits / mc point", i), 10, -0.5, 9.5);
+ fhHitsPerPoint[i].SetDirectory(0);
+ fhHitsPerPoint[i].SetLineWidth(2);
+ fhHitsPerPoint[i].SetOptStat(110);
+ fHistFolder->Add(&fhHitsPerPoint[i]);
+ }
+
+ for (Int_t i = 0; i < fNtrackingStations; i++) {
+
+ double dx = 350;
+ double dy = 350;
+ double dr = sqrt(dx * dx + dy * dy);
+
+ fhEfficiencyR.emplace_back(Form("hEfficiencyR%i", i), Form("Efficiency R: Station %i;R [cm]", i), 100, 0, dr);
+ fhEfficiencyR[i].SetDirectory(0);
+ fhEfficiencyR[i].SetOptStat(1110);
+ fHistFolder->Add(&fhEfficiencyR[i]);
+
+ fhEfficiencyXY.emplace_back(Form("hEfficiencyXY%i", i), Form("Efficiency XY: Station %i;X [cm];Y [cm]", i), 50, -dx,
+ dx, 50, -dy, dy);
+ fhEfficiencyXY[i].SetDirectory(0);
+ fhEfficiencyXY[i].SetOptStat(10);
+ fhEfficiencyXY[i].GetYaxis()->SetTitleOffset(1.4);
+ fHistFolder->Add(&fhEfficiencyXY[i]);
+ }
+
+ fCanvResidual.Clear();
+ fCanvResidual.Divide2D(6);
+
+ fCanvPull.Clear();
+ fCanvPull.Divide2D(6);
+
+ fCanvEfficiencyR.Clear();
+ fCanvEfficiencyR.Divide2D(fNtrackingStations);
+
+ fCanvEfficiencyXY.Clear();
+ fCanvEfficiencyXY.Divide2D(fNtrackingStations);
+
+ fCanvPointsPerHit.Clear();
+ fCanvPointsPerHit.Divide2D(fNtrackingStations);
+
+ fCanvHitsPerPoint.Clear();
+ fCanvHitsPerPoint.Divide2D(fNtrackingStations);
+
+ fOutFolder.Add(&fCanvResidual);
+ fOutFolder.Add(&fCanvPull);
+ fOutFolder.Add(&fCanvEfficiencyR);
+ fOutFolder.Add(&fCanvEfficiencyXY);
+ fOutFolder.Add(&fCanvPointsPerHit);
+ fOutFolder.Add(&fCanvHitsPerPoint);
+
+ // analyse the geometry setup
+ return GeometryQa();
+}
+
+
+// -------------------------------------------------------------------------
+void CbmTrackerInputQaTrd::Exec(Option_t*)
+{
+ if (!fIsTrdInSetup) { return; }
+
+ // update number of processed events
+ fNevents.SetVal(fNevents.GetVal() + 1);
+ LOG(debug) << "Event: " << fNevents.GetVal();
+ ResolutionQa();
+}
+
+
+// -------------------------------------------------------------------------
+void CbmTrackerInputQaTrd::Finish()
+{
+ FairSink* sink = FairRootManager::Instance()->GetSink();
+ if (sink) { sink->WriteObject(&GetQa(), nullptr); }
+}
+
+
+// -------------------------------------------------------------------------
+InitStatus CbmTrackerInputQaTrd::GeometryQa()
+{
+ // check geometry of the TRD modules
+
+ double lastZ = -1;
+ for (int iStation = 0; iStation < fNtrackingStations; iStation++) {
+
+ int ModuleId = fTrdDigiPar->GetModuleId(iStation);
+
+ const CbmTrdParModGeo* pGeo = dynamic_cast<const CbmTrdParModGeo*>(fTrdGeoPar->GetModulePar(ModuleId));
+ if (!pGeo) {
+ LOG(fatal) << " No Geo params for station " << iStation << " module " << ModuleId << " found ";
+ return kERROR;
+ }
+
+ double staZ = pGeo->GetZ();
+
+ // check that the stations are properly ordered in Z
+ if ((iStation > 0) && (staZ <= lastZ)) {
+ LOG(error) << "TRD trackig stations are not properly ordered in Z";
+ return kERROR;
+ }
+ lastZ = staZ;
+
+ //TODO: what are these 3 and 6 here in L1 code? Why are they hardcoed?
+ // if (num == 0 || num == 2 || num == 4) geo.push_back(3);
+ // if (num == 1 || num == 3) geo.push_back(6);
+ }
+
+ return kSUCCESS;
+}
+
+// -------------------------------------------------------------------------
+void CbmTrackerInputQaTrd::ResolutionQa()
+{
+ if (!fDigiManager->IsMatchPresent(ECbmModuleId::kTrd)) return;
+
+ if (!fIsMcPresent) return;
+
+ Int_t nHits = fHits->GetEntriesFast();
+ Int_t nClusters = fClusters->GetEntriesFast();
+ Int_t nDigis = fDigiManager->GetNofDigis(ECbmModuleId::kTrd);
+
+ int nMcEvents = fMcEventList->GetNofEvents();
+
+ // Vector of integers parallel to mc points
+
+ std::vector<std::vector<int>> pointNhits;
+ pointNhits.resize(nMcEvents);
+ for (int iE = 0; iE < nMcEvents; iE++) {
+ int fileId = fMcEventList->GetFileIdByIndex(iE);
+ int eventId = fMcEventList->GetEventIdByIndex(iE);
+ int nPoints = fMcPoints->Size(fileId, eventId);
+ pointNhits[iE].resize(nPoints, 0);
+ }
+
+ for (Int_t iHit = 0; iHit < nHits; iHit++) {
+
+ CbmTrdHit* hit = dynamic_cast<CbmTrdHit*>(fHits->At(iHit));
+ if (!hit) {
+ LOG(error) << "TRD hit N " << iHit << " doesn't exist";
+ return;
+ }
+
+ if ((int) hit->GetClassType() > 1) {
+ // class type 0: trd-1D
+ // class type 1: trd-2D
+ // the return type is (currently) boolean, so this error should never happen
+ LOG(error) << "Unknown detector class type " << hit->GetClassType() << " for TRD hit N " << iHit;
+ return;
+ }
+
+ const int address = hit->GetAddress();
+ int StationIndex = CbmTrdAddress::GetLayerId(address);
+
+ if (StationIndex < 0 || StationIndex >= fNtrackingStations) {
+ LOG(fatal) << "TRD hit layer id " << StationIndex << " is out of range";
+ return;
+ }
+
+ Int_t clusterId = hit->GetRefId();
+ if (clusterId < 0 || clusterId >= nClusters) {
+ LOG(error) << "TRD cluster id " << clusterId << " is out of range";
+ return;
+ }
+
+ CbmTrdCluster* cluster = dynamic_cast<CbmTrdCluster*>(fClusters->At(clusterId));
+ if (!cluster) {
+ LOG(error) << "TRD cluster N " << clusterId << " doesn't exist";
+ return;
+ }
+
+ if (cluster->GetAddress() != address) {
+ LOG(error) << "TRD hit address " << address << " differs from its cluster address " << cluster->GetAddress();
+ return;
+ }
+
+ Int_t nClusterDigis = cluster->GetNofDigis();
+
+ if (nClusterDigis <= 0) {
+ LOG(error) << "hit match for TRD cluster " << clusterId << " has no digis ";
+ return;
+ }
+
+ // custom finder of the digi matches
+
+ CbmMatch myHitMatch;
+ for (Int_t iDigi = 0; iDigi < nClusterDigis; iDigi++) {
+ Int_t digiIdx = cluster->GetDigi(iDigi);
+ if (digiIdx < 0 || digiIdx >= nDigis) {
+ LOG(error) << "TRD cluster: digi index " << digiIdx << " out of range ";
+ return;
+ }
+ const CbmTrdDigi* digi = CbmDigiManager::Instance()->Get<CbmTrdDigi>(digiIdx);
+ if (!digi) {
+ LOG(fatal) << "TRD digi " << iDigi << " not found";
+ return;
+ }
+
+ if (digi->GetAddressModule() != address) {
+ LOG(error) << "TRD hit address " << address << " differs from its digi address " << digi->GetAddressModule();
+ return;
+ }
+ const CbmMatch* match = dynamic_cast<const CbmMatch*>(fDigiManager->GetMatch(ECbmModuleId::kTrd, digiIdx));
+ if (!match) {
+ LOG(fatal) << "TRD digi match " << digiIdx << " not found";
+ return;
+ }
+ myHitMatch.AddLinks(*match);
+ }
+ if (myHitMatch.GetNofLinks() == 0) { continue; }
+
+ const CbmLink& bestLink = myHitMatch.GetMatchedLink();
+
+ { // check if the hit match is correct
+ CbmMatch* hitMatch = dynamic_cast<CbmMatch*>(fHitMatches->At(iHit));
+ if (hitMatch == nullptr) { LOG(error) << "hit match for TRD hit " << iHit << " doesn't exist"; }
+ else {
+ const CbmLink& link = hitMatch->GetMatchedLink();
+ if ((link != bestLink) && (link.GetWeight() != bestLink.GetWeight())) {
+ LOG(error) << "hit match for TRD hit " << iHit << " doesn't correspond to digi matches";
+ }
+ }
+ }
+
+ // how many MC points? fill N hits per mc point
+
+ int nHitPoints = 0;
+ for (int i = 0; i < myHitMatch.GetNofLinks(); i++) {
+ const CbmLink& link = myHitMatch.GetLink(i);
+ if (link.GetIndex() >= 0) { // not a noise digi
+ nHitPoints++;
+ int iE = fMcEventList->GetEventIndex(link);
+ if (iE < 0 || iE >= nMcEvents) {
+ LOG(error) << "link points to a non-existing MC event";
+ return;
+ }
+ if (link.GetIndex() >= (int) pointNhits[iE].size()) {
+ LOG(error) << "link points to a non-existing MC point";
+ return;
+ }
+ pointNhits[iE][link.GetIndex()]++;
+ }
+ }
+
+ fhPointsPerHit[StationIndex].Fill(nHitPoints);
+
+ // take corresponding MC point
+
+ // skip hits from the noise digis
+ if (bestLink.GetIndex() < 0) { continue; }
+
+ CbmTrdPoint* p = dynamic_cast<CbmTrdPoint*>(fMcPoints->Get(bestLink));
+ if (p == nullptr) {
+ LOG(error) << "link points to a non-existing MC point";
+ return;
+ }
+
+ if (p->GetModuleAddress() != (int) CbmTrdAddress::GetModuleAddress(address)) {
+ LOG(error) << "mc point module address differs from the hit module address";
+ return;
+ }
+
+ // check that the nominal station Z is not far from the active material
+ {
+ // the cut of 1 cm is choosen arbitrary and can be changed
+
+ int ModuleId = fTrdDigiPar->GetModuleId(StationIndex);
+
+ const CbmTrdParModGeo* pGeo = dynamic_cast<const CbmTrdParModGeo*>(fTrdGeoPar->GetModulePar(ModuleId));
+ if (!pGeo) {
+ LOG(fatal) << " No Geo params for station " << StationIndex << " module " << ModuleId << " found ";
+ return;
+ }
+
+ double staZ = pGeo->GetZ(); // module->GetZ(); //+ 410;
+ if ((staZ < p->GetZIn() - 1.) || (staZ > p->GetZOut() + 1.)) {
+ LOG(error) << "TRD station " << StationIndex << ": active material Z[" << p->GetZIn() << " cm," << p->GetZOut()
+ << " cm] is too far from the nominal station Z " << staZ << " cm";
+ return;
+ }
+ // the cut of 1 cm is choosen arbitrary and can be changed
+ if (fabs(hit->GetZ() - staZ) > 1.) {
+ LOG(error) << "TRD station " << StationIndex << ": hit Z " << hit->GetZ()
+ << " cm, is too far from the nominal station Z " << staZ << " cm";
+ return;
+ }
+ }
+
+ // residual and pull
+
+ if (nHitPoints != 1) continue; // only check residual for non-mixed hits
+
+ double t0 = fMcEventList->GetEventTime(bestLink);
+ if (t0 < 0) {
+ LOG(error) << "MC event time doesn't exist for a TRD link";
+ return;
+ }
+
+ double x = p->GetX(); // take the "In" point since the time is stored only for this point
+ double y = p->GetY();
+ double z = p->GetZ();
+ double t = t0 + p->GetTime();
+ double px = p->GetPx();
+ double py = p->GetPy();
+ double pz = p->GetPz();
+
+ // skip very slow particles
+ if (fabs(p->GetPzOut()) < 0.01) continue;
+
+ // transport the particle from the MC point to the hit Z
+ double dz = hit->GetZ() - z;
+ x += dz * px / pz;
+ y += dz * py / pz;
+
+ // get particle mass
+ double mass = 0;
+ {
+ CbmLink mcTrackLink = bestLink;
+ mcTrackLink.SetIndex(p->GetTrackID());
+ CbmMCTrack* mcTrack = dynamic_cast<CbmMCTrack*>(fMcTracks->Get(mcTrackLink));
+ if (!mcTrack) {
+ LOG(error) << "MC track " << p->GetTrackID() << " doesn't exist";
+ return;
+ }
+ Int_t pdg = mcTrack->GetPdgCode();
+ if (pdg < 9999999 && ((TParticlePDG*) TDatabasePDG::Instance()->GetParticle(pdg))) {
+ mass = TDatabasePDG::Instance()->GetParticle(pdg)->Mass();
+ }
+ }
+
+ constexpr double speedOfLight = 29.979246; // cm/ns
+ TVector3 mom3;
+ p->Momentum(mom3);
+ t += dz / (pz * speedOfLight) * sqrt(mass * mass + mom3.Mag2());
+
+ double du = hit->GetX() - x;
+ double dv = hit->GetY() - y;
+ double dt = hit->GetTime() - t;
+ double su = hit->GetDx();
+ double sv = hit->GetDy();
+ double st = hit->GetTimeError();
+
+ // residuals
+ if (hit->GetClassType() == 0) {
+ if (StationIndex % 2) {
+ std::swap(du, dv);
+ std::swap(su, sv);
+ }
+ fh1DresidualU.Fill(du);
+ fh1DresidualV.Fill(dv);
+ fh1DresidualT.Fill(hit->GetTime() - t);
+ }
+ else {
+ fh2DresidualX.Fill(du);
+ fh2DresidualY.Fill(dv);
+ fh2DresidualT.Fill(hit->GetTime() - t);
+ }
+
+ // pulls
+ if ((su < 1.e-5) || (sv < 1.e-5) || (st < 1.e-5)) {
+ LOG(error) << "TRD hit errors are not properly set: errX " << hit->GetDx() << " errY " << hit->GetDy() << " errT "
+ << st;
+ return;
+ }
+
+ if (hit->GetClassType() == 0) {
+ fh1DpullU.Fill(du / su);
+ fh1DpullV.Fill(dv / sv);
+ fh1DpullT.Fill(dt / st);
+ }
+ else {
+ fh2DpullX.Fill(du / su);
+ fh2DpullY.Fill(dv / sv);
+ fh2DpullT.Fill(dt / st);
+ }
+
+ } // iHit
+
+ // fill efficiency: N hits per MC point
+
+ for (int iE = 0; iE < nMcEvents; iE++) {
+ int fileId = fMcEventList->GetFileIdByIndex(iE);
+ int eventId = fMcEventList->GetEventIdByIndex(iE);
+ int nPoints = fMcPoints->Size(fileId, eventId);
+ for (int ip = 0; ip < nPoints; ip++) {
+ CbmTrdPoint* p = dynamic_cast<CbmTrdPoint*>(fMcPoints->Get(fileId, eventId, ip));
+ if (p == nullptr) {
+ LOG(error) << "MC point doesn't exist";
+ return;
+ }
+ auto address = p->GetModuleAddress();
+ int StationIndex = CbmTrdAddress::GetLayerId(address);
+ if (StationIndex < 0 || StationIndex >= fNtrackingStations) {
+ LOG(fatal) << "TRD hit layer id " << StationIndex << " for module " << address << " is out of range";
+ return;
+ }
+ fhHitsPerPoint[StationIndex].Fill(pointNhits[iE][ip]);
+ fhEfficiencyR[StationIndex].Fill(sqrt(p->GetX() * p->GetX() + p->GetY() * p->GetY()), (pointNhits[iE][ip] > 0));
+ fhEfficiencyXY[StationIndex].Fill(p->GetX(), p->GetY(), (pointNhits[iE][ip] > 0));
+ }
+ }
+}
+
+
+// -------------------------------------------------------------------------
+TFolder& CbmTrackerInputQaTrd::GetQa()
+{
+ gStyle->SetPaperSize(20, 20);
+
+ for (Int_t i = 0; i < fNtrackingStations; i++) {
+ fCanvHitsPerPoint.cd(i + 1);
+ fhHitsPerPoint[i].DrawCopy("", "");
+ fCanvPointsPerHit.cd(i + 1);
+ fhPointsPerHit[i].DrawCopy("", "");
+
+ fCanvEfficiencyR.cd(i + 1);
+ fhEfficiencyR[i].DrawCopy("colz", "");
+
+ fCanvEfficiencyXY.cd(i + 1);
+ fhEfficiencyXY[i].DrawCopy("colz", "");
+ }
+
+ for (UInt_t i = 0; i < fHistList.size(); i++) {
+ fHistList[i]->Fit("gaus", "Q"); // Q for the quiet mode
+ }
+
+ for (UInt_t i = 0; i < 6; i++) {
+ fCanvResidual.cd(i + 1);
+ fHistList[i]->DrawCopy("", "");
+ fCanvPull.cd(i + 1);
+ fHistList[6 + i]->DrawCopy("", "");
+ }
+
+ return fOutFolder;
+}
diff --git a/reco/L1/qa/CbmTrackerInputQaTrd.h b/reco/L1/qa/CbmTrackerInputQaTrd.h
new file mode 100644
index 0000000000..c7c6f5aad8
--- /dev/null
+++ b/reco/L1/qa/CbmTrackerInputQaTrd.h
@@ -0,0 +1,166 @@
+/* Copyright (C) 2021 Facility for Antiproton and Ion Research in Europe, Darmstadt
+ SPDX-License-Identifier: GPL-3.0-only
+ Authors: Sergey Gorbunov[committer]*/
+
+/// @file CbmTrackerInputQaTrd.h
+/// @author Sergey Gorbunov
+/// @date 02.11.2021
+
+#ifndef CbmTrackerInputQaTrd_H
+#define CbmTrackerInputQaTrd_H
+
+#include "CbmQaCanvas.h"
+#include "CbmQaHist.h"
+
+#include <FairTask.h>
+
+#include <TFolder.h>
+#include <TH1D.h>
+#include <TH1F.h>
+#include <TH1I.h>
+#include <TH2F.h>
+#include <TParameter.h>
+#include <TProfile.h>
+#include <TProfile2D.h>
+
+class CbmDigiManager;
+class CbmMCDataManager;
+class CbmMCEventList;
+class CbmMCDataArray;
+class CbmTimeSlice;
+class CbmTrdParSetGeo;
+class CbmTrdParSetDigi;
+
+class TClonesArray;
+
+///
+/// CbmTrackerInputQaTrd class represents all the CA tracker requirements for the TRD detector.
+/// When this QA test is passed, the tracker must work (at least that is the idea).
+///
+/// The class ensures that the tracker has the correct understanding of the TRD geometry and interfaces
+/// and checks the quality of all tracking-related aspects of the TRD data.
+///
+class CbmTrackerInputQaTrd : public FairTask {
+
+public:
+ /// Constructor
+ CbmTrackerInputQaTrd(const char* name = "TrackerInputQaTrd", Int_t verbose = 1);
+
+ /// Destructor
+ ~CbmTrackerInputQaTrd();
+
+ /// FairTask: Intialisation at begin of run.
+ InitStatus Init() { return ReInit(); }
+
+ /// FairTask: Reinitialisation.
+ InitStatus ReInit();
+
+ /// FairTask: Intialise parameter containers.
+ void SetParContainers();
+
+ /// FairTask: Action at end of run. For this task and all of the subtasks.
+ void Finish();
+
+ /// TTask: Clear all data structures created by a previous execution of a task.
+ void Clear(Option_t* /*option*/ = "") {}
+
+ /// TTask: Process a timeslice
+ void Exec(Option_t*);
+
+ /// Prepare the Qa output and return it as a reference to an internal TFolder.
+ /// The reference is non-const, because FairSink can not write const objects
+ TFolder& GetQa();
+
+private:
+ /// Check the geometry
+ InitStatus GeometryQa();
+
+ /// Analysis of hit uncertainty (pull) distributions
+ void ResolutionQa();
+
+ /// Free the memory etc.
+ void DeInit();
+
+
+ // Setup
+
+ Bool_t fIsTrdInSetup {false};
+ Bool_t fIsMcPresent {false};
+
+ Int_t fNtrackingStations {0};
+
+ CbmTimeSlice* fTimeSlice {nullptr};
+ CbmTrdParSetGeo* fTrdGeoPar {nullptr};
+ CbmTrdParSetDigi* fTrdDigiPar {nullptr};
+ CbmDigiManager* fDigiManager {nullptr};
+
+ /// MC data
+ CbmMCEventList* fMcEventList {nullptr}; // list of MC events connected to the current time slice
+ CbmMCDataManager* fMcManager {nullptr};
+ CbmMCDataArray* fMcTracks {nullptr};
+ CbmMCDataArray* fMcPoints {nullptr};
+
+ /// Data
+ TClonesArray* fClusters {nullptr};
+ TClonesArray* fHits {nullptr};
+ TClonesArray* fHitMatches {nullptr};
+
+ /// Output
+
+
+ TFolder fOutFolder {"TrackerInputQaTrd", "TrackerInputQaTrd"}; /// output folder with histos and canvases
+ TFolder* fHistFolder {nullptr}; /// subfolder for histograms
+
+ TParameter<int> fNevents {"nEvents", 0}; /// number of processed events
+
+ /// Histogram for Residual Distribution
+ CbmQaHist<TH1D> fh1DresidualU {"h1DresidualU", "Trd1D: Residual U;(U_{reco} - U_{MC})(cm)", 100, -.5, .5};
+ CbmQaHist<TH1D> fh1DresidualV {"h1DresidualV", "Trd1D: Residual V;(V_{reco} - V_{MC})(cm)", 100, -10, 10};
+ CbmQaHist<TH1D> fh1DresidualT {"h1DresidualT", "Trd1D: Residual T;(T_{reco} - T_{MC})(ns)", 100, -50, 50};
+
+ CbmQaHist<TH1D> fh2DresidualX {"h2DresidualX", "Trd2D: Residual X;(X_{reco} - X_{MC})(cm)", 100, -5, 5};
+ CbmQaHist<TH1D> fh2DresidualY {"h2DresidualY", "Trd2D: Residual Y;(Y_{reco} - Y_{MC})(cm)", 100, -5, 5};
+ CbmQaHist<TH1D> fh2DresidualT {"h2DresidualT", "Trd2D: Residual T;(T_{reco} - T_{MC})(ns)", 100, -1000, 1000};
+
+ /// Histogram for PULL Distribution
+ CbmQaHist<TH1D> fh1DpullU {"h1DpullU", "Trd1D: Pull U;(U_{reco} - U_{MC}) / #sigmaU_{reco}", 100, -5, 5};
+ CbmQaHist<TH1D> fh1DpullV {"h1DpullV", "Trd1D: Pull V;(V_{reco} - V_{MC}) / #sigmaV_{reco}", 100, -5, 5};
+ CbmQaHist<TH1D> fh1DpullT {"h1DpullT", "Trd1D: Pull T;(T_{reco} - T_{MC}) / #sigmaT_{reco}", 100, -5, 5};
+
+ CbmQaHist<TH1D> fh2DpullX {"h2DpullX", "Trd2D: Pull X;(X_{reco} - X_{MC}) / #sigmaX_{reco}", 100, -5, 5};
+ CbmQaHist<TH1D> fh2DpullY {"h2DpullY", "Trd2D: Pull Y;(Y_{reco} - Y_{MC}) / #sigmaY_{reco}", 100, -5, 5};
+ CbmQaHist<TH1D> fh2DpullT {"h2DpullT", "Trd2D: Pull T;(T_{reco} - T_{MC}) / #sigmaT_{reco}", 100, -5, 5};
+
+ /// List of the above histogramms
+ std::vector<CbmQaHist<TH1D>*> fHistList;
+
+ /// hits purity
+ std::vector<CbmQaHist<TH1I>> fhPointsPerHit;
+
+ /// hits efficiency
+ std::vector<CbmQaHist<TH1I>> fhHitsPerPoint;
+
+ /// hits efficiency
+ std::vector<CbmQaHist<TProfile2D>> fhEfficiencyXY;
+ std::vector<CbmQaHist<TProfile>> fhEfficiencyR;
+
+ /// Canvaces: collection of histogramms
+
+ CbmQaCanvas fCanvResidual {"cResidual", "Residual Distribution", 3 * 500, 2 * 500};
+ CbmQaCanvas fCanvPull {"cPull", "Pull Distribution", 3 * 500, 2 * 500};
+ CbmQaCanvas fCanvEfficiencyXY {"cEfficiencyXY", "Efficiency XY: % reconstructed McPoint", 2 * 500, 2 * 500};
+ CbmQaCanvas fCanvEfficiencyR {"cEfficiencyR", "Efficiency R: % reconstructed McPoint", 2 * 500, 2 * 500};
+ CbmQaCanvas fCanvHitsPerPoint {"cHitsPerMcPoint", "Efficiency: Hits Per McPoint", 2 * 500, 2 * 500};
+ CbmQaCanvas fCanvPointsPerHit {"cMcPointsPerHit", "Purity: McPoints per Hit", 2 * 500, 2 * 500};
+
+private:
+ /// Suppressed copy constructor
+ CbmTrackerInputQaTrd(const CbmTrackerInputQaTrd&) = delete;
+
+ /// Suppressed assignment operator
+ CbmTrackerInputQaTrd& operator=(const CbmTrackerInputQaTrd&) = delete;
+
+ ClassDef(CbmTrackerInputQaTrd, 0);
+};
+
+#endif
--
GitLab