Browse Source

First commit

PeterParfenov 3 months ago
commit
389da999f9

+ 169 - 0
CMakeLists.txt

@@ -0,0 +1,169 @@
+# CMakeLists.txt for FAS package. It creates a library with dictionary and a main program
+cmake_minimum_required(VERSION 3.1 FATAL_ERROR)
+project(ModelAtConverter)
+
+# You need to tell CMake where to find the ROOT installation. This can be done
+# in a number of ways:
+#   - ROOT built with classic configure/make use the provided
+#   $ROOTSYS/etc/cmake/FindROOT.cmake
+#   - ROOT built with CMake. Add in CMAKE_PREFIX_PATH the installation prefix
+#   for ROOT
+
+find_package(Git)
+
+list(APPEND CMAKE_PREFIX_PATH $ENV{ROOTSYS})
+
+#---Locate the ROOT package and defines a number of variables (e.g. ROOT_INCLUDE_DIRS)
+find_package(ROOT REQUIRED COMPONENTS MathCore RIO Hist Tree Net EG)
+
+#---Define useful ROOT functions and macros (e.g. ROOT_GENERATE_DICTIONARY)
+include(${ROOT_USE_FILE})
+
+add_definitions(${ROOT_CXX_FLAGS})
+
+#---Locate MCINI package
+find_library(MCINI NAMES McIniData PATHS $ENV{MCINI}/build)
+IF(MCINI)
+  message(STATUS "mcini is found: ${MCINI}")
+  add_definitions("-D_MCINI_")
+endif()
+IF(NOT MCINI)
+  message(STATUS "mcini is not found. Building without it.")
+endif()
+
+#---Locate PHQMDEvent library
+find_library(PHQMD NAMES PHQMDEvent PATHS $ENV{PHQMD_PATH}/phqmd2root/lib/)
+if(PHQMD)
+  message(STATUS "libPHQMDEvent is found: ${PHQMD}")
+  add_definitions("-D_PHQMD_")
+  add_definitions("-D_HSD_ROOT_")
+endif()
+if(NOT PHQMD)
+  message(STATUS "libPHQMDEvent is not found. Building without it.")
+endif()
+
+#---Locate AnalysisTree library
+find_library(AnalysisTree NAMES AnalysisTreeBase PATHS $ENV{ANALYSISTREE_LIB})
+if (AnalysisTree)
+  message(STATUS "AnalysisTree is found: ${AnalysisTree}")
+endif()
+if (NOT AnalysisTree)
+  message(FATAL_ERROR "Error: AnalysisTree library wasn't found!")
+endif()
+
+set(CMAKE_BUILD_TYPE Debug)
+#set(CMAKE_C_FLAGS "${CMAKE_C_FLAGS} -O0 -Wall -pthread")
+# set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -O3 -Wall -ffast-math")
+set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -O3 -Wall")
+
+set(MATC_INCLUDE_DIRECTORIES
+  ${CMAKE_CURRENT_SOURCE_DIR}
+  # ${CMAKE_CURRENT_SOURCE_DIR}/bin
+  ${CMAKE_CURRENT_SOURCE_DIR}/format
+  ${CMAKE_CURRENT_SOURCE_DIR}/readers
+  ${CMAKE_CURRENT_SOURCE_DIR}/writers
+  ${CMAKE_CURRENT_SOURCE_DIR}/utility
+  ${ROOT_INLCUDE_DIRS}
+  $ENV{ANALYSISTREE_INC}
+)
+
+set(MATC_INCLUDE_LIBRARIES
+  ${ROOT_LIBRARIES} ${AnalysisTree}
+)
+
+set(MATC_LIBRARY_h_files
+  ${CMAKE_CURRENT_SOURCE_DIR}/format/qaParticle.h
+  ${CMAKE_CURRENT_SOURCE_DIR}/format/qaEvent.h
+  ${CMAKE_CURRENT_SOURCE_DIR}/readers/qaReader_smash_root.h
+  ${CMAKE_CURRENT_SOURCE_DIR}/readers/qaReader_mcpico.h
+  ${CMAKE_CURRENT_SOURCE_DIR}/readers/qaReader_manager.h
+  ${CMAKE_CURRENT_SOURCE_DIR}/writers/qaWriter_manager.h
+  ${CMAKE_CURRENT_SOURCE_DIR}/utility/Utility.h
+)
+
+set(MATC_LIBRARY_cxx_files
+  ${CMAKE_CURRENT_SOURCE_DIR}/format/qaParticle.cxx
+  ${CMAKE_CURRENT_SOURCE_DIR}/format/qaEvent.cxx
+  ${CMAKE_CURRENT_SOURCE_DIR}/readers/qaReader_smash_root.cxx
+  ${CMAKE_CURRENT_SOURCE_DIR}/readers/qaReader_mcpico.cxx
+  ${CMAKE_CURRENT_SOURCE_DIR}/readers/qaReader_manager.cxx
+  ${CMAKE_CURRENT_SOURCE_DIR}/writers/qaWriter_manager.cxx
+  ${CMAKE_CURRENT_SOURCE_DIR}/utility/Utility.cxx
+)
+
+if(MCINI)
+  list(APPEND MATC_INCLUDE_DIRECTORIES $ENV{MCINI}/include)
+  list(APPEND MATC_INCLUDE_LIBRARIES ${MCINI})
+  list(APPEND MATC_LIBRARY_h_files ${CMAKE_CURRENT_SOURCE_DIR}/readers/qaReader_mcini.h)
+  list(APPEND MATC_LIBRARY_cxx_files ${CMAKE_CURRENT_SOURCE_DIR}/readers/qaReader_mcini.cxx)
+endif()
+
+if(PHQMD)
+  list(APPEND MATC_INCLUDE_DIRECTORIES $ENV{PHQMD_PATH})
+  list(APPEND MATC_INCLUDE_LIBRARIES ${PHQMD})
+  list(APPEND MATC_LIBRARY_h_files ${CMAKE_CURRENT_SOURCE_DIR}/readers/qaReader_phqmd.h)
+  list(APPEND MATC_LIBRARY_cxx_files ${CMAKE_CURRENT_SOURCE_DIR}/readers/qaReader_phqmd.cxx)
+  list(APPEND MATC_LIBRARY_h_files ${CMAKE_CURRENT_SOURCE_DIR}/readers/qaReader_hsd_root.h)
+  list(APPEND MATC_LIBRARY_cxx_files ${CMAKE_CURRENT_SOURCE_DIR}/readers/qaReader_hsd_root.cxx)
+endif()
+
+include_directories(${MATC_INCLUDE_DIRECTORIES})
+
+set(MATC_LinkDef
+  ${CMAKE_CURRENT_SOURCE_DIR}/matc.LinkDef.h
+)
+
+#---Generate dictionary
+ROOT_GENERATE_DICTIONARY(G__matc
+  ${MATC_LIBRARY_h_files}
+  LINKDEF ${MATC_LinkDef}
+)
+
+#---Compile library
+add_library(matc SHARED ${MATC_LIBRARY_cxx_files} G__matc.cxx)
+target_link_libraries(matc ${MATC_INCLUDE_LIBRARIES})
+
+# Get the current working branch
+execute_process(
+  COMMAND git rev-parse --abbrev-ref HEAD
+  WORKING_DIRECTORY ${CMAKE_SOURCE_DIR}
+  OUTPUT_VARIABLE GIT_BRANCH
+  RESULT_VARIABLE GIT_BRANCH_ERROR_CODE
+  OUTPUT_STRIP_TRAILING_WHITESPACE
+)
+if( GIT_BRANCH_ERROR_CODE )
+  set(GIT_BRANCH 0)
+endif()
+
+# Get the latest abbreviated commit hash of the working branch
+execute_process(
+  COMMAND git rev-parse --short HEAD
+  WORKING_DIRECTORY ${CMAKE_SOURCE_DIR}
+  OUTPUT_VARIABLE GIT_COMMIT_HASH
+  RESULT_VARIABLE GIT_COMMIT_HASH_ERROR_CODE
+  OUTPUT_STRIP_TRAILING_WHITESPACE
+)
+if( GIT_COMMIT_HASH_ERROR_CODE )
+  set(GIT_COMMIT_HASH 0)
+endif()
+
+# Get the current version
+execute_process(
+  COMMAND git describe --tags --match "*"
+  WORKING_DIRECTORY ${CMAKE_SOURCE_DIR}
+  OUTPUT_VARIABLE GIT_TAG_VERSION
+  RESULT_VARIABLE GIT_TAG_VERSION_ERROR_CODE
+  OUTPUT_STRIP_TRAILING_WHITESPACE
+)
+if( GIT_TAG_VERSION_ERROR_CODE )
+  set(GIT_TAG_VERSION 0.0)
+endif()
+
+set_target_properties (matc 
+  PROPERTIES 
+  VERSION ${GIT_TAG_VERSION} SOVERSION ${GIT_BRANCH}-${GIT_COMMIT_HASH}
+)
+
+#---Compile main executable
+# add_executable(qaProcess "${CMAKE_CURRENT_SOURCE_DIR}/bin/main.cpp")
+# target_link_libraries(qaProcess matc)

