The target is a 4-column × 3-row grid of calibrated test zones. +Each zone has a known waveform type and spatial frequency. +All zones share the same Michelson contrast and DC sky level.
+ +| Row | Col 0 | Col 1 | Col 2 | Col 3 |
|---|---|---|---|---|
| 0 | +Sine H f=0.04 c/px |
+ Sine H f=0.08 c/px |
+ Sine H f=0.16 c/px |
+ Sine H f=0.32 c/px |
| 1 | +Square H f=0.04 c/px |
+ Square H f=0.08 c/px |
+ Square H f=0.16 c/px |
+ Square H f=0.32 c/px |
| 2 | +Sine V f=0.08 c/px |
+ Sine 45° f=0.08 c/px |
+ Siemens star (radial chirp) |
+ Slant edge (~5° from vertical) |
The four column frequencies are chosen to align with the app's 4-level +wavelet decomposition band centres:
+| Column | Frequency (c/px) | Wavelength (px) | +Wavelet level | Approx scale |
|---|---|---|---|---|
| 0 | 0.04 | 25 px | 4 (coarsest) | ~16 px |
| 1 | 0.08 | 12 px | 3 | ~8 px |
| 2 | 0.16 | 6 px | 2 | ~4 px |
| 3 | 0.32 | 3 px | 1 (finest / near Nyquist) | ~2 px |
Each zone ramps Michelson contrast linearly from top to bottom: +the top edge has Contrast min and the bottom edge has Contrast max. +Within any horizontal strip across a row, all four columns share the same contrast — +so you can directly compare how different frequencies or waveform types respond at +identical contrast levels. The ramp restarts at each zone row, so row 0 (sine), +row 1 (square), and row 2 (special) all independently cover the full min→max range.
+ +| {_val(va)} | " f"{_val(vb)} | ") + # Weber fraction contrast table + wc_a = sm.get("weber_contrast_a", {}) + wc_b = sm.get("weber_contrast_b", {}) + wc_rows = "" + for ks in sorted(set(list(wc_a.keys()) + list(wc_b.keys()))): + va = wc_a.get(ks) + vb = wc_b.get(ks) + ca, cb = _better_worse_class(va, vb) + wc_rows += (f"|
| {ks} px | " + f"{_val(va, '.4f')} | " + f"{_val(vb, '.4f')} |
| Kernel size | {ra.label} | {rb.label} |
|---|
| Scale | {ra.label} (NC score) | {rb.label} (NC score) | Ratio A/B |
|---|
Side-by-side local σ maps at each kernel size (shared colour scale), each followed by its cross-section profile. The difference map (right) highlights where one filter preserves more local variation.
+Noise-normalised (× noise floor) — shared colour scale is a fair A/B comparison.
+{figs_for("nrm_std_")}| Scale | {ra.label} (NC score) | {rb.label} (NC score) | Ratio A/B |
|---|
|LoG| maps at σ = 1.5, 3, and 6 px (shared colour scale per row), each followed by its cross-section profile. A filter preserving more fine detail shows brighter, more defined boundaries at small σ.
+Noise-normalised (× noise floor) — shared colour scale is a fair A/B comparison.
+{figs_for("nrm_log_")}| Scale | {ra.label} (NC score) | {rb.label} (NC score) | Ratio A/B |
|---|
Reconstructed detail images at levels 2 and 3 (shared colour scale, diverging colourmap), each followed by its cross-section profile. -The difference panel (right) shows where fine structure differs between the two filters.
""" +The difference panel (right) shows where fine structure differs between the two filters. +Noise-normalised (× noise floor) — shared colour scale is a fair A/B comparison.
+{figs_for("nrm_wavelet_")} + +Formula: c = ΔL / L, ' + 'where ΔL = Imax − Imin (local range in the K × K kernel) ' + 'and L = median(kernel) (local background luminance). ' + 'Output is unbounded ≥ 0; a value of 1.0 means the local range equals the background ' + 'luminance, 2.0 means twice, and so on.
' + 'Why median for L: The median represents the background luminance ' + 'the feature is seen against, matching Weber\'s Law. Using the mean would inflate L ' + 'toward bright filaments within the kernel, artificially suppressing contrast values. ' + 'Median is also robust to hot pixels and residual star halos in starless images.
' + 'Scalar metric (table below): 99th percentile of the Weber map ' + 'within the analysis region. Near-zero median pixels over dark sky produce very large ' + 'Weber values; the 99th percentile captures peak structural contrast while ignoring ' + 'isolated dark-floor artefacts.
' + 'Kernel sizes: Small kernels (3 px) respond to sub-pixel-scale ' + 'transitions. Medium kernels (5 px) capture fine filaments. ' + 'Large kernels (9 px) reflect coarser structural contrast such as knots and shell edges.
' + 'Wide dynamic range: Weber contrast is intentionally unbounded. ' + 'Maps are displayed with a square-root colour scale (PowerNorm γ = 0.5) to compress ' + 'the bright end. Selecting a star-free nebula ROI avoids dark-sky pixels that drive ' + 'Weber values very high. Maps use a starless image when one is available.
', + title="Weber fraction contrast")} +| Kernel size | {ra.label} (99th pct c) | {rb.label} (99th pct c) |
|---|
| Scale | {ra.label} (NC score) | {rb.label} (NC score) | Ratio A/B |
|---|
Per-pixel Weber fraction contrast maps (c = ΔL / L, square-root colour +scale, viridis). Brighter regions have higher Weber contrast — the local intensity range +is large relative to the local median luminance. The difference panel (A−B) shows where +one image achieves greater relative contrast. High values over dark-sky regions are +expected; use a nebula ROI for meaningful filter comparison.
+Noise-normalised (× noise floor) — shared colour scale is a fair A/B comparison.
+{figs_for("nrm_weber_")} + +| Scale | {ra.label} (NC score) | {rb.label} (NC score) | Ratio A/B |
|---|
Gradient magnitude maps at σ = 1.5, 3, and 6 px (shared colour scale per row), +each followed by its cross-section profile. +A filter preserving sharper boundaries shows brighter, more defined gradient response.
+Noise-normalised (× noise floor) — shared colour scale is a fair A/B comparison.
+{figs_for("nrm_gradient_")} + +Ratio A/B for every noise-corrected method plotted against its +approximate spatial scale. Scale units differ by method (see 8b–8f methodology +boxes) — use this chart to spot which spatial-scale regime favours which filter, +not to compare absolute ratio values across methods.
""" # ── Section 9: Signal-to-Noise Ratio ───────────────────────────────────── def _plot_snr_pair(self, disp_a, label_a, disp_b, label_b, - vmin: float, vmax: float) -> plt.Figure: + vmin: float, vmax: float) -> plt.Figure | None: """Render both SNR maps side-by-side with a shared plasma color scale.""" panels = [(d, lbl) for d, lbl in [(disp_a, label_a), (disp_b, label_b)] if d is not None] + if not panels: + return None fig, axes = plt.subplots(1, len(panels), figsize=(7 * len(panels), 5)) if len(panels) == 1: axes = [axes] @@ -4148,6 +4469,8 @@ def row(metric, val_a, val_b, fmt=".3f", cr_b = sm_b.get("contrast_ratios_b", {}) if sm_b else {} snr_wav_a = sm_a.get("wavelet_snr_a", {}) snr_wav_b = sm_b.get("wavelet_snr_b", {}) if sm_b else {} + wc_a_s = sm_a.get("weber_contrast_a", {}) + wc_b_s = sm_b.get("weber_contrast_b", {}) if sm_b else {} snr_ma = ra.snr_metrics or {} snr_mb = rb.snr_metrics or {} @@ -4205,6 +4528,7 @@ def row_pm(metric, val_a, val_b, spread_a, spread_b, fmt=".3f", row("Std contrast ratio (15px)", cr_a.get(15), cr_b.get(15)), row("Wavelet SNR level 2", snr_wav_a.get(2), snr_wav_b.get(2)), row("Wavelet SNR level 3", snr_wav_a.get(3), snr_wav_b.get(3)), + row("Weber contrast 99th pct (5px)", wc_a_s.get(5), wc_b_s.get(5), fmt=".4f"), *([row( "Global SNR — starless (σ) ★", (snr_ma.get("starless") or {}).get("snr_global"), diff --git a/synthetic/target_generator.py b/synthetic/target_generator.py new file mode 100644 index 0000000..acb5777 --- /dev/null +++ b/synthetic/target_generator.py @@ -0,0 +1,282 @@ +"""Spatial test target generator for spatial-detail metric calibration. + +Generates a 4×3 grid of calibrated test patterns with known spatial frequencies +and waveform types. Load the clean reference into Image A and the degraded +companion into Image B to calibrate the app's spatial-detail metrics against +known inputs. + +The four column frequencies are deliberately chosen to align with the four +wavelet decomposition band centres (levels 4→1, coarse→fine): + f=0.04 c/px → λ=25 px → wavelet level 4 (~16 px scale) + f=0.08 c/px → λ=12 px → wavelet level 3 (~8 px scale) + f=0.16 c/px → λ=6 px → wavelet level 2 (~4 px scale) + f=0.32 c/px → λ=3 px → wavelet level 1 (~2 px / near Nyquist) +""" + +from __future__ import annotations + +from pathlib import Path + +import numpy as np +from scipy.signal import fftconvolve + +from synthetic.generator import _moffat_2d, _fwhm_to_alpha, _siemens_star + + +# (row, col) → (waveform_type, frequency_cycles_per_pixel) +# Row 0: horizontal sine sweep — cleanest MTF measurement +# Row 1: horizontal square-wave sweep — tests harmonic/ringing behavior +# Row 2: direction variants + special targets +_ZONE_SPEC: dict[tuple[int, int], tuple[str, float]] = { + (0, 0): ("sine_h", 0.04), + (0, 1): ("sine_h", 0.08), + (0, 2): ("sine_h", 0.16), + (0, 3): ("sine_h", 0.32), + (1, 0): ("square_h", 0.04), + (1, 1): ("square_h", 0.08), + (1, 2): ("square_h", 0.16), + (1, 3): ("square_h", 0.32), + (2, 0): ("sine_v", 0.08), + (2, 1): ("sine_45", 0.08), + (2, 2): ("siemens", 0.00), + (2, 3): ("edge", 0.00), +} + +N_ROWS = 3 +N_COLS = 4 + +# Human-readable zone descriptions exported for use in the dialog info panel +ZONE_LABELS: dict[tuple[int, int], str] = { + (0, 0): "Sine H f=0.04 c/px (λ≈25 px, wavelet L4)", + (0, 1): "Sine H f=0.08 c/px (λ≈12 px, wavelet L3)", + (0, 2): "Sine H f=0.16 c/px (λ≈6 px, wavelet L2)", + (0, 3): "Sine H f=0.32 c/px (λ≈3 px, wavelet L1 / near Nyquist)", + (1, 0): "Square H f=0.04 c/px", + (1, 1): "Square H f=0.08 c/px", + (1, 2): "Square H f=0.16 c/px", + (1, 3): "Square H f=0.32 c/px", + (2, 0): "Sine V f=0.08 c/px (vertical direction)", + (2, 1): "Sine 45° f=0.08 c/px (diagonal direction)", + (2, 2): "Siemens star (all-angle radial chirp)", + (2, 3): "Slant edge ~5° (ISO 12233-style ESF input)", +} + + +class SpatialTargetGenerator: + """Generate a calibrated spatial test target FITS image pair. + + Always produces two FITS files: + • clean reference — pure pattern, no PSF, no noise + • degraded companion — same pattern with optional PSF blur and/or noise + Load clean → Image A and degraded → Image B, then run the analysis to see + how each metric responds to known spatial-frequency content. + """ + + def generate(self, params: dict, + preview: bool = False) -> tuple[str, str] | np.ndarray: + """Generate the target pair. + + preview=True → returns degraded ndarray (float32) for live display + preview=False → returns (clean_path, degraded_path) as FITS paths + """ + width = int(params["width"]) + height = int(params["height"]) + contrast_min = float(params.get("contrast_min", params.get("contrast", 0.02))) + contrast_max = float(params.get("contrast_max", params.get("contrast", 0.50))) + sky_adu = float(params["sky_adu"]) + fwhm_px = float(params.get("fwhm_px", 0.0)) + beta = float(params.get("moffat_beta", 4.77)) + rn_adu = float(params.get("read_noise_adu", 10.0)) + + if preview: + scale = min(640 / width, 480 / height, 1.0) + width = max(64, int(width * scale) & ~1) + height = max(64, int(height * scale) & ~1) + + clean = self._build_pattern(width, height, contrast_min, contrast_max, sky_adu) + degraded = clean.copy() + + if params.get("apply_psf") and fwhm_px > 0: + degraded = self._apply_psf(degraded, fwhm_px, beta) + if params.get("apply_noise"): + degraded = self._add_noise(degraded, sky_adu, rn_adu, + seed=int(params.get("seed", 42))) + + clean = np.clip(clean, 0.0, 65535.0).astype(np.float32) + degraded = np.clip(degraded, 0.0, 65535.0).astype(np.float32) + + if preview: + return degraded + + clean_path = self._write_fits(clean, params, is_clean=True) + degraded_path = self._write_fits(degraded, params, is_clean=False) + return (clean_path, degraded_path) + + # ------------------------------------------------------------------ + # Pattern construction + # ------------------------------------------------------------------ + + def _build_pattern(self, W: int, H: int, + contrast_min: float, contrast_max: float, + dc: float) -> np.ndarray: + """Assemble the N_ROWS × N_COLS zone grid as a float64 array.""" + image = np.full((H, W), dc, dtype=np.float64) + zw = W // N_COLS + zh = H // N_ROWS + amplitude_min = dc * contrast_min + amplitude_max = dc * contrast_max + + for (row, col), (waveform, freq) in _ZONE_SPEC.items(): + x0 = col * zw + y0 = row * zh + x1 = x0 + zw if col < N_COLS - 1 else W + y1 = y0 + zh if row < N_ROWS - 1 else H + zW, zH = x1 - x0, y1 - y0 + image[y0:y1, x0:x1] = self._make_zone( + waveform, freq, zW, zH, amplitude_min, amplitude_max, dc) + return image + + def _make_zone(self, waveform: str, freq: float, + W: int, H: int, + amplitude_min: float, amplitude_max: float, + dc: float) -> np.ndarray: + """Return an (H × W) float64 zone with contrast ramping top→bottom. + + Contrast (Michelson) increases linearly from amplitude_min/dc at the + top edge to amplitude_max/dc at the bottom edge of the zone. At any + horizontal row across all N_COLS zones the contrast is identical, + enabling direct cross-frequency comparison at the same contrast level. + """ + y_idx, x_idx = np.mgrid[0:H, 0:W].astype(np.float64) + # Y-ramp: 0 at top of zone, 1 at bottom — independent of absolute image position + y_rel = y_idx / max(H - 1, 1) + amplitude_map = amplitude_min + (amplitude_max - amplitude_min) * y_rel + + if waveform == "sine_h": + pattern = np.sin(2.0 * np.pi * freq * x_idx) + + elif waveform == "square_h": + # Hard square wave — preserve edge sharpness to test Gibbs behavior + pattern = np.where(np.sin(2.0 * np.pi * freq * x_idx) >= 0, 1.0, -1.0) + + elif waveform == "sine_v": + pattern = np.sin(2.0 * np.pi * freq * y_idx) + + elif waveform == "sine_45": + # Diagonal at 45°: spatial frequency along (1,1) direction + d = (x_idx + y_idx) / np.sqrt(2.0) + pattern = np.sin(2.0 * np.pi * freq * d) + + elif waveform == "siemens": + size = min(W, H) + if size % 2 == 0: + size -= 1 # odd so the star centre falls on a pixel + ss = _siemens_star(size) # values 0–~0.5, tapered to 0 at centre/edge + # Centre at 0 and normalise to [-1, 1] so the zone DC = dc + ss_c = ss - ss.mean() + ss_max = np.abs(ss_c).max() + if ss_max > 0: + ss_c /= ss_max + pat = np.zeros((H, W), dtype=np.float64) + sx = (W - size) // 2; sy = (H - size) // 2 + tx = max(0, sx); ty = max(0, sy) + ox = max(0, -sx); oy = max(0, -sy) + cw = min(size - ox, W - tx) + ch = min(size - oy, H - ty) + if cw > 0 and ch > 0: + pat[ty:ty+ch, tx:tx+cw] = ss_c[oy:oy+ch, ox:ox+cw] + pattern = pat + + elif waveform == "edge": + # Slant edge at ~5° from vertical — hard step for ISO 12233-style ESF + angle_rad = np.radians(5.0) + cx, cy = W / 2.0, H / 2.0 + # Signed distance from the slant line (positive = right side) + dist = ((x_idx - cx) * np.cos(angle_rad) + - (y_idx - cy) * np.sin(angle_rad)) + pattern = np.where(dist >= 0, 1.0, -1.0) + + else: + pattern = np.zeros((H, W), dtype=np.float64) + + return dc + amplitude_map * pattern + + # ------------------------------------------------------------------ + # Degradation + # ------------------------------------------------------------------ + + def _apply_psf(self, image: np.ndarray, + fwhm_px: float, beta: float) -> np.ndarray: + """Convolve with an isotropic Moffat PSF (field-centre, no aberrations).""" + alpha = _fwhm_to_alpha(fwhm_px, beta) + stamp = max(5, int(fwhm_px * 6)) + if stamp % 2 == 0: + stamp += 1 + kern = _moffat_2d(stamp, alpha, beta) + kern /= kern.sum() + result = fftconvolve(image, kern, mode="same") + return np.maximum(result, 0.0) + + def _add_noise(self, image: np.ndarray, sky_adu: float, + read_noise_adu: float, seed: int = 42) -> np.ndarray: + """Add Poisson sky shot noise + Gaussian read noise (gain=1 assumed).""" + rng = np.random.default_rng(seed) + sky_e = max(0.1, sky_adu) + poisson = rng.poisson(sky_e, size=image.shape).astype(np.float64) + rdnoise = rng.normal(0.0, read_noise_adu, size=image.shape) + return image + poisson + rdnoise + + # ------------------------------------------------------------------ + # FITS output + # ------------------------------------------------------------------ + + def _write_fits(self, image: np.ndarray, + params: dict, is_clean: bool) -> str: + from astropy.io import fits + + W = image.shape[1] + H = image.shape[0] + contrast_min = float(params.get("contrast_min", params.get("contrast", 0.02))) + contrast_max = float(params.get("contrast_max", params.get("contrast", 0.50))) + sky_adu = float(params["sky_adu"]) + + hdr = fits.Header() + hdr["SIMPLE"] = True + hdr["BITPIX"] = -32 + hdr["NAXIS"] = 2 + hdr["NAXIS1"] = W + hdr["NAXIS2"] = H + hdr["INSTRUME"] = "SpatialTarget" + hdr["EGAIN"] = (1.0, "e-/ADU gain=1 assumed") + hdr["GAIN"] = (1.0, "e-/ADU gain=1 assumed") + # Grid metadata — allows future auto-detection of zone layout in analysis + hdr["TGT_TYPE"] = ("spatial_grid", "generator: spatial test target") + hdr["TGT_ROWS"] = (N_ROWS, "grid rows") + hdr["TGT_COLS"] = (N_COLS, "grid columns") + hdr["TGT_CMIN"] = (contrast_min, "Michelson contrast at zone top (min)") + hdr["TGT_CMAX"] = (contrast_max, "Michelson contrast at zone bottom (max)") + hdr["TGT_SKY"] = (sky_adu, "DC sky level ADU") + hdr["TGT_CLEN"] = (is_clean, "True=clean reference image") + if params.get("apply_psf") and not is_clean: + hdr["TGT_FWHM"] = (float(params.get("fwhm_px", 0.0)), "PSF FWHM px") + hdr["TGT_BETA"] = (float(params.get("moffat_beta", 4.77)), "Moffat beta") + if params.get("apply_noise") and not is_clean: + hdr["TGT_RN"] = (float(params.get("read_noise_adu", 10.0)), "read noise ADU") + hdr["TGT_NSED"] = (int(params.get("seed", 42)), "noise RNG seed") + # Per-zone frequency and waveform type for future per-zone analysis + for (row, col), (waveform, freq) in _ZONE_SPEC.items(): + hdr[f"TGT_{row}{col}F"] = (freq, f"r{row}c{col} freq c/px") + hdr[f"TGT_{row}{col}W"] = (waveform[:8], f"r{row}c{col} waveform") + + hdu = fits.PrimaryHDU(image, header=hdr) + hdul = fits.HDUList([hdu]) + + tag = "clean" if is_clean else "degraded" + stem = f"target_Cmin{contrast_min:.2f}_Cmax{contrast_max:.2f}_SKY{sky_adu:.0f}" + if params.get("apply_psf") and not is_clean: + stem += f"_PSF{params.get('fwhm_px', 0.0):.1f}px" + if params.get("apply_noise") and not is_clean: + stem += f"_RN{params.get('read_noise_adu', 10.0):.0f}" + out_path = Path(params["output_dir"]) / f"{stem}_{tag}.fits" + hdul.writeto(str(out_path), overwrite=True) + return str(out_path) diff --git a/tests/test_analysis/test_edge_analyzer.py b/tests/test_analysis/test_edge_analyzer.py index 676b3ae..bf3636a 100644 --- a/tests/test_analysis/test_edge_analyzer.py +++ b/tests/test_analysis/test_edge_analyzer.py @@ -1,13 +1,38 @@ """Unit tests for analysis/edge_analyzer.py.""" from __future__ import annotations +import numpy as np import pytest +from scipy.ndimage import gaussian_filter from analysis.edge_analyzer import EdgeAnalyzer +from core.models import EDGE_ESF_MIN_MONOTONICITY _RESULT_KEYS = {"edges", "n_edges", "rois_used"} +def _make_clean_edge_roi(angle_deg: float = 30.0, size: int = 60) -> np.ndarray: + """Single straight edge through the box center, background-subtracted + semantics (background ~ 0, signal positive) -- matches real bgsub data.""" + yy, xx = np.mgrid[0:size, 0:size] + c = (size - 1) / 2.0 + theta = np.radians(angle_deg) + d = (xx - c) * np.cos(theta) + (yy - c) * np.sin(theta) + roi = np.where(d > 0, 200.0, 0.0).astype(float) + return gaussian_filter(roi, sigma=1.5) + + +def _make_double_edge_roi(size: int = 60) -> np.ndarray: + """Thin bright stripe crossing the box -- genuinely two edges for any + perpendicular scan direction; no rotation/masking fix can turn this into + a single clean transition.""" + yy, xx = np.mgrid[0:size, 0:size] + c = (size - 1) / 2.0 + stripe = np.abs((xx - c) - 0.3 * (yy - c)) < 5 + roi = np.where(stripe, 200.0, 0.0).astype(float) + return gaussian_filter(roi, sigma=1.5) + + class TestAnalyze: def test_returns_dict(self, astro_image_a): result = EdgeAnalyzer().analyze(astro_image_a) @@ -35,6 +60,154 @@ def test_with_roi(self, astro_image_a): assert isinstance(result, dict) +class TestEsfQuality: + def test_clean_edge_scores_high(self): + ea = EdgeAnalyzer() + roi = _make_clean_edge_roi(angle_deg=30.0) + edge_info = ea._detect_strongest_edge(roi) + _, esf = ea._extract_esf(roi, edge_info) + assert ea._esf_quality(esf) > 0.8 + + def test_double_edge_scores_low(self): + ea = EdgeAnalyzer() + roi = _make_double_edge_roi() + edge_info = ea._detect_strongest_edge(roi) + _, esf = ea._extract_esf(roi, edge_info) + assert ea._esf_quality(esf) < EDGE_ESF_MIN_MONOTONICITY + + def test_perfectly_flat_scores_zero(self): + ea = EdgeAnalyzer() + assert ea._esf_quality(np.full(60, 0.5)) == 0.0 + + @pytest.mark.parametrize("angle_deg", [15.0, 30.0, 60.0, 75.0]) + def test_clean_edge_scores_high_at_various_angles(self, angle_deg): + # Regression: before the disc-mask fix, a clean oblique edge scored as + # low as ~0.0-0.3 (rotate() zero-padding the clipped box corners + # fabricated a second transition), indistinguishable from a genuinely + # bad edge. 45deg is intentionally excluded: a boundary passing exactly + # through both box corners gives Sobel a perfect gradient tie along the + # whole diagonal, so argmax picks a corner-adjacent pixel instead of + # the center -- a narrow, non-representative synthetic degeneracy + # (real edges are never exactly 45.000 deg through both corners) that + # the quality gate correctly flags rather than silently mismeasures. + ea = EdgeAnalyzer() + roi = _make_clean_edge_roi(angle_deg=angle_deg) + edge_info = ea._detect_strongest_edge(roi) + _, esf = ea._extract_esf(roi, edge_info) + assert ea._esf_quality(esf) > 0.8 + + +class TestExtractEsfDiscMask: + def test_no_nan_in_returned_esf(self): + # _extract_esf trims the disc-masked NaN columns internally -- + # callers should never see NaN. + ea = EdgeAnalyzer() + roi = _make_clean_edge_roi(angle_deg=45.0) + edge_info = ea._detect_strongest_edge(roi) + positions, esf = ea._extract_esf(roi, edge_info) + assert esf is not None + assert not np.any(np.isnan(esf)) + assert not np.any(np.isnan(positions)) + + def test_positions_start_at_zero(self): + ea = EdgeAnalyzer() + roi = _make_clean_edge_roi(angle_deg=45.0) + edge_info = ea._detect_strongest_edge(roi) + positions, _ = ea._extract_esf(roi, edge_info) + assert positions[0] == 0.0 + + def test_esf_normalised_to_unit_range(self): + ea = EdgeAnalyzer() + roi = _make_clean_edge_roi(angle_deg=30.0) + edge_info = ea._detect_strongest_edge(roi) + _, esf = ea._extract_esf(roi, edge_info) + assert esf.min() >= -1e-9 + assert esf.max() <= 1.0 + 1e-9 + + +class TestQualityGateAutoDetect: + """Directly control which candidate ROIs _auto_detect_top_rois returns so + the skip/fallback control flow can be tested deterministically, without + depending on real Background2D estimation picking particular peaks.""" + + def test_bad_candidates_skipped_in_favour_of_clean_one(self, astro_image_a, monkeypatch): + bad1 = _make_double_edge_roi() + bad2 = _make_double_edge_roi() + clean = _make_clean_edge_roi(angle_deg=30.0) + candidates = [ + (bad1, (0, 0, 60, 60)), + (bad2, (100, 100, 160, 160)), + (clean, (200, 200, 260, 260)), + ] + monkeypatch.setattr( + EdgeAnalyzer, "_auto_detect_top_rois", + lambda self, bgsub, image, n=9: candidates) + + result = EdgeAnalyzer().analyze(astro_image_a) + assert result["n_edges"] == 1 + assert result["edges"][0]["low_confidence"] is False + assert result["edges"][0]["roi_used"] == (200, 200, 260, 260) + + def test_all_bad_falls_back_to_least_bad_flagged(self, astro_image_a, monkeypatch): + bad1 = _make_double_edge_roi() + bad2 = _make_double_edge_roi() + candidates = [ + (bad1, (0, 0, 60, 60)), + (bad2, (100, 100, 160, 160)), + ] + monkeypatch.setattr( + EdgeAnalyzer, "_auto_detect_top_rois", + lambda self, bgsub, image, n=9: candidates) + + result = EdgeAnalyzer().analyze(astro_image_a) + assert result["n_edges"] == 1 + assert result["edges"][0]["low_confidence"] is True + + def test_edge_entries_have_quality_keys(self, astro_image_a, monkeypatch): + clean = _make_clean_edge_roi(angle_deg=30.0) + monkeypatch.setattr( + EdgeAnalyzer, "_auto_detect_top_rois", + lambda self, bgsub, image, n=9: [(clean, (0, 0, 60, 60))]) + + result = EdgeAnalyzer().analyze(astro_image_a) + assert result["n_edges"] == 1 + entry = result["edges"][0] + assert "esf_quality" in entry + assert "low_confidence" in entry + assert isinstance(entry["esf_quality"], float) + + +class TestQualityGateExplicitRoi: + """A user-drawn or A-matched ROI is fixed -- low quality must still be + flagged, but the edge must not be silently dropped (there's no + alternative candidate to fall back to).""" + + def test_bad_edge_kept_but_flagged_with_explicit_roi(self, tmp_path): + from astropy.io import fits as ap_fits + from core.astro_image import AstroImage + + h = w = 256 + rng = np.random.default_rng(0) + data = rng.normal(1000.0, 20.0, (h, w)).astype(np.float64) + yy, xx = np.mgrid[0:h, 0:w] + stripe = np.abs((xx - 128) - 0.3 * (yy - 128)) < 5 + data[stripe] += 200.0 + data = np.clip(data, 0, 65535).astype(np.float32) + hdr = ap_fits.Header() + hdr["EGAIN"] = 1.0 + hdr["GAIN"] = 1.0 + path = tmp_path / "stripe.fits" + ap_fits.writeto(str(path), data, hdr, overwrite=True) + + img = AstroImage(str(path), label="Stripe") + img.load() + img.estimate_background() + + result = EdgeAnalyzer().analyze(img, roi=(98, 98, 158, 158)) + assert result["n_edges"] == 1 + assert result["edges"][0]["low_confidence"] is True + + class TestAnalyzeCrosshair: def test_degenerate_crosshair_returns_none(self, astro_image_a): # Zero-length crosshair (start == end) diff --git a/tests/test_analysis/test_spatial_detail.py b/tests/test_analysis/test_spatial_detail.py index 5986f8f..2a59626 100644 --- a/tests/test_analysis/test_spatial_detail.py +++ b/tests/test_analysis/test_spatial_detail.py @@ -7,6 +7,7 @@ from analysis.image_filters import SpatialDetailAnalyzer from core.astro_image import AstroImage +from core.models import STD_KERNEL_SIZES, LOG_SIGMAS, WEBER_KERNEL_SIZES, WAVELET_LEVELS class TestAnalyze: @@ -26,6 +27,13 @@ def test_panels_present(self, astro_image_a): result = SpatialDetailAnalyzer().analyze(astro_image_a) assert "panels" in result + def test_original_panel_present_single_image(self, astro_image_a): + result = SpatialDetailAnalyzer().analyze(astro_image_a) + original = result["panels"]["original"] + assert original["a"] is not None + assert original["b"] is None + assert original["diff"] is None + def test_contrast_ratios_are_positive(self, astro_image_a): result = SpatialDetailAnalyzer().analyze(astro_image_a) ratios = result.get("contrast_ratios_a") or [] @@ -55,7 +63,255 @@ def test_minimal_image_no_crash(self, tmp_path): result = SpatialDetailAnalyzer().analyze(img) assert isinstance(result, dict) + def test_weber_contrast_a_present(self, astro_image_a): + result = SpatialDetailAnalyzer().analyze(astro_image_a) + assert "weber_contrast_a" in result + + def test_single_image_weber_contrast_b_empty(self, astro_image_a): + result = SpatialDetailAnalyzer().analyze(astro_image_a) + # Single-image mode: weber_contrast_b is present but empty (same pattern as contrast_ratios_b) + wc_b = result.get("weber_contrast_b") + assert wc_b is None or not wc_b + + def test_weber_contrast_a_positive(self, astro_image_a): + result = SpatialDetailAnalyzer().analyze(astro_image_a) + for v in result.get("weber_contrast_a", {}).values(): + assert v >= 0.0 + def test_with_roi(self, astro_image_a): result = SpatialDetailAnalyzer().analyze(astro_image_a, roi=(50, 50, 450, 450)) assert isinstance(result, dict) + + +# --------------------------------------------------------------------------- +# Noise-corrected multi-scale local contrast (two-image A/B comparison) +# --------------------------------------------------------------------------- + +def _make_nc_test_fits(path, add_texture: bool, seed: int) -> None: + """256x256 FITS with a smooth nebula blob (sigma=25, well above 2*rms after + background subtraction). When add_texture, a fine checkerboard (period 6px, + amplitude 8x sky noise) is added inside the blob so std/LoG/wavelet/Weber/ + gradient all detect meaningfully more local structure than the plain blob.""" + rng = np.random.default_rng(seed) + h, w = 256, 256 + sky_level = 1000.0 + sky_noise = 20.0 + data = rng.normal(sky_level, sky_noise, (h, w)).astype(np.float64) + + yy, xx = np.mgrid[0:h, 0:w] + cx, cy = w / 2, h / 2 + blob_sigma = 25.0 + blob_amp = 600.0 + blob = blob_amp * np.exp(-0.5 * (((xx - cx) ** 2 + (yy - cy) ** 2) / blob_sigma ** 2)) + data += blob + + if add_texture: + texture_amp = 8.0 * sky_noise + checker = (((xx // 3).astype(int) + (yy // 3).astype(int)) % 2 == 0) + blob_mask = blob > (0.3 * blob_amp) + data += np.where(blob_mask & checker, texture_amp, 0.0) + + data = np.clip(data, 0, 65535).astype(np.float32) + hdr = ap_fits.Header() + hdr["EGAIN"] = 1.0 + hdr["GAIN"] = 1.0 + hdr["FOCALLEN"] = 500.0 + hdr["XPIXSZ"] = 3.76 + hdr["EXPTIME"] = 300.0 + hdr["INSTRUME"] = "TestCam" + ap_fits.writeto(str(path), data, hdr, overwrite=True) + + +@pytest.fixture(scope="module") +def nc_image_pair(tmp_path_factory) -> tuple[AstroImage, AstroImage]: + """Image A = blob + fine checkerboard texture; Image B = blob only. + Both share the same blob location/extent, so their nebula masks overlap.""" + out = tmp_path_factory.mktemp("nc_fits") + path_a = out / "nc_a.fits" + path_b = out / "nc_b.fits" + _make_nc_test_fits(path_a, add_texture=True, seed=1) + _make_nc_test_fits(path_b, add_texture=False, seed=2) + img_a = AstroImage(str(path_a), label="A") + img_a.load() + img_a.estimate_background() + img_b = AstroImage(str(path_b), label="B") + img_b.load() + img_b.estimate_background() + return img_a, img_b + + +@pytest.fixture(scope="module") +def nc_result(nc_image_pair) -> dict: + img_a, img_b = nc_image_pair + return SpatialDetailAnalyzer().analyze(img_a, img_b) + + +class TestNoiseCorrectedContrast: + """A/B noise-corrected multi-scale local contrast: shared-nebula-ROI scoring, + per-scale A/B ratios, and noise-normalised display panels.""" + + def test_shared_nebula_pixels_positive(self, nc_result): + assert nc_result["nc_shared_nebula_pixels"] > 0 + + @pytest.mark.parametrize("prefix,scales", [ + ("std", STD_KERNEL_SIZES), + ("log", LOG_SIGMAS), + ("weber", WEBER_KERNEL_SIZES), + ("gm", LOG_SIGMAS), + ]) + def test_nc_score_dict_keys_match_scales(self, nc_result, prefix, scales): + assert set(nc_result[f"{prefix}_nc_score_a"].keys()) == set(scales) + assert set(nc_result[f"{prefix}_nc_score_b"].keys()) == set(scales) + + def test_wavelet_nc_score_keys_match_levels(self, nc_result): + expected = set(range(1, WAVELET_LEVELS + 1)) + assert set(nc_result["wavelet_nc_score_a"].keys()) == expected + assert set(nc_result["wavelet_nc_score_b"].keys()) == expected + + @pytest.mark.parametrize("prefix", ["std", "log", "wavelet", "weber", "gm"]) + def test_nc_noise_floor_positive_or_none(self, nc_result, prefix): + for side in ("a", "b"): + for v in nc_result[f"{prefix}_nc_noise_{side}"].values(): + if v is not None: + assert v > 0.0 + + def test_std_nc_ratio_captures_finer_detail_in_a(self, nc_result): + ratio = nc_result["std_nc_ratio"][min(STD_KERNEL_SIZES)] + assert ratio is not None and ratio > 1.05 + + def test_log_nc_ratio_captures_finer_detail_in_a(self, nc_result): + ratio = nc_result["log_nc_ratio"][min(LOG_SIGMAS)] + assert ratio is not None and ratio > 1.05 + + def test_wavelet_nc_ratio_captures_finer_detail_in_a(self, nc_result): + ratio = nc_result["wavelet_nc_ratio"][2] + assert ratio is not None and ratio > 1.05 + + def test_weber_nc_ratio_captures_finer_detail_in_a(self, nc_result): + ratio = nc_result["weber_nc_ratio"][min(WEBER_KERNEL_SIZES)] + assert ratio is not None and ratio > 1.05 + + def test_gradient_nc_ratio_captures_finer_detail_in_a(self, nc_result): + # Gradient magnitude's peak response scale for a given texture need not be + # the finest sigma (Gaussian smoothing at small sigma can attenuate a very + # fine checkerboard before the derivative is taken) — check the strongest + # response across scales rather than pinning to one bin. + values = [v for v in nc_result["gm_nc_ratio"].values() if v is not None] + assert values and max(values) > 1.05 + + def test_normalized_panels_present_two_image(self, nc_result): + panels = nc_result["panels"] + for ks in STD_KERNEL_SIZES: + assert f"nrm_std_{ks}px" in panels + for sigma in LOG_SIGMAS: + assert f"nrm_log_{sigma}" in panels + assert f"nrm_gradient_{sigma}" in panels + for ks in WEBER_KERNEL_SIZES: + assert f"nrm_weber_{ks}px" in panels + for lvl in (2, 3): + assert f"nrm_wavelet_{lvl}" in panels + + def test_original_panel_present_two_image(self, nc_result): + original = nc_result["panels"]["original"] + assert original["a"] is not None + assert original["b"] is not None + assert original["diff"] is not None + assert original["a"].shape == original["b"].shape + + def test_normalized_panel_values_differ_from_raw(self, nc_result): + panels = nc_result["panels"] + ks = min(STD_KERNEL_SIZES) + raw_a = panels[f"std_{ks}px"]["a"] + nrm_a = panels[f"nrm_std_{ks}px"]["a"] + assert not np.allclose(raw_a, nrm_a) + + def test_nc_ratio_overview_figure_present(self, nc_result): + assert "nc_ratio_overview" in nc_result["figures"] + + def test_with_roi_two_image_no_crash(self, nc_image_pair): + img_a, img_b = nc_image_pair + result = SpatialDetailAnalyzer().analyze( + img_a, img_b, roi=(20, 20, 236, 236)) + assert isinstance(result, dict) + assert "std_nc_ratio" in result + + +class TestNoiseCorrectedContrastSingleImage: + """Single-image mode: every new NC key must be empty, matching the existing + contrast_ratios_b / weber_contrast_b invariant (never None/absent, never + populated with placeholder values on the A side either).""" + + @pytest.fixture(scope="class") + @classmethod + def single_result(cls, astro_image_a): + return SpatialDetailAnalyzer().analyze(astro_image_a) + + @pytest.mark.parametrize("key", [ + "std_nc_score_a", "std_nc_score_b", "std_nc_ratio", + "log_nc_score_a", "log_nc_score_b", "log_nc_ratio", + "wavelet_nc_score_a", "wavelet_nc_score_b", "wavelet_nc_ratio", + "weber_nc_score_a", "weber_nc_score_b", "weber_nc_ratio", + "gm_nc_score_a", "gm_nc_score_b", "gm_nc_ratio", + ]) + def test_nc_key_empty_in_single_image_mode(self, single_result, key): + assert not single_result[key] + + def test_shared_nebula_pixels_zero(self, single_result): + assert single_result["nc_shared_nebula_pixels"] == 0 + + def test_no_normalized_panels_in_single_image_mode(self, single_result): + assert not any(k.startswith("nrm_") for k in single_result["panels"]) + + +class TestNcScoreHelper: + """Direct unit tests of _nc_score's mask-emptiness contract — deterministic, + unlike relying on two real images happening to produce non-overlapping + nebula masks (background-threshold noise can create incidental overlap).""" + + def test_none_mask_neb_shared_returns_none(self): + analyzer = SpatialDetailAnalyzer() + detail = np.ones((20, 20), dtype=np.float32) + bg_mask = np.ones((20, 20), dtype=bool) + score, noise = analyzer._nc_score(detail, None, bg_mask) + assert score is None and noise is None + + def test_empty_shared_nebula_mask_returns_none(self): + analyzer = SpatialDetailAnalyzer() + detail = np.ones((20, 20), dtype=np.float32) + mask_neb_shared = np.zeros((20, 20), dtype=bool) # no shared nebula pixels + bg_mask = np.ones((20, 20), dtype=bool) + score, noise = analyzer._nc_score(detail, mask_neb_shared, bg_mask) + assert score is None and noise is None + + def test_empty_bg_mask_returns_none(self): + analyzer = SpatialDetailAnalyzer() + detail = np.ones((20, 20), dtype=np.float32) + mask_neb_shared = np.ones((20, 20), dtype=bool) + bg_mask = np.zeros((20, 20), dtype=bool) # no background pixels + score, noise = analyzer._nc_score(detail, mask_neb_shared, bg_mask) + assert score is None and noise is None + + def test_zero_noise_floor_returns_none(self): + analyzer = SpatialDetailAnalyzer() + detail = np.zeros((20, 20), dtype=np.float32) + detail[:10, :] = 5.0 # nebula half has signal, background half is exactly 0 + mask_neb_shared = np.zeros((20, 20), dtype=bool) + mask_neb_shared[:10, :] = True + bg_mask = np.zeros((20, 20), dtype=bool) + bg_mask[10:, :] = True + score, noise = analyzer._nc_score(detail, mask_neb_shared, bg_mask) + assert score is None and noise is None + + def test_valid_masks_return_ratio(self): + analyzer = SpatialDetailAnalyzer() + detail = np.zeros((20, 20), dtype=np.float32) + detail[:10, :] = 10.0 + detail[10:, :] = 2.0 + mask_neb_shared = np.zeros((20, 20), dtype=bool) + mask_neb_shared[:10, :] = True + bg_mask = np.zeros((20, 20), dtype=bool) + bg_mask[10:, :] = True + score, noise = analyzer._nc_score(detail, mask_neb_shared, bg_mask) + assert score == pytest.approx(5.0) + assert noise == pytest.approx(2.0) diff --git a/tests/test_core/test_astro_image.py b/tests/test_core/test_astro_image.py index 33bf633..867f67d 100644 --- a/tests/test_core/test_astro_image.py +++ b/tests/test_core/test_astro_image.py @@ -10,8 +10,8 @@ class TestLoad: - def test_data_float64(self, astro_image_a): - assert astro_image_a.data.dtype == np.float64 + def test_data_float32(self, astro_image_a): + assert astro_image_a.data.dtype == np.float32 def test_correct_shape(self, astro_image_a): assert astro_image_a.data.shape == (512, 512) @@ -109,7 +109,7 @@ def test_background_subtracted_shape(self, astro_image_a): def test_background_subtracted_dtype(self, astro_image_a): bs = astro_image_a.background_subtracted() - assert bs.dtype == np.float64 + assert bs.dtype == np.float32 def test_background_subtracted_without_background(self, test_fits_path): img = AstroImage(str(test_fits_path), label="NoBkg") diff --git a/tests/test_report/test_inspector_catalog.py b/tests/test_report/test_inspector_catalog.py new file mode 100644 index 0000000..87f574b --- /dev/null +++ b/tests/test_report/test_inspector_catalog.py @@ -0,0 +1,93 @@ +"""Regression tests for the Report Inspector's .npz panel catalog. + +Covers the bug where _write_inspector_file used a hardcoded, stale panel-name +map (std_15px/std_31px from an old STD_KERNEL_SIZES, no Weber entries at all) +that silently omitted panels SpatialDetailAnalyzer actually computes. The fix +derives the catalog dynamically from panels_a.keys() via _panel_display_name. +""" +from __future__ import annotations + +import json + +import numpy as np +import pytest +from astropy.io import fits as ap_fits + +from analysis.image_filters import SpatialDetailAnalyzer +from core.astro_image import AstroImage +from core.models import AnalysisResult +from report.report_builder import ReportBuilder + + +def _make_two_image_pair(tmp_path_factory) -> tuple[AstroImage, AstroImage]: + out = tmp_path_factory.mktemp("inspector_fits") + rng_a = np.random.default_rng(1) + rng_b = np.random.default_rng(2) + h, w = 256, 256 + for path, rng in ((out / "a.fits", rng_a), (out / "b.fits", rng_b)): + data = rng.normal(1000.0, 20.0, (h, w)).astype(np.float32) + yy, xx = np.mgrid[0:h, 0:w] + blob = 500.0 * np.exp(-0.5 * (((xx - w / 2) ** 2 + (yy - h / 2) ** 2) / 30.0 ** 2)) + data += blob.astype(np.float32) + hdr = ap_fits.Header() + hdr["EGAIN"] = 1.0 + hdr["GAIN"] = 1.0 + ap_fits.writeto(str(path), data, hdr, overwrite=True) + img_a = AstroImage(str(out / "a.fits"), label="A") + img_a.load() + img_a.estimate_background() + img_b = AstroImage(str(out / "b.fits"), label="B") + img_b.load() + img_b.estimate_background() + return img_a, img_b + + +@pytest.fixture(scope="module") +def inspector_npz_path(tmp_path_factory): + img_a, img_b = _make_two_image_pair(tmp_path_factory) + spatial = SpatialDetailAnalyzer().analyze(img_a, img_b) + result_a = AnalysisResult(label="A", spatial_metrics=spatial) + result_b = AnalysisResult(label="B", spatial_metrics=spatial) + + out_path = tmp_path_factory.mktemp("inspector_out") / "report_inspector.npz" + ReportBuilder()._write_inspector_file(out_path, img_a, img_b, result_a, result_b) + return out_path, spatial + + +class TestInspectorCatalogDynamicPanels: + def test_every_computed_panel_is_cataloged(self, inspector_npz_path): + path, spatial = inspector_npz_path + npz = np.load(str(path), allow_pickle=False) + catalog = json.loads(npz["catalog_json"].tobytes().decode("utf-8")) + spatial_entries = catalog["sections"].get("Spatial Detail", []) + cataloged_names = {e["name"] for e in spatial_entries} + + panels = spatial["panels"] + assert len(spatial_entries) == len(panels) + for pkey in panels: + from report.report_builder import _panel_display_name + assert _panel_display_name(pkey) in cataloged_names + + def test_weber_panels_present(self, inspector_npz_path): + # Regression: the old hardcoded _PANEL_IMAGE_SETS never included Weber at all. + _, spatial = inspector_npz_path + weber_keys = [k for k in spatial["panels"] if k.startswith("weber_")] + assert weber_keys, "fixture should have produced Weber panels" + + def test_all_std_kernel_sizes_present(self, inspector_npz_path): + # Regression: the old map referenced std_15px/std_31px from a stale + # STD_KERNEL_SIZES and dropped std_3px/std_10px from the current one. + from core.models import STD_KERNEL_SIZES + _, spatial = inspector_npz_path + for ks in STD_KERNEL_SIZES: + assert f"std_{ks}px" in spatial["panels"] + + def test_noise_normalized_panels_cataloged(self, inspector_npz_path): + path, spatial = inspector_npz_path + npz = np.load(str(path), allow_pickle=False) + catalog = json.loads(npz["catalog_json"].tobytes().decode("utf-8")) + spatial_entries = catalog["sections"].get("Spatial Detail", []) + cataloged_names = {e["name"] for e in spatial_entries} + nrm_keys = [k for k in spatial["panels"] if k.startswith("nrm_")] + assert nrm_keys + assert any("(noise-normalized)" in n for n in cataloged_names) diff --git a/tests/test_report/test_report_helpers.py b/tests/test_report/test_report_helpers.py index fe554d9..2e47abd 100644 --- a/tests/test_report/test_report_helpers.py +++ b/tests/test_report/test_report_helpers.py @@ -16,6 +16,9 @@ _arr_to_b64_png, _fig_to_b64, _psf_stat_test, + _power_ratio_db, + _nc_ratio_rows, + _panel_display_name, ) EM_DASH = "—" @@ -182,3 +185,105 @@ def test_custom_dpi_changes_size(self): lo = _fig_to_b64(self._make_fig(), dpi=72) hi = _fig_to_b64(self._make_fig(), dpi=200) assert len(base64.b64decode(hi)) > len(base64.b64decode(lo)) + + +class TestPowerRatioDb: + def test_known_ratio_in_db(self): + freq = np.linspace(0.0, 0.5, 10) + rp_b = np.full(10, 1.0) + rp_a = rp_b * 2.0 + result = _power_ratio_db(freq, rp_a, freq, rp_b) + assert result is not None + _, ratio_db = result + assert ratio_db == pytest.approx(10.0 * np.log10(2.0), abs=1e-6) + + def test_identical_curves_zero_db(self): + freq = np.linspace(0.0, 0.5, 10) + rp = np.linspace(1.0, 5.0, 10) + _, ratio_db = _power_ratio_db(freq, rp, freq, rp) + assert ratio_db == pytest.approx(0.0, abs=1e-6) + + def test_zero_bins_stay_finite(self): + freq = np.linspace(0.0, 0.5, 5) + rp_a = np.array([0.0, 1.0, 0.0, 2.0, 0.0]) + rp_b = np.array([1.0, 0.0, 0.0, 1.0, 1.0]) + _, ratio_db = _power_ratio_db(freq, rp_a, freq, rp_b) + assert np.all(np.isfinite(ratio_db)) + + def test_none_inputs_return_none(self): + freq = np.linspace(0.0, 0.5, 5) + rp = np.ones(5) + assert _power_ratio_db(None, rp, freq, rp) is None + assert _power_ratio_db(freq, None, freq, rp) is None + assert _power_ratio_db(freq, rp, None, rp) is None + assert _power_ratio_db(freq, rp, freq, None) is None + + def test_mismatched_length_returns_none(self): + freq_a = np.linspace(0.0, 0.5, 5) + freq_b = np.linspace(0.0, 0.5, 8) + rp_a = np.ones(5) + rp_b = np.ones(8) + assert _power_ratio_db(freq_a, rp_a, freq_b, rp_b) is None + + def test_same_length_different_values_returns_none(self): + freq_a = np.linspace(0.0, 0.5, 5) + freq_b = np.linspace(0.0, 0.4, 5) + rp = np.ones(5) + assert _power_ratio_db(freq_a, rp, freq_b, rp) is None + + +class TestNcRatioRows: + def _label(self, scale): + return f"{scale} px" + + def test_known_ratio_formatted(self): + rows = _nc_ratio_rows({3: 2.0}, {3: 1.0}, {3: 2.0}, self._label) + assert "2.000" in rows + assert "1.000" in rows + + def test_missing_entry_renders_em_dash(self): + rows = _nc_ratio_rows({3: 1.5}, {}, {}, self._label) + assert EM_DASH in rows + + def test_none_ratio_renders_em_dash(self): + rows = _nc_ratio_rows({3: 1.5}, {3: None}, {3: None}, self._label) + assert EM_DASH in rows + + def test_empty_inputs_produce_no_rows(self): + rows = _nc_ratio_rows({}, {}, {}, self._label) + assert rows == "" + + def test_scale_label_used(self): + rows = _nc_ratio_rows({5: 1.0}, {5: 1.0}, {5: 1.0}, lambda s: f"σ = {s} px") + assert "σ = 5 px" in rows + + def test_ratio_cell_colored_independently_of_ab_columns(self): + # A/B columns: A < B -> "worse"/"better"; ratio itself (0.5) vs parity (1.0) + # is colored via a SEPARATE _better_worse_class(vr, 1.0) call. + rows = _nc_ratio_rows({3: 1.0}, {3: 2.0}, {3: 0.5}, self._label) + assert "worse" in rows # A's own value vs B's + # the ratio 0.5 < 1.0 is also "worse" here, but computed independently + # (not derived from the A/B column classes) - verify it renders at all + assert "0.500" in rows + + def test_union_of_scale_keys(self): + rows = _nc_ratio_rows({3: 1.0}, {5: 1.0}, {}, self._label) + assert "3 px" in rows + assert "5 px" in rows + + +class TestPanelDisplayName: + @pytest.mark.parametrize("pkey,expected", [ + ("std_3px", "Std Dev 3 px"), + ("std_10px", "Std Dev 10 px"), + ("weber_9px", "Weber 9 px"), + ("log_1.5", "LoG σ 1.5 px"), + ("wavelet_2", "Wavelet level 2"), + ("gradient_3.0", "Gradient σ 3.0 px"), + ("nrm_log_1.5", "LoG σ 1.5 px (noise-normalized)"), + ("nrm_std_3px", "Std Dev 3 px (noise-normalized)"), + ("nrm_gradient_6.0", "Gradient σ 6.0 px (noise-normalized)"), + ("original", "Original Image (ROI)"), + ]) + def test_display_name(self, pkey, expected): + assert _panel_display_name(pkey) == expected