From cae28b9a05b42e0dcb30ff3e27c39d83124f6503 Mon Sep 17 00:00:00 2001 From: ShannonS00 Date: Tue, 10 Mar 2026 16:49:56 +0100 Subject: [PATCH 1/5] Implements a new `PixelResponseNonUniformity` effect in effeects/electronics/noise.py that models PRNU by multiplying each detector pixel by a per pixel gain factor drawn from N(1, prnu_std). The gain map is generated once per detector (keyed by detector ID) and reused across exposures for consistency. `prnu_std` can be a scalar float or a dict keyed by detector ID, allowing different RMS values per detector type (e.g. H2RG vs GeoSnap). --- scopesim/effects/electronic/__init__.py | 3 +- scopesim/effects/electronic/noise.py | 56 +++++++++++++++++++++++++ 2 files changed, 58 insertions(+), 1 deletion(-) diff --git a/scopesim/effects/electronic/__init__.py b/scopesim/effects/electronic/__init__.py index 9ae7c4ff..5d10e0d4 100644 --- a/scopesim/effects/electronic/__init__.py +++ b/scopesim/effects/electronic/__init__.py @@ -9,6 +9,7 @@ - PoorMansHxRGReadoutNoise - simple readout noise for HAWAII detectors - BasicReadoutNoise - readout noise - ShotNoise - realisation of Poissonian photon noise +- PixelResponseNonUniformity - per-pixel gain variation (PRNU) - DarkCurrent - add dark current - LinearityCurve - apply detector (non-)linearity and saturation - ReferencePixelBorder @@ -23,7 +24,7 @@ from .electrons import LinearityCurve, ADConversion, InterPixelCapacitance from .noise import (Bias, PoorMansHxRGReadoutNoise, BasicReadoutNoise, - ShotNoise, DarkCurrent) + ShotNoise, DarkCurrent, PixelResponseNonUniformity) from .exposure import AutoExposure, ExposureIntegration, ExposureOutput from .pixels import ReferencePixelBorder, BinnedImage, UnequalBinnedImage from .dmps import DetectorModePropertiesSetter diff --git a/scopesim/effects/electronic/noise.py b/scopesim/effects/electronic/noise.py index d2ce522a..e356cb60 100644 --- a/scopesim/effects/electronic/noise.py +++ b/scopesim/effects/electronic/noise.py @@ -134,6 +134,62 @@ def plot_hist(self, det, **kwargs): ax.hist(dtcr.data.flatten()) +class PixelResponseNonUniformity(Effect): + """Per-pixel gain variation (PRNU). + + Multiplies each pixel by a gain factor drawn from N(1, prnu_std). + The map is generated once per detector and reused across exposures. + ``prnu_std`` may be a float or a dict keyed by detector ID. + """ + + required_keys: ClassVar[set] = {"prnu_std"} + z_order: ClassVar[tuple[int, ...]] = (805,) + + def __init__(self, **kwargs): + super().__init__(**kwargs) + self.meta["random_seed"] = "!SIM.random.seed" + self.meta.update(kwargs) + check_keys(self.meta, self.required_keys, action="error") + self._gain_maps = {} # keyed by dtcr_id + + def apply_to(self, obj, **kwargs): + if not isinstance(obj, Detector): + return obj + + random_seed = from_currsys(self.meta["random_seed"], self.cmds) + id_key = real_colname("id", obj.meta) + dtcr_id = obj.meta[id_key] if id_key is not None else None + + prnu_std_meta = from_currsys(self.meta["prnu_std"], self.cmds) + if isinstance(prnu_std_meta, dict): + prnu_std = float(from_currsys(prnu_std_meta[dtcr_id], self.cmds)) + elif isinstance(prnu_std_meta, (int, float)): + prnu_std = float(prnu_std_meta) + else: + raise ValueError( + ".meta['prnu_std'] must be a float " + f"or a dict keyed by detector ID, got {type(prnu_std_meta)}" + ) + + shape = obj._hdu.data.shape + if dtcr_id not in self._gain_maps or self._gain_maps[dtcr_id].shape != shape: + seed = (int(abs(hash((int(random_seed), str(dtcr_id))))) % (2**31) + if random_seed is not None else None) + self._gain_maps[dtcr_id] = np.random.default_rng(seed).normal( + loc=1.0, scale=prnu_std, size=shape) + + obj._hdu.data = obj._hdu.data * self._gain_maps[dtcr_id] + return obj + + def plot(self, det_id=None): + if not self._gain_maps: + raise RuntimeError("No gain map yet — run a simulation first.") + key = det_id if det_id in self._gain_maps else next(iter(self._gain_maps)) + fig, ax = figure_factory() + ax.imshow(self._gain_maps[key], origin="lower", aspect="auto") + return fig + + class ShotNoise(Effect): z_order: ClassVar[tuple[int, ...]] = (820,) From acf6daed92625b6f131690a3a875c2ad9439023c Mon Sep 17 00:00:00 2001 From: ShannonS00 Date: Tue, 10 Mar 2026 17:13:57 +0100 Subject: [PATCH 2/5] Added default deviation of 0.01 --- scopesim/effects/electronic/noise.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/scopesim/effects/electronic/noise.py b/scopesim/effects/electronic/noise.py index e356cb60..8a390c09 100644 --- a/scopesim/effects/electronic/noise.py +++ b/scopesim/effects/electronic/noise.py @@ -142,14 +142,14 @@ class PixelResponseNonUniformity(Effect): ``prnu_std`` may be a float or a dict keyed by detector ID. """ - required_keys: ClassVar[set] = {"prnu_std"} + required_keys: ClassVar[set] = set() z_order: ClassVar[tuple[int, ...]] = (805,) def __init__(self, **kwargs): super().__init__(**kwargs) self.meta["random_seed"] = "!SIM.random.seed" + self.meta.setdefault("prnu_std", 0.01) self.meta.update(kwargs) - check_keys(self.meta, self.required_keys, action="error") self._gain_maps = {} # keyed by dtcr_id def apply_to(self, obj, **kwargs): From b44dde5c64a37ab5c649d6ff62bb8f3a785cff70 Mon Sep 17 00:00:00 2001 From: ShannonS00 Date: Wed, 11 Mar 2026 16:10:29 +0100 Subject: [PATCH 3/5] Addressed PR review: Simplify PixelResponseNonUniformity based on PR review - Replace SeedSequence-based per-detector seeding with direct random_seed - Rename seed parameter from random_seed to prnu_seed - Remove default prnu_std (now must be defined in instrument package) - Add 5 unit tests covering multiplicative behaviour, fixed pattern, dict mode - Update all METIS detector YAMLs to use prnu_seed (in IRDB) --- scopesim/effects/electronic/noise.py | 14 +++--- scopesim/tests/tests_effects/test_prnu.py | 54 +++++++++++++++++++++++ 2 files changed, 61 insertions(+), 7 deletions(-) create mode 100644 scopesim/tests/tests_effects/test_prnu.py diff --git a/scopesim/effects/electronic/noise.py b/scopesim/effects/electronic/noise.py index 8a390c09..73f566dd 100644 --- a/scopesim/effects/electronic/noise.py +++ b/scopesim/effects/electronic/noise.py @@ -147,8 +147,6 @@ class PixelResponseNonUniformity(Effect): def __init__(self, **kwargs): super().__init__(**kwargs) - self.meta["random_seed"] = "!SIM.random.seed" - self.meta.setdefault("prnu_std", 0.01) self.meta.update(kwargs) self._gain_maps = {} # keyed by dtcr_id @@ -156,7 +154,7 @@ def apply_to(self, obj, **kwargs): if not isinstance(obj, Detector): return obj - random_seed = from_currsys(self.meta["random_seed"], self.cmds) + random_seed = from_currsys(self.meta.get("prnu_seed"), self.cmds) id_key = real_colname("id", obj.meta) dtcr_id = obj.meta[id_key] if id_key is not None else None @@ -173,9 +171,7 @@ def apply_to(self, obj, **kwargs): shape = obj._hdu.data.shape if dtcr_id not in self._gain_maps or self._gain_maps[dtcr_id].shape != shape: - seed = (int(abs(hash((int(random_seed), str(dtcr_id))))) % (2**31) - if random_seed is not None else None) - self._gain_maps[dtcr_id] = np.random.default_rng(seed).normal( + self._gain_maps[dtcr_id] = np.random.default_rng(random_seed).normal( loc=1.0, scale=prnu_std, size=shape) obj._hdu.data = obj._hdu.data * self._gain_maps[dtcr_id] @@ -185,8 +181,12 @@ def plot(self, det_id=None): if not self._gain_maps: raise RuntimeError("No gain map yet — run a simulation first.") key = det_id if det_id in self._gain_maps else next(iter(self._gain_maps)) + gain_map = self._gain_maps[key] + dev = np.max(np.abs(gain_map - 1.0)) fig, ax = figure_factory() - ax.imshow(self._gain_maps[key], origin="lower", aspect="auto") + im = ax.imshow(gain_map, origin="lower", aspect="auto", + vmin=1 - dev, vmax=1 + dev) + fig.colorbar(im, ax=ax, label="per-pixel gain") return fig diff --git a/scopesim/tests/tests_effects/test_prnu.py b/scopesim/tests/tests_effects/test_prnu.py new file mode 100644 index 00000000..292ca245 --- /dev/null +++ b/scopesim/tests/tests_effects/test_prnu.py @@ -0,0 +1,54 @@ +import numpy as np + +from scopesim.detector import Detector +from scopesim.effects.electronic import PixelResponseNonUniformity +from scopesim.optics.image_plane_utils import header_from_list_of_xy + + +def make_detector(value=1000, size=10): + hdr = header_from_list_of_xy([-size/2, size/2], [-size/2, size/2], 1, "D") + dtcr = Detector(hdr) + dtcr._hdu.data[:] = value + return dtcr + + +class TestApplyTo: + def test_output_std_matches_prnu_std(self): + """Relative std of output should be close to prnu_std.""" + prnu_std = 0.05 + dtcr = make_detector(value=1000, size=100) + PixelResponseNonUniformity(prnu_std=prnu_std, prnu_seed=42).apply_to(dtcr) + rel_std = dtcr._hdu.data.std() / dtcr._hdu.data.mean() + assert abs(rel_std - prnu_std) < 0.01 + + def test_gain_map_is_reused(self): + """Same map should be applied on repeated calls (fixed pattern).""" + dtcr1 = make_detector() + dtcr2 = make_detector() + prnu = PixelResponseNonUniformity(prnu_std=0.01, prnu_seed=42) + prnu.apply_to(dtcr1) + prnu.apply_to(dtcr2) + np.testing.assert_array_equal(dtcr1._hdu.data, dtcr2._hdu.data) + + def test_dict_prnu_std(self): + """Dict mode: different amplitude per detector ID.""" + hdr = header_from_list_of_xy([-5, 5], [-5, 5], 1, "D") + dtcr = Detector(hdr) + dtcr.meta["id"] = "H2RG" + dtcr._hdu.data[:] = 1000 + prnu = PixelResponseNonUniformity( + prnu_std={"H2RG": 0.005, "GeoSnap": 0.020}, prnu_seed=42) + prnu.apply_to(dtcr) + assert dtcr._hdu.data.std() > 0 + + def test_multiplicative_zero_signal(self): + """Zero signal should remain zero — confirms multiplicative (not additive).""" + dtcr = make_detector(value=0) + PixelResponseNonUniformity(prnu_std=0.01, prnu_seed=42).apply_to(dtcr) + assert dtcr._hdu.data.sum() == 0 + + def test_non_detector_passthrough(self): + """Non-Detector objects should be returned unchanged.""" + prnu = PixelResponseNonUniformity(prnu_std=0.01) + result = prnu.apply_to("not a detector") + assert result == "not a detector" From ec11007bf22f13b47142e7b5f1574f6661959b4e Mon Sep 17 00:00:00 2001 From: ShannonS00 Date: Thu, 12 Mar 2026 17:15:38 +0100 Subject: [PATCH 4/5] Fix type error in PixelResponseNonUniformity for invalid prnu_std input and enhance test coverage for error handling and plotting --- scopesim/effects/electronic/noise.py | 2 +- scopesim/tests/tests_effects/test_prnu.py | 17 +++++++++++++++++ 2 files changed, 18 insertions(+), 1 deletion(-) diff --git a/scopesim/effects/electronic/noise.py b/scopesim/effects/electronic/noise.py index 73f566dd..90d95e14 100644 --- a/scopesim/effects/electronic/noise.py +++ b/scopesim/effects/electronic/noise.py @@ -164,7 +164,7 @@ def apply_to(self, obj, **kwargs): elif isinstance(prnu_std_meta, (int, float)): prnu_std = float(prnu_std_meta) else: - raise ValueError( + raise TypeError( ".meta['prnu_std'] must be a float " f"or a dict keyed by detector ID, got {type(prnu_std_meta)}" ) diff --git a/scopesim/tests/tests_effects/test_prnu.py b/scopesim/tests/tests_effects/test_prnu.py index 292ca245..1c5e1eaa 100644 --- a/scopesim/tests/tests_effects/test_prnu.py +++ b/scopesim/tests/tests_effects/test_prnu.py @@ -1,4 +1,6 @@ +import matplotlib.pyplot as plt import numpy as np +import pytest from scopesim.detector import Detector from scopesim.effects.electronic import PixelResponseNonUniformity @@ -52,3 +54,18 @@ def test_non_detector_passthrough(self): prnu = PixelResponseNonUniformity(prnu_std=0.01) result = prnu.apply_to("not a detector") assert result == "not a detector" + + def test_invalid_prnu_std_raises(self): + prnu = PixelResponseNonUniformity(prnu_std="not a number nor a dict") + with pytest.raises(TypeError): + prnu.apply_to(make_detector()) + + def test_plot_raises_before_simulation(self): + prnu = PixelResponseNonUniformity(prnu_std=0.01) + with pytest.raises(RuntimeError): + prnu.plot() + + def test_plot_returns_figure(self): + prnu = PixelResponseNonUniformity(prnu_std=0.01, prnu_seed=42) + prnu.apply_to(make_detector()) + assert isinstance(prnu.plot(), plt.Figure) From be4928b1c4fb5e2deceaad5df668b06095e157b2 Mon Sep 17 00:00:00 2001 From: ShannonS00 Date: Mon, 16 Mar 2026 09:41:57 +0100 Subject: [PATCH 5/5] Updated PixelResponseNonUniformity docstring; updated test to conditionally assert plot creation based on PLOTS flag. (like in test_PSF.py) --- scopesim/effects/electronic/noise.py | 31 ++++++++++++++++++++--- scopesim/tests/tests_effects/test_prnu.py | 6 +++-- 2 files changed, 31 insertions(+), 6 deletions(-) diff --git a/scopesim/effects/electronic/noise.py b/scopesim/effects/electronic/noise.py index 90d95e14..51db1e17 100644 --- a/scopesim/effects/electronic/noise.py +++ b/scopesim/effects/electronic/noise.py @@ -135,11 +135,34 @@ def plot_hist(self, det, **kwargs): class PixelResponseNonUniformity(Effect): - """Per-pixel gain variation (PRNU). + """Pixel Response Non-Uniformity (PRNU). + + Models the fixed pattern of per-pixel gain variations across the detector + arising from manufacturing differences in quantum efficiency. Each pixel is + multiplied by a gain factor drawn from N(1, ``prnu_std``) keyed by detector ID. + The gain map is generated once per detector on first use and reused identically + across all subsequent exposures. + + Parameters + ---------- + prnu_std : float or dict + Standard deviation of the per-pixel gain distribution. + + prnu_seed : int, fixed + + include: "!DET.include_prnu" + + Example + -------- + + - name: prnu + description: Pixel response non-uniformity + class: PixelResponseNonUniformity + kwargs: + prnu_std: 0.001 + prnu_seed: 42 + include: "!DET.include_prnu" - Multiplies each pixel by a gain factor drawn from N(1, prnu_std). - The map is generated once per detector and reused across exposures. - ``prnu_std`` may be a float or a dict keyed by detector ID. """ required_keys: ClassVar[set] = set() diff --git a/scopesim/tests/tests_effects/test_prnu.py b/scopesim/tests/tests_effects/test_prnu.py index 1c5e1eaa..2c7bbcbd 100644 --- a/scopesim/tests/tests_effects/test_prnu.py +++ b/scopesim/tests/tests_effects/test_prnu.py @@ -1,4 +1,4 @@ -import matplotlib.pyplot as plt +from matplotlib import pyplot as plt import numpy as np import pytest @@ -6,6 +6,7 @@ from scopesim.effects.electronic import PixelResponseNonUniformity from scopesim.optics.image_plane_utils import header_from_list_of_xy +PLOTS = False def make_detector(value=1000, size=10): hdr = header_from_list_of_xy([-size/2, size/2], [-size/2, size/2], 1, "D") @@ -68,4 +69,5 @@ def test_plot_raises_before_simulation(self): def test_plot_returns_figure(self): prnu = PixelResponseNonUniformity(prnu_std=0.01, prnu_seed=42) prnu.apply_to(make_detector()) - assert isinstance(prnu.plot(), plt.Figure) + if PLOTS: + assert isinstance(prnu.plot(), plt.Figure)