diff --git a/macro/mcbm/mcbm_digi.C b/macro/mcbm/mcbm_digi.C index 3cc93000e972996b0709dbf243bdbb6a2d7f582e..1ac53fac74654ed3a40dab17fa22ceb17f1fc431 100644 --- a/macro/mcbm/mcbm_digi.C +++ b/macro/mcbm/mcbm_digi.C @@ -68,11 +68,11 @@ void mcbm_digi(Int_t nEvents = 10, // Number of events to p // ----- Digitization run --------------------------------------------- CbmDigitization run; - run.AddInput(inFile, eventRate); + run.AddInput(inFile, cbm::sim::TimeDist::Poisson, eventRate); run.SetOutputFile(outFile, overwrite); run.SetParameterRootFile(parFile); run.SetTimeSliceLength(timeSliceLength); - run.SetEventMode(eventMode); + run.SetMode(eventMode ? cbm::sim::Mode::EventByEvent : cbm::sim::Mode::Timebased); run.Run(nEvents); // ------------------------------------------------------------------------ diff --git a/macro/mvd/qa/mvd_qa3_digitize.C b/macro/mvd/qa/mvd_qa3_digitize.C index 6cb7ef0349f1ec0951a5f0ad3a0527684225b5f3..d2d8816c9cdaa6d1b2948bc058f908bf89f729ec 100644 --- a/macro/mvd/qa/mvd_qa3_digitize.C +++ b/macro/mvd/qa/mvd_qa3_digitize.C @@ -99,7 +99,7 @@ void mvd_qa3_digitize(const char* setup = "sis100_electron") run.AddInput(inFile); run.SetOutputFile(outFile, kTRUE); run.SetParameterRootFile(parFile); - run.SetEventMode(); + run.SetMode(cbm::sim::Mode::EventByEvent); run.SetDigitizer(ECbmModuleId::kMvd, digi); diff --git a/macro/run/run_digi.C b/macro/run/run_digi.C index bede1dc98dffa798b0492370d5bda0c005268f5b..880d1b0f466f0f9b114f56924bf56dc78e69baa4 100644 --- a/macro/run/run_digi.C +++ b/macro/run/run_digi.C @@ -129,16 +129,18 @@ void run_digi(TString inputEvents = "", Int_t nEvents = -1, TString output = "", // ----- Digitization run --------------------------------------------- CbmDigitization run; - run.AddInput(0, evntFile, eventRate); + cbm::sim::Mode mode = (eventRate < 0. ? cbm::sim::Mode::EventByEvent : cbm::sim::Mode::Timebased); + cbm::sim::TimeDist timeDist = cbm::sim::TimeDist::Uniform; + run.AddInput(0, evntFile, timeDist, eventRate); if (!eventMode) { if (!inputSignal.IsNull()) run.EmbedInput(1, signFile, 0); - if (!inputBeam.IsNull()) run.AddInput(2, beamFile, beamRate); + if (!inputBeam.IsNull()) run.AddInput(2, beamFile, timeDist, beamRate); } run.SetOutputFile(outFile, overwrite); run.SetMonitorFile(moniFile); run.SetParameterRootFile(parFile); run.SetTimeSliceLength(tsLength); - run.SetEventMode(eventMode); + run.SetMode(mode); run.SetProduceNoise(kFALSE); if (nEvents < 0) run.Run(); else diff --git a/sim/response/CMakeLists.txt b/sim/response/CMakeLists.txt index dde4f0799b0863782d462edbd99385768a36136d..88992d439d35818f76ef3e927e133e1b5d27d0f2 100644 --- a/sim/response/CMakeLists.txt +++ b/sim/response/CMakeLists.txt @@ -40,3 +40,6 @@ set(INTERFACE_DEPENDENCIES ) generate_cbm_library() + +install(FILES base/Defs.h DESTINATION include/) + diff --git a/sim/response/base/CbmDigitization.cxx b/sim/response/base/CbmDigitization.cxx index 3acf0c1c6e9ad86ba56c77e34aa2fdc78319bc82..9080ba6870b7e02976f5f0cd0ec0e9b1e5e3a34d 100644 --- a/sim/response/base/CbmDigitization.cxx +++ b/sim/response/base/CbmDigitization.cxx @@ -37,28 +37,12 @@ #include <cassert> +using cbm::sim::Mode; +using cbm::sim::TimeDist; + // ----- Constructor ---------------------------------------------------- -CbmDigitization::CbmDigitization() - : TNamed("CbmDigitization", "Digitisation Run") - , fIsInit(kFALSE) - , fEventMode(kFALSE) - , fTimeSliceLength(-1.) - , fProduceNoise(kTRUE) - , fCreateMatches(kTRUE) - , fDigitizers() - , fDaq(new CbmDaq()) - , fSource(new CbmDigitizationSource()) - , fOutFile() - , fParRootFile() - , fParAsciiFiles() - , fMoniFile() - , fOverwriteOutput(kFALSE) - , fGenerateRunInfo(kTRUE) - , fRun(0) -{ - SetDefaultBranches(); -} +CbmDigitization::CbmDigitization() : TNamed("CbmDigitization", "Digitisation Run") { SetDefaultBranches(); } // -------------------------------------------------------------------------- @@ -75,14 +59,14 @@ CbmDigitization::~CbmDigitization() // ----- Add an input file ---------------------------------------------- -void CbmDigitization::AddInput(UInt_t inputId, TString fileName, Double_t eventRate, ECbmTreeAccess mode) +void CbmDigitization::AddInput(UInt_t inputId, TString fileName, TimeDist dist, Double_t eventRate, ECbmTreeAccess mode) { if (gSystem->AccessPathName(fileName)) LOG(fatal) << fName << ": input file " << fileName << " does not exist!"; if (mode != ECbmTreeAccess::kRegular) LOG(fatal) << fName << ": access modes other than kRegular are not yet supported!"; TChain* chain = new TChain("cbmsim"); chain->Add(fileName.Data()); - fSource->AddInput(inputId, chain, eventRate, mode); + fSource->AddInput(inputId, chain, dist, eventRate, mode); } // -------------------------------------------------------------------------- @@ -160,7 +144,7 @@ Int_t CbmDigitization::CreateDefaultDigitizers() if (it->second->GetDigitizer() != nullptr) continue; // --- Skip MVD for time-based mode - if (it->first == ECbmModuleId::kMvd && (!fEventMode)) { + if (it->first == ECbmModuleId::kMvd && fMode == Mode::Timebased) { LOG(info) << "MVD digitizer is not available " << "in time-based mode. "; continue; @@ -362,13 +346,13 @@ void CbmDigitization::Run(Int_t event1, Int_t event2) if (fGenerateRunInfo) LOG(info) << fName << ": Run info will be generated."; // --- Create DAQ - if (fEventMode) fDaq = new CbmDaq(kTRUE); + if (fMode == Mode::EventByEvent) fDaq = new CbmDaq(kTRUE); else fDaq = new CbmDaq(fTimeSliceLength); // --- Register source - if (fEventMode) fSource->SetEventMode(); + fSource->SetMode(fMode); run->SetSource(fSource); @@ -391,7 +375,7 @@ void CbmDigitization::Run(Int_t event1, Int_t event2) CbmDigitizeBase* digitizer = it->second->GetDigitizer(); if (it->second->IsActive() && digitizer != nullptr) { fDaq->SetDigitizer(it->first, digitizer); - if (fEventMode) digitizer->SetEventMode(); + if (fMode == Mode::EventByEvent) digitizer->SetEventMode(); digitizer->SetProduceNoise(fProduceNoise); digitizer->SetCreateMatches(fCreateMatches); run->AddTask(digitizer); @@ -401,7 +385,7 @@ void CbmDigitization::Run(Int_t event1, Int_t event2) // --- In event-by-event mode: also empty events (time-slices) are stored - if (fEventMode) StoreAllTimeSlices(); + if (fMode == Mode::EventByEvent) StoreAllTimeSlices(); // --- Register DAQ diff --git a/sim/response/base/CbmDigitization.h b/sim/response/base/CbmDigitization.h index bea69df2042782525030fdd22ba749393322afad..e2c3cfcd0de0184e3b713403ab50486a06e250c0 100644 --- a/sim/response/base/CbmDigitization.h +++ b/sim/response/base/CbmDigitization.h @@ -25,6 +25,8 @@ #include <map> #include <vector> +#include "Defs.h" + class TGeoManager; class CbmDigitizeBase; @@ -46,8 +48,8 @@ public: ** @param eventRate Rate for events from input file [1/s] ** @param mode Tree access mode (kRegular / kRepeat / kRandom) **/ - void AddInput(UInt_t inputId, TString fileName, Double_t eventRate = -1., - ECbmTreeAccess mode = ECbmTreeAccess::kRegular); + void AddInput(UInt_t inputId, TString fileName, cbm::sim::TimeDist dist = cbm::sim::TimeDist::Poisson, + Double_t eventRate = 0., ECbmTreeAccess mode = ECbmTreeAccess::kRegular); /** @brief Add an input file @@ -58,9 +60,10 @@ public: ** Shortcut for legacy reasons, when only one input file is used. ** This will set the inputId to zero. Repeated use will lead to abort. **/ - void AddInput(TString fileName, Double_t eventRate = -1., ECbmTreeAccess mode = ECbmTreeAccess::kRegular) + void AddInput(TString fileName, cbm::sim::TimeDist dist = cbm::sim::TimeDist::Poisson, Double_t eventRate = 0., + ECbmTreeAccess mode = ECbmTreeAccess::kRegular) { - AddInput(0, fileName, eventRate, mode); + AddInput(0, fileName, dist, eventRate, mode); } /** @brief Add an ASCII parameter file @@ -156,7 +159,7 @@ public: ** In the event-by-event mode, one time slice will be created for ** each input event. There will be no interference between events. **/ - void SetEventMode(Bool_t choice = kTRUE) { fEventMode = choice; } + void SetMode(cbm::sim::Mode mode) { fMode = mode; } /** @brief Set the monitor file name @@ -230,21 +233,22 @@ public: private: - Bool_t fIsInit; - Bool_t fEventMode; - Double_t fTimeSliceLength; - Bool_t fProduceNoise; - Bool_t fCreateMatches; - std::map<ECbmModuleId, CbmDigitizeInfo*> fDigitizers; //! - CbmDaq* fDaq; - CbmDigitizationSource* fSource; ///< Input source - TString fOutFile; ///< Output data (digis) - TString fParRootFile; ///< ROOT parameter file - TList fParAsciiFiles; ///< ASCII parameter files - TString fMoniFile; ///< Resource monitoring information - Bool_t fOverwriteOutput; - Bool_t fGenerateRunInfo; - Int_t fRun; + Bool_t fIsInit = kFALSE; + cbm::sim::Mode fMode = cbm::sim::Mode::Timebased; + Double_t fTimeSliceLength = -1.; + Bool_t fProduceNoise = kTRUE; + ; + Bool_t fCreateMatches = kTRUE; + std::map<ECbmModuleId, CbmDigitizeInfo*> fDigitizers = {}; //! + CbmDaq* fDaq = new CbmDaq(); + CbmDigitizationSource* fSource = new CbmDigitizationSource(); ///< Input source + TString fOutFile = {}; ///< Output data (digis) + TString fParRootFile = {}; ///< ROOT parameter file + TList fParAsciiFiles = {}; ///< ASCII parameter files + TString fMoniFile = {}; ///< Resource monitoring information + Bool_t fOverwriteOutput = kFALSE; + Bool_t fGenerateRunInfo = kTRUE; + Int_t fRun = 0; /** @brief Copy constructor forbidden **/ diff --git a/sim/response/base/CbmDigitizationConfig.cxx b/sim/response/base/CbmDigitizationConfig.cxx index 980771d2327a5ed2cecb66a9567a1a941fd8f7ed..dd4422c4f28df1c292913b4bb45c0220c09cb414 100644 --- a/sim/response/base/CbmDigitizationConfig.cxx +++ b/sim/response/base/CbmDigitizationConfig.cxx @@ -90,7 +90,7 @@ bool CbmDigitizationConfig::SetIO(CbmDigitization& obj, const pt::ptree& moduleT } float rate = pt_input.get<float>("rate", -1.); LOG(info) << "CbmDigitizationConfig: Adding input: " << traFileName; - obj.AddInput(id, path + ".tra.root", rate, treeAccessMode); + obj.AddInput(id, path + ".tra.root", cbm::sim::TimeDist::Poisson, rate, treeAccessMode); } if (parameterSource) { if (!parametersSet) { @@ -134,7 +134,7 @@ bool CbmDigitizationConfig::SetDigitizationParameters(CbmDigitization& obj, cons return false; } } - obj.SetEventMode(eventMode); + obj.SetMode(eventMode ? cbm::sim::Mode::EventByEvent : cbm::sim::Mode::Timebased); if (timeSliceLength) obj.SetTimeSliceLength(timeSliceLength.get()); if (startTime) obj.SetStartTime(startTime.get()); diff --git a/sim/response/base/CbmDigitizationSource.cxx b/sim/response/base/CbmDigitizationSource.cxx index 362a59e4345103fd84a375eceef4dca1b1036500..ca060b36b09cddac612ac4b533398e37a5b8274d 100644 --- a/sim/response/base/CbmDigitizationSource.cxx +++ b/sim/response/base/CbmDigitizationSource.cxx @@ -22,6 +22,9 @@ #include <iomanip> #include <utility> +using cbm::sim::Mode; +using cbm::sim::TimeDist; + // ----- Constructor ----------------------------------------------------- CbmDigitizationSource::CbmDigitizationSource() @@ -38,7 +41,7 @@ CbmDigitizationSource::CbmDigitizationSource() , fCurrentInputId(0) , fCurrentRunId(0) , fFirstCall(kTRUE) - , fEventMode(kFALSE) + , fMode(cbm::sim::Mode::Timebased) , fCurrentInputSet(nullptr) , fSwitchInputSet(kFALSE) { @@ -70,7 +73,7 @@ Bool_t CbmDigitizationSource::ActivateObject(TObject** object, const char* branc // ----- Add a transport input ------------------------------------------- -void CbmDigitizationSource::AddInput(UInt_t inputId, TChain* chain, Double_t rate, ECbmTreeAccess mode) +void CbmDigitizationSource::AddInput(UInt_t inputId, TChain* chain, TimeDist dist, Double_t rate, ECbmTreeAccess mode) { // Catch input ID already being used @@ -80,7 +83,7 @@ void CbmDigitizationSource::AddInput(UInt_t inputId, TChain* chain, Double_t rat } //? Input ID already in use // Create a new inputSet and add the input to it - CbmMCInputSet* inputSet = new CbmMCInputSet(rate); + CbmMCInputSet* inputSet = new CbmMCInputSet(dist, rate); inputSet->AddInput(inputId, chain, mode); // If it is the first input set, it defines the reference branch list. @@ -217,7 +220,7 @@ Bool_t CbmDigitizationSource::Init() LOG(fatal) << "DigitizationSource: No input sets defined!"; return kFALSE; } - if (fEventMode && fInputMap.size() != 1) { + if (fMode == Mode::EventByEvent && fInputMap.size() != 1) { LOG(fatal) << "DigitizationSource: More than one input defined " << "in event-by-event mode!"; return kFALSE; @@ -246,7 +249,7 @@ Bool_t CbmDigitizationSource::Init() ActivateObject(object, "MCEventHeader."); // Set the time of the first event for each input set - if (!fEventMode) { + if (fMode == Mode::Timebased) { for (size_t iSet = 0; iSet < fInputSets.size(); iSet++) { CbmMCInputSet* inputSet = fInputSets.at(iSet); Double_t time = fTimeStart + inputSet->GetDeltaT(); @@ -282,7 +285,7 @@ Int_t CbmDigitizationSource::ReadEvent(UInt_t event) // In the event-by-event mode, get the respective event from the first // input; the event time is fStartTime. - if (fEventMode) return ReadEventByEvent(event); + if (fMode == Mode::EventByEvent) return ReadEventByEvent(event); // If the last used input set was exhausted, switch to a the next one if (fSwitchInputSet) { diff --git a/sim/response/base/CbmDigitizationSource.h b/sim/response/base/CbmDigitizationSource.h index 2a4d3f20e63f38f86cc0e3c08e7707abbc14eb34..17517920ce3984bade26e50a6726805a7d19c527 100644 --- a/sim/response/base/CbmDigitizationSource.h +++ b/sim/response/base/CbmDigitizationSource.h @@ -75,7 +75,8 @@ public: ** @param chain Pointer to input chain ** @param mode Tree access mode (kRegular / kRepeat / kRandom) **/ - void AddInput(UInt_t inputId, TChain* chain, Double_t rate, ECbmTreeAccess mode = ECbmTreeAccess::kRegular); + void AddInput(UInt_t inputId, TChain* chain, cbm::sim::TimeDist dist, Double_t rate, + ECbmTreeAccess mode = ECbmTreeAccess::kRegular); /** @brief Maximal entry number the source can run to @@ -190,7 +191,7 @@ public: ** In the event-by-event mode, only the first input is processed. ** No event start time is generated; the event time is always zero. **/ - void SetEventMode(Bool_t choice = kTRUE) { fEventMode = choice; } + void SetMode(cbm::sim::Mode mode) { fMode = mode; } /** @brief Abstract in base class. No implementation here. @@ -226,7 +227,7 @@ private: Int_t fCurrentInputId; Int_t fCurrentRunId; Bool_t fFirstCall; - Bool_t fEventMode; + cbm::sim::Mode fMode; CbmMCInputSet* fCurrentInputSet; Bool_t fSwitchInputSet; // Flag to switch the input set at next ReadEvent diff --git a/sim/response/base/CbmMCInputSet.cxx b/sim/response/base/CbmMCInputSet.cxx index 3cc0d400ffbaf0e897dbcdbcea13d25b726fa925..74a4be2b68fc1882fc8a771adaa43c7efa2c7839 100644 --- a/sim/response/base/CbmMCInputSet.cxx +++ b/sim/response/base/CbmMCInputSet.cxx @@ -14,24 +14,27 @@ #include <cassert> +using cbm::sim::TimeDist; + // ----- Default constructor --------------------------------------------- -CbmMCInputSet::CbmMCInputSet() : CbmMCInputSet(-1.) {} +CbmMCInputSet::CbmMCInputSet() : CbmMCInputSet(TimeDist::Poisson, -1.) {} // --------------------------------------------------------------------------- // ----- Constructor ----------------------------------------------------- -CbmMCInputSet::CbmMCInputSet(Double_t rate) +CbmMCInputSet::CbmMCInputSet(TimeDist dist, Double_t eventRate) : TObject() - , fRate(rate) + , fTimeDist(dist) + , fEventRate(eventRate) , fInputs() , fInputHandle() , fBranches() , fDeltaDist(nullptr) { - if (rate > 0.) { - Double_t mean = 1.e9 / rate; // mean time between events + if (fTimeDist == TimeDist::Poisson && eventRate > 0.) { + Double_t mean = 1.e9 / eventRate; // mean time between events in ns fDeltaDist = new TF1("DeltaDist", "exp(-x/[0])/[0]", 0., 10. * mean, "NL"); fDeltaDist->SetParameter(0, mean); } @@ -140,8 +143,21 @@ Bool_t CbmMCInputSet::CheckBranchList(CbmMCInput* input) // ----- Time difference to next event ----------------------------------- Double_t CbmMCInputSet::GetDeltaT() { - if (!fDeltaDist) return 0.; - return fDeltaDist->GetRandom(); + + // --- Poisson distribution: delta_t sampled from exponential distribution + if (fTimeDist == TimeDist::Poisson) { + assert(fDeltaDist); + return fDeltaDist->GetRandom(); + } + + // --- Uniform distribution: delta_t is the inverse event rate + else if (fTimeDist == TimeDist::Uniform) { + return 1.e9 / fEventRate; // rate is in 1/s, delta_t in ns + } + + // --- Just in case: catch unknown distribution models + else + return 0.; } // --------------------------------------------------------------------------- diff --git a/sim/response/base/CbmMCInputSet.h b/sim/response/base/CbmMCInputSet.h index ec4960c39cafea2f55d7e6dbff122123e86e27a3..be725502d5f6bbde6e75f45e9cf7a9f8a53955bc 100644 --- a/sim/response/base/CbmMCInputSet.h +++ b/sim/response/base/CbmMCInputSet.h @@ -20,6 +20,8 @@ #include <set> #include <utility> +#include "Defs.h" + /** @class CbmMCInputSet ** @author Volker Friese <v.friese@gsi.de> @@ -47,7 +49,7 @@ public: /** @brief Constructor ** @param rate Event rate [1/s]. Must be positive. **/ - CbmMCInputSet(Double_t rate); + CbmMCInputSet(cbm::sim::TimeDist dist, Double_t rate); /** @brief Destructor **/ @@ -133,7 +135,7 @@ public: /** @brief Event rate ** @value Event rate [1/s] **/ - Double_t GetRate() const { return fRate; } + Double_t GetEventRate() const { return fEventRate; } /** @brief register all input chains to the FairRootManager **/ @@ -141,11 +143,12 @@ public: private: - Double_t fRate; // Event rate [1/s] - std::map<UInt_t, CbmMCInput*> fInputs; // Key is input ID - std::map<UInt_t, CbmMCInput*>::iterator fInputHandle; // Handle for inputs - std::set<TString> fBranches; // List of branch names - TF1* fDeltaDist; // Probability distribution for delta(t) + cbm::sim::TimeDist fTimeDist = cbm::sim::TimeDist::Poisson; // Event time distribution model + Double_t fEventRate; // Event rate [1/s] + std::map<UInt_t, CbmMCInput*> fInputs; // Key is input ID + std::map<UInt_t, CbmMCInput*>::iterator fInputHandle; // Handle for inputs + std::set<TString> fBranches; // List of branch names + TF1* fDeltaDist; // Probability distribution for delta(t) /** @brief Compare an input branch list with the reference branch list diff --git a/sim/response/base/Defs.h b/sim/response/base/Defs.h new file mode 100644 index 0000000000000000000000000000000000000000..0030057a5ba4c3b307b8d657f34d322f33afbcb8 --- /dev/null +++ b/sim/response/base/Defs.h @@ -0,0 +1,37 @@ +/* Copyright (C) 2023 GSI Helmholtzzentrum fuer Schwerionenforschung, Darmstadt + SPDX-License-Identifier: GPL-3.0-only + Authors: Volker Friese [committer] */ + +/** @file CbmMCInput.h + ** @author Volker Friese <v.friese@gsi.de> + ** @date 20.06.2023 + ** + ** Definitions in namespace cbm::sim + **/ + +#ifndef SIM_RESPONSE_BASE_DEFS_H +#define SIM_RESPONSE_BASE_DEFS_H 1 + + +namespace cbm::sim +{ + + /** @enum Simulation mode **/ + enum class Mode + { + Timebased, + EventByEvent + }; + + /** @enum Time distribution model **/ + enum class TimeDist + { + Poisson, + Uniform + }; + + +} // namespace cbm::sim + + +#endif /* SIM_RESPONSE_BASE_DEFS_H */