+ 0 - 0
SetEnv.sh


+ 11 - 0
format/qaEvent.cxx

@@ -0,0 +1,11 @@
+#include <qaEvent.h>
+
+ClassImp(qaEvent);
+
+qaEvent::qaEvent(/* args */)
+{
+}
+
+qaEvent::~qaEvent()
+{
+}

+ 33 - 0
format/qaEvent.h

@@ -0,0 +1,33 @@
+#ifndef QATOOLS_FORMAT_EVENT_H
+#define QATOOLS_FORMAT_EVENT_H
+
+#include <iostream>
+#include <vector>
+
+#include <Rtypes.h>
+
+class qaEvent
+{
+private:
+  Float_t fB;
+  Float_t fPhiRP;
+  Int_t fNparticles;
+
+public:
+  qaEvent(/* args */);
+  virtual ~qaEvent();
+
+  // Setters
+  virtual void SetB(Float_t _a) { fB = _a; }
+  virtual void SetPhiRP(Float_t _a) { fPhiRP = _a; }
+  virtual void SetNparticles(Float_t _a) { fNparticles = _a; }
+
+  // Getters
+  virtual Float_t GetB() const { return fB; }
+  virtual Float_t GetPhiRP() const { return fPhiRP; }
+  virtual Int_t GetNparticles() const { return fNparticles; }
+
+  ClassDef(qaEvent, 0);
+};
+
+#endif

+ 11 - 0
format/qaParticle.cxx

@@ -0,0 +1,11 @@
+#include <qaParticle.h>
+
+ClassImp(qaParticle);
+
+qaParticle::qaParticle(/* args */)
+{
+}
+
+qaParticle::~qaParticle()
+{
+}

+ 58 - 0
format/qaParticle.h

