From 8644a01ff573e75fafaa5b85751ff77e6f8aec4f Mon Sep 17 00:00:00 2001 From: Martin Beyer <martin.beyer@physik.uni-giessen.de> Date: Tue, 7 Nov 2023 10:31:27 +0000 Subject: [PATCH] Rich: Update RichEventDisplay, rich/run macros & logging for HitProducer and Reco --- core/base/CbmHistManager.cxx | 10 +- core/base/CbmHistManager.h | 4 +- core/base/utils/CbmUtils.cxx | 2 + core/detectors/rich/CbmRichEventDisplay.cxx | 151 ++++---- core/detectors/rich/CbmRichEventDisplay.h | 114 +++--- macro/rich/run/.gitignore | 11 + macro/rich/run/batch_job.py | 64 ++-- macro/rich/run/batch_send.py | 11 +- macro/rich/run/run_digi.C | 68 ++-- macro/rich/run/run_many.py | 3 +- macro/rich/run/run_one.py | 72 ++-- macro/rich/run/run_qa.C | 207 ++++++---- macro/rich/run/run_reco.C | 355 +++++++++++------- macro/rich/run/run_transport.C | 104 +++-- reco/detectors/rich/CbmRichHitProducer.cxx | 39 +- reco/detectors/rich/CbmRichHitProducer.h | 1 + reco/detectors/rich/CbmRichReconstruction.cxx | 137 +++++-- reco/detectors/rich/CbmRichReconstruction.h | 16 +- .../CbmRichProjectionProducerAnalytical.cxx | 3 +- .../tracks/CbmRichProjectionProducerBase.h | 6 + .../tracks/CbmRichProjectionProducerTGeo.cxx | 3 +- 21 files changed, 863 insertions(+), 518 deletions(-) create mode 100644 macro/rich/run/.gitignore diff --git a/core/base/CbmHistManager.cxx b/core/base/CbmHistManager.cxx index 86d702c24f..c16b33d9df 100755 --- a/core/base/CbmHistManager.cxx +++ b/core/base/CbmHistManager.cxx @@ -167,11 +167,15 @@ void CbmHistManager::AddTDirectoryObject(TObject* obj) void CbmHistManager::Clear(Option_t*) { - map<string, TNamed*>::iterator it; - for (it = fMap.begin(); it != fMap.end(); it++) { - delete (*it).second; + for (auto& pair : fMap) { + delete pair.second; } fMap.clear(); + + for (auto& canvas : fCanvases) { + delete canvas; + } + fCanvases.clear(); } void CbmHistManager::ShrinkEmptyBinsH1(const string& histName) diff --git a/core/base/CbmHistManager.h b/core/base/CbmHistManager.h index d0eaeb3945..fae5f81a52 100644 --- a/core/base/CbmHistManager.h +++ b/core/base/CbmHistManager.h @@ -446,9 +446,9 @@ public: void AddTDirectoryObject(TObject* obj); /** - * \brief Clear memory. Remove all histograms. + * \brief Clear memory. Remove all histograms and canvases. */ - void Clear(Option_t*); + void Clear(Option_t* = ""); /** * \brief Shrink empty bins in H1. diff --git a/core/base/utils/CbmUtils.cxx b/core/base/utils/CbmUtils.cxx index 7fea819a9f..c260f549dd 100644 --- a/core/base/utils/CbmUtils.cxx +++ b/core/base/utils/CbmUtils.cxx @@ -39,6 +39,8 @@ namespace Cbm SaveCanvasAsImageImpl("eps", c, dir, option); SaveCanvasAsImageImpl("png", c, dir, option); SaveCanvasAsImageImpl("gif", c, dir, option); + SaveCanvasAsImageImpl("pdf", c, dir, option); + SaveCanvasAsImageImpl("svg", c, dir, option); } void SaveCanvasAsImageImpl(const string& imageType, TCanvas* c, const string& dir, const string& option) diff --git a/core/detectors/rich/CbmRichEventDisplay.cxx b/core/detectors/rich/CbmRichEventDisplay.cxx index 3daa5c9967..1bbb34eeb6 100644 --- a/core/detectors/rich/CbmRichEventDisplay.cxx +++ b/core/detectors/rich/CbmRichEventDisplay.cxx @@ -1,6 +1,6 @@ -/* Copyright (C) 2006-2020 UGiessen, JINR-LIT +/* Copyright (C) 2006-2023 UGiessen, JINR-LIT SPDX-License-Identifier: GPL-3.0-only - Authors: Supriya Das, Semen Lebedev [committer], Florian Uhlig */ + Authors: Supriya Das, Semen Lebedev [committer], Martin Beyer, Florian Uhlig */ /** * \file CbmRichEventDisplay.cxx @@ -10,74 +10,49 @@ **/ #include "CbmRichEventDisplay.h" -#include "CbmDrawHist.h" // for DrawH2, SetDefaultDrawStyle -#include "CbmHistManager.h" // for CbmHistManager -#include "CbmRichHit.h" // for CbmRichHit -#include "CbmRichPoint.h" // for CbmRichPoint -#include "CbmRichRing.h" // for CbmRichRing - -#include <FairRootManager.h> // for FairRootManager -#include <FairTask.h> // for FairTask, InitStatus, kSUCCESS -#include <FairTrackParam.h> // for FairTrackParam - -#include <TAxis.h> // for TAxis -#include <TCanvas.h> // for TCanvas -#include <TClonesArray.h> // for TClonesArray -#include <TEllipse.h> // for TEllipse -#include <TGenericClassInfo.h> // for TGenericClassInfo -#include <TH2.h> // for TH2D -#include <TMarker.h> // for TMarker -#include <TVirtualPad.h> // for TVirtualPad, gPad - -#include <iostream> // for string, stringstream, operator<<, bas... -#include <string> // for operator==, basic_string, char_traits - -using namespace std; - -CbmRichEventDisplay::CbmRichEventDisplay() - : FairTask("CbmRichEventDisplay") - , fRichRings(nullptr) - , fRichHits(nullptr) - , fRichPoints(nullptr) - , fRichMatches(nullptr) - , fRichProjections(nullptr) - , fMcTracks(nullptr) - , fOutputDir("") - , fHM(nullptr) - , fEventNum(0) - , fDrawRings(true) - , fDrawHits(true) - , fDrawPoints(true) - , fDrawProjections(true) -{ - SetDefaultDrawStyle(); -} +#include "CbmDigiManager.h" +#include "CbmDrawHist.h" +#include "CbmEvent.h" +#include "CbmHistManager.h" +#include "CbmRichGeoManager.h" +#include "CbmRichHit.h" +#include "CbmRichPoint.h" +#include "CbmRichRing.h" -CbmRichEventDisplay::~CbmRichEventDisplay() {} +#include <FairRootManager.h> +#include <FairTrackParam.h> +#include <TCanvas.h> +#include <TClonesArray.h> +#include <TEllipse.h> +#include <TH2.h> +#include <TMarker.h> + +#include <string> InitStatus CbmRichEventDisplay::Init() { FairRootManager* ioman = FairRootManager::Instance(); - if (nullptr == ioman) { Fatal("CbmRichEventDisplay::Init", "RootManager not instantiated!"); } + if (!ioman) LOG(fatal) << "CbmRichEventDisplay::Init: RootManager not instantiated!"; - fRichHits = (TClonesArray*) ioman->GetObject("RichHit"); - if (nullptr == fRichHits) { Fatal("CbmRichEventDisplay::Init", "No RichHit array!"); } + fDigiManager = CbmDigiManager::Instance(); + fDigiManager->Init(); + if (!fDigiManager->IsPresent(ECbmModuleId::kRich)) LOG(fatal) << "CbmRichEventDisplay::Init: No RichDigi array!"; - fRichRings = (TClonesArray*) ioman->GetObject("RichRing"); - if (nullptr == fRichRings) { Fatal("CbmRichEventDisplay::Init", "No RichRing array!"); } + fCbmEvents = static_cast<TClonesArray*>(ioman->GetObject("CbmEvent")); + if (!fCbmEvents) LOG(fatal) << "CbmRichEventDisplay::Init: No CbmEvent array!"; - fRichPoints = (TClonesArray*) ioman->GetObject("RichPoint"); - if (nullptr == fRichPoints) { Fatal("CbmRichEventDisplay::Init", "No RichPoint array!"); } + fRichPoints = static_cast<TClonesArray*>(ioman->GetObject("RichPoint")); + if (!fRichPoints && fDrawPoints) LOG(fatal) << "CbmRichEventDisplay::Init: No RichPoint array!"; - fRichMatches = (TClonesArray*) ioman->GetObject("RichRingMatch"); - if (nullptr == fRichMatches) { Fatal("CbmRichEventDisplay::Init", "No RichRingMatch array!"); } + fRichHits = static_cast<TClonesArray*>(ioman->GetObject("RichHit")); + if (!fRichHits && fDrawHits) LOG(fatal) << "CbmRichEventDisplay::Init: No RichHit array!"; - fRichProjections = (TClonesArray*) ioman->GetObject("RichProjection"); - if (nullptr == fRichProjections) { Fatal("CbmRichEventDisplay::Init", "No RichProjection array!"); } + fRichRings = static_cast<TClonesArray*>(ioman->GetObject("RichRing")); + if (!fRichRings && fDrawRings) LOG(fatal) << "CbmRichEventDisplay::Init: No RichRing array!"; - fMcTracks = (TClonesArray*) ioman->GetObject("MCTrack"); - if (nullptr == fMcTracks) { Fatal("CbmRichEventDisplay::Init", "No MCTrack array!"); } + fRichProjections = static_cast<TClonesArray*>(ioman->GetObject("RichProjection")); + if (!fRichProjections && fDrawProjections) LOG(fatal) << "CbmRichEventDisplay::Init: No RichProjection array!"; fHM = new CbmHistManager(); @@ -86,14 +61,17 @@ InitStatus CbmRichEventDisplay::Init() void CbmRichEventDisplay::Exec(Option_t* /*opt*/) { - fEventNum++; SetDefaultDrawStyle(); - DrawOneEvent(); + for (Int_t i = 0; i < fCbmEvents->GetEntriesFast(); i++) { + fEventNum++; + CbmEvent* event = static_cast<CbmEvent*>(fCbmEvents->At(i)); + DrawOneEvent(event); + } } -void CbmRichEventDisplay::DrawOneEvent() +void CbmRichEventDisplay::DrawOneEvent(CbmEvent* event) { - stringstream ss; + std::stringstream ss; ss << "rich_event_display_event_" << fEventNum; TCanvas* c = fHM->CreateCanvas(ss.str().c_str(), ss.str().c_str(), 800, 800); c->Divide(1, 2); @@ -103,7 +81,7 @@ void CbmRichEventDisplay::DrawOneEvent() padU->GetYaxis()->SetTitleOffset(0.75); gPad->SetLeftMargin(0.1); gPad->SetRightMargin(0.05); - DrawOnePmtPlane("up"); + DrawOnePmtPlane("up", event); c->cd(2); TH2D* padD = new TH2D("padD", ";x [cm];y [cm]", 1, -120., 120., 1, -210., -120.); @@ -111,17 +89,19 @@ void CbmRichEventDisplay::DrawOneEvent() padD->GetYaxis()->SetTitleOffset(0.75); gPad->SetLeftMargin(0.1); gPad->SetRightMargin(0.05); - DrawOnePmtPlane("down"); + DrawOnePmtPlane("down", event); } -void CbmRichEventDisplay::DrawOnePmtPlane(const string& plane) +void CbmRichEventDisplay::DrawOnePmtPlane(const std::string& plane, CbmEvent* event) { - //Draw Track projections + // Draw Track projections if (fDrawProjections) { - int nofProjections = fRichProjections->GetEntriesFast(); + int nofProjections = event->GetNofData(ECbmDataType::kRichTrackProjection); for (int iP = 0; iP < nofProjections; iP++) { - FairTrackParam* pr = (FairTrackParam*) fRichProjections->At(iP); - if (nullptr == pr) continue; + FairTrackParam* pr = + (FairTrackParam*) fRichProjections->At(event->GetIndex(ECbmDataType::kRichTrackProjection, iP)); + if (!pr) continue; + // FIXME: only draw successful projections if ((plane == "up" && pr->GetY() >= 0.) || (plane == "down" && pr->GetY() < 0.)) { TMarker* m = new TMarker(pr->GetX(), pr->GetY(), 3.); m->SetMarkerSize(0.7); @@ -133,12 +113,11 @@ void CbmRichEventDisplay::DrawOnePmtPlane(const string& plane) // Draw hits if (fDrawHits) { - int nofHits = fRichHits->GetEntriesFast(); + int nofHits = event->GetNofData(ECbmDataType::kRichHit); for (int iH = 0; iH < nofHits; iH++) { - CbmRichHit* hit = (CbmRichHit*) fRichHits->At(iH); - if (nullptr == hit) continue; + CbmRichHit* hit = (CbmRichHit*) fRichHits->At(event->GetIndex(ECbmDataType::kRichHit, iH)); + if (!hit) continue; if ((plane == "up" && hit->GetY() >= 0.) || (plane == "down" && hit->GetY() < 0.)) { - TEllipse* hitDr = new TEllipse(hit->GetX(), hit->GetY(), 0.6); hitDr->SetFillColor(kRed); hitDr->SetLineColor(kRed); @@ -149,31 +128,39 @@ void CbmRichEventDisplay::DrawOnePmtPlane(const string& plane) // Draw rings if (fDrawRings) { - int nofRings = fRichRings->GetEntriesFast(); + int nofRings = event->GetNofData(ECbmDataType::kRichRing); for (int iR = 0; iR < nofRings; iR++) { - CbmRichRing* ring = (CbmRichRing*) fRichRings->At(iR); - if (nullptr == ring) continue; + CbmRichRing* ring = (CbmRichRing*) fRichRings->At(event->GetIndex(ECbmDataType::kRichRing, iR)); + if (!ring) continue; if ((plane == "up" && ring->GetCenterY() >= 0.) || (plane == "down" && ring->GetCenterY() < 0.)) { - DrawCircle(ring); + DrawEllipse(ring); } } } // Draw RICH MC Points + // Slow for large ts if (fDrawPoints) { - int nofPoints = fRichPoints->GetEntriesFast(); + Double_t evStart = event->GetStartTime(); + Double_t evStop = event->GetEndTime(); + int nofPoints = fRichPoints->GetEntriesFast(); for (int iP = 0; iP < nofPoints; iP++) { CbmRichPoint* point = (CbmRichPoint*) fRichPoints->At(iP); - if (nullptr == point) continue; + if (!point) continue; + Double_t pointTime = point->GetTime(); + if (pointTime > evStart && pointTime < evStop) continue; if ((plane == "up" && point->GetY() >= 0.) || (plane == "down" && point->GetY() < 0.)) { - TEllipse* pointDr = new TEllipse(point->GetX(), point->GetY(), 0.4); + TVector3 posPoint(point->GetX(), point->GetY(), point->GetZ()); + TVector3 detPoint; + CbmRichGeoManager::GetInstance().RotatePoint(&posPoint, &detPoint, false); + TEllipse* pointDr = new TEllipse(detPoint.X(), detPoint.Y(), 0.2); pointDr->Draw(); } } } } -void CbmRichEventDisplay::DrawCircle(CbmRichRing* ring) +void CbmRichEventDisplay::DrawEllipse(CbmRichRing* ring) { TEllipse* circle = new TEllipse(ring->GetCenterX(), ring->GetCenterY(), ring->GetRadius()); circle->SetFillStyle(0); @@ -186,6 +173,6 @@ void CbmRichEventDisplay::DrawCircle(CbmRichRing* ring) center->Draw(); } -void CbmRichEventDisplay::Finish() { fHM->SaveCanvasToImage(fOutputDir); } +void CbmRichEventDisplay::Finish() { fHM->SaveCanvasToImage(fOutputDir, fOutputFormat); } ClassImp(CbmRichEventDisplay) diff --git a/core/detectors/rich/CbmRichEventDisplay.h b/core/detectors/rich/CbmRichEventDisplay.h index 779bf39f40..86c770ee09 100644 --- a/core/detectors/rich/CbmRichEventDisplay.h +++ b/core/detectors/rich/CbmRichEventDisplay.h @@ -1,6 +1,6 @@ -/* Copyright (C) 2012-2020 UGiessen/JINR-LIT, Giessen/Dubna +/* Copyright (C) 2012-2023 UGiessen/JINR-LIT, Giessen/Dubna SPDX-License-Identifier: GPL-3.0-only - Authors: Semen Lebedev [committer] */ + Authors: Semen Lebedev [committer], Martin Beyer */ /** * \file CbmRichEventDisplay.h @@ -14,19 +14,17 @@ #ifndef CBM_RICH_EVENT_DISPLAY #define CBM_RICH_EVENT_DISPLAY -#include "CbmHistManager.h" // IWYU pragma: keep needed by RootCling +#include "CbmDigiManager.h" +#include "CbmHistManager.h" -#include "FairTask.h" // for FairTask, InitStatus +#include <FairTask.h> -#include <Rtypes.h> // for THashConsistencyHolder, ClassDef -#include <RtypesCore.h> // for Option_t -#include <TClonesArray.h> // IWYU pragma: keep needed by RootCling +#include <TClonesArray.h> -#include <string> // for string +#include <string> class CbmRichRing; -class TClonesArray; -class TMemberInspector; +class CbmEvent; /** * \class CbmRichEventDisplay @@ -38,34 +36,37 @@ class TMemberInspector; **/ class CbmRichEventDisplay : public FairTask { public: - /** - * \brief Default constructor. - */ - CbmRichEventDisplay(); + /** Default constructor */ + CbmRichEventDisplay() : FairTask("CbmRichEventDisplay") {}; - /** - * \brief Destructor. - */ - virtual ~CbmRichEventDisplay(); + /** Destructor */ + ~CbmRichEventDisplay() = default; - /** - * \brief Inherited from FairTask. - */ - virtual InitStatus Init(); + /** Copy constructor (disabled) */ + CbmRichEventDisplay(const CbmRichEventDisplay&) = delete; - /** - * \brief Inherited from FairTask. - */ - virtual void Exec(Option_t* opt); + /** Assignment operator (disabled) */ + CbmRichEventDisplay& operator=(const CbmRichEventDisplay&) = delete; - /** - * \brief Inherited from FairTask. - */ - virtual void Finish(); + /** Inherited from FairTask */ + InitStatus Init(); + /** Inherited from FairTask */ + void Exec(Option_t* opt); + + /** Inherited from FairTask */ + void Finish(); + + /** Draw rings (Enabled by default) */ void SetDrawRings(bool b) { fDrawRings = b; } + + /** Draw hits (Enabled by default) */ void SetDrawHits(bool b) { fDrawHits = b; } + + /** Draw points (Disabled by default) */ void SetDrawPoints(bool b) { fDrawPoints = b; } + + /** Draw projections (Enabled by default) */ void SetDrawProjections(bool b) { fDrawProjections = b; } /** @@ -74,43 +75,40 @@ public: */ void SetOutputDir(const std::string& dir) { fOutputDir = dir; } -private: - TClonesArray* fRichRings; - TClonesArray* fRichHits; - TClonesArray* fRichPoints; - TClonesArray* fRichMatches; - TClonesArray* fRichProjections; - - TClonesArray* fMcTracks; - - std::string fOutputDir; // output dir for results - CbmHistManager* fHM; + /** + * \brief Set output file format (any combination of png, gif, eps, svg and pdf). + * \param[in] format File format, e.g. png or png,pdf. (Default: png,pdf) + */ + void SetOutputFormat(const std::string& format) { fOutputFormat = format; } - int fEventNum; +private: + int fEventNum {}; - bool fDrawRings; - bool fDrawHits; - bool fDrawPoints; - bool fDrawProjections; + CbmDigiManager* fDigiManager {nullptr}; + CbmHistManager* fHM {nullptr}; - void DrawOneEvent(); + TClonesArray* fCbmEvents {nullptr}; + TClonesArray* fRichPoints {nullptr}; + TClonesArray* fRichHits {nullptr}; + TClonesArray* fRichRings {nullptr}; + TClonesArray* fRichMatches {nullptr}; + TClonesArray* fRichProjections {nullptr}; - void DrawOnePmtPlane(const std::string& plane); + std::string fOutputDir {"RichEventDisplays"}; + std::string fOutputFormat {"png,pdf"}; - void DrawCircle(CbmRichRing* ring); + bool fDrawRings {true}; + bool fDrawHits {true}; + bool fDrawPoints {false}; + bool fDrawProjections {true}; + void DrawOneEvent(CbmEvent* event); - /** - * \brief Copy constructor. - */ - CbmRichEventDisplay(const CbmRichEventDisplay&); + void DrawOnePmtPlane(const std::string& plane, CbmEvent* event); - /** - * \brief Assignment operator. - */ - CbmRichEventDisplay& operator=(const CbmRichEventDisplay&); + void DrawEllipse(CbmRichRing* ring); - ClassDef(CbmRichEventDisplay, 1); + ClassDef(CbmRichEventDisplay, 2); }; #endif diff --git a/macro/rich/run/.gitignore b/macro/rich/run/.gitignore new file mode 100644 index 0000000000..17b28fd1e1 --- /dev/null +++ b/macro/rich/run/.gitignore @@ -0,0 +1,11 @@ +*.root +*.dat +all_*.par +core_dump_* +*.log +*.txt +*.eps +*.pdf +*.png +*.svg +*.gif diff --git a/macro/rich/run/batch_job.py b/macro/rich/run/batch_job.py index 0af067e428..5b2b04ca92 100755 --- a/macro/rich/run/batch_job.py +++ b/macro/rich/run/batch_job.py @@ -6,40 +6,52 @@ import shutil dataDir = sys.argv[1] geoSetup = sys.argv[2] -cbmrootConfigPath = "/lustre/nyx/cbm/users/slebedev/cbm/trunk/build/config.sh" -macroDir = "/lustre/nyx/cbm/users/slebedev/cbm/trunk/cbmroot/macro/rich/" +cbmrootConfigPath = "/lustre/cbm/users/mabeyer/cbmroot/build/config.sh" +macroDir = "/lustre/cbm/users/mabeyer/cbmroot/macro/rich/run" taskId = os.environ.get('SLURM_ARRAY_TASK_ID') jobId = os.environ.get('SLURM_ARRAY_JOB_ID') print("dataDir:" + dataDir) -os.system(("source {}").format(cbmrootConfigPath)) - workDir = dataDir + "/workdir/" + jobId + "_" + taskId + "/" if os.path.exists(workDir): shutil.rmtree(workDir) os.makedirs(workDir) os.chdir(workDir) - -runRecoqaUrqmd = True - -if runRecoqaUrqmd: - nofEvents = 1000 - dataDirNew = dataDir + "/recoqa_urqmd8gev/" - if not os.path.exists(dataDirNew): os.makedirs(dataDirNew) - urqmdFile = "/lustre/nyx/cbm/prod/gen/urqmd/auau/8gev/centr/urqmd.auau.8gev.centr." + str(taskId).zfill(5) + ".root" - mcFile = dataDirNew + "/mc."+ taskId + ".root" - parFile = dataDirNew + "/param."+ taskId + ".root" - digiFile = dataDirNew + "/digi."+ taskId + ".root" - recoFile = dataDirNew + "/reco."+ taskId + ".root" - qaFile = dataDirNew + "/qa."+ taskId + ".root" - geoSimFile = dataDirNew + "/geosim."+ taskId + ".root" - resultDir = dataDirNew + "/results/results_" + geoSetup + "/" - nofElectrons = 5 - nofPositrons = 5 - plutoFile = "" - os.system( ('root -l -b -q {}/run/run_sim.C\(\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",{},{},\\"{}\\",\\"{}\\",{}\)').format(macroDir, urqmdFile, mcFile, parFile, geoSimFile, nofElectrons, nofPositrons, plutoFile, geoSetup, nofEvents) ) - os.system( ('root -l -b -q {}/run/run_digi.C\(\\"{}\\",\\"{}\\",\\"{}\\",{}\)').format(macroDir, mcFile, parFile, digiFile, nofEvents) ) - os.system( ('root -l -b -q {}/run/run_reco.C\(\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",{}\)').format(macroDir, mcFile, parFile, digiFile, recoFile, geoSetup, resultDir, nofEvents) ) - os.system( ('root -l -b -q {}/run/run_qa.C\(\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",{}\)').format(macroDir, mcFile, parFile, digiFile, recoFile, qaFile, geoSetup, resultDir, nofEvents) ) +# Settings +nEvents = 1000 +urqmdFile = "/lustre/cbm/prod/gen/urqmd/auau/10gev/centr/urqmd.auau.10gev.centr." + str(taskId).zfill(5) + ".root" +plutoFile = "" +nElectrons = 5 +nPositrons = 5 +engine = 0 # 0: geant3 1: geant4 + +randomSeed = 0 + +eventRate = -1e7 +tsLength = -1 +eventBuilder = "Ideal" +useMC = 1 +monitor = 0 +resultDir = "" + +# Files +fileName = "data" +traFile = dataDir + fileName + taskId + ".tra.root" +parFile = dataDir + fileName + taskId + ".par.root" +geoFile = dataDir + fileName + taskId + ".geo.root" +digiFile = dataDir + fileName + taskId + ".digi.root" +recoFile = dataDir + fileName + taskId + ".reco.root" +qaFile = dataDir + fileName + taskId + ".qa.root" + + +traCmd=('root -l -b -q {}/run_transport.C\(\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",{},{},{},\\"{}\\",\\"{}\\",{},{}\)').format(macroDir, urqmdFile, traFile, parFile, geoFile, nEvents, nElectrons, nPositrons, plutoFile, geoSetup, engine, randomSeed) +digiCmd=('root -l -b -q {}/run_digi.C\(\\"{}\\",\\"{}\\",\\"{}\\",{},{},{},{},{}\)').format(macroDir, traFile, parFile, digiFile, nEvents, eventRate, tsLength, randomSeed, monitor) +recoCmd=('root -l -b -q {}/run_reco.C\(\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",{},\\"{}\\",\\"{}\\",{},{}\)').format(macroDir, traFile, parFile, digiFile, recoFile, nEvents, geoSetup, eventBuilder, useMC, monitor) +qaCmd=('root -l -b -q {}/run_qa.C\(\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",{},\\"{}\\",\\"{}\\",{}\)').format(macroDir, traFile, parFile, digiFile, recoFile, qaFile, nEvents, geoSetup, resultDir, monitor) + +os.system((". {} -a ; {}").format(cbmrootConfigPath, traCmd)) +os.system((". {} -a ; {}").format(cbmrootConfigPath, digiCmd)) +os.system((". {} -a ; {}").format(cbmrootConfigPath, recoCmd)) +os.system((". {} -a ; {}").format(cbmrootConfigPath, qaCmd)) diff --git a/macro/rich/run/batch_send.py b/macro/rich/run/batch_send.py index 3a49127984..e79c2045cd 100755 --- a/macro/rich/run/batch_send.py +++ b/macro/rich/run/batch_send.py @@ -3,12 +3,12 @@ import os import shutil -nofJobs = 100 +nJobs = 100 timeLimit="08:00:00" -geoSetup="sis100_electron_rich_pcarb_bcarb" -#All data in data dir will be removed +geoSetup="sis100_electron" +# All data in data dir will be removed removeData = False -dataDir = "/lustre/nyx/cbm/users/slebedev/cbm/data/rich_carb/" + geoSetup + "/" +dataDir = "/lustre/cbm/users/mabeyer/cbmdata/test/" logFile = dataDir + "/log/log_slurm-%A_%a.out" errorFile = dataDir + "/error/error_slurm-%A_%a.out" jobName = "RICH" @@ -30,7 +30,6 @@ os.makedirs(os.path.dirname(errorFile)) #-p debug commandStr=('sbatch --job-name={} --time={} --output={} --error={} --array=1-{} batch_job.py {} {}').format(jobName, timeLimit, - logFile, errorFile, nofJobs, dataDir, geoSetup) + logFile, errorFile, nJobs, dataDir, geoSetup) #print(commandStr) os.system(commandStr) - diff --git a/macro/rich/run/run_digi.C b/macro/rich/run/run_digi.C index 20093cb262..c43a07e903 100644 --- a/macro/rich/run/run_digi.C +++ b/macro/rich/run/run_digi.C @@ -1,45 +1,71 @@ -/* Copyright (C) 2018-2020 UGiessen, JINR-LIT +/* Copyright (C) 2018-2023 UGiessen, JINR-LIT SPDX-License-Identifier: GPL-3.0-only - Authors: Semen Lebedev [committer], Andrey Lebedev, Semen Lebedev [committer] */ + Authors: Semen Lebedev [committer], Andrey Lebedev, Martin Beyer */ -void run_digi(const string& traFile = "/Users/slebedev/Development/cbm/data/sim/rich/reco/tra.0.root", - const string& parFile = "/Users/slebedev/Development/cbm/data/sim/rich/reco/par.0.root", - const string& digiFile = "/Users/slebedev/Development/cbm/data/sim/rich/reco/digi.0.root", - int nEvents = 100) +#if !defined(__CLING__) +#include "CbmDigitization.h" +#include "CbmRichDigitizer.h" + +#include <FairLogger.h> + +#include <TRandom.h> +#include <TStopwatch.h> +#include <TTree.h> + +#include <iostream> +#endif + +void run_digi(TString traFile = "data/test.tra.root", TString parFile = "data/test.par.root", + TString digiFile = "data/test.digi.root", Int_t nEvents = -1, Double_t eventRate = -1.e7, + Double_t tsLength = -1, Int_t randomSeed = 0, Bool_t monitor = true) { TTree::SetMaxTreeSize(90000000000); - FairLogger::GetLogger()->SetLogScreenLevel("INFO"); - FairLogger::GetLogger()->SetLogVerbosityLevel("LOW"); - Double_t eventRate = 1.e7; - Double_t timeSliceLength = 1.e4; - Bool_t eventMode = true; - Bool_t overwrite = true; + remove(digiFile); + + gRandom->SetSeed(randomSeed); - remove(digiFile.c_str()); + // ----- Logger settings ---------------------------------------------- + FairLogger::GetLogger()->SetLogScreenLevel("INFO"); + FairLogger::GetLogger()->SetLogVerbosityLevel("LOW"); + // ------------------------------------------------------------------------ + // ----- Timer -------------------------------------------------------- TStopwatch timer; timer.Start(); + // ------------------------------------------------------------------------ + // ----- Digitization run --------------------------------------------- CbmDigitization run; - run.AddInput(traFile.c_str(), eventRate); - run.SetOutputFile(digiFile.c_str(), overwrite); - run.SetParameterRootFile(parFile.c_str()); - run.SetTimeSliceLength(timeSliceLength); - run.SetEventMode(eventMode); + cbm::sim::Mode mode = (eventRate < 0. ? cbm::sim::Mode::EventByEvent : cbm::sim::Mode::Timebased); + cbm::sim::TimeDist timeDist = cbm::sim::TimeDist::Poisson; + run.AddInput(0, traFile, timeDist, eventRate); + run.SetOutputFile(digiFile); + if (monitor) run.SetMonitorFile(" "); + run.SetParameterRootFile(parFile); + run.GenerateRunInfo(kFALSE); + run.SetTimeSliceLength(tsLength); + run.SetMode(mode); run.SetProduceNoise(false); - run.SetMonitorFile(""); + // run.DeactivateAllBut(ECbmModuleId::kRich); + // ----- Rich Digitizer ----------------------------------------------- CbmRichDigitizer* richDigitizer = new CbmRichDigitizer(); run.SetDigitizer(ECbmModuleId::kRich, richDigitizer, "RichPoint", true); + // ------------------------------------------------------------------------ - run.Run(nEvents); + if (nEvents < 0) run.Run(); + else + run.Run(nEvents); + // ------------------------------------------------------------------------ + // ----- Finish ------------------------------------------------------- timer.Stop(); std::cout << std::endl; std::cout << "Macro finished successfully." << std::endl; std::cout << "Digi file is " << digiFile << std::endl; std::cout << "Parameter file is " << parFile << std::endl; std::cout << "Real time " << timer.RealTime() << " s, CPU time " << timer.CpuTime() << " s" << std::endl; - std::cout << "Test passed" << std::endl << "All ok" << std::endl; + std::cout << std::endl << "Test passed" << std::endl << "All ok" << std::endl; + // ------------------------------------------------------------------------ } diff --git a/macro/rich/run/run_many.py b/macro/rich/run/run_many.py index be6fb4e5c2..857d6a4e72 100755 --- a/macro/rich/run/run_many.py +++ b/macro/rich/run/run_many.py @@ -1,6 +1,5 @@ import os -for taskId in range(1,4): +for taskId in range(4): commandXterm = ('xterm -hold -e python3 run_one.py {}&').format(taskId) os.system(commandXterm) - diff --git a/macro/rich/run/run_one.py b/macro/rich/run/run_one.py index 2068bca28c..980a384b92 100755 --- a/macro/rich/run/run_one.py +++ b/macro/rich/run/run_one.py @@ -4,30 +4,48 @@ import os import sys runId = sys.argv[1] -cbmrootConfigPath = "/Users/slebedev/Development/cbm/git/build/config.sh" -macroDir = "/Users/slebedev/Development/cbm/git/cbmroot/macro/rich/" -nofEvents = 100 - -urqmdFile = "/Users/slebedev/Development/cbm/data/urqmd/auau/10gev/centr/urqmd.auau.10gev.centr.00001.root" -dataDir= "/Users/slebedev/Development/cbm/data/sim/rich/reco/" -traFile = dataDir + "/tra."+ runId + ".root" -parFile = dataDir + "/par."+ runId + ".root" -digiFile = dataDir + "/digi."+ runId + ".root" -recoFile = dataDir + "/reco."+ runId + ".root" -qaFile = dataDir + "/qa."+ runId + ".root" -geoSimFile = dataDir + "/geosim."+ runId + ".root" - -# geoSetup = "sis100_electron" -geoSetup = "sis100_electron_" + runId -resultDir = "results_recoqa_" + runId + "/" -targetZ = -44.0 - -traCmd=('root -l -b -q {}/run/run_transport.C\(\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",{},{},\\"{}\\",\\"{}\\",{},{}\)').format(macroDir, urqmdFile, traFile, parFile, geoSimFile, 5, 5, "", geoSetup, nofEvents, targetZ) -digiCmd=('root -l -b -q {}/run/run_digi.C\(\\"{}\\",\\"{}\\",\\"{}\\",{}\)').format(macroDir, traFile, parFile, digiFile, nofEvents) -recoCmd=('root -l -b -q {}/run/run_reco.C\(\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",{}\)').format(macroDir, traFile, parFile, digiFile, recoFile, geoSetup, resultDir, nofEvents) -qaCmd=('root -l -b -q {}/run/run_qa.C\(\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",{}\)').format(macroDir, traFile, parFile, digiFile, recoFile, qaFile, geoSetup, resultDir, nofEvents) - -os.system((". /{} -a; {}").format(cbmrootConfigPath, traCmd)) -os.system((". /{} -a; {}").format(cbmrootConfigPath, digiCmd)) -os.system((". /{} -a; {}").format(cbmrootConfigPath, recoCmd)) -os.system((". /{} -a; {}").format(cbmrootConfigPath, qaCmd)) +cbmrootConfigPath = "/home/martin/cbm/cbmroot/build/config.sh" +macroDir = "/home/martin/cbm/cbmroot/macro/rich/run" + +# Settings +nEvents = 100 +urqmdFile = "default10gev" +plutoFile = "" +nElectrons = 5 +nPositrons = 5 +engine = 0 # 0: geant3 1: geant4 + +geoSetup = "sis100_electron" +randomSeed = 0 + +eventRate = -1e7 +tsLength = -1 +eventBuilder = "Ideal" +useMC = 1 +monitor = 0 +resultDir = "" +# resultDir = "results_recoqa_" + runId + "/" + +# Files +dataDir = "/home/martin/cbm/cbmroot/macro/rich/run/data/" +fileName = "test" +traFile = dataDir + fileName + runId + ".tra.root" +parFile = dataDir + fileName + runId + ".par.root" +geoFile = dataDir + fileName + runId + ".geo.root" +digiFile = dataDir + fileName + runId + ".digi.root" +recoFile = dataDir + fileName + runId + ".reco.root" +qaFile = dataDir + fileName + runId + ".qa.root" +# traFile = "/home/martin/cbm/cbmroot/macro/rich/run/data/test.tra.root" +# parFile = "/home/martin/cbm/cbmroot/macro/rich/run/data/test.par.root" +# digiFile = "/home/martin/cbm/cbmroot/macro/rich/run/data/test.digi.root" + + +traCmd=('root -l -b -q {}/run_transport.C\(\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",{},{},{},\\"{}\\",\\"{}\\",{},{}\)').format(macroDir, urqmdFile, traFile, parFile, geoFile, nEvents, nElectrons, nPositrons, plutoFile, geoSetup, engine, randomSeed) +digiCmd=('root -l -b -q {}/run_digi.C\(\\"{}\\",\\"{}\\",\\"{}\\",{},{},{},{},{}\)').format(macroDir, traFile, parFile, digiFile, nEvents, eventRate, tsLength, randomSeed, monitor) +recoCmd=('root -l -b -q {}/run_reco.C\(\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",{},\\"{}\\",\\"{}\\",{},{}\)').format(macroDir, traFile, parFile, digiFile, recoFile, nEvents, geoSetup, eventBuilder, useMC, monitor) +qaCmd=('root -l -b -q {}/run_qa.C\(\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",\\"{}\\",{},\\"{}\\",\\"{}\\",{}\)').format(macroDir, traFile, parFile, digiFile, recoFile, qaFile, nEvents, geoSetup, resultDir, monitor) + +os.system((". {} -a ; {}").format(cbmrootConfigPath, traCmd)) +os.system((". {} -a ; {}").format(cbmrootConfigPath, digiCmd)) +os.system((". {} -a ; {}").format(cbmrootConfigPath, recoCmd)) +os.system((". {} -a ; {}").format(cbmrootConfigPath, qaCmd)) diff --git a/macro/rich/run/run_qa.C b/macro/rich/run/run_qa.C index 14b95da9e0..a146ff4a92 100644 --- a/macro/rich/run/run_qa.C +++ b/macro/rich/run/run_qa.C @@ -1,26 +1,62 @@ -/* Copyright (C) 2019-2020 UGiessen, JINR-LIT +/* Copyright (C) 2019-2023 UGiessen, JINR-LIT SPDX-License-Identifier: GPL-3.0-only - Authors: Semen Lebedev [committer] */ - -void run_qa(const string& traFile = "/Users/slebedev/Development/cbm/data/sim/rich/reco/tra.0.root", - const string& parFile = "/Users/slebedev/Development/cbm/data/sim/rich/reco/par.0.root", - const string& digiFile = "/Users/slebedev/Development/cbm/data/sim/rich/reco/digi.0.root", - const string& recoFile = "/Users/slebedev/Development/cbm/data/sim/rich/reco/reco.0.root", - const string& qaFile = "/Users/slebedev/Development/cbm/data/sim/rich/reco/qa.0.root", - const string& geoSetup = "sis100_electron", const string& resultDir = "results_recoqa_apr21plus/", - int nEvents = 100) + Authors: Semen Lebedev [committer], Martin Beyer */ + +#if !defined(__CLING__) +#include "CbmL1RichRingQa.h" +#include "CbmLitTrackingQa.h" +#include "CbmMCDataManager.h" +#include "CbmMatchRecoToMC.h" +#include "CbmRichGeoTest.h" +#include "CbmRichMatchRings.h" +#include "CbmRichRecoQa.h" +// #include "CbmRichRingFinderQa.h" +// #include "CbmRichTrackProjQa.h" +#include "CbmSetup.h" + +#include <FairFileSource.h> +#include <FairLogger.h> +#include <FairMonitor.h> +#include <FairParAsciiFileIo.h> +#include <FairParRootFileIo.h> +#include <FairRootFileSink.h> +#include <FairRunAna.h> +#include <FairRuntimeDb.h> + +#include <TGeoManager.h> +#include <TObjString.h> +#include <TROOT.h> +#include <TStopwatch.h> +#include <TTree.h> + +#include <iostream> +#include <string> +#endif + +void run_qa(TString traFile = "data/test.tra.root", TString parFile = "data/test.par.root", + TString digiFile = "data/test.digi.root", TString recoFile = "data/test.reco.root", + TString qaFile = "data/test.qa.root", Int_t nofTimeSlices = -1, TString geoSetup = "sis100_electron", + TString resultDir = "", Bool_t monitor = true) { - TTree::SetMaxTreeSize(90000000000); + + remove(qaFile); + + // ----- Logger settings ---------------------------------------------- FairLogger::GetLogger()->SetLogScreenLevel("INFO"); FairLogger::GetLogger()->SetLogVerbosityLevel("LOW"); - remove(qaFile.c_str()); + // ------------------------------------------------------------------------ + // ----- Environment -------------------------------------------------- TString srcDir = gSystem->Getenv("VMCWORKDIR"); + // ------------------------------------------------------------------------ + // ----- Load the geometry setup ------------------------------------- CbmSetup* geo = CbmSetup::Instance(); - geo->LoadSetup(geoSetup.c_str()); + geo->LoadSetup(geoSetup); + // ------------------------------------------------------------------------ + // ----- Parameter files as input to the runtime database ------------- TList* parFileList = new TList(); TString geoTag; @@ -39,93 +75,124 @@ void run_qa(const string& traFile = "/Users/slebedev/Development/cbm/data/sim/r TObjString* tofBdfFile = new TObjString(srcDir + "/parameters/tof/tof_" + geoTag + ".digibdf.par"); parFileList->Add(tofBdfFile); } + // ------------------------------------------------------------------------ + // ----- Timer -------------------------------------------------------- TStopwatch timer; timer.Start(); - gDebug = 0; + // ------------------------------------------------------------------------ + // ----- FairRunAna --------------------------------------------------- + FairFileSource* inputSource = new FairFileSource(digiFile); + inputSource->AddFriend(traFile); + inputSource->AddFriend(recoFile); - FairRunAna* run = new FairRunAna(); - FairFileSource* inputSource = new FairFileSource(digiFile.c_str()); - inputSource->AddFriend(traFile.c_str()); - inputSource->AddFriend(recoFile.c_str()); + FairRunAna* run = new FairRunAna(); run->SetSource(inputSource); - run->SetOutputFile(qaFile.c_str()); - run->SetGenerateRunInfo(kTRUE); + run->SetGenerateRunInfo(kFALSE); - CbmMCDataManager* mcManager = new CbmMCDataManager("MCManager", 1); - mcManager->AddFile(traFile.c_str()); + FairRootFileSink* sink = new FairRootFileSink(qaFile); + run->SetSink(sink); + + FairMonitor::GetMonitor()->EnableMonitor(monitor); + // ------------------------------------------------------------------------ + + // ----- MCDataManager ------------------------------------------------ + CbmMCDataManager* mcManager = new CbmMCDataManager("MCManager", 0); + mcManager->AddFile(traFile); run->AddTask(mcManager); + // ------------------------------------------------------------------------ - // // RICH reco QA - CbmRichRecoQa* richRecoQa = new CbmRichRecoQa(); - richRecoQa->SetOutputDir(resultDir); - run->AddTask(richRecoQa); + // ----- Match reco to MC --------------------------------------------- + // CbmMatchRecoToMC* match = new CbmMatchRecoToMC(); + // run->AddTask(match); + // ------------------------------------------------------------------------ - // // Reconstruction Qa - CbmLitTrackingQa* trackingQa = new CbmLitTrackingQa(); - trackingQa->SetMinNofPointsSts(4); - trackingQa->SetUseConsecutivePointsInSts(true); - trackingQa->SetMinNofPointsTrd(2); - trackingQa->SetMinNofPointsMuch(10); - trackingQa->SetMinNofPointsTof(1); - trackingQa->SetQuota(0.7); - trackingQa->SetMinNofHitsTrd(2); - trackingQa->SetMinNofHitsMuch(10); - trackingQa->SetVerbose(0); - trackingQa->SetMinNofHitsRich(7); - trackingQa->SetQuotaRich(0.6); - trackingQa->SetOutputDir(resultDir); - trackingQa->SetPRange(12, 0., 6.); - // trackingQa->SetTrdAnnCut(trdAnnCut); - std::vector<std::string> trackCat, richCat; - trackCat.push_back("All"); - trackCat.push_back("Electron"); - richCat.push_back("Electron"); - richCat.push_back("ElectronReference"); - trackingQa->SetTrackCategories(trackCat); - trackingQa->SetRingCategories(richCat); - run->AddTask(trackingQa); - - CbmLitFitQa* fitQa = new CbmLitFitQa(); - fitQa->SetMvdMinNofHits(0); - fitQa->SetStsMinNofHits(6); - fitQa->SetMuchMinNofHits(10); - fitQa->SetTrdMinNofHits(2); - fitQa->SetOutputDir(resultDir); - //run->AddTask(fitQa); - - CbmLitClusteringQa* clusteringQa = new CbmLitClusteringQa(); - clusteringQa->SetOutputDir(resultDir); - //run->AddTask(clusteringQa); - - CbmLitTofQa* tofQa = new CbmLitTofQa(); - tofQa->SetOutputDir(std::string(resultDir)); - //run->AddTask(tofQa); + // ----- RichRingFinderQa --------------------------------------------- + // CbmRichRingFinderQa* finder = new CbmRichRingFinderQa(); + // run->AddTask(finder); + // ------------------------------------------------------------------------ + // ----- CbmRichTrackProjQa ------------------------------------------- + // CbmRichTrackProjQa* proj = new CbmRichTrackProjQa(); + // run->AddTask(proj); + // ------------------------------------------------------------------------ + // ----- RICH Reco QA ------------------------------------------------- + CbmRichRecoQa* richRecoQa = new CbmRichRecoQa(); + richRecoQa->SetOutputDir(resultDir.Data()); + run->AddTask(richRecoQa); + // ------------------------------------------------------------------------ + + // ----- RICH GeoTest ------------------------------------------------- + // CbmRichGeoTest* geoTest = new CbmRichGeoTest(); + // geoTest->SetOutputDir(resultDir.Data()); + // run->AddTask(geoTest); + // ------------------------------------------------------------------------ + + // ----- Lit Tracking Qa ---------------------------------------------- + // CbmLitTrackingQa* trackingQa = new CbmLitTrackingQa(); + // trackingQa->SetMinNofPointsSts(4); + // trackingQa->SetUseConsecutivePointsInSts(true); + // trackingQa->SetMinNofPointsTrd(2); + // trackingQa->SetMinNofPointsMuch(10); + // trackingQa->SetMinNofPointsTof(1); + // trackingQa->SetQuota(0.7); + // trackingQa->SetMinNofHitsTrd(2); + // trackingQa->SetMinNofHitsMuch(10); + // trackingQa->SetVerbose(0); + // trackingQa->SetMinNofHitsRich(7); + // trackingQa->SetQuotaRich(0.6); + // trackingQa->SetPRange(12, 0., 6.); + // // trackingQa->SetTrdAnnCut(trdAnnCut); + // std::vector<std::string> trackCat, richCat; + // trackCat.push_back("All"); + // trackCat.push_back("Electron"); + // richCat.push_back("Electron"); + // richCat.push_back("ElectronReference"); + // trackingQa->SetTrackCategories(trackCat); + // trackingQa->SetRingCategories(richCat); + // trackingQa->SetOutputDir(resultDir.Data()); + // run->AddTask(trackingQa); + // ------------------------------------------------------------------------ + + // ----- Parameter database -------------------------------------------- FairRuntimeDb* rtdb = run->GetRuntimeDb(); FairParRootFileIo* parIo1 = new FairParRootFileIo(); FairParAsciiFileIo* parIo2 = new FairParAsciiFileIo(); - parIo1->open(parFile.c_str(), "UPDATE"); + parIo1->open(parFile, "UPDATE"); rtdb->setFirstInput(parIo1); if (!parFileList->IsEmpty()) { parIo2->open(parFileList, "in"); rtdb->setSecondInput(parIo2); } - run->Init(); + // ------------------------------------------------------------------------ + // ----- Run initialisation ------------------------------------------- + run->Init(); rtdb->setOutput(parIo1); rtdb->saveOutput(); rtdb->print(); + // ------------------------------------------------------------------------ - run->Run(0, nEvents); + // ----- Start run ---------------------------------------------------- + run->Run(0, nofTimeSlices); + // ------------------------------------------------------------------------ + // ----- Finish ------------------------------------------------------- timer.Stop(); - std::cout << std::endl << std::endl; + std::cout << std::endl; + FairMonitor::GetMonitor()->Print(); + std::cout << std::endl; std::cout << "Macro finished succesfully." << std::endl; std::cout << "Output file is " << qaFile << std::endl; std::cout << "Parameter file is " << parFile << std::endl; std::cout << "Real time " << timer.RealTime() << " s, CPU time " << timer.CpuTime() << " s" << std::endl; + if (gROOT->GetVersionInt() >= 60602) { + gGeoManager->GetListOfVolumes()->Delete(); + gGeoManager->GetListOfShapes()->Delete(); + delete gGeoManager; + } std::cout << std::endl << "Test passed" << std::endl << "All ok" << std::endl; + // ------------------------------------------------------------------------ } diff --git a/macro/rich/run/run_reco.C b/macro/rich/run/run_reco.C index e012219487..ffef6ba360 100644 --- a/macro/rich/run/run_reco.C +++ b/macro/rich/run/run_reco.C @@ -1,28 +1,81 @@ -/* Copyright (C) 2009-2020 UGiessen, JINR-LIT +/* Copyright (C) 2009-2023 UGiessen, JINR-LIT SPDX-License-Identifier: GPL-3.0-only - Authors: Semen Lebedev [committer], Semen Lebedev [committer], Andrey Lebedev */ - -void run_reco(const string& traFile = "/Users/slebedev/Development/cbm/data/sim/rich/reco/tra.0.root", - const string& parFile = "/Users/slebedev/Development/cbm/data/sim/rich/reco/par.0.root", - const string& digiFile = "/Users/slebedev/Development/cbm/data/sim/rich/reco/digi.0.root", - const string& recoFile = "/Users/slebedev/Development/cbm/data/sim/rich/reco/reco.0.root", - const string& geoSetup = "sis100_electron", const string& resultDir = "results_recoqa_newqa/", - int nofTimeSlices = 100) + Authors: Semen Lebedev [committer], Andrey Lebedev, Martin Beyer */ + +#if !defined(__CLING__) +#include "CbmBuildEventsFromTracksReal.h" +#include "CbmBuildEventsIdeal.h" +#include "CbmBuildEventsQa.h" +#include "CbmFindPrimaryVertex.h" +#include "CbmKF.h" +#include "CbmL1.h" +#include "CbmL1StsTrackFinder.h" +#include "CbmLitFindGlobalTracks.h" +#include "CbmMCDataManager.h" +#include "CbmMatchRecoToMC.h" +#include "CbmMvdClusterfinder.h" +#include "CbmMvdHitfinder.h" +#include "CbmPVFinderKF.h" +#include "CbmPrimaryVertexFinder.h" +#include "CbmRecoSts.h" +#include "CbmRecoTzero.h" +#include "CbmRichHitProducer.h" +#include "CbmRichReconstruction.h" +#include "CbmSetup.h" +#include "CbmStsFindTracks.h" +#include "CbmStsFindTracksEvents.h" +#include "CbmStsTrackFinder.h" +#include "CbmTaskBuildRawEvents.h" +#include "CbmTofSimpClusterizer.h" +#include "CbmTrackingDetectorInterfaceInit.h" +#include "CbmTrdClusterFinder.h" +#include "CbmTrdHitProducer.h" +#include "CbmTrdSetTracksPidLike.h" + +#include <FairFileSource.h> +#include <FairLogger.h> +#include <FairMonitor.h> +#include <FairParAsciiFileIo.h> +#include <FairParRootFileIo.h> +#include <FairRootFileSink.h> +#include <FairRunAna.h> + +#include <TGeoManager.h> +#include <TObjString.h> +#include <TROOT.h> +#include <TStopwatch.h> +#include <TTree.h> + +#include <iostream> +#endif + +void run_reco(TString traFile = "data/test.tra.root", TString parFile = "data/test.par.root", + TString digiFile = "data/test.digi.root", TString recoFile = "data/test.reco.root", + Int_t nofTimeSlices = -1, TString geoSetup = "sis100_electron", TString sEvBuildRaw = "Ideal", + Bool_t useMC = true, Bool_t monitor = true) { - TTree::SetMaxTreeSize(90000000000); + + remove(recoFile); + + // --- Logger settings ---------------------------------------------------- FairLogger::GetLogger()->SetLogScreenLevel("INFO"); FairLogger::GetLogger()->SetLogVerbosityLevel("LOW"); - remove(recoFile.c_str()); + // ------------------------------------------------------------------------ + // ----- Environment -------------------------------------------------- + TString myName = "run_reco"; TString srcDir = gSystem->Getenv("VMCWORKDIR"); + // ------------------------------------------------------------------------ + // ----- Load the geometry setup -------------------------------------- + std::cout << std::endl; + std::cout << "-I- " << myName << ": Loading setup " << geoSetup << std::endl; CbmSetup* geo = CbmSetup::Instance(); - geo->LoadSetup(geoSetup.c_str()); - - TString sEvBuildRaw = ""; - Bool_t useMC = true; + geo->LoadSetup(geoSetup); + // ------------------------------------------------------------------------ + // ----- Some global switches ----------------------------------------- Bool_t eventBased = (sEvBuildRaw != ""); Bool_t useMvd = geo->IsActive(ECbmModuleId::kMvd); Bool_t useSts = geo->IsActive(ECbmModuleId::kSts); @@ -31,7 +84,12 @@ void run_reco(const string& traFile = "/Users/slebedev/Development/cbm/data/sim Bool_t useTrd = geo->IsActive(ECbmModuleId::kTrd); Bool_t useTof = geo->IsActive(ECbmModuleId::kTof); Bool_t usePsd = geo->IsActive(ECbmModuleId::kPsd); + // Bool_t useFsd = geo->IsActive(ECbmModuleId::kFsd); + // ------------------------------------------------------------------------ + // ----- Parameter files as input to the runtime database ------------- + std::cout << std::endl; + std::cout << "-I- " << myName << ": Defining parameter files " << std::endl; TList* parFileList = new TList(); TString geoTag; @@ -42,6 +100,7 @@ void run_reco(const string& traFile = "/Users/slebedev/Development/cbm/data/sim for (Int_t i(0); i < 4; i++) { trdParFile = new TObjString(srcDir + "/parameters/trd/trd_" + geoTag + "." + npar[i] + ".par"); parFileList->Add(trdParFile); + std::cout << "-I- " << myName << ": Using parameter file " << trdParFile->GetString() << std::endl; } } @@ -49,32 +108,45 @@ void run_reco(const string& traFile = "/Users/slebedev/Development/cbm/data/sim if (CbmSetup::Instance()->GetGeoTag(ECbmModuleId::kTof, geoTag)) { TObjString* tofBdfFile = new TObjString(srcDir + "/parameters/tof/tof_" + geoTag + ".digibdf.par"); parFileList->Add(tofBdfFile); + std::cout << "-I- " << myName << ": Using parameter file " << tofBdfFile->GetString() << std::endl; } + // ------------------------------------------------------------------------ + // ----- Timer -------------------------------------------------------- TStopwatch timer; timer.Start(); + // ------------------------------------------------------------------------ + + // ----- FairRunAna --------------------------------------------------- + FairFileSource* inputSource = new FairFileSource(digiFile); + if (useMC) { inputSource->AddFriend(traFile); } - FairRunAna* run = new FairRunAna(); - FairFileSource* inputSource = new FairFileSource(digiFile.c_str()); - if (useMC) { inputSource->AddFriend(traFile.c_str()); } + FairRunAna* run = new FairRunAna(); run->SetSource(inputSource); - run->SetOutputFile(recoFile.c_str()); - run->SetGenerateRunInfo(true); + run->SetGenerateRunInfo(kFALSE); + + FairRootFileSink* sink = new FairRootFileSink(recoFile); + run->SetSink(sink); + + FairMonitor::GetMonitor()->EnableMonitor(monitor); + // ------------------------------------------------------------------------ + // ----- MCDataManager ------------------------------------------------- if (useMC) { CbmMCDataManager* mcManager = new CbmMCDataManager("MCDataManager", 0); - mcManager->AddFile(traFile.c_str()); + mcManager->AddFile(traFile); run->AddTask(mcManager); + std::cout << "-I- " << myName << ": Added task " << mcManager->GetName() << std::endl; } + // ------------------------------------------------------------------------ + // ----- Raw event building from digis -------------------------------- if (eventBased) { if (sEvBuildRaw.EqualTo("Ideal", TString::ECaseCompare::kIgnoreCase)) { FairTask* evBuildRaw = new CbmBuildEventsIdeal(); run->AddTask(evBuildRaw); - std::cout << "-I- " - << ": Added task " << evBuildRaw->GetName() << std::endl; - eventBased = kTRUE; - } + std::cout << "-I- " << myName << ": Added task " << evBuildRaw->GetName() << std::endl; + } //? Ideal raw event building else if (sEvBuildRaw.EqualTo("Real", TString::ECaseCompare::kIgnoreCase)) { CbmTaskBuildRawEvents* evBuildRaw = new CbmTaskBuildRawEvents(); @@ -85,21 +157,27 @@ void run_reco(const string& traFile = "/Users/slebedev/Development/cbm/data/sim if (!useRich) evBuildRaw->RemoveDetector(kRawEventBuilderDetRich); if (!useMuch) evBuildRaw->RemoveDetector(kRawEventBuilderDetMuch); if (!usePsd) evBuildRaw->RemoveDetector(kRawEventBuilderDetPsd); + // if (!useFsd) evBuildRaw->RemoveDetector(kRawEventBuilderDetFsd); if (!useTof) evBuildRaw->RemoveDetector(kRawEventBuilderDetTof); if (!useTrd) evBuildRaw->RemoveDetector(kRawEventBuilderDetTrd); if (!useSts) { - std::cerr << "-E- " - << ": Sts must be present for raw event " + std::cerr << "-E- " << myName << ": Sts must be present for raw event " << "building using ``Real2019'' option. Terminating macro." << std::endl; return; } // Set STS as reference detector evBuildRaw->SetReferenceDetector(kRawEventBuilderDetSts); - evBuildRaw->SetTsParameters(0.0, 1.e7, 0.0); - // Make Bmon (previous reference detector) a selected detector (with default parameters) evBuildRaw->AddDetector(kRawEventBuilderDetT0); + // Use sliding window seed builder with STS + //evBuildRaw->SetReferenceDetector(kRawEventBuilderDetUndef); + //evBuildRaw->AddSeedTimeFillerToList(kRawEventBuilderDetSts); + //evBuildRaw->SetSlidingWindowSeedFinder(1000, 500, 500); + //evBuildRaw->SetSeedFinderQa(true); // optional QA information for seed finder + + evBuildRaw->SetTsParameters(0.0, 1.e7, 0.0); + // Use CbmMuchDigi instead of CbmMuchBeamtimeDigi evBuildRaw->ChangeMuchBeamtimeDigiFlag(kFALSE); @@ -108,55 +186,55 @@ void run_reco(const string& traFile = "/Users/slebedev/Development/cbm/data/sim evBuildRaw->SetTriggerWindow(ECbmModuleId::kSts, -500, 500); run->AddTask(evBuildRaw); - std::cout << "-I- " - << ": Added task " << evBuildRaw->GetName() << std::endl; - eventBased = kTRUE; + std::cout << "-I- " << myName << ": Added task " << evBuildRaw->GetName() << std::endl; } //? Real raw event building else { - std::cerr << "-E- " - << ": Unknown option " << sEvBuildRaw << " for raw event building! Terminating macro execution." - << std::endl; + std::cerr << "-E- " << myName << ": Unknown option " << sEvBuildRaw + << " for raw event building! Terminating macro execution." << std::endl; return; } - } - + } //? event-based reco + // ------------------------------------------------------------------------ - if (eventBased && useMC) { - CbmBuildEventsQa* evBuildQA = new CbmBuildEventsQa(); - run->AddTask(evBuildQA); - } + // ----- QA for raw event builder ------------------------------------- + // if (eventBased && useMC) { + // CbmBuildEventsQa* evBuildQA = new CbmBuildEventsQa(); + // run->AddTask(evBuildQA); + // std::cout << "-I- " << myName << ": Added task " << evBuildQA->GetName() << std::endl; + // } + // ------------------------------------------------------------------------ + // ----- Local reconstruction in MVD ---------------------------------- if (useMvd) { - CbmMvdClusterfinder* mvdCluster = new CbmMvdClusterfinder("MVD Cluster Finder", 0, 0); run->AddTask(mvdCluster); - std::cout << "-I- : Added task " << mvdCluster->GetName() << std::endl; + std::cout << "-I- " << myName << ": Added task " << mvdCluster->GetName() << std::endl; CbmMvdHitfinder* mvdHit = new CbmMvdHitfinder("MVD Hit Finder", 0, 0); mvdHit->UseClusterfinder(kTRUE); run->AddTask(mvdHit); - std::cout << "-I- " - << ": Added task " << mvdHit->GetName() << std::endl; + std::cout << "-I- " << myName << ": Added task " << mvdHit->GetName() << std::endl; } + // ------------------------------------------------------------------------ - + // ----- Local reconstruction in STS ---------------------------------- if (useSts) { - CbmRecoSts* stsReco = new CbmRecoSts(kCbmRecoTimeslice); - if (eventBased) stsReco->SetMode(kCbmRecoEvent); + CbmRecoSts* stsReco = new CbmRecoSts(ECbmRecoMode::Timeslice); + if (eventBased) stsReco->SetMode(ECbmRecoMode::EventByEvent); run->AddTask(stsReco); - std::cout << "-I- " - << ": Added task " << stsReco->GetName() << std::endl; + std::cout << "-I- " << myName << ": Added task " << stsReco->GetName() << std::endl; } + // ------------------------------------------------------------------------ - + // ----- Local reconstruction in RICH --------------------------------- if (useRich) { CbmRichHitProducer* richHitProd = new CbmRichHitProducer(); run->AddTask(richHitProd); - std::cout << "-I- " - << ": Added task " << richHitProd->GetName() << std::endl; + std::cout << "-I- " << myName << ": Added task " << richHitProd->GetName() << std::endl; } + // ------------------------------------------------------------------------ - + // ----- Local reconstruction in TRD ---------------------------------- if (useTrd) { Double_t triggerThreshold = 0.5e-6; // SIS100 @@ -173,124 +251,128 @@ void run_reco(const string& traFile = "/Users/slebedev/Development/cbm/data/sim // trdCluster->SetUseOnlyEventDigis(kFALSE); run->AddTask(trdCluster); - std::cout << "-I- " - << ": Added task " << trdCluster->GetName() << std::endl; + std::cout << "-I- " << myName << ": Added task " << trdCluster->GetName() << std::endl; CbmTrdHitProducer* trdHit = new CbmTrdHitProducer(); run->AddTask(trdHit); - std::cout << "-I- " - << ": Added task " << trdHit->GetName() << std::endl; + std::cout << "-I- " << myName << ": Added task " << trdHit->GetName() << std::endl; } + // ------------------------------------------------------------------------ + // ----- Local reconstruction in TOF ---------------------------------- if (useTof) { CbmTofSimpClusterizer* tofCluster = new CbmTofSimpClusterizer("TofSimpClusterizer", 0); tofCluster->SetOutputBranchPersistent("TofHit", kTRUE); tofCluster->SetOutputBranchPersistent("TofDigiMatch", kTRUE); run->AddTask(tofCluster); - std::cout << "-I- " - << ": Added task " << tofCluster->GetName() << std::endl; + std::cout << "-I- " << myName << ": Added task " << tofCluster->GetName() << std::endl; } // ------------------------------------------------------------------------ + // TODO: add Fsd + + // ----- Track finding in STS (+ MVD) -------------------------------- if (useMvd || useSts) { - run->AddTask(new CbmTrackingDetectorInterfaceInit()); - CbmKF* kalman = new CbmKF(); + auto pDetIF = new CbmTrackingDetectorInterfaceInit(); + run->AddTask(pDetIF); // Geometry interface initializer for tracker + std::cout << "-I- " << myName << ": Added task " << pDetIF->GetName() << std::endl; + + // Kalman filter + auto kalman = new CbmKF(); run->AddTask(kalman); - CbmL1* l1 = new CbmL1("L1", 0); - - // --- Material budget file names - TString mvdGeoTag; - if (geo->GetGeoTag(ECbmModuleId::kMvd, mvdGeoTag)) { - TString parFile = gSystem->Getenv("VMCWORKDIR"); - parFile += "/parameters/mvd/mvd_matbudget_" + mvdGeoTag + ".root"; - std::cout << "Using material budget file " << parFile << std::endl; - l1->SetMvdMaterialBudgetFileName(parFile.Data()); - } - TString stsGeoTag; - if (geo->GetGeoTag(ECbmModuleId::kSts, stsGeoTag)) { - TString parFile = gSystem->Getenv("VMCWORKDIR"); - parFile += "/parameters/sts/sts_matbudget_" + stsGeoTag + ".root"; - std::cout << "Using material budget file " << parFile << std::endl; - l1->SetStsMaterialBudgetFileName(parFile.Data()); - } + std::cout << "-I- " << myName << ": Added task " << kalman->GetName() << std::endl; + + // L1 tracking + auto l1 = new CbmL1("CA"); + // auto l1 = (useMC) ? new CbmL1("CA", 2, 3) : new CbmL1("CA"); + + // L1 configuration file (optional) + // At the moment, the YAML configuration file defines different parameters for a sequence of track finder + // iterations. The same file should be used in ca::tools::WindowFinder class for hit search window estimation + //l1->SetInputConfigName(TString(gSystem->Getenv("VMCWORKDIR")) + "/reco/L1/L1Algo/L1ConfigExample.yaml"); + run->AddTask(l1); - std::cout << "-I- " - << ": Added task " << l1->GetName() << std::endl; + std::cout << "-I- " << myName << ": Added task " << l1->GetName() << std::endl; CbmStsTrackFinder* stsTrackFinder = new CbmL1StsTrackFinder(); if (eventBased) { FairTask* stsFindTracks = new CbmStsFindTracksEvents(stsTrackFinder, useMvd); run->AddTask(stsFindTracks); - std::cout << "-I- " - << ": Added task " << stsFindTracks->GetName() << std::endl; + std::cout << "-I- " << myName << ": Added task " << stsFindTracks->GetName() << std::endl; } else { FairTask* stsFindTracks = new CbmStsFindTracks(0, stsTrackFinder, useMvd); run->AddTask(stsFindTracks); - std::cout << "-I- " - << ": Added task " << stsFindTracks->GetName() << std::endl; + std::cout << "-I- " << myName << ": Added task " << stsFindTracks->GetName() << std::endl; } } - // ------------------------------------------------------------------------ - - - // --- Global track finding ----------------------------------------- - CbmLitFindGlobalTracks* finder = new CbmLitFindGlobalTracks(); - finder->SetTrackingType("branch"); - finder->SetMergerType("nearest_hit"); - run->AddTask(finder); - std::cout << "-I- : Added task " << finder->GetName() << std::endl; // ---------------------------------------------------------------------- - // --- Particle Id in TRD ----------------------------------------- - if (useTrd) { - CbmTrdSetTracksPidLike* trdLI = new CbmTrdSetTracksPidLike("TRDLikelihood", "TRDLikelihood"); - trdLI->SetUseMCInfo(kTRUE); - trdLI->SetUseMomDependence(kTRUE); - run->AddTask(trdLI); - std::cout << "-I- : Added task " << trdLI->GetName() << std::endl; - } - // ------------------------------------------------------------------------ - - // ----- RICH reconstruction ---------------------------------------- - if (useRich) { - CbmRichReconstruction* richReco = new CbmRichReconstruction(); - run->AddTask(richReco); - std::cout << "-I- : Added task " << richReco->GetName() << std::endl; - } - // ---------------------------------------------------------------------- - - // ==== From here on, the time-based and the event-based reconstruction - // ==== chains differ, since time-based version of primary vertex finding - // ==== and global tracking are not yet available. For time-based - // ==== reconstruction, a track-based event finder is used; no global - // ==== tracks are produced. - if (eventBased) { - - //----- Primary vertex finding ------------------------------------- + // ----- Primary vertex finding ------------------------------------- CbmPrimaryVertexFinder* pvFinder = new CbmPVFinderKF(); CbmFindPrimaryVertex* findVertex = new CbmFindPrimaryVertex(pvFinder); run->AddTask(findVertex); - std::cout << "-I- " - << ": Added task " << findVertex->GetName() << std::endl; + std::cout << "-I- " << myName << ": Added task " << findVertex->GetName() << std::endl; // ---------------------------------------------------------------------- + // --- Global track finding ----------------------------------------- + CbmLitFindGlobalTracks* finder = new CbmLitFindGlobalTracks(); + finder->SetTrackingType("branch"); + finder->SetMergerType("nearest_hit"); + run->AddTask(finder); + std::cout << "-I- " << myName << ": Added task " << finder->GetName() << std::endl; + // ---------------------------------------------------------------------- - } //? event-based reco + // --- Particle Id in TRD ------------------------------------------- + if (useTrd) { + CbmTrdSetTracksPidLike* trdLI = new CbmTrdSetTracksPidLike("TRDLikelihood", "TRDLikelihood"); + trdLI->SetUseMCInfo(kTRUE); + trdLI->SetUseMomDependence(kTRUE); + run->AddTask(trdLI); + std::cout << "-I- " << myName << ": Added task " << trdLI->GetName() << std::endl; + } + // ---------------------------------------------------------------------- - else { + // ----- RICH reconstruction ---------------------------------------- + if (useRich) { + CbmRichReconstruction* richReco = new CbmRichReconstruction(); + run->AddTask(richReco); + std::cout << "-I- " << myName << ": Added task " << richReco->GetName() << std::endl; + } + // ---------------------------------------------------------------------- + + // ----- T0 reconstruction -------------------------------------------- + CbmRecoTzero* recoT0 = new CbmRecoTzero(); + run->AddTask(recoT0); + std::cout << "-I- " << myName << ": Added task " << recoT0->GetName() << std::endl; + // ---------------------------------------------------------------------- + + // ----- Match reco to MC ------------------------------------------- + if (useMC) { + CbmMatchRecoToMC* match = new CbmMatchRecoToMC(); + std::cout << "-I- " << myName << ": Added task " << match->GetName() << std::endl; + run->AddTask(match); + } + // ---------------------------------------------------------------------- + } //? event-based reco + else { + if (useRich) { // Place it here to test rich ring reco on timeslices + CbmRichReconstruction* richReco = new CbmRichReconstruction(); + richReco->SetRunExtrapolation(kFALSE); + richReco->SetRunProjection(kFALSE); + richReco->SetRunTrackAssign(kFALSE); + run->AddTask(richReco); + std::cout << "-I- " << myName << ": Added task " << richReco->GetName() << std::endl; + } // -----Â Event building from STS tracks ----------------------------- run->AddTask(new CbmBuildEventsFromTracksReal()); // ---------------------------------------------------------------------- - + //FIXME: if placing rich reco here: + //FIXME: hits are found in events, but no rings are reconstructed } //? time-based reco - CbmMatchRecoToMC* match = new CbmMatchRecoToMC(); - run->AddTask(match); - - // ----- Parameter database -------------------------------------------- std::cout << std::endl << std::endl; std::cout << "-I- " @@ -298,7 +380,7 @@ void run_reco(const string& traFile = "/Users/slebedev/Development/cbm/data/sim FairRuntimeDb* rtdb = run->GetRuntimeDb(); FairParRootFileIo* parIo1 = new FairParRootFileIo(); FairParAsciiFileIo* parIo2 = new FairParAsciiFileIo(); - parIo1->open(parFile.c_str(), "UPDATE"); + parIo1->open(parFile.Data(), "UPDATE"); rtdb->setFirstInput(parIo1); if (!parFileList->IsEmpty()) { parIo2->open(parFileList, "in"); @@ -306,7 +388,6 @@ void run_reco(const string& traFile = "/Users/slebedev/Development/cbm/data/sim } // ------------------------------------------------------------------------ - // ----- Run initialisation ------------------------------------------- std::cout << std::endl; std::cout << "-I- " @@ -317,7 +398,6 @@ void run_reco(const string& traFile = "/Users/slebedev/Development/cbm/data/sim rtdb->print(); // ------------------------------------------------------------------------ - // ----- Register light ions (d, t, He3, He4) ------------------------- std::cout << std::endl; TString registerLightIonsMacro = gSystem->Getenv("VMCWORKDIR"); @@ -334,11 +414,20 @@ void run_reco(const string& traFile = "/Users/slebedev/Development/cbm/data/sim run->Run(0, nofTimeSlices); // ------------------------------------------------------------------------ + // ----- Finish ------------------------------------------------------- timer.Stop(); - std::cout << std::endl << std::endl; + std::cout << std::endl; + FairMonitor::GetMonitor()->Print(); + std::cout << std::endl; std::cout << "Macro finished succesfully." << std::endl; std::cout << "Output file is " << recoFile << std::endl; std::cout << "Parameter file is " << parFile << std::endl; std::cout << "Real time " << timer.RealTime() << " s, CPU time " << timer.CpuTime() << " s" << std::endl; + if (gROOT->GetVersionInt() >= 60602) { + gGeoManager->GetListOfVolumes()->Delete(); + gGeoManager->GetListOfShapes()->Delete(); + delete gGeoManager; + } std::cout << std::endl << "Test passed" << std::endl << "All ok" << std::endl; + // ------------------------------------------------------------------------ } diff --git a/macro/rich/run/run_transport.C b/macro/rich/run/run_transport.C index af74a747e2..f8753bb3d9 100644 --- a/macro/rich/run/run_transport.C +++ b/macro/rich/run/run_transport.C @@ -1,74 +1,110 @@ -/* Copyright (C) 2021 UGiessen, JINR-LIT +/* Copyright (C) 2021-2023 UGiessen, JINR-LIT SPDX-License-Identifier: GPL-3.0-only - Authors: Semen Lebedev [committer] */ + Authors: Semen Lebedev [committer], Martin Beyer */ -void run_transport(const string& urqmdFile = "/Users/slebedev/Development/cbm/data/urqmd/auau/8gev/centr/" - "urqmd.auau.8gev.centr.00001.root", // if "", no urqmd - const string& traFile = "/Users/slebedev/Development/cbm/data/sim/rich/reco/tra.0.root", - const string& parFile = "/Users/slebedev/Development/cbm/data/sim/rich/reco/par.0.root", - const string& geoFile = "/Users/slebedev/Development/cbm/data/sim/rich/reco/geo.0.root", - int nofElectrons = 5, // number of e- to be generated, if <=0 no e- are embedded into event - int nofPositrons = 5, // number of e+ to be generated, if <=0 no e+ are embedded into event - const string& plutoFile = "", // if "", no pluto particles are embedded into event - const string& geoSetup = "sis100_electron", int nEvents = 2, double targetZ = -44.) -{ +#if !defined(__CLING__) +#include "CbmRich.h" +#include "CbmTransport.h" + +#include <FairBoxGenerator.h> +#include <FairLogger.h> +#include <FairRunSim.h> + +#include <TStopwatch.h> +#include <TTree.h> +#include <iostream> +#endif + +void run_transport(TString urqmdFile = "default10gev", // if "", no urqmd; "default10gev" or "default25gev" + TString traFile = "data/test.tra.root", TString parFile = "data/test.par.root", + TString geoFile = "data/test.geo.root", Int_t nEvents = 1, + Int_t nElectrons = 5, // number of embedded e- into each event, generated by the BoxGenerator + Int_t nPositrons = 5, // number of embedded e+ into each event, generated by the BoxGenerator + TString plutoFile = "", // if "", no pluto particles are embedded into event + TString geoSetup = "sis100_electron", Int_t engine = 0, // geant3: 0, geant4: 1 + Int_t randomSeed = 0) +{ TTree::SetMaxTreeSize(90000000000); + + remove(traFile); + remove(parFile); + remove(geoFile); + + Double_t targetZ = -44.0; + + // ----- Logger settings ---------------------------------------------- FairLogger::GetLogger()->SetLogScreenLevel("INFO"); FairLogger::GetLogger()->SetLogVerbosityLevel("LOW"); + // ------------------------------------------------------------------------ - remove(parFile.c_str()); - remove(traFile.c_str()); - remove(geoFile.c_str()); + // ----- Environment -------------------------------------------------- + TString srcDir = gSystem->Getenv("VMCWORKDIR"); + TString defaultInput10gev = srcDir + "/input/urqmd.auau.10gev.centr.root"; + TString defaultInput25gev = srcDir + "/input/urqmd.auau.25gev.centr.root"; + // ------------------------------------------------------------------------ + // ----- Timer -------------------------------------------------------- TStopwatch timer; timer.Start(); + // ------------------------------------------------------------------------ + // --- Transport run ---------------------------------------------------- CbmTransport run; + run.SetEngine(static_cast<ECbmEngine>(engine)); - if (urqmdFile.length() > 0) { run.AddInput(urqmdFile.c_str()); } - - if (nofElectrons > 0) { - FairBoxGenerator* boxGen1 = new FairBoxGenerator(11, nofElectrons); + if (nElectrons > 0) { + FairBoxGenerator* boxGen1 = new FairBoxGenerator(11, nElectrons); boxGen1->SetPtRange(0., 3.); + // boxGen1->SetPRange(0., 3.); boxGen1->SetPhiRange(0., 360.); boxGen1->SetThetaRange(2.5, 25.); - //boxGen1->SetXYZ(0., 0., targetZ); + // boxGen1->SetXYZ(0., 0., targetZ); // Dont use together with run.SetTarget(..., targetZ), boxGen will be at 2*targetZ boxGen1->SetCosTheta(); boxGen1->Init(); run.AddInput(boxGen1); } - if (nofPositrons > 0) { - FairBoxGenerator* boxGen2 = new FairBoxGenerator(-11, nofPositrons); + if (nPositrons > 0) { + FairBoxGenerator* boxGen2 = new FairBoxGenerator(-11, nPositrons); boxGen2->SetPtRange(0., 3.); + // boxGen2->SetPRange(0., 3.); boxGen2->SetPhiRange(0., 360.); boxGen2->SetThetaRange(2.5, 25.); - // boxGen2->SetXYZ(0., 0., targetZ); + // boxGen2->SetXYZ(0., 0., targetZ); // Dont use together with run.SetTarget(..., targetZ), boxGen will be at 2*targetZ boxGen2->SetCosTheta(); boxGen2->Init(); run.AddInput(boxGen2); } - if (plutoFile.length() > 0) { run.AddInput(plutoFile.c_str(), kPluto); } + if (!plutoFile.IsNull()) { run.AddInput(plutoFile, kPluto); } - run.SetOutFileName(traFile.c_str()); - run.SetParFileName(parFile.c_str()); - run.SetGeoFileName(geoFile.c_str()); - run.LoadSetup(geoSetup.c_str()); + if (!urqmdFile.IsNull()) { + if (urqmdFile == "default10gev") run.AddInput(defaultInput10gev); + else if (urqmdFile == "default25gev") + run.AddInput(defaultInput25gev); + else + run.AddInput(urqmdFile); + } + run.SetOutFileName(traFile); + run.SetParFileName(parFile); + run.SetGeoFileName(geoFile); + run.LoadSetup(geoSetup); run.SetTarget("Gold", 0.025, 2.5, 0, 0, targetZ); run.SetBeamPosition(0., 0., 0.1, 0.1); - //run.SetRandomEventPlane(); - //run.SetEngine(kGeant4); - //run.StoreTrajectories(true); + run.SetRandomEventPlane(); + run.SetRandomSeed(randomSeed); + // run.StoreTrajectories(true); run.Run(nEvents); + // ------------------------------------------------------------------------ - + // ----- Finish ------------------------------------------------------- timer.Stop(); - std::cout << std::endl << std::endl; + std::cout << std::endl; std::cout << "Macro finished successfully." << std::endl; std::cout << "Output file is " << traFile << std::endl; std::cout << "Parameter file is " << parFile << std::endl; std::cout << "Geometry file is " << geoFile << std::endl; std::cout << "Real time " << timer.RealTime() << " s, CPU time " << timer.CpuTime() << "s" << std::endl; - std::cout << "Test passed" << std::endl << "All ok" << std::endl; + std::cout << std::endl << "Test passed" << std::endl << "All ok" << std::endl; + // ------------------------------------------------------------------------ } diff --git a/reco/detectors/rich/CbmRichHitProducer.cxx b/reco/detectors/rich/CbmRichHitProducer.cxx index da37ec3b2c..9595ed8d11 100644 --- a/reco/detectors/rich/CbmRichHitProducer.cxx +++ b/reco/detectors/rich/CbmRichHitProducer.cxx @@ -18,12 +18,11 @@ #include "CbmRichDigiMapManager.h" #include "CbmRichGeoManager.h" #include "CbmRichHit.h" -#include "CbmRichPoint.h" #include <FairRootManager.h> #include <Logger.h> -#include "TClonesArray.h" +#include <TClonesArray.h> #include <TStopwatch.h> #include <iomanip> @@ -47,14 +46,14 @@ InitStatus CbmRichHitProducer::Init() FairRootManager* manager = FairRootManager::Instance(); fCbmEvents = dynamic_cast<TClonesArray*>(manager->GetObject("CbmEvent")); - if (fCbmEvents == nullptr) { LOG(info) << ": CbmEvent NOT found \n \n \n"; } + if (!fCbmEvents) { LOG(info) << ": CbmEvent NOT found \n \n \n"; } else { LOG(info) << ": CbmEvent found \n \n \n"; } fDigiMan = CbmDigiManager::Instance(); fDigiMan->Init(); - if (!fDigiMan->IsPresent(ECbmModuleId::kRich)) { Fatal("CbmRichHitProducer::Init", "No RichDigi array!"); } + if (!fDigiMan->IsPresent(ECbmModuleId::kRich)) { LOG(fatal) << "CbmRichHitProducer::Init: No RichDigi array!"; } fRichHits = new TClonesArray("CbmRichHit"); manager->Register("RichHit", "RICH", fRichHits, IsOutputBranchPersistent("RichHit")); @@ -71,12 +70,13 @@ void CbmRichHitProducer::Exec(Option_t* /*option*/) timer.Start(); Int_t nDigisAll = fDigiMan->GetNofDigis(ECbmModuleId::kRich); Int_t nDigisUsed = 0; + Int_t nHits = 0; Int_t nEvents = 0; Int_t result = 0; fRichHits->Delete(); // Time-slice processing - if (fCbmEvents == nullptr) nDigisUsed = ProcessData(nullptr); + if (!fCbmEvents) nDigisUsed = ProcessData(nullptr); // Event processing else { @@ -89,27 +89,30 @@ void CbmRichHitProducer::Exec(Option_t* /*option*/) } timer.Stop(); + nHits = fRichHits->GetEntriesFast(); stringstream logOut; logOut << setw(20) << left << GetName() << " ["; logOut << fixed << setw(8) << setprecision(1) << right << timer.RealTime() * 1000. << " ms] "; logOut << "TS " << fNofTs; - if (fCbmEvents != nullptr) logOut << ", events " << nEvents; - logOut << ", digis " << nDigisUsed << " / " << nDigisAll; + if (fCbmEvents) logOut << ", events " << nEvents; + logOut << ", digis " << nDigisUsed << " / " << nDigisAll << ", hits " << nHits; LOG(info) << logOut.str(); fNofTs++; + fNofEvents += nEvents; fNofDigisAll += nDigisAll; fNofDigisUsed += nDigisUsed; + fNofHitsAll += nHits; fTime += timer.RealTime(); } Int_t CbmRichHitProducer::ProcessData(CbmEvent* event) { Int_t nDigis = 0; - if (event != NULL) { - Int_t nofDigis = event->GetNofData(ECbmDataType::kRichDigi); + if (event) { + Int_t nofDigis = static_cast<Int_t>(event->GetNofData(ECbmDataType::kRichDigi)); LOG(debug) << GetName() << ": Event mode. Event # " << event->GetNumber() << ", digis: " << nofDigis; for (Int_t iDigi = 0; iDigi < nofDigis; iDigi++) { - Int_t digiIndex = event->GetIndex(ECbmDataType::kRichDigi, iDigi); + Int_t digiIndex = static_cast<Int_t>(event->GetIndex(ECbmDataType::kRichDigi, iDigi)); ProcessDigi(event, digiIndex); } nDigis = nofDigis; @@ -126,7 +129,7 @@ Int_t CbmRichHitProducer::ProcessData(CbmEvent* event) void CbmRichHitProducer::ProcessDigi(CbmEvent* event, Int_t digiIndex) { const CbmRichDigi* digi = fDigiMan->Get<CbmRichDigi>(digiIndex); - if (digi == nullptr) return; + if (!digi) return; if (digi->GetAddress() < 0) return; CbmRichPixelData* data = CbmRichDigiMapManager::GetInstance().GetPixelDataByAddress(digi->GetAddress()); TVector3 posPoint; @@ -142,14 +145,14 @@ void CbmRichHitProducer::AddHit(CbmEvent* event, TVector3& posHit, Double_t time { Int_t nofHits = fRichHits->GetEntriesFast(); new ((*fRichHits)[nofHits]) CbmRichHit(); - CbmRichHit* hit = (CbmRichHit*) fRichHits->At(nofHits); + CbmRichHit* hit = static_cast<CbmRichHit*>(fRichHits->At(nofHits)); hit->SetPosition(posHit); hit->SetDx(fHitError); hit->SetDy(fHitError); hit->SetRefId(index); hit->SetTime(time); - if (event != NULL) { event->AddData(ECbmDataType::kRichHit, nofHits); } + if (event) { event->AddData(ECbmDataType::kRichHit, nofHits); } } void CbmRichHitProducer::Finish() @@ -159,10 +162,18 @@ void CbmRichHitProducer::Finish() LOG(info) << "====================================="; LOG(info) << GetName() << ": Run summary"; LOG(info) << "Time slices : " << fNofTs; - if (fCbmEvents) LOG(info) << "Events : " << fNofEvents; LOG(info) << "Digis / TS : " << fixed << setprecision(2) << Double_t(fNofDigisAll) / Double_t(fNofTs); LOG(info) << "Used digis / TS : " << fixed << setprecision(2) << Double_t(fNofDigisUsed) / Double_t(fNofTs); LOG(info) << "Time / TS : " << fixed << setprecision(2) << 1000. * fTime / Double_t(fNofTs) << " ms"; + if (fCbmEvents) { + LOG(info) << "Events : " << fNofEvents; + LOG(info) << "Events / TS : " << fixed << setprecision(2) << fNofEvents / (Double_t) fNofTs; + if (fNofEvents > 0) { + LOG(info) << "Digis / Ev : " << fixed << setprecision(2) << Double_t(fNofDigisAll) / Double_t(fNofEvents); + LOG(info) << "Used digis / Ev : " << fixed << setprecision(2) << Double_t(fNofDigisUsed) / Double_t(fNofEvents); + LOG(info) << "Time / Ev : " << fixed << setprecision(2) << 1000. * fTime / Double_t(fNofEvents) << " ms"; + } + } LOG(info) << "=====================================\n"; } diff --git a/reco/detectors/rich/CbmRichHitProducer.h b/reco/detectors/rich/CbmRichHitProducer.h index 7934dc25bd..106e657bcb 100644 --- a/reco/detectors/rich/CbmRichHitProducer.h +++ b/reco/detectors/rich/CbmRichHitProducer.h @@ -86,6 +86,7 @@ private: Int_t fNofEvents = 0; // number of events Long64_t fNofDigisAll = 0; // all digis in input Long64_t fNofDigisUsed = 0; // digis used for hit finding + Long64_t fNofHitsAll = 0; // all hits in output Double_t fTime = 0.; // processing time Bool_t fRotationNeeded = kTRUE; diff --git a/reco/detectors/rich/CbmRichReconstruction.cxx b/reco/detectors/rich/CbmRichReconstruction.cxx index a93c762a51..5db3d9d9d2 100644 --- a/reco/detectors/rich/CbmRichReconstruction.cxx +++ b/reco/detectors/rich/CbmRichReconstruction.cxx @@ -49,8 +49,10 @@ #include <iomanip> #include <iostream> -using std::cout; -using std::endl; +using std::fixed; +using std::right; +using std::setprecision; +using std::setw; CbmRichReconstruction::CbmRichReconstruction() : FairTask("CbmRichReconstruction") {} @@ -78,7 +80,7 @@ InitStatus CbmRichReconstruction::Init() fRichTrackParamZ = new TClonesArray("FairTrackParam", 100); manager->Register("RichTrackParamZ", "RICH", fRichTrackParamZ, IsOutputBranchPersistent("RichTrackParamZ")); - fGlobalTracks = (TClonesArray*) manager->GetObject("GlobalTrack"); + fGlobalTracks = static_cast<TClonesArray*>(manager->GetObject("GlobalTrack")); if (fGlobalTracks == nullptr) LOG(fatal) << "CbmRichReconstruction::Init(): No GlobalTrack array."; } @@ -88,7 +90,7 @@ InitStatus CbmRichReconstruction::Init() manager->Register("RichProjection", "RICH", fRichProjections, IsOutputBranchPersistent("RichProjection")); } - fRichHits = (TClonesArray*) manager->GetObject("RichHit"); + fRichHits = static_cast<TClonesArray*>(manager->GetObject("RichHit")); if (fRichHits == nullptr) LOG(fatal) << "CbmRichReconstruction::Init(): No RichHit array."; fRichRings = new TClonesArray("CbmRichRing", 100); @@ -118,43 +120,82 @@ void CbmRichReconstruction::Exec(Option_t* /*opt*/) { TStopwatch timer; timer.Start(); + Int_t nEvents {0}; + Int_t nTrackProj {0}; + Int_t nGlobalTracks {0}; + if (fRichTrackParamZ != nullptr) fRichTrackParamZ->Delete(); if (fRichProjections != nullptr) fRichProjections->Delete(); if (fRichRings != nullptr) fRichRings->Delete(); - fNofTs++; - if (fCbmEvents == nullptr) { - fNofEvents++; ProcessData(nullptr); + nTrackProj += (fRichProjections ? fProjectionProducer->GetSuccessfullProj() : 0); + nGlobalTracks += (fGlobalTracks ? fGlobalTracks->GetEntriesFast() : 0); } else { - int nEvents = fCbmEvents->GetEntriesFast(); + nEvents = fCbmEvents->GetEntriesFast(); fNofEvents += nEvents; for (Int_t iEvent = 0; iEvent < nEvents; iEvent++) { CbmEvent* event = static_cast<CbmEvent*>(fCbmEvents->At(iEvent)); ProcessData(event); + nTrackProj += (fRichProjections ? fProjectionProducer->GetSuccessfullProj() : 0); + nGlobalTracks += (fGlobalTracks ? event->GetNofData(ECbmDataType::kGlobalTrack) : 0); } } + timer.Stop(); - fCalcTime += timer.RealTime(); + std::stringstream logOut; + logOut << setw(20) << left << GetName() << "["; + logOut << fixed << setw(8) << setprecision(1) << right << timer.RealTime() * 1000. << " ms] "; + logOut << "TS " << fNofTs; + if (fCbmEvents) logOut << ", events " << nEvents; + logOut << ", hits " << fRichHits->GetEntriesFast(); + logOut << ", rings " << fRichRings->GetEntriesFast(); + if (fRichProjections) logOut << ", trackProj " << nTrackProj << " / " << nGlobalTracks; + LOG(info) << logOut.str(); + + fNofTs++; + fCalcTime[0] += timer.RealTime(); fTotalNofHits += fRichHits->GetEntriesFast(); fTotalNofRings += fRichRings->GetEntriesFast(); - fTotalNofTrackProj += (fRichProjections ? fRichProjections->GetEntriesFast() : 0); - - LOG(info) << setw(20) << left << GetName() << "[" << fixed << setw(8) << setprecision(1) << right - << timer.RealTime() * 1000. << " ms] " - << "TS " << fNofTs << ", hits " << fRichHits->GetEntriesFast() << ", rings " << fRichRings->GetEntriesFast() - << ", trackProj " << (fRichProjections ? fRichProjections->GetEntriesFast() : 0); + fTotalNofTrackProj += nTrackProj; + fTotalNofGlobalTracks += nGlobalTracks; } void CbmRichReconstruction::ProcessData(CbmEvent* event) { - if (fRunExtrapolation) RunExtrapolation(event); - if (fRunProjection) RunProjection(event); - if (fRunFinder) RunFinder(event); - if (fRunFitter) RunFitter(event); - if (fRunTrackAssign) RunTrackAssign(event); + TStopwatch timer; + if (fRunExtrapolation) { + timer.Start(); + RunExtrapolation(event); + timer.Stop(); + fCalcTime[1] += timer.RealTime(); + } + if (fRunProjection) { + timer.Start(); + RunProjection(event); + timer.Stop(); + fCalcTime[2] += timer.RealTime(); + } + if (fRunFinder) { + timer.Start(); + RunFinder(event); + timer.Stop(); + fCalcTime[3] += timer.RealTime(); + } + if (fRunFitter) { + timer.Start(); + RunFitter(event); + timer.Stop(); + fCalcTime[4] += timer.RealTime(); + } + if (fRunTrackAssign) { + timer.Start(); + RunTrackAssign(event); + timer.Stop(); + fCalcTime[5] += timer.RealTime(); + } } void CbmRichReconstruction::InitExtrapolation() @@ -286,21 +327,55 @@ void CbmRichReconstruction::RunTrackAssign(CbmEvent* event) void CbmRichReconstruction::Finish() { - std::cout << std::endl; LOG(info) << "====================================="; LOG(info) << GetName() << ": Run summary"; LOG(info) << "Time slices : " << fNofTs; - LOG(info) << "Events : " << fNofEvents; - LOG(info) << "Events / TS : " << fixed << setprecision(2) << fNofEvents / (Double_t) fNofTs; - LOG(info) << "Hits / TS : " << fixed << setprecision(2) << fTotalNofHits / (Double_t) fNofTs; - LOG(info) << "Ring / TS : " << fixed << setprecision(2) << fTotalNofRings / (Double_t) fNofTs; - LOG(info) << "TrackProj / TS : " << fixed << setprecision(2) << fTotalNofTrackProj / (Double_t) fNofTs; - LOG(info) << "Time / TS : " << fixed << setprecision(2) << 1000. * fCalcTime / Double_t(fNofTs) << " ms"; - LOG(info) << "Hits / ev : " << fixed << setprecision(2) << fTotalNofHits / (Double_t) fNofEvents; - LOG(info) << "Ring / ev : " << fixed << setprecision(2) << fTotalNofRings / (Double_t) fNofEvents; - LOG(info) << "TrackProj / ev : " << fixed << setprecision(2) << fTotalNofTrackProj / (Double_t) fNofEvents; - LOG(info) << "Time / ev : " << fixed << setprecision(2) << 1000. * fCalcTime / Double_t(fNofEvents) << " ms"; + LOG(info) << "Hits / TS : " << fixed << setprecision(2) << fTotalNofHits / Double_t(fNofTs); + LOG(info) << "Rings / TS : " << fixed << setprecision(2) << fTotalNofRings / Double_t(fNofTs); + if (fRichProjections) + LOG(info) << "TrackProj / TS : " << fixed << setprecision(2) << fTotalNofTrackProj / Double_t(fNofTs); + LOG(info) << "Time / TS : " << fixed << setprecision(2) << 1000. * fCalcTime[0] / Double_t(fNofTs) << " ms"; + if (fCbmEvents) { + LOG(info) << "Events : " << fNofEvents; + LOG(info) << "Events / TS : " << fixed << setprecision(2) << fNofEvents / Double_t(fNofTs); + if (fNofEvents > 0) { + LOG(info) << "Hits / ev : " << fixed << setprecision(2) << fTotalNofHits / Double_t(fNofEvents); + LOG(info) << "Rings / ev : " << fixed << setprecision(2) << fTotalNofRings / Double_t(fNofEvents); + if (fRichProjections) + LOG(info) << "TrackProj / ev : " << fixed << setprecision(2) << fTotalNofTrackProj / Double_t(fNofEvents); + LOG(info) << "Time / ev : " << fixed << setprecision(2) << 1000. * fCalcTime[0] / Double_t(fNofEvents) + << " ms"; + } + } + if (fRichProjections && fTotalNofGlobalTracks > 0) + LOG(info) << "TrackProj / GTr : " << fixed << setprecision(2) + << fTotalNofTrackProj / Double_t(fTotalNofGlobalTracks); + TString eventOrTsStr = fCbmEvents && fNofEvents > 0 ? "ev: " : "TS: "; + Double_t eventOrTsValue = Double_t(fCbmEvents && fNofEvents > 0 ? fNofEvents : fNofTs); + if (fCalcTime[0] != 0.) { + LOG(info) << "===== Time by task (real time) ======"; + if (fTrackExtrapolation) + LOG(info) << "TrackExtrapolation / " << eventOrTsStr << fixed << setprecision(2) << setw(9) << right + << 1000. * fCalcTime[1] / eventOrTsValue << " ms [" << setw(5) << right + << 100 * fCalcTime[1] / fCalcTime[0] << " %]"; + if (fProjectionProducer) + LOG(info) << "TrackProjection / " << eventOrTsStr << fixed << setprecision(2) << setw(9) << right + << 1000. * fCalcTime[2] / eventOrTsValue << " ms [" << setw(5) << right + << 100 * fCalcTime[2] / fCalcTime[0] << " %]"; + if (fRingFinder) + LOG(info) << "RingFinder / " << eventOrTsStr << fixed << setprecision(2) << setw(9) << right + << 1000. * fCalcTime[3] / eventOrTsValue << " ms [" << setw(5) << right + << 100 * fCalcTime[3] / fCalcTime[0] << " %]"; + if (fRingFitter) + LOG(info) << "RingFitter / " << eventOrTsStr << fixed << setprecision(2) << setw(9) << right + << 1000. * fCalcTime[4] / eventOrTsValue << " ms [" << setw(5) << right + << 100 * fCalcTime[4] / fCalcTime[0] << " %]"; + if (fRingTrackAssign) + LOG(info) << "RingTrackAssign / " << eventOrTsStr << fixed << setprecision(2) << setw(9) << right + << 1000. * fCalcTime[5] / eventOrTsValue << " ms [" << setw(5) << right + << 100 * fCalcTime[5] / fCalcTime[0] << " %]"; + } LOG(info) << "=====================================\n"; } diff --git a/reco/detectors/rich/CbmRichReconstruction.h b/reco/detectors/rich/CbmRichReconstruction.h index 70b2b69a0b..8bd66b6e65 100644 --- a/reco/detectors/rich/CbmRichReconstruction.h +++ b/reco/detectors/rich/CbmRichReconstruction.h @@ -14,7 +14,7 @@ #ifndef CBM_RICH_RECONSTRUCTION #define CBM_RICH_RECONSTRUCTION -#include "FairTask.h" +#include <FairTask.h> #include <string> @@ -102,12 +102,14 @@ private: TClonesArray* fGlobalTracks = nullptr; TClonesArray* fCbmEvents = nullptr; - Int_t fNofTs = 0; - Int_t fNofEvents = 0; - Double_t fCalcTime = 0.; - Int_t fTotalNofHits = 0; - Int_t fTotalNofRings = 0; - Int_t fTotalNofTrackProj = 0; + std::array<Double_t, 6> fCalcTime {0., 0., 0., 0., 0., 0.}; + + Int_t fNofTs = 0; + Int_t fNofEvents = 0; + Int_t fTotalNofHits = 0; + Int_t fTotalNofRings = 0; + Int_t fTotalNofTrackProj = 0; + Int_t fTotalNofGlobalTracks = 0; // pointers to the algorithms diff --git a/reco/detectors/rich/tracks/CbmRichProjectionProducerAnalytical.cxx b/reco/detectors/rich/tracks/CbmRichProjectionProducerAnalytical.cxx index 2ad1b42b0c..b68008304f 100644 --- a/reco/detectors/rich/tracks/CbmRichProjectionProducerAnalytical.cxx +++ b/reco/detectors/rich/tracks/CbmRichProjectionProducerAnalytical.cxx @@ -54,7 +54,7 @@ void CbmRichProjectionProducerAnalytical::Init() void CbmRichProjectionProducerAnalytical::DoProjection(CbmEvent* event, TClonesArray* richProj) { - + fnSuccessfullProj = 0; if (richProj == nullptr) { LOG(error) << "CbmRichProjectionProducerAnalytical::DoExtrapolation(): richProj is nullptr."; return; @@ -201,6 +201,7 @@ void CbmRichProjectionProducerAnalytical::DoProjection(CbmEvent* event, TClonesA if (isInsidePmt) { FairTrackParam richtrack(outPos.x(), outPos.y(), outPos.z(), 0., 0., 0., covMat); *(FairTrackParam*) (richProj->At(iT)) = richtrack; + fnSuccessfullProj++; } } } diff --git a/reco/detectors/rich/tracks/CbmRichProjectionProducerBase.h b/reco/detectors/rich/tracks/CbmRichProjectionProducerBase.h index 9e4378e351..278f715a72 100644 --- a/reco/detectors/rich/tracks/CbmRichProjectionProducerBase.h +++ b/reco/detectors/rich/tracks/CbmRichProjectionProducerBase.h @@ -49,6 +49,12 @@ public: **/ virtual void DoProjection(CbmEvent* event, TClonesArray* richProj) = 0; + /** Get number of successful projections */ + int GetSuccessfullProj() const { return fnSuccessfullProj; } + +protected: + int fnSuccessfullProj {}; + private: /** * \brief Copy constructor. diff --git a/reco/detectors/rich/tracks/CbmRichProjectionProducerTGeo.cxx b/reco/detectors/rich/tracks/CbmRichProjectionProducerTGeo.cxx index bf877b0b79..867a2c6820 100644 --- a/reco/detectors/rich/tracks/CbmRichProjectionProducerTGeo.cxx +++ b/reco/detectors/rich/tracks/CbmRichProjectionProducerTGeo.cxx @@ -55,7 +55,7 @@ void CbmRichProjectionProducerTGeo::Init() void CbmRichProjectionProducerTGeo::DoProjection(CbmEvent* event, TClonesArray* richProj) { - + fnSuccessfullProj = 0; if (richProj == nullptr) { LOG(error) << "CbmRichProjectionProducerTGeo::DoExtrapolation(): richProj is nullptr."; return; @@ -138,6 +138,7 @@ void CbmRichProjectionProducerTGeo::DoProjection(CbmEvent* event, TClonesArray* // if(TMath::Abs(yDet) > (fGP.fPmtY-fGP.fPmtWidthY) && TMath::Abs(yDet) < (fGP.fPmtY+fGP.fPmtWidthY)){ FairTrackParam richtrack(xDet, yDet, zDet, 0., 0., 0., covMat); *(FairTrackParam*) (richProj->At(iT)) = richtrack; + fnSuccessfullProj++; // } // } } // j -- GitLab