From 8c13846572fa6f4dc8ce754f7f483791429b72fe Mon Sep 17 00:00:00 2001
From: Dominik Smith <d.smith@gsi.de>
Date: Tue, 21 Sep 2021 16:27:32 +0200
Subject: [PATCH] Cleaned up CbmTrdHitProducerQa and added it to run_qa.C.
---
macro/run/run_qa.C | 2 +-
reco/detectors/trd/qa/CbmTrdHitProducerQa.cxx | 173 +++++-------------
reco/detectors/trd/qa/CbmTrdHitProducerQa.h | 26 ++-
3 files changed, 59 insertions(+), 142 deletions(-)
diff --git a/macro/run/run_qa.C b/macro/run/run_qa.C
index 9104e8c149..b56a8a2007 100644
--- a/macro/run/run_qa.C
+++ b/macro/run/run_qa.C
@@ -172,7 +172,7 @@ void run_qa(TString dataTra = "data/sis100_muon_jpsi_test", TString dataRaw = "d
//run->AddTask(new CbmTrdHitRateQa()); //opens lots of windows
//run->AddTask(new CbmTrdDigitizerPRFQa()); //works put currently doesn't do anything
//run->AddTask(new CbmTrdHitRateFastQa()); //opens lots of windows
- run->AddTask(new CbmTrdHitProducerQa()); //Histograms currently don't appear in output file
+ run->AddTask(new CbmTrdHitProducerQa());
}
// ------------------------------------------------------------------------
diff --git a/reco/detectors/trd/qa/CbmTrdHitProducerQa.cxx b/reco/detectors/trd/qa/CbmTrdHitProducerQa.cxx
index e1fe90ce58..1234f2cfde 100644
--- a/reco/detectors/trd/qa/CbmTrdHitProducerQa.cxx
+++ b/reco/detectors/trd/qa/CbmTrdHitProducerQa.cxx
@@ -9,6 +9,7 @@
#include "CbmTrdHitProducerQa.h"
+#include "CbmDigiManager.h"
#include "CbmMCTrack.h"
#include "CbmMatch.h"
#include "CbmTrdDigi.h"
@@ -39,9 +40,8 @@ CbmTrdHitProducerQa::CbmTrdHitProducerQa() : CbmTrdHitProducerQa("TrdHitProducer
// ---- Standard constructor ------------------------------------------------
CbmTrdHitProducerQa::CbmTrdHitProducerQa(const char* name, const char*)
: FairTask(name)
+ , fOutFolder("TrdHitProducerQA", "TrdHitProducerQA")
, fTrdHitCollection(NULL)
- , fTrdDigiCollection(NULL)
- , fTrdDigiMatchCollection(NULL)
, fTrdPointCollection(NULL)
, fMCTrackArray(NULL)
, fNoTrdStations(-1)
@@ -73,6 +73,20 @@ CbmTrdHitProducerQa::~CbmTrdHitProducerQa() {}
// ---- Initialisation ------------------------------------------------------
InitStatus CbmTrdHitProducerQa::Init()
{
+ fOutFolder.Clear();
+ histFolder = fOutFolder.AddFolder("hist", "Histogramms");
+
+ histFolder->Add(fHitPoolsX);
+ histFolder->Add(fHitPoolsY);
+ histFolder->Add(S1L1eTR15);
+ histFolder->Add(S1L1edEdx15);
+ histFolder->Add(S1L1edE15);
+ histFolder->Add(S1L1edEall);
+ histFolder->Add(S1L1pidE15);
+ histFolder->Add(S1L1pidEall);
+ histFolder->Add(S3L4edEall);
+ histFolder->Add(S3L4pidEall);
+
// Get pointer to the ROOT I/O manager
FairRootManager* rootMgr = FairRootManager::Instance();
if (NULL == rootMgr) {
@@ -81,27 +95,20 @@ InitStatus CbmTrdHitProducerQa::Init()
return kFATAL;
}
- // Get pointer to TRD point array
- fTrdPointCollection = (TClonesArray*) rootMgr->GetObject("TrdPoint");
- if (NULL == fTrdPointCollection) {
+ // Get a pointer to the previous already existing data level
+ fDigiMan = CbmDigiManager::Instance();
+ if (NULL == fDigiMan) {
cout << "-W- CbmTrdHitProducerQa::Init : "
- << "no TRD point array !" << endl;
+ << "no digi manager found !" << endl;
return kERROR;
}
+ fDigiMan->Init();
- // Get pointer to TRD digi array
- fTrdDigiCollection = (TClonesArray*) rootMgr->GetObject("TrdDigi");
- if (NULL == fTrdDigiCollection) {
- cout << "-W- CbmTrdHitProducerQa::Init : "
- << "no TRD digi array !" << endl;
- return kERROR;
- }
-
- // Get pointer to TRD digi array
- fTrdDigiMatchCollection = (TClonesArray*) rootMgr->GetObject("TrdDigiMatch");
- if (NULL == fTrdDigiMatchCollection) {
+ // Get pointer to TRD point array
+ fTrdPointCollection = (TClonesArray*) rootMgr->GetObject("TrdPoint");
+ if (NULL == fTrdPointCollection) {
cout << "-W- CbmTrdHitProducerQa::Init : "
- << "no TRD digi match array !" << endl;
+ << "no TRD point array !" << endl;
return kERROR;
}
@@ -120,6 +127,12 @@ InitStatus CbmTrdHitProducerQa::Init()
return kFATAL;
}
+ // Get pointer to TRD digi array match
+ if (!fDigiMan->IsMatchPresent(ECbmModuleId::kTrd)) {
+ cout << GetName() << ": no TRD match branch in digi manager." << endl;
+ return kERROR;
+ }
+
return kSUCCESS;
}
@@ -129,58 +142,36 @@ InitStatus CbmTrdHitProducerQa::Init()
// ---- Task execution ------------------------------------------------------
void CbmTrdHitProducerQa::Exec(Option_t*)
{
- // Declare variables
- CbmTrdHit* trdHit = NULL;
- // CbmTrdDigi* trdDigi = NULL;
- CbmMatch* trdDigiMatch = NULL;
- CbmTrdPoint* trdPoint = NULL;
- // CbmMCTrack* mctrack = NULL;
-
-
- Float_t hitPoolX = 0;
- Float_t hitPoolY = 0;
- Float_t hitPosX = 0;
- Float_t hitPosY = 0;
// Float_t hitErrX = 0;
// Float_t hitErrY = 0;
-
- Float_t pointPosX = 0;
- Float_t pointPosY = 0;
-
- Int_t plane = 0;
Int_t partID = 0;
-
Float_t momentum;
-
// Event counters
Int_t nTrdHit = fTrdHitCollection->GetEntriesFast();
// Loop over TRD hits
for (Int_t iHit = 0; iHit < nTrdHit; iHit++) {
- trdHit = (CbmTrdHit*) fTrdHitCollection->At(iHit);
+ const CbmTrdHit* trdHit = (CbmTrdHit*) fTrdHitCollection->At(iHit);
if (NULL == trdHit) continue;
- // This will have to change in the future, when the creation of the poin
- // will not be necessarily connected to existence of tyhe point
+ // This will have to change in the future, when the creation of the point
+ // will not be necessarily connected to existence of the point
- trdDigiMatch = (CbmMatch*) fTrdDigiMatchCollection->At(trdHit->GetRefId());
+ const CbmMatch* trdDigiMatch = fDigiMan->GetMatch(ECbmModuleId::kTrd, trdHit->GetRefId());
if (NULL == trdDigiMatch) continue;
- trdPoint = (CbmTrdPoint*) fTrdPointCollection->At(trdDigiMatch->GetMatchedLink().GetIndex());
+ const CbmTrdPoint* trdPoint = (CbmTrdPoint*) fTrdPointCollection->At(trdDigiMatch->GetMatchedLink().GetIndex());
if (NULL == trdPoint) continue;
- plane = trdHit->GetPlaneId();
+ const Int_t plane = trdHit->GetPlaneId();
if (plane == 1) {
-
- partID = (((CbmMCTrack*) fMCTrackArray->At(trdPoint->GetTrackID()))->GetPdgCode());
-
+ partID = (((CbmMCTrack*) fMCTrackArray->At(trdPoint->GetTrackID()))->GetPdgCode());
momentum = TMath::Sqrt((trdPoint->GetPx() * trdPoint->GetPx()) + (trdPoint->GetPy() * trdPoint->GetPy())
+ (trdPoint->GetPz() * trdPoint->GetPz()));
if ((TMath::Abs(partID) == 11) && (momentum > 1.25) && (momentum < 1.75)) {
-
S1L1eTR15->Fill(0); //(trdHit->GetELossTR())*1000000);
S1L1edEdx15->Fill((trdHit->GetELoss()) * 1000000);
S1L1edE15->Fill((trdHit->GetELoss()) * 1000000);
@@ -193,32 +184,26 @@ void CbmTrdHitProducerQa::Exec(Option_t*)
}
if (plane == 12) {
-
- partID = (((CbmMCTrack*) fMCTrackArray->At(trdPoint->GetTrackID()))->GetPdgCode());
-
+ partID = (((CbmMCTrack*) fMCTrackArray->At(trdPoint->GetTrackID()))->GetPdgCode());
momentum = TMath::Sqrt((trdPoint->GetPx() * trdPoint->GetPx()) + (trdPoint->GetPy() * trdPoint->GetPy())
+ (trdPoint->GetPz() * trdPoint->GetPz()));
-
if (TMath::Abs(partID) == 11) { S3L4edEall->Fill((trdHit->GetELoss()) * 1000000); }
if (TMath::Abs(partID) == 211) { S3L4pidEall->Fill((trdHit->GetELoss()) * 1000000); }
}
-
if (plane == 1) {
// compute the Hit pools for X and Y coordinate
- hitPosX = trdHit->GetX();
+ const Float_t hitPosX = trdHit->GetX();
// hitErrX = trdHit->GetDx();
// hitErrX /= 10000; // micrometers to centimeters
- pointPosX = trdPoint->GetX();
- hitPoolX = (hitPosX - pointPosX); ///hitErrX;
+ const Float_t pointPosX = trdPoint->GetX();
+ const Float_t hitPoolX = (hitPosX - pointPosX); ///hitErrX;
-
- hitPosY = trdHit->GetY();
+ const Float_t hitPosY = trdHit->GetY();
// hitErrY = trdHit->GetDy();
// hitErrY /= 10000; // micrometers to centimeters
- pointPosY = trdPoint->GetY();
- hitPoolY = (hitPosY - pointPosY); ///hitErrY;
-
+ const Float_t pointPosY = trdPoint->GetY();
+ const Float_t hitPoolY = (hitPosY - pointPosY); ///hitErrY;
// fill histograms
fHitPoolsX->Fill(hitPoolX);
@@ -226,82 +211,18 @@ void CbmTrdHitProducerQa::Exec(Option_t*)
}
}
}
-
// --------------------------------------------------------------------------
-
// ---- Finish --------------------------------------------------------------
void CbmTrdHitProducerQa::Finish() { WriteHistograms(); }
// --------------------------------------------------------------------------
-
-// ---- Prepare test histograms ---------------------------------------------
-
-/*
-void CbmTrdHitProducerQa::PrepareHistograms()
-{
-
- fHitPoolsX=NULL;
- fHitPoolsY=NULL;
-
- S1L1eTR15=NULL;
- S1L1edEdx15=NULL;
- S1L1edE15=NULL;
- S1L1edEall=NULL;
- S1L1pidE15=NULL;
- S1L1pidEall=NULL;
-
- S3L4eTR15=NULL;
- S3L4edEdx15=NULL;
- S3L4edE15=NULL;
- S3L4edEall=NULL;
- S3L4pidE15=NULL;
- S3L4pidEall=NULL;
-
- fHitPoolsX = new TH1F("fHitPoolsX", "", 500, -50, 50);
- fHitPoolsY = new TH1F("fHitPoolsY", "", 500, -50, 50);
-
- S1L1eTR15 = new TH1F("S1L1eTR15","TR of e- for first layer ",600,0.,60.);
- S1L1edEdx15 = new TH1F("S1L1edEdx15","dEdx of e- for first layer ",600,0.,60.);
- S1L1edE15 = new TH1F("S1L1edE15","dEdx+TR of e- for first layer ",600,0.,60.);
- S1L1edEall = new TH1F("S1L1edEall","dEdx+TR of e- for first layer ",600,0.,60.);
- S1L1pidE15 = new TH1F("S1L1pidE15","dEdx+TR of pi- for first layer ",600,0.,60.);
- S1L1pidEall = new TH1F("S1L1pidEall","dEdx+TR of pi- for first layer ",600,0.,60.);
- S3L4edEall = new TH1F("S3L4edEall","dEdx+TR of e- for layer 12",600,0.,60.);
- S3L4pidEall = new TH1F("S3L4pidEall","dEdx+TR of pi- for layer 12",600,0.,60.);
-
-}
-*/
-// --------------------------------------------------------------------------
-
-
// ---- Write test histograms ------------------------------------------------
-
void CbmTrdHitProducerQa::WriteHistograms()
{
- gDirectory->mkdir("TrdHitProducerQA");
- gDirectory->cd("TrdHitProducerQA");
-
- if (fHitPoolsX) fHitPoolsX->Write();
- if (fHitPoolsY) fHitPoolsY->Write();
-
- if (S1L1eTR15) S1L1eTR15->Write();
- if (S1L1edEdx15) S1L1edEdx15->Write();
- if (S1L1edE15) S1L1edE15->Write();
- if (S1L1edEall) S1L1edEall->Write();
- if (S1L1pidE15) S1L1pidE15->Write();
- if (S1L1pidEall) S1L1pidEall->Write();
-
- if (S3L4eTR15) S3L4eTR15->Write();
- if (S3L4edEdx15) S3L4edEdx15->Write();
- if (S3L4edE15) S3L4edE15->Write();
- if (S3L4edEall) S3L4edEall->Write();
- if (S3L4pidE15) S3L4pidE15->Write();
- if (S3L4pidEall) S3L4pidEall->Write();
-
- gDirectory->cd("..");
+ FairSink* sink = FairRootManager::Instance()->GetSink();
+ sink->WriteObject(&fOutFolder, nullptr);
}
-
// --------------------------------------------------------------------------
ClassImp(CbmTrdHitProducerQa);
diff --git a/reco/detectors/trd/qa/CbmTrdHitProducerQa.h b/reco/detectors/trd/qa/CbmTrdHitProducerQa.h
index c8838185ed..6d542a568c 100644
--- a/reco/detectors/trd/qa/CbmTrdHitProducerQa.h
+++ b/reco/detectors/trd/qa/CbmTrdHitProducerQa.h
@@ -23,6 +23,9 @@
#include "FairTask.h"
+#include <TFolder.h>
+
+class CbmDigiManager;
class TClonesArray;
class TH1F;
@@ -46,12 +49,14 @@ public:
/* Finish at the end of each event */
virtual void Finish();
-
private:
+ TFolder fOutFolder; /// output folder with histos and canvases
+ TFolder* histFolder; /// subfolder for histograms
+
+ CbmDigiManager* fDigiMan = nullptr;
+
/* Data branches*/
TClonesArray* fTrdHitCollection;
- TClonesArray* fTrdDigiCollection;
- TClonesArray* fTrdDigiMatchCollection;
TClonesArray* fTrdPointCollection;
TClonesArray* fMCTrackArray;
@@ -61,6 +66,8 @@ private:
/** Number of layers per station **/
Int_t fNoTrdPerStation;
+ /* Write test histograms */
+ void WriteHistograms();
/* Test histograms*/
TH1F* fHitPoolsX; // = ((Hit - Point) / HitError) in X
@@ -80,18 +87,7 @@ private:
TH1F* S3L4pidE15; //
TH1F* S3L4pidEall; //
-private:
- CbmTrdHitProducerQa(const CbmTrdHitProducerQa&);
- CbmTrdHitProducerQa& operator=(const CbmTrdHitProducerQa&);
-
- //public:
-
- /* Write test histograms */
- void WriteHistograms();
-
-
- ClassDef(CbmTrdHitProducerQa, 2)
+ ClassDef(CbmTrdHitProducerQa, 3)
};
-
#endif
--
GitLab