@@ -0,0 +1,58 @@
+#ifndef QATOOLS_FORMAT_PARTICLE_H
+#define QATOOLS_FORMAT_PARTICLE_H
+
+#include <iostream>
+
+#include <Rtypes.h>
+#include <TVector3.h>
+
+class qaParticle
+{
+private:
+  TVector3 fMom;
+  TVector3 fR;
+  Float_t fEn;
+  Float_t fT;
+  Int_t fPdg;
+  Int_t fCharge;
+
+public:
+  qaParticle(/* args */);
+  virtual ~qaParticle();
+
+  // Setters
+  virtual void SetPxPyPz(Float_t _px, Float_t _py, Float_t _pz) { fMom.SetXYZ(_px, _py, _pz); }
+  virtual void SetXYZ(Float_t _x, Float_t _y, Float_t _z) { fR.SetXYZ(_x, _y, _z); }
+  virtual void SetPx(Float_t _a) { fMom.SetX(_a); }
+  virtual void SetPy(Float_t _a) { fMom.SetY(_a); }
+  virtual void SetPz(Float_t _a) { fMom.SetZ(_a); }
+  virtual void SetX(Float_t _a) { fR.SetX(_a); }
+  virtual void SetY(Float_t _a) { fR.SetY(_a); }
+  virtual void SetZ(Float_t _a) { fR.SetZ(_a); }
+  virtual void SetPdg(Int_t _a) { fPdg = _a; }
+  virtual void SetEnergy(Float_t _a) { fEn = _a; }
+  virtual void SetTime(Float_t _a) { fT = _a; }
+  virtual void SetCharge(Int_t _a) { fCharge = _a; }
+
+  // Getters
+  virtual TVector3 GetMom() { return fMom; }
+  virtual TVector3 GetR() { return fR; }
+  virtual Float_t GetPx() { return fMom.X(); }
+  virtual Float_t GetPy() { return fMom.Y(); }
+  virtual Float_t GetPz() { return fMom.Z(); }
+  virtual Float_t GetPt() { return fMom.Perp(); }
+  virtual Float_t GetP() { return fMom.Mag(); }
+  virtual Float_t GetEta() { return fMom.Eta(); }
+  virtual Float_t GetPhi() { return fMom.Phi(); }
+  virtual Float_t GetT() { return fT; }
+  virtual Float_t GetX() { return fR.X(); }
+  virtual Float_t GetY() { return fR.Y(); }
+  virtual Float_t GetZ() { return fR.Z(); }
+  virtual Int_t GetPdg() const { return fPdg; }
+  virtual Float_t GetEnergy() const { return fEn; }
+  virtual Int_t GetCharge() const { return fCharge; }
+
+  ClassDef(qaParticle, 0);
+};
+
+#endif

+ 28 - 0
matc.LinkDef.h

@@ -0,0 +1,28 @@
+#ifdef __CINT__
+
+#pragma link off all global;
+#pragma link off all classes;
+#pragma link off all functions;
+//#pragma link off all extern;
+
+#pragma link C++ nestedclass;
+#pragma link C++ nestedtypedef;
+
+#pragma link C++ class qaParticle+;
+#pragma link C++ class qaEvent+;
+#pragma link C++ class qaReader_smash_root+;
+#pragma link C++ class qaReader_mcpico+;
+#ifdef _MCINI_
+#pragma link C++ class qaReader_mcini+;
+#endif
+#ifdef _PHQMD_
+#pragma link C++ class qaReader_phqmd+;
+#endif
+#ifdef _HSD_ROOT_
+#pragma link C++ class qaReader_hsd_root+;
+#endif
+#pragma link C++ class qaReader_manager+;
+#pragma link C++ class qaWriter_manager+;
+#pragma link C++ class qaUtility+;
+
+#endif

+ 136 - 0
readers/qaReader_hsd_root.cxx

@@ -0,0 +1,136 @@
+#include <qaReader_hsd_root.h>
+
+ClassImp(qaReader_hsd_root);
+
+qaReader_hsd_root::qaReader_hsd_root(/* args */) : is_init(false), fCurrentEvent(-1), fEvent(nullptr), fParticle(nullptr)
+{
+}
+
+qaReader_hsd_root::~qaReader_hsd_root()
+{
+}
+
+Bool_t qaReader_hsd_root::ChainCheck()
+{
+  if (!is_init)
+  {
+    return false;
+  }
+
+  if (fCurrentEvent == -1)
+  {
+    return false;
+  }
+
+  if (!fChainPHSD->GetEntry(fCurrentEvent))
+  {
+    return false;
+  }
+  if (!fChainFrag->GetEntry(fCurrentEvent))
+  {
+    return false;
+  }
+
+  return true;
+}
+
+void qaReader_hsd_root::SetChain(const TString &inputFileName)
+{
+  fChainPHSD = (TChain *)qaUtility::GetInstance()->initChain(inputFileName, fChainPHSDName);
+
+  fChainPHSD->SetBranchAddress("event", &fEvent);
+
+  fChainFrag = (TChain *)qaUtility::GetInstance()->initChain(inputFileName, fChainFragName);
+
+  fChainFrag->SetBranchAddress("eventFrag", &fEventFrag);
+
+  is_init = kTRUE;
+}
+
+qaEvent *qaReader_hsd_root::ReadEvent(Long64_t iev)
+{
+  fCurrentEvent = iev;
+
+  if (!ChainCheck())
+  {
+    return nullptr;
+  }
+
+  qaEvent *event = new qaEvent();
+
+  event->SetB(fEvent->GetB());
+  event->SetPhiRP(0.); // in PHQMD axis orientation is rotated by 180 deg. In HSD/PHSD it's ok.
+  //event->SetPhiRP(TMath::Pi()); // in PHQMD axis orientation is rotated by 180 deg. In HSD/PHSD it's ok.
+  event->SetNparticles(fEvent->GetNpart() + fEventFrag->GetNfragMST()); // number of particles from PHSD and number of MST fragments
+
+  return event;
+}
+
+qaParticle *qaReader_hsd_root::ReadParticle(Int_t ipart)
+{
+  if (!ChainCheck())
+  {
+    return nullptr;
+  }
+
+  if (ipart >= fEvent->GetNpart() + fEventFrag->GetNfragMST())
+  {
+    return nullptr;
+  }
+
+  if (ipart < 0)
+  {
+    return nullptr;
+  }
+
+  bool is_particle = false;
+  bool is_fragment = false;
+
+  if (ipart < fEvent->GetNpart())
+  {
+    is_particle = true;
+  }
+  if (ipart >= fEvent->GetNpart())
+  {
+    is_fragment = true;
+  }
+
+  if (!is_particle && !is_fragment)
+  {
+    return nullptr;
+  }
+
+  if (is_particle && is_fragment)
+  {
+    return nullptr;
+  }
+
+  qaParticle *particle = new qaParticle();
+
+  if (is_particle)
+  {
+    fParticle = fEvent->GetParticle(ipart);
+    if (fParticle->IsInMST()) // Check whether the particle from PHSD was used to create MST fragment
+    {
+      return nullptr;
+    }
+    particle->SetEnergy(fParticle->E());
+    particle->SetPdg(fParticle->GetPdg());
+    particle->SetPxPyPz(fParticle->Px(), fParticle->Py(), fParticle->Pz());
+    particle->SetTime(0.);
+    particle->SetXYZ(0., 0., 0.);
+    particle->SetCharge(fParticle->GetCharge());
+  }
+  if (is_fragment)
+  {
+    fFragmentMST = fEventFrag->GetFragmentMST(ipart - fEvent->GetNpart());
+    particle->SetEnergy(fFragmentMST->E());
+    particle->SetPdg(fFragmentMST->GetPdg());
+    particle->SetPxPyPz(fFragmentMST->Px(), fFragmentMST->Py(), fFragmentMST->Pz());
+    particle->SetTime(0.);
+    particle->SetXYZ(0., 0., 0.);
+    particle->SetCharge(fFragmentMST->GetZ());
+  }
+
+  return particle;
+}

