From f00a1b690530319f6211fc05f3d1da7aa1e3ed68 Mon Sep 17 00:00:00 2001 From: "copilot-swe-agent[bot]" <198982749+Copilot@users.noreply.github.com> Date: Sun, 17 May 2026 15:04:36 +0000 Subject: [PATCH 01/11] Enable DRICH generation-stage QE in npsim Agent-Logs-Url: https://github.com/eic/npsim/sessions/f21a9097-2f8b-4486-92e3-771e522e8a8f Co-authored-by: wdconinc <4656391+wdconinc@users.noreply.github.com> --- src/dd4pod/python/npsim.py | 21 +++++++++++++++++- .../OpticalPhotonEfficiencyStackingAction.h | 22 +++++++++++++++++++ 2 files changed, 42 insertions(+), 1 deletion(-) diff --git a/src/dd4pod/python/npsim.py b/src/dd4pod/python/npsim.py index 2c35e91..d902f31 100755 --- a/src/dd4pod/python/npsim.py +++ b/src/dd4pod/python/npsim.py @@ -69,8 +69,27 @@ def setupCerenkov(kernel): RUNNER.action.mapActions['PFRICH'] = 'Geant4OpticalTrackerAction' RUNNER.action.mapActions['DIRC'] = 'Geant4OpticalTrackerAction' - # Use the optical photon efficiency stacking action for hpDIRC + # Use the optical photon efficiency stacking action + drich_wavelengths = [ # [nm] + 315, 325, 340, 350, 370, 400, 450, 500, 550, 600, 650, 700, 750, 800, 850, 900, 1000 + ] + drich_efficiency = [ # [%] + 0.00, 0.04, 0.10, 0.20, 0.30, 0.35, 0.40, 0.38, 0.35, 0.27, 0.20, 0.15, 0.12, 0.08, 0.06, 0.04, 0.00 + ] RUNNER.action.stack = [ + *[ + { + "name": "OpticalPhotonEfficiencyStackingAction", + "parameter": { + "LambdaMin": f"{drich_wavelengths[0]}*nm", + "LambdaMax": f"{drich_wavelengths[-1]}*nm", + "Wavelengths": drich_wavelengths, + "LogicalVolume": f"DRICH_pss_sec{isec}", + "Efficiency": drich_efficiency, + } + } + for isec in range(6) + ], { "name": "OpticalPhotonEfficiencyStackingAction", "parameter": { diff --git a/src/plugins/include/npdet/OpticalPhotonEfficiencyStackingAction.h b/src/plugins/include/npdet/OpticalPhotonEfficiencyStackingAction.h index e902f4a..122bedb 100644 --- a/src/plugins/include/npdet/OpticalPhotonEfficiencyStackingAction.h +++ b/src/plugins/include/npdet/OpticalPhotonEfficiencyStackingAction.h @@ -17,6 +17,7 @@ #include "DDG4/Geant4StackingAction.h" #include "G4OpticalPhoton.hh" #include "G4Track.hh" +#include /// Namespace for the AIDA detector description toolkit namespace dd4hep { @@ -31,6 +32,7 @@ namespace dd4hep { : Geant4StackingAction(c, n) { declareProperty("LambdaMin", m_lambda_min); declareProperty("LambdaMax", m_lambda_max); + declareProperty("Wavelengths", m_wavelengths); declareProperty("Efficiency", m_efficiency); declareProperty("LogicalVolume", m_logical_volume); }; @@ -80,6 +82,25 @@ namespace dd4hep { // Single constant value over lambda range efficiency = m_efficiency.front(); + } else if (m_wavelengths.size() == m_efficiency.size() && m_wavelengths.size() > 1) { + // Linear interpolation on non-uniform wavelength grid + auto lambda_nm = lambda / CLHEP::nm; + if (lambda_nm <= m_wavelengths.front()) { + efficiency = m_efficiency.front(); + } else if (lambda_nm >= m_wavelengths.back()) { + efficiency = m_efficiency.back(); + } else { + auto upper = std::upper_bound(m_wavelengths.begin(), m_wavelengths.end(), lambda_nm); + auto i = std::distance(m_wavelengths.begin(), upper) - 1; + double a_lambda = m_wavelengths[i]; + double b_lambda = m_wavelengths[i+1]; + double t = (lambda_nm - a_lambda) / (b_lambda - a_lambda); + double a = m_efficiency[i]; + double b = m_efficiency[i+1]; + efficiency = a + t * (b - a); + printout(VERBOSE, name(), "a = %f, b = %f, t = %f", a, b, t); + printout(VERBOSE, name(), "efficiency %f", efficiency); + } } else { // Linear interpolation on lambda grid double lambda_step = (m_lambda_max - m_lambda_min) / (m_efficiency.size() - 1); @@ -120,6 +141,7 @@ namespace dd4hep { }; private: double m_lambda_min{0.}, m_lambda_max{0.}; + std::vector m_wavelengths; std::vector m_efficiency; std::string m_logical_volume; std::size_t m_total_photons{0}, m_killed_photons{0}; From d7fbedcfc36114f280f921d7bb99a0bf74a30ade Mon Sep 17 00:00:00 2001 From: "copilot-swe-agent[bot]" <198982749+Copilot@users.noreply.github.com> Date: Sun, 17 May 2026 15:29:37 +0000 Subject: [PATCH 02/11] Address PR feedback on QE interpolation and volume matching Agent-Logs-Url: https://github.com/eic/npsim/sessions/3adbb508-8d8e-415a-a9be-edc3274a6e86 Co-authored-by: wdconinc <4656391+wdconinc@users.noreply.github.com> --- src/dd4pod/python/npsim.py | 35 +++--- .../OpticalPhotonEfficiencyStackingAction.h | 100 ++++++++++-------- 2 files changed, 80 insertions(+), 55 deletions(-) diff --git a/src/dd4pod/python/npsim.py b/src/dd4pod/python/npsim.py index d902f31..fba97f8 100755 --- a/src/dd4pod/python/npsim.py +++ b/src/dd4pod/python/npsim.py @@ -70,26 +70,35 @@ def setupCerenkov(kernel): RUNNER.action.mapActions['DIRC'] = 'Geant4OpticalTrackerAction' # Use the optical photon efficiency stacking action - drich_wavelengths = [ # [nm] + drich_lambda_values = [ # [nm] 315, 325, 340, 350, 370, 400, 450, 500, 550, 600, 650, 700, 750, 800, 850, 900, 1000 ] drich_efficiency = [ # [%] 0.00, 0.04, 0.10, 0.20, 0.30, 0.35, 0.40, 0.38, 0.35, 0.27, 0.20, 0.15, 0.12, 0.08, 0.06, 0.04, 0.00 ] + pfrich_lambda_values = [ # [nm] + 315, 325, 340, 350, 370, 400, 450, 500, 550, 600, 650, 700, 750, 800, 850, 900, 1000 + ] + pfrich_efficiency = [ # [%] + 0.00, 0.04, 0.10, 0.20, 0.30, 0.35, 0.40, 0.38, 0.35, 0.27, 0.20, 0.15, 0.12, 0.08, 0.06, 0.04, 0.00 + ] RUNNER.action.stack = [ - *[ - { - "name": "OpticalPhotonEfficiencyStackingAction", - "parameter": { - "LambdaMin": f"{drich_wavelengths[0]}*nm", - "LambdaMax": f"{drich_wavelengths[-1]}*nm", - "Wavelengths": drich_wavelengths, - "LogicalVolume": f"DRICH_pss_sec{isec}", - "Efficiency": drich_efficiency, - } + { + "name": "OpticalPhotonEfficiencyStackingAction", + "parameter": { + "LambdaValues": drich_lambda_values, + "LogicalVolume": "DRICH_pss_sec[0-9]+", + "Efficiency": drich_efficiency, + } + }, + { + "name": "OpticalPhotonEfficiencyStackingAction", + "parameter": { + "LambdaValues": pfrich_lambda_values, + "LogicalVolume": "PFRICH-photocathode-[0-9]+", + "Efficiency": pfrich_efficiency, } - for isec in range(6) - ], + }, { "name": "OpticalPhotonEfficiencyStackingAction", "parameter": { diff --git a/src/plugins/include/npdet/OpticalPhotonEfficiencyStackingAction.h b/src/plugins/include/npdet/OpticalPhotonEfficiencyStackingAction.h index 122bedb..fd3bc26 100644 --- a/src/plugins/include/npdet/OpticalPhotonEfficiencyStackingAction.h +++ b/src/plugins/include/npdet/OpticalPhotonEfficiencyStackingAction.h @@ -18,6 +18,10 @@ #include "G4OpticalPhoton.hh" #include "G4Track.hh" #include +#include +#include +#include +#include /// Namespace for the AIDA detector description toolkit namespace dd4hep { @@ -32,7 +36,7 @@ namespace dd4hep { : Geant4StackingAction(c, n) { declareProperty("LambdaMin", m_lambda_min); declareProperty("LambdaMax", m_lambda_max); - declareProperty("Wavelengths", m_wavelengths); + declareProperty("LambdaValues", m_lambda_values); declareProperty("Efficiency", m_efficiency); declareProperty("LogicalVolume", m_logical_volume); }; @@ -54,6 +58,7 @@ namespace dd4hep { virtual void prepare(G4StackManager*) override { }; /// Return TrackClassification with enum G4ClassificationOfNewTrack or NoTrackClassification virtual TrackClassification classifyNewTrack(G4StackManager*, const G4Track* aTrack) override { + std::call_once(m_init_once, [this]() { this->initialize(); }); // Only apply to optical photons if (aTrack->GetDefinition() == G4OpticalPhoton::OpticalPhotonDefinition()) { auto* pv = aTrack->GetVolume(); @@ -61,57 +66,25 @@ namespace dd4hep { auto* lv = pv->GetLogicalVolume(); printout(VERBOSE, name(), "photon in pv %s lv %s", pv->GetName().c_str(), lv->GetName().c_str()); - // Only apply to specified logical volume - if (lv->GetName() == m_logical_volume) { + // Only apply to specified logical volume regex + if (std::regex_match(lv->GetName().c_str(), m_logical_volume_regex)) { double mom = aTrack->GetMomentum().mag(); double lambda = CLHEP::hbarc * CLHEP::twopi / mom; + double lambda_nm = lambda / CLHEP::nm; printout(VERBOSE, name(), "with mom = %f eV, lambda = %f nm", mom / CLHEP::eV, lambda / CLHEP::nm); m_total_photons++; - if (m_lambda_min < lambda && lambda < m_lambda_max) { - double efficiency{0.}; + auto lambda_min = !m_interp_lambda_values.empty() ? m_interp_lambda_values.front() : m_lambda_min / CLHEP::nm; + auto lambda_max = !m_interp_lambda_values.empty() ? m_interp_lambda_values.back() : m_lambda_max / CLHEP::nm; + if (lambda_min < lambda_nm && lambda_nm < lambda_max) { + double efficiency = interpolate(lambda_nm); if (m_efficiency.size() == 0) { // No efficiency specified, assume zero efficiency = 0.; // which means kill ++m_killed_photons; return TrackClassification(fKill); - - } else if (m_efficiency.size() == 1) { - // Single constant value over lambda range - efficiency = m_efficiency.front(); - - } else if (m_wavelengths.size() == m_efficiency.size() && m_wavelengths.size() > 1) { - // Linear interpolation on non-uniform wavelength grid - auto lambda_nm = lambda / CLHEP::nm; - if (lambda_nm <= m_wavelengths.front()) { - efficiency = m_efficiency.front(); - } else if (lambda_nm >= m_wavelengths.back()) { - efficiency = m_efficiency.back(); - } else { - auto upper = std::upper_bound(m_wavelengths.begin(), m_wavelengths.end(), lambda_nm); - auto i = std::distance(m_wavelengths.begin(), upper) - 1; - double a_lambda = m_wavelengths[i]; - double b_lambda = m_wavelengths[i+1]; - double t = (lambda_nm - a_lambda) / (b_lambda - a_lambda); - double a = m_efficiency[i]; - double b = m_efficiency[i+1]; - efficiency = a + t * (b - a); - printout(VERBOSE, name(), "a = %f, b = %f, t = %f", a, b, t); - printout(VERBOSE, name(), "efficiency %f", efficiency); - } - } else { - // Linear interpolation on lambda grid - double lambda_step = (m_lambda_max - m_lambda_min) / (m_efficiency.size() - 1); - double div = (lambda - m_lambda_min) / lambda_step; - auto i = std::llround(std::floor(div)); - double t = div - i; - double a = m_efficiency[i]; - double b = m_efficiency[i+1]; - efficiency = a + t * (b - a); - printout(VERBOSE, name(), "a = %f, b = %f, t = %f", a, b, t); - printout(VERBOSE, name(), "efficiency %f", efficiency); } // Edge cases @@ -131,7 +104,7 @@ namespace dd4hep { return TrackClassification(fKill); } } else { - printout(VERBOSE, name(), "outside lambda range [%f,%f] nm", m_lambda_min / CLHEP::nm, m_lambda_max / CLHEP::nm); + printout(VERBOSE, name(), "outside lambda range [%f,%f] nm", lambda_min, lambda_max); } } else { printout(VERBOSE, name(), "not in volume %s", m_logical_volume.c_str()); @@ -140,10 +113,53 @@ namespace dd4hep { return TrackClassification(); }; private: + void initialize() { + try { + m_logical_volume_regex = std::regex(m_logical_volume); + } catch (const std::regex_error&) { + throw std::runtime_error("Invalid LogicalVolume regex"); + } + if (m_efficiency.size() > 1) { + if (m_lambda_values.size() == m_efficiency.size()) { + m_interp_lambda_values = m_lambda_values; + } else { + auto lambda_min_nm = m_lambda_min / CLHEP::nm; + auto lambda_max_nm = m_lambda_max / CLHEP::nm; + m_interp_lambda_values.reserve(m_efficiency.size()); + auto lambda_step = (lambda_max_nm - lambda_min_nm) / (m_efficiency.size() - 1); + for (std::size_t i = 0; i < m_efficiency.size(); ++i) { + m_interp_lambda_values.push_back(lambda_min_nm + i * lambda_step); + } + } + } else if (m_lambda_values.size() > 1) { + m_interp_lambda_values = m_lambda_values; + } + } + double interpolate(double lambda_nm) const { + if (m_efficiency.size() == 1) { + return m_efficiency.front(); + } + if (m_efficiency.size() < 2 || m_interp_lambda_values.size() < 2) { + return 0.0; + } + auto upper = std::upper_bound(m_interp_lambda_values.begin(), m_interp_lambda_values.end(), lambda_nm); + auto i = std::distance(m_interp_lambda_values.begin(), upper) - 1; + double a_lambda = m_interp_lambda_values[i]; + double b_lambda = m_interp_lambda_values[i+1]; + double t = (lambda_nm - a_lambda) / (b_lambda - a_lambda); + double a = m_efficiency[i]; + double b = m_efficiency[i+1]; + printout(VERBOSE, name(), "a = %f, b = %f, t = %f", a, b, t); + printout(VERBOSE, name(), "efficiency %f", a + t * (b - a)); + return a + t * (b - a); + } double m_lambda_min{0.}, m_lambda_max{0.}; - std::vector m_wavelengths; + std::vector m_lambda_values; + std::vector m_interp_lambda_values; std::vector m_efficiency; std::string m_logical_volume; + std::regex m_logical_volume_regex{}; + std::once_flag m_init_once{}; std::size_t m_total_photons{0}, m_killed_photons{0}; }; } // End namespace sim From f795cf65187e647b72ab22027cae0fe6731bfc6c Mon Sep 17 00:00:00 2001 From: "copilot-swe-agent[bot]" <198982749+Copilot@users.noreply.github.com> Date: Sun, 17 May 2026 16:50:49 +0000 Subject: [PATCH 03/11] Require LambdaValues in dd4hep units Agent-Logs-Url: https://github.com/eic/npsim/sessions/1a0bce3d-a9bd-4818-9cc7-9f6d4ea19785 Co-authored-by: wdconinc <4656391+wdconinc@users.noreply.github.com> --- src/dd4pod/python/npsim.py | 9 ++++--- .../OpticalPhotonEfficiencyStackingAction.h | 25 +++++++++---------- 2 files changed, 17 insertions(+), 17 deletions(-) diff --git a/src/dd4pod/python/npsim.py b/src/dd4pod/python/npsim.py index fba97f8..43c0fda 100755 --- a/src/dd4pod/python/npsim.py +++ b/src/dd4pod/python/npsim.py @@ -70,15 +70,16 @@ def setupCerenkov(kernel): RUNNER.action.mapActions['DIRC'] = 'Geant4OpticalTrackerAction' # Use the optical photon efficiency stacking action - drich_lambda_values = [ # [nm] + nm = 1e-6 # dd4hep length unit + drich_lambda_values = [l * nm for l in [ # [nm] 315, 325, 340, 350, 370, 400, 450, 500, 550, 600, 650, 700, 750, 800, 850, 900, 1000 - ] + ]] drich_efficiency = [ # [%] 0.00, 0.04, 0.10, 0.20, 0.30, 0.35, 0.40, 0.38, 0.35, 0.27, 0.20, 0.15, 0.12, 0.08, 0.06, 0.04, 0.00 ] - pfrich_lambda_values = [ # [nm] + pfrich_lambda_values = [l * nm for l in [ # [nm] 315, 325, 340, 350, 370, 400, 450, 500, 550, 600, 650, 700, 750, 800, 850, 900, 1000 - ] + ]] pfrich_efficiency = [ # [%] 0.00, 0.04, 0.10, 0.20, 0.30, 0.35, 0.40, 0.38, 0.35, 0.27, 0.20, 0.15, 0.12, 0.08, 0.06, 0.04, 0.00 ] diff --git a/src/plugins/include/npdet/OpticalPhotonEfficiencyStackingAction.h b/src/plugins/include/npdet/OpticalPhotonEfficiencyStackingAction.h index fd3bc26..0bb3c7b 100644 --- a/src/plugins/include/npdet/OpticalPhotonEfficiencyStackingAction.h +++ b/src/plugins/include/npdet/OpticalPhotonEfficiencyStackingAction.h @@ -70,15 +70,14 @@ namespace dd4hep { if (std::regex_match(lv->GetName().c_str(), m_logical_volume_regex)) { double mom = aTrack->GetMomentum().mag(); double lambda = CLHEP::hbarc * CLHEP::twopi / mom; - double lambda_nm = lambda / CLHEP::nm; printout(VERBOSE, name(), "with mom = %f eV, lambda = %f nm", mom / CLHEP::eV, lambda / CLHEP::nm); m_total_photons++; - auto lambda_min = !m_interp_lambda_values.empty() ? m_interp_lambda_values.front() : m_lambda_min / CLHEP::nm; - auto lambda_max = !m_interp_lambda_values.empty() ? m_interp_lambda_values.back() : m_lambda_max / CLHEP::nm; - if (lambda_min < lambda_nm && lambda_nm < lambda_max) { - double efficiency = interpolate(lambda_nm); + auto lambda_min = !m_interp_lambda_values.empty() ? m_interp_lambda_values.front() : m_lambda_min; + auto lambda_max = !m_interp_lambda_values.empty() ? m_interp_lambda_values.back() : m_lambda_max; + if (lambda_min < lambda && lambda < lambda_max) { + double efficiency = interpolate(lambda); if (m_efficiency.size() == 0) { // No efficiency specified, assume zero efficiency = 0.; @@ -104,7 +103,7 @@ namespace dd4hep { return TrackClassification(fKill); } } else { - printout(VERBOSE, name(), "outside lambda range [%f,%f] nm", lambda_min, lambda_max); + printout(VERBOSE, name(), "outside lambda range [%f,%f] nm", lambda_min / CLHEP::nm, lambda_max / CLHEP::nm); } } else { printout(VERBOSE, name(), "not in volume %s", m_logical_volume.c_str()); @@ -123,30 +122,30 @@ namespace dd4hep { if (m_lambda_values.size() == m_efficiency.size()) { m_interp_lambda_values = m_lambda_values; } else { - auto lambda_min_nm = m_lambda_min / CLHEP::nm; - auto lambda_max_nm = m_lambda_max / CLHEP::nm; + auto lambda_min = m_lambda_min; + auto lambda_max = m_lambda_max; m_interp_lambda_values.reserve(m_efficiency.size()); - auto lambda_step = (lambda_max_nm - lambda_min_nm) / (m_efficiency.size() - 1); + auto lambda_step = (lambda_max - lambda_min) / (m_efficiency.size() - 1); for (std::size_t i = 0; i < m_efficiency.size(); ++i) { - m_interp_lambda_values.push_back(lambda_min_nm + i * lambda_step); + m_interp_lambda_values.push_back(lambda_min + i * lambda_step); } } } else if (m_lambda_values.size() > 1) { m_interp_lambda_values = m_lambda_values; } } - double interpolate(double lambda_nm) const { + double interpolate(double lambda) const { if (m_efficiency.size() == 1) { return m_efficiency.front(); } if (m_efficiency.size() < 2 || m_interp_lambda_values.size() < 2) { return 0.0; } - auto upper = std::upper_bound(m_interp_lambda_values.begin(), m_interp_lambda_values.end(), lambda_nm); + auto upper = std::upper_bound(m_interp_lambda_values.begin(), m_interp_lambda_values.end(), lambda); auto i = std::distance(m_interp_lambda_values.begin(), upper) - 1; double a_lambda = m_interp_lambda_values[i]; double b_lambda = m_interp_lambda_values[i+1]; - double t = (lambda_nm - a_lambda) / (b_lambda - a_lambda); + double t = (lambda - a_lambda) / (b_lambda - a_lambda); double a = m_efficiency[i]; double b = m_efficiency[i+1]; printout(VERBOSE, name(), "a = %f, b = %f, t = %f", a, b, t); From 8eadf60ec97912a026e5fb0c6d9023f1d100a267 Mon Sep 17 00:00:00 2001 From: "copilot-swe-agent[bot]" <198982749+Copilot@users.noreply.github.com> Date: Sun, 17 May 2026 16:53:52 +0000 Subject: [PATCH 04/11] Normalize RICH efficiency tables to fractions Agent-Logs-Url: https://github.com/eic/npsim/sessions/ba4e2433-66ae-4dbf-8578-1179788a6269 Co-authored-by: wdconinc <4656391+wdconinc@users.noreply.github.com> --- src/dd4pod/python/npsim.py | 12 ++++++------ 1 file changed, 6 insertions(+), 6 deletions(-) diff --git a/src/dd4pod/python/npsim.py b/src/dd4pod/python/npsim.py index 43c0fda..4122760 100755 --- a/src/dd4pod/python/npsim.py +++ b/src/dd4pod/python/npsim.py @@ -74,15 +74,15 @@ def setupCerenkov(kernel): drich_lambda_values = [l * nm for l in [ # [nm] 315, 325, 340, 350, 370, 400, 450, 500, 550, 600, 650, 700, 750, 800, 850, 900, 1000 ]] - drich_efficiency = [ # [%] - 0.00, 0.04, 0.10, 0.20, 0.30, 0.35, 0.40, 0.38, 0.35, 0.27, 0.20, 0.15, 0.12, 0.08, 0.06, 0.04, 0.00 - ] + drich_efficiency = [e/100. for e in [ # [%] + 0.00, 4.00, 10.0, 20.0, 30.0, 35.0, 40.0, 38.0, 35.0, 27.0, 20.0, 15.0, 12.0, 8.00, 6.00, 4.00, 0.00 + ]] pfrich_lambda_values = [l * nm for l in [ # [nm] 315, 325, 340, 350, 370, 400, 450, 500, 550, 600, 650, 700, 750, 800, 850, 900, 1000 ]] - pfrich_efficiency = [ # [%] - 0.00, 0.04, 0.10, 0.20, 0.30, 0.35, 0.40, 0.38, 0.35, 0.27, 0.20, 0.15, 0.12, 0.08, 0.06, 0.04, 0.00 - ] + pfrich_efficiency = [e/100. for e in [ # [%] + 0.00, 4.00, 10.0, 20.0, 30.0, 35.0, 40.0, 38.0, 35.0, 27.0, 20.0, 15.0, 12.0, 8.00, 6.00, 4.00, 0.00 + ]] RUNNER.action.stack = [ { "name": "OpticalPhotonEfficiencyStackingAction", From 3134a1ed898e7520db6b724efc01924dd6340260 Mon Sep 17 00:00:00 2001 From: "copilot-swe-agent[bot]" <198982749+Copilot@users.noreply.github.com> Date: Sun, 17 May 2026 17:05:15 +0000 Subject: [PATCH 05/11] Use std::lerp in stacking interpolation Agent-Logs-Url: https://github.com/eic/npsim/sessions/45003b90-638f-41c1-82a2-67a63baea600 Co-authored-by: wdconinc <4656391+wdconinc@users.noreply.github.com> --- .../npdet/OpticalPhotonEfficiencyStackingAction.h | 12 +++++++++--- 1 file changed, 9 insertions(+), 3 deletions(-) diff --git a/src/plugins/include/npdet/OpticalPhotonEfficiencyStackingAction.h b/src/plugins/include/npdet/OpticalPhotonEfficiencyStackingAction.h index 0bb3c7b..7947922 100644 --- a/src/plugins/include/npdet/OpticalPhotonEfficiencyStackingAction.h +++ b/src/plugins/include/npdet/OpticalPhotonEfficiencyStackingAction.h @@ -18,6 +18,7 @@ #include "G4OpticalPhoton.hh" #include "G4Track.hh" #include +#include #include #include #include @@ -134,7 +135,7 @@ namespace dd4hep { m_interp_lambda_values = m_lambda_values; } } - double interpolate(double lambda) const { + inline double interpolate(double lambda) const { if (m_efficiency.size() == 1) { return m_efficiency.front(); } @@ -149,8 +150,13 @@ namespace dd4hep { double a = m_efficiency[i]; double b = m_efficiency[i+1]; printout(VERBOSE, name(), "a = %f, b = %f, t = %f", a, b, t); - printout(VERBOSE, name(), "efficiency %f", a + t * (b - a)); - return a + t * (b - a); +#if defined(__cpp_lib_interpolate) && __cpp_lib_interpolate >= 201902L + double efficiency = std::lerp(a, b, t); +#else + double efficiency = a + t * (b - a); +#endif + printout(VERBOSE, name(), "efficiency %f", efficiency); + return efficiency; } double m_lambda_min{0.}, m_lambda_max{0.}; std::vector m_lambda_values; From 2ed77b230a540b011db3fd7b664937db0ded924f Mon Sep 17 00:00:00 2001 From: "copilot-swe-agent[bot]" <198982749+Copilot@users.noreply.github.com> Date: Sun, 17 May 2026 17:10:45 +0000 Subject: [PATCH 06/11] Remove redundant zero-efficiency branch Agent-Logs-Url: https://github.com/eic/npsim/sessions/cadf4be9-5ff6-4152-b088-8895b02fa8a8 Co-authored-by: wdconinc <4656391+wdconinc@users.noreply.github.com> --- .../include/npdet/OpticalPhotonEfficiencyStackingAction.h | 8 -------- 1 file changed, 8 deletions(-) diff --git a/src/plugins/include/npdet/OpticalPhotonEfficiencyStackingAction.h b/src/plugins/include/npdet/OpticalPhotonEfficiencyStackingAction.h index 7947922..eb894e7 100644 --- a/src/plugins/include/npdet/OpticalPhotonEfficiencyStackingAction.h +++ b/src/plugins/include/npdet/OpticalPhotonEfficiencyStackingAction.h @@ -79,14 +79,6 @@ namespace dd4hep { auto lambda_max = !m_interp_lambda_values.empty() ? m_interp_lambda_values.back() : m_lambda_max; if (lambda_min < lambda && lambda < lambda_max) { double efficiency = interpolate(lambda); - if (m_efficiency.size() == 0) { - // No efficiency specified, assume zero - efficiency = 0.; - // which means kill - ++m_killed_photons; - return TrackClassification(fKill); - } - // Edge cases if (efficiency == 0.0) { ++m_killed_photons; From 99c993596a69d9a8a6c715d7732ba814234e2dcd Mon Sep 17 00:00:00 2001 From: Wouter Deconinck Date: Sun, 17 May 2026 14:47:12 -0500 Subject: [PATCH 07/11] Kill photons outside QE range and promote suppression summary to WARNING - Photons outside the specified lambda range now get killed (efficiency=0 outside the QE table boundaries) instead of being kept. This is physically correct: a photocathode with 0% QE at its wavelength boundaries has 0% QE outside that range too. Previously, UV Cherenkov photons (below 315 nm) were passing unfiltered, dramatically reducing the expected hit count reduction. - Change the destructor summary printout from DEBUG to WARNING so the 'Suppressed N of M photons in lv ...' message appears in CI runs (which use -v WARNING). This provides in-CI verification that the stacking action is active and killing photons at the expected rate. Co-authored-by: Copilot <223556219+Copilot@users.noreply.github.com> --- .../npdet/OpticalPhotonEfficiencyStackingAction.h | 12 ++++-------- 1 file changed, 4 insertions(+), 8 deletions(-) diff --git a/src/plugins/include/npdet/OpticalPhotonEfficiencyStackingAction.h b/src/plugins/include/npdet/OpticalPhotonEfficiencyStackingAction.h index eb894e7..d947171 100644 --- a/src/plugins/include/npdet/OpticalPhotonEfficiencyStackingAction.h +++ b/src/plugins/include/npdet/OpticalPhotonEfficiencyStackingAction.h @@ -43,15 +43,8 @@ namespace dd4hep { }; /// Default destructor virtual ~OpticalPhotonEfficiencyStackingAction() { - printout(DEBUG, name(), "Suppressed %d of %d photons in lv %s", + printout(WARNING, name(), "Suppressed %zu of %zu photons in lv %s", m_killed_photons, m_total_photons, m_logical_volume.c_str()); - printout(DEBUG, name(), "lambda range: [%f,%f] nm", - m_lambda_min / CLHEP::nm, m_lambda_max / CLHEP::nm); - std::ostringstream oss_efficiency; - std::copy(m_efficiency.begin(), m_efficiency.end(), - std::ostream_iterator(oss_efficiency, " ")); - std::string str_efficiency = oss_efficiency.str(); - printout(DEBUG, name(), "efficiency: %s", str_efficiency.c_str()); }; /// New-stage callback virtual void newStage(G4StackManager*) override { }; @@ -96,7 +89,10 @@ namespace dd4hep { return TrackClassification(fKill); } } else { + // Outside the QE-specified wavelength range: efficiency is 0 printout(VERBOSE, name(), "outside lambda range [%f,%f] nm", lambda_min / CLHEP::nm, lambda_max / CLHEP::nm); + ++m_killed_photons; + return TrackClassification(fKill); } } else { printout(VERBOSE, name(), "not in volume %s", m_logical_volume.c_str()); From 60136b57996699403bbf1e474f22d48a406588c0 Mon Sep 17 00:00:00 2001 From: Wouter Deconinck Date: Sun, 17 May 2026 16:30:25 -0500 Subject: [PATCH 08/11] Apply QE stacking action to Cherenkov radiator volumes, not sensor volumes The stacking action should filter photons when they are created in the radiator volumes, not when they arrive at the sensor surface. This is the correct place to apply generation-stage QE: the photon wavelength is set at creation, and any photon that would not pass the sensor QE can be killed immediately, saving simulation time. DRICH radiator volumes (both gas ring Cherenkov and proximity-focus aerogel): DRICH_gas and DRICH_aerogel -> regex DRICH_(gas|aerogel) PFRICH radiator volumes (gas volume and individual aerogel tiles): PFRICH_GasVol and PFRICH-aerogel-{ir}-{ia} -> regex PFRICH(_GasVol|-aerogel-[0-9]+-[0-9]+) Co-authored-by: Copilot <223556219+Copilot@users.noreply.github.com> --- src/dd4pod/python/npsim.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/dd4pod/python/npsim.py b/src/dd4pod/python/npsim.py index 4122760..ac3adb9 100755 --- a/src/dd4pod/python/npsim.py +++ b/src/dd4pod/python/npsim.py @@ -88,7 +88,7 @@ def setupCerenkov(kernel): "name": "OpticalPhotonEfficiencyStackingAction", "parameter": { "LambdaValues": drich_lambda_values, - "LogicalVolume": "DRICH_pss_sec[0-9]+", + "LogicalVolume": "DRICH_(gas|aerogel)", "Efficiency": drich_efficiency, } }, @@ -96,7 +96,7 @@ def setupCerenkov(kernel): "name": "OpticalPhotonEfficiencyStackingAction", "parameter": { "LambdaValues": pfrich_lambda_values, - "LogicalVolume": "PFRICH-photocathode-[0-9]+", + "LogicalVolume": "PFRICH(_GasVol|-aerogel-[0-9]+-[0-9]+)", "Efficiency": pfrich_efficiency, } }, From 68039cfcb49c8b948e8db1b634d94e7c9d9ceec8 Mon Sep 17 00:00:00 2001 From: Wouter Deconinck Date: Sun, 17 May 2026 16:31:48 -0500 Subject: [PATCH 09/11] Print percentage of suppressed photons in stacking action summary Co-authored-by: Copilot <223556219+Copilot@users.noreply.github.com> --- .../include/npdet/OpticalPhotonEfficiencyStackingAction.h | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/src/plugins/include/npdet/OpticalPhotonEfficiencyStackingAction.h b/src/plugins/include/npdet/OpticalPhotonEfficiencyStackingAction.h index d947171..4a57f6d 100644 --- a/src/plugins/include/npdet/OpticalPhotonEfficiencyStackingAction.h +++ b/src/plugins/include/npdet/OpticalPhotonEfficiencyStackingAction.h @@ -43,8 +43,9 @@ namespace dd4hep { }; /// Default destructor virtual ~OpticalPhotonEfficiencyStackingAction() { - printout(WARNING, name(), "Suppressed %zu of %zu photons in lv %s", - m_killed_photons, m_total_photons, m_logical_volume.c_str()); + double pct = m_total_photons > 0 ? 100.0 * m_killed_photons / m_total_photons : 0.0; + printout(WARNING, name(), "Suppressed %zu of %zu photons (%.1f%%) in lv %s", + m_killed_photons, m_total_photons, pct, m_logical_volume.c_str()); }; /// New-stage callback virtual void newStage(G4StackManager*) override { }; From cd0d00a68e67926185c519bb339f9b1b80cbdd2a Mon Sep 17 00:00:00 2001 From: Wouter Deconinck Date: Sun, 17 May 2026 17:13:36 -0500 Subject: [PATCH 10/11] fix: import DDG4 and use g4units for nm definition --- src/dd4pod/python/npsim.py | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/src/dd4pod/python/npsim.py b/src/dd4pod/python/npsim.py index ac3adb9..c5609a0 100755 --- a/src/dd4pod/python/npsim.py +++ b/src/dd4pod/python/npsim.py @@ -12,7 +12,8 @@ import sys from DDSim.DD4hepSimulation import DD4hepSimulation - +import DDG4 +from g4units import keV if __name__ == "__main__": logging.basicConfig(format='%(name)-16s %(levelname)s %(message)s', level=logging.INFO, stream=sys.stdout) @@ -70,7 +71,7 @@ def setupCerenkov(kernel): RUNNER.action.mapActions['DIRC'] = 'Geant4OpticalTrackerAction' # Use the optical photon efficiency stacking action - nm = 1e-6 # dd4hep length unit + nm = 1e-6 * mm # dd4hep length unit drich_lambda_values = [l * nm for l in [ # [nm] 315, 325, 340, 350, 370, 400, 450, 500, 550, 600, 650, 700, 750, 800, 850, 900, 1000 ]] From 6379d10566c0751cb3f0352121b6479447c8876a Mon Sep 17 00:00:00 2001 From: Wouter Deconinck Date: Sun, 17 May 2026 17:17:10 -0500 Subject: [PATCH 11/11] fix: from g4units import mm not keV Updated import statements for consistency and clarity. --- src/dd4pod/python/npsim.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/dd4pod/python/npsim.py b/src/dd4pod/python/npsim.py index c5609a0..c3e8d77 100755 --- a/src/dd4pod/python/npsim.py +++ b/src/dd4pod/python/npsim.py @@ -12,8 +12,8 @@ import sys from DDSim.DD4hepSimulation import DD4hepSimulation -import DDG4 -from g4units import keV + +from g4units import mm if __name__ == "__main__": logging.basicConfig(format='%(name)-16s %(levelname)s %(message)s', level=logging.INFO, stream=sys.stdout)