From ff5df8d93f0f90174e89ee79d6a4e3a94cc634ff Mon Sep 17 00:00:00 2001
From: Markus Frank <Markus.Frank@cern.ch>
Date: Mon, 12 Jun 2023 13:59:51 +0200
Subject: [PATCH] Move DDDigi edm4hep I/O to new frame writing. Improve DDG4
edm4hep I/O
---
.../include/DDDigi/DigiContainerProcessor.h | 16 +-
DDDigi/io/DigiEdm4hepInput.h | 27 +--
DDDigi/io/DigiEdm4hepOutput.cpp | 156 ++++++++++--------
DDDigi/plugins/DigiDepositDropKilled.cpp | 58 +++----
DDDigi/plugins/DigiDepositEnergyCut.cpp | 32 ++--
DDDigi/plugins/DigiDepositNoiseOnSignal.cpp | 36 ++--
DDDigi/plugins/DigiDepositRecalibEnergy.cpp | 44 ++---
DDDigi/plugins/DigiDepositSmearEnergy.cpp | 116 ++++++-------
.../DigiDepositSmearPositionResolution.cpp | 62 +++----
.../plugins/DigiDepositSmearPositionTrack.cpp | 96 +++++------
DDDigi/plugins/DigiDepositSmearTime.cpp | 48 +++---
.../plugins/DigiDepositWeightedPosition.cpp | 48 +++---
DDDigi/plugins/DigiDepositZeroSuppress.cpp | 36 ++--
DDDigi/plugins/DigiHitHistoryDrop.cpp | 94 +++++------
DDDigi/plugins/DigiIPCreate.cpp | 30 ++--
DDDigi/plugins/DigiSimpleADCResponse.cpp | 46 +++---
DDG4/edm4hep/Geant4Output2EDM4hep.cpp | 27 +--
17 files changed, 495 insertions(+), 477 deletions(-)
diff --git a/DDDigi/include/DDDigi/DigiContainerProcessor.h b/DDDigi/include/DDDigi/DigiContainerProcessor.h
index c5d3d6ff1..9875baa17 100644
--- a/DDDigi/include/DDDigi/DigiContainerProcessor.h
+++ b/DDDigi/include/DDDigi/DigiContainerProcessor.h
@@ -167,10 +167,18 @@ namespace dd4hep {
/// Main functional callback adapter
virtual void execute(context_t& context, work_t& work, const predicate_t& predicate) const;
};
-
-#define DEPOSIT_PROCESSOR_BIND_HANDLERS(X) { using namespace std::placeholders; \
- this->m_handleVector = std::bind( &X<DepositVector>, this, _1, _2, _3, _4); \
- this->m_handleMapping = std::bind( &X<DepositMapping>, this, _1, _2, _3, _4); }
+
+#define DEPOSIT_PROCESSOR_BIND_HANDLERS(X) \
+ this->m_handleVector = std::bind( &X<DepositVector>, this, \
+ std::placeholders::_1, \
+ std::placeholders::_2, \
+ std::placeholders::_3, \
+ std::placeholders::_4); \
+ this->m_handleMapping = std::bind( &X<DepositMapping>, this, \
+ std::placeholders::_1, \
+ std::placeholders::_2, \
+ std::placeholders::_3, \
+ std::placeholders::_4)
/// Worker class act on containers in an event identified by input masks and container name
/**
diff --git a/DDDigi/io/DigiEdm4hepInput.h b/DDDigi/io/DigiEdm4hepInput.h
index e80105cbf..ad1a8078e 100644
--- a/DDDigi/io/DigiEdm4hepInput.h
+++ b/DDDigi/io/DigiEdm4hepInput.h
@@ -74,28 +74,29 @@ namespace dd4hep {
/// Generic conversion function for hits
template <typename HIT_TYPE, typename EDM4HEP_COLLECTION_TYPE>
void hits_from_edm4hep(DigiContext& context,
- DataSegment& segment,
- Key::mask_type mask,
- const std::string& nam,
- const EDM4HEP_COLLECTION_TYPE* collection) const;
+ DataSegment& segment,
+ Key::mask_type mask,
+ const std::string& nam,
+ const EDM4HEP_COLLECTION_TYPE* collection) const;
/// Generic conversion function for MC particles
void parts_from_edm4hep(DigiContext& context,
- DataSegment& segment,
- int mask,
- const std::string& nam,
- const edm4hep::MCParticleCollection* collection) const;
+ DataSegment& segment,
+ int mask,
+ const std::string& nam,
+ const edm4hep::MCParticleCollection* collection) const;
/// Generic conversion function for event parameter settings
void params_from_edm4hep(DigiContext& context,
- DataSegment& segment,
- int mask,
- const std::string& nam,
- const podio::Frame& frame,
- const edm4hep::EventHeaderCollection* collection) const;
+ DataSegment& segment,
+ int mask,
+ const std::string& nam,
+ const podio::Frame& frame,
+ const edm4hep::EventHeaderCollection* collection) const;
/// Callback to handle single branch
virtual void operator()(DigiContext& context, work_t& work) const;
+
/// Event action callback
virtual void execute(DigiContext& context) const override;
};
diff --git a/DDDigi/io/DigiEdm4hepOutput.cpp b/DDDigi/io/DigiEdm4hepOutput.cpp
index 82d1653f5..17fc4a632 100644
--- a/DDDigi/io/DigiEdm4hepOutput.cpp
+++ b/DDDigi/io/DigiEdm4hepOutput.cpp
@@ -20,8 +20,9 @@
#include "DigiIO.h"
/// edm4hep include files
-#include <podio/EventStore.h>
-#include <podio/ROOTWriter.h>
+#include <podio/CollectionBase.h>
+#include <podio/ROOTFrameWriter.h>
+#include <podio/Frame.h>
#include <edm4hep/SimTrackerHit.h>
#include <edm4hep/MCParticleCollection.h>
#include <edm4hep/TrackerHitCollection.h>
@@ -45,33 +46,36 @@ namespace dd4hep {
*/
class DigiEdm4hepOutput::internals_t {
public:
- DigiEdm4hepOutput* m_parent { nullptr };
- /// Reference to podio store
- std::unique_ptr<podio::EventStore> m_store { };
+ using particlecollection_t = std::pair<std::string,std::unique_ptr<edm4hep::MCParticleCollection> >;
+ using headercollection_t = std::pair<std::string,std::unique_ptr<edm4hep::EventHeaderCollection> >;
+ DigiEdm4hepOutput* m_parent { nullptr };
/// Reference to podio writer
- std::unique_ptr<podio::ROOTWriter> m_file { };
+ std::unique_ptr<podio::ROOTFrameWriter> m_writer { };
+ /// Reference to podio store
+ podio::Frame m_frame { };
/// edm4hep event header collection
- edm4hep::EventHeaderCollection* m_header { nullptr };
+ headercollection_t m_header { };
/// MC particle collection
- edm4hep::MCParticleCollection* m_particles { nullptr };
- /// Collection of all edm4hep object collections
- std::map<std::string, podio::CollectionBase*> m_collections;
+ particlecollection_t m_particles { };
+ /// Collection of all edm4hep tracker object collections
+ std::map<std::string, std::unique_ptr<edm4hep::TrackerHitCollection> > m_tracker_collections;
+ /// Collection of all edm4hep calorimeter object collections
+ std::map<std::string, std::unique_ptr<edm4hep::CalorimeterHitCollection> > m_calo_collections;
+ std::string m_section_name { };
/// 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(DigiEdm4hepOutput* parent);
/// Default destructor
~internals_t();
+ /// Clear local data content
+ void clear();
/// Commit data at end of filling procedure
void commit();
/// Open new output stream
@@ -88,58 +92,80 @@ namespace dd4hep {
/// Default constructor
DigiEdm4hepOutput::internals_t::internals_t(DigiEdm4hepOutput* parent) : m_parent(parent)
{
- m_store = std::make_unique<podio::EventStore>();
}
/// Default destructor
DigiEdm4hepOutput::internals_t::~internals_t() {
- if ( m_file ) close();
- m_store.reset();
- }
-
- template <typename T> T* DigiEdm4hepOutput::internals_t::register_collection(const std::string& nam, T* coll) {
- m_collections.emplace(nam, coll);
- m_store->registerCollection(nam, coll);
- m_parent->debug("+++ created collection %s <%s>", nam.c_str(), coll->getTypeName().data());
- return coll;
+ if ( m_writer ) close();
+ m_tracker_collections.clear();
+ m_calo_collections.clear();
+ m_particles.second.reset();
+ m_header.second.reset();
}
/// Create all collections according to the parent setup
void DigiEdm4hepOutput::internals_t::create_collections() {
- if ( nullptr == m_header ) {
- m_header = register_collection("EventHeader", new edm4hep::EventHeaderCollection());
+ if ( !m_header.second ) {
+ m_header = std::make_pair("EventHeader", std::make_unique<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());
+ m_particles = std::make_pair(nam, std::make_unique<edm4hep::MCParticleCollection>());
}
else if ( typ == "TrackerHits" ) {
- register_collection(nam, new edm4hep::TrackerHitCollection());
+ m_tracker_collections.emplace(nam, std::make_unique<edm4hep::TrackerHitCollection>());
}
else if ( typ == "CalorimeterHits" ) {
- register_collection(nam, new edm4hep::CalorimeterHitCollection());
+ m_calo_collections.emplace(nam, std::make_unique<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* DigiEdm4hepOutput::internals_t::get_collection(const T& cont) {
- auto iter = m_collections.find(cont.name);
- if ( iter == m_collections.end() ) {
- m_parent->except("Error");
+ switch(cont.data_type) {
+ case SegmentEntry::TRACKER_HITS: {
+ auto iter = m_tracker_collections.find(cont.name);
+ if ( iter == m_tracker_collections.end() )
+ m_parent->except("Error");
+ return iter->second.get();
+ }
+ case SegmentEntry::CALORIMETER_HITS: {
+ auto iter = m_calo_collections.find(cont.name);
+ if ( iter == m_calo_collections.end() )
+ m_parent->except("Error");
+ return iter->second.get();
}
- return iter->second;
+ default:
+ return nullptr;
+ }
+ };
+
+ /// Clear local data content
+ void DigiEdm4hepOutput::internals_t::clear() {
+ m_header.second->clear();
+ m_particles.second->clear();
+ for( const auto& c : m_tracker_collections )
+ c.second->clear();
+ for( const auto& c : m_calo_collections )
+ c.second->clear();
}
/// Commit data at end of filling procedure
void DigiEdm4hepOutput::internals_t::commit() {
- if ( m_file ) {
- m_file->writeEvent();
- m_store->clearCollections();
+ if ( m_writer ) {
+ m_frame.put( std::move(*m_header.second), m_header.first);
+ m_frame.put( std::move(*m_particles.second), m_particles.first);
+ for( const auto& c : m_tracker_collections )
+ m_frame.put( std::move(*c.second), c.first);
+ for( const auto& c : m_calo_collections )
+ m_frame.put( std::move(*c.second), c.first);
+
+ m_writer->writeFrame(m_frame, m_section_name);
+ clear();
return;
}
m_parent->except("+++ Failed to write output file. [Stream is not open]");
@@ -147,25 +173,23 @@ namespace dd4hep {
/// Open new output stream
void DigiEdm4hepOutput::internals_t::open() {
- if ( m_file ) {
+ if ( m_writer ) {
close();
}
- m_file.reset();
+ clear();
+ m_writer.reset();
std::string fname = m_parent->next_stream_name();
- m_file = std::make_unique<podio::ROOTWriter>(fname, m_store.get());
+ m_writer = std::make_unique<podio::ROOTFrameWriter>(fname);
m_parent->info("+++ Opened EDM4HEP output file %s", fname.c_str());
- for( const auto& c : m_collections ) {
- m_file->registerForWrite(c.first);
- }
}
/// Commit data to disk and close output stream
void DigiEdm4hepOutput::internals_t::close() {
m_parent->info("+++ Closing EDM4HEP output file.");
- if ( m_file ) {
- m_file->finish();
+ if ( m_writer ) {
+ m_writer->finish();
}
- m_file.reset();
+ m_writer.reset();
}
/// Standard constructor
@@ -189,9 +213,9 @@ namespace dd4hep {
auto* act = dynamic_cast<DigiEdm4hepOutputProcessor*>(c.second);
if ( act ) { // This is not nice! Need to think about something better.
act->internals = this->internals;
- continue;
+ continue;
}
- except("Error: Invalid processor type for EDM4HEP output: %s", c.second->c_name());
+ except("Error: Invalid processor type for EDM4HEP output: %s", c.second->c_name());
}
m_parallel = false;
internals->create_collections();
@@ -199,7 +223,7 @@ namespace dd4hep {
/// Check for valid output stream
bool DigiEdm4hepOutput::have_output() const {
- return internals->m_file.get() != nullptr;
+ return internals->m_writer.get() != nullptr;
}
/// Open new output stream
@@ -227,24 +251,22 @@ namespace dd4hep {
}
void DigiEdm4hepOutputProcessor::convert_particles(DigiContext& ctxt,
- const ParticleMapping& cont) const
+ const ParticleMapping& cont) const
{
- auto* coll = internals->m_particles;
- std::size_t start = coll->size();
- data_io<edm4hep_input>::_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().data());
+ auto& parts = internals->m_particles.second;
+ data_io<edm4hep_input>::_to_edm4hep(cont, parts.get());
+ info("%s+++ %-24s added %6ld entries from mask: %04X to %s",
+ ctxt.event->id(), cont.name.c_str(), parts->size(), cont.key.mask(),
+ parts->getTypeName().data());
}
template <typename T> void
DigiEdm4hepOutputProcessor::convert_depos(const T& cont,
- const predicate_t& predicate,
- edm4hep::TrackerHitCollection* collection) const
+ const predicate_t& predicate,
+ edm4hep::TrackerHitCollection* collection) const
{
std::array<float,6> covMat = {0., 0., m_pointResoutionRPhi*m_pointResoutionRPhi,
- 0., 0., m_pointResoutionZ*m_pointResoutionZ
+ 0., 0., m_pointResoutionZ*m_pointResoutionZ
};
for ( const auto& depo : cont ) {
if ( predicate(depo) ) {
@@ -255,8 +277,8 @@ namespace dd4hep {
template <typename T> void
DigiEdm4hepOutputProcessor::convert_depos(const T& cont,
- const predicate_t& predicate,
- edm4hep::CalorimeterHitCollection* collection) const
+ const predicate_t& predicate,
+ edm4hep::CalorimeterHitCollection* collection) const
{
for ( const auto& depo : cont ) {
if ( predicate(depo) ) {
@@ -267,8 +289,8 @@ namespace dd4hep {
template <typename T> void
DigiEdm4hepOutputProcessor::convert_deposits(DigiContext& ctxt,
- const T& cont,
- const predicate_t& predicate) const
+ const T& cont,
+ const predicate_t& predicate) const
{
podio::CollectionBase* coll = internals->get_collection(cont);
std::size_t start = coll->size();
@@ -292,9 +314,9 @@ namespace dd4hep {
}
void DigiEdm4hepOutputProcessor::convert_history(DigiContext& ctxt,
- const DepositsHistory& cont,
- work_t& work,
- const predicate_t& predicate) const
+ 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),
diff --git a/DDDigi/plugins/DigiDepositDropKilled.cpp b/DDDigi/plugins/DigiDepositDropKilled.cpp
index e9c115ce7..7928cd57e 100644
--- a/DDDigi/plugins/DigiDepositDropKilled.cpp
+++ b/DDDigi/plugins/DigiDepositDropKilled.cpp
@@ -34,35 +34,35 @@ namespace dd4hep {
/// Main functional callback
virtual void execute(DigiContext& context, work_t& work, const predicate_t&) const override final {
- std::size_t killed = 0, total = 0, i = 0;
- if ( auto* v = work.get_input<DepositVector>() ) {
- total = v->size();
- for( auto iter = v->begin(); iter != v->end(); ++iter, ++i ) {
- if ( v->at(i).flag&EnergyDeposit::KILLED ) {
- v->remove(iter);
- iter = v->begin() + (--i);
- }
- }
- killed = total - v->size();
- }
- else if ( auto* m = work.get_input<DepositMapping>() ) {
- CellID last_cell = ~0x0LL;
- total = m->size();
- for( auto iter = m->begin(); iter != m->end(); ++iter ) {
- if ( iter->second.flag&EnergyDeposit::KILLED ) {
- m->remove(iter);
- iter = (last_cell != ~0x0LL) ? m->data.find(last_cell) : m->begin();
- continue;
- }
- last_cell = iter->first;
- }
- killed = total - m->size();
- }
- else {
- except("Request to handle unknown data type: %s", work.input_type_name().c_str());
- }
- if ( m_monitor ) m_monitor->count_shift(total, killed);
- info("%s+++ Killed %6ld out of %6ld deposit entries from container: %s",context.event->id(), killed, total);
+ std::size_t killed = 0, total = 0, i = 0;
+ if ( auto* v = work.get_input<DepositVector>() ) {
+ total = v->size();
+ for( auto iter = v->begin(); iter != v->end(); ++iter, ++i ) {
+ if ( v->at(i).flag&EnergyDeposit::KILLED ) {
+ v->remove(iter);
+ iter = v->begin() + (--i);
+ }
+ }
+ killed = total - v->size();
+ }
+ else if ( auto* m = work.get_input<DepositMapping>() ) {
+ CellID last_cell = ~0x0ULL;
+ total = m->size();
+ for( auto iter = m->begin(); iter != m->end(); ++iter ) {
+ if ( iter->second.flag&EnergyDeposit::KILLED ) {
+ m->remove(iter);
+ iter = (last_cell != ~0x0ULL) ? m->data.find(last_cell) : m->begin();
+ continue;
+ }
+ last_cell = iter->first;
+ }
+ killed = total - m->size();
+ }
+ else {
+ except("Request to handle unknown data type: %s", work.input_type_name().c_str());
+ }
+ if ( m_monitor ) m_monitor->count_shift(total, killed);
+ info("%s+++ Killed %6ld out of %6ld deposit entries from container: %s",context.event->id(), killed, total);
}
};
} // End namespace digi
diff --git a/DDDigi/plugins/DigiDepositEnergyCut.cpp b/DDDigi/plugins/DigiDepositEnergyCut.cpp
index c695f9015..91b676554 100644
--- a/DDDigi/plugins/DigiDepositEnergyCut.cpp
+++ b/DDDigi/plugins/DigiDepositEnergyCut.cpp
@@ -43,27 +43,27 @@ namespace dd4hep {
/// Create deposit mapping with updates on same cellIDs
template <typename T> void
cut_energy(context_t& context, T& cont, work_t& /* work */, const predicate_t& predicate) const {
- std::size_t dropped = 0UL;
- for( auto& dep : cont ) {
- if ( predicate(dep) ) {
- EnergyDeposit& depo = dep.second;
- if ( depo.deposit < m_cutoff ) {
- depo.flag |= EnergyDeposit::KILLED;
- ++dropped;
- }
- }
- }
- if ( m_monitor ) m_monitor->count_shift(cont.size(), dropped);
- info("%s+++ %-32s dropped %6ld out of %6ld entries from mask: %04X",
- context.event->id(), cont.name.c_str(), dropped, cont.size(), cont.key.mask());
+ std::size_t dropped = 0UL;
+ for( auto& dep : cont ) {
+ if ( predicate(dep) ) {
+ EnergyDeposit& depo = dep.second;
+ if ( depo.deposit < m_cutoff ) {
+ depo.flag |= EnergyDeposit::KILLED;
+ ++dropped;
+ }
+ }
+ }
+ if ( m_monitor ) m_monitor->count_shift(cont.size(), dropped);
+ info("%s+++ %-32s dropped %6ld out of %6ld entries from mask: %04X",
+ context.event->id(), cont.name.c_str(), dropped, cont.size(), cont.key.mask());
}
/// Standard constructor
DigiDepositEnergyCut(const DigiKernel& krnl, const std::string& nam)
- : DigiDepositsProcessor(krnl, nam)
+ : DigiDepositsProcessor(krnl, nam)
{
- declareProperty("deposit_cutoff", m_cutoff);
- DEPOSIT_PROCESSOR_BIND_HANDLERS(DigiDepositEnergyCut::cut_energy)
+ declareProperty("deposit_cutoff", m_cutoff);
+ DEPOSIT_PROCESSOR_BIND_HANDLERS(DigiDepositEnergyCut::cut_energy);
}
};
} // End namespace digi
diff --git a/DDDigi/plugins/DigiDepositNoiseOnSignal.cpp b/DDDigi/plugins/DigiDepositNoiseOnSignal.cpp
index be3f3e7a2..19143bc9b 100644
--- a/DDDigi/plugins/DigiDepositNoiseOnSignal.cpp
+++ b/DDDigi/plugins/DigiDepositNoiseOnSignal.cpp
@@ -42,29 +42,29 @@ namespace dd4hep {
/// Create deposit mapping with updates on same cellIDs
template <typename T> void
create_noise(DigiContext& context, T& cont, work_t& /* work */, const predicate_t& predicate) const {
- auto& random = context.randomGenerator();
- std::size_t updated = 0UL;
- for( auto& dep : cont ) {
- if ( predicate(dep) ) {
- int flag = EnergyDeposit::DEPOSIT_NOISE;
- double delta_E = random.gaussian(m_mean, m_sigma);
- if ( m_monitor ) m_monitor->energy_shift(dep, delta_E);
- dep.second.deposit += delta_E;
- dep.second.flag |= flag;
- ++updated;
- }
- }
- info("%s+++ %-32s Noise on signal: %6ld entries, updated %6ld entries. mask: %04X",
- context.event->id(), cont.name.c_str(), cont.size(), updated, cont.key.mask());
+ auto& random = context.randomGenerator();
+ std::size_t updated = 0UL;
+ for( auto& dep : cont ) {
+ if ( predicate(dep) ) {
+ int flag = EnergyDeposit::DEPOSIT_NOISE;
+ double delta_E = random.gaussian(m_mean, m_sigma);
+ if ( m_monitor ) m_monitor->energy_shift(dep, delta_E);
+ dep.second.deposit += delta_E;
+ dep.second.flag |= flag;
+ ++updated;
+ }
+ }
+ info("%s+++ %-32s Noise on signal: %6ld entries, updated %6ld entries. mask: %04X",
+ context.event->id(), cont.name.c_str(), cont.size(), updated, cont.key.mask());
}
/// Standard constructor
DigiDepositNoiseOnSignal(const DigiKernel& krnl, const std::string& nam)
- : DigiDepositsProcessor(krnl, nam)
+ : DigiDepositsProcessor(krnl, nam)
{
- declareProperty("mean", m_mean);
- declareProperty("sigma", m_sigma);
- DEPOSIT_PROCESSOR_BIND_HANDLERS(DigiDepositNoiseOnSignal::create_noise)
+ declareProperty("mean", m_mean);
+ declareProperty("sigma", m_sigma);
+ DEPOSIT_PROCESSOR_BIND_HANDLERS(DigiDepositNoiseOnSignal::create_noise);
}
};
} // End namespace digi
diff --git a/DDDigi/plugins/DigiDepositRecalibEnergy.cpp b/DDDigi/plugins/DigiDepositRecalibEnergy.cpp
index 8d61d2323..c3315a632 100644
--- a/DDDigi/plugins/DigiDepositRecalibEnergy.cpp
+++ b/DDDigi/plugins/DigiDepositRecalibEnergy.cpp
@@ -42,34 +42,34 @@ namespace dd4hep {
public:
/// Standard constructor
DigiDepositRecalibEnergy(const DigiKernel& krnl, const std::string& nam)
- : DigiDepositsProcessor(krnl, nam)
+ : DigiDepositsProcessor(krnl, nam)
{
- declareProperty("e0", m_e0);
- declareProperty("parameters" , m_polynomial);
- DEPOSIT_PROCESSOR_BIND_HANDLERS(DigiDepositRecalibEnergy::recalibrate_deposit_energy)
+ declareProperty("e0", m_e0);
+ declareProperty("parameters" , m_polynomial);
+ DEPOSIT_PROCESSOR_BIND_HANDLERS(DigiDepositRecalibEnergy::recalibrate_deposit_energy);
}
/// Create deposit mapping with updates on same cellIDs
template <typename T> void
recalibrate_deposit_energy(DigiContext& context, T& cont, work_t& /* work */, const predicate_t& predicate) const {
- std::size_t updated = 0UL;
- for( auto& dep : cont ) {
- if ( predicate(dep) ) {
- double deposit = dep.second.deposit;
- double delta_E = deposit - m_e0;
- double fac = 1e0;
- for( double param : m_polynomial ) {
- deposit += param * fac;
- fac *= delta_E;
- }
- if ( m_monitor ) m_monitor->energy_shift(dep, delta_E);
- dep.second.deposit = deposit;
- dep.second.flag |= EnergyDeposit::RECALIBRATED;
- ++updated;
- }
- }
- info("%s+++ %-32s Recalibrating: updated %6ld out of %6ld entries from mask: %04X",
- context.event->id(), cont.name.c_str(), updated, cont.size(), cont.key.mask());
+ std::size_t updated = 0UL;
+ for( auto& dep : cont ) {
+ if ( predicate(dep) ) {
+ double deposit = dep.second.deposit;
+ double delta_E = deposit - m_e0;
+ double fac = 1e0;
+ for( double param : m_polynomial ) {
+ deposit += param * fac;
+ fac *= delta_E;
+ }
+ if ( m_monitor ) m_monitor->energy_shift(dep, delta_E);
+ dep.second.deposit = deposit;
+ dep.second.flag |= EnergyDeposit::RECALIBRATED;
+ ++updated;
+ }
+ }
+ info("%s+++ %-32s Recalibrating: updated %6ld out of %6ld entries from mask: %04X",
+ context.event->id(), cont.name.c_str(), updated, cont.size(), cont.key.mask());
}
};
} // End namespace digi
diff --git a/DDDigi/plugins/DigiDepositSmearEnergy.cpp b/DDDigi/plugins/DigiDepositSmearEnergy.cpp
index 66a781707..299f61334 100644
--- a/DDDigi/plugins/DigiDepositSmearEnergy.cpp
+++ b/DDDigi/plugins/DigiDepositSmearEnergy.cpp
@@ -87,69 +87,69 @@ namespace dd4hep {
public:
/// Standard constructor
DigiDepositSmearEnergy(const DigiKernel& krnl, const std::string& nam)
- : DigiDepositsProcessor(krnl, nam)
+ : DigiDepositsProcessor(krnl, nam)
{
- declareProperty("intrinsic_fluctuation", m_intrinsic_fluctuation = 0e0);
- declareProperty("instrumentation_resolution", m_instrumentation_resolution = 0e0);
- declareProperty("systematic_resolution", m_systematic_resolution = 0e0);
- declareProperty("pair_ionisation_energy", m_pair_ionization_energy = 3.6*dd4hep::eV);
- declareProperty("ionization_fluctuation", m_ionization_fluctuation = false);
- declareProperty("modify_energy", m_modify_energy = true);
- DEPOSIT_PROCESSOR_BIND_HANDLERS(DigiDepositSmearEnergy::smear)
+ declareProperty("intrinsic_fluctuation", m_intrinsic_fluctuation = 0e0);
+ declareProperty("instrumentation_resolution", m_instrumentation_resolution = 0e0);
+ declareProperty("systematic_resolution", m_systematic_resolution = 0e0);
+ declareProperty("pair_ionisation_energy", m_pair_ionization_energy = 3.6*dd4hep::eV);
+ declareProperty("ionization_fluctuation", m_ionization_fluctuation = false);
+ declareProperty("modify_energy", m_modify_energy = true);
+ DEPOSIT_PROCESSOR_BIND_HANDLERS(DigiDepositSmearEnergy::smear);
}
/// Create deposit mapping with updates on same cellIDs
template <typename T> void
smear(DigiContext& context, T& cont, work_t& /* work */, const predicate_t& predicate) const {
- auto& random = context.randomGenerator();
- std::size_t updated = 0UL;
+ auto& random = context.randomGenerator();
+ std::size_t updated = 0UL;
- for( auto& dep : cont ) {
- if ( predicate(dep) ) {
- CellID cell = dep.first;
- EnergyDeposit& depo = dep.second;
- double deposit = depo.deposit;
- double delta_E = 0e0;
- double energy = deposit / dd4hep::GeV; // E in units of GeV
- double sigma_E_systematic = m_systematic_resolution * energy;
- double sigma_E_intrin_fluct = m_intrinsic_fluctuation * std::sqrt(energy);
- double sigma_E_instrument = m_instrumentation_resolution / dd4hep::GeV;
- double delta_ion = 0e0, num_pairs = 0e0;
- constexpr static double eps = std::numeric_limits<double>::epsilon();
- if ( sigma_E_systematic > eps ) {
- delta_E += sigma_E_systematic * random.gaussian(0e0, 1e0);
- }
- if ( sigma_E_intrin_fluct > eps ) {
- delta_E += sigma_E_intrin_fluct * random.gaussian(0e0, 1e0);
- }
- if ( sigma_E_instrument > eps ) {
- delta_E += sigma_E_instrument * random.gaussian(0e0, 1e0);
- }
- if ( m_ionization_fluctuation ) {
- num_pairs = energy / (m_pair_ionization_energy/dd4hep::GeV);
- delta_ion = energy * (random.poisson(num_pairs)/num_pairs);
- delta_E += delta_ion;
- }
- if ( dd4hep::isActivePrintLevel(outputLevel()) ) {
- print("%s+++ %016lX [GeV] E:%9.2e [%9.2e %9.2e] intrin_fluct:%9.2e systematic:%9.2e instrument:%9.2e ioni:%9.2e/%.0f",
- context.event->id(), cell, energy, deposit/dd4hep::GeV, delta_E,
- sigma_E_intrin_fluct, sigma_E_systematic, sigma_E_instrument, delta_ion, num_pairs);
- }
- /// delta_E is in GeV
- delta_E *= dd4hep::GeV;
- depo.depositError = delta_E;
- if ( m_monitor ) {
- m_monitor->energy_shift(dep, delta_E);
- }
- if ( m_modify_energy ) {
- depo.deposit = deposit + delta_E;
- depo.flag |= EnergyDeposit::ENERGY_SMEARED;
- }
- ++updated;
- }
- }
- info("%s+++ %-32s Smear energy: updated %6ld out of %6ld entries from mask: %04X",
- context.event->id(), cont.name.c_str(), updated, cont.size(), cont.key.mask());
+ for( auto& dep : cont ) {
+ if ( predicate(dep) ) {
+ CellID cell = dep.first;
+ EnergyDeposit& depo = dep.second;
+ double deposit = depo.deposit;
+ double delta_E = 0e0;
+ double energy = deposit / dd4hep::GeV; // E in units of GeV
+ double sigma_E_systematic = m_systematic_resolution * energy;
+ double sigma_E_intrin_fluct = m_intrinsic_fluctuation * std::sqrt(energy);
+ double sigma_E_instrument = m_instrumentation_resolution / dd4hep::GeV;
+ double delta_ion = 0e0, num_pairs = 0e0;
+ constexpr static double eps = std::numeric_limits<double>::epsilon();
+ if ( sigma_E_systematic > eps ) {
+ delta_E += sigma_E_systematic * random.gaussian(0e0, 1e0);
+ }
+ if ( sigma_E_intrin_fluct > eps ) {
+ delta_E += sigma_E_intrin_fluct * random.gaussian(0e0, 1e0);
+ }
+ if ( sigma_E_instrument > eps ) {
+ delta_E += sigma_E_instrument * random.gaussian(0e0, 1e0);
+ }
+ if ( m_ionization_fluctuation ) {
+ num_pairs = energy / (m_pair_ionization_energy/dd4hep::GeV);
+ delta_ion = energy * (random.poisson(num_pairs)/num_pairs);
+ delta_E += delta_ion;
+ }
+ if ( dd4hep::isActivePrintLevel(outputLevel()) ) {
+ print("%s+++ %016lX [GeV] E:%9.2e [%9.2e %9.2e] intrin_fluct:%9.2e systematic:%9.2e instrument:%9.2e ioni:%9.2e/%.0f",
+ context.event->id(), cell, energy, deposit/dd4hep::GeV, delta_E,
+ sigma_E_intrin_fluct, sigma_E_systematic, sigma_E_instrument, delta_ion, num_pairs);
+ }
+ /// delta_E is in GeV
+ delta_E *= dd4hep::GeV;
+ depo.depositError = delta_E;
+ if ( m_monitor ) {
+ m_monitor->energy_shift(dep, delta_E);
+ }
+ if ( m_modify_energy ) {
+ depo.deposit = deposit + delta_E;
+ depo.flag |= EnergyDeposit::ENERGY_SMEARED;
+ }
+ ++updated;
+ }
+ }
+ info("%s+++ %-32s Smear energy: updated %6ld out of %6ld entries from mask: %04X",
+ context.event->id(), cont.name.c_str(), updated, cont.size(), cont.key.mask());
}
};
@@ -164,9 +164,9 @@ namespace dd4hep {
public:
/// Standard constructor
DigiDepositSetEnergyError(const DigiKernel& krnl, const std::string& nam)
- : DigiDepositSmearEnergy(krnl, nam)
+ : DigiDepositSmearEnergy(krnl, nam)
{
- m_modify_energy = false;
+ m_modify_energy = false;
}
};
} // End namespace digi
diff --git a/DDDigi/plugins/DigiDepositSmearPositionResolution.cpp b/DDDigi/plugins/DigiDepositSmearPositionResolution.cpp
index 3946c5304..dc689da76 100644
--- a/DDDigi/plugins/DigiDepositSmearPositionResolution.cpp
+++ b/DDDigi/plugins/DigiDepositSmearPositionResolution.cpp
@@ -48,47 +48,47 @@ namespace dd4hep {
public:
/// Standard constructor
DigiDepositSmearPositionResolution(const DigiKernel& krnl, const std::string& nam)
- : DigiDepositsProcessor(krnl, nam)
+ : DigiDepositsProcessor(krnl, nam)
{
- declareProperty("resolution_u", m_resolution_u);
- declareProperty("resolution_v", m_resolution_v);
- DEPOSIT_PROCESSOR_BIND_HANDLERS(DigiDepositSmearPositionResolution::smear)
+ declareProperty("resolution_u", m_resolution_u);
+ declareProperty("resolution_v", m_resolution_v);
+ DEPOSIT_PROCESSOR_BIND_HANDLERS(DigiDepositSmearPositionResolution::smear);
}
/// Create deposit mapping with updates on same cellIDs
template <typename T> void
smear(DigiContext& context, T& cont, work_t& /* work */, const predicate_t& predicate) const {
- VolumeManager volMgr = m_kernel.detectorDescription().volumeManager();
- auto& random = context.randomGenerator();
- std::size_t updated = 0UL;
+ VolumeManager volMgr = m_kernel.detectorDescription().volumeManager();
+ auto& random = context.randomGenerator();
+ std::size_t updated = 0UL;
- for( auto& dep : cont ) {
- if ( predicate(dep) ) {
- CellID cell = dep.first;
- EnergyDeposit& depo = dep.second;
- auto* ctxt = volMgr.lookupContext(cell);
- Position local_pos = ctxt->worldToLocal(depo.position);
- double delta_u = m_resolution_u * random.gaussian();
- double delta_v = m_resolution_v * random.gaussian();
- Position delta_pos(delta_u, delta_v, 0e0);
- Position oldpos = depo.position;
- Position newpos = ctxt->localToWorld(local_pos + delta_pos);
+ for( auto& dep : cont ) {
+ if ( predicate(dep) ) {
+ CellID cell = dep.first;
+ EnergyDeposit& depo = dep.second;
+ auto* ctxt = volMgr.lookupContext(cell);
+ Position local_pos = ctxt->worldToLocal(depo.position);
+ double delta_u = m_resolution_u * random.gaussian();
+ double delta_v = m_resolution_v * random.gaussian();
+ Position delta_pos(delta_u, delta_v, 0e0);
+ Position oldpos = depo.position;
+ Position newpos = ctxt->localToWorld(local_pos + delta_pos);
- if ( m_local_monitor ) m_local_monitor->position_shift(dep, local_pos, delta_pos);
- delta_pos = newpos - oldpos;
- if ( m_monitor ) m_monitor->position_shift(dep, oldpos, delta_pos);
+ if ( m_local_monitor ) m_local_monitor->position_shift(dep, local_pos, delta_pos);
+ delta_pos = newpos - oldpos;
+ if ( m_monitor ) m_monitor->position_shift(dep, oldpos, delta_pos);
- info("+++ %016lX Smeared position [%9.2e %9.2e %9.2e] by [%9.2e %9.2e %9.2e] to [%9.2e %9.2e %9.2e] [mm]",
- oldpos.X(), oldpos.Y(), oldpos.Z(), delta_pos.X(), delta_pos.Y(), delta_pos.Z(),
- newpos.X(), newpos.Y(), newpos.Z());
+ info("+++ %016lX Smeared position [%9.2e %9.2e %9.2e] by [%9.2e %9.2e %9.2e] to [%9.2e %9.2e %9.2e] [mm]",
+ oldpos.X(), oldpos.Y(), oldpos.Z(), delta_pos.X(), delta_pos.Y(), delta_pos.Z(),
+ newpos.X(), newpos.Y(), newpos.Z());
- depo.position = newpos;
- depo.flag |= EnergyDeposit::POSITION_SMEARED;
- ++updated;
- }
- }
- info("%s+++ %-32s Smear position(resolution): updated %6ld out of %6ld entries from mask: %04X",
- context.event->id(), cont.name.c_str(), updated, cont.size(), cont.key.mask());
+ depo.position = newpos;
+ depo.flag |= EnergyDeposit::POSITION_SMEARED;
+ ++updated;
+ }
+ }
+ info("%s+++ %-32s Smear position(resolution): updated %6ld out of %6ld entries from mask: %04X",
+ context.event->id(), cont.name.c_str(), updated, cont.size(), cont.key.mask());
}
};
} // End namespace digi
diff --git a/DDDigi/plugins/DigiDepositSmearPositionTrack.cpp b/DDDigi/plugins/DigiDepositSmearPositionTrack.cpp
index 3c321809a..a05b49d9b 100644
--- a/DDDigi/plugins/DigiDepositSmearPositionTrack.cpp
+++ b/DDDigi/plugins/DigiDepositSmearPositionTrack.cpp
@@ -49,12 +49,12 @@ namespace dd4hep {
public:
/// Standard constructor
DigiDepositSmearPositionTrack(const DigiKernel& krnl, const std::string& nam)
- : DigiDepositsProcessor(krnl, nam)
+ : DigiDepositsProcessor(krnl, nam)
{
- declareProperty("resolution_u", m_resolution_u);
- declareProperty("resolution_v", m_resolution_v);
- m_kernel.register_initialize(std::bind(&DigiDepositSmearPositionTrack::initialize,this));
- DEPOSIT_PROCESSOR_BIND_HANDLERS(DigiDepositSmearPositionTrack::smear)
+ declareProperty("resolution_u", m_resolution_u);
+ declareProperty("resolution_v", m_resolution_v);
+ m_kernel.register_initialize(std::bind(&DigiDepositSmearPositionTrack::initialize,this));
+ DEPOSIT_PROCESSOR_BIND_HANDLERS(DigiDepositSmearPositionTrack::smear);
}
/// Processor initialization
@@ -64,49 +64,49 @@ namespace dd4hep {
/// Create deposit mapping with updates on same cellIDs
template <typename T> void
smear(DigiContext& context, T& cont, work_t& /* work */, const predicate_t& predicate) const {
- constexpr double eps = detail::numeric_epsilon;
- auto& random = context.randomGenerator();
- const auto& ev = *(context.event);
- std::size_t updated = 0UL;
-
- VolumeManager volMgr = m_kernel.detectorDescription().volumeManager();
- for( auto& dep : cont ) {
- if ( predicate(dep) ) {
- CellID cell = dep.first;
- EnergyDeposit& depo = dep.second;
- auto* ctxt = volMgr.lookupContext(cell);
- Direction part_momentum = depo.history.average_particle_momentum(ev);
- Position local_pos = ctxt->worldToLocal(depo.position);
- Position local_dir = ctxt->worldToLocal(part_momentum).unit();
- double cos_u = local_dir.Dot(Position(1,0,0));
- double sin_u = std::sqrt(1e0 - cos_u*cos_u);
- double tan_u = sin_u/(std::abs(cos_u)>eps ? cos_u : eps);
- double delta_u = tan_u * m_resolution_u * random.gaussian();
-
- double cos_v = local_dir.Dot(Position(0,1,0));
- double sin_v = std::sqrt(1e0 - cos_v*cos_v);
- double tan_v = sin_v/(std::abs(cos_v)>eps ? cos_v : eps);
- double delta_v = tan_v * m_resolution_v * random.gaussian();
-
- Position delta_pos(delta_u, delta_v, 0e0);
- Position oldpos = depo.position;
- Position newpos = ctxt->localToWorld(local_pos + delta_pos);
-
- if ( m_local_monitor ) m_local_monitor->position_shift(dep, local_pos, delta_pos);
- delta_pos = newpos - oldpos;
- if ( m_monitor ) m_monitor->position_shift(dep, oldpos, delta_pos);
-
- info("%s+++ %016lX Smeared position [%9.2e %9.2e %9.2e] by [%9.2e %9.2e %9.2e] to [%9.2e %9.2e %9.2e] [mm]",
- ev.id(), oldpos.X(), oldpos.Y(), oldpos.Z(), delta_pos.X(), delta_pos.Y(), delta_pos.Z(),
- newpos.X(), newpos.Y(), newpos.Z());
-
- depo.position = newpos;
- depo.flag |= EnergyDeposit::POSITION_SMEARED;
- ++updated;
- }
- }
- info("%s+++ %-32s Smear position(track): updated %6ld out of %6ld entries from mask: %04X",
- context.event->id(), cont.name.c_str(), updated, cont.size(), cont.key.mask());
+ constexpr double eps = detail::numeric_epsilon;
+ auto& random = context.randomGenerator();
+ const auto& ev = *(context.event);
+ std::size_t updated = 0UL;
+
+ VolumeManager volMgr = m_kernel.detectorDescription().volumeManager();
+ for( auto& dep : cont ) {
+ if ( predicate(dep) ) {
+ CellID cell = dep.first;
+ EnergyDeposit& depo = dep.second;
+ auto* ctxt = volMgr.lookupContext(cell);
+ Direction part_momentum = depo.history.average_particle_momentum(ev);
+ Position local_pos = ctxt->worldToLocal(depo.position);
+ Position local_dir = ctxt->worldToLocal(part_momentum).unit();
+ double cos_u = local_dir.Dot(Position(1,0,0));
+ double sin_u = std::sqrt(1e0 - cos_u*cos_u);
+ double tan_u = sin_u/(std::abs(cos_u)>eps ? cos_u : eps);
+ double delta_u = tan_u * m_resolution_u * random.gaussian();
+
+ double cos_v = local_dir.Dot(Position(0,1,0));
+ double sin_v = std::sqrt(1e0 - cos_v*cos_v);
+ double tan_v = sin_v/(std::abs(cos_v)>eps ? cos_v : eps);
+ double delta_v = tan_v * m_resolution_v * random.gaussian();
+
+ Position delta_pos(delta_u, delta_v, 0e0);
+ Position oldpos = depo.position;
+ Position newpos = ctxt->localToWorld(local_pos + delta_pos);
+
+ if ( m_local_monitor ) m_local_monitor->position_shift(dep, local_pos, delta_pos);
+ delta_pos = newpos - oldpos;
+ if ( m_monitor ) m_monitor->position_shift(dep, oldpos, delta_pos);
+
+ info("%s+++ %016lX Smeared position [%9.2e %9.2e %9.2e] by [%9.2e %9.2e %9.2e] to [%9.2e %9.2e %9.2e] [mm]",
+ ev.id(), oldpos.X(), oldpos.Y(), oldpos.Z(), delta_pos.X(), delta_pos.Y(), delta_pos.Z(),
+ newpos.X(), newpos.Y(), newpos.Z());
+
+ depo.position = newpos;
+ depo.flag |= EnergyDeposit::POSITION_SMEARED;
+ ++updated;
+ }
+ }
+ info("%s+++ %-32s Smear position(track): updated %6ld out of %6ld entries from mask: %04X",
+ context.event->id(), cont.name.c_str(), updated, cont.size(), cont.key.mask());
}
};
} // End namespace digi
diff --git a/DDDigi/plugins/DigiDepositSmearTime.cpp b/DDDigi/plugins/DigiDepositSmearTime.cpp
index ff7821296..a6b619b00 100644
--- a/DDDigi/plugins/DigiDepositSmearTime.cpp
+++ b/DDDigi/plugins/DigiDepositSmearTime.cpp
@@ -44,35 +44,35 @@ namespace dd4hep {
/// Create deposit mapping with updates on same cellIDs
template <typename T> void
smear(DigiContext& context, T& cont, work_t& /* work */, const predicate_t& predicate) const {
- auto& random = context.randomGenerator();
- std::size_t killed = 0UL;
- std::size_t updated = 0UL;
- for( auto& dep : cont ) {
- if ( predicate(dep) ) {
- int flag = EnergyDeposit::TIME_SMEARED;
- double delta_T = m_resolution_time * random.gaussian();
- if ( delta_T < m_window_time.first || delta_T > m_window_time.second ) {
- flag |= EnergyDeposit::KILLED;
- ++killed;
- }
- if ( m_monitor ) m_monitor->time_shift(dep, delta_T);
- dep.second.time += delta_T;
- dep.second.flag |= flag;
- ++updated;
- }
- }
- if ( m_monitor ) m_monitor->count_shift(cont.size(), -killed);
- info("%s+++ %-32s Smeared time resolution: %6ld entries, updated %6ld killed %6ld entries from mask: %04X",
- context.event->id(), cont.name.c_str(), cont.size(), updated, killed, cont.key.mask());
+ auto& random = context.randomGenerator();
+ std::size_t killed = 0UL;
+ std::size_t updated = 0UL;
+ for( auto& dep : cont ) {
+ if ( predicate(dep) ) {
+ int flag = EnergyDeposit::TIME_SMEARED;
+ double delta_T = m_resolution_time * random.gaussian();
+ if ( delta_T < m_window_time.first || delta_T > m_window_time.second ) {
+ flag |= EnergyDeposit::KILLED;
+ ++killed;
+ }
+ if ( m_monitor ) m_monitor->time_shift(dep, delta_T);
+ dep.second.time += delta_T;
+ dep.second.flag |= flag;
+ ++updated;
+ }
+ }
+ if ( m_monitor ) m_monitor->count_shift(cont.size(), -killed);
+ info("%s+++ %-32s Smeared time resolution: %6ld entries, updated %6ld killed %6ld entries from mask: %04X",
+ context.event->id(), cont.name.c_str(), cont.size(), updated, killed, cont.key.mask());
}
/// Standard constructor
DigiDepositSmearTime(const DigiKernel& krnl, const std::string& nam)
- : DigiDepositsProcessor(krnl, nam)
+ : DigiDepositsProcessor(krnl, nam)
{
- declareProperty("resolution_time", m_resolution_time);
- declareProperty("window_time", m_window_time);
- DEPOSIT_PROCESSOR_BIND_HANDLERS(DigiDepositSmearTime::smear)
+ declareProperty("resolution_time", m_resolution_time);
+ declareProperty("window_time", m_window_time);
+ DEPOSIT_PROCESSOR_BIND_HANDLERS(DigiDepositSmearTime::smear);
}
};
} // End namespace digi
diff --git a/DDDigi/plugins/DigiDepositWeightedPosition.cpp b/DDDigi/plugins/DigiDepositWeightedPosition.cpp
index a0cc68750..ffa2ec940 100644
--- a/DDDigi/plugins/DigiDepositWeightedPosition.cpp
+++ b/DDDigi/plugins/DigiDepositWeightedPosition.cpp
@@ -43,35 +43,35 @@ namespace dd4hep {
/// Create deposit mapping with updates on same cellIDs
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, cont.data_type );
- std::size_t dropped = 0UL, updated = 0UL, added = 0UL;
- for( const auto& dep : cont ) {
- if ( predicate(dep) ) {
- const EnergyDeposit& depo = dep.second;
- if ( depo.deposit >= m_cutoff ) {
- CellID cell = dep.first;
- auto iter = m.data.find(cell);
- if ( iter == m.data.end() )
- m.data.emplace(dep.first, depo), ++added;
- else
- iter->second.update_deposit_weighted(depo), ++updated;
- continue;
- }
- ++dropped;
- }
- }
- info("%s+++ %-32s added %6ld updated %6ld dropped %6ld entries (now: %6ld) from mask: %04X to mask: %04X",
- context.event->id(), cont.name.c_str(), added, updated, dropped, m.size(), cont.key.mask(), m.key.mask());
- work.environ.output.data.put(m.key, std::move(m));
+ Key key(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) ) {
+ const EnergyDeposit& depo = dep.second;
+ if ( depo.deposit >= m_cutoff ) {
+ CellID cell = dep.first;
+ auto iter = m.data.find(cell);
+ if ( iter == m.data.end() )
+ m.data.emplace(dep.first, depo), ++added;
+ else
+ iter->second.update_deposit_weighted(depo), ++updated;
+ continue;
+ }
+ ++dropped;
+ }
+ }
+ info("%s+++ %-32s added %6ld updated %6ld dropped %6ld entries (now: %6ld) from mask: %04X to mask: %04X",
+ context.event->id(), cont.name.c_str(), added, updated, dropped, m.size(), cont.key.mask(), m.key.mask());
+ work.environ.output.data.put(m.key, std::move(m));
}
/// Standard constructor
DigiDepositWeightedPosition(const DigiKernel& krnl, const std::string& nam)
- : DigiDepositsProcessor(krnl, nam)
+ : DigiDepositsProcessor(krnl, nam)
{
- declareProperty("deposit_cutoff", m_cutoff);
- DEPOSIT_PROCESSOR_BIND_HANDLERS(DigiDepositWeightedPosition::create_deposits)
+ declareProperty("deposit_cutoff", m_cutoff);
+ DEPOSIT_PROCESSOR_BIND_HANDLERS(DigiDepositWeightedPosition::create_deposits);
}
};
} // End namespace digi
diff --git a/DDDigi/plugins/DigiDepositZeroSuppress.cpp b/DDDigi/plugins/DigiDepositZeroSuppress.cpp
index 4f245f43a..b8208ebc7 100644
--- a/DDDigi/plugins/DigiDepositZeroSuppress.cpp
+++ b/DDDigi/plugins/DigiDepositZeroSuppress.cpp
@@ -41,28 +41,28 @@ namespace dd4hep {
/// Create deposit mapping with updates on same cellIDs
template <typename T> void
handle_deposits(DigiContext& context, T& cont, work_t& /* work */, const predicate_t& predicate) const {
- std::size_t killed = 0UL;
- std::size_t handled = 0UL;
- for( auto& dep : cont ) {
- if ( predicate(dep) ) {
- int flag = EnergyDeposit::ZERO_SUPPRESSED;
- if ( dep.second.deposit * dd4hep::GeV < m_energy_threshold ) {
- flag |= EnergyDeposit::KILLED;
- ++killed;
- }
- dep.second.flag |= flag;
- ++handled;
- }
- }
- info("%s+++ %-32s Zero suppression: entries: %6ld handled: %6ld killed %6ld entries from mask: %04X",
- context.event->id(), cont.name.c_str(), cont.size(), handled, killed, cont.key.mask());
+ std::size_t killed = 0UL;
+ std::size_t handled = 0UL;
+ for( auto& dep : cont ) {
+ if ( predicate(dep) ) {
+ int flag = EnergyDeposit::ZERO_SUPPRESSED;
+ if ( dep.second.deposit * dd4hep::GeV < m_energy_threshold ) {
+ flag |= EnergyDeposit::KILLED;
+ ++killed;
+ }
+ dep.second.flag |= flag;
+ ++handled;
+ }
+ }
+ info("%s+++ %-32s Zero suppression: entries: %6ld handled: %6ld killed %6ld entries from mask: %04X",
+ context.event->id(), cont.name.c_str(), cont.size(), handled, killed, cont.key.mask());
}
/// Standard constructor
DigiDepositZeroSuppress(const DigiKernel& krnl, const std::string& nam)
- : DigiDepositsProcessor(krnl, nam)
+ : DigiDepositsProcessor(krnl, nam)
{
- declareProperty("threshold", m_energy_threshold);
- DEPOSIT_PROCESSOR_BIND_HANDLERS(DigiDepositZeroSuppress::handle_deposits)
+ declareProperty("threshold", m_energy_threshold);
+ DEPOSIT_PROCESSOR_BIND_HANDLERS(DigiDepositZeroSuppress::handle_deposits);
}
};
} // End namespace digi
diff --git a/DDDigi/plugins/DigiHitHistoryDrop.cpp b/DDDigi/plugins/DigiHitHistoryDrop.cpp
index eca3c87f7..445d2639a 100644
--- a/DDDigi/plugins/DigiHitHistoryDrop.cpp
+++ b/DDDigi/plugins/DigiHitHistoryDrop.cpp
@@ -46,69 +46,69 @@ namespace dd4hep {
/// Default destructor
virtual ~DigiHitHistoryDrop() {
- InstanceCount::decrement(this);
+ InstanceCount::decrement(this);
}
public:
/// Standard constructor
DigiHitHistoryDrop(const DigiKernel& krnl, const std::string& nam)
- : DigiEventAction(krnl, nam)
+ : DigiEventAction(krnl, nam)
{
- declareProperty("masks", m_masks);
- declareProperty("input", m_input = "inputs");
- InstanceCount::increment(this);
+ declareProperty("masks", m_masks);
+ declareProperty("input", m_input = "inputs");
+ InstanceCount::increment(this);
}
template <typename T>
std::pair<std::size_t, std::size_t> drop_history(T& cont) const {
- std::size_t num_drop_hit = 0;
- std::size_t num_drop_particle = 0;
- for( auto& dep : cont ) {
- auto ret = dep.second.history.drop();
- num_drop_hit += ret.first;
- num_drop_particle += ret.second;
- }
- return std::make_pair(num_drop_hit,num_drop_particle);
+ std::size_t num_drop_hit = 0;
+ std::size_t num_drop_particle = 0;
+ for( auto& dep : cont ) {
+ auto ret = dep.second.history.drop();
+ num_drop_hit += ret.first;
+ num_drop_particle += ret.second;
+ }
+ return std::make_pair(num_drop_hit,num_drop_particle);
}
std::pair<std::size_t, std::size_t> drop_history(DetectorHistory& cont) const {
- std::size_t num_drop_hit = 0;
- std::size_t num_drop_particle = 0;
- for( auto& dep : cont ) {
- auto ret = dep.second.drop();
- num_drop_hit += ret.first;
- num_drop_particle += ret.second;
- }
- return std::make_pair(num_drop_hit,num_drop_particle);
+ std::size_t num_drop_hit = 0;
+ std::size_t num_drop_particle = 0;
+ for( auto& dep : cont ) {
+ auto ret = dep.second.drop();
+ num_drop_hit += ret.first;
+ num_drop_particle += ret.second;
+ }
+ return std::make_pair(num_drop_hit,num_drop_particle);
}
/// Main functional callback
virtual void execute(DigiContext& context) const final {
- auto& inputs = context.event->get_segment(m_input);
- std::size_t num_drop_hit = 0;
- std::size_t num_drop_particle = 0;
- for( auto& i : inputs ) {
- Key key(i.first);
- auto im = std::find(m_masks.begin(), m_masks.end(), key.mask());
- if( im != m_masks.end() ) {
- if ( DepositMapping* m = std::any_cast<DepositMapping>(&i.second) ) {
- auto ret = drop_history(*m);
- num_drop_hit += ret.first;
- num_drop_particle += ret.second;
- }
- else if ( DepositVector* v = std::any_cast<DepositVector>(&i.second) ) {
- auto ret = drop_history(*v);
- num_drop_hit += ret.first;
- num_drop_particle += ret.second;
- }
- else if( DetectorHistory* h = std::any_cast<DetectorHistory>(&i.second) ) {
- auto [nhit, npart] = drop_history(*h);
- num_drop_hit += nhit;
- num_drop_particle += npart;
- }
- }
- }
- info("%s+++ Dropped history of %6ld hits %6ld particles",
- context.event->id(), num_drop_hit, num_drop_particle);
+ auto& inputs = context.event->get_segment(m_input);
+ std::size_t num_drop_hit = 0;
+ std::size_t num_drop_particle = 0;
+ for( auto& i : inputs ) {
+ Key key(i.first);
+ auto im = std::find(m_masks.begin(), m_masks.end(), key.mask());
+ if( im != m_masks.end() ) {
+ if ( DepositMapping* m = std::any_cast<DepositMapping>(&i.second) ) {
+ auto ret = drop_history(*m);
+ num_drop_hit += ret.first;
+ num_drop_particle += ret.second;
+ }
+ else if ( DepositVector* v = std::any_cast<DepositVector>(&i.second) ) {
+ auto ret = drop_history(*v);
+ num_drop_hit += ret.first;
+ num_drop_particle += ret.second;
+ }
+ else if( DetectorHistory* h = std::any_cast<DetectorHistory>(&i.second) ) {
+ auto [nhit, npart] = drop_history(*h);
+ num_drop_hit += nhit;
+ num_drop_particle += npart;
+ }
+ }
+ }
+ info("%s+++ Dropped history of %6ld hits %6ld particles",
+ context.event->id(), num_drop_hit, num_drop_particle);
}
};
} // End namespace digi
diff --git a/DDDigi/plugins/DigiIPCreate.cpp b/DDDigi/plugins/DigiIPCreate.cpp
index 70f973d15..2ebd4a7fc 100644
--- a/DDDigi/plugins/DigiIPCreate.cpp
+++ b/DDDigi/plugins/DigiIPCreate.cpp
@@ -31,26 +31,26 @@ namespace dd4hep {
public:
/// Standard constructor
DigiIPCreate(const DigiKernel& krnl, const std::string& nam)
- : DigiEventAction(krnl, nam)
+ : DigiEventAction(krnl, nam)
{
- declareProperty("offset_ip", m_offset_ip);
- declareProperty("sigma_ip", m_sigma_ip);
- declareProperty("interaction_point", m_interaction_point);
+ declareProperty("offset_ip", m_offset_ip);
+ declareProperty("sigma_ip", m_sigma_ip);
+ declareProperty("interaction_point", m_interaction_point);
}
/// Main functional callback
virtual void execute(DigiContext& context) const override final {
- auto& rndm = context.randomGenerator();
- double theta = rndm.uniform(0.0, M_PI); // theta = [0,pi]
- double phi = rndm.uniform(0.0, 2.0*M_PI); // phi = [0,2*pi]
- double radius = rndm.uniform(0.0, 1.0); // radius = [0,1]
- radius = std::sqrt( -std::log(radius) ); // radius in a gaussian distribution
-
- double st = std::sin(theta);
- double xx = radius * st * std::cos(phi) * m_sigma_ip.X() + m_offset_ip.X();
- double yy = radius * st * std::sin(phi) * m_sigma_ip.Y() + m_offset_ip.Y();
- double zz = radius * std::cos(theta) * m_sigma_ip.Z() + m_offset_ip.Z();
- m_interaction_point = Position(xx, yy, zz);
+ auto& rndm = context.randomGenerator();
+ double theta = rndm.uniform(0.0, M_PI); // theta = [0,pi]
+ double phi = rndm.uniform(0.0, 2.0*M_PI); // phi = [0,2*pi]
+ double radius = rndm.uniform(0.0, 1.0); // radius = [0,1]
+ radius = std::sqrt( -std::log(radius) ); // radius in a gaussian distribution
+
+ double st = std::sin(theta);
+ double xx = radius * st * std::cos(phi) * m_sigma_ip.X() + m_offset_ip.X();
+ double yy = radius * st * std::sin(phi) * m_sigma_ip.Y() + m_offset_ip.Y();
+ double zz = radius * std::cos(theta) * m_sigma_ip.Z() + m_offset_ip.Z();
+ m_interaction_point = Position(xx, yy, zz);
}
};
} // End namespace digi
diff --git a/DDDigi/plugins/DigiSimpleADCResponse.cpp b/DDDigi/plugins/DigiSimpleADCResponse.cpp
index d4ce388c4..7c09d0734 100644
--- a/DDDigi/plugins/DigiSimpleADCResponse.cpp
+++ b/DDDigi/plugins/DigiSimpleADCResponse.cpp
@@ -33,37 +33,37 @@ namespace dd4hep {
double m_adc_offset { 0e0 };
std::size_t m_adc_resolution { 1024 };
- public:
+ public:
/// Standard constructor
DigiSimpleADCResponse(const DigiKernel& krnl, const std::string& nam)
- : DigiDepositsProcessor(krnl, nam)
+ : DigiDepositsProcessor(krnl, nam)
{
- declareProperty("saturation", m_signal_saturation);
- declareProperty("adc_resolution", m_adc_resolution);
- declareProperty("response_postfix", m_response_postfix);
- declareProperty("history_postfix", m_history_postfix);
- DEPOSIT_PROCESSOR_BIND_HANDLERS(DigiSimpleADCResponse::emulate_adc)
+ declareProperty("saturation", m_signal_saturation);
+ declareProperty("adc_resolution", m_adc_resolution);
+ declareProperty("response_postfix", m_response_postfix);
+ declareProperty("history_postfix", m_history_postfix);
+ DEPOSIT_PROCESSOR_BIND_HANDLERS(DigiSimpleADCResponse::emulate_adc);
}
/// Create container with ADC counts and register it to the output segment
template <typename T>
void emulate_adc(DigiContext& context, const T& input, work_t& work, const predicate_t& predicate) const {
- const char* tag = context.event->id();
- std::string postfix = predicate.segmentation ? "."+predicate.segmentation->identifier(predicate.id) : std::string();
- std::string response_name = input.name + postfix + m_response_postfix;
- DetectorResponse response(response_name, work.environ.output.mask);
- for( const auto& dep : input ) {
- if ( predicate(dep) ) {
- CellID cell = dep.first;
- const auto& depo = dep.second;
- double offset_ene = depo.deposit-m_adc_offset;
- ADCValue::value_t adc_count = std::round(((offset_ene) * m_adc_resolution) / m_signal_saturation);
- adc_count = std::min(adc_count, ADCValue::value_t(m_adc_resolution));
- response.emplace(cell, {adc_count, ADCValue::address_t(cell)});
- }
- }
- info("%s+++ %-32s %6ld ADC values. Input: %-32s %6ld deposits", tag,
- response_name.c_str(), response.size(), input.name.c_str(), input.size());
+ const char* tag = context.event->id();
+ std::string postfix = predicate.segmentation ? "."+predicate.segmentation->identifier(predicate.id) : std::string();
+ std::string response_name = input.name + postfix + m_response_postfix;
+ DetectorResponse response(response_name, work.environ.output.mask);
+ for( const auto& dep : input ) {
+ if ( predicate(dep) ) {
+ CellID cell = dep.first;
+ const auto& depo = dep.second;
+ double offset_ene = depo.deposit-m_adc_offset;
+ ADCValue::value_t adc_count = std::round(((offset_ene) * m_adc_resolution) / m_signal_saturation);
+ adc_count = std::min(adc_count, ADCValue::value_t(m_adc_resolution));
+ response.emplace(cell, {adc_count, ADCValue::address_t(cell)});
+ }
+ }
+ info("%s+++ %-32s %6ld ADC values. Input: %-32s %6ld deposits", tag,
+ response_name.c_str(), response.size(), input.name.c_str(), input.size());
work.environ.output.data.put(response.key, std::move(response));
}
};
diff --git a/DDG4/edm4hep/Geant4Output2EDM4hep.cpp b/DDG4/edm4hep/Geant4Output2EDM4hep.cpp
index 6a71a20f5..02ca17daa 100644
--- a/DDG4/edm4hep/Geant4Output2EDM4hep.cpp
+++ b/DDG4/edm4hep/Geant4Output2EDM4hep.cpp
@@ -235,13 +235,9 @@ void Geant4Output2EDM4hep::commit( OutputContext<G4Event>& /* ctxt */) {
G4AutoLock protection_lock(&action_mutex);
m_frame.put( std::move(m_particles), "MCParticles");
for (auto it = m_trackerHits.begin(); it != m_trackerHits.end(); ++it) {
- it->second.getDatamodelRegistryIndex();
- it->second.size();
m_frame.put( std::move(it->second), it->first);
}
for (auto& [colName, calorimeterHits] : m_calorimeterHits) {
- calorimeterHits.first.getDatamodelRegistryIndex();
- calorimeterHits.second.getDatamodelRegistryIndex();
m_frame.put( std::move(calorimeterHits.first), colName);
m_frame.put( std::move(calorimeterHits.second), colName + "Contributions");
}
@@ -276,24 +272,11 @@ void Geant4Output2EDM4hep::saveRun(const G4Run* run) {
void Geant4Output2EDM4hep::begin(const G4Event* event) {
/// Create event frame object
- G4HCofThisEvent* hce = event->GetHCofThisEvent();
m_eventNo = event->GetEventID();
m_frame = {};
m_particles = {};
m_trackerHits.clear();
m_calorimeterHits.clear();
- if ( hce ) {
- int nCol = hce->GetNumberOfCollections();
- for (int i = 0; i < nCol; ++i) {
- Geant4HitCollection* coll = dynamic_cast<Geant4HitCollection*>(hce->GetHC(i));
- if( typeid( Geant4Tracker::Hit ) == coll->type().type() )
- m_trackerHits[coll->GetName()] = {};
- else if( typeid( Geant4Calorimeter::Hit ) == coll->type().type() )
- m_calorimeterHits[coll->GetName()] = {};
- else
- except("begin(Event): Unknown container type!");
- }
- }
}
/// Data conversion interface for MC particles to EDM4hep format
@@ -458,8 +441,10 @@ void Geant4Output2EDM4hep::saveCollection(OutputContext<G4Event>& /*ctxt*/, G4VH
debug("+++ Saving EDM4hep collection %s with %d entries.", colName.c_str(), int(nhits));
//-------------------------------------------------------------------
if( typeid( Geant4Tracker::Hit ) == coll->type().type() ){
+ // Create the hit container even if there are no entries!
+ auto& hits = m_trackerHits[colName] = {};
for(unsigned i=0 ; i < nhits ; ++i){
- auto sth = m_trackerHits[colName]->create();
+ auto sth = hits->create();
const Geant4Tracker::Hit* hit = coll->hit(i);
const Geant4Tracker::Hit::Contribution& t = hit->truth;
int trackID = pm->particleID(t.trackID);
@@ -487,8 +472,10 @@ void Geant4Output2EDM4hep::saveCollection(OutputContext<G4Event>& /*ctxt*/, G4VH
else if( typeid( Geant4Calorimeter::Hit ) == coll->type().type() ){
Geant4Sensitive* sd = coll->sensitive();
int hit_creation_mode = sd->hitCreationMode();
+ // Create the hit container even if there are no entries!
+ auto& hits = m_calorimeterHits[colName] = {};
for(unsigned i=0 ; i < nhits ; ++i){
- auto sch = m_calorimeterHits[colName].first->create();
+ auto sch = hits.first->create();
const Geant4Calorimeter::Hit* hit = coll->hit(i);
const auto& pos = hit->position;
sch.setCellID( hit->cellID );
@@ -499,7 +486,7 @@ void Geant4Output2EDM4hep::saveCollection(OutputContext<G4Event>& /*ctxt*/, G4VH
// now add the individual step contributions
for(auto ci=hit->truth.begin(); ci != hit->truth.end(); ++ci){
- auto sCaloHitCont = m_calorimeterHits[colName].second->create();
+ auto sCaloHitCont = hits.second->create();
sch.addToContributions( sCaloHitCont );
const Geant4HitData::Contribution& c = *ci;
--
GitLab