+ 48 - 0
readers/qaReader_hsd_root.h

@@ -0,0 +1,48 @@
+#ifndef QATOOLS_READERS_HSD_ROOT_H
+#define QATOOLS_READERS_HSD_ROOT_H
+
+#include <Rtypes.h>
+#include <TChain.h>
+#include <TLorentzVector.h>
+
+#include <qaEvent.h>
+#include <qaParticle.h>
+#include <Utility.h>
+
+#include <qaReader_manager.h>
+
+#include <phqmd2root/src/libPHQMDEvent.h>
+
+class qaReader_hsd_root : virtual public qaReader_manager
+{
+private:
+  TChain *fChainPHSD, *fChainFrag;
+
+  const char *fChainPHSDName = "PHQMDtree";
+  const char *fChainFragName = "FRAGtree";
+
+  bool is_init;
+  Long64_t fCurrentEvent;
+
+  Event *fEvent;
+  EventClust *fEventFrag;
+  Particle *fParticle;
+  Baryon *fBaryonMST;
+  Fragment *fFragmentMST;
+  TLorentzVector fMomentum;
+
+  Bool_t ChainCheck();
+
+public:
+  qaReader_hsd_root(/* args */);
+  virtual ~qaReader_hsd_root();
+
+  virtual void SetChain(const TString &inputFileName);
+  virtual Long64_t GetEntries() { return fChainFrag->GetEntries(); }
+  virtual qaEvent *ReadEvent(Long64_t iev);
+  virtual qaParticle* ReadParticle(Int_t ipart);
+
+  ClassDef(qaReader_hsd_root, 0);
+};
+
+#endif

+ 7 - 0
readers/qaReader_manager.cxx

@@ -0,0 +1,7 @@
+#include <qaReader_manager.h>
+
+ClassImp(qaReader_manager);
+
+qaReader_manager::~qaReader_manager()
+{
+}

+ 26 - 0
readers/qaReader_manager.h

@@ -0,0 +1,26 @@
+#ifndef QATOOLS_READERS_MANAGER_H
+#define QATOOLS_READERS_MANAGER_H
+
+#include <Rtypes.h>
+#include <TString.h>
+
+#include <qaEvent.h>
+#include <qaParticle.h>
+#include <Utility.h>
+
+class qaReader_manager
+{
+private:
+  /* data */
+public:
+  virtual ~qaReader_manager();
+
+  virtual void SetChain(const TString &inputFileName) = 0;
+  virtual Long64_t GetEntries() = 0;
+  virtual qaEvent *ReadEvent(Long64_t iev) = 0;
+  virtual qaParticle* ReadParticle(Int_t ipart) = 0;
+
+  ClassDef(qaReader_manager,0);
+};
+
+#endif

+ 85 - 0
readers/qaReader_mcini.cxx

@@ -0,0 +1,85 @@
+#include <qaReader_mcini.h>
+
+ClassImp(qaReader_mcini);
+
+qaReader_mcini::qaReader_mcini(/* args */) : is_init(false), fCurrentEvent(-1), fEvent(nullptr), fParticle(nullptr)
+{
+}
+
+qaReader_mcini::~qaReader_mcini()
+{
+}
+
+Bool_t qaReader_mcini::ChainCheck()
+{
+  if (!is_init)
+  {
+    return false;
+  }
+
+  if (fCurrentEvent == -1)
+  {
+    return false;
+  }
+
+  if (!fChain->GetEntry(fCurrentEvent))
+  {
+    return false;
+  }
+
+  return true;
+}
+
+void qaReader_mcini::SetChain(const TString &inputFileName)
+{
+  fChain = (TChain *)qaUtility::GetInstance()->initChain(inputFileName, fChainName);
+
+  fChain->SetBranchAddress("event", &fEvent);
+
+  is_init = kTRUE;
+}
+
+qaEvent *qaReader_mcini::ReadEvent(Long64_t iev)
+{
+  fCurrentEvent = iev;
+
+  if (!ChainCheck())
+  {
+    return nullptr;
+  }
+
+  qaEvent *event = new qaEvent();
+
+  event->SetB(fEvent->GetB());
+  event->SetPhiRP(fEvent->GetPhi());
+  event->SetNparticles(fEvent->GetNpa());
+
+  return event;
+}
+
+qaParticle *qaReader_mcini::ReadParticle(Int_t ipart)
+{
+  if (!ChainCheck())
+  {
+    return nullptr;
+  }
+
+  if (ipart >= fEvent->GetNpa())
+  {
+    return nullptr;
+  }
+
+  fParticle = fEvent->GetParticle(ipart);
+  fMomentum = fParticle->GetMomentum();
+
+  qaParticle *particle = new qaParticle();
+
+  particle->SetEnergy(fMomentum.E());
+  particle->SetPdg(fParticle->GetPdg());
+  particle->SetPxPyPz(fMomentum.Px(), fMomentum.Py(), fMomentum.Pz());
+  particle->SetTime(0.);
+  particle->SetXYZ(0., 0., 0.);
+  particle->SetCharge(qaUtility::GetInstance()->GetCharge(fParticle->GetPdg()));
+
+  return particle;
+}

