From df0b9aff769ee02288a2535fec89259d3954df4f Mon Sep 17 00:00:00 2001 From: mhkim Date: Mon, 22 Sep 2025 17:34:00 -0500 Subject: [PATCH 01/16] First commit for PulseDigi --- src/algorithms/digi/PulseDigi.cc | 84 ++++++++++++++++++++++++++++++++ 1 file changed, 84 insertions(+) create mode 100644 src/algorithms/digi/PulseDigi.cc diff --git a/src/algorithms/digi/PulseDigi.cc b/src/algorithms/digi/PulseDigi.cc new file mode 100644 index 0000000000..85bc4245c4 --- /dev/null +++ b/src/algorithms/digi/PulseDigi.cc @@ -0,0 +1,84 @@ +// SPDX-License-Identifier: LGPL-3.0-or-later +// Copyright (C) 2025 Minho Kim +// + +#include +#include +#include +#include + +#include "PulseDigi.h" + +namespace eicrecon { + + void PulseDigi::init() {} + + void PulseDigi::process(const PulseDigi::Input& input, const PulseDigi::Output& output) const { + const auto [in_pulses] = input; + auto [out_digi_hits, out_digi_samples] = output; + + for (const auto& pulse : *in_pulses) { + bool sample_boundary = false; + double amplitude_max = -std::numeric_limits::infinity(); + std::vector amplitudes_neighboring(2); + double toa = 0; + double tot = 0; + bool totInProgress = false; + bool totComplete = false; + + auto out_digi_hit = out_digi_hits->create(); + + for (std::size_t i = 0; i < pulse.getAmplitude().size(); i++) { + double t = pulse.getTime() + i * pulse.getInterval(); + sample_boundary = is_sample_boundary(t, m_cfg.sample_period); + amplitude = pulse.getAmplitude()[i]; + amplitudes_neighboring[i % 2] = amplitude; + amplitude_max = std::max(amplitude_max, amplitude); + + if(!totInProgress && amplitude > m_cfg.threshold) { + toa = get_crossing_time(m_cfg.threshold, + t, t - pulse.getInterval(), + amplitudes_neiboring[i % 2], amplitudes_neighboring[1 - i % 2]); + totInProgress = true; + } + + if(totInProgress && !totComplete && amplitude < m_cfg.threshold) { + totComplete = true; + totInProgress = false; + tot = (((m_cfg.threshold - amps_neighbor[i % 2]) * pulse.getInterval()) / (amps_neighbor[i % 2] - amps_neighbor[1 - (i % 2)])) + t - toa; + } + + if(sample_boundary){ + auto out_digi_sample = out_digi_samples->create(); + + out_digi_sample.setADC(amplitude_max); + out_digi_sample.setTimeOfArrival(toa); + out_digi_sample.setTimeOverThreshold(tot); + out_digi_sample.setTOTInProgress(totInProgress); + out_digi_sample.setTOTComplete(totComplete); + out_digi_hit.addToSamples(out_digi_hit); + + amplitude_max = -std::numeric_limits::infinity(); + toa = 0; + tot = 0; + } + } + + out_digi_hit.setCellID(pulse.getCellID()); + out_digi_hit.setSamplePahse(); + out_digi_hit.setTimeStap(); + } + } // PulseDigi:process + + bool PulseDigi::is_sample_boundary(double t, double period) const { + double n = std::llround(t / period); + return std::abs(t - n * period) < 1.0e-5; + } + + double PulseDigi::get_crossing_time(double thres, double t1, double t2, double amp1, double amp2) const { + double numerator = (thres - amp2) * (t2 - t1); + double denomenator = amp2 - amp1; + double added = t2; + return (numerator / denominator) + added; + } +} // namespace eicrecon From 43b2761db01d45a308c73c465973da9163fc7e38 Mon Sep 17 00:00:00 2001 From: mhkim Date: Mon, 22 Sep 2025 17:37:54 -0500 Subject: [PATCH 02/16] code-formatter applied --- src/algorithms/digi/PulseDigi.cc | 142 ++++++++++++++++--------------- 1 file changed, 72 insertions(+), 70 deletions(-) diff --git a/src/algorithms/digi/PulseDigi.cc b/src/algorithms/digi/PulseDigi.cc index 85bc4245c4..c88e65f23e 100644 --- a/src/algorithms/digi/PulseDigi.cc +++ b/src/algorithms/digi/PulseDigi.cc @@ -11,74 +11,76 @@ namespace eicrecon { - void PulseDigi::init() {} - - void PulseDigi::process(const PulseDigi::Input& input, const PulseDigi::Output& output) const { - const auto [in_pulses] = input; - auto [out_digi_hits, out_digi_samples] = output; - - for (const auto& pulse : *in_pulses) { - bool sample_boundary = false; - double amplitude_max = -std::numeric_limits::infinity(); - std::vector amplitudes_neighboring(2); - double toa = 0; - double tot = 0; - bool totInProgress = false; - bool totComplete = false; - - auto out_digi_hit = out_digi_hits->create(); - - for (std::size_t i = 0; i < pulse.getAmplitude().size(); i++) { - double t = pulse.getTime() + i * pulse.getInterval(); - sample_boundary = is_sample_boundary(t, m_cfg.sample_period); - amplitude = pulse.getAmplitude()[i]; - amplitudes_neighboring[i % 2] = amplitude; - amplitude_max = std::max(amplitude_max, amplitude); - - if(!totInProgress && amplitude > m_cfg.threshold) { - toa = get_crossing_time(m_cfg.threshold, - t, t - pulse.getInterval(), - amplitudes_neiboring[i % 2], amplitudes_neighboring[1 - i % 2]); - totInProgress = true; - } - - if(totInProgress && !totComplete && amplitude < m_cfg.threshold) { - totComplete = true; - totInProgress = false; - tot = (((m_cfg.threshold - amps_neighbor[i % 2]) * pulse.getInterval()) / (amps_neighbor[i % 2] - amps_neighbor[1 - (i % 2)])) + t - toa; - } - - if(sample_boundary){ - auto out_digi_sample = out_digi_samples->create(); - - out_digi_sample.setADC(amplitude_max); - out_digi_sample.setTimeOfArrival(toa); - out_digi_sample.setTimeOverThreshold(tot); - out_digi_sample.setTOTInProgress(totInProgress); - out_digi_sample.setTOTComplete(totComplete); - out_digi_hit.addToSamples(out_digi_hit); - - amplitude_max = -std::numeric_limits::infinity(); - toa = 0; - tot = 0; - } - } - - out_digi_hit.setCellID(pulse.getCellID()); - out_digi_hit.setSamplePahse(); - out_digi_hit.setTimeStap(); - } - } // PulseDigi:process - - bool PulseDigi::is_sample_boundary(double t, double period) const { - double n = std::llround(t / period); - return std::abs(t - n * period) < 1.0e-5; - } - - double PulseDigi::get_crossing_time(double thres, double t1, double t2, double amp1, double amp2) const { - double numerator = (thres - amp2) * (t2 - t1); - double denomenator = amp2 - amp1; - double added = t2; - return (numerator / denominator) + added; - } +void PulseDigi::init() {} + +void PulseDigi::process(const PulseDigi::Input& input, const PulseDigi::Output& output) const { + const auto [in_pulses] = input; + auto [out_digi_hits, out_digi_samples] = output; + + for (const auto& pulse : *in_pulses) { + bool sample_boundary = false; + double amplitude_max = -std::numeric_limits::infinity(); + std::vector amplitudes_neighboring(2); + double toa = 0; + double tot = 0; + bool totInProgress = false; + bool totComplete = false; + + auto out_digi_hit = out_digi_hits->create(); + + for (std::size_t i = 0; i < pulse.getAmplitude().size(); i++) { + double t = pulse.getTime() + i * pulse.getInterval(); + sample_boundary = is_sample_boundary(t, m_cfg.sample_period); + amplitude = pulse.getAmplitude()[i]; + amplitudes_neighboring[i % 2] = amplitude; + amplitude_max = std::max(amplitude_max, amplitude); + + if (!totInProgress && amplitude > m_cfg.threshold) { + toa = get_crossing_time(m_cfg.threshold, t, t - pulse.getInterval(), + amplitudes_neiboring[i % 2], amplitudes_neighboring[1 - i % 2]); + totInProgress = true; + } + + if (totInProgress && !totComplete && amplitude < m_cfg.threshold) { + totComplete = true; + totInProgress = false; + tot = (((m_cfg.threshold - amps_neighbor[i % 2]) * pulse.getInterval()) / + (amps_neighbor[i % 2] - amps_neighbor[1 - (i % 2)])) + + t - toa; + } + + if (sample_boundary) { + auto out_digi_sample = out_digi_samples->create(); + + out_digi_sample.setADC(amplitude_max); + out_digi_sample.setTimeOfArrival(toa); + out_digi_sample.setTimeOverThreshold(tot); + out_digi_sample.setTOTInProgress(totInProgress); + out_digi_sample.setTOTComplete(totComplete); + out_digi_hit.addToSamples(out_digi_hit); + + amplitude_max = -std::numeric_limits::infinity(); + toa = 0; + tot = 0; + } + } + + out_digi_hit.setCellID(pulse.getCellID()); + out_digi_hit.setSamplePahse(); + out_digi_hit.setTimeStap(); + } +} // PulseDigi:process + +bool PulseDigi::is_sample_boundary(double t, double period) const { + double n = std::llround(t / period); + return std::abs(t - n * period) < 1.0e-5; +} + +double PulseDigi::get_crossing_time(double thres, double t1, double t2, double amp1, + double amp2) const { + double numerator = (thres - amp2) * (t2 - t1); + double denomenator = amp2 - amp1; + double added = t2; + return (numerator / denominator) + added; +} } // namespace eicrecon From 11e9eeaaeaf41f15e4371fc2c9a034b36bee6f76 Mon Sep 17 00:00:00 2001 From: mhkim Date: Mon, 22 Sep 2025 17:41:26 -0500 Subject: [PATCH 03/16] Fix typo --- src/algorithms/digi/PulseDigi.cc | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/algorithms/digi/PulseDigi.cc b/src/algorithms/digi/PulseDigi.cc index c88e65f23e..9ba2ede81b 100644 --- a/src/algorithms/digi/PulseDigi.cc +++ b/src/algorithms/digi/PulseDigi.cc @@ -44,8 +44,8 @@ void PulseDigi::process(const PulseDigi::Input& input, const PulseDigi::Output& if (totInProgress && !totComplete && amplitude < m_cfg.threshold) { totComplete = true; totInProgress = false; - tot = (((m_cfg.threshold - amps_neighbor[i % 2]) * pulse.getInterval()) / - (amps_neighbor[i % 2] - amps_neighbor[1 - (i % 2)])) + + tot = get_crossing_time(((m_cfg.threshold - amps_neighbor[i % 2]) * pulse.getInterval()) / + (amps_neighbor[i % 2] - amps_neighbor[1 - (i % 2)])) + t - toa; } From 6202d549db6265453e6ad72a9f2e9bd401f41727 Mon Sep 17 00:00:00 2001 From: mhkim Date: Sun, 28 Sep 2025 19:10:56 -0500 Subject: [PATCH 04/16] More optimized to run efficiently --- src/algorithms/digi/PulseDigi.cc | 148 ++++++++++++++++--------------- 1 file changed, 76 insertions(+), 72 deletions(-) diff --git a/src/algorithms/digi/PulseDigi.cc b/src/algorithms/digi/PulseDigi.cc index 9ba2ede81b..7241dcbbb8 100644 --- a/src/algorithms/digi/PulseDigi.cc +++ b/src/algorithms/digi/PulseDigi.cc @@ -9,78 +9,82 @@ #include "PulseDigi.h" +class HGCROCRawSample { + public: + HGCROCRawSample(std::size_t n_meas) + : meas_types(n_meas, 0), + meas_values(n_meas, 0) {} + + void addAmplitude(std::size_t idx, float amp) { meas_values[idx] = std::max(meas_values[idx], amp); } + void addUpCrossTime(std::size_t idx, double t) { + meas_types[idx] = 1; + meas_values[idx] = t; + } + void addDownCrossTime(std::size_t idx, double t) { + meas_types[idx] = 2; + meas_values[idx] = t; + } + + std::vector getMeasTypes() { return meas_types; } + std::vector getMeasValues() { return meas_values; } + + private: + std::vector meas_types; + std::vector meas_values; +}; + namespace eicrecon { -void PulseDigi::init() {} - -void PulseDigi::process(const PulseDigi::Input& input, const PulseDigi::Output& output) const { - const auto [in_pulses] = input; - auto [out_digi_hits, out_digi_samples] = output; - - for (const auto& pulse : *in_pulses) { - bool sample_boundary = false; - double amplitude_max = -std::numeric_limits::infinity(); - std::vector amplitudes_neighboring(2); - double toa = 0; - double tot = 0; - bool totInProgress = false; - bool totComplete = false; - - auto out_digi_hit = out_digi_hits->create(); - - for (std::size_t i = 0; i < pulse.getAmplitude().size(); i++) { - double t = pulse.getTime() + i * pulse.getInterval(); - sample_boundary = is_sample_boundary(t, m_cfg.sample_period); - amplitude = pulse.getAmplitude()[i]; - amplitudes_neighboring[i % 2] = amplitude; - amplitude_max = std::max(amplitude_max, amplitude); - - if (!totInProgress && amplitude > m_cfg.threshold) { - toa = get_crossing_time(m_cfg.threshold, t, t - pulse.getInterval(), - amplitudes_neiboring[i % 2], amplitudes_neighboring[1 - i % 2]); - totInProgress = true; - } - - if (totInProgress && !totComplete && amplitude < m_cfg.threshold) { - totComplete = true; - totInProgress = false; - tot = get_crossing_time(((m_cfg.threshold - amps_neighbor[i % 2]) * pulse.getInterval()) / - (amps_neighbor[i % 2] - amps_neighbor[1 - (i % 2)])) + - t - toa; - } - - if (sample_boundary) { - auto out_digi_sample = out_digi_samples->create(); - - out_digi_sample.setADC(amplitude_max); - out_digi_sample.setTimeOfArrival(toa); - out_digi_sample.setTimeOverThreshold(tot); - out_digi_sample.setTOTInProgress(totInProgress); - out_digi_sample.setTOTComplete(totComplete); - out_digi_hit.addToSamples(out_digi_hit); - - amplitude_max = -std::numeric_limits::infinity(); - toa = 0; - tot = 0; - } - } - - out_digi_hit.setCellID(pulse.getCellID()); - out_digi_hit.setSamplePahse(); - out_digi_hit.setTimeStap(); - } -} // PulseDigi:process - -bool PulseDigi::is_sample_boundary(double t, double period) const { - double n = std::llround(t / period); - return std::abs(t - n * period) < 1.0e-5; -} - -double PulseDigi::get_crossing_time(double thres, double t1, double t2, double amp1, - double amp2) const { - double numerator = (thres - amp2) * (t2 - t1); - double denomenator = amp2 - amp1; - double added = t2; - return (numerator / denominator) + added; -} + void PulseDigi::init() {} + + void PulseDigi::process(const PulseDigi::Input& input, const PulseDigi::Output& output) const { + const auto [in_pulses] = input; + auto [out_digi_hits, out_digi_samples] = output; + + for (const auto& pulse : *in_pulses) { + float pulse_t = pulse.getTime(); + float pulse_dt = pulse.getInterval(); + std::size_t n_amplitudes = pulse.getAmplitude().size();; + + const std::size_t idx_begin = static_cast(std::floor(pulse_t / m_cfg.time_window)); + const std::size_t idx_end = static_cast(std::floor((pulse_t + (n_amplitudes - 1) * pulse_dt) / m_cfg.time_window)); + + HGCROCRawSample raw_sample(idx_end - idx_begin + 1); + + bool tot_progress = false; + bool tot_complete = false; + + for (std::size_t i = 0; i < n_amplitudes; i++) { + double t = pulse_t + i * pulse_dt; + const std::size_t idx = static_cast(std::floor(t / m_cfg.time_window)) - idx_begin; + + if(!tot_in_progess){ + raw_sample.addAmplitude(idx, pulse.getAmplitude()[i]); + + if(pulse.getAmplitude()[i] > m_cfg.threshold) { + raw_sample.addUpCrossTime(idx, get_crossing_time(m_cfg.threshold, t, pulse_dt, + pulse.getAmplitude()[i], pulse.getAmplitude()[i-1])); + tot_progress = true; + tot_complete = false; + } + } + + if(tot_in_progress && !tot_complete && pulse.getAmplitude()[i] < m_cfg.threshold){ + pulse_info.addDownCrossTime(idx, get_crossing_time(m_cfg.threshold, t, pulse_dt, + pulse.getAmplitude()[i], pulse.getAmplitude()[i-1])); + tot_complete = true; + tot_in_progess = false; + } + } + } + } // PulseDigi:process + + + double PulseDigi::get_crossing_time(double thres, double t1, double t2, double amp1, + double amp2) const { + double numerator = (thres - amp2) * (t2 - t1); + double denomenator = amp2 - amp1; + double added = t2; + return (numerator / denominator) + added; + } } // namespace eicrecon From 89199f6f37ffacebfeb49f8d5a2872d0e8ad2358 Mon Sep 17 00:00:00 2001 From: mhkim Date: Sun, 28 Sep 2025 19:12:38 -0500 Subject: [PATCH 05/16] clang-format applied --- src/algorithms/digi/PulseDigi.cc | 150 ++++++++++++++++--------------- 1 file changed, 77 insertions(+), 73 deletions(-) diff --git a/src/algorithms/digi/PulseDigi.cc b/src/algorithms/digi/PulseDigi.cc index 7241dcbbb8..54d33e096c 100644 --- a/src/algorithms/digi/PulseDigi.cc +++ b/src/algorithms/digi/PulseDigi.cc @@ -10,81 +10,85 @@ #include "PulseDigi.h" class HGCROCRawSample { - public: - HGCROCRawSample(std::size_t n_meas) - : meas_types(n_meas, 0), - meas_values(n_meas, 0) {} - - void addAmplitude(std::size_t idx, float amp) { meas_values[idx] = std::max(meas_values[idx], amp); } - void addUpCrossTime(std::size_t idx, double t) { - meas_types[idx] = 1; - meas_values[idx] = t; - } - void addDownCrossTime(std::size_t idx, double t) { - meas_types[idx] = 2; - meas_values[idx] = t; - } - - std::vector getMeasTypes() { return meas_types; } - std::vector getMeasValues() { return meas_values; } - - private: - std::vector meas_types; - std::vector meas_values; +public: + HGCROCRawSample(std::size_t n_meas) : meas_types(n_meas, 0), meas_values(n_meas, 0) {} + + void addAmplitude(std::size_t idx, float amp) { + meas_values[idx] = std::max(meas_values[idx], amp); + } + void addUpCrossTime(std::size_t idx, double t) { + meas_types[idx] = 1; + meas_values[idx] = t; + } + void addDownCrossTime(std::size_t idx, double t) { + meas_types[idx] = 2; + meas_values[idx] = t; + } + + std::vector getMeasTypes() { return meas_types; } + std::vector getMeasValues() { return meas_values; } + +private: + std::vector meas_types; + std::vector meas_values; }; namespace eicrecon { - void PulseDigi::init() {} - - void PulseDigi::process(const PulseDigi::Input& input, const PulseDigi::Output& output) const { - const auto [in_pulses] = input; - auto [out_digi_hits, out_digi_samples] = output; - - for (const auto& pulse : *in_pulses) { - float pulse_t = pulse.getTime(); - float pulse_dt = pulse.getInterval(); - std::size_t n_amplitudes = pulse.getAmplitude().size();; - - const std::size_t idx_begin = static_cast(std::floor(pulse_t / m_cfg.time_window)); - const std::size_t idx_end = static_cast(std::floor((pulse_t + (n_amplitudes - 1) * pulse_dt) / m_cfg.time_window)); - - HGCROCRawSample raw_sample(idx_end - idx_begin + 1); - - bool tot_progress = false; - bool tot_complete = false; - - for (std::size_t i = 0; i < n_amplitudes; i++) { - double t = pulse_t + i * pulse_dt; - const std::size_t idx = static_cast(std::floor(t / m_cfg.time_window)) - idx_begin; - - if(!tot_in_progess){ - raw_sample.addAmplitude(idx, pulse.getAmplitude()[i]); - - if(pulse.getAmplitude()[i] > m_cfg.threshold) { - raw_sample.addUpCrossTime(idx, get_crossing_time(m_cfg.threshold, t, pulse_dt, - pulse.getAmplitude()[i], pulse.getAmplitude()[i-1])); - tot_progress = true; - tot_complete = false; - } - } - - if(tot_in_progress && !tot_complete && pulse.getAmplitude()[i] < m_cfg.threshold){ - pulse_info.addDownCrossTime(idx, get_crossing_time(m_cfg.threshold, t, pulse_dt, - pulse.getAmplitude()[i], pulse.getAmplitude()[i-1])); - tot_complete = true; - tot_in_progess = false; - } - } - } - } // PulseDigi:process - - - double PulseDigi::get_crossing_time(double thres, double t1, double t2, double amp1, - double amp2) const { - double numerator = (thres - amp2) * (t2 - t1); - double denomenator = amp2 - amp1; - double added = t2; - return (numerator / denominator) + added; - } +void PulseDigi::init() {} + +void PulseDigi::process(const PulseDigi::Input& input, const PulseDigi::Output& output) const { + const auto [in_pulses] = input; + auto [out_digi_hits, out_digi_samples] = output; + + for (const auto& pulse : *in_pulses) { + float pulse_t = pulse.getTime(); + float pulse_dt = pulse.getInterval(); + std::size_t n_amplitudes = pulse.getAmplitude().size(); + ; + + const std::size_t idx_begin = static_cast(std::floor(pulse_t / m_cfg.time_window)); + const std::size_t idx_end = static_cast( + std::floor((pulse_t + (n_amplitudes - 1) * pulse_dt) / m_cfg.time_window)); + + HGCROCRawSample raw_sample(idx_end - idx_begin + 1); + + bool tot_progress = false; + bool tot_complete = false; + + for (std::size_t i = 0; i < n_amplitudes; i++) { + double t = pulse_t + i * pulse_dt; + const std::size_t idx = + static_cast(std::floor(t / m_cfg.time_window)) - idx_begin; + + if (!tot_in_progess) { + raw_sample.addAmplitude(idx, pulse.getAmplitude()[i]); + + if (pulse.getAmplitude()[i] > m_cfg.threshold) { + raw_sample.addUpCrossTime(idx, get_crossing_time(m_cfg.threshold, t, pulse_dt, + pulse.getAmplitude()[i], + pulse.getAmplitude()[i - 1])); + tot_progress = true; + tot_complete = false; + } + } + + if (tot_in_progress && !tot_complete && pulse.getAmplitude()[i] < m_cfg.threshold) { + pulse_info.addDownCrossTime(idx, get_crossing_time(m_cfg.threshold, t, pulse_dt, + pulse.getAmplitude()[i], + pulse.getAmplitude()[i - 1])); + tot_complete = true; + tot_in_progess = false; + } + } + } +} // PulseDigi:process + +double PulseDigi::get_crossing_time(double thres, double t1, double t2, double amp1, + double amp2) const { + double numerator = (thres - amp2) * (t2 - t1); + double denomenator = amp2 - amp1; + double added = t2; + return (numerator / denominator) + added; +} } // namespace eicrecon From 5e57416db5554b1291dc28414b3cfb19e899f3bd Mon Sep 17 00:00:00 2001 From: mhkim Date: Sun, 28 Sep 2025 19:15:35 -0500 Subject: [PATCH 06/16] Some typos were fixed --- src/algorithms/digi/PulseDigi.cc | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/src/algorithms/digi/PulseDigi.cc b/src/algorithms/digi/PulseDigi.cc index 54d33e096c..0aef1fe5c3 100644 --- a/src/algorithms/digi/PulseDigi.cc +++ b/src/algorithms/digi/PulseDigi.cc @@ -61,7 +61,7 @@ void PulseDigi::process(const PulseDigi::Input& input, const PulseDigi::Output& const std::size_t idx = static_cast(std::floor(t / m_cfg.time_window)) - idx_begin; - if (!tot_in_progess) { + if (!tot_progess) { raw_sample.addAmplitude(idx, pulse.getAmplitude()[i]); if (pulse.getAmplitude()[i] > m_cfg.threshold) { @@ -73,12 +73,12 @@ void PulseDigi::process(const PulseDigi::Input& input, const PulseDigi::Output& } } - if (tot_in_progress && !tot_complete && pulse.getAmplitude()[i] < m_cfg.threshold) { + if (tot_progress && !tot_complete && pulse.getAmplitude()[i] < m_cfg.threshold) { pulse_info.addDownCrossTime(idx, get_crossing_time(m_cfg.threshold, t, pulse_dt, pulse.getAmplitude()[i], pulse.getAmplitude()[i - 1])); tot_complete = true; - tot_in_progess = false; + tot_progess = false; } } } From 766dc84762028853689835ec47953837113d067c Mon Sep 17 00:00:00 2001 From: mhkim Date: Mon, 29 Sep 2025 11:04:47 -0500 Subject: [PATCH 07/16] Add some explanations --- src/algorithms/digi/PulseDigi.cc | 8 +++++++- 1 file changed, 7 insertions(+), 1 deletion(-) diff --git a/src/algorithms/digi/PulseDigi.cc b/src/algorithms/digi/PulseDigi.cc index 0aef1fe5c3..cdecdf37d1 100644 --- a/src/algorithms/digi/PulseDigi.cc +++ b/src/algorithms/digi/PulseDigi.cc @@ -29,6 +29,9 @@ class HGCROCRawSample { std::vector getMeasValues() { return meas_values; } private: + // meas_type 0: Pulse crossed above the threshold. The crossing time is pushed back onto meas_values. + // meas_type 1: Pulse crossed below the threshold. The crossing time is pushed back onto meas_values. + // meas_type 2: Pulse didn't cross the threshold. The maximum amplitude is pushed back onto meas_values. std::vector meas_types; std::vector meas_values; }; @@ -45,8 +48,8 @@ void PulseDigi::process(const PulseDigi::Input& input, const PulseDigi::Output& float pulse_t = pulse.getTime(); float pulse_dt = pulse.getInterval(); std::size_t n_amplitudes = pulse.getAmplitude().size(); - ; + // Estimate the number of samples. const std::size_t idx_begin = static_cast(std::floor(pulse_t / m_cfg.time_window)); const std::size_t idx_end = static_cast( std::floor((pulse_t + (n_amplitudes - 1) * pulse_dt) / m_cfg.time_window)); @@ -61,9 +64,11 @@ void PulseDigi::process(const PulseDigi::Input& input, const PulseDigi::Output& const std::size_t idx = static_cast(std::floor(t / m_cfg.time_window)) - idx_begin; + // Pulse didn't cross the threshold. if (!tot_progess) { raw_sample.addAmplitude(idx, pulse.getAmplitude()[i]); + // Pulse crossed above the threshold. if (pulse.getAmplitude()[i] > m_cfg.threshold) { raw_sample.addUpCrossTime(idx, get_crossing_time(m_cfg.threshold, t, pulse_dt, pulse.getAmplitude()[i], @@ -73,6 +78,7 @@ void PulseDigi::process(const PulseDigi::Input& input, const PulseDigi::Output& } } + // Pulse crossed below the threshold. if (tot_progress && !tot_complete && pulse.getAmplitude()[i] < m_cfg.threshold) { pulse_info.addDownCrossTime(idx, get_crossing_time(m_cfg.threshold, t, pulse_dt, pulse.getAmplitude()[i], From 542443433821e0445a7c68688d6a704a349d9c47 Mon Sep 17 00:00:00 2001 From: mhkim Date: Mon, 29 Sep 2025 11:11:10 -0500 Subject: [PATCH 08/16] Combined two member functions --- src/algorithms/digi/PulseDigi.cc | 12 ++++-------- 1 file changed, 4 insertions(+), 8 deletions(-) diff --git a/src/algorithms/digi/PulseDigi.cc b/src/algorithms/digi/PulseDigi.cc index cdecdf37d1..872f7acc74 100644 --- a/src/algorithms/digi/PulseDigi.cc +++ b/src/algorithms/digi/PulseDigi.cc @@ -16,12 +16,8 @@ class HGCROCRawSample { void addAmplitude(std::size_t idx, float amp) { meas_values[idx] = std::max(meas_values[idx], amp); } - void addUpCrossTime(std::size_t idx, double t) { - meas_types[idx] = 1; - meas_values[idx] = t; - } - void addDownCrossTime(std::size_t idx, double t) { - meas_types[idx] = 2; + void addCrossTime(uint8_t type, std::size_t idx, double t) { + meas_types[idx] = type; meas_values[idx] = t; } @@ -70,7 +66,7 @@ void PulseDigi::process(const PulseDigi::Input& input, const PulseDigi::Output& // Pulse crossed above the threshold. if (pulse.getAmplitude()[i] > m_cfg.threshold) { - raw_sample.addUpCrossTime(idx, get_crossing_time(m_cfg.threshold, t, pulse_dt, + raw_sample.addCrossTime(1, idx, get_crossing_time(m_cfg.threshold, t, pulse_dt, pulse.getAmplitude()[i], pulse.getAmplitude()[i - 1])); tot_progress = true; @@ -80,7 +76,7 @@ void PulseDigi::process(const PulseDigi::Input& input, const PulseDigi::Output& // Pulse crossed below the threshold. if (tot_progress && !tot_complete && pulse.getAmplitude()[i] < m_cfg.threshold) { - pulse_info.addDownCrossTime(idx, get_crossing_time(m_cfg.threshold, t, pulse_dt, + pulse_info.addCrossTime(2, idx, get_crossing_time(m_cfg.threshold, t, pulse_dt, pulse.getAmplitude()[i], pulse.getAmplitude()[i - 1])); tot_complete = true; From 4eb18fd4ebddf8f3c0ffc838154ab5f47922865e Mon Sep 17 00:00:00 2001 From: mhkim Date: Mon, 29 Sep 2025 11:26:04 -0500 Subject: [PATCH 09/16] Added const --- src/algorithms/digi/PulseDigi.cc | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/algorithms/digi/PulseDigi.cc b/src/algorithms/digi/PulseDigi.cc index 872f7acc74..341a551707 100644 --- a/src/algorithms/digi/PulseDigi.cc +++ b/src/algorithms/digi/PulseDigi.cc @@ -21,8 +21,8 @@ class HGCROCRawSample { meas_values[idx] = t; } - std::vector getMeasTypes() { return meas_types; } - std::vector getMeasValues() { return meas_values; } + const std::vector& getMeasTypes() const { return meas_types; } + const std::vector& getMeasValues() const { return meas_values; } private: // meas_type 0: Pulse crossed above the threshold. The crossing time is pushed back onto meas_values. From 5d336d44040b8f03826053236bfc14027ab859a5 Mon Sep 17 00:00:00 2001 From: mhkim Date: Mon, 6 Oct 2025 19:20:29 -0500 Subject: [PATCH 10/16] Amplitudes are measured using tick --- src/algorithms/digi/PulseDigi.cc | 117 +++++++++++++++++++------------ 1 file changed, 72 insertions(+), 45 deletions(-) diff --git a/src/algorithms/digi/PulseDigi.cc b/src/algorithms/digi/PulseDigi.cc index 341a551707..7dc2fc68e7 100644 --- a/src/algorithms/digi/PulseDigi.cc +++ b/src/algorithms/digi/PulseDigi.cc @@ -2,34 +2,52 @@ // Copyright (C) 2025 Minho Kim // +#include +#include +#include +#include +#include #include #include +#include +#include #include +#include +#include +#include +#include #include #include "PulseDigi.h" class HGCROCRawSample { public: - HGCROCRawSample(std::size_t n_meas) : meas_types(n_meas, 0), meas_values(n_meas, 0) {} + HGCROCRawSample(std::size_t n_meas) + : meas_types(n_meas, 0), amplitudes(n_meas, 0), TOAs(n_meas, 0), TOTs(n_meas, 0) {} - void addAmplitude(std::size_t idx, float amp) { - meas_values[idx] = std::max(meas_values[idx], amp); + void addAmplitude(std::size_t idx, double amp) { amplitudes[idx] = amp; } + void addTOA(std::size_t idx, double toa) { + meas_types[idx] = 1; + TOAs[idx] = toa; } - void addCrossTime(uint8_t type, std::size_t idx, double t) { - meas_types[idx] = type; - meas_values[idx] = t; + void addTOT(std::size_t idx, double tot) { + meas_types[idx] = 2; + TOTs[idx] = tot; } const std::vector& getMeasTypes() const { return meas_types; } - const std::vector& getMeasValues() const { return meas_values; } + const std::vector& getAmplitudes() const { return amplitudes; } + const std::vector& getTOAs() const { return TOAs; } + const std::vector& getTOTs() const { return TOTs; } private: - // meas_type 0: Pulse crossed above the threshold. The crossing time is pushed back onto meas_values. - // meas_type 1: Pulse crossed below the threshold. The crossing time is pushed back onto meas_values. - // meas_type 2: Pulse didn't cross the threshold. The maximum amplitude is pushed back onto meas_values. + // meas_type 1: Pulse crossed below the threshold. + // meas_type 2: Pulse didn't cross the threshold. + // meas_type 0: Neither of the above cases. std::vector meas_types; - std::vector meas_values; + std::vector amplitudes; + std::vector TOAs; + std::vector TOTs; }; namespace eicrecon { @@ -37,50 +55,59 @@ namespace eicrecon { void PulseDigi::init() {} void PulseDigi::process(const PulseDigi::Input& input, const PulseDigi::Output& output) const { - const auto [in_pulses] = input; - auto [out_digi_hits, out_digi_samples] = output; + const auto [in_pulses] = input; + auto [out_digi_hits] = output; for (const auto& pulse : *in_pulses) { - float pulse_t = pulse.getTime(); - float pulse_dt = pulse.getInterval(); - std::size_t n_amplitudes = pulse.getAmplitude().size(); + double pulse_t = pulse.getTime(); + double pulse_dt = pulse.getInterval(); + std::size_t n_amps = pulse.getAmplitude().size(); // Estimate the number of samples. - const std::size_t idx_begin = static_cast(std::floor(pulse_t / m_cfg.time_window)); - const std::size_t idx_end = static_cast( - std::floor((pulse_t + (n_amplitudes - 1) * pulse_dt) / m_cfg.time_window)); + const std::size_t timeIdx_begin = + static_cast(std::floor(pulse_t / m_cfg.time_window)); + const std::size_t timeIdx_end = static_cast( + std::floor((pulse_t + (n_amps - 1) * pulse_dt) / m_cfg.time_window)); - HGCROCRawSample raw_sample(idx_end - idx_begin + 1); + HGCROCRawSample raw_sample(timeIdx_end - timeIdx_begin + 1); - bool tot_progress = false; - bool tot_complete = false; + int sample_tick = std::llround(m_cfg.sample_period / pulse_dt); + int adc_counter = + std::llround((timeIdx_begin * m_cfg.sample_period + m_cfg.adc_phase - pulse_t) / pulse_dt); - for (std::size_t i = 0; i < n_amplitudes; i++) { + bool tot_progress = false; + bool tot_complete = false; + std::size_t toaIdx = 0; + + for (std::size_t i = 0; i < n_amps; i++) { double t = pulse_t + i * pulse_dt; - const std::size_t idx = - static_cast(std::floor(t / m_cfg.time_window)) - idx_begin; - - // Pulse didn't cross the threshold. - if (!tot_progess) { - raw_sample.addAmplitude(idx, pulse.getAmplitude()[i]); - - // Pulse crossed above the threshold. - if (pulse.getAmplitude()[i] > m_cfg.threshold) { - raw_sample.addCrossTime(1, idx, get_crossing_time(m_cfg.threshold, t, pulse_dt, - pulse.getAmplitude()[i], - pulse.getAmplitude()[i - 1])); - tot_progress = true; - tot_complete = false; - } + const std::size_t sampleIdx = + static_cast(std::floor(t / m_cfg.sample_period)) - timeIdx_begin; + + adc_counter++; + // Measure amplitudes + if (adc_counter == sample_tick) { + raw_sample.addAmplitude(sampleIdx, pulse.getAmplitude()[i]); + adc_counter = 0; + } + + // Measure crossing point for TOA + if (!tot_progress && pulse.getAmplitude()[i] > m_cfg.toa_thres) { + toaIdx = sampleIdx; + raw_sample.addTOA(sampleIdx, + get_crossing_time(m_cfg.threshold, t, pulse_dt, pulse.getAmplitude()[i], + pulse.getAmplitude()[i - 1])); + tot_progress = true; + tot_complete = false; } - // Pulse crossed below the threshold. + // Measure crossing point for TOT if (tot_progress && !tot_complete && pulse.getAmplitude()[i] < m_cfg.threshold) { - pulse_info.addCrossTime(2, idx, get_crossing_time(m_cfg.threshold, t, pulse_dt, - pulse.getAmplitude()[i], - pulse.getAmplitude()[i - 1])); - tot_complete = true; - tot_progess = false; + raw_sample.addTOT(toaIdx, + get_crossing_time(m_cfg.threshold, t, pulse_dt, pulse.getAmplitude()[i], + pulse.getAmplitude()[i - 1])); + tot_complete = true; + tot_progess = false; } } } @@ -89,7 +116,7 @@ void PulseDigi::process(const PulseDigi::Input& input, const PulseDigi::Output& double PulseDigi::get_crossing_time(double thres, double t1, double t2, double amp1, double amp2) const { double numerator = (thres - amp2) * (t2 - t1); - double denomenator = amp2 - amp1; + double denominator = amp2 - amp1; double added = t2; return (numerator / denominator) + added; } From 4d621f2c686a46e2d972db32d9cf920165ce4a40 Mon Sep 17 00:00:00 2001 From: mhkim Date: Mon, 6 Oct 2025 19:46:19 -0500 Subject: [PATCH 11/16] add -> set --- src/algorithms/digi/PulseDigi.cc | 12 ++++++------ 1 file changed, 6 insertions(+), 6 deletions(-) diff --git a/src/algorithms/digi/PulseDigi.cc b/src/algorithms/digi/PulseDigi.cc index 7dc2fc68e7..67a885dfd2 100644 --- a/src/algorithms/digi/PulseDigi.cc +++ b/src/algorithms/digi/PulseDigi.cc @@ -25,12 +25,12 @@ class HGCROCRawSample { HGCROCRawSample(std::size_t n_meas) : meas_types(n_meas, 0), amplitudes(n_meas, 0), TOAs(n_meas, 0), TOTs(n_meas, 0) {} - void addAmplitude(std::size_t idx, double amp) { amplitudes[idx] = amp; } - void addTOA(std::size_t idx, double toa) { + void setAmplitude(std::size_t idx, double amp) { amplitudes[idx] = amp; } + void setTOA(std::size_t idx, double toa) { meas_types[idx] = 1; TOAs[idx] = toa; } - void addTOT(std::size_t idx, double tot) { + void setTOT(std::size_t idx, double tot) { meas_types[idx] = 2; TOTs[idx] = tot; } @@ -87,14 +87,14 @@ void PulseDigi::process(const PulseDigi::Input& input, const PulseDigi::Output& adc_counter++; // Measure amplitudes if (adc_counter == sample_tick) { - raw_sample.addAmplitude(sampleIdx, pulse.getAmplitude()[i]); + raw_sample.setAmplitude(sampleIdx, pulse.getAmplitude()[i]); adc_counter = 0; } // Measure crossing point for TOA if (!tot_progress && pulse.getAmplitude()[i] > m_cfg.toa_thres) { toaIdx = sampleIdx; - raw_sample.addTOA(sampleIdx, + raw_sample.setTOA(sampleIdx, get_crossing_time(m_cfg.threshold, t, pulse_dt, pulse.getAmplitude()[i], pulse.getAmplitude()[i - 1])); tot_progress = true; @@ -103,7 +103,7 @@ void PulseDigi::process(const PulseDigi::Input& input, const PulseDigi::Output& // Measure crossing point for TOT if (tot_progress && !tot_complete && pulse.getAmplitude()[i] < m_cfg.threshold) { - raw_sample.addTOT(toaIdx, + raw_sample.setTOT(toaIdx, get_crossing_time(m_cfg.threshold, t, pulse_dt, pulse.getAmplitude()[i], pulse.getAmplitude()[i - 1])); tot_complete = true; From ecf9fd7e5c8fdbd08d841211a2c8e43ba1ae2d9f Mon Sep 17 00:00:00 2001 From: mhkim Date: Thu, 30 Oct 2025 15:59:38 -0500 Subject: [PATCH 12/16] Update for PR --- src/algorithms/digi/PulseDigi.cc | 210 ++++++++++++++----------- src/algorithms/digi/PulseDigi.h | 40 +++++ src/algorithms/digi/PulseDigiConfig.h | 27 ++++ src/factories/digi/PulseDigi_factory.h | 48 ++++++ 4 files changed, 230 insertions(+), 95 deletions(-) create mode 100644 src/algorithms/digi/PulseDigi.h create mode 100644 src/algorithms/digi/PulseDigiConfig.h create mode 100644 src/factories/digi/PulseDigi_factory.h diff --git a/src/algorithms/digi/PulseDigi.cc b/src/algorithms/digi/PulseDigi.cc index 67a885dfd2..ff280382b9 100644 --- a/src/algorithms/digi/PulseDigi.cc +++ b/src/algorithms/digi/PulseDigi.cc @@ -21,103 +21,123 @@ #include "PulseDigi.h" class HGCROCRawSample { -public: - HGCROCRawSample(std::size_t n_meas) - : meas_types(n_meas, 0), amplitudes(n_meas, 0), TOAs(n_meas, 0), TOTs(n_meas, 0) {} - - void setAmplitude(std::size_t idx, double amp) { amplitudes[idx] = amp; } - void setTOA(std::size_t idx, double toa) { - meas_types[idx] = 1; - TOAs[idx] = toa; - } - void setTOT(std::size_t idx, double tot) { - meas_types[idx] = 2; - TOTs[idx] = tot; - } - - const std::vector& getMeasTypes() const { return meas_types; } - const std::vector& getAmplitudes() const { return amplitudes; } - const std::vector& getTOAs() const { return TOAs; } - const std::vector& getTOTs() const { return TOTs; } - -private: - // meas_type 1: Pulse crossed below the threshold. - // meas_type 2: Pulse didn't cross the threshold. - // meas_type 0: Neither of the above cases. - std::vector meas_types; - std::vector amplitudes; - std::vector TOAs; - std::vector TOTs; + public: + struct RawEntry { + double adc{0}; + double toa{0}; + double tot{0}; + bool totProgress{false}; + bool totComplete{false}; + }; + + HGCROCRawSample(std::size_t n_samp) { rawEntries.resize(n_samp); } + + void setAmplitude(std::size_t idx, double amp) { rawEntries[idx].adc = amp; } + void setTOA(std::size_t idx, double toa) { rawEntries[idx].toa = toa; } + void setTOT(std::size_t idx, double cross_t) { rawEntries[idx].tot = cross_t - rawEntries[idx].toa; } + void setTotProgress(std::size_t idx) { rawEntries[idx].totProgress = true; } + void setTotComplete(std::size_t idx) { rawEntries[idx].totComplete = true; } + + const std::vector& getEntries() const { return rawEntries; } + + private: + std::vector rawEntries; }; namespace eicrecon { -void PulseDigi::init() {} - -void PulseDigi::process(const PulseDigi::Input& input, const PulseDigi::Output& output) const { - const auto [in_pulses] = input; - auto [out_digi_hits] = output; - - for (const auto& pulse : *in_pulses) { - double pulse_t = pulse.getTime(); - double pulse_dt = pulse.getInterval(); - std::size_t n_amps = pulse.getAmplitude().size(); - - // Estimate the number of samples. - const std::size_t timeIdx_begin = - static_cast(std::floor(pulse_t / m_cfg.time_window)); - const std::size_t timeIdx_end = static_cast( - std::floor((pulse_t + (n_amps - 1) * pulse_dt) / m_cfg.time_window)); - - HGCROCRawSample raw_sample(timeIdx_end - timeIdx_begin + 1); - - int sample_tick = std::llround(m_cfg.sample_period / pulse_dt); - int adc_counter = - std::llround((timeIdx_begin * m_cfg.sample_period + m_cfg.adc_phase - pulse_t) / pulse_dt); - - bool tot_progress = false; - bool tot_complete = false; - std::size_t toaIdx = 0; - - for (std::size_t i = 0; i < n_amps; i++) { - double t = pulse_t + i * pulse_dt; - const std::size_t sampleIdx = - static_cast(std::floor(t / m_cfg.sample_period)) - timeIdx_begin; - - adc_counter++; - // Measure amplitudes - if (adc_counter == sample_tick) { - raw_sample.setAmplitude(sampleIdx, pulse.getAmplitude()[i]); - adc_counter = 0; - } - - // Measure crossing point for TOA - if (!tot_progress && pulse.getAmplitude()[i] > m_cfg.toa_thres) { - toaIdx = sampleIdx; - raw_sample.setTOA(sampleIdx, - get_crossing_time(m_cfg.threshold, t, pulse_dt, pulse.getAmplitude()[i], - pulse.getAmplitude()[i - 1])); - tot_progress = true; - tot_complete = false; - } - - // Measure crossing point for TOT - if (tot_progress && !tot_complete && pulse.getAmplitude()[i] < m_cfg.threshold) { - raw_sample.setTOT(toaIdx, - get_crossing_time(m_cfg.threshold, t, pulse_dt, pulse.getAmplitude()[i], - pulse.getAmplitude()[i - 1])); - tot_complete = true; - tot_progess = false; - } - } - } -} // PulseDigi:process - -double PulseDigi::get_crossing_time(double thres, double t1, double t2, double amp1, - double amp2) const { - double numerator = (thres - amp2) * (t2 - t1); - double denominator = amp2 - amp1; - double added = t2; - return (numerator / denominator) + added; -} + void PulseDigi::init() {} + + void PulseDigi::process(const PulseDigi::Input& input, const PulseDigi::Output& output) const { + const auto [in_pulses] = input; + auto [out_digi_hits] = output; + + for (const auto& pulse : *in_pulses) { + double pulse_t = pulse.getTime(); + double pulse_dt = pulse.getInterval(); + std::size_t n_amps = pulse.getAmplitude().size(); + + // Estimate the number of samples. + const std::size_t timeIdx_begin = + static_cast(std::floor(pulse_t / m_cfg.time_window)); + const std::size_t timeIdx_end = static_cast( + std::floor((pulse_t + (n_amps - 1) * pulse_dt) / m_cfg.time_window)); + + HGCROCRawSample raw_sample(timeIdx_end - timeIdx_begin + 1); + + // For ADC, amplitude is measured with a fixed phase. + // This was reproduced by sample_tick and adc_counter as follows. + // Amplitude is measured whenever adc_counter reaches sample_tick. + int sample_tick = std::llround(m_cfg.time_window / pulse_dt); + int adc_counter = + std::llround((timeIdx_begin * m_cfg.time_window + m_cfg.adc_phase - pulse_t) / pulse_dt); + + bool tot_progress = false; + bool tot_complete = false; + std::size_t toaIdx = 0; + + for (std::size_t i = 0; i < n_amps; i++) { + double t = pulse_t + i * pulse_dt; + const std::size_t sampleIdx = + static_cast(std::floor(t / m_cfg.time_window)) - timeIdx_begin; + + adc_counter++; + + // Measure amplitudes for ADC + if (adc_counter == sample_tick) { + raw_sample.setAmplitude(sampleIdx, pulse.getAmplitude()[i]); + adc_counter = 0; + if(tot_progress) raw_sample.setTotProgress(sampleIdx); + } + + // Measure up-crossing time for TOA + if (!tot_progress && pulse.getAmplitude()[i] > m_cfg.toa_thres) { + toaIdx = sampleIdx; + raw_sample.setTOA(sampleIdx, + get_crossing_time(m_cfg.toa_thres, pulse_dt, t, pulse.getAmplitude()[i], + pulse.getAmplitude()[i - 1])); + tot_progress = true; + tot_complete = false; + raw_sample.setTotProgress(sampleIdx); + } + + // Measure down-crossing time for TOT + if (tot_progress && !tot_complete && pulse.getAmplitude()[i] < m_cfg.tot_thres) { + raw_sample.setTOT(toaIdx, + get_crossing_time(m_cfg.tot_thres, pulse_dt, t, pulse.getAmplitude()[i], + pulse.getAmplitude()[i - 1])); + tot_progress = false; + tot_complete = true; + raw_sample.setTotComplete(sampleIdx); + } + } + + // Fill HGCROCSamples and RawHGCROCHit + auto out_digi_hit = out_digi_hits->create(); + out_digi_hit.setCellID(pulse.getCellID()); + + const auto& entries = raw_sample.getEntries(); + + for (const auto& entry : entries) { + edm4eic::HGCROCSample sample; + auto adc = std::max(std::llround(entry.adc / m_cfg.dyRangeADC * m_cfg.capADC), 0LL); + sample.ADC = adc > m_cfg.capADC ? m_cfg.capADC : adc; + auto toa = std::max(std::llround(entry.toa / m_cfg.dyRangeTOA * m_cfg.capTOA), 0LL); + sample.timeOfArrival = toa > m_cfg.capTOA ? m_cfg.capTOA : toa; + auto tot = std::max(std::llround(entry.tot / m_cfg.dyRangeTOT * m_cfg.capTOT), 0LL); + sample.timeOverThreshold = tot > m_cfg.capTOT ? m_cfg.capTOT : tot; + sample.TOTInProgress = entry.totProgress; + sample.TOTComplete = entry.totComplete; + out_digi_hit.addToSamples(sample); + } + } + } // PulseDigi:process + + double PulseDigi::get_crossing_time(double thres, double dt, double t, double amp1, + double amp2) const { + double numerator = (thres - amp2) * dt; + double denominator = amp2 - amp1; + double added = t; + return (numerator / denominator) + added; + } } // namespace eicrecon diff --git a/src/algorithms/digi/PulseDigi.h b/src/algorithms/digi/PulseDigi.h new file mode 100644 index 0000000000..7f9ba3a6be --- /dev/null +++ b/src/algorithms/digi/PulseDigi.h @@ -0,0 +1,40 @@ +// SPDX-License-Identifier: LGPL-3.0-or-later +// Copyright (C) 2025 Minho Kim + +#pragma once + +#include +#include +#include +#include +#include +#include +#include +#include + +#include "algorithms/digi/PulseDigiConfig.h" +#include "algorithms/interfaces/WithPodConfig.h" + +namespace eicrecon { + +using PulseDigiAlgorithm = algorithms::Algorithm< + algorithms::Input, + algorithms::Output>; + +class PulseDigi : public PulseDigiAlgorithm, public WithPodConfig { + +public: + PulseDigi(std::string_view name) + : PulseDigiAlgorithm{name, {"InputPulses"}, {"OutputRawHGCROCHits"}, {}} {} + virtual void init() final; + void process(const Input&, const Output&) const; + +private: + edm4eic::HGCROCSample sample; + +private: + double get_crossing_time(double thres, double dt, double t, + double amp1, double amp2) const; +}; + +} // namespace eicrecon diff --git a/src/algorithms/digi/PulseDigiConfig.h b/src/algorithms/digi/PulseDigiConfig.h new file mode 100644 index 0000000000..00ad64530a --- /dev/null +++ b/src/algorithms/digi/PulseDigiConfig.h @@ -0,0 +1,27 @@ +// SPDX-License-Identifier: LGPL-3.0-or-later +// Copyright (C) 2025 Minho Kim + +#pragma once + +#include + +namespace eicrecon { + +struct PulseDigiConfig { + + // Variables for HGCROC measurement + double time_window{25 * edm4eic::unit::ns}; + double adc_phase{10 * edm4eic::unit::ns}; + double toa_thres{1}; + double tot_thres{1}; + + // Variables for digitization + unsigned int capADC{1024}; + double dyRangeADC{1}; + unsigned int capTOA{1024}; + double dyRangeTOA{1 * edm4eic::unit::ns}; + unsigned int capTOT{1024}; + double dyRangeTOT{1 * edm4eic::unit::ns}; +}; + +} // namespace eicrecon diff --git a/src/factories/digi/PulseDigi_factory.h b/src/factories/digi/PulseDigi_factory.h new file mode 100644 index 0000000000..4511d99fd0 --- /dev/null +++ b/src/factories/digi/PulseDigi_factory.h @@ -0,0 +1,48 @@ +// SPDX-License-Identifier: LGPL-3.0-or-later +// Copyright (C) 2025 Minho Kim + +#pragma once + +#include "algorithms/digi/PulseDigi.h" +#include "services/algorithms_init/AlgorithmsInit_service.h" +#include "extensions/jana/JOmniFactory.h" + +namespace eicrecon { + +class PulseDigi_factory : public JOmniFactory { + +public: + using AlgoT = eicrecon::PulseDigi; + +private: + std::unique_ptr m_algo; + + PodioInput m_pulse_input{this}; + PodioOutput m_digi_output{this}; + + ParameterRef m_time_window{this, "time_window", config().time_window}; + ParameterRef m_adc_phase{this, "adc_phase", config().adc_phase}; + ParameterRef m_toa_thres{this, "toa_thres", config().toa_thres}; + ParameterRef m_tot_thres{this, "tot_thres", config().tot_thres}; + ParameterRef m_capADC{this, "capADC", config().capADC}; + ParameterRef m_dyRangeADC{this, "dyRangeADC", config().dyRangeADC}; + ParameterRef m_capTOA{this, "capTOA", config().capTOA}; + ParameterRef m_dyRangeTOA{this, "dyRangeTOA", config().dyRangeTOA}; + ParameterRef m_capTOT{this, "capTOT", config().capTOT}; + ParameterRef m_dyRangeTOT{this, "dyRangeTOT", config().dyRangeTOT}; + + Service m_algorithmsInit{this}; + +public: + void Configure() { + m_algo = std::make_unique(GetPrefix()); + m_algo->level(static_cast(logger()->level())); + m_algo->applyConfig(config()); + m_algo->init(); + } + + void Process(int32_t /* run_number */, uint64_t /* event_number */) { + m_algo->process({m_pulse_input()}, {m_digi_output().get()}); + } +}; +} // namespace eicrecon From 594f668728377b5208e9a36325434d9fe6293182 Mon Sep 17 00:00:00 2001 From: mhkim Date: Thu, 30 Oct 2025 16:11:55 -0500 Subject: [PATCH 13/16] code formatter has been applied --- src/algorithms/digi/PulseDigi.cc | 233 +++++++++++++------------- src/algorithms/digi/PulseDigi.h | 15 +- src/algorithms/digi/PulseDigiConfig.h | 24 +-- 3 files changed, 139 insertions(+), 133 deletions(-) diff --git a/src/algorithms/digi/PulseDigi.cc b/src/algorithms/digi/PulseDigi.cc index ff280382b9..dcbc158db3 100644 --- a/src/algorithms/digi/PulseDigi.cc +++ b/src/algorithms/digi/PulseDigi.cc @@ -21,123 +21,126 @@ #include "PulseDigi.h" class HGCROCRawSample { - public: - struct RawEntry { - double adc{0}; - double toa{0}; - double tot{0}; - bool totProgress{false}; - bool totComplete{false}; - }; - - HGCROCRawSample(std::size_t n_samp) { rawEntries.resize(n_samp); } - - void setAmplitude(std::size_t idx, double amp) { rawEntries[idx].adc = amp; } - void setTOA(std::size_t idx, double toa) { rawEntries[idx].toa = toa; } - void setTOT(std::size_t idx, double cross_t) { rawEntries[idx].tot = cross_t - rawEntries[idx].toa; } - void setTotProgress(std::size_t idx) { rawEntries[idx].totProgress = true; } - void setTotComplete(std::size_t idx) { rawEntries[idx].totComplete = true; } - - const std::vector& getEntries() const { return rawEntries; } - - private: - std::vector rawEntries; +public: + struct RawEntry { + double adc{0}; + double toa{0}; + double tot{0}; + bool totProgress{false}; + bool totComplete{false}; + }; + + HGCROCRawSample(std::size_t n_samp) { rawEntries.resize(n_samp); } + + void setAmplitude(std::size_t idx, double amp) { rawEntries[idx].adc = amp; } + void setTOA(std::size_t idx, double toa) { rawEntries[idx].toa = toa; } + void setTOT(std::size_t idx, double cross_t) { + rawEntries[idx].tot = cross_t - rawEntries[idx].toa; + } + void setTotProgress(std::size_t idx) { rawEntries[idx].totProgress = true; } + void setTotComplete(std::size_t idx) { rawEntries[idx].totComplete = true; } + + const std::vector& getEntries() const { return rawEntries; } + +private: + std::vector rawEntries; }; namespace eicrecon { - void PulseDigi::init() {} - - void PulseDigi::process(const PulseDigi::Input& input, const PulseDigi::Output& output) const { - const auto [in_pulses] = input; - auto [out_digi_hits] = output; - - for (const auto& pulse : *in_pulses) { - double pulse_t = pulse.getTime(); - double pulse_dt = pulse.getInterval(); - std::size_t n_amps = pulse.getAmplitude().size(); - - // Estimate the number of samples. - const std::size_t timeIdx_begin = - static_cast(std::floor(pulse_t / m_cfg.time_window)); - const std::size_t timeIdx_end = static_cast( - std::floor((pulse_t + (n_amps - 1) * pulse_dt) / m_cfg.time_window)); - - HGCROCRawSample raw_sample(timeIdx_end - timeIdx_begin + 1); - - // For ADC, amplitude is measured with a fixed phase. - // This was reproduced by sample_tick and adc_counter as follows. - // Amplitude is measured whenever adc_counter reaches sample_tick. - int sample_tick = std::llround(m_cfg.time_window / pulse_dt); - int adc_counter = - std::llround((timeIdx_begin * m_cfg.time_window + m_cfg.adc_phase - pulse_t) / pulse_dt); - - bool tot_progress = false; - bool tot_complete = false; - std::size_t toaIdx = 0; - - for (std::size_t i = 0; i < n_amps; i++) { - double t = pulse_t + i * pulse_dt; - const std::size_t sampleIdx = - static_cast(std::floor(t / m_cfg.time_window)) - timeIdx_begin; - - adc_counter++; - - // Measure amplitudes for ADC - if (adc_counter == sample_tick) { - raw_sample.setAmplitude(sampleIdx, pulse.getAmplitude()[i]); - adc_counter = 0; - if(tot_progress) raw_sample.setTotProgress(sampleIdx); - } - - // Measure up-crossing time for TOA - if (!tot_progress && pulse.getAmplitude()[i] > m_cfg.toa_thres) { - toaIdx = sampleIdx; - raw_sample.setTOA(sampleIdx, - get_crossing_time(m_cfg.toa_thres, pulse_dt, t, pulse.getAmplitude()[i], - pulse.getAmplitude()[i - 1])); - tot_progress = true; - tot_complete = false; - raw_sample.setTotProgress(sampleIdx); - } - - // Measure down-crossing time for TOT - if (tot_progress && !tot_complete && pulse.getAmplitude()[i] < m_cfg.tot_thres) { - raw_sample.setTOT(toaIdx, - get_crossing_time(m_cfg.tot_thres, pulse_dt, t, pulse.getAmplitude()[i], - pulse.getAmplitude()[i - 1])); - tot_progress = false; - tot_complete = true; - raw_sample.setTotComplete(sampleIdx); - } - } - - // Fill HGCROCSamples and RawHGCROCHit - auto out_digi_hit = out_digi_hits->create(); - out_digi_hit.setCellID(pulse.getCellID()); - - const auto& entries = raw_sample.getEntries(); - - for (const auto& entry : entries) { - edm4eic::HGCROCSample sample; - auto adc = std::max(std::llround(entry.adc / m_cfg.dyRangeADC * m_cfg.capADC), 0LL); - sample.ADC = adc > m_cfg.capADC ? m_cfg.capADC : adc; - auto toa = std::max(std::llround(entry.toa / m_cfg.dyRangeTOA * m_cfg.capTOA), 0LL); - sample.timeOfArrival = toa > m_cfg.capTOA ? m_cfg.capTOA : toa; - auto tot = std::max(std::llround(entry.tot / m_cfg.dyRangeTOT * m_cfg.capTOT), 0LL); - sample.timeOverThreshold = tot > m_cfg.capTOT ? m_cfg.capTOT : tot; - sample.TOTInProgress = entry.totProgress; - sample.TOTComplete = entry.totComplete; - out_digi_hit.addToSamples(sample); - } - } - } // PulseDigi:process - - double PulseDigi::get_crossing_time(double thres, double dt, double t, double amp1, - double amp2) const { - double numerator = (thres - amp2) * dt; - double denominator = amp2 - amp1; - double added = t; - return (numerator / denominator) + added; - } +void PulseDigi::init() {} + +void PulseDigi::process(const PulseDigi::Input& input, const PulseDigi::Output& output) const { + const auto [in_pulses] = input; + auto [out_digi_hits] = output; + + for (const auto& pulse : *in_pulses) { + double pulse_t = pulse.getTime(); + double pulse_dt = pulse.getInterval(); + std::size_t n_amps = pulse.getAmplitude().size(); + + // Estimate the number of samples. + const std::size_t timeIdx_begin = + static_cast(std::floor(pulse_t / m_cfg.time_window)); + const std::size_t timeIdx_end = static_cast( + std::floor((pulse_t + (n_amps - 1) * pulse_dt) / m_cfg.time_window)); + + HGCROCRawSample raw_sample(timeIdx_end - timeIdx_begin + 1); + + // For ADC, amplitude is measured with a fixed phase. + // This was reproduced by sample_tick and adc_counter as follows. + // Amplitude is measured whenever adc_counter reaches sample_tick. + int sample_tick = std::llround(m_cfg.time_window / pulse_dt); + int adc_counter = + std::llround((timeIdx_begin * m_cfg.time_window + m_cfg.adc_phase - pulse_t) / pulse_dt); + + bool tot_progress = false; + bool tot_complete = false; + std::size_t toaIdx = 0; + + for (std::size_t i = 0; i < n_amps; i++) { + double t = pulse_t + i * pulse_dt; + const std::size_t sampleIdx = + static_cast(std::floor(t / m_cfg.time_window)) - timeIdx_begin; + + adc_counter++; + + // Measure amplitudes for ADC + if (adc_counter == sample_tick) { + raw_sample.setAmplitude(sampleIdx, pulse.getAmplitude()[i]); + adc_counter = 0; + if (tot_progress) + raw_sample.setTotProgress(sampleIdx); + } + + // Measure up-crossing time for TOA + if (!tot_progress && pulse.getAmplitude()[i] > m_cfg.toa_thres) { + toaIdx = sampleIdx; + raw_sample.setTOA(sampleIdx, + get_crossing_time(m_cfg.toa_thres, pulse_dt, t, pulse.getAmplitude()[i], + pulse.getAmplitude()[i - 1])); + tot_progress = true; + tot_complete = false; + raw_sample.setTotProgress(sampleIdx); + } + + // Measure down-crossing time for TOT + if (tot_progress && !tot_complete && pulse.getAmplitude()[i] < m_cfg.tot_thres) { + raw_sample.setTOT(toaIdx, + get_crossing_time(m_cfg.tot_thres, pulse_dt, t, pulse.getAmplitude()[i], + pulse.getAmplitude()[i - 1])); + tot_progress = false; + tot_complete = true; + raw_sample.setTotComplete(sampleIdx); + } + } + + // Fill HGCROCSamples and RawHGCROCHit + auto out_digi_hit = out_digi_hits->create(); + out_digi_hit.setCellID(pulse.getCellID()); + + const auto& entries = raw_sample.getEntries(); + + for (const auto& entry : entries) { + edm4eic::HGCROCSample sample; + auto adc = std::max(std::llround(entry.adc / m_cfg.dyRangeADC * m_cfg.capADC), 0LL); + sample.ADC = adc > m_cfg.capADC ? m_cfg.capADC : adc; + auto toa = std::max(std::llround(entry.toa / m_cfg.dyRangeTOA * m_cfg.capTOA), 0LL); + sample.timeOfArrival = toa > m_cfg.capTOA ? m_cfg.capTOA : toa; + auto tot = std::max(std::llround(entry.tot / m_cfg.dyRangeTOT * m_cfg.capTOT), 0LL); + sample.timeOverThreshold = tot > m_cfg.capTOT ? m_cfg.capTOT : tot; + sample.TOTInProgress = entry.totProgress; + sample.TOTComplete = entry.totComplete; + out_digi_hit.addToSamples(sample); + } + } +} // PulseDigi:process + +double PulseDigi::get_crossing_time(double thres, double dt, double t, double amp1, + double amp2) const { + double numerator = (thres - amp2) * dt; + double denominator = amp2 - amp1; + double added = t; + return (numerator / denominator) + added; +} } // namespace eicrecon diff --git a/src/algorithms/digi/PulseDigi.h b/src/algorithms/digi/PulseDigi.h index 7f9ba3a6be..cacd22c651 100644 --- a/src/algorithms/digi/PulseDigi.h +++ b/src/algorithms/digi/PulseDigi.h @@ -17,15 +17,19 @@ namespace eicrecon { -using PulseDigiAlgorithm = algorithms::Algorithm< - algorithms::Input, - algorithms::Output>; +using PulseDigiAlgorithm = + algorithms::Algorithm, + algorithms::Output>; class PulseDigi : public PulseDigiAlgorithm, public WithPodConfig { public: PulseDigi(std::string_view name) - : PulseDigiAlgorithm{name, {"InputPulses"}, {"OutputRawHGCROCHits"}, {}} {} + : PulseDigiAlgorithm{name, + {"InputPulses"}, + {"OutputDigiHits"}, + {"ADC, TOA, and TOT are measured referring to the" + "working principle of the HGCROC."}} {} virtual void init() final; void process(const Input&, const Output&) const; @@ -33,8 +37,7 @@ class PulseDigi : public PulseDigiAlgorithm, public WithPodConfig Date: Sat, 22 Nov 2025 19:14:12 -0600 Subject: [PATCH 14/16] Update src/factories/digi/PulseDigi_factory.h Co-authored-by: Dmitry Kalinkin --- src/factories/digi/PulseDigi_factory.h | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/src/factories/digi/PulseDigi_factory.h b/src/factories/digi/PulseDigi_factory.h index 4511d99fd0..c1f407c6cf 100644 --- a/src/factories/digi/PulseDigi_factory.h +++ b/src/factories/digi/PulseDigi_factory.h @@ -20,10 +20,10 @@ class PulseDigi_factory : public JOmniFactory m_pulse_input{this}; PodioOutput m_digi_output{this}; - ParameterRef m_time_window{this, "time_window", config().time_window}; - ParameterRef m_adc_phase{this, "adc_phase", config().adc_phase}; - ParameterRef m_toa_thres{this, "toa_thres", config().toa_thres}; - ParameterRef m_tot_thres{this, "tot_thres", config().tot_thres}; + ParameterRef m_time_window{this, "timeWindow", config().time_window}; + ParameterRef m_adc_phase{this, "adcPhase", config().adc_phase}; + ParameterRef m_toa_thres{this, "toaThres", config().toa_thres}; + ParameterRef m_tot_thres{this, "totThres", config().tot_thres}; ParameterRef m_capADC{this, "capADC", config().capADC}; ParameterRef m_dyRangeADC{this, "dyRangeADC", config().dyRangeADC}; ParameterRef m_capTOA{this, "capTOA", config().capTOA}; From 62efe5d0be8504cf48e904fa4c5df0d235d88ffe Mon Sep 17 00:00:00 2001 From: mhkim Date: Sat, 22 Nov 2025 19:38:09 -0600 Subject: [PATCH 15/16] Removed an unused variable --- src/algorithms/digi/PulseDigi.h | 3 --- 1 file changed, 3 deletions(-) diff --git a/src/algorithms/digi/PulseDigi.h b/src/algorithms/digi/PulseDigi.h index cacd22c651..462e952521 100644 --- a/src/algorithms/digi/PulseDigi.h +++ b/src/algorithms/digi/PulseDigi.h @@ -33,9 +33,6 @@ class PulseDigi : public PulseDigiAlgorithm, public WithPodConfig Date: Sat, 22 Nov 2025 22:41:21 -0600 Subject: [PATCH 16/16] Considered EDM4EIC_VERSION_MIN --- CMakeLists.txt | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/CMakeLists.txt b/CMakeLists.txt index e8f2fc140c..aa24f335c3 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -58,7 +58,7 @@ endif() set(Acts_VERSION_MIN 34.0.0) set(algorithms_VERSION_MIN 1.0.0) set(DD4hep_VERSION_MIN 1.21) -set(EDM4EIC_VERSION_MIN 8.0) +set(EDM4EIC_VERSION_MIN 8.3.0) set(EDM4HEP_VERSION_MIN 0.7.1) set(Eigen3_VERSION_MIN 3.3) set(FastJet_VERSION_MIN 3)