diff --git a/DDCore/python/dd4hep_base.py b/DDCore/python/dd4hep_base.py index e2606b70951b169be1382266372f0bbad1b047b3..fcee673e9023d6d79867ec9eb309a2eb5fbf5668 100644 --- a/DDCore/python/dd4hep_base.py +++ b/DDCore/python/dd4hep_base.py @@ -326,6 +326,10 @@ class Logger: "Logger constructor" self.name = name + def setPrintLevel(self, level): + "Adjust printout level of dd4hep" + dd4hep.setPrintLevel(level) + def always(self, msg): "Call dd4hep printout function with level ALWAYS" dd4hep.always(self.name, msg) diff --git a/DDDigi/CMakeLists.txt b/DDDigi/CMakeLists.txt index 151afca05c8e8c4349091be4c56e95ef9c6f961c..f8bd4ce3c025db7ef78235cfd8ae78484db54c7f 100644 --- a/DDDigi/CMakeLists.txt +++ b/DDDigi/CMakeLists.txt @@ -45,24 +45,44 @@ dd4hep_add_plugin(DDDigiPlugins GENERATED G__DDDigi.cxx USES DD4hep::DDDigi ) -#--------------------------- Plugin library to read input from DDG4 -------------- -if (DD4HEP_USE_GEANT4) +# +#---------------- I/O Plugin library to read/write DDG4 and edm4hep files --------- +set(DDDigiIO_USES "DD4hep::DDDigi;DD4hep::DDCore") +set(DDDigiIO_SOURCES "io/DigiIO.cpp") +set(DDDigiIO_GENERATED) +set(DDDigiIO_DEFINITIONS) +if(DD4HEP_USE_GEANT4) dd4hep_add_dictionary(G__DDDigi_DDG4_IO - SOURCES ../DDCore/include/ROOT/Warnings.h ddg4/IO.cpp + SOURCES ../DDCore/include/ROOT/Warnings.h io/DDG4IO.cpp LINKDEF ../DDCore/include/ROOT/LinkDef.h USES DD4hep::DDG4 DD4hep::DDCore ) - dd4hep_add_plugin(DDDigi_DDG4_IO - SOURCES ddg4/*.cpp - GENERATED G__DDDigi_DDG4_IO.cxx - USES DD4hep::DDDigi DD4hep::DDG4 DD4hep::DDCore - ) - set_target_properties(DDDigi_DDG4_IO PROPERTIES VERSION ${DD4hep_VERSION} SOVERSION ${DD4hep_SOVERSION}) - install(TARGETS DDDigi_DDG4_IO EXPORT DD4hep ARCHIVE DESTINATION lib LIBRARY DESTINATION lib) + list(APPEND DDDigiIO_DEFINITIONS DD4HEP_USE_DDG4=1) + list(APPEND DDDigiIO_GENERATED G__DDDigi_DDG4_IO.cxx) + list(APPEND DDDigiIO_SOURCES "io/DDG4IO.cpp;io/DigiDDG4Input.cpp") + list(APPEND DDDigiIO_USES DD4hep::DDG4) else() dd4hep_print( "|++> Geant4 not used. DDDigi will not be able to read DDG4 output.") endif() # +if(DD4HEP_USE_EDM4HEP) + list(APPEND DDDigiIO_SOURCES io/Digi2edm4hep.cpp) + list(APPEND DDDigiIO_DEFINITIONS DD4HEP_USE_EDM4HEP=1) + list(APPEND DDDigiIO_USES "EDM4HEP::edm4hep;EDM4HEP::edm4hepDict;podio::podio;podio::podioDict;podio::podioRootIO") +else() + dd4hep_print( "|++> EDM4HEP not used. DDDigi will not be able to write EDM4HEP output.") +endif() +# +list(REMOVE_DUPLICATES DDDigiIO_USES) +dd4hep_add_plugin(DDDigi_IO + SOURCES ${DDDigiIO_SOURCES} + GENERATED ${DDDigiIO_GENERATED} + USES ${DDDigiIO_USES} + ) + target_compile_definitions(DDDigi_IO PRIVATE ${DDDigiIO_DEFINITIONS}) + set_target_properties(DDDigi_IO PROPERTIES VERSION ${DD4hep_VERSION} SOVERSION ${DD4hep_SOVERSION}) + install(TARGETS DDDigi_IO EXPORT DD4hep ARCHIVE DESTINATION lib LIBRARY DESTINATION lib) +# #---Package installation procedure(s) ---------------------------------------------- set_target_properties(DDDigi DDDigiPlugins PROPERTIES VERSION ${DD4hep_VERSION} SOVERSION ${DD4hep_SOVERSION}) # diff --git a/DDDigi/ddg4/DigiDDG4Input.cpp b/DDDigi/ddg4/DigiDDG4Input.cpp deleted file mode 100644 index ca56883f364838f96c19260649b3a6f7c50fcc0c..0000000000000000000000000000000000000000 --- a/DDDigi/ddg4/DigiDDG4Input.cpp +++ /dev/null @@ -1,191 +0,0 @@ -//========================================================================== -// AIDA Detector description implementation -//-------------------------------------------------------------------------- -// Copyright (C) Organisation europeenne pour la Recherche nucleaire (CERN) -// All rights reserved. -// -// For the licensing terms see $DD4hepINSTALL/LICENSE. -// For the list of contributors see $DD4hepINSTALL/doc/CREDITS. -// -// Author : M.Frank -// -//========================================================================== - -// Framework include files -#include <DD4hep/DD4hepUnits.h> - -#include <DDG4/Geant4Data.h> -#include <DDG4/Geant4Particle.h> - -#include <DDDigi/DigiData.h> -#include <DDDigi/DigiROOTInput.h> -#include <DDDigi/DigiFactories.h> - -// ROOT include files -#include <TBranch.h> - -// C/C++ include files -#include <stdexcept> -#include <any> - -/// Namespace for the AIDA detector description toolkit -namespace dd4hep { - - /// Namespace for the Digitization part of the AIDA detector description toolkit - namespace digi { - - class DigiDDG4ROOT : public DigiROOTInput { - public: - static constexpr double epsilon = std::numeric_limits<double>::epsilon(); - - /// Property to generate extra history records - bool m_keep_raw { false }; - - mutable TClass* m_trackerHitClass { nullptr }; - mutable TClass* m_caloHitClass { nullptr }; - mutable TClass* m_particlesClass { nullptr }; - - template <typename T> union input_data { - const void* raw; - std::vector<T*>* items; - input_data(const void* p) { this->raw = p; } - }; - - TClass* get_class_pointers(TClass* c , const char* nam) const { - if ( strstr(nam, "vector<dd4hep::sim::Geant4Particle*") ) - return m_particlesClass = c; - if ( strstr(nam, "vector<dd4hep::sim::Geant4Tracker::Hit*") ) - return m_trackerHitClass = c; - if ( strstr(nam, "vector<dd4hep::sim::Geant4Calorimeter::Hit*") ) - return m_caloHitClass = c; - except("Unknown data class: %s Cannot be handled", nam); - return c; - } - - void add_particle_history(const sim::Geant4Calorimeter::Hit* hit, Key key, History& hist) const { - for( const auto& truth : hit->truth ) { - key.set_item(truth.trackID); - hist.particles.emplace_back(key, truth.deposit); - } - } - void add_particle_history(const sim::Geant4Tracker::Hit* hit, Key key, History& hist) const { - key.set_item(hit->truth.trackID); - hist.particles.emplace_back(key, hit->truth.deposit); - } - - public: - /// Initializing constructor - DigiDDG4ROOT(const DigiKernel& krnl, const std::string& nam) - : DigiROOTInput(krnl, nam) - { - declareProperty("keep_raw", m_keep_raw); - } - - template <typename T> - void conv_hits(DigiContext& context, DataSegment& segment, - const std::string& tag, Key::mask_type mask, const char* name, void* ptr) const - { - using wrap_t = std::shared_ptr<sim::Geant4HitData>; - std::size_t len = 0; - std::string nam = name; - std::vector<wrap_t> geant4_hits; - DepositVector out(nam, mask); - if ( ptr ) { - input_data<T> data(ptr); - for(size_t i=0; i < data.items->size(); ++i) { - auto* p = (*data.items)[i]; - if ( p->energyDeposit > epsilon ) { - Key history_key; - CellID cell = p->cellID; - wrap_t raw(p); - EnergyDeposit dep { }; - Position pos = p->position; - pos *= 1./dd4hep::mm; - - dep.flag = p->flag; - dep.deposit = p->energyDeposit; - dep.position = pos; - - history_key.set_mask(mask); - history_key.set_item(out.size()); - history_key.set_segment(segment.id); - dep.history.hits.emplace_back(history_key, dep.deposit); - add_particle_history(p, history_key, dep.history); - out.emplace(cell, std::move(dep)); - if ( m_keep_raw ) { - geant4_hits.emplace_back(std::move(raw)); - } - } - } - len = data.items->size(); - data.items->clear(); - } - info("%s++ %-24s Converted %6ld DDG4 %-14s hits to %6ld cell deposits", - context.event->id(), name, len, tag.c_str(), out.size()); - Key depo_key(out.name, mask); - segment.put(depo_key, std::move(out)); - if ( m_keep_raw ) { - Key src_key(nam+".ddg4", mask); - segment.emplace_any(src_key, std::make_any<std::vector<wrap_t> >(std::move(geant4_hits))); - } - } - - void conv_particles(DigiContext& context, DataSegment& segment, - int mask, const char* name, void* ptr) const - { - using wrap_t = std::shared_ptr<sim::Geant4Particle>; - ParticleMapping out(name, mask); - std::map<Key, wrap_t> source; - Key part_key(name, mask); - if ( ptr ) { - auto* items = (std::vector<sim::Geant4Particle*>*)ptr; - Key key; - key.set_segment(segment.id); - key.set_mask(mask); - for( auto* p : *items ) { - Particle part; - p->mask = mask; - part.start_position = Position(p->vsx, p->vsy, p->vsz); - part.end_position = Position(p->vex, p->vey, p->vez); - part.momentum = Direction(p->psx, p->psy, p->psz); - part.pdgID = p->pdgID; - part.charge = p->charge; - part.mass = p->mass; - part.history = key.set_item(part_key.item()); - key.set_item(out.size()); - out.push(key, std::move(part)); - source[key] = wrap_t(p); - } - items->clear(); - } - debug("%s++ Converted %ld DDG4 particles", context.event->id(), out.size()); - std::string nam = name; - segment.put(part_key, std::move(out)); - if ( m_keep_raw ) { - Key src_key(nam+".ddg4", mask); - segment.emplace_any(src_key, std::make_any<std::map<Key,wrap_t> >(std::move(source))); - } - } - - /// Callback to handle single branch - virtual void operator()(DigiContext& context, work_t& work) const override { - TBranch& br = work.branch; - TClass& cl = work.cl; - TClass* c = get_class_pointers(&cl, br.GetClassName()); - void** addr = (void**)br.GetAddress(); - const char* nam = br.GetName(); - int mask = work.input_key.mask(); - auto& seg = work.input_segment; - if ( c == m_caloHitClass ) - conv_hits<sim::Geant4Calorimeter::Hit>(context, seg, "calorimeter", mask, nam, *addr); - else if ( c == m_trackerHitClass ) - conv_hits<sim::Geant4Tracker::Hit>(context, seg, "tracker", mask, nam, *addr); - else if ( c == m_particlesClass ) - conv_particles(context, seg, mask, nam, *addr); - else - except("Unknown data type encountered in branch: %s", nam); - } - }; - } // End namespace digi -} // End namespace dd4hep -DECLARE_DIGIACTION_NS(dd4hep::digi,DigiDDG4ROOT) diff --git a/DDDigi/include/DDDigi/DigiContainerProcessor.h b/DDDigi/include/DDDigi/DigiContainerProcessor.h index 78ccab793bd5ec0afbeea758e4ded6db9e50701a..9f6bf290235a974e227c2ef18e6d015004346c64 100644 --- a/DDDigi/include/DDDigi/DigiContainerProcessor.h +++ b/DDDigi/include/DDDigi/DigiContainerProcessor.h @@ -136,7 +136,7 @@ namespace dd4hep { /// Default destructor virtual ~DigiContainerProcessor(); /// Adopt monitoring action - void adopt_monitor(DigiDepositMonitor* monitor); + virtual void adopt_monitor(DigiDepositMonitor* monitor); /// Main functional callback adapter virtual void execute(context_t& context, work_t& work, const predicate_t& predicate) const; }; @@ -214,7 +214,7 @@ namespace dd4hep { /// Standard constructor DigiContainerSequence(const kernel_t& kernel, const std::string& name); /// Set the default predicate - void set_predicate(const predicate_t& predicate); + virtual void set_predicate(const predicate_t& predicate); /// Adopt new parallel worker virtual void adopt_processor(DigiContainerProcessor* action); /// Main functional callback adapter @@ -251,9 +251,9 @@ namespace dd4hep { using self_t = DigiContainerSequenceAction; /// Work definition structure: Argument structure for client calls struct work_t { - env_t& environ; - work_items_t& input_items; - const self_t& parent; + env_t& environ; + work_items_t& input_items; + const self_t& parent; }; using worker_t = DigiParallelWorker<processor_t, work_t, std::size_t, self_t&>; @@ -290,9 +290,13 @@ namespace dd4hep { /// Default destructor virtual ~DigiContainerSequenceAction(); + /// Initialization callback virtual void initialize(); + /// Finalization callback + virtual void finalize(); + /// Get hold of the registered processor for a given container worker_t* need_registered_worker(Key item_key, bool exc=true) const; @@ -300,11 +304,11 @@ namespace dd4hep { /// Standard constructor DigiContainerSequenceAction(const kernel_t& kernel, const std::string& name); /// Set the default predicate - void set_predicate(const predicate_t& predicate); + virtual void set_predicate(const predicate_t& predicate); /// Adopt new parallel worker acting on one single container - void adopt_processor(DigiContainerProcessor* action, const std::string& container); + virtual void adopt_processor(DigiContainerProcessor* action, const std::string& container); /// Adopt new parallel worker acting on multiple containers - void adopt_processor(DigiContainerProcessor* action, const std::vector<std::string>& containers); + virtual void adopt_processor(DigiContainerProcessor* action, const std::vector<std::string>& containers); /// Main functional callback if specific work is known virtual void execute(context_t& context) const override; }; @@ -377,9 +381,9 @@ namespace dd4hep { return this->m_input_masks; } /// Set the default predicate - void set_predicate(const predicate_t& predicate); + virtual void set_predicate(const predicate_t& predicate); /// Adopt new parallel worker - void adopt_processor(DigiContainerProcessor* action, const std::vector<std::string>& containers); + virtual void adopt_processor(DigiContainerProcessor* action, const std::vector<std::string>& containers); /// Main functional callback virtual void execute(context_t& context) const; }; diff --git a/DDDigi/include/DDDigi/DigiData.h b/DDDigi/include/DDDigi/DigiData.h index 3ca4981feabddd4d18614826215a192fca592338..5212f621d56293796109816519829acf56cd6c07 100644 --- a/DDDigi/include/DDDigi/DigiData.h +++ b/DDDigi/include/DDDigi/DigiData.h @@ -87,6 +87,8 @@ namespace dd4hep { explicit Key(const char* item, mask_type mask); /// Initializing constructor with key generation using hash algorithm explicit Key(const std::string& item, mask_type mask); + /// Initializing constructor with key generation using hash algorithm + explicit Key(const std::string& item, segment_type segment, mask_type mask); /// Assignment operator Key& operator = (const Key::key_type) = delete; /// Assignment operator @@ -102,6 +104,8 @@ namespace dd4hep { /// Generate key using hash algorithm Key& set(const std::string& name, int mask); + /// Generate key using hash algorithm + Key& set(const std::string& name, int segment, int mask); /// Set key mask Key& set_mask(mask_type m); @@ -171,6 +175,11 @@ namespace dd4hep { this->set(item, mask); } + /// Initializing constructor with key generation using hash algorithm + inline Key::Key(const std::string& item, segment_type segmnt, mask_type msk) { + this->set(item, segmnt, msk); + } + /// Move assignment operator inline Key& Key::operator = (Key&& copy) { this->key = copy.key; @@ -269,11 +278,20 @@ namespace dd4hep { */ class SegmentEntry { public: + enum data_type_t { + UNKNOWN = 0, + PARTICLES = 1 << 1, + TRACKER_HITS = 1 << 2, + CALORIMETER_HITS = 1 << 3, + HISTORY = 1 << 4, + DETECTOR_RESPONSE = 1 << 5, + }; std::string name { }; Key key { }; + data_type_t data_type { UNKNOWN }; public: /// Initializing constructor - SegmentEntry(const std::string& name, Key::mask_type mask); + SegmentEntry(const std::string& name, Key::mask_type mask, data_type_t typ); /// Default constructor SegmentEntry() = default; /// Disable move constructor @@ -289,8 +307,8 @@ namespace dd4hep { }; /// Initializing constructor - inline SegmentEntry::SegmentEntry(const std::string& nam, Key::mask_type msk) - : name(nam) + inline SegmentEntry::SegmentEntry(const std::string& nam, Key::mask_type msk, data_type_t typ) + : name(nam), data_type(typ) { key.set(nam, msk); } @@ -308,14 +326,13 @@ namespace dd4hep { Position end_position { }; Direction momentum { }; double mass { 0e0 }; + int time { 0 }; int pdgID { 0 }; char charge { 0 }; /// Source contributing - Key history; + std::any source { }; //! not persistent public: - /// Initializing constructor - Particle(Key history); /// Default constructor Particle() = default; /// Disable move constructor @@ -328,14 +345,10 @@ namespace dd4hep { Particle& operator=(Particle&& copy) = default; /// Disable copy assignment Particle& operator=(const Particle& copy) = default; + /// Move particle + void move_position(const Position& delta); }; - /// Initializing constructor - inline Particle::Particle(Key h) - : history(h) - { - } - /// Particle mapping definition for digitization /** * @@ -403,7 +416,7 @@ namespace dd4hep { /// Initializing constructor inline ParticleMapping::ParticleMapping(const std::string& nam, Key::mask_type msk) - : SegmentEntry(nam, msk) + : SegmentEntry(nam, msk, SegmentEntry::PARTICLES) { } @@ -478,13 +491,13 @@ namespace dd4hep { class EnergyDeposit { public: enum { - KILLED = 1 << 0, - ENERGY_SMEARED = 1 << 1, - POSITION_SMEARED = 1 << 2, - TIME_SMEARED = 1 << 3, - ZERO_SUPPRESSED = 1 << 4, - DEPOSIT_NOISE = 1 << 5, - RECALIBRATED = 1 << 6 + KILLED = 1 << 0, + ENERGY_SMEARED = 1 << 1, + POSITION_SMEARED = 1 << 2, + TIME_SMEARED = 1 << 3, + ZERO_SUPPRESSED = 1 << 4, + DEPOSIT_NOISE = 1 << 5, + RECALIBRATED = 1 << 6, }; /// Hit position @@ -498,7 +511,7 @@ namespace dd4hep { /// Proper creation time of the deposit with rescpect to beam crossing double time { 0 }; /// Optional flag for user masks - long flag { 0 }; + uint64_t flag { 0UL }; /// Source mask of this deposit Key::mask_type mask { 0 }; /// Deposit history @@ -524,6 +537,37 @@ namespace dd4hep { }; + /// Base class to select energy deposits + /** + * + * \author M.Frank + * \version 1.0 + * \ingroup DD4HEP_DIGITIZATION + */ + template <typename USERDATA=int> class DepositPredicate { + public: + using userdata_t = USERDATA; + /// User data block + userdata_t data; + /// Standard constructor + DepositPredicate(const userdata_t& data=userdata_t()); + /// Default destructor + virtual ~DepositPredicate() = default; + /// Selector function + template <typename ARG> bool operator()(ARG argument) const; + }; + + /// Standard constructor + template <typename USERDATA> inline + DepositPredicate<USERDATA>::DepositPredicate(const userdata_t& d) + : data(d) + { + } + + struct EnergyCut { + double cutoff; + }; + /// Energy deposit vector definition for digitization /** * @@ -541,11 +585,11 @@ namespace dd4hep { using const_iterator = container_t::const_iterator; protected: - container_t data { }; + container_t data { }; public: /// Initializing constructor - DepositVector(const std::string& name, Key::mask_type mask); + DepositVector(const std::string& name, Key::mask_type mask, data_type_t typ); /// Default constructor DepositVector() = default; /// Disable move constructor @@ -596,8 +640,8 @@ namespace dd4hep { }; /// Initializing constructor - inline DepositVector::DepositVector(const std::string& nam, Key::mask_type msk) - : SegmentEntry(nam, msk) + inline DepositVector::DepositVector(const std::string& nam, Key::mask_type msk, data_type_t typ) + : SegmentEntry(nam, msk, typ) { } @@ -622,11 +666,11 @@ namespace dd4hep { using iterator = container_t::iterator; using const_iterator = container_t::const_iterator; - container_t data { }; + container_t data { }; public: /// Initializing constructor - DepositMapping(const std::string& name, Key::mask_type mask); + DepositMapping(const std::string& name, Key::mask_type mask, data_type_t typ); /// Default constructor DepositMapping() = default; /// Disable move constructor @@ -673,8 +717,8 @@ namespace dd4hep { }; /// Initializing constructor - inline DepositMapping::DepositMapping(const std::string& nam, Key::mask_type msk) - : SegmentEntry(nam, msk) + inline DepositMapping::DepositMapping(const std::string& nam, Key::mask_type msk, data_type_t typ) + : SegmentEntry(nam, msk, typ) { } @@ -743,7 +787,7 @@ namespace dd4hep { /// Initializing constructor inline DetectorResponse::DetectorResponse(const std::string& nam, Key::mask_type msk) - : SegmentEntry(nam, msk) + : SegmentEntry(nam, msk, SegmentEntry::DETECTOR_RESPONSE) { } @@ -812,7 +856,7 @@ namespace dd4hep { /// Initializing constructor inline DetectorHistory::DetectorHistory(const std::string& nam, Key::mask_type msk) - : SegmentEntry(nam, msk) + : SegmentEntry(nam, msk, SegmentEntry::HISTORY) { } @@ -926,31 +970,45 @@ namespace dd4hep { }; /// Access data as reference by key. If not existing, an exception is thrown - template<typename T> inline T& DataSegment::get(Key key) { - if ( T* ptr = std::any_cast<T>(this->get_item(key, true)) ) - return *ptr; - throw std::runtime_error(this->invalid_cast(key, typeid(T))); + template<typename DATA> inline DATA& DataSegment::get(Key key) { + if ( DATA* ptr = std::any_cast<DATA>(this->get_item(key, true)) ) + return *ptr; + throw std::runtime_error(this->invalid_cast(key, typeid(DATA))); } /// Access data as reference by key. If not existing, an exception is thrown - template<typename T> inline const T& DataSegment::get(Key key) const { - if ( const T* ptr = std::any_cast<T>(this->get_item(key, true)) ) - return *ptr; - throw std::runtime_error(this->invalid_cast(key, typeid(T))); + template<typename DATA> inline const DATA& DataSegment::get(Key key) const { + if ( const DATA* ptr = std::any_cast<DATA>(this->get_item(key, true)) ) + return *ptr; + throw std::runtime_error(this->invalid_cast(key, typeid(DATA))); } /// Access data as pointers by key. If not existing, nullptr is returned - template<typename T> inline T* DataSegment::pointer(Key key) { - if ( T* ptr = std::any_cast<T>(this->get_item(key, false)) ) - return ptr; + template<typename DATA> inline DATA* DataSegment::pointer(Key key) { + if ( DATA* ptr = std::any_cast<DATA>(this->get_item(key, false)) ) + return ptr; return nullptr; } /// Access data as pointers by key. If not existing, nullptr is returned - template<typename T> inline const T* DataSegment::pointer(Key key) const { - if ( const T* ptr = std::any_cast<T>(this->get_item(key, false)) ) - return ptr; + template<typename DATA> inline const DATA* DataSegment::pointer(Key key) const { + if ( const DATA* ptr = std::any_cast<DATA>(this->get_item(key, false)) ) + return ptr; return nullptr; } + /// Move data items other than std::any to the data segment + template <typename DATA> inline bool DataSegment::put(Key key, DATA&& value) { + key.set_segment(this->id); + value.key.set_segment(this->id); + std::any item = std::make_any<DATA>(std::move(value)); + return this->emplace_any(key, std::move(item)); + } + + /// Helper to place data to data segment + template <typename KEY, typename DATA> + bool put_data(DataSegment& segment, KEY key, DATA& data) { + return segment.emplace_any(key, std::any(data)); + } + /// User event data for DDDigi /** * diff --git a/DDDigi/include/DDDigi/DigiOutputAction.h b/DDDigi/include/DDDigi/DigiOutputAction.h index ae67ca80b934d5cb334bfe952710eba8c5c6ff1a..fb2013ceb0cbba8922d1cfda4564fc9899074826 100644 --- a/DDDigi/include/DDDigi/DigiOutputAction.h +++ b/DDDigi/include/DDDigi/DigiOutputAction.h @@ -52,6 +52,7 @@ namespace dd4hep { public: /// Standard constructor DigiOutputAction(const kernel_t& kernel, const std::string& nam); + /// Default destructor virtual ~DigiOutputAction(); diff --git a/DDDigi/ddg4/IO.cpp b/DDDigi/io/DDG4IO.cpp similarity index 100% rename from DDDigi/ddg4/IO.cpp rename to DDDigi/io/DDG4IO.cpp diff --git a/DDDigi/io/Digi2edm4hep.cpp b/DDDigi/io/Digi2edm4hep.cpp new file mode 100644 index 0000000000000000000000000000000000000000..a32a64d78f139f3883611367b94698157f6e6768 --- /dev/null +++ b/DDDigi/io/Digi2edm4hep.cpp @@ -0,0 +1,480 @@ +//========================================================================== +// AIDA Detector description implementation +//-------------------------------------------------------------------------- +// Copyright (C) Organisation europeenne pour la Recherche nucleaire (CERN) +// All rights reserved. +// +// For the licensing terms see $DD4hepINSTALL/LICENSE. +// For the list of contributors see $DD4hepINSTALL/doc/CREDITS. +// +// Author : M.Frank +// +//========================================================================== +#ifndef DIGI_DIGI2EDM4HEP_H +#define DIGI_DIGI2EDM4HEP_H + +// Framework include files +#include <DDDigi/DigiContainerProcessor.h> + +/// C/C++ include files + +namespace edm4hep { + class TrackerHitCollection; + class CalorimeterHitCollection; +} + +/// Namespace for the AIDA detector description toolkit +namespace dd4hep { + + /// Namespace for the Digitization part of the AIDA detector description toolkit + namespace digi { + + /// Event action to support edm4hep output format from DDDigi + /** + * Supported output containers types are: + * - edm4hep::MCParticles aka "MCParticles" + * - edm4hep::CalorimeterHitCollection aka "CalorimeterHits" + * - edm4hep::TrackerHitCollection aka "TracketHits" + * + * \author M.Frank + * \version 1.0 + * \ingroup DD4HEP_DIGITIZATION + */ + class Digi2edm4hepAction : public DigiContainerSequenceAction { + public: + class internals_t; + + protected: + /// Property: Container names to be loaded + std::string m_output; + /// Property: Processor type to manage containers + std::string m_processor_type; + /// Property: Container / data type mapping + std::map<std::string, std::string> m_containers { }; + /// Reference to internals + std::shared_ptr<internals_t> internals; + + protected: + /// Define standard assignments and constructors + DDDIGI_DEFINE_ACTION_CONSTRUCTORS(Digi2edm4hepAction); + /// Default destructor + virtual ~Digi2edm4hepAction(); + + public: + /// Standard constructor + Digi2edm4hepAction(const kernel_t& kernel, const std::string& nam); + + /// Initialization callback + virtual void initialize(); + + /// Finalization callback + virtual void finalize(); + + /// Adopt new parallel worker + virtual void adopt_processor(DigiContainerProcessor* action, + const std::string& container) override final; + + /// Adopt new parallel worker acting on multiple containers + virtual void adopt_processor(DigiContainerProcessor* action, + const std::vector<std::string>& containers); + + /// Callback to store the run information + void beginRun(); + + /// Callback to store the run information + void endRun(); + + /// Callback to store the Geant4 run information + void saveRun(); + + /// Main functional callback + virtual void execute(context_t& context) const; + }; + + /// Actor to select energy deposits according to the supplied segmentation + /** Actor to select energy deposits according to the supplied segmentation + * + * The selected deposits are placed in the output container + * supplied by the arguments. + * + * \author M.Frank + * \version 1.0 + * \ingroup DD4HEP_DIGITIZATION + */ + class Digi2edm4hepProcessor : public DigiContainerProcessor { + friend class Digi2edm4hepAction; + + protected: + std::shared_ptr<Digi2edm4hepAction::internals_t> internals; + float pointResoRPhi = 0.004; + float pointResoZ = 0.004; + + public: + /// Standard constructor + Digi2edm4hepProcessor(const DigiKernel& krnl, const std::string& nam); + + /// Standard destructor + virtual ~Digi2edm4hepProcessor() = default; + + template <typename T> void + convert_depos(const T& cont, const predicate_t& predicate, edm4hep::TrackerHitCollection* collection) const; + + template <typename T> void + convert_depos(const T& cont, const predicate_t& predicate, edm4hep::CalorimeterHitCollection* collection) const; + + template <typename T> void + convert_deposits(DigiContext& context, const T& cont, const predicate_t& predicate) const; + + void convert_history(DigiContext& context, const DepositsHistory& cont, work_t& work, const predicate_t& predicate) const; + + void convert_particles(DigiContext& context, const ParticleMapping& cont) const; + + /// Main functional callback + virtual void execute(DigiContext& context, work_t& work, const predicate_t& predicate) const override final; + }; + } // End namespace digi +} // End namespace dd4hep +#endif // DIGI_DIGI2EDM4HEP_H + +//========================================================================== +// AIDA Detector description implementation +//-------------------------------------------------------------------------- +// Copyright (C) Organisation europeenne pour la Recherche nucleaire (CERN) +// All rights reserved. +// +// For the licensing terms see $DD4hepINSTALL/LICENSE. +// For the list of contributors see $DD4hepINSTALL/doc/CREDITS. +// +// Author : M.Frank +// +//========================================================================== + +/// Framework include files +#include <DD4hep/InstanceCount.h> +#include <DDDigi/DigiContext.h> +#include <DDDigi/DigiPlugins.h> +#include <DDDigi/DigiKernel.h> +#include "DigiIO.h" + +/// edm4hep include files +#include <podio/EventStore.h> +#include <podio/ROOTWriter.h> +#include <edm4hep/SimTrackerHit.h> +#include <edm4hep/MCParticleCollection.h> +#include <edm4hep/TrackerHitCollection.h> +#include <edm4hep/EventHeaderCollection.h> +#include <edm4hep/CalorimeterHitCollection.h> +#include <edm4hep/CaloHitContributionCollection.h> + + +/// Namespace for the AIDA detector description toolkit +namespace dd4hep { + + /// Namespace for the Digitization part of the AIDA detector description toolkit + namespace digi { + + class Digi2edm4hepAction::internals_t { + public: + Digi2edm4hepAction* m_parent { nullptr }; + std::unique_ptr<podio::EventStore> m_store { }; + std::unique_ptr<podio::ROOTWriter> m_file { }; + edm4hep::EventHeaderCollection* m_header { nullptr }; + /// MC particle collection + edm4hep::MCParticleCollection* m_particles { nullptr }; + /// Collection of all edm4hep object collections + std::map<std::string, podio::CollectionBase*> m_collections; + + /// Total numbe rof events to be processed + long num_events { -1 }; + /// Running event counter + long event_count { 0 }; + + private: + /// Helper to register single collection + template <typename T> T* register_collection(const std::string& name, T* collection); + + public: + /// Default constructor + internals_t() = default; + /// Default destructor + ~internals_t() = default; + + /// Commit data at end of filling procedure + void commit(); + /// Commit data to disk and close output stream + void close(); + + /// Create all collections according to the parent setup (locked) + void create_collections(); + /// Access named collection: throws exception ifd the collection is not present (unlocked!) + template <typename T> podio::CollectionBase* get_collection(const T&); + }; + + template <typename T> T* Digi2edm4hepAction::internals_t::register_collection(const std::string& nam, T* coll) { + m_collections.emplace(nam, coll); + m_store->registerCollection(nam, coll); + m_file->registerForWrite(nam); + m_parent->debug("+++ created collection %s <%s>", nam.c_str(), coll->getTypeName().c_str()); + return coll; + } + + /// Create all collections according to the parent setup + void Digi2edm4hepAction::internals_t::create_collections() { + if ( nullptr == m_header ) { + std::string fname = m_parent->m_output; + m_store = std::make_unique<podio::EventStore>(); + m_file = std::make_unique<podio::ROOTWriter>(fname, m_store.get()); + m_header = register_collection("EventHeader", new edm4hep::EventHeaderCollection()); + for( auto& cont : m_parent->m_containers ) { + const std::string& nam = cont.first; + const std::string& typ = cont.second; + if ( typ == "MCParticles" ) { + m_particles = register_collection(nam, new edm4hep::MCParticleCollection()); + } + else if ( typ == "TrackerHits" ) { + register_collection(nam, new edm4hep::TrackerHitCollection()); + } + else if ( typ == "CalorimeterHits" ) { + register_collection(nam, new edm4hep::CalorimeterHitCollection()); + } + } + m_parent->info("+++ Will save %ld events to %s", num_events, m_parent->m_output.c_str()); + } + } + + /// Access named collection: throws exception ifd the collection is not present + template <typename T> + podio::CollectionBase* Digi2edm4hepAction::internals_t::get_collection(const T& cont) { + auto iter = m_collections.find(cont.name); + if ( iter == m_collections.end() ) { + m_parent->except("Error"); + } + return iter->second; + } + + /// Commit data at end of filling procedure + void Digi2edm4hepAction::internals_t::commit() { + if ( m_file ) { + m_file->writeEvent(); + m_store->clearCollections(); + return; + } + m_parent->except("+++ Failed to write output file. [Stream is not open]"); + } + + /// Commit data to disk and close output stream + void Digi2edm4hepAction::internals_t::close() { + if ( m_file ) { + m_file->finish(); + } + m_file.reset(); + m_store.reset(); + } + + /// Standard constructor + Digi2edm4hepAction::Digi2edm4hepAction(const DigiKernel& krnl, const std::string& nam) + : DigiContainerSequenceAction(krnl, nam) + { + internals = std::make_shared<internals_t>(); + declareProperty("output", m_output); + declareProperty("containers", m_containers); + declareProperty("processor_type", m_processor_type = "Digi2edm4hepProcessor"); + internals->m_parent = this; + InstanceCount::increment(this); + } + + /// Default destructor + Digi2edm4hepAction::~Digi2edm4hepAction() {{ + std::lock_guard<std::mutex> lock(m_kernel.global_io_lock()); + internals->close(); + } + internals.reset(); + InstanceCount::decrement(this); + } + + /// Initialization callback + void Digi2edm4hepAction::initialize() { + if ( m_containers.empty() ) { + warning("+++ No input containers given for attenuation action -- no action taken"); + return; + } + internals->num_events = m_kernel.property("numEvents").value<long>(); + for ( const auto& c : m_containers ) { + Key key(c.first, 0, 0); + auto it = m_registered_processors.find(key); + if ( it == m_registered_processors.end() ) { + std::string nam = name() + ".E4H." + c.first; + auto* conv = createAction<DigiContainerProcessor>(m_processor_type, m_kernel, nam); + if ( !conv ) { + except("+++ Failed to create edm4hep processor: %s of type: %s", + nam.c_str(), m_processor_type.c_str()); + } + conv->property("OutputLevel").set(int(outputLevel())); + adopt_processor(conv, c.first); + conv->release(); // Release processor **after** adoption. + } + } + std::lock_guard<std::mutex> lock(m_kernel.global_io_lock()); + m_parallel = false; + internals->create_collections(); + this->DigiContainerSequenceAction::initialize(); + } + + /// Finalization callback + void Digi2edm4hepAction::finalize() { + internals->close(); + } + + /// Adopt new parallel worker + void Digi2edm4hepAction::adopt_processor(DigiContainerProcessor* action, + const std::string& container) + { + std::size_t idx = container.find('/'); + if ( idx != std::string::npos ) { + std::string nam = container.substr(0, idx); + std::string typ = container.substr(idx+1); + auto* act = dynamic_cast<Digi2edm4hepProcessor*>(action); + if ( act ) { // This is not nice! Need to think about something better. + act->internals = this->internals; + } + this->DigiContainerSequenceAction::adopt_processor(action, nam); + m_containers.emplace(nam, typ); + return; + } + except("+++ Invalid container specification: %s. %s", + container.c_str(), "Specify container as tuple: \"<name>/<type>\" !"); + } + + /// Adopt new parallel worker acting on multiple containers + void Digi2edm4hepAction::adopt_processor(DigiContainerProcessor* action, + const std::vector<std::string>& containers) + { + DigiContainerSequenceAction::adopt_processor(action, containers); + } + + /// Main functional callback + void Digi2edm4hepAction::execute(DigiContext& context) const { + std::lock_guard<std::mutex> lock(context.global_io_lock()); + this->DigiContainerSequenceAction::execute(context); + this->internals->commit(); + if ( ++internals->event_count == internals->num_events ) { + internals->close(); + } + } + + /// Callback to store the run information + void Digi2edm4hepAction::beginRun() { + saveRun(); + } + + /// Callback to store the run information + void Digi2edm4hepAction::endRun() { + // saveRun(run); + } + + /// Callback to store the Geant4 run information + void Digi2edm4hepAction::saveRun() { + warning("saveRun(): RunHeader not implemented in EDM4hep, nothing written ..."); + } + + + /// Standard constructor + Digi2edm4hepProcessor::Digi2edm4hepProcessor(const DigiKernel& krnl, const std::string& nam) + : DigiContainerProcessor(krnl, nam) + { + } + + void Digi2edm4hepProcessor::convert_particles(DigiContext& ctxt, + const ParticleMapping& cont) const + { + auto* coll = internals->m_particles; + std::size_t start = coll->size(); + digi_io()._to_edm4hep(cont, coll); + std::size_t end = internals->m_particles->size(); + info("%s+++ %-24s added %6ld/%6ld entries from mask: %04X to %s", + ctxt.event->id(), cont.name.c_str(), end-start, end, cont.key.mask(), + coll->getTypeName().c_str()); + } + + template <typename T> void + Digi2edm4hepProcessor::convert_depos(const T& cont, + const predicate_t& predicate, + edm4hep::TrackerHitCollection* collection) const + { + std::array<float,6> covMat = {0., 0., pointResoRPhi*pointResoRPhi, 0., 0., pointResoZ*pointResoZ}; + for ( const auto& depo : cont ) { + if ( predicate(depo) ) { + digi_io::_to_edm4hep(depo, covMat, collection); + } + } + } + + template <typename T> void + Digi2edm4hepProcessor::convert_depos(const T& cont, + const predicate_t& predicate, + edm4hep::CalorimeterHitCollection* collection) const + { + for ( const auto& depo : cont ) { + if ( predicate(depo) ) { + digi_io::_to_edm4hep(depo, collection); + } + } + } + + template <typename T> void + Digi2edm4hepProcessor::convert_deposits(DigiContext& ctxt, + const T& cont, + const predicate_t& predicate) const + { + podio::CollectionBase* coll = internals->get_collection(cont); + std::size_t start = coll->size(); + if ( !cont.empty() ) { + switch(cont.data_type) { + case SegmentEntry::TRACKER_HITS: + convert_depos(cont, predicate, static_cast<edm4hep::TrackerHitCollection*>(coll)); + break; + case SegmentEntry::CALORIMETER_HITS: + convert_depos(cont, predicate, static_cast<edm4hep::CalorimeterHitCollection*>(coll)); + break; + default: + except("Error: Unknown energy deposit type: %d", int(cont.data_type)); + break; + } + } + std::size_t end = coll->size(); + info("%s+++ %-24s added %6ld/%6ld entries from mask: %04X to %s", + ctxt.event->id(), cont.name.c_str(), end-start, end, cont.key.mask(), + coll->getTypeName().c_str()); + } + + void Digi2edm4hepProcessor::convert_history(DigiContext& ctxt, + const DepositsHistory& cont, + work_t& work, + const predicate_t& predicate) const + { + info("%s+++ %-32s Segment: %d Predicate:%s Conversion to edm4hep not implemented!", + ctxt.event->id(), cont.name.c_str(), int(work.input.segment->id), + typeName(typeid(predicate)).c_str()); + } + + /// Main functional callback + void Digi2edm4hepProcessor::execute(DigiContext& ctxt, work_t& work, const predicate_t& predicate) const { + if ( const auto* p = work.get_input<ParticleMapping>() ) + convert_particles(ctxt, *p); + else if ( const auto* m = work.get_input<DepositMapping>() ) + convert_deposits(ctxt, *m, predicate); + else if ( const auto* v = work.get_input<DepositVector>() ) + convert_deposits(ctxt, *v, predicate); + else if ( const auto* h = work.get_input<DepositsHistory>() ) + convert_history(ctxt, *h, work, predicate); + else + except("Request to handle unknown data type: %s", work.input_type_name().c_str()); + } + + } // End namespace digi +} // End namespace dd4hep + +/// Factory instantiation: +#include <DDDigi/DigiFactories.h> +DECLARE_DIGIACTION_NS(dd4hep::digi,Digi2edm4hepAction) +DECLARE_DIGIACTION_NS(dd4hep::digi,Digi2edm4hepProcessor) diff --git a/DDDigi/io/Digi2edm4hepEx.cpp b/DDDigi/io/Digi2edm4hepEx.cpp new file mode 100644 index 0000000000000000000000000000000000000000..5e55aace0e87eb4f90e38e51c3517c2dd9614754 --- /dev/null +++ b/DDDigi/io/Digi2edm4hepEx.cpp @@ -0,0 +1,98 @@ +//========================================================================== +// AIDA Detector description implementation +//-------------------------------------------------------------------------- +// Copyright (C) Organisation europeenne pour la Recherche nucleaire (CERN) +// All rights reserved. +// +// For the licensing terms see $DD4hepINSTALL/LICENSE. +// For the list of contributors see $DD4hepINSTALL/doc/CREDITS. +// +// Author : M.Frank +// +//========================================================================== + +// Framework include files +#include <DDDigi/DigiContext.h> +#include <DDDigi/DigiContainerProcessor.h> + +// edm4hep include files +#include <podio/EventStore.h> +#include <podio/ROOTWriter.h> +#include <edm4hep/MCParticleCollection.h> +#include <edm4hep/EventHeaderCollection.h> +#include <edm4hep/SimTrackerHitCollection.h> +#include <edm4hep/SimCalorimeterHitCollection.h> +#include <edm4hep/CaloHitContributionCollection.h> + +/// Namespace for the AIDA detector description toolkit +namespace dd4hep { + + /// Namespace for the Digitization part of the AIDA detector description toolkit + namespace digi { + + /// Actor to select energy deposits according to the supplied segmentation + /** Actor to select energy deposits according to the supplied segmentation + * + * The selected deposits are placed in the output container + * supplied by the arguments. + * + * \author M.Frank + * \version 1.0 + * \ingroup DD4HEP_DIGITIZATION + */ + class Digi2edm4hep : public DigiContainerProcessor { + + std::unordered_map<std::string, podio::CollectionBase*> m_collections; + public: + /// Standard constructor + Digi2edm4hep(const DigiKernel& krnl, const std::string& nam) + : DigiContainerProcessor(krnl, nam) + { + //declareProperty("mean", m_mean); + //declareProperty("sigma", m_sigma); + } + + /// Standard destructor + virtual ~Digi2edm4hep() = default; + + template <typename T> void + convert_deposits(const char* tag, const T& cont, work_t& work, const predicate_t& predicate) const { + Key key(cont.name, work.environ.output.mask); + DepositMapping m(cont.name, work.environ.output.mask, cont.data_type); + std::size_t start = m.size(); + for( const auto& dep : cont ) { + if ( predicate(dep) ) { + m.data.emplace(dep.first, EnergyDeposit()); + } + } + std::size_t end = m.size(); + work.environ.output.data.put(m.key, std::move(m)); + info("%s+++ %-32s added %6ld entries (now: %6ld) from mask: %04X to mask: %04X", + tag, cont.name.c_str(), end-start, end, cont.key.mask(), m.key.mask()); + } + void convert_history(const char* tag, const DepositsHistory& cont, work_t& work, const predicate_t& predicate) const { + } + void convert_particles(const char* tag, const ParticleMapping& cont, work_t& work, const predicate_t& predicate) const { + //auto mcpc = static_cast<edm4hep::MCParticleCollection*>(m_collections["MCParticles"]); + + } + /// Main functional callback + virtual void execute(DigiContext& context, work_t& work, const predicate_t& predicate) const override final { + if ( const auto* m = work.get_input<DepositMapping>() ) + convert_deposits(context.event->id(), *m, work, predicate); + else if ( const auto* v = work.get_input<DepositVector>() ) + convert_deposits(context.event->id(), *v, work, predicate); + else if ( const auto* h = work.get_input<DepositsHistory>() ) + convert_history(context.event->id(), *h, work, predicate); + else if ( const auto* p = work.get_input<ParticleMapping>() ) + convert_particles(context.event->id(), *p, work, predicate); + else + except("Request to handle unknown data type: %s", work.input_type_name().c_str()); + } + }; + } // End namespace digi +} // End namespace dd4hep + +/// Factory instantiation: +#include <DDDigi/DigiFactories.h> +DECLARE_DIGIACTION_NS(dd4hep::digi,Digi2edm4hep) diff --git a/DDDigi/io/DigiDDG4Input.cpp b/DDDigi/io/DigiDDG4Input.cpp new file mode 100644 index 0000000000000000000000000000000000000000..1dda13b5b6bf330571c4b4629e84c17e4aff058b --- /dev/null +++ b/DDDigi/io/DigiDDG4Input.cpp @@ -0,0 +1,126 @@ +//========================================================================== +// AIDA Detector description implementation +//-------------------------------------------------------------------------- +// Copyright (C) Organisation europeenne pour la Recherche nucleaire (CERN) +// All rights reserved. +// +// For the licensing terms see $DD4hepINSTALL/LICENSE. +// For the list of contributors see $DD4hepINSTALL/doc/CREDITS. +// +// Author : M.Frank +// +//========================================================================== + +// Framework include files +#include <DDDigi/DigiROOTInput.h> +#include <DDDigi/DigiFactories.h> +#include "DigiIO.h" + +#include <DDG4/Geant4Data.h> + +// ROOT include files +#include <TBranch.h> +#include <TClass.h> + +// C/C++ include files + +/// Namespace for the AIDA detector description toolkit +namespace dd4hep { + + /// Namespace for the Digitization part of the AIDA detector description toolkit + namespace digi { + + using namespace std::placeholders; + + class DigiDDG4ROOT : public DigiROOTInput { + public: + static constexpr double epsilon = std::numeric_limits<double>::epsilon(); + + /// Property to generate extra history records + bool m_keep_raw { true }; + mutable TClass* m_trackerHitClass { nullptr }; + mutable TClass* m_caloHitClass { nullptr }; + mutable TClass* m_particlesClass { nullptr }; + + TClass* get_class(TClass* c) const { + if ( c == m_particlesClass || c == m_trackerHitClass || c == m_caloHitClass ) + return c; + const char* nam = c->GetName(); + if ( strstr(nam, "vector<dd4hep::sim::Geant4Particle*") ) + return m_particlesClass = c; + else if ( strstr(nam, "vector<dd4hep::sim::Geant4Tracker::Hit*") ) + return m_trackerHitClass = c; + else if ( strstr(nam, "vector<dd4hep::sim::Geant4Calorimeter::Hit*") ) + return m_caloHitClass = c; + except("Unknown data class: %s Cannot be handled", nam); + return c; + } + + public: + /// Initializing constructor + DigiDDG4ROOT(const DigiKernel& krnl, const std::string& nam) + : DigiROOTInput(krnl, nam) + { + declareProperty("keep_raw", m_keep_raw); + } + + template <typename T> + void conv_hits(DigiContext& context, DataSegment& segment, + const std::string& tag, + Key::mask_type mask, + const char* nam, + void* ptr) const + { + DepositVector out(nam, mask, SegmentEntry::UNKNOWN); + std::map<CellID, std::shared_ptr<T> > hits; + std::size_t len = 0; + if ( ptr ) { + input_data<T> data(ptr); + const DepositPredicate<EnergyCut> predicate ({ this->epsilon }); + len = data.size(); + data_io<ddg4_input>::_to_digi_if(data.get(), hits, predicate); + data_io<ddg4_input>::_to_digi(Key(nam, segment.id, mask), hits, out); + data.clear(); + } + info("%s+++ %-24s Converted %6ld DDG4 %-14s hits to %6ld cell deposits", + context.event->id(), nam, len, tag.c_str(), out.size()); + put_data(segment, Key(out.name, mask), out); + if ( m_keep_raw ) { + put_data(segment, Key(std::string(nam)+".ddg4", mask, segment.id), hits); + } + } + + void conv_particles(DigiContext& context, DataSegment& segment, + int mask, const std::string& nam, void* ptr) const + { + ParticleMapping particles(nam, mask); + if ( ptr ) { + input_data<sim::Geant4Particle> data(ptr); + data_io<ddg4_input>::_to_digi(Key(nam, segment.id, mask), data.get(), particles); + data.clear(); + } + debug("%s+++ Converted %ld DDG4 particles", context.event->id(), particles.size()); + put_data(segment, Key(nam, mask), particles); + } + + /// Callback to handle single branch + virtual void operator()(DigiContext& context, work_t& work) const override { + TBranch& br = work.branch; + int msk = work.input_key.mask(); + auto& seg = work.input_segment; + void** add = (void**)br.GetAddress(); + TClass* cls = get_class(&work.cl); + const char* nam = br.GetName(); + if ( cls == m_caloHitClass ) + conv_hits<sim::Geant4Calorimeter::Hit>(context, seg, "calorimeter", msk, nam, *add); + else if ( cls == m_trackerHitClass ) + conv_hits<sim::Geant4Tracker::Hit>(context, seg, "tracker", msk, nam, *add); + else if ( cls == m_particlesClass ) + conv_particles(context, seg, msk, nam, *add); + else + except("Unknown data type encountered in branch: %s", nam); + } + }; + } // End namespace digi +} // End namespace dd4hep +DECLARE_DIGIACTION_NS(dd4hep::digi,DigiDDG4ROOT) diff --git a/DDDigi/io/DigiIO.cpp b/DDDigi/io/DigiIO.cpp new file mode 100644 index 0000000000000000000000000000000000000000..b54490689995ad439756cfe92e457f7b908ab5f5 --- /dev/null +++ b/DDDigi/io/DigiIO.cpp @@ -0,0 +1,430 @@ +//========================================================================== +// AIDA Detector description implementation +//-------------------------------------------------------------------------- +// Copyright (C) Organisation europeenne pour la Recherche nucleaire (CERN) +// All rights reserved. +// +// For the licensing terms see $DD4hepINSTALL/LICENSE. +// For the list of contributors see $DD4hepINSTALL/doc/CREDITS. +// +// Author : M.Frank +// +//========================================================================== + +/// Framework include files +#include <DD4hep/Primitives.h> +#include <DD4hep/Printout.h> +#include "DigiIO.h" + +// ========================================================================= +// EDM4HEP specific stuff +// ========================================================================= + +#ifdef DD4HEP_USE_EDM4HEP + +/// edm4hep include files +#include <edm4hep/SimTrackerHit.h> +#include <edm4hep/MCParticle.h> +#include <edm4hep/MCParticleCollection.h> + +/// Namespace for the AIDA detector description toolkit +namespace dd4hep { + + /// Namespace for the Digitization part of the AIDA detector description toolkit + namespace digi { + + struct bla { + class my_part; + typedef my_part particle_type; + }; + + + edm4hep::Vector3d _toVectorD(const dd4hep::Position& ep) { + return { ep.x(), ep.y(), ep.z() }; + } + + edm4hep::Vector3f _toVectorF(const dd4hep::Position& ep) { + return { float(ep.x()), float(ep.y()), float(ep.z()) }; + } + namespace { + template <typename DATA> bool internal_can_handle(const DATA&, const std::type_info&) { + return true; + } + template <> bool internal_can_handle(const ParticleMapping::value_type& data, const std::type_info& info) { + return (data.second.source.type() == info); + } + } + + template <typename T> template <typename DATA> + bool data_io<T>::_can_handle(const DATA& data) { + return internal_can_handle(data, typeid(pwrap_t)); + } + + template <typename T> template <typename CONT> + void data_io<T>::_pre_create(CONT* coll, std::size_t n) { + /// We have to pre-create all objects to be able to fill the parent-daughter relationships + for ( std::size_t i=0; i<n; ++i ) { + coll->create(); + } + } + + template <typename T> template <typename CONT> + std::vector<const typename data_io<T>::particle_t*> + data_io<T>::_to_vector(const CONT& cont) { + std::vector<const particle_t*> vec; + vec.reserve(cont.size()); + for ( const auto& part : cont ) { + const auto& p = part.second; + if ( p.source.type() == typeid(pwrap_t) ) { + const auto* ptr = std::any_cast<pwrap_t>(&p.source); + vec.emplace_back(ptr->get()); + } + } + if ( cont.size() != vec.size() ) { + except("data_io","_to_vector: Containers of mixed origin are not supported!"); + } + return vec; + } + + template <typename T> template <typename FIRST, typename SECOND> + void data_io<T>::_to_edm4hep(const FIRST&, SECOND) { + except("data_io::_to_edm4hep","(%s&, %s): Implementation not present!", + typeName(typeid(FIRST)).c_str(), typeName(typeid(SECOND)).c_str()); + } + + /// Set all properties of the MutableMCParticle + template <> template <> + void data_io<digi_input>::_to_edm4hep(const particle_t& p, + edm4hep::MutableMCParticle mcp) { + mcp.setPDG(p.pdgID); + mcp.setTime(p.time); + mcp.setMass(p.mass); + mcp.setCharge(3.0*p.charge); + mcp.setVertex( _toVectorD(p.start_position) ); + mcp.setEndpoint( _toVectorD(p.end_position) ); + mcp.setMomentum( _toVectorF(p.momentum) ); + mcp.setMomentumAtEndpoint( _toVectorF(p.momentum) ); + } + + template <> template <> + void data_io<digi_input>::_to_edm4hep(const std::vector<const particle_t*>& cont, + edm4hep::MCParticleCollection* coll) { + std::size_t i, n = cont.size(); + _pre_create(coll, n); + /// Convert particle body + for ( i=0; i<n; ++i) { + _to_edm4hep(*cont[i], coll->at(i)); + } + } + + /// Set all properties of the MutableMCParticle + template <> template <> + void data_io<edm4hep_input>::_to_edm4hep(const particle_t& p, + edm4hep::MutableMCParticle mcp) + { + mcp.setPDG( p.getPDG() ); + mcp.setMomentum( p.getMomentum() ); + mcp.setMomentumAtEndpoint( p.getMomentumAtEndpoint() ); + mcp.setVertex( p.getVertex() ); + mcp.setEndpoint( p.getEndpoint() ); + mcp.setTime( p.getTime() ); + mcp.setMass( p.getMass() ); + mcp.setCharge( p.getCharge() ); + mcp.setGeneratorStatus( p.getGeneratorStatus() ); + mcp.setSimulatorStatus( p.getSimulatorStatus() ); + mcp.setSpin(p.getSpin()); + mcp.setColorFlow(p.getColorFlow()); + } + + template <> template <> + void data_io<edm4hep_input>::_to_edm4hep(const std::vector<const particle_t*>& cont, + edm4hep::MCParticleCollection* coll) + { + std::size_t i, n = cont.size(); + _pre_create(coll, n); + /// Convert particle body + for ( i=0; i<n; ++i) { + const particle_t* p = cont[i]; + auto mcp = coll->at(i); + _to_edm4hep(*p, mcp); +#if 0 + /// Relationships are already resolved and kept in order: Just copy indices + for (std::size_t idau = 0; idau < p->daughters_size(); ++idau) { + mcp.addToDaughters(coll->at(idau)); + } + for (auto ipar : p->parents) { + mcp.addToParents(coll->at(ipar)); + } +#endif + } + } + } // End namespace digi +} // End namespace dd4hep + +// ========================================================================= +// DDG4 specific stuff +// ========================================================================= +#if defined(DD4HEP_USE_DDG4) + +#include <DDG4/Geant4Data.h> +#include <DDG4/Geant4Particle.h> + +/// Namespace for the AIDA detector description toolkit +namespace dd4hep { + + /// Namespace for the Digitization part of the AIDA detector description toolkit + namespace digi { + + using PropertyMask = dd4hep::detail::ReferenceBitMask<int>; + + template <> template <> + bool DepositPredicate<EnergyCut>::operator()(sim::Geant4Tracker::Hit* h) const { + return h->energyDeposit > data.cutoff; + } + + template <> template <> + bool DepositPredicate<EnergyCut>::operator()(sim::Geant4Calorimeter::Hit* h) const { + return h->energyDeposit > data.cutoff; + } + + void add_particle_history(const sim::Geant4Calorimeter::Hit* hit, Key key, History& hist) { + for( const auto& truth : hit->truth ) { + key.set_item(truth.trackID); + hist.particles.emplace_back(key, truth.deposit); + } + } + + void add_particle_history(const sim::Geant4Tracker::Hit* hit, Key key, History& hist) { + key.set_item(hit->truth.trackID); + hist.particles.emplace_back(key, hit->truth.deposit); + } + + /// Set all properties of the MutableMCParticle + template <> template <> + void data_io<ddg4_input>::_to_edm4hep(const particle_t& p, + edm4hep::MutableMCParticle mcp) + { + auto status = p.status; + const PropertyMask mask(status); + mcp.setPDG(p.pdgID); + + mcp.setMomentum( _toVectorF( { p.psx, p.psy, p.psz } ) ); + mcp.setMomentumAtEndpoint( _toVectorF( {p.pex, p.pey, p.pez} ) ); + mcp.setVertex( _toVectorD( { p.vsx, p.vsy, p.vsz } ) ); + mcp.setEndpoint( _toVectorD( { p.vex, p.vey, p.vez } ) ); + + mcp.setTime(p.time); + mcp.setMass(p.mass); + mcp.setCharge(3.0*float(p.charge)); + + // Set generator status + mcp.setGeneratorStatus(0); + if( p.genStatus ) { + mcp.setGeneratorStatus( p.genStatus ) ; + } else { + if ( mask.isSet(sim::G4PARTICLE_GEN_STABLE) ) mcp.setGeneratorStatus(1); + else if ( mask.isSet(sim::G4PARTICLE_GEN_DECAYED) ) mcp.setGeneratorStatus(2); + else if ( mask.isSet(sim::G4PARTICLE_GEN_DOCUMENTATION) ) mcp.setGeneratorStatus(3); + else if ( mask.isSet(sim::G4PARTICLE_GEN_BEAM) ) mcp.setGeneratorStatus(4); + else if ( mask.isSet(sim::G4PARTICLE_GEN_OTHER) ) mcp.setGeneratorStatus(9); + } + + // Set simulation status + mcp.setCreatedInSimulation( mask.isSet(sim::G4PARTICLE_SIM_CREATED) ); + mcp.setBackscatter( mask.isSet(sim::G4PARTICLE_SIM_BACKSCATTER) ); + mcp.setVertexIsNotEndpointOfParent( mask.isSet(sim::G4PARTICLE_SIM_PARENT_RADIATED) ); + mcp.setDecayedInTracker( mask.isSet(sim::G4PARTICLE_SIM_DECAY_TRACKER) ); + mcp.setDecayedInCalorimeter( mask.isSet(sim::G4PARTICLE_SIM_DECAY_CALO) ); + mcp.setHasLeftDetector( mask.isSet(sim::G4PARTICLE_SIM_LEFT_DETECTOR) ); + mcp.setStopped( mask.isSet(sim::G4PARTICLE_SIM_STOPPED) ); + mcp.setOverlay( false ); + + //fg: if simstatus !=0 we have to set the generator status to 0: + if( mcp.isCreatedInSimulation() ) + mcp.setGeneratorStatus( 0 ); + + mcp.setSpin(p.spin); + mcp.setColorFlow(p.colorFlow); + } + + template <> template <> + void data_io<ddg4_input>::_to_edm4hep(const std::vector<const particle_t*>& cont, + edm4hep::MCParticleCollection* coll) + { + std::size_t i, n = cont.size(); + _pre_create(coll, n); + /// Convert particle body + for ( i=0; i<n; ++i) { + const particle_t* p = cont[i]; + auto mcp = coll->at(i); + _to_edm4hep(*p, mcp); + /// Relationships are already resolved and kept in order: Just copy indices + for (auto idau : p->daughters) + mcp.addToDaughters(coll->at(idau)); + for (auto ipar : p->parents) + mcp.addToParents(coll->at(ipar)); + } + } + + template <> template <> + void data_io<ddg4_input>::_to_digi_if(const std::vector<sim::Geant4Tracker::Hit*>& data, + std::map<CellID, std::shared_ptr<sim::Geant4Tracker::Hit> >& hits, + const DepositPredicate<EnergyCut>& predicate) { + for( auto* p : data ) { + std::shared_ptr<sim::Geant4Tracker::Hit> ptr(p); + if ( predicate(p) ) { + CellID cell = ptr->cellID; + hits.emplace(cell, std::move(ptr)); + } + } + } + + template <> template <> + void data_io<ddg4_input>::_to_digi_if(const std::vector<sim::Geant4Calorimeter::Hit*>& data, + std::map<CellID, std::shared_ptr<sim::Geant4Calorimeter::Hit> >& hits, + const DepositPredicate<EnergyCut>& predicate) { + for( auto* p : data ) { + std::shared_ptr<sim::Geant4Calorimeter::Hit> ptr(p); + if ( predicate(p) ) { + CellID cell = ptr->cellID; + hits.emplace(cell, std::move(ptr)); + } + } + } + + template <> template <> + void data_io<ddg4_input>::_to_digi(Key key, + const std::vector<sim::Geant4Particle*>& input, + ParticleMapping& particles) + { + Key mkey = key; + for( auto* part_ptr : input ) { + std::shared_ptr<sim::Geant4Particle> p(part_ptr); + Particle part; + part.start_position = Position(p->vsx, p->vsy, p->vsz); + part.end_position = Position(p->vex, p->vey, p->vez); + part.momentum = Direction(p->psx,p->psy, p->psz); + part.pdgID = p->pdgID; + part.charge = p->charge; + part.mass = p->mass; + part.time = p->time; + mkey.set_item(particles.size()); + part.source = std::make_any<std::shared_ptr<sim::Geant4Particle> >(std::move(p)); + particles.push(mkey, std::move(part)); + } + } + + template <typename T> + static void cnv_to_digi(Key key, + const std::pair<const CellID, std::shared_ptr<T> >& depo, + DepositVector& out) { + Key history_key; + EnergyDeposit dep { }; + const auto* h = depo.second.get(); + Position pos = h->position; + pos *= 1./dd4hep::mm; + + dep.flag = h->flag; + dep.deposit = h->energyDeposit; + dep.position = pos; + + history_key.set_mask(key.mask()); + history_key.set_item(out.size()); + history_key.set_segment(key.segment()); + dep.history.hits.emplace_back(history_key, dep.deposit); + add_particle_history(h, history_key, dep.history); + out.emplace(depo.first, std::move(dep)); + } + + template <> template <> + void data_io<ddg4_input>::_to_digi(Key key, + const std::map<CellID, std::shared_ptr<sim::Geant4Calorimeter::Hit> >& hits, + DepositVector& out) { + out.data_type = SegmentEntry::CALORIMETER_HITS; + for( const auto& p : hits ) + cnv_to_digi(key, p, out); + } + + template <> template <> + void data_io<ddg4_input>::_to_digi(Key key, + const std::map<CellID, std::shared_ptr<sim::Geant4Tracker::Hit> >& hits, + DepositVector& out) { + out.data_type = SegmentEntry::TRACKER_HITS; + for( const auto& p : hits ) + cnv_to_digi(key, p, out); + } + + } // End namespace digi +} // End namespace dd4hep +#endif // DD4HEP_USE_DDG4 + +/// Namespace for the AIDA detector description toolkit +namespace dd4hep { + + /// Namespace for the Digitization part of the AIDA detector description toolkit + namespace digi { + + template <> + void digi_io::_to_edm4hep(const ParticleMapping& cont, + edm4hep::MCParticleCollection* coll) + { + if ( cont.empty() ) { + return; + } + if ( data_io<edm4hep_input>::_can_handle(*cont.begin()) ) { + data_io<edm4hep_input> io; + auto vec = io._to_vector(cont); + if ( !vec.empty() ) { + io._to_edm4hep(vec, coll); + } + return; + } + else if ( data_io<ddg4_input>::_can_handle(*cont.begin()) ) { + data_io<ddg4_input> io; + auto vec = io._to_vector(cont); + if ( !vec.empty() ) { + io._to_edm4hep(vec, coll); + } + return; + } + // Catch-all: convert what we have at hands + data_io<digi_input> io; + auto vec = io._to_vector(cont); + if ( !vec.empty() ) { + io._to_edm4hep(vec, coll); + } + } + + template <> + void digi_io::_to_edm4hep(const std::pair<const CellID, EnergyDeposit>& dep, + const std::array<float, 6>& covMat, + edm4hep::TrackerHitCollection* collection) + { + const EnergyDeposit& de = dep.second; + auto hit = collection->create(); + hit.setType(0 /* edm4hep::SIMTRACKERHIT */); + hit.setTime(de.time); + hit.setCovMatrix(covMat); + hit.setCellID(dep.first); + hit.setEDep(de.deposit); + hit.setEDepError(de.deposit/1e0); + hit.setEdx(de.deposit/de.length); + hit.setPosition( _toVectorD(de.position) ); + } + + template <> + void digi_io::_to_edm4hep(const std::pair<const CellID, EnergyDeposit>& dep, + edm4hep::CalorimeterHitCollection* collection) + { + const EnergyDeposit& de = dep.second; + auto hit = collection->create(); + hit.setType(0 /* edm4hep::SIMCALORIMETERHIT */); + hit.setTime(de.time); + hit.setCellID(dep.first); + hit.setEnergy(de.deposit); + hit.setEnergyError(de.deposit/10e0); + hit.setPosition( _toVectorF(de.position) ); + } + } // End namespace digi +} // End namespace dd4hep +#endif // DD4HEP_USE_EDM4HEP diff --git a/DDDigi/io/DigiIO.h b/DDDigi/io/DigiIO.h new file mode 100644 index 0000000000000000000000000000000000000000..84cd054ccfaa9734a332605e335d4ec91886a2e2 --- /dev/null +++ b/DDDigi/io/DigiIO.h @@ -0,0 +1,131 @@ +//========================================================================== +// AIDA Detector description implementation +//-------------------------------------------------------------------------- +// Copyright (C) Organisation europeenne pour la Recherche nucleaire (CERN) +// All rights reserved. +// +// For the licensing terms see $DD4hepINSTALL/LICENSE. +// For the list of contributors see $DD4hepINSTALL/doc/CREDITS. +// +// Author : M.Frank +// +//========================================================================== +#ifndef DDDIGI_DIGIIO_H +#define DDDIGI_DIGIIO_H + +// Framework include files +#include <DD4hep/Printout.h> +#include <DD4hep/DD4hepUnits.h> +#include <DDDigi/DigiData.h> + +// C/C++ include files +#include <any> +#include <array> +#include <stdexcept> + +/// Namespace for the AIDA detector description toolkit +namespace dd4hep { + + /// Namespace for the Digitization part of the AIDA detector description toolkit + namespace digi { + + template <typename T> union input_data { + const void* m_raw; + std::vector<T*>* m_items; + input_data(const void* p) { this->m_raw = p; } + void clear() { if ( this->m_items ) this->m_items->clear(); } + std::size_t size() { return (this->m_items) ? this->m_items->size() : 0UL; } + std::vector<T*>& get() { + if ( this->m_items ) return *(this->m_items); + throw std::runtime_error("input_data: Invalid data!"); + } + }; + } // End namespace digi +} // End namespace dd4hep + + +#ifdef DD4HEP_USE_EDM4HEP + +/// edm4hep include files +#include <edm4hep/MCParticleCollection.h> +#include <edm4hep/TrackerHitCollection.h> +#include <edm4hep/SimTrackerHitCollection.h> +#include <edm4hep/CalorimeterHitCollection.h> +#include <edm4hep/SimCalorimeterHitCollection.h> + +/// Namespace for the AIDA detector description toolkit +namespace dd4hep { + + /// Namespace for the Geant4 based simulation part of the AIDA detector description toolkit + namespace sim { + // Forward declarations + class Geant4Particle; + } + + /// Namespace for the Digitization part of the AIDA detector description toolkit + namespace digi { + + edm4hep::Vector3d _toVectorD(const Position& ep); + edm4hep::Vector3f _toVectorF(const Position& ep); + + struct digi_io { + public: + digi_io() = default; + template <typename FIRST_TYPE, typename OUTPUT_TYPE> static + void _to_edm4hep(const FIRST_TYPE& first, OUTPUT_TYPE output); + + template <typename FIRST_TYPE, typename SECOND_TYPE, typename OUTPUT_TYPE> static + void _to_edm4hep(const FIRST_TYPE& first, const SECOND_TYPE& second, OUTPUT_TYPE output); + }; + + struct ddg4_input { + typedef sim::Geant4Particle particle_type; + struct input_trackerhit_type {}; + struct input_calorimeterhit_type {}; + }; + + struct edm4hep_input { + typedef edm4hep::MutableMCParticle particle_type; + struct input_trackerhit_type {}; + struct input_calorimeterhit_type {}; + }; + + struct digi_input { + typedef Particle particle_type; + struct input_trackerhit_type {}; + struct input_calorimeterhit_type {}; + }; + + template <typename T> struct data_io { + using particle_t = typename T::particle_type; + using trackerhit_t = typename T::input_trackerhit_type; + using calorimeterhit_t = typename T::input_calorimeterhit_type; + using pwrap_t = std::shared_ptr<particle_t>; + using twrap_t = std::shared_ptr<trackerhit_t>; + using cwrap_t = std::shared_ptr<calorimeterhit_t>; + + data_io() = default; + + template <typename DATA> + static bool _can_handle(const DATA& data); + + template <typename CONT> static + std::vector<const particle_t*> _to_vector(const CONT& cont); + + template <typename CONT> static + void _pre_create(CONT* coll, std::size_t n); + + template <typename FIRST, typename SECOND> static + void _to_edm4hep(const FIRST& cont, SECOND coll); + + template <typename FIRST, typename SECOND, typename THIRD> static + void _to_digi(FIRST first, const SECOND& second, THIRD& third); + + template <typename FIRST, typename SECOND, typename PREDICATE> static + void _to_digi_if(const FIRST& first, SECOND& second, const PREDICATE& pred); + }; + } // End namespace digi +} // End namespace dd4hep +#endif + +#endif // DDDIGI_DIGIIO_H diff --git a/DDDigi/plugins/DigiDepositMapCreator.cpp b/DDDigi/plugins/DigiDepositMapCreator.cpp index 19fcf5180ef4f582df3cf190a8b7afc6b9e42eff..72bd2d316146be5ae66e454312e1290f6231c7e1 100644 --- a/DDDigi/plugins/DigiDepositMapCreator.cpp +++ b/DDDigi/plugins/DigiDepositMapCreator.cpp @@ -38,7 +38,7 @@ namespace dd4hep { template <typename T> void create_deposits(const char* tag, const T& cont, work_t& work, const predicate_t& predicate) const { Key key(cont.name, work.environ.output.mask); - DepositMapping m(cont.name, work.environ.output.mask); + DepositMapping m(cont.name, work.environ.output.mask, cont.data_type); std::size_t start = m.size(); for( const auto& dep : cont ) { if ( predicate(dep) ) { @@ -62,6 +62,6 @@ namespace dd4hep { }; } // End namespace digi } // End namespace dd4hep - +/// Factory instantiation: #include <DDDigi/DigiFactories.h> DECLARE_DIGIACTION_NS(dd4hep::digi,DigiDepositMapCreator) diff --git a/DDDigi/plugins/DigiDepositNoiseOnSignal.cpp b/DDDigi/plugins/DigiDepositNoiseOnSignal.cpp index 62d63043505d3c258f9419dc035c13dd25aa474e..be3f3e7a25eacee73289fcc5a8b07de351d02886 100644 --- a/DDDigi/plugins/DigiDepositNoiseOnSignal.cpp +++ b/DDDigi/plugins/DigiDepositNoiseOnSignal.cpp @@ -69,6 +69,6 @@ namespace dd4hep { }; } // End namespace digi } // End namespace dd4hep - +/// Factory instantiation: #include <DDDigi/DigiFactories.h> DECLARE_DIGIACTION_NS(dd4hep::digi,DigiDepositNoiseOnSignal) diff --git a/DDDigi/plugins/DigiDepositSmearEnergy.cpp b/DDDigi/plugins/DigiDepositSmearEnergy.cpp index d7d97aa6a96f9511dd0f771e18564a61bc3652cf..25e8e6b4290d6b2294b54194ec49738059f1dc09 100644 --- a/DDDigi/plugins/DigiDepositSmearEnergy.cpp +++ b/DDDigi/plugins/DigiDepositSmearEnergy.cpp @@ -55,6 +55,8 @@ namespace dd4hep { * Sampling resolution: sigma(E)/E ~ 0.09 * sqrt(1000*dE/E) * dE = energy loss of a single charged particle in one sampling layer * + * (See all: https://handwiki.org/wiki/Physics:Energy_resolution_in_calorimeters) + * * * \author M.Frank * \version 1.0 diff --git a/DDDigi/plugins/DigiDepositWeightedPosition.cpp b/DDDigi/plugins/DigiDepositWeightedPosition.cpp index 6fdd9bc6fc1ed0d9464ad510d32fa0eb9e935de2..a0cc68750005b21311a74c639c4d0d5e3684f7b8 100644 --- a/DDDigi/plugins/DigiDepositWeightedPosition.cpp +++ b/DDDigi/plugins/DigiDepositWeightedPosition.cpp @@ -44,7 +44,7 @@ namespace dd4hep { template <typename T> void create_deposits(context_t& context, const T& cont, work_t& work, const predicate_t& predicate) const { Key key(cont.name, work.environ.output.mask); - DepositMapping m(cont.name, work.environ.output.mask); + DepositMapping m(cont.name, work.environ.output.mask, cont.data_type ); std::size_t dropped = 0UL, updated = 0UL, added = 0UL; for( const auto& dep : cont ) { if ( predicate(dep) ) { diff --git a/DDDigi/plugins/DigiIPMover.cpp b/DDDigi/plugins/DigiIPMover.cpp index 804089bf03ba8a329f039d869b7152aaa83a4606..60748053243f84d6e89e0614c12e97d7a7791b15 100644 --- a/DDDigi/plugins/DigiIPMover.cpp +++ b/DDDigi/plugins/DigiIPMover.cpp @@ -59,11 +59,9 @@ namespace dd4hep { move_particles(const char* tag, PARTICLES& cont, const Position& delta) const { info("%s+++ %-32s [%6ld] IP: x:%7.3f y:%7.3f z:%7.3f", tag, cont.name.c_str(), cont.size(), delta.X(), delta.Y(), delta.Z()); - for( auto& p : cont ) { - auto& part = p.second; - part.end_position += delta; - part.start_position += delta; - } + /// Move particle + for( auto& p : cont ) + p.second.move_position(delta); return cont.size(); } @@ -84,6 +82,6 @@ namespace dd4hep { }; } // End namespace digi } // End namespace dd4hep - +/// Factory instantiation: #include <DDDigi/DigiFactories.h> DECLARE_DIGIACTION_NS(dd4hep::digi,DigiIPMover) diff --git a/DDDigi/plugins/DigiResegment.cpp b/DDDigi/plugins/DigiResegment.cpp index 6a4d3df9b1e2fdd4df1ee9a2931fdec6a1f2acbb..07f37deb32f691a3e0595f0c1672efd454b6dd92 100644 --- a/DDDigi/plugins/DigiResegment.cpp +++ b/DDDigi/plugins/DigiResegment.cpp @@ -106,7 +106,7 @@ namespace dd4hep { template <typename T> void resegment_deposits(const T& cont, work_t& work, const predicate_t& predicate) const { Key key(cont.name, work.environ.output.mask); - DepositVector m(cont.name, key.mask()); + DepositVector m(cont.name, key.mask(), cont.data_type); std::size_t start = m.size(); for( const auto& dep : cont ) { if( predicate(dep) ) { diff --git a/DDDigi/plugins/DigiSegmentDepositExtractor.cpp b/DDDigi/plugins/DigiSegmentDepositExtractor.cpp index 17ea0c3e767b5723d139d37c4be2eb65df8c2c69..af8d2e63245b591b28f0b4545b8adc351ea3ad70 100644 --- a/DDDigi/plugins/DigiSegmentDepositExtractor.cpp +++ b/DDDigi/plugins/DigiSegmentDepositExtractor.cpp @@ -38,14 +38,14 @@ namespace dd4hep { : DigiContainerProcessor(kernel, nam) {} template <typename T> void copy_deposits(const T& cont, work_t& work, const predicate_t& predicate) const { - DepositVector deposits(cont.name, work.environ.output.mask); - for( const auto& dep : cont ) { - if( predicate(dep) ) { - CellID cell = dep.first; - EnergyDeposit depo = dep.second; - deposits.emplace(cell, std::move(depo)); - } - } + DepositVector deposits(cont.name, work.environ.output.mask, cont.data_type); + for( const auto& dep : cont ) { + if( predicate(dep) ) { + CellID cell = dep.first; + EnergyDeposit depo = dep.second; + deposits.emplace(cell, std::move(depo)); + } + } work.environ.output.data.put(deposits.key, std::move(deposits)); } /// Main functional callback diff --git a/DDDigi/python/DDDigiDict.C b/DDDigi/python/DDDigiDict.C index 23a654e59c5ffb5b2b0916655059f093369e3ab9..d5b8067d001113458079cd1cc8c5a5f46fffef37 100644 --- a/DDDigi/python/DDDigiDict.C +++ b/DDDigi/python/DDDigiDict.C @@ -128,6 +128,7 @@ namespace dd4hep { static DigiEventAction* toEventAction(DigiAction* a) { return cst<DigiEventAction>(a); } static DigiDepositMonitor* toDepositMonitor(DigiAction* a) { return cst<DigiDepositMonitor>(a); } static DigiContainerProcessor* toContainerProcessor(DigiAction* a) { return cst<DigiContainerProcessor>(a); } + static DigiContainerSequenceAction* toContainerSequenceAction(DigiAction* a) { return cst<DigiContainerSequenceAction>(a); } /// Access DigiAction property static PropertyResult getProperty(DigiAction* action, const std::string& name) { diff --git a/DDDigi/python/dddigi.py b/DDDigi/python/dddigi.py index da0b31f0aab8bcac513f44de3e1de24f69dd6a70..4d35ef7b553f5136e29f1e8094bd8595666d2ab3 100644 --- a/DDDigi/python/dddigi.py +++ b/DDDigi/python/dddigi.py @@ -42,10 +42,10 @@ def loadDDDigi(): logger.info('DDDigi.py: Successfully loaded DDDigi plugin library libDDDigiPlugins!') # # import with ROOT the I/O module to read DDG4 output - result = gSystem.Load("libDDDigi_DDG4_IO") + result = gSystem.Load("libDDDigi_IO") if result < 0: - raise Exception('DDDigi.py: Failed to load the DDG4 IO library libDDDigi_DDG4_IO: ' + gSystem.GetErrorStr()) - logger.info('DDDigi.py: Successfully loaded DDG4 IO plugin library libDDDigi_DDG4_IO!') + raise Exception('DDDigi.py: Failed to load the DDDigi IO library libDDDigi_IO: ' + gSystem.GetErrorStr()) + logger.info('DDDigi.py: Successfully loaded DDDigi IO plugin library libDDDigi_IO!') # # import the main dd4hep module from ROOT from ROOT import dd4hep as module @@ -216,7 +216,8 @@ def _adopt_event_action(self, action): def _adopt_container_processor(self, action, processor_argument): " Helper to convert DigiActions objects to DigiEventAction " - attr = getattr(_get_action(self), 'adopt_processor') + parent = Interface.toContainerSequenceAction(_get_action(self)) + attr = getattr(parent, 'adopt_processor') proc = Interface.toContainerProcessor(_get_action(action)) attr(proc, processor_argument) # --------------------------------------------------------------------------- @@ -325,7 +326,8 @@ _props('ActionHandle', add_set_property=_add_new_set_property, add_list_property=_add_new_list_property, add_vector_property=_add_new_vector_property, - add_mapped_property=_add_new_mapped_property) + add_mapped_property=_add_new_mapped_property, + adopt_container_processor=_adopt_container_processor) _props('DigiSynchronize', adopt=_adopt_event_action, adopt_action=_adopt_sequence_action) _props('DigiActionSequence', adopt=_adopt_event_action, adopt_action=_adopt_sequence_action) _props('DigiParallelActionSequence', adopt_action=_adopt_sequence_action) @@ -333,6 +335,7 @@ _props('DigiSequentialActionSequence', adopt_action=_adopt_sequence_action) _props('DigiContainerSequenceAction', adopt_container_processor=_adopt_container_processor) _props('DigiMultiContainerProcessor', adopt_processor=_adopt_processor) _props('DigiSegmentSplitter', adopt_segment_processor=_adopt_segment_processor) + # --------------------------------------------------------------------------- diff --git a/DDDigi/python/digitize.py b/DDDigi/python/digitize.py index 86673200e73ff0a2a99700e15541100204d24f5d..a1098575cc06eb9a8e23aab7da473724a3856023 100644 --- a/DDDigi/python/digitize.py +++ b/DDDigi/python/digitize.py @@ -32,12 +32,18 @@ class Digitize(dd4hep.Logger, dd4hep.CommandLine): """ def __init__(self, kernel=None): dd4hep.Logger.__init__(self, 'dddigi') + dd4hep.CommandLine.__init__(self) self._kernel = kernel self._main_processor = None self._input_processor = None self._event_processor = None self._parallel = True + self._dddigi = dddigi self.description = self._kernel.detectorDescription() + if self.output_level: + lvl = int(self.output_level) + self.setPrintLevel(lvl) + self._kernel.OutputLevel = lvl """ Access the worker kernel object. @@ -154,7 +160,7 @@ class Digitize(dd4hep.Logger, dd4hep.CommandLine): self.info('+++ List of sensitive detectors:') dets = self.activeDetectors() for d in dets: - self.info('+++ %-32s ---> type:%-12s' % (d['name'], d['type'],)) + self.always('+++ %-32s ---> type:%-12s' % (d['name'], d['type'],)) """ Configure ROOT output for the event digitization diff --git a/DDDigi/src/DigiAttenuator.cpp b/DDDigi/src/DigiAttenuator.cpp index 1d63f324e54388dcbff2a74209486ae531b06e98..62845b03b0e4c7225c1fe7474f5dece37043e2d0 100644 --- a/DDDigi/src/DigiAttenuator.cpp +++ b/DDDigi/src/DigiAttenuator.cpp @@ -126,6 +126,7 @@ void DigiAttenuatorSequence::initialize() { att->property("factor").set(factor); att->property("OutputLevel").set(int(outputLevel())); adopt_processor(att, c.first); + att->release(); // Release processor **after** adoption. } this->DigiContainerSequenceAction::initialize(); } diff --git a/DDDigi/src/DigiContainerCombine.cpp b/DDDigi/src/DigiContainerCombine.cpp index 9608facdd1a2a7ac69ce766456b619a03b084c59..ff7d3dbb6adfaf4077f8087c4aace7300f1c8c31 100644 --- a/DDDigi/src/DigiContainerCombine.cpp +++ b/DDDigi/src/DigiContainerCombine.cpp @@ -86,7 +86,10 @@ public: std::size_t cnt = 0; const auto& nam = input.name; Key::mask_type mask = input.key.mask(); - + if ( output.data_type == SegmentEntry::UNKNOWN ) + output.data_type = input.data_type; + else if ( output.data_type != input.data_type ) + combine->except("+++ Digitization does not allow to mix data of different type!"); if ( combine->m_erase_combined ) { cnt = output.merge(std::move(input)); } @@ -101,7 +104,7 @@ public: /// Generic deposit merger: implicitly assume identical item types are mapped sequentially void merge(const std::string& nam, size_t start, int thr) { Key key = keys[start]; - DepositVector out(nam, combine->m_deposit_mask); + DepositVector out(nam, combine->m_deposit_mask, SegmentEntry::UNKNOWN); for( std::size_t j = start; j < keys.size(); ++j ) { if ( keys[j].item() == key.item() ) { if ( DepositMapping* m = std::any_cast<DepositMapping>(work[j]) ) diff --git a/DDDigi/src/DigiContainerProcessor.cpp b/DDDigi/src/DigiContainerProcessor.cpp index 9dfe9021508b88f97a2bc8d01d92390abcf9f857..d41894678d92fabe8f3b8d2ff1831cbeaf5d0722 100644 --- a/DDDigi/src/DigiContainerProcessor.cpp +++ b/DDDigi/src/DigiContainerProcessor.cpp @@ -143,6 +143,11 @@ DigiContainerSequence::~DigiContainerSequence() { InstanceCount::decrement(this); } +/// Set the default predicate +void DigiContainerSequence::set_predicate(const predicate_t& predicate) { + m_worker_predicate = predicate; +} + /// Adopt new parallel worker void DigiContainerSequence::adopt_processor(DigiContainerProcessor* action) { if ( !action ) { @@ -176,6 +181,7 @@ DigiContainerSequenceAction::DigiContainerSequenceAction(const kernel_t& krnl, c declareProperty("output_mask", m_output_mask); declareProperty("output_segment", m_output_segment); m_kernel.register_initialize(std::bind(&DigiContainerSequenceAction::initialize,this)); + m_kernel.register_terminate(std::bind(&DigiContainerSequenceAction::finalize,this)); InstanceCount::increment(this); } @@ -194,6 +200,11 @@ void dd4hep::digi::DigiContainerSequenceAction::initialize() { } } +/// Finalization callback +void dd4hep::digi::DigiContainerSequenceAction::finalize() { + m_workers.clear(); +} + /// Set the default predicate void DigiContainerSequenceAction::set_predicate(const predicate_t& predicate) { m_worker_predicate = predicate; @@ -203,7 +214,7 @@ void DigiContainerSequenceAction::set_predicate(const predicate_t& predicate) void DigiContainerSequenceAction::adopt_processor(DigiContainerProcessor* action, const std::string& container) { - Key key(container, 0x0); + Key key(container.c_str(), 0x0); auto it = m_registered_processors.find(key); if ( it != m_registered_processors.end() ) { if ( action != it->second ) { @@ -219,7 +230,7 @@ void DigiContainerSequenceAction::adopt_processor(DigiContainerProcessor* action action->addRef(); m_registered_processors.emplace(key, action); info("+++ Adding processor: %s for container: [%08X] %s", - action->c_name(), key.item(), container.c_str()); + action->c_name(), key.item(), ('"'+container+'"').c_str()); } /// Adopt new parallel worker acting on multiple containers @@ -234,8 +245,9 @@ void DigiContainerSequenceAction::adopt_processor(DigiContainerProcessor* action /// Get hold of the registered processor for a given container DigiContainerSequenceAction::worker_t* -DigiContainerSequenceAction::need_registered_worker(Key item_key, bool exc) const { - item_key.set_mask(0); +DigiContainerSequenceAction::need_registered_worker(Key key, bool exc) const { + Key item_key; + item_key.set_item(key.item()); auto it = m_registered_workers.find(item_key); if ( it != m_registered_workers.end() ) { return it->second; @@ -249,14 +261,14 @@ DigiContainerSequenceAction::need_registered_worker(Key item_key, bool exc) co /// Main functional callback if specific work is known void DigiContainerSequenceAction::execute(context_t& context) const { std::vector<ParallelWorker*> event_workers; + work_items_t items; auto& event = *context.event; auto& input = event.get_segment(m_input_segment); auto& output = event.get_segment(m_output_segment); - output_t out { m_output_mask, output }; - env_t env { context, m_properties, out }; - work_items_t items; - work_t arg { env, items, *this }; + output_t out { m_output_mask, output }; + env_t env { context, m_properties, out }; work_item_t itm { nullptr, { }, nullptr }; + work_t arg { env, items, *this }; arg.input_items.resize(m_workers.size(), itm); event_workers.reserve(input.size()); diff --git a/DDDigi/src/DigiData.cpp b/DDDigi/src/DigiData.cpp index c3ddff9f7daf9619c3fb76e888b7e02af67ee411..f89f3560d61cf11d6b14197fe1958318ace6a589 100644 --- a/DDDigi/src/DigiData.cpp +++ b/DDDigi/src/DigiData.cpp @@ -41,6 +41,7 @@ std::string dd4hep::digi::digiTypeName(const std::type_info& info) { if ( idx != std::string::npos ) typ = str_replace(typ, typ.substr(idx), ">"); typ = str_replace(str_replace(typ,"std::",""),"dd4hep::",""); typ = str_replace(str_replace(typ,"sim::",""),"digi::",""); + typ = str_replace(str_replace(typ,", less<Key>",""),", less<long long>",""); return typ; } @@ -49,14 +50,20 @@ std::string dd4hep::digi::digiTypeName(const std::any& data) { } Key& Key::set(const std::string& name, int mask) { + return this->set(name, 0x0, mask); +} + +/// Generate key using hash algorithm +Key& Key::set(const std::string& name, int segment, int mask) { auto& _k = keys(); if ( name.empty() ) { std::lock_guard<std::mutex> lock(_k.lock); except("DDDigi::Key", "+++ No key name was specified -- this is illegal!"); } this->key = 0; - this->set_mask(Key::mask_type(0xFFFF&mask)); this->set_item(detail::hash32(name)); + this->set_mask(Key::mask_type(0xFFFF&mask)); + this->set_segment(Key::segment_type(0xFF&segment)); std::lock_guard<std::mutex> lock(_k.lock); _k.map[this->item()] = name; return *this; @@ -301,6 +308,12 @@ void DepositMapping::remove(iterator position) { data.erase(position); } +/// Move particle +void Particle::move_position(const Position& delta) { + this->start_position += delta; + this->end_position += delta; +} + /// Merge new deposit map onto existing map std::size_t ParticleMapping::insert(const ParticleMapping& updates) { std::size_t update_size = updates.size(); @@ -407,14 +420,6 @@ bool DataSegment::emplace_any(Key key, std::any&& item) { return ret; } -/// Move data items other than std::any to the data segment -template <typename DATA> bool DataSegment::put(Key key, DATA&& value) { - key.set_segment(this->id); - value.key.set_segment(this->id); - std::any item = std::make_any<DATA>(std::move(value)); - return this->emplace_any(key, std::move(item)); -} - template bool DataSegment::put(Key key, DepositVector&& data); template bool DataSegment::put(Key key, DepositMapping&& data); template bool DataSegment::put(Key key, ParticleMapping&& data); diff --git a/DDDigi/src/DigiKernel.cpp b/DDDigi/src/DigiKernel.cpp index 5a9f2807a34229ed1d10c8e900b616481d117b62..92ea50c06c9cad726663ed6a94035f17f0e1afa0 100644 --- a/DDDigi/src/DigiKernel.cpp +++ b/DDDigi/src/DigiKernel.cpp @@ -39,11 +39,7 @@ namespace tbb { struct global_control { enum { max_allowed_parallelism = -1 }; #include <memory> #include <chrono> -using namespace dd4hep; using namespace dd4hep::digi; -namespace { - static std::mutex kernel_mutex; -} /// DigiKernel herlp class: Container of instance variabled /* @@ -120,8 +116,12 @@ public: Internals() = default; /// Default destructor ~Internals() = default; + + static std::mutex kernel_mutex; }; +std::mutex DigiKernel::Internals::kernel_mutex {}; + /// DigiKernel herlp class: TBB wrapper to execute DigiAction instances /* * @@ -205,7 +205,7 @@ DigiKernel::DigiKernel(Detector& description_ref) /// Default destructor DigiKernel::~DigiKernel() { - std::lock_guard<std::mutex> lock(kernel_mutex); + std::lock_guard<std::mutex> lock(Internals::kernel_mutex); internals->tbb_init.reset(); detail::releasePtr(internals->monitor_handler); detail::releasePtr(internals->output_action); @@ -221,7 +221,7 @@ DigiKernel::~DigiKernel() { DigiKernel& DigiKernel::instance(Detector& description) { static dd4hep::dd4hep_ptr<DigiKernel> s_main_instance(0); if ( 0 == s_main_instance.get() ) { - std::lock_guard<std::mutex> lock(kernel_mutex); + std::lock_guard<std::mutex> lock(Internals::kernel_mutex); if ( 0 == s_main_instance.get() ) { // Need to check again! s_main_instance.adopt(new DigiKernel(description)); } @@ -367,18 +367,27 @@ void DigiKernel::register_monitor(DigiAction* action, TNamed* object) const /// Submit a bunch of actions to be executed in parallel void DigiKernel::submit (DigiContext& context, ParallelCall*const algorithms[], std::size_t count, void* data, bool parallel) const { - (void)parallel; // Silence compiler warning when not using TBB const char* tag = context.event->id(); #ifdef DD4HEP_USE_TBB bool para = parallel && (internals->tbb_init && internals->num_threads > 0); if ( para ) { tbb::task_group que; info("%s+++ Executing chunk of %3ld execution entries in parallel", tag, count); - for( std::size_t i=0; i<count; ++i) - que.run( Wrapper<ParallelCall,void*>(algorithms[i], data) ); - que.wait(); + try { + for( std::size_t i=0; i<count && !internals->stop; ++i) + que.run( Wrapper<ParallelCall,void*>(algorithms[i], data) ); + que.wait(); + } + catch(const std::exception& e) { + std::exception_ptr eptr = std::current_exception(); + internals->stop = true; + error("%s+++ C++ exception. STOP event loop. [%s]", tag, e.what()); + std::rethrow_exception(eptr); + } return; } +#else + (void)parallel; // Silence compiler warning when not using TBB #endif info("%s+++ Executing chunk of %3ld execution entries sequentially", tag, count); for( std::size_t i=0; i<count; ++i) @@ -438,22 +447,29 @@ int DigiKernel::run() { info("+++ Total number of events: %d",internals->numEvents); #ifdef DD4HEP_USE_TBB if ( !internals->tbb_init && internals->num_threads > 0 ) { - using ctrl_t = tbb::global_control; - if ( 0 == internals->num_threads ) { - internals->num_threads = ctrl_t::max_allowed_parallelism; - } - info("+++ Number of TBB threads: %d",internals->num_threads); - info("+++ Number of parallel events: %d",internals->maxEventsParallel); - internals->tbb_init = std::make_unique<ctrl_t>(ctrl_t::max_allowed_parallelism,internals->num_threads+1); - if ( internals->maxEventsParallel >= 0 ) { - int todo_evt = internals->events_todo; - int num_proc = std::min(todo_evt,internals->maxEventsParallel); - tbb::task_group main_group; - for(int i=0; i < num_proc; ++i) - main_group.run(Processor(*this)); - main_group.wait(); - } - debug("+++ All event processing threads Synchronized --- Done!"); + using ctrl_t = tbb::global_control; + if ( 0 == internals->num_threads ) { + internals->num_threads = ctrl_t::max_allowed_parallelism; + } + info("+++ Number of TBB threads: %d",internals->num_threads); + info("+++ Number of parallel events: %d",internals->maxEventsParallel); + internals->tbb_init = std::make_unique<ctrl_t>(ctrl_t::max_allowed_parallelism,internals->num_threads+1); + if ( internals->maxEventsParallel >= 0 ) { + int todo_evt = internals->events_todo; + int num_proc = std::min(todo_evt,internals->maxEventsParallel); + tbb::task_group main_group; + try { + for(int i=0; i < num_proc; ++i) + main_group.run(Processor(*this)); + main_group.wait(); + } + catch(const std::exception& e) { + internals->stop = true; + error("run: +++ C++ exception. Event loop stop. [%s]", e.what()); + main_group.wait(); + } + } + debug("+++ All event processing threads Synchronized --- Done!"); } else #endif @@ -464,6 +480,7 @@ int DigiKernel::run() { ++internals->events_submitted; } } + std::chrono::duration<double> duration = std::chrono::system_clock::now() - start; double sec = std::chrono::duration_cast<std::chrono::seconds>(duration).count(); info("+++ %d Events out of %d processed. " diff --git a/DDDigi/src/DigiStoreDump.cpp b/DDDigi/src/DigiStoreDump.cpp index d8a0b37a3693e7fda6dcf61c71ee2260ba813595..84c4744ee1e7150ff091170fea24e94ffe93141c 100644 --- a/DDDigi/src/DigiStoreDump.cpp +++ b/DDDigi/src/DigiStoreDump.cpp @@ -55,15 +55,22 @@ void DigiStoreDump::initialize() { template <typename T> std::string DigiStoreDump::data_header(Key key, const std::string& tag, const T& container) const { return this->format("%04X %04X %08X %-32s: %6ld %-12s [%s]", key.segment(), key.mask(), key.item(), - Key::key_name(key).c_str(), container.size(), tag.c_str(), + ('"'+Key::key_name(key)+'"').c_str(), container.size(), tag.c_str(), digiTypeName(typeid(T)).c_str()); } template <> std::string DigiStoreDump::data_header(Key key, const std::string& tag, const std::type_info& info) const { + std::string typ = digiTypeName(info); + if ( tag.empty() ) { + return this->format("%04X %04X %08X %-32s: %6s %s", + key.segment(), key.mask(), key.item(), + ('"'+Key::key_name(key)+'"').c_str(), "", + typ.c_str()); + } return this->format("%04X %04X %08X %-32s: %-12s %s", key.segment(), key.mask(), key.item(), - Key::key_name(key).c_str(), tag.c_str(), - digiTypeName(info).c_str()); + ('"'+Key::key_name(key)+'"'), tag.c_str(), + typ.c_str()); } template <> std::string DigiStoreDump::data_header(Key key, const std::string& tag, const std::any& data) const { @@ -78,9 +85,7 @@ DigiStoreDump::dump_history(DigiContext& context, { std::stringstream str; auto& ev = *(context.event); - const auto& item = data.second; Key particle_key = data.first; - Key history_key = item.history; std::vector<std::string> records; History::hist_entry_t hist { particle_key, 1.0 }; const Particle& par = hist.get_particle(ev); @@ -89,9 +94,8 @@ DigiStoreDump::dump_history(DigiContext& context, str << Key::key_name(key) << "[" << seq_no << "]:"; std::string line = - format("|---- %-18s Segment:%04X Mask:%04X Item:%04X History: %s Segment:%04X Mask:%04X Item:%04X %016lX", + format("|---- %-18s Segment:%04X Mask:%04X Item:%04X %016lX", str.str().c_str(), particle_key.segment(), particle_key.mask(), particle_key.item(), - Key::key_name(history_key).c_str(), history_key.segment(), history_key.mask(), history_key.item(), long(&data.second)); records.emplace_back(line); line = format("| PDG:%6d Charge:%-2d Mass:%7.2f v:%7.5g %7.5g %7.5g p:%12g %12g %12g %016lX", @@ -295,29 +299,35 @@ void DigiStoreDump::dump_headers(const std::string& tag, const DigiEvent& event, const DataSegment& segment) const { + int first = 1; std::string str; std::vector<std::string> records; - info("%s+--- %-12s segment: %ld entries", event.id(), tag.c_str(), segment.size()); std::lock_guard<std::mutex> lock(segment.lock); + records.push_back(format("+--- %-12s segment: %ld entries", tag.c_str(), segment.size())); + records.push_back(format("| Segt Mask Item-id Item-name")); for ( const auto& entry : segment ) { Key key {entry.first}; const std::any& data = entry.second; if ( const auto* mapping = std::any_cast<DepositMapping>(&data) ) - str = data_header(key, "deposits", *mapping); + str = "| " + data_header(key, "deposits", *mapping); else if ( const auto* vector = std::any_cast<DepositVector>(&data) ) - str = data_header(key, "deposits", *vector); + str = "| " + data_header(key, "deposits", *vector); else if ( const auto* parts = std::any_cast<ParticleMapping>(&data) ) - str = data_header(key, "particles", *parts); + str = "| " + data_header(key, "particles", *parts); else if ( const auto* adcs = std::any_cast<DetectorResponse>(&data) ) - str = data_header(key, "ADC values", *adcs); + str = "| " + data_header(key, "ADC values", *adcs); else if ( const auto* hist = std::any_cast<DetectorHistory>(&data) ) - str = data_header(key, "histories", *hist); + str = "| " + data_header(key, "histories", *hist); else if ( data.type() == typeid(void) ) - str = data_header(key, "void data", data); + str = "| " + data_header(key, "void data", data); else - str = data_header(key, "", data); + str = "| " + data_header(key, "", data); + if ( first ) { + first = false; + } records.push_back(str); } + if ( records.size() == 2 ) records.pop_back(); std::lock_guard<std::mutex> record_lock(m_kernel.global_output_lock()); for(const auto& s : records) info("%s|---- %s", event.id(), s.c_str()); @@ -335,3 +345,4 @@ void DigiStoreDump::execute(DigiContext& context) const { dump_headers(seg, *event, event->get_segment(segment)); } } + diff --git a/examples/DDDigi/CMakeLists.txt b/examples/DDDigi/CMakeLists.txt index 0fb7b52cf50c105370a4913a921e0751c9681697..76dac9f350fe5f5ec99edf47c48722475327c778 100644 --- a/examples/DDDigi/CMakeLists.txt +++ b/examples/DDDigi/CMakeLists.txt @@ -48,6 +48,7 @@ dd4hep_add_test_reg(DDDigi_framework dd4hep_add_test_reg(DDDigi_colored_noise COMMAND "${CMAKE_INSTALL_PREFIX}/bin/run_test_DDDigi.sh" EXEC_ARGS geoPluginRun -ui -plugin DD4hep_FalphaNoise -shots 1000000 -variance 1 -alpha 0.5 + DEPENDS DDDigi_framework REGEX_PASS "FalphaNoise INFO Distribution RMS 1.0" REGEX_FAIL "Error;ERROR;Exception" ) @@ -77,6 +78,14 @@ if (DD4HEP_USE_GEANT4) REGEX_PASS "\\+\\+\\+ 5 Events out of 5 processed" REGEX_FAIL "Error;ERROR;Exception" ) + # Test DDDigi exception while processing + dd4hep_add_test_reg(DDDigi_test_processing_exception + COMMAND "${CMAKE_INSTALL_PREFIX}/bin/run_test_DDDigi.sh" + EXEC_ARGS ${Python_EXECUTABLE} ${CMAKE_INSTALL_PREFIX}/examples/DDDigi/scripts/TestInput.py + -num_events 1000 + DEPENDS DDDigi_generate_data + REGEX_PASS "\\+\\+ Terminate Digi and delete associated actions." + ) # Test signal attenuation for spillover dd4hep_add_test_reg(DDDigi_test_attenuate COMMAND "${CMAKE_INSTALL_PREFIX}/bin/run_test_DDDigi.sh" @@ -197,6 +206,15 @@ if (DD4HEP_USE_GEANT4) REGEX_PASS "\\+\\+\\+ 5 Events out of 5 processed" REGEX_FAIL "Error;ERROR;Exception" ) + # Test edm4hep write (needs to be expanded) + dd4hep_add_test_reg(DDDigi_test_edm4hep_write_1 + COMMAND "${CMAKE_INSTALL_PREFIX}/bin/run_test_DDDigi.sh" + EXEC_ARGS ${Python_EXECUTABLE} ${CMAKE_INSTALL_PREFIX}/examples/DDDigi/scripts/TestWriteEdm4hep.py + -num_events 5 -num_threads 10 -events_parallel 4 + DEPENDS DDDigi_generate_data + REGEX_PASS "\\+\\+\\+ 5 Events out of 5 processed." + REGEX_FAIL "Error;ERROR;Exception" + ) # endif() # diff --git a/examples/DDDigi/scripts/DigiTest.py b/examples/DDDigi/scripts/DigiTest.py index 2f144b30dc19a4d7e208fa85292f4ca361ff208b..34db86b5161fecfa195f88380e0c61a268c04738 100644 --- a/examples/DDDigi/scripts/DigiTest.py +++ b/examples/DDDigi/scripts/DigiTest.py @@ -18,9 +18,9 @@ from dddigi import DEBUG, INFO, WARNING, ERROR # noqa: F401 logging.basicConfig(format='%(levelname)s: %(message)s', level=logging.INFO) logger = logging.getLogger(__name__) -attenuation = {'MiniTel1Hits': 50 * units.ns, - 'MiniTel2Hits': 50 * units.ns, - 'MiniTel3Hits': 50 * units.ns, +attenuation = {'Minitel1Hits': 50 * units.ns, + 'Minitel2Hits': 50 * units.ns, + 'Minitel3Hits': 50 * units.ns, } @@ -111,12 +111,20 @@ class Test(dddigi.Digitize): def run_checked(self, num_events=5, num_threads=5, parallel=3): result = "FAILED" + if self.num_events: + num_events = int(self.num_events) + if self.num_threads: + num_threads = int(self.num_threads) + if self.events_parallel: + parallel = int(self.events_parallel) evt_done = self.run(num_events=num_events, num_threads=num_threads, parallel=parallel) if evt_done == num_events: result = "PASSED" - self.always('%s Test finished after processing %d events.' % (result, evt_done,)) + self.always('%s Test finished after processing %d events. [%d parallel threads, %d parallel events]' \ + % (result, evt_done, num_threads, parallel, )) self.always('Test done. Exiting') - + self.kernel().terminate() + return evt_done # ========================================================================================================== def test_setup_1(digi, print_level=WARNING, parallel=True): diff --git a/examples/DDDigi/scripts/TestSpillover.py b/examples/DDDigi/scripts/TestSpillover.py index f25109a7e0a61a095b21f47d3c2809a685f32930..6986b143313ce11a33f5268e5d4bae39b4f03054 100644 --- a/examples/DDDigi/scripts/TestSpillover.py +++ b/examples/DDDigi/scripts/TestSpillover.py @@ -15,12 +15,23 @@ from g4units import ns def run(): import DigiTest digi = DigiTest.Test(geometry=None) - attenuation = digi.attenuation + keep_raw = False rdr_output = DigiTest.DEBUG + erase_combined = True + + if digi.keep_raw: + keep_raw = True + if digi.output_level: + rdr_output = int(digi.output_level) + if digi.erase_combined: + erase_combined = int(erase_combined) + + attenuation = digi.attenuation input_action = digi.input_action('DigiParallelActionSequence/READER') # ======================================================================================================== input_action.adopt_action('DigiDDG4ROOT/SignalReader', mask=0x0, + keep_raw = keep_raw, input=[digi.next_input()], OutputLevel=rdr_output) digi.info('Created input.signal') @@ -30,6 +41,7 @@ def run(): spillover = input_action.adopt_action('DigiSequentialActionSequence/Spillover-25') evtreader = spillover.adopt_action('DigiDDG4ROOT/Reader-25ns', mask=0x1, + keep_raw = keep_raw, input=[digi.next_input()], OutputLevel=rdr_output) attenuate = spillover.adopt_action('DigiAttenuatorSequence/Att-25ns', @@ -46,6 +58,7 @@ def run(): spillover = input_action.adopt_action('DigiSequentialActionSequence/Spillover-50') evtreader = spillover.adopt_action('DigiDDG4ROOT/Reader-50ns', mask=0x2, + keep_raw = keep_raw, input=[digi.next_input()], OutputLevel=rdr_output) attenuate = spillover.adopt_action('DigiAttenuatorSequence/Att-50ns', @@ -60,6 +73,7 @@ def run(): spillover = input_action.adopt_action('DigiSequentialActionSequence/Spillover-75') evtreader = spillover.adopt_action('DigiDDG4ROOT/Reader-75ns', mask=0x3, + keep_raw = keep_raw, input=[digi.next_input()], OutputLevel=rdr_output) attenuate = spillover.adopt_action('DigiAttenuatorSequence/Att-75ns', @@ -76,6 +90,7 @@ def run(): spillover = input_action.adopt_action('DigiSequentialActionSequence/Spillover+25') evtreader = spillover.adopt_action('DigiDDG4ROOT/Reader+25ns', mask=0x4, + keep_raw = keep_raw, input=[digi.next_input()], OutputLevel=rdr_output) attenuate = spillover.adopt_action('DigiAttenuatorSequence/Att+25ns', @@ -90,6 +105,7 @@ def run(): spillover = input_action.adopt_action('DigiSequentialActionSequence/Spillover+50') evtreader = spillover.adopt_action('DigiDDG4ROOT/Reader+50ns', mask=0x5, + keep_raw = keep_raw, input=[digi.next_input()], OutputLevel=rdr_output) attenuate = spillover.adopt_action('DigiAttenuatorSequence/Att+50ns', @@ -104,6 +120,7 @@ def run(): spillover = input_action.adopt_action('DigiSequentialActionSequence/Spillover+75') evtreader = spillover.adopt_action('DigiDDG4ROOT/Reader+75ns', mask=0x6, + keep_raw = keep_raw, input=[digi.next_input()], OutputLevel=rdr_output) attenuate = spillover.adopt_action('DigiAttenuatorSequence/Att+75ns', @@ -121,7 +138,7 @@ def run(): input_segment='inputs', output_mask=0xFEED, output_segment='deposits') - combine.erase_combined = True # Not thread-safe! only do in SequentialActionSequence + combine.erase_combined = erase_combined evtdump = event.adopt_action('DigiStoreDump/StoreDump') digi.check_creation([combine, evtdump]) digi.info('Created event.dump') diff --git a/examples/DDDigi/scripts/TestWriteEdm4hep.py b/examples/DDDigi/scripts/TestWriteEdm4hep.py new file mode 100644 index 0000000000000000000000000000000000000000..596f7a14279690c4c98a8d2d8d857c838f9e4213 --- /dev/null +++ b/examples/DDDigi/scripts/TestWriteEdm4hep.py @@ -0,0 +1,38 @@ +# ========================================================================== +# AIDA Detector description implementation +# -------------------------------------------------------------------------- +# Copyright (C) Organisation europeenne pour la Recherche nucleaire (CERN) +# All rights reserved. +# +# For the licensing terms see $DD4hepINSTALL/LICENSE. +# For the list of contributors see $DD4hepINSTALL/doc/CREDITS. +# +# ========================================================================== +from __future__ import absolute_import + + +# --------------------------------------------------------------------------- +def run(): + import DigiTest + digi = DigiTest.Test(geometry=None) + read = digi.input_action('DigiDDG4ROOT/SignalReader', mask=0x0, input=[digi.next_input()]) + dump = digi.event_action('DigiStoreDump/StoreDump', parallel=False) + writ = digi.output_action('Digi2edm4hepAction/Writer', + parallel=True, + input_mask=0x0, + input_segment='input', + output='test_edm4hep.root') + proc = digi.create_action('Digi2edm4hepProcessor/edm4hep') + hit_type = 'TrackerHits' + if digi.hit_type: + hit_type = digi.hit_type + cont = [c + '/' + hit_type for c in digi.containers()] + writ.adopt_container_processor(proc, cont) + writ.adopt_container_processor(proc, 'MCParticles/MCParticles') + digi.check_creation([read, dump]) + digi.run_checked(num_events=10, num_threads=10, parallel=3) + + +# --------------------------------------------------------------------------- +if __name__ == '__main__': + run()