+ 47 - 0
readers/qaReader_mcini.h

@@ -0,0 +1,47 @@
+#ifndef QATOOLS_READERS_MCINI_H
+#define QATOOLS_READERS_MCINI_H
+
+#include <Rtypes.h>
+#include <TChain.h>
+#include <TLorentzVector.h>
+
+#include <qaEvent.h>
+#include <qaParticle.h>
+#include <Utility.h>
+
+#include <qaReader_manager.h>
+
+#include <URun.h>
+#include <UEvent.h>
+#include <UParticle.h>
+#include <EventInitialState.h>
+
+class qaReader_mcini : virtual public qaReader_manager
+{
+private:
+  TChain *fChain;
+
+  const char *fChainName = "events";
+
+  bool is_init;
+  Long64_t fCurrentEvent;
+
+  UEvent *fEvent;
+  UParticle *fParticle;
+  TLorentzVector fMomentum;
+
+  Bool_t ChainCheck();
+
+public:
+  qaReader_mcini(/* args */);
+  virtual ~qaReader_mcini();
+
+  virtual void SetChain(const TString &inputFileName);
+  virtual Long64_t GetEntries() { return fChain->GetEntries(); }
+  virtual qaEvent *ReadEvent(Long64_t iev);
+  virtual qaParticle* ReadParticle(Int_t ipart);
+
+  ClassDef(qaReader_mcini, 0);
+};
+
+#endif

+ 95 - 0
readers/qaReader_mcpico.cxx

@@ -0,0 +1,95 @@
+#include <qaReader_mcpico.h>
+
+ClassImp(qaReader_mcpico);
+
+qaReader_mcpico::qaReader_mcpico(/* args */) : is_init(false), fCurrentEvent(-1)
+{
+}
+
+qaReader_mcpico::~qaReader_mcpico()
+{
+}
+
+Bool_t qaReader_mcpico::ChainCheck()
+{
+  if (!is_init)
+  {
+    return false;
+  }
+
+  if (fCurrentEvent == -1)
+  {
+    return false;
+  }
+
+  if (!fChain->GetEntry(fCurrentEvent))
+  {
+    return false;
+  }
+
+  return true;
+}
+
+void qaReader_mcpico::SetChain(const TString &inputFileName)
+{
+  fChain = (TChain *)qaUtility::GetInstance()->initChain(inputFileName, fChainName);
+
+  fChain->SetBranchAddress("bimp", &bimp);
+  fChain->SetBranchAddress("phi2", &phi2);
+  fChain->SetBranchAddress("phi3", &phi3);
+  fChain->SetBranchAddress("ecc2", &ecc2);
+  fChain->SetBranchAddress("ecc3", &ecc3);
+  fChain->SetBranchAddress("npart", &npart);
+  fChain->SetBranchAddress("nh", &nh);
+  fChain->SetBranchAddress("momx", momx);
+  fChain->SetBranchAddress("momy", momy);
+  fChain->SetBranchAddress("momz", momz);
+  fChain->SetBranchAddress("ene", ene);
+  fChain->SetBranchAddress("hid", hid);
+  fChain->SetBranchAddress("pdg", pdg);
+  fChain->SetBranchAddress("charge", charge);
+
+  is_init = kTRUE;
+}
+
+qaEvent *qaReader_mcpico::ReadEvent(Long64_t iev)
+{
+  fCurrentEvent = iev;
+
+  if (!ChainCheck())
+  {
+    return nullptr;
+  }
+
+  qaEvent *event = new qaEvent();
+
+  event->SetB(bimp);
+  event->SetPhiRP(phi2);
+  event->SetNparticles(nh);
+
+  return event;
+}
+
+qaParticle *qaReader_mcpico::ReadParticle(Int_t ipart)
+{
+  if (!ChainCheck())
+  {
+    return nullptr;
+  }
+
+  if (ipart >= nh)
+  {
+    return nullptr;
+  }
+
+  qaParticle *particle = new qaParticle();
+
+  particle->SetEnergy(ene[ipart]);
+  particle->SetPdg(pdg[ipart]);
+  particle->SetPxPyPz(momx[ipart], momy[ipart], momz[ipart]);
+  particle->SetTime(0.);
+  particle->SetXYZ(0., 0., 0.);
+  particle->SetCharge(qaUtility::GetInstance()->GetCharge(pdg[ipart]));
+
+  return particle;
+}

+ 54 - 0
readers/qaReader_mcpico.h

