From 9a55db44f10b0a5361746f7b481a7909f7b51eea Mon Sep 17 00:00:00 2001
From: Oleksii Lubynets <o.lubynets@gsi.de>
Date: Mon, 1 Aug 2022 18:22:18 +0000
Subject: [PATCH] Event matching for time based mode
---
.../analysis_tree_converter/CMakeLists.txt | 2 +
.../CbmAnalysisTreeInterfaceLinkDef.h | 1 +
.../CbmConverterManager.cxx | 35 +++-
.../CbmConverterManager.h | 21 +-
.../CbmConverterTask.h | 22 +-
.../CbmMatchEvents.cxx | 60 ++++++
.../analysis_tree_converter/CbmMatchEvents.h | 39 ++++
.../CbmPsdModulesConverter.cxx | 18 +-
.../CbmPsdModulesConverter.h | 7 +-
.../CbmRecEventHeaderConverter.cxx | 46 ++++-
.../CbmRecEventHeaderConverter.h | 6 +-
.../CbmRichRingsConverter.cxx | 32 +--
.../CbmRichRingsConverter.h | 4 +-
.../CbmSimEventHeaderConverter.cxx | 54 +++--
.../CbmSimEventHeaderConverter.h | 7 +-
.../CbmSimTracksConverter.cxx | 189 +++++++++---------
.../CbmSimTracksConverter.h | 49 ++---
.../CbmStsTracksConverter.cxx | 163 +++++++++------
.../CbmStsTracksConverter.h | 16 +-
.../CbmTofHitsConverter.cxx | 54 +++--
.../CbmTofHitsConverter.h | 14 +-
.../CbmTrdTracksConverter.cxx | 36 ++--
.../CbmTrdTracksConverter.h | 11 +-
external/InstallAnalysisTree.cmake | 2 +-
external/InstallAnalysisTreeQA.cmake | 2 +-
.../run_analysis_tree_maker.C | 12 +-
.../analysis/common/qa/run_analysistree_qa.C | 8 +
27 files changed, 605 insertions(+), 305 deletions(-)
create mode 100644 analysis/common/analysis_tree_converter/CbmMatchEvents.cxx
create mode 100644 analysis/common/analysis_tree_converter/CbmMatchEvents.h
diff --git a/analysis/common/analysis_tree_converter/CMakeLists.txt b/analysis/common/analysis_tree_converter/CMakeLists.txt
index de5f664286..1fee3b8512 100644
--- a/analysis/common/analysis_tree_converter/CMakeLists.txt
+++ b/analysis/common/analysis_tree_converter/CMakeLists.txt
@@ -11,6 +11,7 @@ set(SRCS
CbmTofHitsConverter.cxx
CbmTrdTracksConverter.cxx
CbmRichRingsConverter.cxx
+ CbmMatchEvents.cxx
)
@@ -26,6 +27,7 @@ Set(INCLUDE_DIRECTORIES
${CBMROOT_SOURCE_DIR}/core/data/tof
${CBMROOT_SOURCE_DIR}/core/data/psd
${CBMROOT_SOURCE_DIR}/core/data/trd
+ ${CBMROOT_SOURCE_DIR}/core/data/rich
${CBMROOT_SOURCE_DIR}/sim/transport/generators
${CBMROOT_SOURCE_DIR}/sim/transport/steer
diff --git a/analysis/common/analysis_tree_converter/CbmAnalysisTreeInterfaceLinkDef.h b/analysis/common/analysis_tree_converter/CbmAnalysisTreeInterfaceLinkDef.h
index d09b7a7516..eac02d0865 100644
--- a/analysis/common/analysis_tree_converter/CbmAnalysisTreeInterfaceLinkDef.h
+++ b/analysis/common/analysis_tree_converter/CbmAnalysisTreeInterfaceLinkDef.h
@@ -18,6 +18,7 @@
#pragma link C++ class CbmPsdModulesConverter + ;
#pragma link C++ class CbmTrdTracksConverter + ;
#pragma link C++ class CbmRichRingsConverter + ;
+//#pragma link C++ class CbmMatchEvents + ;
//#pragma link C++ class CbmMuchTracksConverter+;
diff --git a/analysis/common/analysis_tree_converter/CbmConverterManager.cxx b/analysis/common/analysis_tree_converter/CbmConverterManager.cxx
index 15a98a75d7..5a923cf6fc 100644
--- a/analysis/common/analysis_tree_converter/CbmConverterManager.cxx
+++ b/analysis/common/analysis_tree_converter/CbmConverterManager.cxx
@@ -5,7 +5,9 @@
#include "CbmConverterManager.h"
#include "CbmConverterTask.h"
+#include "CbmEvent.h"
+#include "TClonesArray.h"
#include "TGeoBBox.h"
#include "TGeoManager.h"
@@ -14,34 +16,55 @@
#include "AnalysisTree/DataHeader.hpp"
#include "AnalysisTree/TaskManager.hpp"
-ClassImp(CbmConverterManager)
+ClassImp(CbmConverterManager);
- InitStatus CbmConverterManager::Init()
+InitStatus CbmConverterManager::Init()
{
task_manager_->Init();
FillDataHeader();
+ InitEvent();
return kSUCCESS;
}
-void CbmConverterManager::Exec(Option_t* /*opt*/)
+void CbmConverterManager::AddTask(CbmConverterTask* task)
+{
+ tasks_.emplace_back(task);
+ task_manager_->AddTask(reinterpret_cast<AnalysisTree::Task*>(task));
+}
+
+void CbmConverterManager::ProcessData(CbmEvent* event)
{
index_map_.clear();
for (auto* task : tasks_) {
task->SetIndexesMap(&index_map_);
- task->Exec();
+ task->ProcessData(event);
index_map_.insert(std::make_pair(task->GetOutputBranchName(), task->GetOutIndexesMap()));
}
task_manager_->FillOutput();
}
+void CbmConverterManager::Exec(Option_t* /*opt*/)
+{
+ if (events_ != nullptr) {
+ auto n_events = events_->GetEntriesFast();
+ for (int i_event = 0; i_event < n_events; ++i_event) {
+ auto* event = (CbmEvent*) events_->At(i_event);
+ ProcessData(event);
+ }
+ }
+ else {
+ LOG(info) << "Event based mode\n";
+ ProcessData(nullptr);
+ }
+}
+
+
void CbmConverterManager::Finish()
{
TDirectory* curr = gDirectory; // TODO check why this is needed
TFile* currentFile = gFile;
- task_manager_->GetChain()->Write();
-
task_manager_->Finish();
gFile = currentFile;
diff --git a/analysis/common/analysis_tree_converter/CbmConverterManager.h b/analysis/common/analysis_tree_converter/CbmConverterManager.h
index 117c8ad45e..4934c7d524 100644
--- a/analysis/common/analysis_tree_converter/CbmConverterManager.h
+++ b/analysis/common/analysis_tree_converter/CbmConverterManager.h
@@ -18,6 +18,7 @@ namespace AnalysisTree
} // namespace AnalysisTree
class CbmConverterTask;
+class CbmEvent;
class CbmConverterManager : public FairTask {
@@ -29,19 +30,26 @@ public:
void Exec(Option_t* opt) override;
void Finish() override;
- void AddTask(CbmConverterTask* task)
- {
- tasks_.emplace_back(task);
- task_manager_->AddTask(reinterpret_cast<AnalysisTree::Task*>(task));
- }
+ void AddTask(CbmConverterTask* task);
void SetSystem(const std::string& system) { system_ = system; }
void SetBeamMomentum(float beam_mom) { beam_mom_ = beam_mom; }
- void SetOutputName(std::string file, std::string tree = "rTree") { task_manager_->SetOutputName(file, tree); }
+ void SetOutputName(std::string file, std::string tree = "rTree")
+ {
+ task_manager_->SetOutputName(std::move(file), std::move(tree));
+ }
+
+ void InitEvent()
+ {
+ auto* ioman = FairRootManager::Instance();
+ events_ = (TClonesArray*) ioman->GetObject("CbmEvent");
+ }
+
private:
void FillDataHeader();
+ void ProcessData(CbmEvent* event);
AnalysisTree::TaskManager* task_manager_ {AnalysisTree::TaskManager::GetInstance()};
@@ -51,6 +59,7 @@ private:
std::vector<CbmConverterTask*> tasks_ {};
std::map<std::string, std::map<int, int>> index_map_ {}; ///< map CbmRoot to AT of indexes for a given branch
+ TClonesArray* events_ {nullptr};
ClassDefOverride(CbmConverterManager, 1)
};
diff --git a/analysis/common/analysis_tree_converter/CbmConverterTask.h b/analysis/common/analysis_tree_converter/CbmConverterTask.h
index 707bf0fd5d..f2d4e24847 100644
--- a/analysis/common/analysis_tree_converter/CbmConverterTask.h
+++ b/analysis/common/analysis_tree_converter/CbmConverterTask.h
@@ -5,17 +5,23 @@
#ifndef ANALYSIS_TREE_CONVERTERTASK_H_
#define ANALYSIS_TREE_CONVERTERTASK_H_
+#include <CbmLink.h>
+
+#include <FairLogger.h>
+
#include <map>
#include <string>
#include "AnalysisTree/Task.hpp"
class FairRootManager;
+class CbmEvent;
class CbmConverterTask : public AnalysisTree::Task {
+ using MapType = std::map<int, int>;
public:
- CbmConverterTask() = delete;
+ CbmConverterTask() = default;
explicit CbmConverterTask(std::string out_branch_name, std::string match_to = "")
{
out_branch_ = std::move(out_branch_name);
@@ -24,17 +30,21 @@ public:
~CbmConverterTask() override = default;
- const std::map<int, int>& GetOutIndexesMap() const { return out_indexes_map_; }
+ virtual void ProcessData(CbmEvent* event) = 0;
+
+ void Exec() final { throw std::runtime_error("Should not be used!"); };
- void SetIndexesMap(std::map<std::string, std::map<int, int>>* indexes_map) { indexes_map_ = indexes_map; }
+ const MapType& GetOutIndexesMap() const { return out_indexes_map_; }
+
+ void SetIndexesMap(std::map<std::string, MapType>* indexes_map) { indexes_map_ = indexes_map; }
const std::string& GetOutputBranchName() const { return out_branch_; }
+
protected:
- std::map<int, int> out_indexes_map_ {}; ///< CbmRoot to AnalysisTree indexes
- ///< map for output branch
+ MapType out_indexes_map_ {}; ///< CbmRoot to AnalysisTree indexes map for output branch
std::string out_branch_ {};
- std::map<std::string, std::map<int, int>>* indexes_map_ {}; ///< CbmRoot to AnalysisTree indexes map for branches
+ std::map<std::string, MapType>* indexes_map_ {}; ///< CbmRoot to AnalysisTree indexes map for branches
///< from other tasks
std::string match_to_ {}; ///< AT branch to match
};
diff --git a/analysis/common/analysis_tree_converter/CbmMatchEvents.cxx b/analysis/common/analysis_tree_converter/CbmMatchEvents.cxx
new file mode 100644
index 0000000000..b08cc56ffa
--- /dev/null
+++ b/analysis/common/analysis_tree_converter/CbmMatchEvents.cxx
@@ -0,0 +1,60 @@
+/* Copyright (C) 2020-2021 Physikalisches Institut, Eberhard Karls Universität Tuebingen, Tuebingen
+ SPDX-License-Identifier: GPL-3.0-only
+ Authors: Viktor Klochkov [committer] */
+
+#include "CbmMatchEvents.h"
+
+#include "CbmEvent.h"
+#include "CbmMCDataManager.h"
+#include "CbmTrackMatchNew.h"
+
+#include "TClonesArray.h"
+
+#include "AnalysisTree/TaskManager.hpp"
+
+//ClassImp(CbmMatchEvents);
+
+void CbmMatchEvents::Init()
+{
+ auto* ioman = FairRootManager::Instance();
+ cbm_sts_match_ = (TClonesArray*) ioman->GetObject("StsTrackMatch");
+}
+
+void CbmMatchEvents::ProcessData(CbmEvent* event)
+{
+ if (!event) { throw std::runtime_error("No event to match"); }
+
+ LOG(info) << "Event: " << event->GetNumber() << " t_start " << event->GetStartTime()
+ << " Msts = " << event->GetNofStsTracks();
+
+ count_map_.clear();
+
+ const int n_sts_tracks = event->GetNofStsTracks();
+ for (short i_track = 0; i_track < n_sts_tracks; ++i_track) {
+ const auto track_index = event->GetStsTrackIndex(i_track);
+
+ auto* match = (CbmTrackMatchNew*) cbm_sts_match_->At(track_index);
+ if (match->GetNofLinks() == 0) { continue; }
+
+ const auto& link = match->GetMatchedLink();
+ auto file = link.GetFile();
+ auto entry = link.GetEntry();
+
+ count_map_[{file, entry}]++;
+ }
+ // TODO remove later if not needed
+ // using pair_type = decltype(count_map_)::value_type;
+ // auto matched_event = std::max_element(count_map_.begin(), count_map_.end(),
+ // [](const pair_type& a, const pair_type& b) -> bool{ return a.second < b.second; } );
+
+ auto match_event = new CbmMatch();
+
+ int i {0};
+ for (const auto& match : count_map_) {
+ auto weight = float(match.second) / n_sts_tracks;
+ match_event->AddLink(weight, i++, match.first.entry, match.first.file);
+ LOG(info) << " matched to " << match.first.entry << " " << match.second << " weight = " << weight;
+ }
+
+ if (match_event->GetNofLinks() > 0) event->SetMatch(match_event);
+}
diff --git a/analysis/common/analysis_tree_converter/CbmMatchEvents.h b/analysis/common/analysis_tree_converter/CbmMatchEvents.h
new file mode 100644
index 0000000000..91445eeb63
--- /dev/null
+++ b/analysis/common/analysis_tree_converter/CbmMatchEvents.h
@@ -0,0 +1,39 @@
+/* Copyright (C) 2020-2021 Physikalisches Institut, Eberhard Karls Universität Tuebingen, Tuebingen
+ SPDX-License-Identifier: GPL-3.0-only
+ Authors: Viktor Klochkov [committer] */
+
+#ifndef CBMROOT_ANALYSIS_COMMON_ANALYSIS_TREE_CONVERTER_CBMMATCHEVENTS_H_
+#define CBMROOT_ANALYSIS_COMMON_ANALYSIS_TREE_CONVERTER_CBMMATCHEVENTS_H_
+
+#include "CbmConverterTask.h"
+
+class CbmMCDataManager;
+class CbmMCDataArray;
+
+class CbmMatchEvents final : public CbmConverterTask {
+public:
+ explicit CbmMatchEvents() = default;
+
+ ~CbmMatchEvents() final = default;
+
+ void Init() final;
+
+ void ProcessData(CbmEvent* event) final;
+
+ void Finish() final {}
+
+ struct EventId {
+ EventId(int f, int e) : file(f), entry(e) {}
+ int file {0};
+ int entry {0};
+ bool operator<(const EventId& other) const { return this->file < other.file || this->entry < other.entry; }
+ };
+
+private:
+ TClonesArray* cbm_sts_match_ {nullptr}; ///< non-owning pointer
+ std::map<EventId, int> count_map_ {}; ///!
+
+ // ClassDef(CbmMatchEvents, 1);
+};
+
+#endif //CBMROOT_ANALYSIS_COMMON_ANALYSIS_TREE_CONVERTER_CBMMATCHEVENTS_H_
diff --git a/analysis/common/analysis_tree_converter/CbmPsdModulesConverter.cxx b/analysis/common/analysis_tree_converter/CbmPsdModulesConverter.cxx
index 357ab54885..10240fff5c 100644
--- a/analysis/common/analysis_tree_converter/CbmPsdModulesConverter.cxx
+++ b/analysis/common/analysis_tree_converter/CbmPsdModulesConverter.cxx
@@ -4,6 +4,8 @@
#include "CbmPsdModulesConverter.h"
+#include "CbmDefs.h"
+#include "CbmEvent.h"
#include "CbmPsdHit.h"
#include "FairRootManager.h"
@@ -16,9 +18,9 @@
#include "AnalysisTree/Detector.hpp"
-ClassImp(CbmPsdModulesConverter)
+ClassImp(CbmPsdModulesConverter);
- void CbmPsdModulesConverter::Init()
+void CbmPsdModulesConverter::Init()
{
assert(!out_branch_.empty());
auto* ioman = FairRootManager::Instance();
@@ -28,11 +30,11 @@ ClassImp(CbmPsdModulesConverter)
AnalysisTree::BranchConfig psd_branch(out_branch_, AnalysisTree::DetType::kModule);
auto* man = AnalysisTree::TaskManager::GetInstance();
- man->AddBranch(out_branch_, psd_modules_, psd_branch);
+ man->AddBranch(psd_modules_, psd_branch);
}
-void CbmPsdModulesConverter::Exec()
+void CbmPsdModulesConverter::ProcessData(CbmEvent* event)
{
assert(cbm_psd_hits_);
psd_modules_->ClearChannels();
@@ -45,13 +47,19 @@ void CbmPsdModulesConverter::Exec()
const auto& branch = config->GetBranchConfig(out_branch_);
const int n_psd_modules = data_header->GetModulePositions(0).GetNumberOfChannels();
+
psd_modules_->Reserve(n_psd_modules);
for (int i = 0; i < n_psd_modules; ++i) {
auto& module = psd_modules_->AddChannel(branch);
module.SetSignal(0.f);
}
- const int nPsdHits = cbm_psd_hits_->GetEntriesFast();
+ const int nPsdHits = event ? event->GetNofData(ECbmDataType::kPsdHit) : cbm_psd_hits_->GetEntriesFast();
+ if (nPsdHits <= 0) {
+ LOG(warn) << "No PSD hits!";
+ return;
+ }
+
for (int i = 0; i < nPsdHits; ++i) {
hit = (CbmPsdHit*) cbm_psd_hits_->At(i);
if (hit == nullptr) continue;
diff --git a/analysis/common/analysis_tree_converter/CbmPsdModulesConverter.h b/analysis/common/analysis_tree_converter/CbmPsdModulesConverter.h
index 6b0725e60c..d089537d32 100644
--- a/analysis/common/analysis_tree_converter/CbmPsdModulesConverter.h
+++ b/analysis/common/analysis_tree_converter/CbmPsdModulesConverter.h
@@ -14,15 +14,12 @@ class TClonesArray;
class CbmPsdModulesConverter final : public CbmConverterTask {
public:
explicit CbmPsdModulesConverter(std::string out_branch_name, std::string match_to = "")
- : CbmConverterTask(std::move(out_branch_name), std::move(match_to))
- {
- out_branch_ = "PsdModules";
- };
+ : CbmConverterTask(std::move(out_branch_name), std::move(match_to)) {};
~CbmPsdModulesConverter() final;
void Init() final;
- void Exec() final;
+ void ProcessData(CbmEvent* event) final;
void Finish() final;
private:
diff --git a/analysis/common/analysis_tree_converter/CbmRecEventHeaderConverter.cxx b/analysis/common/analysis_tree_converter/CbmRecEventHeaderConverter.cxx
index e776304406..44c12f519c 100644
--- a/analysis/common/analysis_tree_converter/CbmRecEventHeaderConverter.cxx
+++ b/analysis/common/analysis_tree_converter/CbmRecEventHeaderConverter.cxx
@@ -17,9 +17,9 @@
#include "cassert"
#include "iostream"
-ClassImp(CbmRecEventHeaderConverter)
+ClassImp(CbmRecEventHeaderConverter);
- void CbmRecEventHeaderConverter::Init()
+void CbmRecEventHeaderConverter::Init()
{
assert(!out_branch_.empty());
auto* ioman = FairRootManager::Instance();
@@ -35,33 +35,59 @@ ClassImp(CbmRecEventHeaderConverter)
RecEventHeaderBranch.AddField<float>("Epsd", "GeV, full energy deposit in PSD");
RecEventHeaderBranch.AddField<int>("M", "total multiplicity in STS(+MVD)");
RecEventHeaderBranch.AddField<int>("evt_id", "event identifier");
+ RecEventHeaderBranch.AddField<float>("start_time", "Start time of the event, ns");
+ RecEventHeaderBranch.AddField<float>("end_time", "End time of the event, ns");
+ RecEventHeaderBranch.AddField<float>("match_weight", "");
auto* man = AnalysisTree::TaskManager::GetInstance();
- man->AddBranch(out_branch_, rec_event_header_, RecEventHeaderBranch);
+ man->AddBranch(rec_event_header_, RecEventHeaderBranch);
rec_event_header_->Init(RecEventHeaderBranch);
}
-void CbmRecEventHeaderConverter::Exec()
+void CbmRecEventHeaderConverter::ProcessData(CbmEvent* event)
{
-
auto* out_config_ = AnalysisTree::TaskManager::GetInstance()->GetConfig();
const auto& branch = out_config_->GetBranchConfig(out_branch_);
+ if (event) { cbm_prim_vertex_ = event->GetVertex(); }
+
if (!cbm_prim_vertex_) { throw std::runtime_error("No fPrimVtx"); }
rec_event_header_->SetVertexPosition3({cbm_prim_vertex_->GetX(), cbm_prim_vertex_->GetY(), cbm_prim_vertex_->GetZ()});
- rec_event_header_->SetField(int(cbm_sts_tracks_->GetEntries()), branch.GetFieldId("M"));
- rec_event_header_->SetField(GetPsdEnergy(), branch.GetFieldId("Epsd"));
- rec_event_header_->SetField(int(cbm_header_->GetEventID()), branch.GetFieldId("evt_id"));
rec_event_header_->SetField(float(cbm_prim_vertex_->GetChi2() / cbm_prim_vertex_->GetNDF()),
branch.GetFieldId("vtx_chi2"));
+
+ const int n_sts_tracks = event ? event->GetNofStsTracks() : cbm_sts_tracks_->GetEntries();
+ rec_event_header_->SetField(n_sts_tracks, branch.GetFieldId("M"));
+ rec_event_header_->SetField(GetPsdEnergy(event), branch.GetFieldId("Epsd"));
+
+ int evt_id;
+ float match_weight, start_time, end_time;
+ if (event) {
+ evt_id = event->GetUniqueID();
+ start_time = event->GetStartTime();
+ end_time = event->GetEndTime();
+ if (event->GetMatch()) match_weight = float(event->GetMatch()->GetMatchedLink().GetWeight());
+ else
+ match_weight = 0.;
+ }
+ else {
+ evt_id = cbm_header_->GetEventID();
+ start_time = cbm_header_->GetT();
+ end_time = cbm_header_->GetT();
+ match_weight = 1.;
+ }
+ rec_event_header_->SetField(evt_id, branch.GetFieldId("evt_id"));
+ rec_event_header_->SetField(start_time, branch.GetFieldId("start_time"));
+ rec_event_header_->SetField(end_time, branch.GetFieldId("end_time"));
+ rec_event_header_->SetField(match_weight, branch.GetFieldId("match_weight"));
}
-float CbmRecEventHeaderConverter::GetPsdEnergy()
+float CbmRecEventHeaderConverter::GetPsdEnergy(CbmEvent* event)
{
//TODO avoid duplicating the code with PsdModulesConverter
float psd_energy {0.f};
- const int nPsdHits = cbm_psd_hits_->GetEntriesFast();
+ const int nPsdHits = event ? event->GetNofData(ECbmDataType::kPsdHit) : cbm_psd_hits_->GetEntriesFast();
for (int i = 0; i < nPsdHits; ++i) {
auto* hit = (CbmPsdHit*) cbm_psd_hits_->At(i);
if (hit == nullptr) continue;
diff --git a/analysis/common/analysis_tree_converter/CbmRecEventHeaderConverter.h b/analysis/common/analysis_tree_converter/CbmRecEventHeaderConverter.h
index 3f542c74c2..bf1e5d20e3 100644
--- a/analysis/common/analysis_tree_converter/CbmRecEventHeaderConverter.h
+++ b/analysis/common/analysis_tree_converter/CbmRecEventHeaderConverter.h
@@ -6,6 +6,8 @@
#define ANALYSIS_TREE_RECEVENTHEADERCONVERTER_H_
#include "CbmConverterTask.h"
+#include "CbmDefs.h"
+#include "CbmEvent.h"
#include "AnalysisTree/EventHeader.hpp"
@@ -19,11 +21,11 @@ public:
~CbmRecEventHeaderConverter() final = default;
void Init() final;
- void Exec() final;
+ void ProcessData(CbmEvent* event) final;
void Finish() final { delete rec_event_header_; };
private:
- float GetPsdEnergy();
+ float GetPsdEnergy(CbmEvent* event);
AnalysisTree::EventHeader* rec_event_header_ {nullptr};
diff --git a/analysis/common/analysis_tree_converter/CbmRichRingsConverter.cxx b/analysis/common/analysis_tree_converter/CbmRichRingsConverter.cxx
index ba703d790a..1e9dcc737d 100644
--- a/analysis/common/analysis_tree_converter/CbmRichRingsConverter.cxx
+++ b/analysis/common/analysis_tree_converter/CbmRichRingsConverter.cxx
@@ -4,8 +4,10 @@
#include "CbmRichRingsConverter.h"
+#include "CbmDefs.h"
+#include "CbmEvent.h"
#include <CbmGlobalTrack.h>
-#include <rich/CbmRichRing.h>
+#include <CbmRichRing.h>
#include <FairMCPoint.h>
#include <FairRootManager.h>
@@ -38,13 +40,12 @@ void CbmRichRingsConverter::Init()
rich_branch.AddField<float>("radial_angle", "(0||1||2)*pi +- atan( abs((+-100-y)/-x) )");
auto* man = AnalysisTree::TaskManager::GetInstance();
- man->AddBranch(out_branch_, rich_rings_, rich_branch);
+ man->AddBranch(rich_rings_, rich_branch);
man->AddMatching(match_to_, out_branch_, vtx_tracks_2_rich_);
}
-void CbmRichRingsConverter::FillRichRings()
+void CbmRichRingsConverter::ProcessData(CbmEvent* event)
{
-
assert(cbm_rich_rings_);
rich_rings_->ClearChannels();
vtx_tracks_2_rich_->Clear();
@@ -65,22 +66,27 @@ void CbmRichRingsConverter::FillRichRings()
if (it == indexes_map_->end()) { throw std::runtime_error(match_to_ + " is not found to match with TOF hits"); }
auto rec_tracks_map = it->second;
- rich_rings_->Reserve(cbm_global_tracks_->GetEntriesFast());
+ const int n_tracks = event ? event->GetNofData(ECbmDataType::kGlobalTrack) : cbm_global_tracks_->GetEntriesFast();
+ if (n_tracks <= 0) {
+ LOG(warn) << "No global tracks!";
+ return;
+ }
+ rich_rings_->Reserve(n_tracks);
- for (Int_t igt = 0; igt < cbm_global_tracks_->GetEntriesFast(); igt++) {
- const auto* global_track = static_cast<const CbmGlobalTrack*>(cbm_global_tracks_->At(igt));
+ for (Int_t igt = 0; igt < n_tracks; igt++) {
+ const auto trackIndex = event ? event->GetIndex(ECbmDataType::kGlobalTrack, igt) : igt;
+ const auto* global_track = static_cast<const CbmGlobalTrack*>(cbm_global_tracks_->At(trackIndex));
Int_t i_rich = global_track->GetRichRingIndex();
if (i_rich < 0) continue;
auto rich_ring = static_cast<CbmRichRing*>(cbm_rich_rings_->At(i_rich));
- // auto rich_match = static_cast<CbmTrackMatchNew*>( cbm_rich_rings_->At(i_rich));
- // Int_t itrdMC = (rich_match ? rich_match->GetMatchedLink().GetIndex() : -1);
- // richProj = static_cast<FairTrackParam*>(fRichProjection->At(i));
+ if (rich_ring == nullptr) {
+ LOG(warn) << "No RICH ring!";
+ continue;
+ }
auto& ring = rich_rings_->AddChannel(branch);
-
ring.SetPosition(rich_ring->GetCenterX(), rich_ring->GetCenterY(), 0.f);
-
ring.SetField(int(rich_ring->GetNofHits()), i_n_hits);
ring.SetField(int(rich_ring->GetNofHitsOnRing()), i_n_hits + 1);
ring.SetField(float(rich_ring->GetAaxis()), i_axis);
@@ -101,8 +107,6 @@ void CbmRichRingsConverter::FillRichRings()
}
}
-void CbmRichRingsConverter::Exec() { FillRichRings(); }
-
CbmRichRingsConverter::~CbmRichRingsConverter()
{
delete rich_rings_;
diff --git a/analysis/common/analysis_tree_converter/CbmRichRingsConverter.h b/analysis/common/analysis_tree_converter/CbmRichRingsConverter.h
index 6acf1809f3..5e6c22dd27 100644
--- a/analysis/common/analysis_tree_converter/CbmRichRingsConverter.h
+++ b/analysis/common/analysis_tree_converter/CbmRichRingsConverter.h
@@ -24,12 +24,10 @@ public:
~CbmRichRingsConverter() final;
void Init() final;
- void Exec() final;
+ void ProcessData(CbmEvent* event) final;
void Finish() final {};
private:
- void FillRichRings();
-
TClonesArray* cbm_global_tracks_ {nullptr};
TClonesArray* cbm_rich_rings_ {nullptr};
diff --git a/analysis/common/analysis_tree_converter/CbmSimEventHeaderConverter.cxx b/analysis/common/analysis_tree_converter/CbmSimEventHeaderConverter.cxx
index c5fd77fb55..a62f2d2254 100644
--- a/analysis/common/analysis_tree_converter/CbmSimEventHeaderConverter.cxx
+++ b/analysis/common/analysis_tree_converter/CbmSimEventHeaderConverter.cxx
@@ -4,6 +4,11 @@
#include "CbmSimEventHeaderConverter.h"
+#include "CbmEvent.h"
+#include "CbmMCDataArray.h"
+#include "CbmMCDataManager.h"
+#include "CbmMCDataObject.h"
+
#include "FairMCEventHeader.h"
#include "FairRootManager.h"
@@ -14,40 +19,63 @@
#include "cassert"
#include "iostream"
-ClassImp(CbmSimEventHeaderConverter)
+ClassImp(CbmSimEventHeaderConverter);
- void CbmSimEventHeaderConverter::Init()
+void CbmSimEventHeaderConverter::Init()
{
assert(!out_branch_.empty());
auto* ioman = FairRootManager::Instance();
assert(ioman != nullptr);
- cbm_header_ = (FairMCEventHeader*) ioman->GetObject("MCEventHeader.");
+
+ cbm_mc_manager_ = dynamic_cast<CbmMCDataManager*>(ioman->GetObject("MCDataManager"));
+ cbm_header_obj_ = cbm_mc_manager_->GetObject("MCEventHeader.");
// ***** SimEventHeader *******
AnalysisTree::BranchConfig SimEventHeaderBranch("SimEventHeader", AnalysisTree::DetType::kEventHeader);
SimEventHeaderBranch.AddField<float>("psi_RP", "reaction plane orientation");
SimEventHeaderBranch.AddField<float>("b", "fm, impact parameter");
- SimEventHeaderBranch.AddField<float>("time", "ns, event time");
+ SimEventHeaderBranch.AddFields<float>({"start_time", "end_time"}, "ns (?), event time");
SimEventHeaderBranch.AddField<int>("run_id", "run identifier");
SimEventHeaderBranch.AddField<int>("event_id", "event identifier");
auto* man = AnalysisTree::TaskManager::GetInstance();
- man->AddBranch(out_branch_, sim_event_header_, SimEventHeaderBranch);
+ man->AddBranch(sim_event_header_, SimEventHeaderBranch);
sim_event_header_->Init(SimEventHeaderBranch);
}
-void CbmSimEventHeaderConverter::Exec()
+void CbmSimEventHeaderConverter::ProcessData(CbmEvent* event)
{
- if (!cbm_header_) { throw std::runtime_error("CbmSimEventHeaderConverter::Exec - ERROR! No fHeader!"); }
+ FairMCEventHeader* cbm_header {nullptr};
+
+ if (event) {
+ if (!event->GetMatch()) {
+ LOG(error) << "No match to SimEvent!";
+ return;
+ }
+ const auto& link = event->GetMatch()->GetMatchedLink();
+ cbm_header = (FairMCEventHeader*) (cbm_header_obj_->Get(link));
+ }
+ else {
+ cbm_header = (FairMCEventHeader*) (cbm_header_obj_->Get(0, FairRootManager::Instance()->GetEntryNr()));
+ }
+
+ LOG(info) << "MCEvent " << cbm_header->GetEventID() << " " << cbm_header->GetT();
+
+ if (!cbm_header) { throw std::runtime_error("CbmSimEventHeaderConverter::Exec - ERROR! No fHeader!"); }
auto* out_config_ = AnalysisTree::TaskManager::GetInstance()->GetConfig();
const auto& branch = out_config_->GetBranchConfig(out_branch_);
- TVector3 pos {cbm_header_->GetX(), cbm_header_->GetY(), cbm_header_->GetZ()};
+ TVector3 pos {cbm_header->GetX(), cbm_header->GetY(), cbm_header->GetZ()};
sim_event_header_->SetVertexPosition3(pos);
- sim_event_header_->SetField(float(cbm_header_->GetRotZ()), branch.GetFieldId("psi_RP"));
- sim_event_header_->SetField(float(cbm_header_->GetB()), branch.GetFieldId("b"));
- sim_event_header_->SetField(int(cbm_header_->GetEventID()), branch.GetFieldId("event_id"));
- sim_event_header_->SetField(int(cbm_header_->GetRunID()), branch.GetFieldId("run_id"));
- sim_event_header_->SetField(float(cbm_header_->GetT()), branch.GetFieldId("time"));
+ sim_event_header_->SetField(float(cbm_header->GetRotZ()), branch.GetFieldId("psi_RP"));
+ sim_event_header_->SetField(float(cbm_header->GetB()), branch.GetFieldId("b"));
+ sim_event_header_->SetField(int(cbm_header->GetEventID()), branch.GetFieldId("event_id"));
+ sim_event_header_->SetField(int(cbm_header->GetRunID()), branch.GetFieldId("run_id"));
+
+ if (event) {
+ LOG(info) << "TIME: " << event->GetStartTime() << " " << event->GetEndTime();
+ sim_event_header_->SetField(float(event->GetStartTime()), branch.GetFieldId("start_time"));
+ sim_event_header_->SetField(float(event->GetEndTime()), branch.GetFieldId("end_time"));
+ }
}
diff --git a/analysis/common/analysis_tree_converter/CbmSimEventHeaderConverter.h b/analysis/common/analysis_tree_converter/CbmSimEventHeaderConverter.h
index 35bdcee506..f994f0e663 100644
--- a/analysis/common/analysis_tree_converter/CbmSimEventHeaderConverter.h
+++ b/analysis/common/analysis_tree_converter/CbmSimEventHeaderConverter.h
@@ -11,6 +11,9 @@
class FairMCEventHeader;
class CbmVertex;
+class CbmMCDataManager;
+class CbmMCDataArray;
+class CbmMCDataObject;
class CbmSimEventHeaderConverter final : public CbmConverterTask {
public:
@@ -18,12 +21,14 @@ public:
~CbmSimEventHeaderConverter() final = default;
void Init() final;
- void Exec() final;
+ void ProcessData(CbmEvent* event) final;
void Finish() final { delete sim_event_header_; };
private:
AnalysisTree::EventHeader* sim_event_header_ {nullptr};
FairMCEventHeader* cbm_header_ {nullptr};
+ CbmMCDataManager* cbm_mc_manager_ {nullptr};
+ CbmMCDataObject* cbm_header_obj_ {nullptr};
ClassDef(CbmSimEventHeaderConverter, 1)
};
diff --git a/analysis/common/analysis_tree_converter/CbmSimTracksConverter.cxx b/analysis/common/analysis_tree_converter/CbmSimTracksConverter.cxx
index f175627ab4..37d8d6dcb9 100644
--- a/analysis/common/analysis_tree_converter/CbmSimTracksConverter.cxx
+++ b/analysis/common/analysis_tree_converter/CbmSimTracksConverter.cxx
@@ -4,6 +4,10 @@
#include "CbmSimTracksConverter.h"
+#include "CbmDefs.h"
+#include "CbmEvent.h"
+#include "CbmMCDataArray.h"
+#include "CbmMCDataManager.h"
#include "CbmMCTrack.h"
#include "FairLogger.h"
@@ -25,17 +29,42 @@
#include "UParticle.h"
#include "URun.h"
-ClassImp(CbmSimTracksConverter)
+ClassImp(CbmSimTracksConverter);
- void CbmSimTracksConverter::Init()
+void CbmSimTracksConverter::InitUnigen()
{
+ unigen_file_ = TFile::Open(unigen_file_name_.c_str(), "READ");
+ unigen_file_->Print();
+ if (unigen_file_->IsOpen()) {
+ unigen_tree_ = unigen_file_->Get<TTree>("events");
+ if (unigen_tree_) use_unigen_ = kTRUE;
+ unigen_tree_->SetBranchAddress("event", &unigen_event_);
+ URun* run = unigen_file_->Get<URun>("run");
+ if (run == nullptr) {
+ LOG(error) << "CbmSimTracksConverter: No run description in urqmd file!";
+ delete unigen_file_;
+ unigen_file_ = nullptr;
+ use_unigen_ = kFALSE;
+ }
+ else {
+ Double_t mProt = 0.938272;
+ Double_t pTarg = run->GetPTarg(); // target momentum per nucleon
+ Double_t pProj = run->GetPProj(); // projectile momentum per nucleon
+ Double_t eTarg = TMath::Sqrt(pProj * pProj + mProt * mProt);
+ beta_cm_ = pTarg / eTarg;
+ }
+ }
+}
+void CbmSimTracksConverter::Init()
+{
assert(!out_branch_.empty());
auto* ioman = FairRootManager::Instance();
+ // cbm_mc_tracks_ = (TClonesArray*) ioman->GetObject("MCTrack");
+ cbm_header_ = (FairMCEventHeader*) ioman->GetObject("MCEventHeader.");
-
- cbm_mc_tracks_ = (TClonesArray*) ioman->GetObject("MCTrack");
- cbm_header_ = (FairMCEventHeader*) ioman->GetObject(fEventHeaderBranch);
+ cbm_mc_manager_ = dynamic_cast<CbmMCDataManager*>(ioman->GetObject("MCDataManager"));
+ cbm_mc_tracks_new_ = cbm_mc_manager_->InitBranch("MCTrack");
AnalysisTree::BranchConfig sim_particles_branch(out_branch_, AnalysisTree::DetType::kParticle);
sim_particles_branch.AddField<int>("mother_id", "id of mother particle, -1 for primaries");
@@ -43,84 +72,73 @@ ClassImp(CbmSimTracksConverter)
sim_particles_branch.AddField<int>("geant_process_id", "");
sim_particles_branch.AddFields<int>({"n_hits_mvd", "n_hits_sts", "n_hits_trd"}, "Number of hits in the detector");
- //TDirectory* currentDir = gDirectory;
- if (fUnigenFile.Length() > 0) {
- fFile = new TFile(fUnigenFile, "READ");
- fFile->Print();
- if (fFile->IsOpen()) {
- fTree = fFile->Get<TTree>("events");
- if (fTree) fUseUnigen = kTRUE;
- fTree->SetBranchAddress("event", &fUnigenEvent);
- URun* run = fFile->Get<URun>("run");
- if (run == nullptr) {
- LOG(error) << "CbmSimTracksConverter: No run description in urqmd file!";
- delete fFile;
- fFile = nullptr;
- fUseUnigen = kFALSE;
- }
- else {
- Double_t mProt = 0.938272;
- Double_t pTarg = run->GetPTarg(); // target momentum per nucleon
- Double_t pProj = run->GetPProj(); // projectile momentum per nucleon
- Double_t eTarg = TMath::Sqrt(pProj * pProj + mProt * mProt);
- fBetaCM = pTarg / eTarg;
- }
- }
- }
+ if (!unigen_file_name_.empty()) { InitUnigen(); }
else {
- LOG(info) << "lack of unigen file" << fUnigenFile;
+ LOG(info) << "lack of unigen file" << unigen_file_name_;
}
- if (fUseUnigen) {
- sim_particles_branch.AddField<float>("xfreez", "x freezout coordinate fm/c");
- sim_particles_branch.AddField<float>("yfreez", "y freezout coordinate fm/c");
- sim_particles_branch.AddField<float>("zfreez", "z freezout coordinate fm/c");
- sim_particles_branch.AddField<float>("tfreez", "t freezout coordinate fm/c");
- }
+ sim_particles_branch.AddFields<float>({"start_x", "start_y", "start_z"}, "Start position, cm");
+ sim_particles_branch.AddField<float>("start_t", "t freezout coordinate fm/c");
+
auto* man = AnalysisTree::TaskManager::GetInstance();
- man->AddBranch(out_branch_, sim_tracks_, sim_particles_branch);
+ man->AddBranch(sim_tracks_, sim_particles_branch);
}
-void CbmSimTracksConverter::Exec()
+void CbmSimTracksConverter::ProcessData(CbmEvent* event)
{
- assert(cbm_mc_tracks_);
+ assert(cbm_mc_tracks_new_);
out_indexes_map_.clear();
- if (fUseUnigen) {
- fTree->GetEntry(fEntry++);
- const Double_t unigen_phi = fUnigenEvent->GetPhi();
+ float delta_phi {0.f};
+ if (use_unigen_) {
+ unigen_tree_->GetEntry(entry_++);
+ const Double_t unigen_phi = unigen_event_->GetPhi();
const Double_t mc_phi = cbm_header_->GetRotZ();
- fDeltaPhi = mc_phi - unigen_phi;
+ delta_phi = mc_phi - unigen_phi;
}
+
sim_tracks_->ClearChannels();
auto* out_config_ = AnalysisTree::TaskManager::GetInstance()->GetConfig();
const auto& branch = out_config_->GetBranchConfig(out_branch_);
- const int nMcTracks = cbm_mc_tracks_->GetEntriesFast();
+ int file_id {0}, event_id {0};
+ if (event) {
+ auto match = event->GetMatch();
+ if (!match) return;
+ file_id = event->GetMatch()->GetMatchedLink().GetFile();
+ event_id = event->GetMatch()->GetMatchedLink().GetEntry();
+ }
+ else {
+ event_id = FairRootManager::Instance()->GetEntryNr();
+ }
+
+ LOG(info) << "Writing MC-tracks from event # " << event_id << " file " << file_id;
+
+ const int nMcTracks = cbm_mc_tracks_new_->Size(file_id, event_id);
+
+ if (nMcTracks <= 0) {
+ LOG(warn) << "No MC tracks!";
+ return;
+ }
const int imother_id = branch.GetFieldId("mother_id");
const int igeant_id = branch.GetFieldId("geant_process_id");
const int in_hits = branch.GetFieldId("n_hits_mvd");
const int icbm_id = branch.GetFieldId("cbmroot_id");
+ const int istart_x = branch.GetFieldId("start_x");
sim_tracks_->Reserve(nMcTracks);
- int freezX = 0, freezY = 0, freezZ = 0, freezT = 0;
- if (fUseUnigen) {
- freezX = out_config_->GetBranchConfig(sim_tracks_->GetId()).GetFieldId("xfreez");
- freezY = out_config_->GetBranchConfig(sim_tracks_->GetId()).GetFieldId("yfreez");
- freezZ = out_config_->GetBranchConfig(sim_tracks_->GetId()).GetFieldId("zfreez");
- freezT = out_config_->GetBranchConfig(sim_tracks_->GetId()).GetFieldId("tfreez");
- }
const Double_t nsTofmc = 1. / (0.3356 * 1E-15);
for (int iMcTrack = 0; iMcTrack < nMcTracks; ++iMcTrack) {
- const auto* mctrack = (CbmMCTrack*) cbm_mc_tracks_->At(iMcTrack);
+ const auto trackIndex = iMcTrack; //event ? event->GetIndex(ECbmDataType::kMCTrack, iMcTrack) : iMcTrack;
+ const auto* mctrack = (CbmMCTrack*) cbm_mc_tracks_new_->Get(file_id, event_id, trackIndex);
if (mctrack->GetPdgCode() == 50000050) { //Cherenkov
continue;
}
auto& track = sim_tracks_->AddChannel(branch);
- out_indexes_map_.insert(std::make_pair(iMcTrack, track.GetId()));
+ out_indexes_map_.insert(std::make_pair(trackIndex, track.GetId()));
track.SetMomentum(mctrack->GetPx(), mctrack->GetPy(), mctrack->GetPz());
track.SetMass(float(mctrack->GetMass()));
@@ -131,6 +149,31 @@ void CbmSimTracksConverter::Exec()
track.SetField(int(mctrack->GetNPoints(ECbmModuleId::kTrd)), in_hits + 2);
track.SetField(int(mctrack->GetUniqueID()), icbm_id);
+ if (mctrack->GetMotherId() >= 0) { // secondary
+ track.SetField(float(mctrack->GetStartX() - cbm_header_->GetX()), istart_x);
+ track.SetField(float(mctrack->GetStartY() - cbm_header_->GetY()), istart_x + 1);
+ track.SetField(float(mctrack->GetStartZ() - cbm_header_->GetZ()), istart_x + 2);
+ track.SetField(float(nsTofmc * (mctrack->GetStartT() - cbm_header_->GetT())), istart_x + 3);
+ }
+ else { // primary
+ if (use_unigen_ && trackIndex < unigen_event_->GetNpa()) {
+ UParticle* p = unigen_event_->GetParticle(trackIndex);
+ TLorentzVector boostedX = p->GetPosition();
+ boostedX.Boost(0, 0, -beta_cm_);
+ boostedX.RotateZ(delta_phi);
+ track.SetField(float(boostedX.X() * 1e-13), istart_x);
+ track.SetField(float(boostedX.Y() * 1e-13), istart_x + 1);
+ track.SetField(float(boostedX.Z() * 1e-13), istart_x + 2);
+ track.SetField(float(boostedX.T()), istart_x + 3);
+ }
+ else {
+ track.SetField(0.f, istart_x);
+ track.SetField(0.f, istart_x + 1);
+ track.SetField(0.f, istart_x + 2);
+ track.SetField(0.f, istart_x + 3);
+ }
+ }
+
// mother id should < than track id, so we can do it here
if (mctrack->GetMotherId() == -1) { track.SetField(int(-1), imother_id); }
else {
@@ -141,44 +184,6 @@ void CbmSimTracksConverter::Exec()
track.SetField(int(p->second), imother_id);
}
}
- if (fUseUnigen) {
- if (mctrack->GetMotherId() == -1) { // is primary
- if (iMcTrack >= fUnigenEvent->GetNpa()) { // skip embedded (?) tracks
- track.SetField(float(gRandom->Gaus(0, 200)), freezX);
- track.SetField(float(gRandom->Gaus(0, 200)), freezY);
- track.SetField(float(gRandom->Gaus(0, 200)), freezZ);
- track.SetField(float(gRandom->Gaus(0, 200)), freezT);
- }
- else {
- UParticle* p = fUnigenEvent->GetParticle(iMcTrack);
- TLorentzVector boostedX = p->GetPosition();
- boostedX.Boost(0, 0, -fBetaCM);
- boostedX.RotateZ(fDeltaPhi);
- track.SetField(float(boostedX.X()), freezX);
- track.SetField(float(boostedX.Y()), freezY);
- track.SetField(float(boostedX.Z()), freezZ);
- track.SetField(float(boostedX.T()), freezT);
- }
- }
- else {
- Double_t x = mctrack->GetStartX() - cbm_header_->GetX();
- Double_t y = mctrack->GetStartY() - cbm_header_->GetY();
- Double_t z = mctrack->GetStartZ() - cbm_header_->GetZ();
- Double_t t = mctrack->GetStartT() - cbm_header_->GetT();
-
- track.SetField(float(x * 1.E+13), freezX);
- track.SetField(float(y * 1.E+13), freezY);
- track.SetField(float(z * 1.E+13), freezZ);
- track.SetField(float(t * nsTofmc), freezT);
- }
- }
}
}
-CbmSimTracksConverter::~CbmSimTracksConverter()
-{
- delete sim_tracks_;
- if (fFile) {
- // fFile->Close();
- // delete fFile;
- }
-};
+CbmSimTracksConverter::~CbmSimTracksConverter() { delete sim_tracks_; };
diff --git a/analysis/common/analysis_tree_converter/CbmSimTracksConverter.h b/analysis/common/analysis_tree_converter/CbmSimTracksConverter.h
index 6dba3b8c06..1f27dbf76b 100644
--- a/analysis/common/analysis_tree_converter/CbmSimTracksConverter.h
+++ b/analysis/common/analysis_tree_converter/CbmSimTracksConverter.h
@@ -16,43 +16,44 @@ class UEvent;
class TFile;
class TTree;
class FairMCEventHeader;
+class CbmMCDataManager;
+class CbmMCDataArray;
class CbmSimTracksConverter final : public CbmConverterTask {
-
public:
- explicit CbmSimTracksConverter(std::string out_branch_name, std::string match_to = "", std::string unigen_name = "",
- std::string mc_eventheader = "MCEventHeader.")
- : CbmConverterTask(std::move(out_branch_name), std::move(match_to))
- , fUnigenEvent(nullptr)
- , fFile(nullptr)
- , fTree(nullptr)
- , fUnigenFile(unigen_name)
- , fEventHeaderBranch(mc_eventheader)
- , fEntry(0)
- , fUseUnigen(kFALSE)
- , fDeltaPhi(0)
- , fBetaCM(0) {};
+ explicit CbmSimTracksConverter(std::string out_branch_name, std::string match_to = "")
+ : CbmConverterTask(std::move(out_branch_name), std::move(match_to)) {};
~CbmSimTracksConverter() final;
+ void SetUnigenInfo(const std::string& unigen_name)
+ {
+ use_unigen_ = true;
+ unigen_file_name_ = unigen_name;
+ }
+
void Init() final;
- void Exec() final;
+ void ProcessData(CbmEvent* event) final;
void Finish() final {};
private:
AnalysisTree::Particles* sim_tracks_ {nullptr};
FairMCEventHeader* cbm_header_ {nullptr};
- TClonesArray* cbm_mc_tracks_ {nullptr};
- UEvent* fUnigenEvent;
- TFile* fFile;
- TTree* fTree;
- TString fUnigenFile;
- TString fEventHeaderBranch;
- Int_t fEntry;
- Bool_t fUseUnigen;
- Double_t fDeltaPhi;
- Double_t fBetaCM; // CM velocity in the lab frame
+ // TClonesArray* cbm_mc_tracks_ {nullptr};
+
+ CbmMCDataManager* cbm_mc_manager_ {nullptr};
+ CbmMCDataArray* cbm_mc_tracks_new_ {nullptr};
+
+ void InitUnigen();
+
+ UEvent* unigen_event_ {nullptr};
+ TFile* unigen_file_ {nullptr};
+ TTree* unigen_tree_ {nullptr};
+ std::string unigen_file_name_;
+ Int_t entry_ {0};
+ Double_t beta_cm_ {0}; ///< CM velocity in the lab frame
+ Bool_t use_unigen_ {false};
ClassDef(CbmSimTracksConverter, 1)
};
diff --git a/analysis/common/analysis_tree_converter/CbmStsTracksConverter.cxx b/analysis/common/analysis_tree_converter/CbmStsTracksConverter.cxx
index 1244ac5918..b45f40dbb7 100644
--- a/analysis/common/analysis_tree_converter/CbmStsTracksConverter.cxx
+++ b/analysis/common/analysis_tree_converter/CbmStsTracksConverter.cxx
@@ -4,6 +4,9 @@
#include "CbmStsTracksConverter.h"
+#include "CbmDefs.h"
+#include "CbmEvent.h"
+#include "CbmMCDataManager.h"
#include "CbmMCTrack.h"
#include "CbmStsTrack.h"
#include "CbmTrackMatchNew.h"
@@ -24,17 +27,17 @@
#include "L1Field.h"
-ClassImp(CbmStsTracksConverter)
+ClassImp(CbmStsTracksConverter);
- void CbmStsTracksConverter::Exec()
+void CbmStsTracksConverter::ProcessData(CbmEvent* event)
{
- assert(cbm_sts_tracks_ != nullptr && cbm_mc_tracks_ != nullptr);
+ assert(cbm_sts_tracks_ != nullptr && cbm_mc_tracks_new_ != nullptr);
vtx_tracks_2_sim_->Clear();
out_indexes_map_.clear();
- ReadVertexTracks();
- MapTracks();
+ ReadVertexTracks(event);
+ MapTracks(event);
}
CbmStsTracksConverter::~CbmStsTracksConverter()
@@ -49,8 +52,11 @@ void CbmStsTracksConverter::InitInput()
cbm_prim_vertex_ = (CbmVertex*) ioman->GetObject("PrimaryVertex.");
cbm_sts_tracks_ = (TClonesArray*) ioman->GetObject("StsTrack");
- cbm_mc_tracks_ = (TClonesArray*) ioman->GetObject("MCTrack");
- cbm_sts_match_ = (TClonesArray*) ioman->GetObject("StsTrackMatch");
+ // cbm_mc_tracks_ = (TClonesArray*) ioman->GetObject("MCTrack");
+ cbm_sts_match_ = (TClonesArray*) ioman->GetObject("StsTrackMatch");
+
+ cbm_mc_manager_ = dynamic_cast<CbmMCDataManager*>(ioman->GetObject("MCDataManager"));
+ cbm_mc_tracks_new_ = cbm_mc_manager_->InitBranch("MCTrack");
}
void CbmStsTracksConverter::Init()
@@ -66,6 +72,7 @@ void CbmStsTracksConverter::Init()
vtx_tracks_config.AddField<int>("q", "charge");
vtx_tracks_config.AddField<int>("nhits", "number of hits (total MVD+STS)");
vtx_tracks_config.AddField<int>("nhits_mvd", "number of hits in MVD");
+ vtx_tracks_config.AddField<float>("match_weight", "true over all hits ratio for a matched MC-track");
if (is_write_kfinfo_) {
vtx_tracks_config.AddFields<float>({"x", "y", "z", "tx", "ty", "qp"}, "track parameters");
@@ -75,6 +82,8 @@ void CbmStsTracksConverter::Init()
"cov11", "cov12", "cov13", "cov14", "cov15"},
"covarience matrix");
+ vtx_tracks_config.AddField<float>("dE_over_dx");
+
vtx_tracks_config.AddField<int>("mother_pdg", "PDG code of mother particle");
vtx_tracks_config.AddField<int>("mc_pdg", "MC-true PDG code");
vtx_tracks_config.AddField<bool>("pass_cuts", "ask Oleksii");
@@ -87,7 +96,7 @@ void CbmStsTracksConverter::Init()
ipasscuts_ = vtx_tracks_config.GetFieldId("pass_cuts");
}
auto* man = AnalysisTree::TaskManager::GetInstance();
- man->AddBranch(out_branch_, vtx_tracks_, vtx_tracks_config);
+ man->AddBranch(vtx_tracks_, vtx_tracks_config);
man->AddMatching(out_branch_, match_to_, vtx_tracks_2_sim_);
}
@@ -122,7 +131,7 @@ float CbmStsTracksConverter::ExtrapolateToVertex(CbmStsTrack* sts_track, Analysi
return chi2_to_vtx[0];
}
-void CbmStsTracksConverter::ReadVertexTracks()
+void CbmStsTracksConverter::ReadVertexTracks(CbmEvent* event)
{
assert(cbm_prim_vertex_ && cbm_sts_tracks_);
@@ -137,20 +146,26 @@ void CbmStsTracksConverter::ReadVertexTracks()
const int inhits_mvd = branch.GetFieldId("nhits_mvd");
const int idcax = branch.GetFieldId("dcax");
const int ivtx_chi2 = branch.GetFieldId("vtx_chi2");
+ const int ide_dx = branch.GetFieldId("dE_over_dx");
- const int n_sts_tracks = cbm_sts_tracks_->GetEntries();
+ const int n_sts_tracks = event ? event->GetNofStsTracks() : cbm_sts_tracks_->GetEntries();
+ if (n_sts_tracks <= 0) {
+ LOG(warn) << "No STS tracks!";
+ return;
+ }
vtx_tracks_->Reserve(n_sts_tracks);
for (short i_track = 0; i_track < n_sts_tracks; ++i_track) {
- auto* sts_track = (CbmStsTrack*) cbm_sts_tracks_->At(i_track);
+ const int track_index = event ? event->GetStsTrackIndex(i_track) : i_track;
+ auto* sts_track = (CbmStsTrack*) cbm_sts_tracks_->At(track_index);
if (!sts_track) { throw std::runtime_error("empty track!"); }
auto& track = vtx_tracks_->AddChannel(branch);
- int pdg = GetMcPid((CbmTrackMatchNew*) cbm_sts_match_->At(i_track), track);
+ // int pdg = GetMcPid((CbmTrackMatchNew*) cbm_sts_match_->At(track_index), track);
bool is_good_track = IsGoodCovMatrix(sts_track);
- float chi2_vertex = ExtrapolateToVertex(sts_track, track, pdg);
+ float chi2_vertex = ExtrapolateToVertex(sts_track, track, 211);
const FairTrackParam* trackParamFirst = sts_track->GetParamFirst();
TVector3 momRec;
trackParamFirst->Momentum(momRec);
@@ -166,8 +181,9 @@ void CbmStsTracksConverter::ReadVertexTracks()
track.SetField(float(trackParamFirst->GetZ() - cbm_prim_vertex_->GetZ()), idcax + 2);
track.SetField(int(sts_track->GetNofMvdHits()), inhits_mvd);
track.SetField(float(chi2_vertex), ivtx_chi2);
+ track.SetField(float(sts_track->GetELoss()), ide_dx);
- out_indexes_map_.insert(std::make_pair(i_track, track.GetId()));
+ out_indexes_map_.insert(std::make_pair(track_index, track.GetId()));
if (is_write_kfinfo_) { WriteKFInfo(track, sts_track, is_good_track); }
}
@@ -196,7 +212,6 @@ void CbmStsTracksConverter::WriteKFInfo(AnalysisTree::Track& track, const CbmSts
bool CbmStsTracksConverter::IsGoodCovMatrix(const CbmStsTrack* sts_track) const
{
-
if (!is_reproduce_cbmkfpf_) return true;
assert(sts_track);
@@ -229,52 +244,7 @@ bool CbmStsTracksConverter::IsGoodCovMatrix(const CbmStsTrack* sts_track) const
return ok;
}
-int CbmStsTracksConverter::GetMcPid(const CbmTrackMatchNew* match, AnalysisTree::Track& track) const
-{
-
- if (!is_write_kfinfo_) { return -2; }
- //-----------PID as in MZ's
- //CbmKFPF----------------------------------------------------------
- Int_t nMCTracks = cbm_mc_tracks_->GetEntriesFast();
-
- int pdgCode = -2;
- if (match->GetNofLinks() > 0) { // there is at least one matched MC-track
- Float_t bestWeight = 0.f;
- Float_t totalWeight = 0.f;
- Int_t mcTrackId = -1;
-
- for (int iLink = 0; iLink < match->GetNofLinks(); iLink++) {
- totalWeight += match->GetLink(iLink).GetWeight();
- if (match->GetLink(iLink).GetWeight() > bestWeight) {
- bestWeight = match->GetLink(iLink).GetWeight();
- mcTrackId = match->GetLink(iLink).GetIndex();
- }
- }
-
- if (!((bestWeight / totalWeight < 0.7) || (mcTrackId >= nMCTracks || mcTrackId < 0))) {
- auto* mctrack = static_cast<CbmMCTrack*>(cbm_mc_tracks_->At(mcTrackId));
-
- if (!(TMath::Abs(mctrack->GetPdgCode()) == 11 || TMath::Abs(mctrack->GetPdgCode()) == 13
- || TMath::Abs(mctrack->GetPdgCode()) == 211 || TMath::Abs(mctrack->GetPdgCode()) == 321
- || TMath::Abs(mctrack->GetPdgCode()) == 2212 || TMath::Abs(mctrack->GetPdgCode()) == 1000010020
- || TMath::Abs(mctrack->GetPdgCode()) == 1000010030 || TMath::Abs(mctrack->GetPdgCode()) == 1000020030
- || TMath::Abs(mctrack->GetPdgCode()) == 1000020040)) {
- pdgCode = -1;
- }
- else {
- pdgCode = mctrack->GetPdgCode();
- }
- if (mctrack->GetMotherId() > -1) {
- track.SetField(int(((CbmMCTrack*) cbm_mc_tracks_->At(mctrack->GetMotherId()))->GetPdgCode()), imother_pdg_);
- }
- }
- }
- track.SetField(pdgCode, imc_pdg_);
-
- return pdgCode;
-}
-
-void CbmStsTracksConverter::MapTracks()
+void CbmStsTracksConverter::MapTracks(CbmEvent* event)
{
const auto it = indexes_map_->find(match_to_);
if (it == indexes_map_->end()) { throw std::runtime_error(match_to_ + " is not found to match with vertex tracks"); }
@@ -283,21 +253,86 @@ void CbmStsTracksConverter::MapTracks()
if (out_indexes_map_.empty() || sim_tracks_map.empty()) return;
CbmTrackMatchNew* match {nullptr};
+ int file_id {0}, event_id {0};
+ if (event) {
+ file_id = event->GetMatch()->GetMatchedLink().GetFile();
+ event_id = event->GetMatch()->GetMatchedLink().GetEntry();
+ }
+ else {
+ event_id = FairRootManager::Instance()->GetEntryNr();
+ }
+ auto weight_id = config_->GetBranchConfig(out_branch_).GetFieldId("match_weight");
+
for (const auto& track_id : out_indexes_map_) {
const int cbm_id = track_id.first;
const int out_id = track_id.second;
- match = (CbmTrackMatchNew*) cbm_sts_match_->At(cbm_id);
+ auto& track = vtx_tracks_->Channel(out_id);
+
+ match = (CbmTrackMatchNew*) cbm_sts_match_->At(cbm_id);
if (match->GetNofLinks() > 0) { // there is at least one matched MC-track
- const int mc_track_id = match->GetMatchedLink().GetIndex();
- if (match->GetTrueOverAllHitsRatio() < 0.5) // not enough hits to be a match
+ const auto& link = match->GetMatchedLink();
+
+ if (link.GetFile() != file_id || link.GetEntry() != event_id) { // match from different event
+ LOG(warn) << "MC track from different event";
continue;
+ }
+
+ const int mc_track_id = link.GetIndex();
auto p = sim_tracks_map.find(mc_track_id);
if (p == sim_tracks_map.end()) // match is not found
continue;
+ // LOG(info) << match->GetTrueOverAllHitsRatio();
+
+ track.SetField(float(match->GetTrueOverAllHitsRatio()), weight_id);
vtx_tracks_2_sim_->AddMatch(out_id, p->second);
}
}
}
+
+//int CbmStsTracksConverter::GetMcPid(const CbmTrackMatchNew* match, AnalysisTree::Track& track) const
+//{
+//
+// if (!is_write_kfinfo_) { return -2; }
+// //-----------PID as in MZ's
+// //CbmKFPF----------------------------------------------------------
+// Int_t nMCTracks = cbm_mc_tracks_->GetEntriesFast();
+//
+// int pdgCode = -2;
+// if (match->GetNofLinks() > 0) { // there is at least one matched MC-track
+// Float_t bestWeight = 0.f;
+// Float_t totalWeight = 0.f;
+// Int_t mcTrackId = -1;
+//
+// for (int iLink = 0; iLink < match->GetNofLinks(); iLink++) {
+// totalWeight += match->GetLink(iLink).GetWeight();
+// if (match->GetLink(iLink).GetWeight() > bestWeight) {
+// bestWeight = match->GetLink(iLink).GetWeight();
+// mcTrackId = match->GetLink(iLink).GetIndex();
+// }
+// }
+//
+// if (!((bestWeight / totalWeight < 0.7) || (mcTrackId >= nMCTracks || mcTrackId < 0))) {
+// auto* mctrack = static_cast<CbmMCTrack*>(cbm_mc_tracks_->At(mcTrackId));
+//
+// if (!(TMath::Abs(mctrack->GetPdgCode()) == 11 || TMath::Abs(mctrack->GetPdgCode()) == 13
+// || TMath::Abs(mctrack->GetPdgCode()) == 211 || TMath::Abs(mctrack->GetPdgCode()) == 321
+// || TMath::Abs(mctrack->GetPdgCode()) == 2212 || TMath::Abs(mctrack->GetPdgCode()) == 1000010020
+// || TMath::Abs(mctrack->GetPdgCode()) == 1000010030 || TMath::Abs(mctrack->GetPdgCode()) == 1000020030
+// || TMath::Abs(mctrack->GetPdgCode()) == 1000020040)) {
+// pdgCode = -1;
+// }
+// else {
+// pdgCode = mctrack->GetPdgCode();
+// }
+// if (mctrack->GetMotherId() > -1) {
+// track.SetField(int(((CbmMCTrack*) cbm_mc_tracks_->At(mctrack->GetMotherId()))->GetPdgCode()), imother_pdg_);
+// }
+// }
+// }
+// track.SetField(pdgCode, imc_pdg_);
+//
+// return pdgCode;
+//}
diff --git a/analysis/common/analysis_tree_converter/CbmStsTracksConverter.h b/analysis/common/analysis_tree_converter/CbmStsTracksConverter.h
index f02ab3942a..98b5400bf4 100644
--- a/analysis/common/analysis_tree_converter/CbmStsTracksConverter.h
+++ b/analysis/common/analysis_tree_converter/CbmStsTracksConverter.h
@@ -6,6 +6,7 @@
#define ANALYSIS_TREE_STSTRACKSCONVERTER_H_
#include "CbmConverterTask.h"
+#include <CbmMCDataManager.h>
#include "AnalysisTree/Detector.hpp"
@@ -13,6 +14,8 @@ class TClonesArray;
class CbmVertex;
class CbmStsTrack;
class CbmTrackMatchNew;
+class CbmMCDataManager;
+class CbmMCDataArray;
namespace AnalysisTree
{
@@ -30,30 +33,33 @@ public:
~CbmStsTracksConverter() final;
void Init() final;
- void Exec() final;
+ void ProcessData(CbmEvent* event) final;
void Finish() final {}
void SetIsWriteKFInfo(bool is = true) { is_write_kfinfo_ = is; }
void SetIsReproduceCbmKFPF(bool is = true) { is_reproduce_cbmkfpf_ = is; }
private:
- void ReadVertexTracks();
- void MapTracks();
+ void ReadVertexTracks(CbmEvent* event);
+ void MapTracks(CbmEvent* event);
void InitInput();
float ExtrapolateToVertex(CbmStsTrack* sts_track, AnalysisTree::Track& track, int pdg);
void WriteKFInfo(AnalysisTree::Track& track, const CbmStsTrack* sts_track, bool is_good_track) const;
bool IsGoodCovMatrix(const CbmStsTrack* sts_track) const;
- int GetMcPid(const CbmTrackMatchNew* match, AnalysisTree::Track& track) const;
+ // int GetMcPid(const CbmTrackMatchNew* match, AnalysisTree::Track& track) const;
AnalysisTree::TrackDetector* vtx_tracks_ {nullptr}; ///< raw pointers are needed for TTree::Branch
AnalysisTree::Matching* vtx_tracks_2_sim_ {nullptr}; ///< raw pointers are needed for TTree::Branch
CbmVertex* cbm_prim_vertex_ {nullptr}; ///< non-owning pointer
- TClonesArray* cbm_mc_tracks_ {nullptr}; ///< non-owning pointer
+ // TClonesArray* cbm_mc_tracks_ {nullptr}; ///< non-owning pointer
TClonesArray* cbm_sts_tracks_ {nullptr}; ///< non-owning pointer
TClonesArray* cbm_sts_match_ {nullptr}; ///< non-owning pointer
+ CbmMCDataManager* cbm_mc_manager_ {nullptr};
+ CbmMCDataArray* cbm_mc_tracks_new_ {nullptr};
+
bool is_write_kfinfo_ {true};
bool is_reproduce_cbmkfpf_ {true};
diff --git a/analysis/common/analysis_tree_converter/CbmTofHitsConverter.cxx b/analysis/common/analysis_tree_converter/CbmTofHitsConverter.cxx
index 766cfe94bf..2909165515 100644
--- a/analysis/common/analysis_tree_converter/CbmTofHitsConverter.cxx
+++ b/analysis/common/analysis_tree_converter/CbmTofHitsConverter.cxx
@@ -4,8 +4,11 @@
#include "CbmTofHitsConverter.h"
+#include "CbmDefs.h"
+#include "CbmEvent.h"
+#include "CbmMCDataArray.h"
+#include "CbmMCDataManager.h"
#include <CbmGlobalTrack.h>
-#include <CbmMCTrack.h>
#include <CbmTofHit.h>
#include <CbmTrackMatchNew.h>
@@ -16,7 +19,6 @@
#include <AnalysisTree/TaskManager.hpp>
#include <cassert>
-#include <vector>
#include "AnalysisTree/Matching.hpp"
@@ -30,8 +32,12 @@ void CbmTofHitsConverter::Init()
cbm_tof_hits_ = (TClonesArray*) ioman->GetObject("TofHit");
cbm_global_tracks_ = (TClonesArray*) ioman->GetObject("GlobalTrack");
cbm_tof_match_ = (TClonesArray*) ioman->GetObject("TofHitMatch");
- cbm_tof_points_ = (TClonesArray*) ioman->GetObject("TofPoint");
- cbm_mc_tracks_ = (TClonesArray*) ioman->GetObject("MCTrack");
+ // cbm_tof_points_ = (TClonesArray*) ioman->GetObject("TofPoint");
+ // cbm_mc_tracks_ = (TClonesArray*) ioman->GetObject("MCTrack");
+
+ cbm_mc_manager_ = dynamic_cast<CbmMCDataManager*>(ioman->GetObject("MCDataManager"));
+ cbm_mc_tracks_new_ = cbm_mc_manager_->InitBranch("MCTrack");
+ cbm_tof_points_new_ = cbm_mc_manager_->InitBranch("TofPoint");
AnalysisTree::BranchConfig tof_branch(out_branch_, AnalysisTree::DetType::kHit);
tof_branch.AddField<float>("mass2", "Mass squared");
@@ -42,7 +48,7 @@ void CbmTofHitsConverter::Init()
tof_branch.AddField<int>("mc_pdg", "MC-true PDG code of particle with highest contribution to TOF hit");
auto* man = AnalysisTree::TaskManager::GetInstance();
- man->AddBranch(out_branch_, tof_hits_, tof_branch);
+ man->AddBranch(tof_hits_, tof_branch);
man->AddMatching(match_to_, out_branch_, vtx_tracks_2_tof_);
man->AddMatching(out_branch_, mc_tracks_, tof_hits_2_mc_tracks_);
}
@@ -60,8 +66,7 @@ void CbmTofHitsConverter::ExtrapolateStraightLine(FairTrackParam* params, float
params->SetPosition({x, y, z});
}
-
-void CbmTofHitsConverter::FillTofHits()
+void CbmTofHitsConverter::ProcessData(CbmEvent* event)
{
assert(cbm_tof_hits_);
tof_hits_->ClearChannels();
@@ -80,10 +85,29 @@ void CbmTofHitsConverter::FillTofHits()
auto rec_tracks_map = GetMatchMap(match_to_);
auto sim_tracks_map = GetMatchMap(mc_tracks_);
- tof_hits_->Reserve(cbm_global_tracks_->GetEntriesFast());
+ int file_id {0}, event_id {0};
+ if (event) {
+ auto match = event->GetMatch();
+ if (!match) return;
+ file_id = event->GetMatch()->GetMatchedLink().GetFile();
+ event_id = event->GetMatch()->GetMatchedLink().GetEntry();
+ }
+ else {
+ event_id = FairRootManager::Instance()->GetEntryNr();
+ }
+
+ const int n_tracks = event ? event->GetNofData(ECbmDataType::kGlobalTrack) : cbm_global_tracks_->GetEntriesFast();
+ const int n_tof_hits = event ? event->GetNofData(ECbmDataType::kTofHit) : cbm_tof_hits_->GetEntriesFast();
- for (Int_t igt = 0; igt < cbm_global_tracks_->GetEntriesFast(); igt++) {
- const auto* globalTrack = static_cast<const CbmGlobalTrack*>(cbm_global_tracks_->At(igt));
+ if (n_tracks <= 0) {
+ LOG(warn) << "No global tracks!";
+ return;
+ }
+ tof_hits_->Reserve(n_tracks);
+
+ for (Int_t igt = 0; igt < n_tracks; igt++) {
+ const auto trackIndex = event ? event->GetIndex(ECbmDataType::kGlobalTrack, igt) : igt;
+ const auto* globalTrack = static_cast<const CbmGlobalTrack*>(cbm_global_tracks_->At(trackIndex));
const Int_t tofHitIndex = globalTrack->GetTofHitIndex();
if (tofHitIndex < 0) continue;
@@ -125,7 +149,14 @@ void CbmTofHitsConverter::FillTofHits()
const auto* tofMatch = dynamic_cast<CbmMatch*>(cbm_tof_match_->At(tofHitIndex));
if (tofMatch && tofMatch->GetNofLinks() > 0) {
- const auto* tofPoint = dynamic_cast<FairMCPoint*>(cbm_tof_points_->At(tofMatch->GetMatchedLink().GetIndex()));
+ const auto& link = tofMatch->GetMatchedLink();
+ if (link.GetFile() != file_id || link.GetEntry() != event_id) { // match from different event
+ LOG(warn) << "match from different event";
+ // continue;
+ }
+ const auto* tofPoint = dynamic_cast<FairMCPoint*>(cbm_tof_points_new_->Get(link));
+ // const auto* tofPoint = dynamic_cast<FairMCPoint*>(cbm_tof_points_->At(link.GetIndex()));
+
if (!tofPoint) { throw std::runtime_error("no TOF point"); }
Int_t mc_track_id = tofPoint->GetTrackID();
@@ -139,7 +170,6 @@ void CbmTofHitsConverter::FillTofHits()
}
}
-void CbmTofHitsConverter::Exec() { FillTofHits(); }
CbmTofHitsConverter::~CbmTofHitsConverter()
{
diff --git a/analysis/common/analysis_tree_converter/CbmTofHitsConverter.h b/analysis/common/analysis_tree_converter/CbmTofHitsConverter.h
index 83af254910..cd7dd733e3 100644
--- a/analysis/common/analysis_tree_converter/CbmTofHitsConverter.h
+++ b/analysis/common/analysis_tree_converter/CbmTofHitsConverter.h
@@ -11,6 +11,8 @@
class TClonesArray;
class FairTrackParam;
+class CbmMCDataManager;
+class CbmMCDataArray;
namespace AnalysisTree
{
@@ -27,14 +29,13 @@ public:
~CbmTofHitsConverter() final;
void Init() final;
- void Exec() final;
+ void ProcessData(CbmEvent* event) final;
void Finish() final {}
private:
- void FillTofHits();
static void ExtrapolateStraightLine(FairTrackParam* params, float z);
- const std::map<int, int>& GetMatchMap(std::string name) const
+ const std::map<int, int>& GetMatchMap(const std::string& name) const
{
const auto& it = indexes_map_->find(name);
if (it == indexes_map_->end()) { throw std::runtime_error(name + " is not found to match with TOF hits"); }
@@ -44,9 +45,12 @@ private:
TClonesArray* cbm_global_tracks_ {nullptr};
TClonesArray* cbm_tof_hits_ {nullptr};
- TClonesArray* cbm_tof_points_ {nullptr};
+ // TClonesArray* cbm_tof_points_ {nullptr};
TClonesArray* cbm_tof_match_ {nullptr};
- TClonesArray* cbm_mc_tracks_ {nullptr};
+ // TClonesArray* cbm_mc_tracks_ {nullptr};
+ CbmMCDataManager* cbm_mc_manager_ {nullptr};
+ CbmMCDataArray* cbm_mc_tracks_new_ {nullptr};
+ CbmMCDataArray* cbm_tof_points_new_ {nullptr};
AnalysisTree::HitDetector* tof_hits_ {nullptr};
AnalysisTree::Matching* vtx_tracks_2_tof_ {nullptr};
diff --git a/analysis/common/analysis_tree_converter/CbmTrdTracksConverter.cxx b/analysis/common/analysis_tree_converter/CbmTrdTracksConverter.cxx
index 43fb1d2e07..a29ea1c8a6 100644
--- a/analysis/common/analysis_tree_converter/CbmTrdTracksConverter.cxx
+++ b/analysis/common/analysis_tree_converter/CbmTrdTracksConverter.cxx
@@ -4,8 +4,9 @@
#include "CbmTrdTracksConverter.h"
+#include "CbmDefs.h"
+#include "CbmEvent.h"
#include <CbmGlobalTrack.h>
-#include <CbmTrackMatchNew.h>
#include <CbmTrdHit.h>
#include <CbmTrdTrack.h>
@@ -16,15 +17,14 @@
#include <AnalysisTree/TaskManager.hpp>
#include <cassert>
-#include <vector>
#include "AnalysisTree/Matching.hpp"
-ClassImp(CbmTrdTracksConverter)
+ClassImp(CbmTrdTracksConverter);
- void CbmTrdTracksConverter::Init()
-{
+void CbmTrdTracksConverter::Init()
+{
assert(!out_branch_.empty());
auto* ioman = FairRootManager::Instance();
@@ -42,13 +42,12 @@ ClassImp(CbmTrdTracksConverter)
trd_branch.AddField<int>("n_hits", "Number of hits");
auto* man = AnalysisTree::TaskManager::GetInstance();
- man->AddBranch(out_branch_, trd_tracks_, trd_branch);
+ man->AddBranch(trd_tracks_, trd_branch);
man->AddMatching(match_to_, out_branch_, vtx_tracks_2_trd_);
}
-void CbmTrdTracksConverter::FillTrdTracks()
+void CbmTrdTracksConverter::ProcessData(CbmEvent* event)
{
-
assert(cbm_trd_tracks_);
trd_tracks_->ClearChannels();
vtx_tracks_2_trd_->Clear();
@@ -66,16 +65,27 @@ void CbmTrdTracksConverter::FillTrdTracks()
if (it == indexes_map_->end()) { throw std::runtime_error(match_to_ + " is not found to match with TOF hits"); }
auto rec_tracks_map = it->second;
- trd_tracks_->Reserve(cbm_global_tracks_->GetEntriesFast());
+ const int n_tracks = event ? event->GetNofData(ECbmDataType::kGlobalTrack) : cbm_global_tracks_->GetEntriesFast();
+ if (n_tracks <= 0) {
+ LOG(warn) << "No global tracks!";
+ return;
+ }
+
+ trd_tracks_->Reserve(n_tracks);
- for (Int_t igt = 0; igt < cbm_global_tracks_->GetEntriesFast(); igt++) {
- const auto* global_track = static_cast<const CbmGlobalTrack*>(cbm_global_tracks_->At(igt));
+ for (Int_t igt = 0; igt < n_tracks; igt++) {
+ const auto trackIndex = event ? event->GetIndex(ECbmDataType::kGlobalTrack, igt) : igt;
+
+ const auto* global_track = static_cast<const CbmGlobalTrack*>(cbm_global_tracks_->At(trackIndex));
Int_t itrd = global_track->GetTrdTrackIndex();
if (itrd < 0) continue;
auto trd_track = static_cast<CbmTrdTrack*>(cbm_trd_tracks_->At(itrd));
-
+ if (trd_track == nullptr) {
+ LOG(warn) << "No TRD track!";
+ continue;
+ }
auto& track = trd_tracks_->AddChannel(branch);
TVector3 mom, mom_last;
trd_track->GetParamFirst()->Momentum(mom);
@@ -116,8 +126,6 @@ void CbmTrdTracksConverter::FillTrdTracks()
}
}
-void CbmTrdTracksConverter::Exec() { FillTrdTracks(); }
-
CbmTrdTracksConverter::~CbmTrdTracksConverter()
{
delete trd_tracks_;
diff --git a/analysis/common/analysis_tree_converter/CbmTrdTracksConverter.h b/analysis/common/analysis_tree_converter/CbmTrdTracksConverter.h
index 7bfd5c8e5b..4a54d052cb 100644
--- a/analysis/common/analysis_tree_converter/CbmTrdTracksConverter.h
+++ b/analysis/common/analysis_tree_converter/CbmTrdTracksConverter.h
@@ -18,13 +18,6 @@ namespace AnalysisTree
class CbmTrdTracksConverter final : public CbmConverterTask {
- enum kInBranches
- {
- kTrdTrack = 0,
- kGlobalTrack,
- kTrdHit
- };
-
public:
explicit CbmTrdTracksConverter(std::string out_branch_name, std::string match_to = "")
: CbmConverterTask(std::move(out_branch_name), std::move(match_to))
@@ -34,12 +27,10 @@ public:
~CbmTrdTracksConverter() final;
void Init() final;
- void Exec() final;
+ void ProcessData(CbmEvent* event) final;
void Finish() final {};
private:
- void FillTrdTracks();
-
TClonesArray* cbm_global_tracks_ {nullptr};
TClonesArray* cbm_trd_tracks_ {nullptr};
TClonesArray* cbm_trd_hits_ {nullptr};
diff --git a/external/InstallAnalysisTree.cmake b/external/InstallAnalysisTree.cmake
index 0adf61ee62..43d3e6c20b 100644
--- a/external/InstallAnalysisTree.cmake
+++ b/external/InstallAnalysisTree.cmake
@@ -1,4 +1,4 @@
-set(ANALYSISTREE_VERSION dedb33e2227dfc348aa0fbc862a48f606b09858d) #v2.1.2
+set(ANALYSISTREE_VERSION c338c2b7b47cf47c0d6f2a53871374d2c416bb8a) #v2.2.6
set(ANALYSISTREE_SRC_URL "https://github.com/HeavyIonAnalysis/AnalysisTree.git")
set(ANALYSISTREE_DESTDIR "${CMAKE_BINARY_DIR}/external/ANALYSISTREE-prefix")
diff --git a/external/InstallAnalysisTreeQA.cmake b/external/InstallAnalysisTreeQA.cmake
index eff500d4c1..7328409c89 100644
--- a/external/InstallAnalysisTreeQA.cmake
+++ b/external/InstallAnalysisTreeQA.cmake
@@ -1,4 +1,4 @@
-set(ANALYSISTREEQA_VERSION 76e75a57826c5f4c9ee558328f0e4ab369b9bc78) #v2.0.1
+set(ANALYSISTREEQA_VERSION e3feeef9c50303a6f17a271ff0bcf1cd62b48a4e) #v2.1.0
set(ANALYSISTREEQA_SRC_URL "https://github.com/HeavyIonAnalysis/AnalysisTreeQA.git")
set(ANALYSISTREEQA_DESTDIR "${CMAKE_BINARY_DIR}/external/ANALYSISTREEQA-prefix")
diff --git a/macro/analysis/common/analysis_tree_converter/run_analysis_tree_maker.C b/macro/analysis/common/analysis_tree_converter/run_analysis_tree_maker.C
index 8c61c205db..c20e169781 100644
--- a/macro/analysis/common/analysis_tree_converter/run_analysis_tree_maker.C
+++ b/macro/analysis/common/analysis_tree_converter/run_analysis_tree_maker.C
@@ -3,7 +3,7 @@
Authors: Viktor Klochkov [committer] */
void run_analysis_tree_maker(TString dataSet = "../../../run/test", TString setupName = "sis100_electron",
- TString unigenFile = "")
+ TString unigenFile = "", bool is_event_base = true)
{
const std::string system = "Au+Au"; // TODO can we read it automatically?
const float beam_mom = 12.;
@@ -69,7 +69,7 @@ void run_analysis_tree_maker(TString dataSet = "../../../run/test", TString setu
// ------------------------------------------------------------------------
// ----- Mc Data Manager ------------------------------------------------
- auto* mcManager = new CbmMCDataManager("MCManager", 1);
+ auto* mcManager = new CbmMCDataManager("MCManager", is_event_base);
mcManager->AddFile(traFile);
run->AddTask(mcManager);
// ------------------------------------------------------------------------
@@ -112,7 +112,8 @@ void run_analysis_tree_maker(TString dataSet = "../../../run/test", TString setu
man->SetBeamMomentum(beam_mom);
man->SetOutputName(outFile, "rTree");
- // man->SetOutTreeName("aTree");
+
+ if (!is_event_base) { man->AddTask(new CbmMatchEvents()); }
man->AddTask(new CbmSimEventHeaderConverter("SimEventHeader"));
man->AddTask(new CbmRecEventHeaderConverter("RecEventHeader"));
@@ -126,8 +127,7 @@ void run_analysis_tree_maker(TString dataSet = "../../../run/test", TString setu
man->AddTask(new CbmRichRingsConverter("RichRings", "VtxTracks"));
man->AddTask(new CbmTofHitsConverter("TofHits", "VtxTracks"));
man->AddTask(new CbmTrdTracksConverter("TrdTracks", "VtxTracks"));
- man->AddTask(new CbmPsdModulesConverter("PsdModules"));
-
+ if (is_event_base) { man->AddTask(new CbmPsdModulesConverter("PsdModules")); }
run->AddTask(man);
// ----- Parameter database --------------------------------------------
@@ -179,5 +179,5 @@ void run_analysis_tree_maker(TString dataSet = "../../../run/test", TString setu
// Generate_CTest_Dependency_File(depFile);
// ------------------------------------------------------------------------
- RemoveGeoManager();
+ // RemoveGeoManager();
}
diff --git a/macro/analysis/common/qa/run_analysistree_qa.C b/macro/analysis/common/qa/run_analysistree_qa.C
index 80c1e33a88..ebbe34db0b 100644
--- a/macro/analysis/common/qa/run_analysistree_qa.C
+++ b/macro/analysis/common/qa/run_analysistree_qa.C
@@ -189,6 +189,14 @@ void SimEventHeaderQA(QA::Task& task)
void RecEventHeaderQA(QA::Task& task)
{
+ task.AddH1({"Match weight", {rec_event_header, "match_weight"}, {QA::gNbins, 0, 1.01}});
+
+ task.AddH1({"Start time", {rec_event_header, "start_time"}, {QA::gNbins, 1e6, 1e6}});
+ task.AddH1({"End time", {rec_event_header, "end_time"}, {QA::gNbins, 1e6, 1e6}});
+
+ task.AddH2({"M_{tracks}", {rec_event_header, "M"}, {800, 0, 800}},
+ {"Match weight", {rec_event_header, "match_weight"}, {QA::gNbins, 0, 1.01}});
+
task.AddH1({"x_{vertex} (cm)", {rec_event_header, "vtx_x"}, {QA::gNbins, -1, 1}});
task.AddH1({"y_{vertex} (cm)", {rec_event_header, "vtx_y"}, {QA::gNbins, -1, 1}});
task.AddH1({"z_{vertex} (cm)", {rec_event_header, "vtx_z"}, {QA::gNbins, -1, 1}});
--
GitLab