@@ -0,0 +1,54 @@
+#ifndef QATOOLS_READERS_MCPICO_H
+#define QATOOLS_READERS_MCPICO_H
+
+#include <Rtypes.h>
+#include <TChain.h>
+
+#include <qaEvent.h>
+#include <qaParticle.h>
+#include <Utility.h>
+
+#include <qaReader_manager.h>
+
+#define MAX_TRACKS 15000
+
+class qaReader_mcpico : virtual public qaReader_manager
+{
+private:
+  TChain *fChain;
+
+  const char *fChainName = "mctree";
+
+  bool is_init;
+  Long64_t fCurrentEvent;
+
+   Float_t         bimp;
+   Float_t         phi2;
+   Float_t         phi3;
+   Float_t         ecc2;
+   Float_t         ecc3;
+   Int_t           npart;
+   Int_t           nh;
+   Float_t         momx[MAX_TRACKS];   //[nh]
+   Float_t         momy[MAX_TRACKS];   //[nh]
+   Float_t         momz[MAX_TRACKS];   //[nh]
+   Float_t         ene[MAX_TRACKS];   //[nh]
+   Int_t           hid[MAX_TRACKS];   //[nh]
+   Int_t           pdg[MAX_TRACKS];   //[nh]
+   Short_t         charge[MAX_TRACKS];   //[nh]
+
+  Bool_t ChainCheck();
+
+public:
+  qaReader_mcpico(/* args */);
+  virtual ~qaReader_mcpico();
+
+  virtual void SetChain(const TString &inputFileName);
+  virtual Long64_t GetEntries() { return fChain->GetEntries(); }
+  virtual qaEvent *ReadEvent(Long64_t iev);
+  virtual qaParticle* ReadParticle(Int_t ipart);
+
+  ClassDef(qaReader_mcpico, 0);
+};
+
+#endif

+ 136 - 0
readers/qaReader_phqmd.cxx

@@ -0,0 +1,136 @@
+#include <qaReader_phqmd.h>
+
+ClassImp(qaReader_phqmd);
+
+qaReader_phqmd::qaReader_phqmd(/* args */) : is_init(false), fCurrentEvent(-1), fEvent(nullptr), fParticle(nullptr)
+{
+}
+
+qaReader_phqmd::~qaReader_phqmd()
+{
+}
+
+Bool_t qaReader_phqmd::ChainCheck()
+{
+  if (!is_init)
+  {
+    return false;
+  }
+
+  if (fCurrentEvent == -1)
+  {
+    return false;
+  }
+
+  if (!fChainPHSD->GetEntry(fCurrentEvent))
+  {
+    return false;
+  }
+  if (!fChainFrag->GetEntry(fCurrentEvent))
+  {
+    return false;
+  }
+
+  return true;
+}
+
+void qaReader_phqmd::SetChain(const TString &inputFileName)
+{
+  fChainPHSD = (TChain *)qaUtility::GetInstance()->initChain(inputFileName, fChainPHSDName);
+
+  fChainPHSD->SetBranchAddress("event", &fEvent);
+
+  fChainFrag = (TChain *)qaUtility::GetInstance()->initChain(inputFileName, fChainFragName);
+
+  fChainFrag->SetBranchAddress("eventFrag", &fEventFrag);
+
+  is_init = kTRUE;
+}
+
+qaEvent *qaReader_phqmd::ReadEvent(Long64_t iev)
+{
+  fCurrentEvent = iev;
+
+  if (!ChainCheck())
+  {
+    return nullptr;
+  }
+
+  qaEvent *event = new qaEvent();
+
+  event->SetB(fEvent->GetB());
+  // event->SetPhiRP(0.); // in PHQMD axis orientation is rotated by 180 deg. In HSD/PHSD it's ok.
+  event->SetPhiRP(TMath::Pi()); // in PHQMD axis orientation is rotated by 180 deg. In HSD/PHSD it's ok.
+  event->SetNparticles(fEvent->GetNpart() + fEventFrag->GetNfragMST()); // number of particles from PHSD and number of MST fragments
+
+  return event;
+}
+
+qaParticle *qaReader_phqmd::ReadParticle(Int_t ipart)
+{
+  if (!ChainCheck())
+  {
+    return nullptr;
+  }
+
+  if (ipart >= fEvent->GetNpart() + fEventFrag->GetNfragMST())
+  {
+    return nullptr;
+  }
+
+  if (ipart < 0)
+  {
+    return nullptr;
+  }
+
+  bool is_particle = false;
+  bool is_fragment = false;
+
+  if (ipart < fEvent->GetNpart())
+  {
+    is_particle = true;
+  }
+  if (ipart >= fEvent->GetNpart())
+  {
+    is_fragment = true;
+  }
+
+  if (!is_particle && !is_fragment)
+  {
+    return nullptr;
+  }
+
+  if (is_particle && is_fragment)
+  {
+    return nullptr;
+  }
+
+  qaParticle *particle = new qaParticle();
+
+  if (is_particle)
+  {
+    fParticle = fEvent->GetParticle(ipart);
+    if (fParticle->IsInMST()) // Check whether the particle from PHSD was used to create MST fragment
+    {
+      return nullptr;
+    }
+    particle->SetEnergy(fParticle->E());
+    particle->SetPdg(fParticle->GetPdg());
+    particle->SetPxPyPz(fParticle->Px(), fParticle->Py(), fParticle->Pz());
+    particle->SetTime(0.);
+    particle->SetXYZ(0., 0., 0.);
+    particle->SetCharge(fParticle->GetCharge());
+  }
+  if (is_fragment)
+  {
+    fFragmentMST = fEventFrag->GetFragmentMST(ipart - fEvent->GetNpart());
+    particle->SetEnergy(fFragmentMST->E());
+    particle->SetPdg(fFragmentMST->GetPdg());
+    particle->SetPxPyPz(fFragmentMST->Px(), fFragmentMST->Py(), fFragmentMST->Pz());
+    particle->SetTime(0.);
+    particle->SetXYZ(0., 0., 0.);
+    particle->SetCharge(fFragmentMST->GetZ());
+  }
+
+  return particle;
+}

+ 48 - 0
readers/qaReader_phqmd.h

@@ -0,0 +1,48 @@
+#ifndef QATOOLS_READERS_PHQMD_H
+#define QATOOLS_READERS_PHQMD_H
+
+#include <Rtypes.h>
+#include <TChain.h>
+#include <TLorentzVector.h>
+
+#include <qaEvent.h>
+#include <qaParticle.h>
+#include <Utility.h>
+
+#include <qaReader_manager.h>
+
+#include <phqmd2root/src/libPHQMDEvent.h>
+
+class qaReader_phqmd : virtual public qaReader_manager
+{
+private:
+  TChain *fChainPHSD, *fChainFrag;
+
+  const char *fChainPHSDName = "PHQMDtree";
+  const char *fChainFragName = "FRAGtree";
+
+  bool is_init;
+  Long64_t fCurrentEvent;
+
+  Event *fEvent;
+  EventClust *fEventFrag;
+  Particle *fParticle;
+  Baryon *fBaryonMST;
+  Fragment *fFragmentMST;
+  TLorentzVector fMomentum;
+
+  Bool_t ChainCheck();
+
+public:
+  qaReader_phqmd(/* args */);
+  virtual ~qaReader_phqmd();
+
+  virtual void SetChain(const TString &inputFileName);
+  virtual Long64_t GetEntries() { return fChainFrag->GetEntries(); }
+  virtual qaEvent *ReadEvent(Long64_t iev);
+  virtual qaParticle* ReadParticle(Int_t ipart);
+
+  ClassDef(qaReader_phqmd, 0);
+};
+
+#endif

+ 92 - 0
readers/qaReader_smash_root.cxx

@@ -0,0 +1,92 @@
+#include <qaReader_smash_root.h>
+
+ClassImp(qaReader_smash_root);
+
+qaReader_smash_root::qaReader_smash_root(/* args */) : is_init(false), fCurrentEvent(-1)
+{
+}
+
+qaReader_smash_root::~qaReader_smash_root()
+{
+}
+
+Bool_t qaReader_smash_root::ChainCheck()
+{
+  if (!is_init)
+  {
+    return false;
+  }
+
+  if (fCurrentEvent == -1)
+  {
+    return false;
+  }
+
+  if (!fChain->GetEntry(fCurrentEvent))
+  {
+    return false;
+  }
+
+  return true;
+}
+
+void qaReader_smash_root::SetChain(const TString &inputFileName)
+{
+  fChain = (TChain *)qaUtility::GetInstance()->initChain(inputFileName, fChainName);
+
+  fChain->SetBranchAddress("npart", &npart);
+  fChain->SetBranchAddress("impact_b", &impact_b);
+  fChain->SetBranchAddress("pdgcode", pdgcode);
+  fChain->SetBranchAddress("p0", p0);
+  fChain->SetBranchAddress("px", px);
+  fChain->SetBranchAddress("py", py);
+  fChain->SetBranchAddress("pz", pz);
+  fChain->SetBranchAddress("t", t);
+  fChain->SetBranchAddress("x", x);
+  fChain->SetBranchAddress("y", y);
+  fChain->SetBranchAddress("z", z);
+
+  is_init = kTRUE;
+}
+
+qaEvent *qaReader_smash_root::ReadEvent(Long64_t iev)
+{
+  fCurrentEvent = iev;
+
+  if (!ChainCheck())
+  {
+    return nullptr;
+  }
+
+  qaEvent *event = new qaEvent();
+
+  event->SetB(impact_b);
+  event->SetPhiRP(0.);
+  event->SetNparticles(npart);
+
+  return event;
+}
+
+qaParticle *qaReader_smash_root::ReadParticle(Int_t ipart)
+{
+  if (!ChainCheck())
+  {
+    return nullptr;
+  }
+
+  if (ipart >= npart)
+  {
+    return nullptr;
+  }
+
+  qaParticle *particle = new qaParticle();
+
+  particle->SetEnergy(p0[ipart]);
+  particle->SetPdg(pdgcode[ipart]);
+  particle->SetPxPyPz(px[ipart], py[ipart], pz[ipart]);
+  particle->SetTime(t[ipart]);
+  particle->SetXYZ(x[ipart], y[ipart], z[ipart]);
+  particle->SetCharge(qaUtility::GetInstance()->GetCharge(pdgcode[ipart]));
+
+  return particle;
+}

+ 53 - 0
readers/qaReader_smash_root.h

@@ -0,0 +1,53 @@
+#ifndef QATOOLS_READERS_SMASH_ROOT_H
+#define QATOOLS_READERS_SMASH_ROOT_H
+
+#include <Rtypes.h>
+#include <TChain.h>
+
+#include <qaEvent.h>
+#include <qaParticle.h>
+#include <Utility.h>
+
+#include <qaReader_manager.h>
+
+#define MAX_TRACKS 15000
+
+class qaReader_smash_root : virtual public qaReader_manager
+{
+private:
+  TChain *fChain;
+
+  const char *fChainName = "particles";
+
+  bool is_init;
+  Long64_t fCurrentEvent;
+
+  Int_t npart;
+  Double_t impact_b;
+  Int_t ev;
+  Int_t tcounter;
+  Int_t pdgcode[MAX_TRACKS]; //[npart]
+  Double_t p0[MAX_TRACKS];   //[npart]
+  Double_t px[MAX_TRACKS];   //[npart]
+  Double_t py[MAX_TRACKS];   //[npart]
+  Double_t pz[MAX_TRACKS];   //[npart]
+  Double_t t[MAX_TRACKS];    //[npart]
+  Double_t x[MAX_TRACKS];    //[npart]
+  Double_t y[MAX_TRACKS];    //[npart]
+  Double_t z[MAX_TRACKS];    //[npart]
+
+  Bool_t ChainCheck();
+
+public:
+  qaReader_smash_root(/* args */);
+  virtual ~qaReader_smash_root();
+
+  virtual void SetChain(const TString &inputFileName);
+  virtual Long64_t GetEntries() { return fChain->GetEntries(); }
+  virtual qaEvent *ReadEvent(Long64_t iev);
+  virtual qaParticle* ReadParticle(Int_t ipart);
+
+  ClassDef(qaReader_smash_root, 0);
+};
+
+#endif

+ 82 - 0
utility/Utility.cxx

@@ -0,0 +1,82 @@
+#include <Utility.h>
+#include <TDatabasePDG.h>
+#include <TParticlePDG.h>
+
+ClassImp(qaUtility);
+
+qaUtility *qaUtility::fUtility = nullptr;
+;
+
+qaUtility::qaUtility() : Nevents(-1),
+                         debug(0),
+                         format("mctree")
+{
+}
+
+qaUtility::~qaUtility()
+{
+}
+
+qaUtility *qaUtility::GetInstance()
+{
+  if (fUtility == nullptr)
+  {
+    fUtility = new qaUtility();
+  }
+  return fUtility;
+}
+
+TChain *qaUtility::initChain(const TString &inputFileName, const char *chainName)
+{
+  TChain *chain = new TChain(chainName);
+  std::ifstream file(inputFileName.Data());
+  std::string line;
+  while (std::getline(file, line))
+  {
+    chain->Add(line.c_str());
+  }
+
+  return chain;
+}
+
+std::vector<Float_t> qaUtility::ParseVector(std::string _input)
+{
+  std::vector<Float_t> vB;
+
+  std::istringstream iss(_input);
+
+  std::copy(std::istream_iterator<Float_t>(iss),
+            std::istream_iterator<Float_t>(),
+            std::back_inserter(vB));
+
+  return vB;
+}
+
+Double_t qaUtility::GetCharge(Int_t pdg)
+{
+  auto particle = (TParticlePDG *)TDatabasePDG::Instance()->GetParticle(pdg);
+  if (!particle)
+    return error_code;
+  return particle->Charge() / 3.;
+}
+
+Int_t qaUtility::GetPdgId(Int_t pdg)
+{
+  if (pdg == vpdg.at(1))
+    return 1;
+  if (pdg == vpdg.at(2))
+    return 2;
+  if (pdg == vpdg.at(3))
+    return 3;
+  if (pdg == vpdg.at(4))
+    return 4;
+  if (pdg == vpdg.at(6))
+    return 6;
+  if (pdg == vpdg.at(7))
+    return 7;
+  if (pdg == vpdg.at(8))
+    return 8;
+  if (pdg == vpdg.at(9))
+    return 9;
+  return -1;
+}

+ 55 - 0
utility/Utility.h

@@ -0,0 +1,55 @@
+#ifndef QATOOLS_UTILITY_H
+#define QATOOLS_UTILITY_H
+
+#include <iostream>
+#include <fstream>
+#include <vector>
+#include <algorithm>
+#include <iterator>
+#include <cassert>
+#include <sstream>
+#include <string>
+
+#include <Rtypes.h>
+#include <TString.h>
+#include <TChain.h>
+#include <TEnv.h>
+
+#include <qaEvent.h>
+#include <qaParticle.h>
+
+class qaUtility
+{
+protected:
+  qaUtility();
+  virtual ~qaUtility();
+  static qaUtility *fUtility;
+
+public:
+  qaUtility(qaUtility &other) = delete;
+  void operator=(const qaUtility &) = delete;
+
+  static qaUtility *GetInstance();
+
+  const Double_t error_code = -999.;
+  const Int_t npid = 10;
+  const std::vector<Int_t> vpdg = {0, 211, 321, 2212, 2112, 0, -211, -321, -2212, -2112};
+  const std::vector<Double_t> mpdg = {error_code, 0.13957, 0.493677, 0.938272, 0.88277299, error_code, 0.13957, 0.493677, 0.938272, 0.88277299};
+  const Int_t maxCentBins = 20;
+
+  Int_t Nevents;
+  Int_t debug;
+  std::string format;
+
+  TChain *initChain(const TString &inputFileName, const char *chainName);
+
+  std::vector<Float_t> ParseVector(std::string _input);
+  Bool_t initCentrality();
+
+  Double_t GetCharge(Int_t pdg);
+  Int_t GetPdgId(Int_t pdg);
+
+  ClassDef(qaUtility, 0);
+}; // class qaUtility
+
+#endif

+ 7 - 0
writers/qaWriter_manager.cxx

@@ -0,0 +1,7 @@
+#include <qaWriter_manager.h>
+
+ClassImp(qaWriter_manager);
+
+qaWriter_manager::~qaWriter_manager()
+{
+}

+ 26 - 0
writers/qaWriter_manager.h

@@ -0,0 +1,26 @@
+#ifndef QATOOLS_WRITERS_MANAGER_H
+#define QATOOLS_WRITERS_MANAGER_H
+
+#include <Rtypes.h>
+#include <TString.h>
+
+#include <qaEvent.h>
+#include <qaParticle.h>
+#include <Utility.h>
+
+class qaWriter_manager
+{
+private:
+  /* data */
+public:
+  virtual ~qaWriter_manager();
+
+  virtual void Init(std::string filename, std::string treename) = 0;
+  virtual void WriteEvent(qaEvent* event) = 0;
+  virtual void WriteParticle(qaParticle* particle) = 0;
+  virtual void WriteTree() = 0;
+
+  ClassDef(qaWriter_manager,0);
+};
+
+#endif