From 9752f221db6a7700c8e3614f975502d97d5174a1 Mon Sep 17 00:00:00 2001 From: Brent <52629076+brentmantooth@users.noreply.github.com> Date: Mon, 13 Jul 2026 23:33:15 -0400 Subject: [PATCH 1/8] Add Section 8a diff-distribution violin plot for spatial detail comparison Single scalar noise-corrected scores were hard to trust as "the right metric" for filter comparison. Each A-B diff map (std, LoG, gradient, wavelet, Weber, plus the raw normalised-image diff) now also yields a subsampled pixel-value distribution, split by the same shared nebula/background masks already used for the noise-corrected scores, rendered as one combined violin+IQR-box figure alongside the existing per-scale maps and tables. Co-Authored-By: Claude Sonnet 5 --- analysis/image_filters.py | 124 ++++++++++++++++++++++++++++------- core/models.py | 1 + report/report_builder.py | 134 ++++++++++++++++++++++++++++++++++++++ 3 files changed, 236 insertions(+), 23 deletions(-) diff --git a/analysis/image_filters.py b/analysis/image_filters.py index 9415f5f..d5338f1 100644 --- a/analysis/image_filters.py +++ b/analysis/image_filters.py @@ -18,7 +18,7 @@ from core.models import (STD_KERNEL_SIZES, LOG_SIGMAS, WAVELET_NAME, WAVELET_LEVELS, WEBER_KERNEL_SIZES, XS_LINE_ALPHA, SECTION8_BORDER_CROP_FRACTION, SECTION8_ANALYSIS_CMAP, - XS_SNR_REGION_WIDTH) + XS_SNR_REGION_WIDTH, SECTION8_DIFF_DIST_MAX_SAMPLES) MAX_DIM_FOR_STD = 2048 # downsample to this before generic_filter (performance) _DISPLAY_SMOOTH_SIGMA = 1.0 # applied to maps before plotting; does NOT affect metrics @@ -58,6 +58,7 @@ def analyze(self, image_a: AstroImage, image_b: AstroImage | None = None, "weber_contrast_a": {}, "weber_contrast_b": {}, "panels": {}, + "diff_dist": {}, "nc_shared_nebula_pixels": 0, "std_nc_score_a": {}, "std_nc_score_b": {}, "std_nc_noise_a": {}, "std_nc_noise_b": {}, "std_nc_ratio": {}, @@ -120,26 +121,36 @@ def _clip01(v): return max(0.0, min(1.0, v)) result["display_roi"] = display_roi + # Shared nebula/background regions for noise-corrected A/B scoring and diff + # distributions: pixels BOTH images independently classify the same way. + # None in single-image mode. + mask_neb_shared = mask_bg_shared = None + if mask_neb_b is not None: + h_s = min(mask_neb_a.shape[0], mask_neb_b.shape[0]) + w_s = min(mask_neb_a.shape[1], mask_neb_b.shape[1]) + mask_neb_shared = mask_neb_a[:h_s, :w_s] & mask_neb_b[:h_s, :w_s] + mask_bg_shared = mask_bg_a[:h_s, :w_s] & mask_bg_b[:h_s, :w_s] + result["nc_shared_nebula_pixels"] = ( + int(np.count_nonzero(mask_neb_shared)) if mask_neb_shared is not None else 0 + ) + + # Fixed seed so the diff-distribution subsampling below is reproducible + # across report generations for the same input images. + diff_dist_rng = np.random.default_rng(42) + # Export the exact preprocessed array every std/LoG/wavelet/Weber/gradient # map is computed from (mean-normalised, ROI-cropped if a ROI was given) as # its own panel, so the Report Inspector can show the source image content # for a map alongside the map itself, pixel-aligned to the same crop. + original_diff = (analysis_a - analysis_b).astype(np.float32) if analysis_b is not None else None result["panels"]["original"] = { "a": analysis_a.astype(np.float32), "b": analysis_b.astype(np.float32) if analysis_b is not None else None, - "diff": (analysis_a - analysis_b).astype(np.float32) if analysis_b is not None else None, + "diff": original_diff, } - - # Shared nebula ROI for noise-corrected A/B scoring: pixels BOTH images - # independently classify as nebula. None in single-image mode. - mask_neb_shared = None - if mask_neb_b is not None: - h_s = min(mask_neb_a.shape[0], mask_neb_b.shape[0]) - w_s = min(mask_neb_a.shape[1], mask_neb_b.shape[1]) - mask_neb_shared = mask_neb_a[:h_s, :w_s] & mask_neb_b[:h_s, :w_s] - result["nc_shared_nebula_pixels"] = ( - int(np.count_nonzero(mask_neb_shared)) if mask_neb_shared is not None else 0 - ) + if original_diff is not None: + result["diff_dist"]["original"] = self._diff_distribution( + original_diff, mask_neb_shared, mask_bg_shared, diff_dist_rng) _label_b = image_b.label if image_b is not None else None @@ -156,6 +167,7 @@ def _clip01(v): return max(0.0, min(1.0, v)) display_roi=display_roi, crosshair=crosshair_roi, mask_neb_shared=mask_neb_shared, + mask_bg_shared=mask_bg_shared, diff_dist_rng=diff_dist_rng, ) _f_log = _ex.submit(self._log_analysis, analysis_a, analysis_b, log_sigmas, @@ -163,6 +175,7 @@ def _clip01(v): return max(0.0, min(1.0, v)) display_roi=display_roi, crosshair=crosshair_roi, mask_neb_shared=mask_neb_shared, mask_bg_a=mask_bg_a, mask_bg_b=mask_bg_b, + mask_bg_shared=mask_bg_shared, diff_dist_rng=diff_dist_rng, ) _f_wav = _ex.submit(self._wavelet_analysis, analysis_a, analysis_b, wavelet, levels, @@ -170,6 +183,7 @@ def _clip01(v): return max(0.0, min(1.0, v)) display_roi=display_roi, crosshair=crosshair_roi, mask_neb_shared=mask_neb_shared, mask_bg_a=mask_bg_a, mask_bg_b=mask_bg_b, + mask_bg_shared=mask_bg_shared, diff_dist_rng=diff_dist_rng, ) _f_web = _ex.submit(self._weber_analysis, analysis_a, analysis_b, @@ -177,6 +191,7 @@ def _clip01(v): return max(0.0, min(1.0, v)) image_a.label, _label_b, display_roi=display_roi, mask_neb_shared=mask_neb_shared, mask_bg_a=mask_bg_a, mask_bg_b=mask_bg_b, + mask_bg_shared=mask_bg_shared, diff_dist_rng=diff_dist_rng, ) _f_grad = _ex.submit(self._gradient_analysis, analysis_a, analysis_b, log_sigmas, @@ -184,6 +199,7 @@ def _clip01(v): return max(0.0, min(1.0, v)) display_roi=display_roi, crosshair=crosshair_roi, mask_neb_shared=mask_neb_shared, mask_bg_a=mask_bg_a, mask_bg_b=mask_bg_b, + mask_bg_shared=mask_bg_shared, diff_dist_rng=diff_dist_rng, ) std_b64, std_partial = _f_std.result() log_b64, log_partial = _f_log.result() @@ -209,6 +225,11 @@ def _clip01(v): return max(0.0, min(1.0, v)) result["panels"].update(wav_partial["panels"]) result["panels"].update(web_partial["panels"]) result["panels"].update(grad_partial["panels"]) + result["diff_dist"].update(std_partial["diff_dist"]) + result["diff_dist"].update(log_partial["diff_dist"]) + result["diff_dist"].update(wav_partial["diff_dist"]) + result["diff_dist"].update(web_partial["diff_dist"]) + result["diff_dist"].update(grad_partial["diff_dist"]) # Merge noise-corrected scores/noise-floors and compute A/B ratios centrally. for prefix, partial in (("std", std_partial), ("log", log_partial), @@ -337,13 +358,15 @@ def _std_analysis(self, norm_a, norm_b, kernel_sizes, label_a, label_b, display_roi=None, crosshair=None, - mask_neb_shared=None) -> tuple[dict, dict]: + mask_neb_shared=None, + mask_bg_shared=None, diff_dist_rng=None) -> tuple[dict, dict]: figures = {} partial: dict = { "contrast_ratios_a": {}, "contrast_ratios_b": {}, "std_nc_score_a": {}, "std_nc_score_b": {}, "std_nc_noise_a": {}, "std_nc_noise_b": {}, "panels": {}, + "diff_dist": {}, } single = norm_b is None for ks in kernel_sizes: @@ -366,11 +389,15 @@ def _std_analysis(self, norm_a, norm_b, partial["std_nc_score_b"][ks] = nc_b partial["std_nc_noise_b"][ks] = noise_b + diff = (std_a - std_b).astype(np.float32) if not single else None partial["panels"][f"std_{ks}px"] = { "a": std_a.astype(np.float32), "b": std_b.astype(np.float32) if std_b is not None else None, - "diff": (std_a - std_b).astype(np.float32) if std_b is not None else None, + "diff": diff, } + if diff is not None: + partial["diff_dist"][f"std_{ks}px"] = self._diff_distribution( + diff, mask_neb_shared, mask_bg_shared, diff_dist_rng) if not single and noise_a and noise_b: partial["panels"][f"nrm_std_{ks}px"] = { "a": (std_a / noise_a).astype(np.float32), @@ -476,6 +503,33 @@ def _nc_score(self, detail_map: np.ndarray, return None, None return float(bn.median(neb_vals)) / noise_floor, noise_floor + @staticmethod + def _diff_distribution(diff_map: np.ndarray, + mask_neb_shared: np.ndarray | None, + mask_bg_shared: np.ndarray | None, + rng: np.random.Generator) -> dict: + """Random-subsampled A-B diff pixel populations for nebula vs background. + + Returns {"nebula": ndarray, "background": ndarray} (float32, signed diff + values, up to SECTION8_DIFF_DIST_MAX_SAMPLES each). Either array is empty + if the corresponding mask is unavailable (single-image mode) or selects + zero pixels. + """ + out = {"nebula": np.empty(0, dtype=np.float32), + "background": np.empty(0, dtype=np.float32)} + if mask_neb_shared is None or mask_bg_shared is None: + return out + h = min(diff_map.shape[0], mask_neb_shared.shape[0], mask_bg_shared.shape[0]) + w = min(diff_map.shape[1], mask_neb_shared.shape[1], mask_bg_shared.shape[1]) + cropped = diff_map[:h, :w] + for key, mask in (("nebula", mask_neb_shared), ("background", mask_bg_shared)): + vals = cropped[mask[:h, :w]] + if vals.size > SECTION8_DIFF_DIST_MAX_SAMPLES: + idx = rng.choice(vals.size, SECTION8_DIFF_DIST_MAX_SAMPLES, replace=False) + vals = vals[idx] + out[key] = vals.astype(np.float32) + return out + @staticmethod def _compute_nc_ratios(score_a: dict, score_b: dict) -> dict: """Per-scale A/B ratio of noise-corrected scores; {} if either side is empty @@ -496,12 +550,14 @@ def _log_analysis(self, norm_a, norm_b, sigmas, label_a, label_b, display_roi=None, crosshair=None, - mask_neb_shared=None, mask_bg_a=None, mask_bg_b=None) -> tuple[dict, dict]: + mask_neb_shared=None, mask_bg_a=None, mask_bg_b=None, + mask_bg_shared=None, diff_dist_rng=None) -> tuple[dict, dict]: figures = {} partial: dict = { "log_nc_score_a": {}, "log_nc_score_b": {}, "log_nc_noise_a": {}, "log_nc_noise_b": {}, "panels": {}, + "diff_dist": {}, } single = norm_b is None for sigma in sigmas: @@ -517,11 +573,15 @@ def _log_analysis(self, norm_a, norm_b, sigmas, partial["log_nc_score_b"][sigma] = nc_b partial["log_nc_noise_b"][sigma] = noise_b + diff = (log_a - log_b).astype(np.float32) if log_b is not None else None partial["panels"][f"log_{sigma}"] = { "a": log_a.astype(np.float32), "b": log_b.astype(np.float32) if log_b is not None else None, - "diff": (log_a - log_b).astype(np.float32) if log_b is not None else None, + "diff": diff, } + if diff is not None: + partial["diff_dist"][f"log_{sigma}"] = self._diff_distribution( + diff, mask_neb_shared, mask_bg_shared, diff_dist_rng) if not single and noise_a and noise_b: partial["panels"][f"nrm_log_{sigma}"] = { "a": (log_a / noise_a).astype(np.float32), @@ -576,7 +636,8 @@ def _gradient_analysis(self, norm_a, norm_b, sigmas, label_a, label_b, display_roi=None, crosshair=None, - mask_neb_shared=None, mask_bg_a=None, mask_bg_b=None) -> tuple[dict, dict]: + mask_neb_shared=None, mask_bg_a=None, mask_bg_b=None, + mask_bg_shared=None, diff_dist_rng=None) -> tuple[dict, dict]: """G = |gradient| at Gaussian scale sigma (first spatial derivative magnitude). Reuses the LOG_SIGMAS scale set so gradient and |LoG| are directly comparable at identical spatial scales. Structured identically to _log_analysis.""" @@ -585,6 +646,7 @@ def _gradient_analysis(self, norm_a, norm_b, sigmas, "gm_nc_score_a": {}, "gm_nc_score_b": {}, "gm_nc_noise_a": {}, "gm_nc_noise_b": {}, "panels": {}, + "diff_dist": {}, } single = norm_b is None for sigma in sigmas: @@ -600,11 +662,15 @@ def _gradient_analysis(self, norm_a, norm_b, sigmas, partial["gm_nc_score_b"][sigma] = nc_b partial["gm_nc_noise_b"][sigma] = noise_b + diff = (gm_a - gm_b).astype(np.float32) if gm_b is not None else None partial["panels"][f"gradient_{sigma}"] = { "a": gm_a.astype(np.float32), "b": gm_b.astype(np.float32) if gm_b is not None else None, - "diff": (gm_a - gm_b).astype(np.float32) if gm_b is not None else None, + "diff": diff, } + if diff is not None: + partial["diff_dist"][f"gradient_{sigma}"] = self._diff_distribution( + diff, mask_neb_shared, mask_bg_shared, diff_dist_rng) if not single and noise_a and noise_b: partial["panels"][f"nrm_gradient_{sigma}"] = { "a": (gm_a / noise_a).astype(np.float32), @@ -659,7 +725,8 @@ def _wavelet_analysis(self, norm_a, norm_b, wavelet, levels, label_a, label_b, display_roi=None, crosshair=None, - mask_neb_shared=None, mask_bg_a=None, mask_bg_b=None) -> tuple[dict, dict]: + mask_neb_shared=None, mask_bg_a=None, mask_bg_b=None, + mask_bg_shared=None, diff_dist_rng=None) -> tuple[dict, dict]: figures = {} partial: dict = { "sigma_noise_a": None, "sigma_noise_b": None, @@ -667,6 +734,7 @@ def _wavelet_analysis(self, norm_a, norm_b, wavelet, levels, "wavelet_nc_score_a": {}, "wavelet_nc_score_b": {}, "wavelet_nc_noise_a": {}, "wavelet_nc_noise_b": {}, "panels": {}, + "diff_dist": {}, } single = norm_b is None @@ -712,11 +780,15 @@ def _wavelet_analysis(self, norm_a, norm_b, wavelet, levels, continue # display/panels only for levels 2-3, unchanged from prior behaviour display_level = human_level + diff = (rec_a - rec_b).astype(np.float32) if rec_b is not None else None partial["panels"][f"wavelet_{display_level}"] = { "a": rec_a.astype(np.float32), "b": rec_b.astype(np.float32) if rec_b is not None else None, - "diff": (rec_a - rec_b).astype(np.float32) if rec_b is not None else None, + "diff": diff, } + if diff is not None: + partial["diff_dist"][f"wavelet_{display_level}"] = self._diff_distribution( + diff, mask_neb_shared, mask_bg_shared, diff_dist_rng) if not single and noise_a and noise_b: partial["panels"][f"nrm_wavelet_{display_level}"] = { "a": (rec_a / noise_a).astype(np.float32), @@ -800,7 +872,8 @@ def _reconstruct_level(self, coeffs, target_coeff_idx: int, def _weber_analysis(self, norm_a, norm_b, kernel_sizes, label_a, label_b, display_roi=None, - mask_neb_shared=None, mask_bg_a=None, mask_bg_b=None) -> tuple[dict, dict]: + mask_neb_shared=None, mask_bg_a=None, mask_bg_b=None, + mask_bg_shared=None, diff_dist_rng=None) -> tuple[dict, dict]: figures = {} partial: dict = { "weber_contrast_a": {}, @@ -808,6 +881,7 @@ def _weber_analysis(self, norm_a, norm_b, kernel_sizes, "weber_nc_score_a": {}, "weber_nc_score_b": {}, "weber_nc_noise_a": {}, "weber_nc_noise_b": {}, "panels": {}, + "diff_dist": {}, } single = norm_b is None for ks in kernel_sizes: @@ -828,11 +902,15 @@ def _weber_analysis(self, norm_a, norm_b, kernel_sizes, partial["weber_nc_score_b"][ks] = nc_b partial["weber_nc_noise_b"][ks] = noise_b + diff = (wc_a - wc_b).astype(np.float32) if wc_b is not None else None partial["panels"][f"weber_{ks}px"] = { "a": wc_a.astype(np.float32), "b": wc_b.astype(np.float32) if wc_b is not None else None, - "diff": (wc_a - wc_b).astype(np.float32) if wc_b is not None else None, + "diff": diff, } + if diff is not None: + partial["diff_dist"][f"weber_{ks}px"] = self._diff_distribution( + diff, mask_neb_shared, mask_bg_shared, diff_dist_rng) if not single and noise_a and noise_b: partial["panels"][f"nrm_weber_{ks}px"] = { "a": (wc_a / noise_a).astype(np.float32), diff --git a/core/models.py b/core/models.py index 8c03695..0974772 100644 --- a/core/models.py +++ b/core/models.py @@ -40,6 +40,7 @@ EPSF_MAX_STARS = 600 # maximum candidate stars passed to EPSFBuilder; limits computation time SECTION8_BORDER_CROP_FRACTION = 0.05 # fraction of each image dimension cropped from perimeter in Section 8 display maps SECTION8_ANALYSIS_CMAP = "viridis" # colormap for Section 8 A/B analysis map panels (std, LoG, wavelet) +SECTION8_DIFF_DIST_MAX_SAMPLES = 64000 # per masked population, per scale — caps violin/KDE cost on full-res diff maps PSF_SPATIAL_MAP_SIZE = 150 # px; long-axis resolution of FWHM / eccentricity spatial maps PSF_SPATIAL_MAP_SMOOTH_SIGMA = 5.0 # Gaussian smoothing sigma (px) applied to spatial maps before display diff --git a/report/report_builder.py b/report/report_builder.py index 96bdd27..6ca4d86 100644 --- a/report/report_builder.py +++ b/report/report_builder.py @@ -428,6 +428,133 @@ def _draw_boxwhisker(ax, vals_list): return img_html, caption_html +# Display order and labels for _spatial_diff_distributions_figure. Keys must match +# analysis/image_filters.py's partial["diff_dist"] keys exactly (same keys used for +# partial["panels"], see SpatialDetailAnalyzer._std_analysis/_log_analysis/etc.). +_SPATIAL_DIFF_DIST_ROWS = [ + ("original", "Original (normalised image)"), + ("std_3px", "Local σ — 3 px"), + ("std_5px", "Local σ — 5 px"), + ("std_10px", "Local σ — 10 px"), + ("log_1.5", "|LoG| — σ=1.5 px"), + ("log_3.0", "|LoG| — σ=3 px"), + ("log_6.0", "|LoG| — σ=6 px"), + ("gradient_1.5", "Gradient |G| — σ=1.5 px"), + ("gradient_3.0", "Gradient |G| — σ=3 px"), + ("gradient_6.0", "Gradient |G| — σ=6 px"), + ("wavelet_2", "Wavelet — level 2"), + ("wavelet_3", "Wavelet — level 3"), + ("weber_3px", "Weber contrast — 3 px"), + ("weber_5px", "Weber contrast — 5 px"), + ("weber_9px", "Weber contrast — 9 px"), +] + + +def _spatial_diff_distributions_figure(diff_dist: dict) -> tuple[str, str]: + """Combined figure: one row per Section 8 calculation, each row showing a + nebula-region violin and a background-region violin of the A-B diff pixel + values, with an IQR box-plot overlay (median magenta, IQR cyan — same styling + as the ePSF section's _psf_distributions_figure). Always violin+box, never + strip/swarm — diff populations here are always high-N. + + Returns (img_html, caption_html), or ("", "") if no row has enough data + (including single-image mode, where diff_dist is empty). + """ + import seaborn as sns + import pandas as pd + + rows = [(k, label) for k, label in _SPATIAL_DIFF_DIST_ROWS if k in diff_dist] + has_data = any( + diff_dist[k]["nebula"].size >= 3 and diff_dist[k]["background"].size >= 3 + for k, _ in rows + ) + if not rows or not has_data: + return "", "" + + fig, axes = plt.subplots(len(rows), 1, figsize=(7, 1.3 * len(rows) + 1)) + fig.subplots_adjust(hspace=0.65, left=0.22, right=0.97, top=0.97, bottom=0.04) + palette = {"Nebula": "steelblue", "Background": "tomato"} + order = ["Nebula", "Background"] + + def _draw_boxwhisker(ax, vals_list): + for i, vals in enumerate(vals_list): + ax.boxplot( + [vals], positions=[i], vert=False, + widths=0.45, zorder=5, + patch_artist=True, + manage_ticks=False, + boxprops=dict(facecolor="none", edgecolor="#00e5ff", linewidth=1.5, alpha=0.9), + medianprops=dict(color="magenta", linewidth=2.0, alpha=0.9), + whiskerprops=dict(color="#00e5ff", linewidth=1.5, alpha=0.9), + capprops=dict(color="#00e5ff", linewidth=1.5, alpha=0.9), + flierprops=dict(marker="", visible=False), + ) + + for ax, (key, title) in zip(np.atleast_1d(axes), rows): + neb = diff_dist[key]["nebula"] + bg = diff_dist[key]["background"] + if neb.size < 3 or bg.size < 3: + ax.set_visible(False) + continue + + combined = np.concatenate([neb, bg]) + df = pd.DataFrame({ + "value": combined, + "region": (["Nebula"] * neb.size) + (["Background"] * bg.size), + }) + sns.violinplot(data=df, x="value", y="region", order=order, + palette=palette, inner=None, linewidth=0.8, ax=ax) + _draw_boxwhisker(ax, [neb, bg]) + + # Some detail maps (e.g. Weber contrast, unbounded near dark-sky pixels — + # see 8e methodology) have rare extreme-outlier diffs that stretch the axis + # so far the IQR box becomes an invisible sliver. Clip the *view* to the + # 1st-99th percentile of this row's own pooled data; the boxplot/violin + # values themselves are unaffected, only what's visible is cropped. + lo, hi = np.percentile(combined, [1.0, 99.0]) + if hi > lo: + pad = 0.05 * (hi - lo) + ax.set_xlim(lo - pad, hi + pad) + + ax.axvline(0.0, color="red", linestyle="--", linewidth=1.0, alpha=0.8, zorder=3) + ax.set_title(title, fontsize=8, loc="left", pad=2) + ax.set_xlabel("", fontsize=7) + ax.set_ylabel("", fontsize=7) + ax.tick_params(axis="x", labelsize=7) + ax.tick_params(axis="y", labelsize=7, pad=1) + ax.spines[["top", "right"]].set_visible(False) + + img_html = _img_tag(fig, "spatial_diff_distributions") + + caption_html = ( + '

' + "Pixel-wise A−B difference distributions. " + "Each row is one Section 8 calculation, shown as a " + "violin plot (kernel density estimate of the diff pixel values, " + "randomly subsampled for display) with an IQR box-plot overlay: " + "a cyan box spanning Q1–Q3, " + "a magenta centre line at the " + "median. Nebula = pixels both " + "images classify as nebula; Background " + "= pixels both images classify as background sky — the same shared masks used " + "for the noise-corrected scores in 8b–8f, here shown as full distributions " + "rather than a single median ratio. " + "A red dashed line marks zero " + "(A = B, no difference). " + "How to read it: a Background row centred at zero with a narrow IQR is " + "the expected noise floor; a Nebula row with a similarly narrow, zero-centred " + "distribution means the two filters agree at that scale. A Nebula median shifted " + "away from zero, or an IQR visibly wider than the Background row's, indicates a " + "real structural difference between the filters at that scale rather than noise. " + "Each row's x-axis is independently clipped to its own 1st–99th percentile " + "range so the IQR box stays visible; rows with rare extreme-outlier diffs (e.g. " + "Weber contrast near dark-sky pixels, see 8e) may have a small fraction of the " + "violin's tail extend beyond the visible axis." + "

" + ) + return img_html, caption_html + + def _epsf_stars_cell(psf_metrics: dict) -> str: """Format the number of stars used to build the ePSF.""" n = psf_metrics.get("epsf_n_stars") @@ -3873,6 +4000,12 @@ def _section_spatial(self, ra: AnalysisResult, rb: AnalysisResult) -> str: sm = ra.spatial_metrics or {} figs = sm.get("figures", {}) + dist_img, dist_caption = _spatial_diff_distributions_figure(sm.get("diff_dist", {})) + dist_html = ( + "

8a. Difference Distribution Overview (A−B)

" + dist_img + dist_caption + if dist_img else "" + ) + cr_a = sm.get("contrast_ratios_a", {}) cr_b = sm.get("contrast_ratios_b", {}) @@ -4053,6 +4186,7 @@ def paired_figs_for(img_prefix: str, xs_prefix: str) -> str: 'across different filter bandwidths. Images are shown side-by-side with a shared ' 'colour scale; the third panel shows the difference A−B.', title="Spatial detail maps overview")} +{dist_html} {nc_methodology_box} {nc_empty_note} From 5240b0e533b6d001a4444976fd2dc03e30aeebe7 Mon Sep 17 00:00:00 2001 From: Brent <52629076+brentmantooth@users.noreply.github.com> Date: Tue, 14 Jul 2026 21:37:10 -0400 Subject: [PATCH 2/8] Replace Section 8 diff metric with log-ratio; add mask illustration and correlation plots Section 8 (Spatial Detail) previously summarised A-vs-B agreement as a plain per-pixel difference, which conflates images at different absolute brightness scales and buries the signal this report exists to surface: which filter preserves more structure/contrast. - Replace the A-B diff with log10(|A|/|B|) across every metric family (std, LoG, gradient, wavelet, Weber, original) - panels, diff-panel colour maps, the Section 8a violin plots, and the Report Inspector's exported label all now reflect the ratio metric. - Add a mask illustration figure (translucent nebula/background overlay on Image A) explaining how the shared masks used by the violin and correlation plots are detected. - Add a new Section 8g: per-scale correlation scatter plots (A vs B, 1:1 reference line) for all 14 metric/scale combinations, with full unclipped axis ranges plus dual linear fits (whole-population and tail-restricted) so divergence in the high-detail regime is visible and quantified separately from bulk agreement. Co-Authored-By: Claude Sonnet 5 --- CLAUDE.md | 15 ++ analysis/image_filters.py | 281 +++++++++++++++++++-- core/models.py | 5 +- report/report_builder.py | 111 ++++++-- tests/test_analysis/test_spatial_detail.py | 73 ++++++ 5 files changed, 442 insertions(+), 43 deletions(-) diff --git a/CLAUDE.md b/CLAUDE.md index 8f74ecf..8a047c3 100644 --- a/CLAUDE.md +++ b/CLAUDE.md @@ -59,6 +59,9 @@ synthetic/ | `_extract_cutout(data, xc, yc, radius)` | `gui/halo_dialog.py` | 2r×2r patch centred on star, zero-padded at image edges | | `_annular_rdf(log_data, xc, yc, radius)` | `gui/halo_dialog.py` | 1-px annular mean/std in log10 space; mirrors `HaloAnalyzer._annular_stats` | | `_power_ratio_db(freq_a, rp_a, freq_b, rp_b)` | `report_builder.py` | 10·log10 dB ratio between two radial power curves; returns `None` on missing data or misaligned frequency bins | +| `_log_ratio_map(a, b)` | `analysis/image_filters.py` | Per-pixel `log10(\|A\|/\|B\|)` map with percentile-based epsilon floor and defensive shape crop — the Section 8 replacement for plain `A − B` diff | +| `_plot_mask_illustration(base, mask_neb, mask_bg)` | `analysis/image_filters.py` | Translucent steelblue/tomato mask overlay on a grayscale base image | +| `_plot_metric_correlation(map_a, map_b, mask_neb, mask_bg, ...)` | `analysis/image_filters.py` | 1×2 masked-region scatter (A vs B) with a 1:1 line plus overall + tail-restricted dual linear fits | --- @@ -239,6 +242,16 @@ When adding a new A-vs-B ratio curve to a report figure (precedent: `_power_rati guaranteed whenever an auto-selected ROI is involved (`_extract_roi` in `analysis/power_spectrum.py` computes `N` independently per image when no explicit ROI is set). Degrade gracefully — return `None` / skip the curve — rather than crash. +- **Epsilon-flooring a per-pixel ratio map needs a percentile, not a raw minimum.** + `_power_ratio_db`'s `positive.min() * 0.01` floor is fine for small 1-D arrays + (frequency bins, ~10²–10³ samples) but fragile at per-pixel map scale (10⁵–10⁷ + samples): the minimum order statistic over that many samples can be pathologically + tiny and let one spurious pixel dominate the ratio's dynamic range. `_log_ratio_map` + (`analysis/image_filters.py`) instead floors both operands at a low percentile + (`SECTION8_LOGRATIO_EPS_PERCENTILE`, default 1st) of the pooled positive `|A|,|B|` + values — same tool as the existing display-clipping precedent + (`_plot_side_by_side`'s `np.percentile(arr, 0.5)`), applied to the epsilon floor + instead of just the color scale. --- @@ -337,6 +350,8 @@ pytest tests/ --cov=analysis,core,synthetic,report --cov-report=html | Adding a float64 cast in analysis code | Don't. All image data is float32 after `AstroImage.load()`. The only float64 exception is `astroalign` in `gui/analysis_thread.py`. Redundant float64 casts waste memory and defeat the float32 performance gains. | | Mixed float32/float64 arithmetic silently widens to float64 | NumPy upcasts when operands differ (e.g. `float32_array - float64_scalar`). If photutils ever returns a float64 background model, `background_subtracted()` will silently return float64. Guard by adding `.astype(np.float32)` at the end of `background_subtracted()` in `astro_image.py` if this is observed. | | `_section_snr` crashed when SNR metric is unchecked | `_plot_snr_pair` (`report_builder.py`) built a `panels` list filtered to non-`None` entries but never checked whether it was empty before calling `plt.subplots(1, len(panels), ...)` — 0 columns raised `ValueError: Number of columns must be a positive integer, not 0`. Hit whenever SNR is unchecked while another metric (e.g. Power Spectrum) is run. Fixed with an early `if not panels: return None` guard, matching `_plot_radial_overlay`/`_plot_radial_ratio_db`; both call sites already pipe the result through `_img_tag`, which turns `None` into `""`. | +| Whole-population regression slope hides tail-specific divergence | A single OLS fit over an entire A-vs-B scatter is dominated by leverage (distance from the mean), not point count — it can show near-zero correlation even when a sparse high-value tail diverges sharply (or vice versa), because it blends "the bulk agrees" and "the tail diverges" into one ambiguous number. When the signal of interest lives specifically in the tail (Section 8g: "which filter shows more detail in the brightest/most-structured pixels"), fit two lines — an overall fit and a second fit restricted to the top N% by combined magnitude (`_plot_metric_correlation`, `SECTION8_SCATTER_TAIL_PERCENTILE`) — and report slope + R² for each. High tail R² with slope far from 1 is a real, systematic effect; low tail R² means the apparent divergence is mostly noise. | +| New Section 8 panel key doesn't need Report Inspector code changes | `gui/report_inspector.py` is fully generic — driven entirely by a companion `_inspector.npz` (raw float32/uint8 arrays) plus an embedded `catalog_json` built in `report_builder.py::_write_inspector_file`. `_panel_display_name`/`_panel_concept` dynamically parse any `panels` dict key prefix, so a new `SpatialDetailAnalyzer` panel family auto-appears in the inspector with zero inspector-side changes. A genuinely new *visual type* is a different story: the inspector only knows how to `imshow` 2D/RGB arrays (side-by-side or slider-reveal), so scatter-style plots (Section 8g's `corr_*` figures) must stay static-HTML-only unless new inspector canvas code is written. | --- diff --git a/analysis/image_filters.py b/analysis/image_filters.py index d5338f1..9377c3a 100644 --- a/analysis/image_filters.py +++ b/analysis/image_filters.py @@ -10,7 +10,9 @@ import matplotlib.colors as mcolors matplotlib.use("Agg") import matplotlib.pyplot as plt +from matplotlib.patches import Patch from scipy.ndimage import generic_filter, gaussian_filter, gaussian_laplace, gaussian_gradient_magnitude, map_coordinates, zoom, maximum_filter, minimum_filter, median_filter +from scipy.stats import linregress import pywt from core.astro_image import AstroImage @@ -18,7 +20,9 @@ from core.models import (STD_KERNEL_SIZES, LOG_SIGMAS, WAVELET_NAME, WAVELET_LEVELS, WEBER_KERNEL_SIZES, XS_LINE_ALPHA, SECTION8_BORDER_CROP_FRACTION, SECTION8_ANALYSIS_CMAP, - XS_SNR_REGION_WIDTH, SECTION8_DIFF_DIST_MAX_SAMPLES) + XS_SNR_REGION_WIDTH, SECTION8_DIFF_DIST_MAX_SAMPLES, + SECTION8_LOGRATIO_EPS_PERCENTILE, SECTION8_SCATTER_MAX_SAMPLES, + SECTION8_SCATTER_TAIL_PERCENTILE) MAX_DIM_FOR_STD = 2048 # downsample to this before generic_filter (performance) _DISPLAY_SMOOTH_SIGMA = 1.0 # applied to maps before plotting; does NOT affect metrics @@ -142,7 +146,7 @@ def _clip01(v): return max(0.0, min(1.0, v)) # map is computed from (mean-normalised, ROI-cropped if a ROI was given) as # its own panel, so the Report Inspector can show the source image content # for a map alongside the map itself, pixel-aligned to the same crop. - original_diff = (analysis_a - analysis_b).astype(np.float32) if analysis_b is not None else None + original_diff = self._log_ratio_map(analysis_a, analysis_b) if analysis_b is not None else None result["panels"]["original"] = { "a": analysis_a.astype(np.float32), "b": analysis_b.astype(np.float32) if analysis_b is not None else None, @@ -152,6 +156,13 @@ def _clip01(v): return max(0.0, min(1.0, v)) result["diff_dist"]["original"] = self._diff_distribution( original_diff, mask_neb_shared, mask_bg_shared, diff_dist_rng) + # Mask illustration: only meaningful in two-image mode, mirroring the + # violin plots' own empty-in-single-image-mode behaviour. + if mask_neb_shared is not None: + mask_fig = self._plot_mask_illustration( + result["panels"]["original"]["a"], mask_neb_shared, mask_bg_shared) + figures["mask_illustration"] = fig_to_b64(mask_fig, dpi=150) + _label_b = image_b.label if image_b is not None else None # 1-5. Local std, LoG, wavelet, Weber, gradient — all read norm_a/norm_b with no @@ -389,7 +400,7 @@ def _std_analysis(self, norm_a, norm_b, partial["std_nc_score_b"][ks] = nc_b partial["std_nc_noise_b"][ks] = noise_b - diff = (std_a - std_b).astype(np.float32) if not single else None + diff = self._log_ratio_map(std_a, std_b) if not single else None partial["panels"][f"std_{ks}px"] = { "a": std_a.astype(np.float32), "b": std_b.astype(np.float32) if std_b is not None else None, @@ -398,6 +409,12 @@ def _std_analysis(self, norm_a, norm_b, if diff is not None: partial["diff_dist"][f"std_{ks}px"] = self._diff_distribution( diff, mask_neb_shared, mask_bg_shared, diff_dist_rng) + if not single: + corr_fig = self._plot_metric_correlation( + std_a, std_b, mask_neb_shared, mask_bg_shared, + label_a, label_b, f"Local σ (kernel {ks}px)", diff_dist_rng) + if corr_fig is not None: + figures[f"corr_std_{ks}px"] = corr_fig if not single and noise_a and noise_b: partial["panels"][f"nrm_std_{ks}px"] = { "a": (std_a / noise_a).astype(np.float32), @@ -411,7 +428,7 @@ def _std_analysis(self, norm_a, norm_b, self._crop_border(std_b, SECTION8_BORDER_CROP_FRACTION), f"Local σ — kernel {ks}px — {label_a}", f"Local σ — kernel {ks}px — {label_b}", - diff_title=f"Diff (A−B), kernel {ks}px", + diff_title=f"Log ratio (A/B), kernel {ks}px", cmap=SECTION8_ANALYSIS_CMAP, nonlinear_norm=True, display_roi=None, @@ -431,7 +448,7 @@ def _std_analysis(self, norm_a, norm_b, self._crop_border(std_b / noise_b, SECTION8_BORDER_CROP_FRACTION), f"Local σ (× noise floor) — kernel {ks}px — {label_a}", f"Local σ (× noise floor) — kernel {ks}px — {label_b}", - diff_title=f"Diff (A−B), noise-normalised, kernel {ks}px", + diff_title=f"Log ratio (A/B), noise-normalised, kernel {ks}px", cmap=SECTION8_ANALYSIS_CMAP, display_roi=None, ) @@ -503,17 +520,52 @@ def _nc_score(self, detail_map: np.ndarray, return None, None return float(bn.median(neb_vals)) / noise_floor, noise_floor + @staticmethod + def _log_ratio_map(a: np.ndarray, b: np.ndarray) -> np.ndarray: + """Per-pixel log10(|A|/|B|) map, replacing plain A-B difference. + + Sign is discarded (via abs()) before the ratio so the result is always + well-defined even for map families that can go negative (original + background-subtracted flux, wavelet band-pass reconstructions) — this + turns the comparison into "which image shows more structure/contrast at + this scale", not "which image is signed-brighter". Non-negative families + (std, |LoG|, gradient magnitude, Weber contrast) are unaffected by the + abs() since they're already >= 0. + + Epsilon-floors both operands using a low percentile (not a raw minimum — + a raw minimum over a multi-megapixel array can be pathologically tiny and + let a single spurious pixel dominate the log-ratio's dynamic range) of the + pooled positive values from both inputs, mirroring report_builder.py's + _power_ratio_db epsilon pattern but adapted for per-pixel map sizes. + + A and B are defensively cropped to their common shape first — two-image + analysis can reach here with mismatched shapes when astroalign + registration fails and analysis proceeds unaligned. + """ + abs_a, abs_b = np.abs(a), np.abs(b) + h = min(abs_a.shape[0], abs_b.shape[0]) + w = min(abs_a.shape[1], abs_b.shape[1]) + abs_a, abs_b = abs_a[:h, :w], abs_b[:h, :w] + + positive = np.concatenate([abs_a[abs_a > 0].ravel(), abs_b[abs_b > 0].ravel()]) + eps = max(float(np.percentile(positive, SECTION8_LOGRATIO_EPS_PERCENTILE)), 1e-12) \ + if positive.size > 0 else 1e-12 + + ratio = np.clip(abs_a, eps, None) / np.clip(abs_b, eps, None) + return np.log10(ratio).astype(np.float32) + @staticmethod def _diff_distribution(diff_map: np.ndarray, mask_neb_shared: np.ndarray | None, mask_bg_shared: np.ndarray | None, rng: np.random.Generator) -> dict: - """Random-subsampled A-B diff pixel populations for nebula vs background. + """Random-subsampled log10(|A|/|B|) ratio pixel populations for nebula vs + background. - Returns {"nebula": ndarray, "background": ndarray} (float32, signed diff - values, up to SECTION8_DIFF_DIST_MAX_SAMPLES each). Either array is empty - if the corresponding mask is unavailable (single-image mode) or selects - zero pixels. + Returns {"nebula": ndarray, "background": ndarray} (float32, signed + log-ratio values — 0 means A=B — up to SECTION8_DIFF_DIST_MAX_SAMPLES + each). Either array is empty if the corresponding mask is unavailable + (single-image mode) or selects zero pixels. """ out = {"nebula": np.empty(0, dtype=np.float32), "background": np.empty(0, dtype=np.float32)} @@ -573,7 +625,7 @@ def _log_analysis(self, norm_a, norm_b, sigmas, partial["log_nc_score_b"][sigma] = nc_b partial["log_nc_noise_b"][sigma] = noise_b - diff = (log_a - log_b).astype(np.float32) if log_b is not None else None + diff = self._log_ratio_map(log_a, log_b) if log_b is not None else None partial["panels"][f"log_{sigma}"] = { "a": log_a.astype(np.float32), "b": log_b.astype(np.float32) if log_b is not None else None, @@ -582,6 +634,12 @@ def _log_analysis(self, norm_a, norm_b, sigmas, if diff is not None: partial["diff_dist"][f"log_{sigma}"] = self._diff_distribution( diff, mask_neb_shared, mask_bg_shared, diff_dist_rng) + if not single: + corr_fig = self._plot_metric_correlation( + log_a, log_b, mask_neb_shared, mask_bg_shared, + label_a, label_b, f"|LoG| (σ={sigma}px)", diff_dist_rng) + if corr_fig is not None: + figures[f"corr_log_{sigma}"] = corr_fig if not single and noise_a and noise_b: partial["panels"][f"nrm_log_{sigma}"] = { "a": (log_a / noise_a).astype(np.float32), @@ -595,7 +653,7 @@ def _log_analysis(self, norm_a, norm_b, sigmas, self._crop_border(log_b, SECTION8_BORDER_CROP_FRACTION), f"|LoG| σ={sigma}px — {label_a}", f"|LoG| σ={sigma}px — {label_b}", - diff_title=f"LoG diff (A−B), σ={sigma}px", + diff_title=f"|LoG| log-ratio (A/B), σ={sigma}px", cmap=SECTION8_ANALYSIS_CMAP, nonlinear_norm=True, display_roi=None, @@ -615,7 +673,7 @@ def _log_analysis(self, norm_a, norm_b, sigmas, self._crop_border(log_b / noise_b, SECTION8_BORDER_CROP_FRACTION), f"|LoG| (× noise floor) σ={sigma}px — {label_a}", f"|LoG| (× noise floor) σ={sigma}px — {label_b}", - diff_title=f"Diff (A−B), noise-normalised, σ={sigma}px", + diff_title=f"Log ratio (A/B), noise-normalised, σ={sigma}px", cmap=SECTION8_ANALYSIS_CMAP, display_roi=None, ) @@ -662,7 +720,7 @@ def _gradient_analysis(self, norm_a, norm_b, sigmas, partial["gm_nc_score_b"][sigma] = nc_b partial["gm_nc_noise_b"][sigma] = noise_b - diff = (gm_a - gm_b).astype(np.float32) if gm_b is not None else None + diff = self._log_ratio_map(gm_a, gm_b) if gm_b is not None else None partial["panels"][f"gradient_{sigma}"] = { "a": gm_a.astype(np.float32), "b": gm_b.astype(np.float32) if gm_b is not None else None, @@ -671,6 +729,12 @@ def _gradient_analysis(self, norm_a, norm_b, sigmas, if diff is not None: partial["diff_dist"][f"gradient_{sigma}"] = self._diff_distribution( diff, mask_neb_shared, mask_bg_shared, diff_dist_rng) + if not single: + corr_fig = self._plot_metric_correlation( + gm_a, gm_b, mask_neb_shared, mask_bg_shared, + label_a, label_b, f"Gradient |G| (σ={sigma}px)", diff_dist_rng) + if corr_fig is not None: + figures[f"corr_gradient_{sigma}"] = corr_fig if not single and noise_a and noise_b: partial["panels"][f"nrm_gradient_{sigma}"] = { "a": (gm_a / noise_a).astype(np.float32), @@ -684,7 +748,7 @@ def _gradient_analysis(self, norm_a, norm_b, sigmas, self._crop_border(gm_b, SECTION8_BORDER_CROP_FRACTION), f"Gradient |G| σ={sigma}px — {label_a}", f"Gradient |G| σ={sigma}px — {label_b}", - diff_title=f"Gradient diff (A−B), σ={sigma}px", + diff_title=f"Gradient log-ratio (A/B), σ={sigma}px", cmap=SECTION8_ANALYSIS_CMAP, nonlinear_norm=True, display_roi=None, @@ -704,7 +768,7 @@ def _gradient_analysis(self, norm_a, norm_b, sigmas, self._crop_border(gm_b / noise_b, SECTION8_BORDER_CROP_FRACTION), f"Gradient (× noise floor) σ={sigma}px — {label_a}", f"Gradient (× noise floor) σ={sigma}px — {label_b}", - diff_title=f"Diff (A−B), noise-normalised, σ={sigma}px", + diff_title=f"Log ratio (A/B), noise-normalised, σ={sigma}px", cmap=SECTION8_ANALYSIS_CMAP, display_roi=None, ) @@ -780,7 +844,7 @@ def _wavelet_analysis(self, norm_a, norm_b, wavelet, levels, continue # display/panels only for levels 2-3, unchanged from prior behaviour display_level = human_level - diff = (rec_a - rec_b).astype(np.float32) if rec_b is not None else None + diff = self._log_ratio_map(rec_a, rec_b) if rec_b is not None else None partial["panels"][f"wavelet_{display_level}"] = { "a": rec_a.astype(np.float32), "b": rec_b.astype(np.float32) if rec_b is not None else None, @@ -789,6 +853,15 @@ def _wavelet_analysis(self, norm_a, norm_b, wavelet, levels, if diff is not None: partial["diff_dist"][f"wavelet_{display_level}"] = self._diff_distribution( diff, mask_neb_shared, mask_bg_shared, diff_dist_rng) + if not single: + # Raw signed reconstructions (not abs()) — complementary to the + # sign-discarding log-ratio map, shows whether band-pass detail + # flips sign between the two filters at a given pixel. + corr_fig = self._plot_metric_correlation( + rec_a, rec_b, mask_neb_shared, mask_bg_shared, + label_a, label_b, f"Wavelet level {display_level}", diff_dist_rng) + if corr_fig is not None: + figures[f"corr_wavelet_{display_level}"] = corr_fig if not single and noise_a and noise_b: partial["panels"][f"nrm_wavelet_{display_level}"] = { "a": (rec_a / noise_a).astype(np.float32), @@ -801,7 +874,7 @@ def _wavelet_analysis(self, norm_a, norm_b, wavelet, levels, self._crop_border(rec_b, SECTION8_BORDER_CROP_FRACTION), f"Wavelet level {display_level} — {label_a}", f"Wavelet level {display_level} — {label_b}", - diff_title=f"Level {display_level} diff (A−B)", + diff_title=f"Level {display_level} log-ratio (A/B)", cmap=SECTION8_ANALYSIS_CMAP, display_roi=None, ) @@ -819,7 +892,7 @@ def _wavelet_analysis(self, norm_a, norm_b, wavelet, levels, self._crop_border(rec_b / noise_b, SECTION8_BORDER_CROP_FRACTION), f"Wavelet level {display_level} (× noise floor) — {label_a}", f"Wavelet level {display_level} (× noise floor) — {label_b}", - diff_title=f"Level {display_level} diff (A−B), noise-normalised", + diff_title=f"Level {display_level} log-ratio (A/B), noise-normalised", cmap=SECTION8_ANALYSIS_CMAP, display_roi=None, ) @@ -902,7 +975,7 @@ def _weber_analysis(self, norm_a, norm_b, kernel_sizes, partial["weber_nc_score_b"][ks] = nc_b partial["weber_nc_noise_b"][ks] = noise_b - diff = (wc_a - wc_b).astype(np.float32) if wc_b is not None else None + diff = self._log_ratio_map(wc_a, wc_b) if wc_b is not None else None partial["panels"][f"weber_{ks}px"] = { "a": wc_a.astype(np.float32), "b": wc_b.astype(np.float32) if wc_b is not None else None, @@ -911,6 +984,12 @@ def _weber_analysis(self, norm_a, norm_b, kernel_sizes, if diff is not None: partial["diff_dist"][f"weber_{ks}px"] = self._diff_distribution( diff, mask_neb_shared, mask_bg_shared, diff_dist_rng) + if not single: + corr_fig = self._plot_metric_correlation( + wc_a, wc_b, mask_neb_shared, mask_bg_shared, + label_a, label_b, f"Weber contrast (kernel {ks}px)", diff_dist_rng) + if corr_fig is not None: + figures[f"corr_weber_{ks}px"] = corr_fig if not single and noise_a and noise_b: partial["panels"][f"nrm_weber_{ks}px"] = { "a": (wc_a / noise_a).astype(np.float32), @@ -924,7 +1003,7 @@ def _weber_analysis(self, norm_a, norm_b, kernel_sizes, self._crop_border(wc_b, SECTION8_BORDER_CROP_FRACTION), f"Weber contrast — kernel {ks}px — {label_a}", f"Weber contrast — kernel {ks}px — {label_b}", - diff_title=f"Weber diff (A−B), kernel {ks}px", + diff_title=f"Weber log-ratio (A/B), kernel {ks}px", cmap=SECTION8_ANALYSIS_CMAP, nonlinear_norm=True, # unbounded output; PowerNorm compresses dynamic range display_roi=None, # _crop_border already applied @@ -944,7 +1023,7 @@ def _weber_analysis(self, norm_a, norm_b, kernel_sizes, self._crop_border(wc_b / noise_b, SECTION8_BORDER_CROP_FRACTION), f"Weber (× noise floor) — kernel {ks}px — {label_a}", f"Weber (× noise floor) — kernel {ks}px — {label_b}", - diff_title=f"Diff (A−B), noise-normalised, kernel {ks}px", + diff_title=f"Log ratio (A/B), noise-normalised, kernel {ks}px", cmap=SECTION8_ANALYSIS_CMAP, display_roi=None, ) @@ -1009,10 +1088,8 @@ def _plot_side_by_side(self, arr_a: np.ndarray, arr_b: np.ndarray, # Sqrt (PowerNorm gamma=0.5) compresses bright stars, reveals faint nebula norm = mcolors.PowerNorm(gamma=0.5, vmin=vmin, vmax=vmax) if nonlinear_norm else None - # Difference panel (computed before any possible shape mismatch) - h_min = min(arr_a.shape[0], arr_b.shape[0]) - w_min = min(arr_a.shape[1], arr_b.shape[1]) - diff = arr_a[:h_min, :w_min] - arr_b[:h_min, :w_min] + # Log-ratio panel (helper handles any shape mismatch defensively) + diff = self._log_ratio_map(arr_a, arr_b) # Symmetric about zero so the "bwr" midpoint (white) always means no difference. d_max = float(np.percentile(np.abs(diff), 99.5)) or 1.0 @@ -1074,6 +1151,48 @@ def _plot_single(self, arr_a: np.ndarray, title_a: str, fig.colorbar(im, ax=ax, fraction=0.046, pad=0.04) return fig + def _plot_mask_illustration(self, base: np.ndarray, mask_neb: np.ndarray, + mask_bg: np.ndarray, alpha: float = 0.45) -> plt.Figure: + """Image A shown with a translucent nebula/background mask overlay. + + Uses the same steelblue/tomato color convention as the Section 8a violin + plots (report_builder.py's palette = {"Nebula": "steelblue", + "Background": "tomato"}) so the two figures read as one visual language. + Unclassified pixels (neither mask) are left plain grayscale. mask_neb/ + mask_bg are expected to already be the two-image shared classification + (mask_neb_a & mask_neb_b) — the exact masks that feed the violin plots + and correlation scatter plots — so this figure depicts what those plots + are actually gated on. + """ + gray = self._stretch_for_display(base) + h = min(gray.shape[0], mask_neb.shape[0], mask_bg.shape[0]) + w = min(gray.shape[1], mask_neb.shape[1], mask_bg.shape[1]) + gray = gray[:h, :w] + mask_neb = mask_neb[:h, :w] + mask_bg = mask_bg[:h, :w] + + rgb = np.stack([gray, gray, gray], axis=-1) + neb_color = np.array(mcolors.to_rgb("steelblue")) + bg_color = np.array(mcolors.to_rgb("tomato")) + rgb[mask_neb] = (1 - alpha) * rgb[mask_neb] + alpha * neb_color + rgb[mask_bg] = (1 - alpha) * rgb[mask_bg] + alpha * bg_color + + aspect_ratio = h / max(w, 1) + panel_w = 8.0 + fig, ax = plt.subplots(figsize=(panel_w, panel_w * aspect_ratio + 1.0)) + ax.imshow(rgb, origin="upper", interpolation="nearest", aspect="equal") + ax.axis("off") + ax.set_title("Nebula / background mask regions (shared A∩B classification), " + "shown on Image A", fontsize=10) + legend_handles = [ + Patch(facecolor="steelblue", edgecolor="none", alpha=0.7, label="Nebula"), + Patch(facecolor="tomato", edgecolor="none", alpha=0.7, label="Background"), + Patch(facecolor="0.5", edgecolor="none", label="Unclassified"), + ] + ax.legend(handles=legend_handles, loc="lower right", fontsize=8, framealpha=0.8) + fig.tight_layout() + return fig + def _plot_snr_bars(self, snr_a: dict, snr_b: dict, label_a: str, label_b: str | None, levels: int) -> plt.Figure: fig, ax = plt.subplots(figsize=(7, 4)) @@ -1184,6 +1303,116 @@ def _plot_cross_section(pos: np.ndarray, prof_a: np.ndarray, prof_b: np.ndarray, ax1.set_title(title, fontsize=10) return fig + @staticmethod + def _plot_metric_correlation(map_a: np.ndarray, map_b: np.ndarray, + mask_neb_shared: np.ndarray, mask_bg_shared: np.ndarray, + label_a: str, label_b: str, metric_title: str, + rng: np.random.Generator, + max_samples: int = SECTION8_SCATTER_MAX_SAMPLES) -> plt.Figure | None: + """1x2 correlation scatter (Nebula | Background): metric value in A (y) + vs. metric value in B (x), with a dashed 1:1 reference line plus two + linear fits — an overall OLS fit over all points, and a second fit + restricted to the top (100 - SECTION8_SCATTER_TAIL_PERCENTILE)% of + points by combined magnitude ("tail fit"). The tail fit isolates the + high-detail regime's slope from the bulk's, since a whole-population + slope is dominated by leverage rather than point count and conflates + the two regimes into one ambiguous number. + + Axis limits reflect the FULL pooled population range per subplot (not + percentile-clipped like the Section 8a violin plots) so upper-tail + divergence from the 1:1 line — the signal this plot exists to surface — + stays visible. Point clouds are randomly subsampled (up to max_samples) + for render cost only; the axis range and both fits are always computed + from the full, unsampled population. Returns None if both subplots have + too few points. + """ + h = min(map_a.shape[0], map_b.shape[0], mask_neb_shared.shape[0], mask_bg_shared.shape[0]) + w = min(map_a.shape[1], map_b.shape[1], mask_neb_shared.shape[1], mask_bg_shared.shape[1]) + map_a = map_a[:h, :w] + map_b = map_b[:h, :w] + mask_neb_shared = mask_neb_shared[:h, :w] + mask_bg_shared = mask_bg_shared[:h, :w] + + # Dark-mode-aware reference-line color (project convention). + is_dark = matplotlib.rcParams.get("figure.facecolor", "white") not in ("white", "#ffffff", 1.0) + orig_color = "white" if is_dark else "black" + + fig, axes = plt.subplots(1, 2, figsize=(9, 4.5)) + any_data = False + for ax, region_name, mask in zip(axes, ("Nebula", "Background"), + (mask_neb_shared, mask_bg_shared)): + a_vals = map_a[mask] + b_vals = map_b[mask] + n = a_vals.size + if n < 3: + ax.set_visible(False) + continue + any_data = True + + lo = float(min(a_vals.min(), b_vals.min())) + hi = float(max(a_vals.max(), b_vals.max())) + pad = 0.1 * (hi - lo) if hi > lo else 1.0 + lo -= pad + hi += pad + + if n > max_samples: + idx = rng.choice(n, max_samples, replace=False) + a_plot, b_plot = a_vals[idx], b_vals[idx] + else: + a_plot, b_plot = a_vals, b_vals + + # Linear fits computed on the FULL population (not the render subsample) + # for statistical accuracy. Two fits, not one: an overall OLS fit over + # all points is dominated by leverage, not point count, so it conflates + # "the bulk agrees near 1:1" with "the tail diverges" into one ambiguous + # number. A second fit restricted to the top (100 - + # SECTION8_SCATTER_TAIL_PERCENTILE)% of points by combined magnitude + # (A + B, a proxy for "how far out along the diagonal") isolates the + # high-detail regime specifically — that slope is the more direct + # answer to "does one filter show disproportionately more detail/ + # contrast in the brightest/most-structured regions". + overall_fit = linregress(b_vals, a_vals) + m_all, c_all, r_all = overall_fit.slope, overall_fit.intercept, overall_fit.rvalue + + combined = a_vals + b_vals + tail_thresh = np.percentile(combined, SECTION8_SCATTER_TAIL_PERCENTILE) + tail_mask = combined >= tail_thresh + tail_fit = None + if np.count_nonzero(tail_mask) >= 3: + tf = linregress(b_vals[tail_mask], a_vals[tail_mask]) + tail_fit = (tf.slope, tf.intercept, tf.rvalue) + + ax.scatter(b_plot, a_plot, alpha=0.55, s=8, color="steelblue", + zorder=3, edgecolors="none", rasterized=True) + ax.plot([lo, hi], [lo, hi], color=orig_color, linestyle="--", + linewidth=1.2, zorder=4, label="Slope = 1 (A = B)") + fit_x = np.array([lo, hi]) + ax.plot(fit_x, m_all * fit_x + c_all, color="darkorange", linestyle="-", + linewidth=1.6, zorder=5, + label=f"Fit: y={m_all:.2f}x{c_all:+.2f} (R²={r_all ** 2:.2f})") + if tail_fit is not None: + m_t, c_t, r_t = tail_fit + tail_pct = 100 - SECTION8_SCATTER_TAIL_PERCENTILE + ax.plot(fit_x, m_t * fit_x + c_t, color="crimson", linestyle=":", + linewidth=1.8, zorder=6, + label=f"Tail fit (top {tail_pct:.0f}%): " + f"y={m_t:.2f}x{c_t:+.2f} (R²={r_t ** 2:.2f})") + ax.set_xlim(lo, hi) + ax.set_ylim(lo, hi) + ax.set_aspect("equal") + ax.set_xlabel(f"{metric_title} — {label_b} (x)", fontsize=8) + ax.set_ylabel(f"{metric_title} — {label_a} (y)", fontsize=8) + ax.set_title(f"{region_name} (n={n})", fontsize=9) + ax.tick_params(labelsize=7) + ax.legend(fontsize=6.5, loc="best", labelspacing=0.3) + ax.grid(True, alpha=0.3) + + if not any_data: + plt.close(fig) + return None + fig.tight_layout() + return fig + @staticmethod def _crop_border(arr: np.ndarray, fraction: float) -> np.ndarray: n = max(1, int(min(arr.shape[0], arr.shape[1]) * fraction)) diff --git a/core/models.py b/core/models.py index 0974772..31fd51d 100644 --- a/core/models.py +++ b/core/models.py @@ -40,7 +40,10 @@ EPSF_MAX_STARS = 600 # maximum candidate stars passed to EPSFBuilder; limits computation time SECTION8_BORDER_CROP_FRACTION = 0.05 # fraction of each image dimension cropped from perimeter in Section 8 display maps SECTION8_ANALYSIS_CMAP = "viridis" # colormap for Section 8 A/B analysis map panels (std, LoG, wavelet) -SECTION8_DIFF_DIST_MAX_SAMPLES = 64000 # per masked population, per scale — caps violin/KDE cost on full-res diff maps +SECTION8_DIFF_DIST_MAX_SAMPLES = 1000000 # per masked population, per scale — caps violin/KDE cost on full-res diff maps +SECTION8_LOGRATIO_EPS_PERCENTILE = 1.0 # percentile of pooled positive |A|,|B| values used as the epsilon floor in log10(|A|/|B|) +SECTION8_SCATTER_MAX_SAMPLES = 50000 # per masked population, per scale — caps render cost of Section 8g correlation scatter plots +SECTION8_SCATTER_TAIL_PERCENTILE = 80 # percentile of (A+B) combined magnitude above which points are fit separately as the "tail" regression in Section 8g PSF_SPATIAL_MAP_SIZE = 150 # px; long-axis resolution of FWHM / eccentricity spatial maps PSF_SPATIAL_MAP_SMOOTH_SIGMA = 5.0 # Gaussian smoothing sigma (px) applied to spatial maps before display diff --git a/report/report_builder.py b/report/report_builder.py index 6ca4d86..f520462 100644 --- a/report/report_builder.py +++ b/report/report_builder.py @@ -449,6 +449,10 @@ def _draw_boxwhisker(ax, vals_list): ("weber_9px", "Weber contrast — 9 px"), ] +# Same order/labels as above minus "original" (no kernel scale) — used for the +# Section 8g per-scale correlation scatter plots. +_SPATIAL_CORR_ROWS = [(k, label) for k, label in _SPATIAL_DIFF_DIST_ROWS if k != "original"] + def _spatial_diff_distributions_figure(diff_dist: dict) -> tuple[str, str]: """Combined figure: one row per Section 8 calculation, each row showing a @@ -507,8 +511,8 @@ def _draw_boxwhisker(ax, vals_list): _draw_boxwhisker(ax, [neb, bg]) # Some detail maps (e.g. Weber contrast, unbounded near dark-sky pixels — - # see 8e methodology) have rare extreme-outlier diffs that stretch the axis - # so far the IQR box becomes an invisible sliver. Clip the *view* to the + # see 8e methodology) have rare extreme-outlier log-ratios that stretch the + # axis so far the IQR box becomes an invisible sliver. Clip the *view* to the # 1st-99th percentile of this row's own pooled data; the boxplot/violin # values themselves are unaffected, only what's visible is cropped. lo, hi = np.percentile(combined, [1.0, 99.0]) @@ -528,28 +532,33 @@ def _draw_boxwhisker(ax, vals_list): caption_html = ( '

' - "Pixel-wise A−B difference distributions. " + "Pixel-wise log₁₀(A / B) ratio distributions. " "Each row is one Section 8 calculation, shown as a " - "violin plot (kernel density estimate of the diff pixel values, " + "violin plot (kernel density estimate of the log-ratio pixel values, " "randomly subsampled for display) with an IQR box-plot overlay: " "a cyan box spanning Q1–Q3, " "a magenta centre line at the " "median. Nebula = pixels both " "images classify as nebula; Background " - "= pixels both images classify as background sky — the same shared masks used " - "for the noise-corrected scores in 8b–8f, here shown as full distributions " - "rather than a single median ratio. " + "= pixels both images classify as background sky — the same shared masks " + "illustrated above and used for the noise-corrected scores in 8b–8f, here " + "shown as full distributions rather than a single median ratio. " "A red dashed line marks zero " - "(A = B, no difference). " + "(A = B, equal). Units are log₁₀ — ±0.3 ≈ a " + "2× difference, ±1.0 ≈ a 10× difference. For the " + "Original and Wavelet rows specifically (the only two families that can go " + "negative), sign was discarded before the ratio — these rows compare the " + "magnitude of structure, not signed brightness. " "How to read it: a Background row centred at zero with a narrow IQR is " "the expected noise floor; a Nebula row with a similarly narrow, zero-centred " "distribution means the two filters agree at that scale. A Nebula median shifted " "away from zero, or an IQR visibly wider than the Background row's, indicates a " "real structural difference between the filters at that scale rather than noise. " "Each row's x-axis is independently clipped to its own 1st–99th percentile " - "range so the IQR box stays visible; rows with rare extreme-outlier diffs (e.g. " + "range so the IQR box stays visible; rows with rare extreme-outlier log-ratios (e.g. " "Weber contrast near dark-sky pixels, see 8e) may have a small fraction of the " - "violin's tail extend beyond the visible axis." + "violin's tail extend beyond the visible axis — see 8g for the full, unclipped " + "upper-tail behaviour." "

" ) return img_html, caption_html @@ -1099,7 +1108,7 @@ def _add_options_entry(section: str, name: str, options: dict, if pa_panel.get("diff") is not None: npz_diff = f"sp_{pkey}_diff" _add(npz_diff, pa_panel["diff"]) - sp_opts["Diff (A−B)"] = npz_diff + sp_opts["Log Ratio (A/B)"] = npz_diff _add_options_entry("Spatial Detail", img_set_name, sp_opts, concept=_panel_concept(pkey)) @@ -4000,9 +4009,30 @@ def _section_spatial(self, ra: AnalysisResult, rb: AnalysisResult) -> str: sm = ra.spatial_metrics or {} figs = sm.get("figures", {}) + mask_fig = figs.get("mask_illustration") dist_img, dist_caption = _spatial_diff_distributions_figure(sm.get("diff_dist", {})) + mask_html = ( + "

Nebula / Background Mask Regions

" + + _hires_img_tag(mask_fig, "Mask illustration") + + '

' + "How the masks are detected. Each image is independently classified " + "per-pixel from its own background-subtracted flux and noise RMS: a pixel is " + "Nebula if its background-subtracted " + "value exceeds 2× the image's RMS, and " + "Background if it's below 0.5× the " + "RMS (a top-5% percentile fallback applies if no pixel clears the nebula threshold). " + "Pixels between the two thresholds are left unclassified (plain grayscale). " + "How they're used below. The log-ratio distributions and per-scale " + "correlation plots (8g) use the two-image intersection of these masks — " + "a pixel counts as Nebula only if both images classify it as Nebula (same for " + "Background). This is deliberately conservative: it excludes pixels where one " + f"image's registration, PSF, or local noise disagrees with the other's. Shown on " + f"{ra.label}, the same array as the "Original" row below." + "

" + if mask_fig else "" + ) dist_html = ( - "

8a. Difference Distribution Overview (A−B)

" + dist_img + dist_caption + "

8a. Log-Ratio Distribution & Mask Overview

" + mask_html + dist_img + dist_caption if dist_img else "" ) @@ -4112,6 +4142,47 @@ def paired_figs_for(img_prefix: str, xs_prefix: str) -> str: out += _hires_img_tag(figs[xs_key], xs_key) + "\n" return out + def corr_figs_html() -> str: + """Emit each present correlation scatter figure in _SPATIAL_CORR_ROWS + order (not alphabetic — e.g. std_10px must not sort before std_3px).""" + out = "" + for key, label in _SPATIAL_CORR_ROWS: + fig = figs.get(f"corr_{key}") + if fig: + out += f"

{label}

" + _hires_img_tag(fig, f"corr_{key}") + "\n" + return out + + _corr_figs = corr_figs_html() + corr_section_html = ( + "

8g. Per-Pixel Metric Correlation (A vs. B)

" + + _info_box('Each row plots the raw (pre-ratio) metric value at every shared pixel: ' + f'{ra.label} on the y-axis vs. {rb.label} on the x-axis, split into a ' + 'Nebula panel and a Background panel using the shared masks explained in 8a. ' + 'The black dashed diagonal is the 1:1 line (perfect agreement, ' + 'A = B). Axis limits span the full data range for that panel — ' + 'unlike the 8a violin plots, they are not clipped to a percentile — so ' + 'the behaviour of the upper tail is always visible. ' + '

Two fitted lines are drawn on every panel: a solid ' + 'orange overall fit through every ' + 'point, and a dotted crimson tail fit ' + 'restricted to the brightest/most-structured points (the top fraction by ' + 'combined A+B magnitude — see the legend for the exact percentage). A ' + 'whole-population slope is dominated by leverage (distance from the mean), not ' + 'point count, so it blends "the bulk agrees near 1:1" and "the tail diverges" ' + 'into one ambiguous number. The tail fit\'s slope is the more direct ' + 'answer to which filter shows disproportionately more detail/contrast in the ' + 'brightest, most-structured pixels — a tail slope far from 1 with a ' + 'high R² is strong, consistent evidence of a real difference; a tail slope near ' + '1 (even if the point cloud looks scattered) means the apparent divergence is ' + 'mostly noise, not a systematic effect.

' + 'Point clouds are randomly subsampled for rendering; the axis range, both fits, ' + 'and each panel\'s point count (n) are always computed from the full, unsampled ' + 'population.', + title="Per-pixel correlation methodology") + + _corr_figs + if _corr_figs else "" + ) + has_crosshair = sm.get("crosshair") is not None xs_note = _info_box( 'ℹ Cross-section profiles below are extracted along ' @@ -4184,7 +4255,12 @@ def paired_figs_for(img_prefix: str, xs_prefix: str) -> str: {_info_box('All maps below are computed on mean-signal-normalised data ' '(each image divided by its own mean signal), making them dimensionless and comparable ' 'across different filter bandwidths. Images are shown side-by-side with a shared ' - 'colour scale; the third panel shows the difference A−B.', + 'colour scale; the third panel shows the per-pixel log-ratio ' + 'log10(A / B): red (positive) means A is stronger at that pixel, blue ' + '(negative) means B is stronger, white means the two are equal. ±0.3 ≈ a 2× ' + 'difference, ±1.0 ≈ a 10× difference. For the Original and Wavelet maps — the ' + 'only two families whose values can be negative — the ratio uses |A| and |B|, so it ' + 'compares the magnitude of structure rather than signed brightness.', title="Spatial detail maps overview")} {dist_html} {nc_methodology_box} @@ -4213,7 +4289,7 @@ def paired_figs_for(img_prefix: str, xs_prefix: str) -> str: {xs_note}{paired_figs_for("std_", "xs_std_")}

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.

+The log-ratio 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_")}

8c. Laplacian of Gaussian (LoG) Maps

@@ -4257,7 +4333,8 @@ def paired_figs_for(img_prefix: str, xs_prefix: str) -> str: {paired_figs_for("wavelet_level", "xs_wavelet_level")}

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 log-ratio panel (right) shows where fine structure differs between the two filters +(sign discarded — |A|/|B| — since wavelet reconstructions can be negative).

Noise-normalised (× noise floor) — shared colour scale is a fair A/B comparison.

{figs_for("nrm_wavelet_")} @@ -4295,7 +4372,7 @@ def paired_figs_for(img_prefix: str, xs_prefix: str) -> str: {figs_for("weber_")}

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 +is large relative to the local median luminance. The log-ratio 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.

@@ -4326,6 +4403,8 @@ def paired_figs_for(img_prefix: str, xs_prefix: str) -> str:

Noise-normalised (× noise floor) — shared colour scale is a fair A/B comparison.

{figs_for("nrm_gradient_")} +{corr_section_html} +

8h. Noise-Corrected Contrast — Cross-Method Overview

{_hires_img_tag(figs.get("nc_ratio_overview"), "NC ratio overview")}

Ratio A/B for every noise-corrected method plotted against its diff --git a/tests/test_analysis/test_spatial_detail.py b/tests/test_analysis/test_spatial_detail.py index 2a59626..ec0f99a 100644 --- a/tests/test_analysis/test_spatial_detail.py +++ b/tests/test_analysis/test_spatial_detail.py @@ -315,3 +315,76 @@ def test_valid_masks_return_ratio(self): score, noise = analyzer._nc_score(detail, mask_neb_shared, bg_mask) assert score == pytest.approx(5.0) assert noise == pytest.approx(2.0) + + +class TestLogRatioHelper: + """Direct unit tests of _log_ratio_map's epsilon-floor and sign-discard contract.""" + + def test_equal_inputs_give_zero(self): + a = np.full((10, 10), 5.0, dtype=np.float32) + b = np.full((10, 10), 5.0, dtype=np.float32) + result = SpatialDetailAnalyzer._log_ratio_map(a, b) + assert np.allclose(result, 0.0, atol=1e-5) + + def test_double_ratio_gives_log10_two(self): + b = np.full((10, 10), 3.0, dtype=np.float32) + a = 2.0 * b + result = SpatialDetailAnalyzer._log_ratio_map(a, b) + assert np.allclose(result, np.log10(2.0), atol=1e-5) + + def test_opposite_sign_equal_magnitude_gives_zero(self): + a = np.full((10, 10), -5.0, dtype=np.float32) + b = np.full((10, 10), 5.0, dtype=np.float32) + result = SpatialDetailAnalyzer._log_ratio_map(a, b) + assert np.allclose(result, 0.0, atol=1e-5) + + def test_near_zero_pixel_amid_nonzero_data_stays_finite(self): + a = np.full((10, 10), 5.0, dtype=np.float32) + b = np.full((10, 10), 5.0, dtype=np.float32) + b[0, 0] = 0.0 + result = SpatialDetailAnalyzer._log_ratio_map(a, b) + assert np.all(np.isfinite(result)) + + def test_mismatched_shapes_crop_to_common_size(self): + a = np.full((10, 10), 5.0, dtype=np.float32) + b = np.full((8, 9), 5.0, dtype=np.float32) + result = SpatialDetailAnalyzer._log_ratio_map(a, b) + assert result.shape == (8, 9) + + def test_all_zero_inputs_stay_finite(self): + a = np.zeros((10, 10), dtype=np.float32) + b = np.zeros((10, 10), dtype=np.float32) + result = SpatialDetailAnalyzer._log_ratio_map(a, b) + assert np.all(np.isfinite(result)) + + +class TestMaskIllustrationFigure: + def test_absent_in_single_image_mode(self, astro_image_a): + result = SpatialDetailAnalyzer().analyze(astro_image_a) + assert "mask_illustration" not in result["figures"] + + def test_present_in_two_image_mode(self, nc_result): + assert "mask_illustration" in nc_result["figures"] + assert isinstance(nc_result["figures"]["mask_illustration"], str) + assert len(nc_result["figures"]["mask_illustration"]) > 0 + + +class TestCorrelationScatterFigures: + _CORR_KEYS = ( + [f"corr_std_{ks}px" for ks in STD_KERNEL_SIZES] + + [f"corr_log_{s}" for s in LOG_SIGMAS] + + [f"corr_gradient_{s}" for s in LOG_SIGMAS] + + ["corr_wavelet_2", "corr_wavelet_3"] + + [f"corr_weber_{ks}px" for ks in WEBER_KERNEL_SIZES] + ) + + @pytest.mark.parametrize("key", _CORR_KEYS) + def test_present_in_two_image_mode(self, nc_result, key): + assert key in nc_result["figures"] + assert isinstance(nc_result["figures"][key], str) + assert len(nc_result["figures"][key]) > 0 + + @pytest.mark.parametrize("key", _CORR_KEYS) + def test_absent_in_single_image_mode(self, astro_image_a, key): + result = SpatialDetailAnalyzer().analyze(astro_image_a) + assert key not in result["figures"] From 1c7ecc17b06e3f5f3d699bff5cc943203de8763b Mon Sep 17 00:00:00 2001 From: Brent <52629076+brentmantooth@users.noreply.github.com> Date: Tue, 14 Jul 2026 22:43:29 -0400 Subject: [PATCH 3/8] Merge Section 8 cross-section figures into a single 2x2 grid per group Each spatial-detail metric group (e.g. "Local sigma, kernel 3px") used to render as a tall 3x1 stack (Image A, Image B, log-ratio diff) plus a separate standalone cross-section image glued on below it when a crosshair was set - visually awkward and treating the cross-section as a bolt-on. - _plot_side_by_side now renders a 2x2 grid (A | B on top, log-ratio diff | cross-section on bottom) via a new xs_data parameter, always at the same geometry whether or not a crosshair is set (blank bottom-right panel when it isn't). - _plot_cross_section is replaced by _draw_cross_section, which draws into an existing Axes instead of creating its own Figure - one image per group instead of two. - All 5 metric families (std, LoG, gradient, wavelet, Weber) sample the cross-section before building the figure, and the noise-normalised variant now gets its own real cross-section (sampled from the noise-normalised arrays) instead of no cross-section at all. - Weber contrast gains crosshair/cross-section support for the first time, matching the other four families. - report_builder.py: removed the now-unused paired_figs_for/xs_figs_for helpers, relocated the cross-section methodology note to apply globally, and reworded captions across 8b-8f for the new embedded-panel layout. Co-Authored-By: Claude Sonnet 5 --- CLAUDE.md | 2 +- analysis/image_filters.py | 176 ++++++++++++++------- report/report_builder.py | 98 ++++++------ tests/test_analysis/test_spatial_detail.py | 44 ++++++ 4 files changed, 210 insertions(+), 110 deletions(-) diff --git a/CLAUDE.md b/CLAUDE.md index 8a047c3..52c64b7 100644 --- a/CLAUDE.md +++ b/CLAUDE.md @@ -230,7 +230,7 @@ When adding a new A-vs-B ratio curve to a report figure (precedent: `_power_rati exception, no obviously-wrong output, just a subtly incorrect number. - **Don't add the ratio via `ax.twinx()`** onto the existing absolute-value plot unless both axes are the same kind of quantity (linear-vs-linear, as in - `analysis/image_filters.py::_plot_cross_section`'s A−B difference line). A linear, + `analysis/image_filters.py::_draw_cross_section`'s A−B difference line). A linear, zero-centered ratio next to a log-scale absolute axis has no principled vertical alignment between the two scales — matplotlib's independent autoscaling invents a relationship that isn't in the data. Build a separate, dedicated figure/panel instead. diff --git a/analysis/image_filters.py b/analysis/image_filters.py index 9377c3a..70f10c1 100644 --- a/analysis/image_filters.py +++ b/analysis/image_filters.py @@ -201,6 +201,7 @@ def _clip01(v): return max(0.0, min(1.0, v)) weber_kernel_sizes, image_a.label, _label_b, display_roi=display_roi, + crosshair=crosshair_roi, mask_neb_shared=mask_neb_shared, mask_bg_a=mask_bg_a, mask_bg_b=mask_bg_b, mask_bg_shared=mask_bg_shared, diff_dist_rng=diff_dist_rng, ) @@ -422,6 +423,13 @@ def _std_analysis(self, norm_a, norm_b, "diff": None, } + xs_raw = None + if crosshair is not None and not single: + pos, pa = self._sample_line(std_a, **crosshair) + _, pb = self._sample_line(std_b, **crosshair) + xs_raw = (pos, pa, pb, label_a, label_b, + f"Cross-section — Local σ, kernel {ks}px") + if not single: fig = self._plot_side_by_side( self._crop_border(std_a, SECTION8_BORDER_CROP_FRACTION), @@ -432,6 +440,7 @@ def _std_analysis(self, norm_a, norm_b, cmap=SECTION8_ANALYSIS_CMAP, nonlinear_norm=True, display_roi=None, + xs_data=xs_raw, ) else: fig = self._plot_single( @@ -443,6 +452,12 @@ def _std_analysis(self, norm_a, norm_b, figures[f"std_{ks}px"] = fig if not single and noise_a and noise_b: + xs_nrm = None + if crosshair is not None: + pos_n, pa_n = self._sample_line(std_a / noise_a, **crosshair) + _, pb_n = self._sample_line(std_b / noise_b, **crosshair) + xs_nrm = (pos_n, pa_n, pb_n, label_a, label_b, + f"Cross-section — Local σ (× noise floor), kernel {ks}px") figures[f"nrm_std_{ks}px"] = self._plot_side_by_side( self._crop_border(std_a / noise_a, SECTION8_BORDER_CROP_FRACTION), self._crop_border(std_b / noise_b, SECTION8_BORDER_CROP_FRACTION), @@ -451,16 +466,9 @@ def _std_analysis(self, norm_a, norm_b, diff_title=f"Log ratio (A/B), noise-normalised, kernel {ks}px", cmap=SECTION8_ANALYSIS_CMAP, display_roi=None, + xs_data=xs_nrm, ) - if crosshair is not None: - pos, pa = self._sample_line(std_a, **crosshair) - if not single: - _, pb = self._sample_line(std_b, **crosshair) - figures[f"xs_std_{ks}px"] = self._plot_cross_section( - pos, pa, pb, label_a, label_b, - f"Cross-section — Local σ, kernel {ks}px") - return figs_to_b64(figures, dpi=150), partial def _compute_std_map(self, norm: np.ndarray, kernel_size: int) -> np.ndarray: @@ -647,6 +655,13 @@ def _log_analysis(self, norm_a, norm_b, sigmas, "diff": None, } + xs_raw = None + if crosshair is not None and not single: + pos, pa = self._sample_line(log_a, **crosshair) + _, pb = self._sample_line(log_b, **crosshair) + xs_raw = (pos, pa, pb, label_a, label_b, + f"Cross-section — |LoG|, σ={sigma}px") + if not single: fig = self._plot_side_by_side( self._crop_border(log_a, SECTION8_BORDER_CROP_FRACTION), @@ -657,6 +672,7 @@ def _log_analysis(self, norm_a, norm_b, sigmas, cmap=SECTION8_ANALYSIS_CMAP, nonlinear_norm=True, display_roi=None, + xs_data=xs_raw, ) else: fig = self._plot_single( @@ -668,6 +684,12 @@ def _log_analysis(self, norm_a, norm_b, sigmas, figures[f"log_sigma{sigma}"] = fig if not single and noise_a and noise_b: + xs_nrm = None + if crosshair is not None: + pos_n, pa_n = self._sample_line(log_a / noise_a, **crosshair) + _, pb_n = self._sample_line(log_b / noise_b, **crosshair) + xs_nrm = (pos_n, pa_n, pb_n, label_a, label_b, + f"Cross-section — |LoG| (× noise floor), σ={sigma}px") figures[f"nrm_log_{sigma}"] = self._plot_side_by_side( self._crop_border(log_a / noise_a, SECTION8_BORDER_CROP_FRACTION), self._crop_border(log_b / noise_b, SECTION8_BORDER_CROP_FRACTION), @@ -676,14 +698,8 @@ def _log_analysis(self, norm_a, norm_b, sigmas, diff_title=f"Log ratio (A/B), noise-normalised, σ={sigma}px", cmap=SECTION8_ANALYSIS_CMAP, display_roi=None, + xs_data=xs_nrm, ) - - if crosshair is not None and not single: - pos, pa = self._sample_line(log_a, **crosshair) - _, pb = self._sample_line(log_b, **crosshair) - figures[f"xs_log_sigma{sigma}"] = self._plot_cross_section( - pos, pa, pb, label_a, label_b, - f"Cross-section — |LoG|, σ={sigma}px") return figs_to_b64(figures, dpi=150), partial # ------------------------------------------------------------------ @@ -742,6 +758,13 @@ def _gradient_analysis(self, norm_a, norm_b, sigmas, "diff": None, } + xs_raw = None + if crosshair is not None and not single: + pos, pa = self._sample_line(gm_a, **crosshair) + _, pb = self._sample_line(gm_b, **crosshair) + xs_raw = (pos, pa, pb, label_a, label_b, + f"Cross-section — Gradient, σ={sigma}px") + if not single: fig = self._plot_side_by_side( self._crop_border(gm_a, SECTION8_BORDER_CROP_FRACTION), @@ -752,6 +775,7 @@ def _gradient_analysis(self, norm_a, norm_b, sigmas, cmap=SECTION8_ANALYSIS_CMAP, nonlinear_norm=True, display_roi=None, + xs_data=xs_raw, ) else: fig = self._plot_single( @@ -763,6 +787,12 @@ def _gradient_analysis(self, norm_a, norm_b, sigmas, figures[f"gradient_{sigma}"] = fig if not single and noise_a and noise_b: + xs_nrm = None + if crosshair is not None: + pos_n, pa_n = self._sample_line(gm_a / noise_a, **crosshair) + _, pb_n = self._sample_line(gm_b / noise_b, **crosshair) + xs_nrm = (pos_n, pa_n, pb_n, label_a, label_b, + f"Cross-section — Gradient (× noise floor), σ={sigma}px") figures[f"nrm_gradient_{sigma}"] = self._plot_side_by_side( self._crop_border(gm_a / noise_a, SECTION8_BORDER_CROP_FRACTION), self._crop_border(gm_b / noise_b, SECTION8_BORDER_CROP_FRACTION), @@ -771,14 +801,8 @@ def _gradient_analysis(self, norm_a, norm_b, sigmas, diff_title=f"Log ratio (A/B), noise-normalised, σ={sigma}px", cmap=SECTION8_ANALYSIS_CMAP, display_roi=None, + xs_data=xs_nrm, ) - - if crosshair is not None and not single: - pos, pa = self._sample_line(gm_a, **crosshair) - _, pb = self._sample_line(gm_b, **crosshair) - figures[f"xs_gradient_{sigma}"] = self._plot_cross_section( - pos, pa, pb, label_a, label_b, - f"Cross-section — Gradient, σ={sigma}px") return figs_to_b64(figures, dpi=150), partial # ------------------------------------------------------------------ @@ -868,6 +892,13 @@ def _wavelet_analysis(self, norm_a, norm_b, wavelet, levels, "b": (rec_b / noise_b).astype(np.float32), "diff": None, } + xs_raw = None + if crosshair is not None and not single: + pos, pa = self._sample_line(rec_a, **crosshair) + _, pb = self._sample_line(rec_b, **crosshair) + xs_raw = (pos, pa, pb, label_a, label_b, + f"Cross-section — Wavelet level {display_level}") + if not single: fig = self._plot_side_by_side( self._crop_border(rec_a, SECTION8_BORDER_CROP_FRACTION), @@ -877,6 +908,7 @@ def _wavelet_analysis(self, norm_a, norm_b, wavelet, levels, diff_title=f"Level {display_level} log-ratio (A/B)", cmap=SECTION8_ANALYSIS_CMAP, display_roi=None, + xs_data=xs_raw, ) else: fig = self._plot_single( @@ -887,6 +919,12 @@ def _wavelet_analysis(self, norm_a, norm_b, wavelet, levels, figures[f"wavelet_level{display_level}"] = fig if not single and noise_a and noise_b: + xs_nrm = None + if crosshair is not None: + pos_n, pa_n = self._sample_line(rec_a / noise_a, **crosshair) + _, pb_n = self._sample_line(rec_b / noise_b, **crosshair) + xs_nrm = (pos_n, pa_n, pb_n, label_a, label_b, + f"Cross-section — Wavelet level {display_level} (× noise floor)") figures[f"nrm_wavelet_{display_level}"] = self._plot_side_by_side( self._crop_border(rec_a / noise_a, SECTION8_BORDER_CROP_FRACTION), self._crop_border(rec_b / noise_b, SECTION8_BORDER_CROP_FRACTION), @@ -895,15 +933,9 @@ def _wavelet_analysis(self, norm_a, norm_b, wavelet, levels, diff_title=f"Level {display_level} log-ratio (A/B), noise-normalised", cmap=SECTION8_ANALYSIS_CMAP, display_roi=None, + xs_data=xs_nrm, ) - if crosshair is not None and not single: - pos, pa = self._sample_line(rec_a, **crosshair) - _, pb = self._sample_line(rec_b, **crosshair) - figures[f"xs_wavelet_level{display_level}"] = self._plot_cross_section( - pos, pa, pb, label_a, label_b, - f"Cross-section — Wavelet level {display_level}") - return figs_to_b64(figures, dpi=150), partial def _estimate_noise(self, coeffs) -> float: @@ -945,6 +977,7 @@ def _reconstruct_level(self, coeffs, target_coeff_idx: int, def _weber_analysis(self, norm_a, norm_b, kernel_sizes, label_a, label_b, display_roi=None, + crosshair=None, mask_neb_shared=None, mask_bg_a=None, mask_bg_b=None, mask_bg_shared=None, diff_dist_rng=None) -> tuple[dict, dict]: figures = {} @@ -997,6 +1030,13 @@ def _weber_analysis(self, norm_a, norm_b, kernel_sizes, "diff": None, } + xs_raw = None + if crosshair is not None and not single: + pos, pa = self._sample_line(wc_a, **crosshair) + _, pb = self._sample_line(wc_b, **crosshair) + xs_raw = (pos, pa, pb, label_a, label_b, + f"Cross-section — Weber contrast, kernel {ks}px") + if not single: fig = self._plot_side_by_side( self._crop_border(wc_a, SECTION8_BORDER_CROP_FRACTION), @@ -1007,6 +1047,7 @@ def _weber_analysis(self, norm_a, norm_b, kernel_sizes, cmap=SECTION8_ANALYSIS_CMAP, nonlinear_norm=True, # unbounded output; PowerNorm compresses dynamic range display_roi=None, # _crop_border already applied + xs_data=xs_raw, ) else: fig = self._plot_single( @@ -1018,6 +1059,12 @@ def _weber_analysis(self, norm_a, norm_b, kernel_sizes, figures[f"weber_{ks}px"] = fig if not single and noise_a and noise_b: + xs_nrm = None + if crosshair is not None: + pos_n, pa_n = self._sample_line(wc_a / noise_a, **crosshair) + _, pb_n = self._sample_line(wc_b / noise_b, **crosshair) + xs_nrm = (pos_n, pa_n, pb_n, label_a, label_b, + f"Cross-section — Weber contrast (× noise floor), kernel {ks}px") figures[f"nrm_weber_{ks}px"] = self._plot_side_by_side( self._crop_border(wc_a / noise_a, SECTION8_BORDER_CROP_FRACTION), self._crop_border(wc_b / noise_b, SECTION8_BORDER_CROP_FRACTION), @@ -1026,6 +1073,7 @@ def _weber_analysis(self, norm_a, norm_b, kernel_sizes, diff_title=f"Log ratio (A/B), noise-normalised, kernel {ks}px", cmap=SECTION8_ANALYSIS_CMAP, display_roi=None, + xs_data=xs_nrm, ) return figs_to_b64(figures, dpi=150), partial @@ -1066,7 +1114,10 @@ def _plot_side_by_side(self, arr_a: np.ndarray, arr_b: np.ndarray, cmap: str = "viridis", nonlinear_norm: bool = False, display_roi=None, - smooth_display: bool = True) -> plt.Figure: + smooth_display: bool = True, + xs_data: tuple | None = None) -> plt.Figure: + """xs_data, if given, is (pos, prof_a, prof_b, label_a, label_b, xs_title) + for the embedded cross-section panel; None leaves that panel blank.""" # Crop to bright-feature ROI if available if display_roi is not None: r0, r1, c0, c1 = display_roi @@ -1095,16 +1146,18 @@ def _plot_side_by_side(self, arr_a: np.ndarray, arr_b: np.ndarray, d_max = float(np.percentile(np.abs(diff), 99.5)) or 1.0 dvmin, dvmax = -d_max, d_max - # 3×1 column layout — A, B, then diff stacked vertically. - # Use 1:1 pixel aspect; size figure based on the cropped array dimensions. + # 2×2 grid: A|B on top, log-ratio diff | cross-section on bottom. A/B/diff + # share the source array's pixel aspect ratio (aspect="equal" imshow); the + # cross-section panel is a line plot with no such constraint but occupies + # an equal-size grid cell so the other three panels stay geometrically + # identical whether or not a crosshair (and thus xs_data) is set. h, w = arr_a.shape[:2] aspect_ratio = h / max(w, 1) - panel_w = 10.0 + panel_w = 5.0 # half the old single-column width — 2 columns now share it panel_h = panel_w * aspect_ratio - fig_h = panel_h * 3 + 2.0 # 3 panels + headroom for colorbars/titles - fig, axes = plt.subplots(3, 1, figsize=(panel_w, fig_h), + fig, axes = plt.subplots(2, 2, figsize=(panel_w * 2, panel_h * 2 + 1.5), constrained_layout=True) - ax_a, ax_b, ax_diff = axes + (ax_a, ax_b), (ax_diff, ax_xs) = axes for ax, arr, title in zip([ax_a, ax_b], [arr_a, arr_b], [title_a, title_b]): im = ax.imshow(arr, origin="upper", cmap=cmap, @@ -1123,6 +1176,13 @@ def _plot_side_by_side(self, arr_a: np.ndarray, arr_b: np.ndarray, ax_diff.axis("off") fig.colorbar(im_diff, ax=ax_diff, fraction=0.046, pad=0.04) + if xs_data is not None: + pos, prof_a, prof_b, xs_label_a, xs_label_b, xs_title = xs_data + self._draw_cross_section(ax_xs, pos, prof_a, prof_b, + xs_label_a, xs_label_b, xs_title) + else: + ax_xs.axis("off") + return fig def _plot_single(self, arr_a: np.ndarray, title_a: str, @@ -1278,30 +1338,34 @@ def _sample_line(arr: np.ndarray, x0: float, y0: float, return positions, values @staticmethod - def _plot_cross_section(pos: np.ndarray, prof_a: np.ndarray, prof_b: np.ndarray, - label_a: str, label_b: str, title: str) -> plt.Figure: + def _draw_cross_section(ax, pos: np.ndarray, prof_a: np.ndarray, prof_b: np.ndarray, + label_a: str, label_b: str, title: str) -> None: + """Draw a cross-section profile (+ linear A−B diff on a twinx — CLAUDE.md's + sanctioned linear-vs-linear twinx precedent, unrelated to the log-ratio + metric used elsewhere in Section 8) into an existing Axes: the bottom-right + quadrant of _plot_side_by_side's 2×2 grid. Fonts/linewidths are tuned down + from the metric's old standalone-figure sizes since this panel is now + roughly a quarter the area.""" # Images with slightly different pixel dimensions produce different-length profiles n = min(len(pos), len(prof_a), len(prof_b)) pos, prof_a, prof_b = pos[:n], prof_a[:n], prof_b[:n] - fig, ax1 = plt.subplots(figsize=(9, 4), constrained_layout=True) - ax1.plot(pos, prof_a, color="steelblue", linewidth=1.5, alpha=XS_LINE_ALPHA, label=label_a) - ax1.plot(pos, prof_b, color="tomato", linewidth=1.5, alpha=XS_LINE_ALPHA, label=label_b) - ax1.set_xlabel("Position along line (px)") - ax1.set_ylabel("Map value") - ax1.legend(loc="upper left", fontsize=8) - ax1.grid(True, alpha=0.3) - - ax2 = ax1.twinx() - n = min(len(prof_a), len(prof_b)) - diff = prof_a[:n] - prof_b[:n] - ax2.plot(pos[:n], diff, color="#2ca02c", linewidth=1.2, + ax.plot(pos, prof_a, color="steelblue", linewidth=1.1, alpha=XS_LINE_ALPHA, label=label_a) + ax.plot(pos, prof_b, color="tomato", linewidth=1.1, alpha=XS_LINE_ALPHA, label=label_b) + ax.set_xlabel("Position along line (px)", fontsize=8) + ax.set_ylabel("Map value", fontsize=8) + ax.tick_params(labelsize=7) + ax.legend(loc="upper left", fontsize=6.5, labelspacing=0.3) + ax.grid(True, alpha=0.3) + + ax2 = ax.twinx() + diff = prof_a - prof_b + ax2.plot(pos, diff, color="#2ca02c", linewidth=0.9, linestyle="--", alpha=0.85, label="A−B") - ax2.axhline(0, color="#2ca02c", linewidth=0.8, alpha=0.3) # zero-crossing reference - ax2.set_ylabel("Difference (A−B)", color="#2ca02c") - ax2.tick_params(axis="y", labelcolor="#2ca02c") - ax2.legend(loc="upper right", fontsize=8) - ax1.set_title(title, fontsize=10) - return fig + ax2.axhline(0, color="#2ca02c", linewidth=0.6, alpha=0.3) # zero-crossing reference + ax2.set_ylabel("Difference (A−B)", color="#2ca02c", fontsize=8) + ax2.tick_params(axis="y", labelcolor="#2ca02c", labelsize=7) + ax2.legend(loc="upper right", fontsize=6.5, labelspacing=0.3) + ax.set_title(title, fontsize=9) @staticmethod def _plot_metric_correlation(map_a: np.ndarray, map_b: np.ndarray, diff --git a/report/report_builder.py b/report/report_builder.py index f520462..e358457 100644 --- a/report/report_builder.py +++ b/report/report_builder.py @@ -4124,24 +4124,6 @@ def figs_for(prefix): out += _hires_img_tag(figs[key], key) + "\n" return out - def xs_figs_for(prefix: str) -> str: - out = "" - for key in sorted(figs): - if key.startswith(prefix): - out += _hires_img_tag(figs[key], key) + "\n" - return out - - def paired_figs_for(img_prefix: str, xs_prefix: str) -> str: - """Emit each image immediately followed by its matching cross-section.""" - out = "" - for img_key in sorted(k for k in figs if k.startswith(img_prefix)): - suffix = img_key[len(img_prefix):] - out += _hires_img_tag(figs[img_key], img_key) + "\n" - xs_key = xs_prefix + suffix - if xs_key in figs: - out += _hires_img_tag(figs[xs_key], xs_key) + "\n" - return out - def corr_figs_html() -> str: """Emit each present correlation scatter figure in _SPATIAL_CORR_ROWS order (not alphabetic — e.g. std_10px must not sort before std_3px).""" @@ -4185,9 +4167,11 @@ def corr_figs_html() -> str: has_crosshair = sm.get("crosshair") is not None xs_note = _info_box( - 'ℹ Cross-section profiles below are extracted along ' - 'the line selected in the viewer. Left axis: both images ' - '(steelblue = A, tomato = B). Right axis (green dashed): difference A−B.', + 'ℹ When a cross-section line is set in the viewer, its profile is embedded ' + 'as the bottom-right panel of each map-pair figure below (8b–8f): Image A ' + '(top-left), Image B (top-right), log-ratio map (bottom-left), cross-section ' + 'profile (bottom-right). Left axis: both images (steelblue = A, tomato = B). ' + 'Right axis (green dashed): difference A−B.', title="Cross-section profiles", open=True, ) if has_crosshair else "" @@ -4241,7 +4225,8 @@ def corr_figs_html() -> str: f'thin filaments), Level 3 ≈ 8 px (medium structures — emission knots, shell edges), ' f'Level 4 ≈ 16 px (broader features). A higher SNR at Level 2 indicates the filter ' f'preserves sub-arcsecond detail better; Level 3 reflects medium-scale structure. ' - f'Cross-section profiles show how detail amplitude varies spatially along the selected line.', + f'When a cross-section line is set, its profile (embedded in each map figure below) ' + f'shows how detail amplitude varies spatially along the selected line.', title="Wavelet decomposition", ) @@ -4254,17 +4239,20 @@ def corr_figs_html() -> str: {smooth_note} {_info_box('All maps below are computed on mean-signal-normalised data ' '(each image divided by its own mean signal), making them dimensionless and comparable ' - 'across different filter bandwidths. Images are shown side-by-side with a shared ' - 'colour scale; the third panel shows the per-pixel log-ratio ' - 'log10(A / B): red (positive) means A is stronger at that pixel, blue ' - '(negative) means B is stronger, white means the two are equal. ±0.3 ≈ a 2× ' - 'difference, ±1.0 ≈ a 10× difference. For the Original and Wavelet maps — the ' - 'only two families whose values can be negative — the ratio uses |A| and |B|, so it ' - 'compares the magnitude of structure rather than signed brightness.', + 'across different filter bandwidths. Each figure is a 2×2 grid: Image A ' + '(top-left) and Image B (top-right) share a colour scale; the bottom-left panel shows ' + 'the per-pixel log-ratio log10(A / B): red (positive) means A is stronger at ' + 'that pixel, blue (negative) means B is stronger, white means the two are equal. ±0.3 ≈ ' + 'a 2× difference, ±1.0 ≈ a 10× difference. For the Original and Wavelet maps ' + '— the only two families whose values can be negative — the ratio uses |A| and |B|, so ' + 'it compares the magnitude of structure rather than signed brightness. The ' + 'bottom-right panel shows the cross-section profile along the selected line when one is ' + 'set, otherwise left blank.', title="Spatial detail maps overview")} {dist_html} {nc_methodology_box} {nc_empty_note} +{xs_note}

8b. Local Standard Deviation Maps

{_info_box('Measures how much pixel values vary within a neighbourhood. ' @@ -4275,8 +4263,8 @@ def corr_figs_html() -> str: 'window. Brighter regions contain more local variation — typically nebula filaments, ' 'star halos, or noise. A filter with higher std values in targeted emission regions ' 'preserves more structure; higher std in blank sky regions indicates more photon noise. ' - 'The cross-section profiles below each map pair show how local detail amplitude varies ' - 'along the selected line.', + 'When a cross-section line is set, its profile is embedded in the bottom-right panel of ' + 'each map figure below, showing how local detail amplitude varies along the selected line.', title="Local standard deviation")} @@ -4286,10 +4274,10 @@ def corr_figs_html() -> str: {std_nc_rows}
Kernel size{ra.label}{rb.label}
Scale{ra.label} (NC score){rb.label} (NC score)Ratio A/B
-{xs_note}{paired_figs_for("std_", "xs_std_")} -

Side-by-side local σ maps at each kernel size (shared colour scale), -each followed by its cross-section profile. -The log-ratio map (right) highlights where one filter preserves more local variation.

+{figs_for("std_")} +

Local σ maps at each kernel size (shared colour scale): Image A (top-left), +Image B (top-right), log-ratio map (bottom-left) highlighting where one filter preserves more +local variation, and — when a cross-section line is set — its profile (bottom-right).

Noise-normalised (× noise floor) — shared colour scale is a fair A/B comparison.

{figs_for("nrm_std_")}

8c. Laplacian of Gaussian (LoG) Maps

@@ -4301,17 +4289,18 @@ def corr_figs_html() -> str: 'then computing the Laplacian (second spatial derivative), which peaks at intensity ' 'boundaries. |LoG| is shown so bright-to-dark and dark-to-bright edges are treated ' 'equally. Compare maps at each σ: a sharper or higher-contrast filter will show ' - 'brighter LoG response at small σ values. Cross-section profiles reveal subtle ' - 'differences in edge sharpness along the selected line.', + 'brighter LoG response at small σ values. When a cross-section line is set, its profile ' + 'reveals subtle differences in edge sharpness along the selected line.', title="Laplacian of Gaussian (LoG)")} {log_nc_rows}
Scale{ra.label} (NC score){rb.label} (NC score)Ratio A/B
-{paired_figs_for("log_", "xs_log_")} -

|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 σ.

+{figs_for("log_")} +

|LoG| maps at σ = 1.5, 3, and 6 px (shared colour scale per figure): +Image A (top-left), Image B (top-right), log-ratio map (bottom-left), and — when a +cross-section line is set — its profile (bottom-right). 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_")}

8d. Wavelet Decomposition

@@ -4330,11 +4319,12 @@ def corr_figs_html() -> str: {wavelet_nc_rows} -{paired_figs_for("wavelet_level", "xs_wavelet_level")} +{figs_for("wavelet_level")}

Reconstructed detail images at levels 2 and 3 (shared colour scale, -diverging colourmap), each followed by its cross-section profile. -The log-ratio panel (right) shows where fine structure differs between the two filters -(sign discarded — |A|/|B| — since wavelet reconstructions can be negative).

+diverging colourmap): Image A (top-left), Image B (top-right), log-ratio panel (bottom-left) +showing where fine structure differs between the two filters (sign discarded — |A|/|B| — +since wavelet reconstructions can be negative), and — when a cross-section line is set — +its profile (bottom-right).

Noise-normalised (× noise floor) — shared colour scale is a fair A/B comparison.

{figs_for("nrm_wavelet_")} @@ -4371,10 +4361,11 @@ def corr_figs_html() -> str: {figs_for("weber_")}

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 log-ratio 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.

+scale, viridis): Image A (top-left), Image B (top-right), log-ratio panel (bottom-left) +showing where one image achieves greater relative contrast, and — when a cross-section +line is set — its profile (bottom-right). Brighter regions have higher Weber contrast — +the local intensity range is large relative to the local median luminance. 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_")} @@ -4396,10 +4387,11 @@ def corr_figs_html() -> str: Scale{ra.label} (NC score){rb.label} (NC score)Ratio A/B {gm_nc_rows} -{paired_figs_for("gradient_", "xs_gradient_")} -

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.

+{figs_for("gradient_")} +

Gradient magnitude maps at σ = 1.5, 3, and 6 px (shared colour scale per +figure): Image A (top-left), Image B (top-right), log-ratio map (bottom-left), and — when a +cross-section line is set — its profile (bottom-right). 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_")} diff --git a/tests/test_analysis/test_spatial_detail.py b/tests/test_analysis/test_spatial_detail.py index ec0f99a..b0fde89 100644 --- a/tests/test_analysis/test_spatial_detail.py +++ b/tests/test_analysis/test_spatial_detail.py @@ -388,3 +388,47 @@ def test_present_in_two_image_mode(self, nc_result, key): def test_absent_in_single_image_mode(self, astro_image_a, key): result = SpatialDetailAnalyzer().analyze(astro_image_a) assert key not in result["figures"] + + +@pytest.fixture(scope="module") +def nc_result_with_crosshair(nc_image_pair) -> dict: + img_a, img_b = nc_image_pair + crosshair = {"x0": 0.1, "y0": 0.5, "x1": 0.9, "y1": 0.5} + return SpatialDetailAnalyzer().analyze(img_a, img_b, crosshair=crosshair) + + +class TestCrosshairEmbeddedCrossSections: + """Section 8 cross-sections are embedded panels inside the 2x2 map figures, + not separate figures — no standalone xs_* keys should ever appear, and all + 5 families (including Weber, newly crosshair-capable) must still produce + their usual figure keys with or without a crosshair, without crashing.""" + + def test_no_crash_with_crosshair(self, nc_result_with_crosshair): + assert isinstance(nc_result_with_crosshair, dict) + + def test_no_standalone_xs_keys_with_crosshair(self, nc_result_with_crosshair): + figs = nc_result_with_crosshair["figures"] + assert not any(k.startswith(("xs_std_", "xs_log_", "xs_wavelet_", + "xs_gradient_", "xs_weber_")) for k in figs) + + def test_no_standalone_xs_keys_without_crosshair(self, nc_result): + figs = nc_result["figures"] + assert not any(k.startswith(("xs_std_", "xs_log_", "xs_wavelet_", + "xs_gradient_", "xs_weber_")) for k in figs) + + @pytest.mark.parametrize("key_prefix,scales", [ + ("std_", STD_KERNEL_SIZES), ("weber_", WEBER_KERNEL_SIZES), + ]) + def test_family_figures_present_with_crosshair(self, nc_result_with_crosshair, + key_prefix, scales): + figs = nc_result_with_crosshair["figures"] + for scale in scales: + assert f"{key_prefix}{scale}px" in figs + + def test_weber_no_crash_without_crosshair(self, nc_result): + # Regression: weber previously had no crosshair param at all. + assert all(f"weber_{ks}px" in nc_result["figures"] for ks in WEBER_KERNEL_SIZES) + + def test_weber_nrm_figures_present_with_crosshair(self, nc_result_with_crosshair): + figs = nc_result_with_crosshair["figures"] + assert any(k.startswith("nrm_weber_") for k in figs) From 298217b3e2bbff55b939322d9413a322af584504 Mon Sep 17 00:00:00 2001 From: Brent <52629076+brentmantooth@users.noreply.github.com> Date: Tue, 14 Jul 2026 22:47:57 -0400 Subject: [PATCH 4/8] Add timeout-minutes to CI/release jobs to fail fast on runner hangs A release build silently hung on the ubuntu-latest apt-get step for the full default 6-hour job timeout before being canceled. Bound job and apt-get step durations so a future hang fails in minutes, not hours. Co-Authored-By: Claude Sonnet 5 --- .github/workflows/ci.yml | 3 +++ .github/workflows/release.yml | 2 ++ 2 files changed, 5 insertions(+) diff --git a/.github/workflows/ci.yml b/.github/workflows/ci.yml index 34b5ceb..977da0d 100644 --- a/.github/workflows/ci.yml +++ b/.github/workflows/ci.yml @@ -13,6 +13,7 @@ jobs: os: [windows-latest, ubuntu-latest, macos-latest] runs-on: ${{ matrix.os }} + timeout-minutes: 20 steps: - uses: actions/checkout@v4 @@ -44,6 +45,7 @@ jobs: os: [windows-latest, ubuntu-latest, macos-latest] runs-on: ${{ matrix.os }} + timeout-minutes: 45 steps: - uses: actions/checkout@v4 @@ -57,6 +59,7 @@ jobs: # to import PyQt6 during the analysis phase. - name: Install Linux Qt system libraries if: runner.os == 'Linux' + timeout-minutes: 5 run: | sudo apt-get update -qq sudo apt-get install -y libgl1 libegl1 libxcb-cursor0 libxkbcommon-x11-0 diff --git a/.github/workflows/release.yml b/.github/workflows/release.yml index be16023..589c54a 100644 --- a/.github/workflows/release.yml +++ b/.github/workflows/release.yml @@ -16,6 +16,7 @@ jobs: os: [windows-latest, ubuntu-latest, macos-latest] runs-on: ${{ matrix.os }} + timeout-minutes: 45 steps: - uses: actions/checkout@v4 @@ -27,6 +28,7 @@ jobs: - name: Install Linux Qt system libraries if: runner.os == 'Linux' + timeout-minutes: 5 run: | sudo apt-get update -qq sudo apt-get install -y libgl1 libegl1 libxcb-cursor0 libxkbcommon-x11-0 From 72ab11373557b5a495a4002c4d82aa95d46be2be Mon Sep 17 00:00:00 2001 From: Brent <52629076+brentmantooth@users.noreply.github.com> Date: Wed, 15 Jul 2026 19:24:08 -0400 Subject: [PATCH 5/8] Simplify Section 8 cross-sections, interleave correlation plots, dedupe background calc Cross-section panels no longer plot a distracting A-B diff trace on a twinx axis; per-scale correlation scatters now sit directly next to their matching map figure (8b-8f) instead of in one late "8g" dump, fixing a latent alphabetic mis-ordering (std_10px before std_3px) along the way. Also confirmed the nebula/background mask already uses full-image background stats, not the ROI, and closed a real inefficiency found while checking: estimate_background() is now idempotent and precomputed once per image object before parallel analyzer dispatch, instead of being recomputed by every one of the 6 analyzers that touch the same AstroImage instance. Co-Authored-By: Claude Sonnet 5 --- CLAUDE.md | 37 +++++++++-- analysis/image_filters.py | 13 +--- core/astro_image.py | 2 + gui/analysis_thread.py | 14 ++++ report/report_builder.py | 131 +++++++++++++++++++++----------------- 5 files changed, 122 insertions(+), 75 deletions(-) diff --git a/CLAUDE.md b/CLAUDE.md index 52c64b7..6ee2888 100644 --- a/CLAUDE.md +++ b/CLAUDE.md @@ -62,6 +62,7 @@ synthetic/ | `_log_ratio_map(a, b)` | `analysis/image_filters.py` | Per-pixel `log10(\|A\|/\|B\|)` map with percentile-based epsilon floor and defensive shape crop — the Section 8 replacement for plain `A − B` diff | | `_plot_mask_illustration(base, mask_neb, mask_bg)` | `analysis/image_filters.py` | Translucent steelblue/tomato mask overlay on a grayscale base image | | `_plot_metric_correlation(map_a, map_b, mask_neb, mask_bg, ...)` | `analysis/image_filters.py` | 1×2 masked-region scatter (A vs B) with a 1:1 line plus overall + tail-restricted dual linear fits | +| `_family_figs_with_corr(rows, map_key_fn)` | `report_builder.py` | Emits a Section 8 family's map figure immediately followed by its `corr_*` correlation scatter, one scale at a time, in numeric order (`_SPATIAL_CORR_ROWS`) — the pattern to follow when adding any new per-scale Section 8 figure pair | --- @@ -229,11 +230,13 @@ When adding a new A-vs-B ratio curve to a report figure (precedent: `_power_rati `20 * np.log10(ratio)`. Using the wrong constant is silently off by 2× in dB — no exception, no obviously-wrong output, just a subtly incorrect number. - **Don't add the ratio via `ax.twinx()`** onto the existing absolute-value plot unless - both axes are the same kind of quantity (linear-vs-linear, as in - `analysis/image_filters.py::_draw_cross_section`'s A−B difference line). A linear, - zero-centered ratio next to a log-scale absolute axis has no principled vertical - alignment between the two scales — matplotlib's independent autoscaling invents a - relationship that isn't in the data. Build a separate, dedicated figure/panel instead. + both axes are the same kind of quantity (linear-vs-linear). A linear, zero-centered + ratio next to a log-scale absolute axis has no principled vertical alignment between + the two scales — matplotlib's independent autoscaling invents a relationship that + isn't in the data. Build a separate, dedicated figure/panel instead. (The codebase's + prior linear-vs-linear precedent, `_draw_cross_section`'s A−B difference line, was + removed as unnecessary clutter — there is currently no `ax.twinx()` usage anywhere + in the codebase, so treat this as a rule to apply fresh, not an existing pattern to copy.) - **Guard bin alignment before dividing two arrays from different analyses.** Two per-image radial/frequency arrays are only safely divisible bin-for-bin when they share the same shape *and* values (`freq_a.shape == freq_b.shape and @@ -253,6 +256,27 @@ When adding a new A-vs-B ratio curve to a report figure (precedent: `_power_rati (`_plot_side_by_side`'s `np.percentile(arr, 0.5)`), applied to the epsilon floor instead of just the color scale. +### Background estimation — compute once via the pre-pass, never redundantly + +`AstroImage.estimate_background()` (`core/astro_image.py`) is idempotent: it returns +immediately if `self.background is not None`, since `self.data` is only ever set once, +during `load()`. Every analyzer (`SNRAnalyzer`, `PSFAnalyzer`, `HaloAnalyzer`, +`EdgeAnalyzer`, `PowerSpectrumAnalyzer`, `SpatialDetailAnalyzer`) still calls +`estimate_background()` unconditionally at the top of its `analyze()` — that's +intentional and does not need to change; the idempotency guard just makes each of +those calls a cheap no-op once the object's background has already been computed. + +`gui/analysis_thread.py::_execute()` runs a pre-pass — after alignment, before task +dispatch — that calls `estimate_background()` once per distinct `AstroImage` object +(`img_a`, `img_b`, `self._starless_a`, `self._starless_b`) via a small +`ThreadPoolExecutor`. This exists because multiple analyzers share the same image +object and can run concurrently under `parallel=True`; without the pre-pass each one +would independently trigger a full `Background2D` computation (expensive) and race to +write `self.background` / `self.background_rms` on the same object. When adding a new +analyzer that needs background stats, just call `image.estimate_background()` as +normal at the top of `analyze()` — do not add another pre-pass call site; the existing +one in `_execute()` already covers every image object the thread constructs. + --- ## Collaboration Rules @@ -351,7 +375,8 @@ pytest tests/ --cov=analysis,core,synthetic,report --cov-report=html | Mixed float32/float64 arithmetic silently widens to float64 | NumPy upcasts when operands differ (e.g. `float32_array - float64_scalar`). If photutils ever returns a float64 background model, `background_subtracted()` will silently return float64. Guard by adding `.astype(np.float32)` at the end of `background_subtracted()` in `astro_image.py` if this is observed. | | `_section_snr` crashed when SNR metric is unchecked | `_plot_snr_pair` (`report_builder.py`) built a `panels` list filtered to non-`None` entries but never checked whether it was empty before calling `plt.subplots(1, len(panels), ...)` — 0 columns raised `ValueError: Number of columns must be a positive integer, not 0`. Hit whenever SNR is unchecked while another metric (e.g. Power Spectrum) is run. Fixed with an early `if not panels: return None` guard, matching `_plot_radial_overlay`/`_plot_radial_ratio_db`; both call sites already pipe the result through `_img_tag`, which turns `None` into `""`. | | Whole-population regression slope hides tail-specific divergence | A single OLS fit over an entire A-vs-B scatter is dominated by leverage (distance from the mean), not point count — it can show near-zero correlation even when a sparse high-value tail diverges sharply (or vice versa), because it blends "the bulk agrees" and "the tail diverges" into one ambiguous number. When the signal of interest lives specifically in the tail (Section 8g: "which filter shows more detail in the brightest/most-structured pixels"), fit two lines — an overall fit and a second fit restricted to the top N% by combined magnitude (`_plot_metric_correlation`, `SECTION8_SCATTER_TAIL_PERCENTILE`) — and report slope + R² for each. High tail R² with slope far from 1 is a real, systematic effect; low tail R² means the apparent divergence is mostly noise. | -| New Section 8 panel key doesn't need Report Inspector code changes | `gui/report_inspector.py` is fully generic — driven entirely by a companion `_inspector.npz` (raw float32/uint8 arrays) plus an embedded `catalog_json` built in `report_builder.py::_write_inspector_file`. `_panel_display_name`/`_panel_concept` dynamically parse any `panels` dict key prefix, so a new `SpatialDetailAnalyzer` panel family auto-appears in the inspector with zero inspector-side changes. A genuinely new *visual type* is a different story: the inspector only knows how to `imshow` 2D/RGB arrays (side-by-side or slider-reveal), so scatter-style plots (Section 8g's `corr_*` figures) must stay static-HTML-only unless new inspector canvas code is written. | +| New Section 8 panel key doesn't need Report Inspector code changes | `gui/report_inspector.py` is fully generic — driven entirely by a companion `_inspector.npz` (raw float32/uint8 arrays) plus an embedded `catalog_json` built in `report_builder.py::_write_inspector_file`. `_panel_display_name`/`_panel_concept` dynamically parse any `panels` dict key prefix, so a new `SpatialDetailAnalyzer` panel family auto-appears in the inspector with zero inspector-side changes. A genuinely new *visual type* is a different story: the inspector only knows how to `imshow` 2D/RGB arrays (side-by-side or slider-reveal), so scatter-style plots (Section 8's `corr_*` correlation figures, interleaved into 8b–8f right after each map figure via `_family_figs_with_corr`) must stay static-HTML-only unless new inspector canvas code is written. | +| Renumbering a Section 8 subsection misses caption cross-references | Section 8's sub-heading letters (8a–8g) are referenced by literal string in caption/info-box text scattered throughout `_section_spatial` — not just in the `

` tags (e.g. "see 8g for…", "(8b–8f, 8g)"). After adding, removing, or renumbering a subsection, `grep` the function for every old *and* new heading letter — HTML renders a stale cross-reference without error, it just silently misdirects the reader to the wrong subsection. | --- diff --git a/analysis/image_filters.py b/analysis/image_filters.py index 70f10c1..ce81827 100644 --- a/analysis/image_filters.py +++ b/analysis/image_filters.py @@ -1340,9 +1340,7 @@ def _sample_line(arr: np.ndarray, x0: float, y0: float, @staticmethod def _draw_cross_section(ax, pos: np.ndarray, prof_a: np.ndarray, prof_b: np.ndarray, label_a: str, label_b: str, title: str) -> None: - """Draw a cross-section profile (+ linear A−B diff on a twinx — CLAUDE.md's - sanctioned linear-vs-linear twinx precedent, unrelated to the log-ratio - metric used elsewhere in Section 8) into an existing Axes: the bottom-right + """Draw a cross-section profile into an existing Axes: the bottom-right quadrant of _plot_side_by_side's 2×2 grid. Fonts/linewidths are tuned down from the metric's old standalone-figure sizes since this panel is now roughly a quarter the area.""" @@ -1356,15 +1354,6 @@ def _draw_cross_section(ax, pos: np.ndarray, prof_a: np.ndarray, prof_b: np.ndar ax.tick_params(labelsize=7) ax.legend(loc="upper left", fontsize=6.5, labelspacing=0.3) ax.grid(True, alpha=0.3) - - ax2 = ax.twinx() - diff = prof_a - prof_b - ax2.plot(pos, diff, color="#2ca02c", linewidth=0.9, - linestyle="--", alpha=0.85, label="A−B") - ax2.axhline(0, color="#2ca02c", linewidth=0.6, alpha=0.3) # zero-crossing reference - ax2.set_ylabel("Difference (A−B)", color="#2ca02c", fontsize=8) - ax2.tick_params(axis="y", labelcolor="#2ca02c", labelsize=7) - ax2.legend(loc="upper right", fontsize=6.5, labelspacing=0.3) ax.set_title(title, fontsize=9) @staticmethod diff --git a/core/astro_image.py b/core/astro_image.py index add122c..e293209 100644 --- a/core/astro_image.py +++ b/core/astro_image.py @@ -272,6 +272,8 @@ def _extract_metadata(self) -> None: def estimate_background(self, box_size: int = 64) -> None: if self.data is None: raise RuntimeError("Image not loaded") + if self.background is not None: + return # already computed for this instance's data; self.data never changes post-load with warnings.catch_warnings(): warnings.simplefilter("ignore") self.background = Background2D( diff --git a/gui/analysis_thread.py b/gui/analysis_thread.py index 7289ba1..c8f9c56 100644 --- a/gui/analysis_thread.py +++ b/gui/analysis_thread.py @@ -109,6 +109,20 @@ def _execute(self) -> None: else: self.progress.emit(2, "Single-image mode — skipping alignment…") + # Pre-compute background once per distinct image object, ahead of the parallel + # dispatch below. Every analyzer (SNR, PSF, Halo, Edge, PowerSpectrum, Spatial) + # calls estimate_background() unconditionally on whichever image it's given; + # with parallel=True several of them can share the same img_a/img_b/starless + # instance and race to recompute the same Background2D concurrently. Computing + # it here first makes every downstream call a no-op (see the idempotency guard + # in AstroImage.estimate_background) instead of a redundant recomputation. + _bg_images = {id(im): im for im in + (img_a, img_b, self._starless_a, self._starless_b) + if im is not None} + if _bg_images: + with concurrent.futures.ThreadPoolExecutor(max_workers=len(_bg_images)) as ex: + list(ex.map(lambda im: im.estimate_background(), _bg_images.values())) + # Build ordered task list: (metric_key, display_label, callable) # Each callable is a zero-arg function that writes into result_a / result_b. tasks: list[tuple[str, str, Callable[[], None]]] = [] diff --git a/report/report_builder.py b/report/report_builder.py index e358457..f40a1ed 100644 --- a/report/report_builder.py +++ b/report/report_builder.py @@ -449,8 +449,9 @@ def _draw_boxwhisker(ax, vals_list): ("weber_9px", "Weber contrast — 9 px"), ] -# Same order/labels as above minus "original" (no kernel scale) — used for the -# Section 8g per-scale correlation scatter plots. +# Same order/labels as above minus "original" (no kernel scale) — used to order the +# per-scale correlation scatter plots interleaved into Section 8b-8f, each one placed +# immediately after its corresponding map figure. _SPATIAL_CORR_ROWS = [(k, label) for k, label in _SPATIAL_DIFF_DIST_ROWS if k != "original"] @@ -557,7 +558,8 @@ def _draw_boxwhisker(ax, vals_list): "Each row's x-axis is independently clipped to its own 1st–99th percentile " "range so the IQR box stays visible; rows with rare extreme-outlier log-ratios (e.g. " "Weber contrast near dark-sky pixels, see 8e) may have a small fraction of the " - "violin's tail extend beyond the visible axis — see 8g for the full, unclipped " + "violin's tail extend beyond the visible axis — see the per-pixel correlation " + "scatter next to each map figure below (8b–8f) for the full, unclipped " "upper-tail behaviour." "

" ) @@ -677,7 +679,7 @@ def _nc_ratio_rows(score_a: dict, score_b: dict, ratio: dict, scale_label, val_f " Contrast (Weber's law, ΔL/L)" ' Local range vs. local median background — the most literal "Contrast" metric in ' ' this section' - ' Noise-corrected (NC) score (8b–8f, 8h)same scale as parent method' + ' Noise-corrected (NC) score (8b–8f, 8g)same scale as parent method' ' SNR of Detail' " Whether a detail-map response in the nebula is real structure or just this image's " ' own noise floor at that scale' @@ -4022,8 +4024,9 @@ def _section_spatial(self, ra: AnalysisResult, rb: AnalysisResult) -> str: "Background if it's below 0.5× the " "RMS (a top-5% percentile fallback applies if no pixel clears the nebula threshold). " "Pixels between the two thresholds are left unclassified (plain grayscale). " - "How they're used below. The log-ratio distributions and per-scale " - "correlation plots (8g) use the two-image intersection of these masks — " + "How they're used below. The log-ratio distributions and the per-scale " + "correlation plots embedded in 8b–8f use the two-image intersection " + "of these masks — " "a pixel counts as Nebula only if both images classify it as Nebula (same for " "Background). This is deliberately conservative: it excludes pixels where one " f"image's registration, PSF, or local noise disagrees with the other's. Shown on " @@ -4105,7 +4108,7 @@ def _section_spatial(self, ra: AnalysisResult, rb: AnalysisResult) -> str: 'after accounting for each image\'s own noise level. Scale units differ by ' 'method (kernel px for std/Weber, Gaussian σ px for LoG/gradient, ≈2level ' 'px for wavelet) and ratios should not be compared numerically across methods — ' - 'see 8h for a cross-method overview. See also Section 7 for the frequency-domain ' + 'see 8g for a cross-method overview. See also Section 7 for the frequency-domain ' 'view of this same question. Maps below are also shown in noise-normalised form ' '(map ÷ noise floor) so that a shared colour scale is a fair visual comparison ' 'between A and B, even when their absolute noise levels differ.', @@ -4124,45 +4127,56 @@ def figs_for(prefix): out += _hires_img_tag(figs[key], key) + "\n" return out - def corr_figs_html() -> str: - """Emit each present correlation scatter figure in _SPATIAL_CORR_ROWS - order (not alphabetic — e.g. std_10px must not sort before std_3px).""" + def _family_figs_with_corr(rows, map_key_fn) -> str: + """Emit each row's raw map figure immediately followed by its + per-pixel correlation scatter (when present), in _SPATIAL_CORR_ROWS + (numeric-scale) order — not figs_for's alphabetic key sort, which + would put e.g. std_10px before std_3px.""" out = "" - for key, label in _SPATIAL_CORR_ROWS: - fig = figs.get(f"corr_{key}") + for key, _label in rows: + map_key = map_key_fn(key) + fig = figs.get(map_key) if fig: - out += f"

{label}

" + _hires_img_tag(fig, f"corr_{key}") + "\n" + out += _hires_img_tag(fig, map_key) + "\n" + corr_fig = figs.get(f"corr_{key}") + if corr_fig: + out += _hires_img_tag(corr_fig, f"corr_{key}") + "\n" return out - _corr_figs = corr_figs_html() - corr_section_html = ( - "

8g. Per-Pixel Metric Correlation (A vs. B)

" - + _info_box('Each row plots the raw (pre-ratio) metric value at every shared pixel: ' - f'{ra.label} on the y-axis vs. {rb.label} on the x-axis, split into a ' - 'Nebula panel and a Background panel using the shared masks explained in 8a. ' - 'The black dashed diagonal is the 1:1 line (perfect agreement, ' - 'A = B). Axis limits span the full data range for that panel — ' - 'unlike the 8a violin plots, they are not clipped to a percentile — so ' - 'the behaviour of the upper tail is always visible. ' - '

Two fitted lines are drawn on every panel: a solid ' - 'orange overall fit through every ' - 'point, and a dotted crimson tail fit ' - 'restricted to the brightest/most-structured points (the top fraction by ' - 'combined A+B magnitude — see the legend for the exact percentage). A ' - 'whole-population slope is dominated by leverage (distance from the mean), not ' - 'point count, so it blends "the bulk agrees near 1:1" and "the tail diverges" ' - 'into one ambiguous number. The tail fit\'s slope is the more direct ' - 'answer to which filter shows disproportionately more detail/contrast in the ' - 'brightest, most-structured pixels — a tail slope far from 1 with a ' - 'high R² is strong, consistent evidence of a real difference; a tail slope near ' - '1 (even if the point cloud looks scattered) means the apparent divergence is ' - 'mostly noise, not a systematic effect.

' - 'Point clouds are randomly subsampled for rendering; the axis range, both fits, ' - 'and each panel\'s point count (n) are always computed from the full, unsampled ' - 'population.', - title="Per-pixel correlation methodology") - + _corr_figs - if _corr_figs else "" + _std_rows = [(k, l) for k, l in _SPATIAL_CORR_ROWS if k.startswith("std_")] + _log_rows = [(k, l) for k, l in _SPATIAL_CORR_ROWS if k.startswith("log_")] + _gradient_rows = [(k, l) for k, l in _SPATIAL_CORR_ROWS if k.startswith("gradient_")] + _wavelet_rows = [(k, l) for k, l in _SPATIAL_CORR_ROWS if k.startswith("wavelet_")] + _weber_rows = [(k, l) for k, l in _SPATIAL_CORR_ROWS if k.startswith("weber_")] + + _has_any_corr = any(f"corr_{k}" in figs for k, _l in _SPATIAL_CORR_ROWS) + corr_methodology_box = ( + _info_box('Each per-pixel correlation scatter (embedded next to its map figure ' + 'below) plots the raw (pre-ratio) metric value at every shared pixel: ' + f'{ra.label} on the y-axis vs. {rb.label} on the x-axis, split into a ' + 'Nebula panel and a Background panel using the shared masks explained in 8a. ' + 'The black dashed diagonal is the 1:1 line (perfect agreement, ' + 'A = B). Axis limits span the full data range for that panel — ' + 'unlike the 8a violin plots, they are not clipped to a percentile — so ' + 'the behaviour of the upper tail is always visible. ' + '

Two fitted lines are drawn on every panel: a solid ' + 'orange overall fit through every ' + 'point, and a dotted crimson tail fit ' + 'restricted to the brightest/most-structured points (the top fraction by ' + 'combined A+B magnitude — see the legend for the exact percentage). A ' + 'whole-population slope is dominated by leverage (distance from the mean), not ' + 'point count, so it blends "the bulk agrees near 1:1" and "the tail diverges" ' + 'into one ambiguous number. The tail fit\'s slope is the more direct ' + 'answer to which filter shows disproportionately more detail/contrast in the ' + 'brightest, most-structured pixels — a tail slope far from 1 with a ' + 'high R² is strong, consistent evidence of a real difference; a tail slope near ' + '1 (even if the point cloud looks scattered) means the apparent divergence is ' + 'mostly noise, not a systematic effect.

' + 'Point clouds are randomly subsampled for rendering; the axis range, both fits, ' + 'and each panel\'s point count (n) are always computed from the full, unsampled ' + 'population.', + title="Per-pixel correlation methodology") + if _has_any_corr else "" ) has_crosshair = sm.get("crosshair") is not None @@ -4170,8 +4184,7 @@ def corr_figs_html() -> str: 'ℹ When a cross-section line is set in the viewer, its profile is embedded ' 'as the bottom-right panel of each map-pair figure below (8b–8f): Image A ' '(top-left), Image B (top-right), log-ratio map (bottom-left), cross-section ' - 'profile (bottom-right). Left axis: both images (steelblue = A, tomato = B). ' - 'Right axis (green dashed): difference A−B.', + 'profile (bottom-right, steelblue = A, tomato = B).', title="Cross-section profiles", open=True, ) if has_crosshair else "" @@ -4252,6 +4265,7 @@ def corr_figs_html() -> str: {dist_html} {nc_methodology_box} {nc_empty_note} +{corr_methodology_box} {xs_note}

8b. Local Standard Deviation Maps

@@ -4274,10 +4288,11 @@ def corr_figs_html() -> str: Scale{ra.label} (NC score){rb.label} (NC score)Ratio A/B {std_nc_rows} -{figs_for("std_")} +{_family_figs_with_corr(_std_rows, lambda k: k)}

Local σ maps at each kernel size (shared colour scale): Image A (top-left), Image B (top-right), log-ratio map (bottom-left) highlighting where one filter preserves more -local variation, and — when a cross-section line is set — its profile (bottom-right).

+local variation, and — when a cross-section line is set — its profile (bottom-right). Each map +is immediately followed by its per-pixel correlation scatter (see methodology above).

Noise-normalised (× noise floor) — shared colour scale is a fair A/B comparison.

{figs_for("nrm_std_")}

8c. Laplacian of Gaussian (LoG) Maps

@@ -4296,11 +4311,12 @@ def corr_figs_html() -> str: Scale{ra.label} (NC score){rb.label} (NC score)Ratio A/B {log_nc_rows} -{figs_for("log_")} +{_family_figs_with_corr(_log_rows, lambda k: "log_sigma" + k.split("_", 1)[1])}

|LoG| maps at σ = 1.5, 3, and 6 px (shared colour scale per figure): Image A (top-left), Image B (top-right), log-ratio map (bottom-left), and — when a cross-section line is set — its profile (bottom-right). A filter preserving more fine -detail shows brighter, more defined boundaries at small σ.

+detail shows brighter, more defined boundaries at small σ. Each map is immediately +followed by its per-pixel correlation scatter (see methodology above).

Noise-normalised (× noise floor) — shared colour scale is a fair A/B comparison.

{figs_for("nrm_log_")}

8d. Wavelet Decomposition

@@ -4319,12 +4335,13 @@ def corr_figs_html() -> str: {wavelet_nc_rows} -{figs_for("wavelet_level")} +{_family_figs_with_corr(_wavelet_rows, lambda k: "wavelet_level" + k.split("_", 1)[1])}

Reconstructed detail images at levels 2 and 3 (shared colour scale, diverging colourmap): Image A (top-left), Image B (top-right), log-ratio panel (bottom-left) showing where fine structure differs between the two filters (sign discarded — |A|/|B| — since wavelet reconstructions can be negative), and — when a cross-section line is set — -its profile (bottom-right).

+its profile (bottom-right). Each map is immediately followed by its per-pixel correlation +scatter (see methodology above).

Noise-normalised (× noise floor) — shared colour scale is a fair A/B comparison.

{figs_for("nrm_wavelet_")} @@ -4359,13 +4376,14 @@ def corr_figs_html() -> str: Scale{ra.label} (NC score){rb.label} (NC score)Ratio A/B {weber_nc_rows} -{figs_for("weber_")} +{_family_figs_with_corr(_weber_rows, lambda k: k)}

Per-pixel Weber fraction contrast maps (c = ΔL / L, square-root colour scale, viridis): Image A (top-left), Image B (top-right), log-ratio panel (bottom-left) showing where one image achieves greater relative contrast, and — when a cross-section line is set — its profile (bottom-right). Brighter regions have higher Weber contrast — the local intensity range is large relative to the local median luminance. High values -over dark-sky regions are expected; use a nebula ROI for meaningful filter comparison.

+over dark-sky regions are expected; use a nebula ROI for meaningful filter comparison. +Each map is immediately followed by its per-pixel correlation scatter (see methodology above).

Noise-normalised (× noise floor) — shared colour scale is a fair A/B comparison.

{figs_for("nrm_weber_")} @@ -4387,17 +4405,16 @@ def corr_figs_html() -> str: Scale{ra.label} (NC score){rb.label} (NC score)Ratio A/B {gm_nc_rows} -{figs_for("gradient_")} +{_family_figs_with_corr(_gradient_rows, lambda k: k)}

Gradient magnitude maps at σ = 1.5, 3, and 6 px (shared colour scale per figure): Image A (top-left), Image B (top-right), log-ratio map (bottom-left), and — when a cross-section line is set — its profile (bottom-right). A filter preserving sharper -boundaries shows brighter, more defined gradient response.

+boundaries shows brighter, more defined gradient response. Each map is immediately followed +by its per-pixel correlation scatter (see methodology above).

Noise-normalised (× noise floor) — shared colour scale is a fair A/B comparison.

{figs_for("nrm_gradient_")} -{corr_section_html} - -

8h. Noise-Corrected Contrast — Cross-Method Overview

+

8g. Noise-Corrected Contrast — Cross-Method Overview

{_hires_img_tag(figs.get("nc_ratio_overview"), "NC ratio overview")}

Ratio A/B for every noise-corrected method plotted against its approximate spatial scale. Scale units differ by method (see 8b–8f methodology From 8d9219f2bccd1b29b2cf651ef4df66426acaccfb Mon Sep 17 00:00:00 2001 From: Brent <52629076+brentmantooth@users.noreply.github.com> Date: Wed, 15 Jul 2026 20:18:24 -0400 Subject: [PATCH 6/8] Clear and warn on a stale ROI that no longer fits the loaded images Section 8 was crashing with "index -1 is out of bounds for axis 0 with size 0" whenever a previously-drawn ROI outlived the image pair it was drawn on: MainWindow._roi is never cleared on a new image load, so a stale, now out-of-bounds ROI silently slices to a zero-size array (NumPy doesn't raise on an out-of-range slice), which only fails much later and far from the cause, deep inside np.percentile. Power Spectrum and Edge Detection have the identical unguarded ROI slice, so validate once at the actual boundary (MainWindow._on_run, the only path that constructs AnalysisThread) instead of patching every analyzer: an out-of-bounds ROI is now cleared with an explanatory warning rather than silently corrupting downstream arrays. Co-Authored-By: Claude Sonnet 5 --- CLAUDE.md | 1 + gui/main_window.py | 29 +++++++++++++++++++++++++++++ 2 files changed, 30 insertions(+) diff --git a/CLAUDE.md b/CLAUDE.md index 6ee2888..3ebda62 100644 --- a/CLAUDE.md +++ b/CLAUDE.md @@ -377,6 +377,7 @@ pytest tests/ --cov=analysis,core,synthetic,report --cov-report=html | Whole-population regression slope hides tail-specific divergence | A single OLS fit over an entire A-vs-B scatter is dominated by leverage (distance from the mean), not point count — it can show near-zero correlation even when a sparse high-value tail diverges sharply (or vice versa), because it blends "the bulk agrees" and "the tail diverges" into one ambiguous number. When the signal of interest lives specifically in the tail (Section 8g: "which filter shows more detail in the brightest/most-structured pixels"), fit two lines — an overall fit and a second fit restricted to the top N% by combined magnitude (`_plot_metric_correlation`, `SECTION8_SCATTER_TAIL_PERCENTILE`) — and report slope + R² for each. High tail R² with slope far from 1 is a real, systematic effect; low tail R² means the apparent divergence is mostly noise. | | New Section 8 panel key doesn't need Report Inspector code changes | `gui/report_inspector.py` is fully generic — driven entirely by a companion `_inspector.npz` (raw float32/uint8 arrays) plus an embedded `catalog_json` built in `report_builder.py::_write_inspector_file`. `_panel_display_name`/`_panel_concept` dynamically parse any `panels` dict key prefix, so a new `SpatialDetailAnalyzer` panel family auto-appears in the inspector with zero inspector-side changes. A genuinely new *visual type* is a different story: the inspector only knows how to `imshow` 2D/RGB arrays (side-by-side or slider-reveal), so scatter-style plots (Section 8's `corr_*` correlation figures, interleaved into 8b–8f right after each map figure via `_family_figs_with_corr`) must stay static-HTML-only unless new inspector canvas code is written. | | Renumbering a Section 8 subsection misses caption cross-references | Section 8's sub-heading letters (8a–8g) are referenced by literal string in caption/info-box text scattered throughout `_section_spatial` — not just in the `

` tags (e.g. "see 8g for…", "(8b–8f, 8g)"). After adding, removing, or renumbering a subsection, `grep` the function for every old *and* new heading letter — HTML renders a stale cross-reference without error, it just silently misdirects the reader to the wrong subsection. | +| Stale ROI crashes Section 8 with "index -1 is out of bounds for axis 0 with size 0" | `MainWindow._roi` (`gui/main_window.py`) is never cleared when a new image is loaded into either panel. If the user draws an ROI on one image pair, then loads a smaller replacement pair without clearing it, the stale coordinates go out of bounds for the new image. NumPy doesn't raise on an out-of-range slice — `norm_a[ry0:ry1, rx0:rx1]` silently returns a zero-size array — so the crash surfaces much later and far from the real cause: `SpatialDetailAnalyzer._plot_mask_illustration → _stretch_for_display → np.percentile(empty_array, ...)`. The same unguarded `bgsub[y0:y1, x0:x1]` pattern exists in `power_spectrum.py::_extract_roi` and `edge_analyzer.py::analyze`, so a stale ROI can corrupt those sections too (with different, equally misleading errors) if they happen to run. Fixed at the single real boundary — `MainWindow._on_run()`, which is the only path that constructs `AnalysisThread` — by validating `self._roi` against every loaded image's `data.shape` right before `settings["roi"]` is set; an out-of-bounds ROI is cleared (falls back to auto-detect/full-image) with a `QMessageBox` explaining why, rather than patching each analyzer's slice individually. | --- diff --git a/gui/main_window.py b/gui/main_window.py index 7f653f9..5b5cf7f 100644 --- a/gui/main_window.py +++ b/gui/main_window.py @@ -296,6 +296,35 @@ def _on_run(self, settings: dict) -> None: self._control.set_run_enabled(True) return + # A previously-drawn ROI is never auto-cleared when a new image is loaded, so it + # can silently go out of bounds for a smaller replacement image. An out-of-range + # ROI slices to a zero-size array in every ROI-aware analyzer (Section 8, Power + # Spectrum, Edge Detection) with no exception at the slice itself — it fails much + # later with a confusing "index -1 is out of bounds for axis 0 with size 0" (or + # similar) deep inside a percentile/zoom/regression call. Validate once here so + # every analyzer downstream always receives either None or an ROI that fits. + if self._roi is not None: + rx0, ry0, rx1, ry1 = self._roi + roi_fits = True + for img in (img_a, img_b): + if img is None: + continue + h, w = img.data.shape[:2] + if not (0 <= rx0 < rx1 <= w and 0 <= ry0 < ry1 <= h): + roi_fits = False + break + if not roi_fits: + self._roi = None + self._control.set_roi(None) + QMessageBox.warning( + self, "ROI no longer fits loaded images", + "The previously selected ROI no longer fits within the currently " + "loaded image(s) and has been cleared. Analysis will proceed using " + "the full image (or an auto-detected region) instead.\n\n" + "Use \"Select ROI…\" to redraw a region if you want to restrict " + "analysis to a specific area." + ) + # Merge ROI and crosshair from window state into settings settings["roi"] = self._roi settings["crosshair"] = self._crosshair From bd10854ead2bffbd2f63bc9e72052748a44553f9 Mon Sep 17 00:00:00 2001 From: Brent <52629076+brentmantooth@users.noreply.github.com> Date: Wed, 15 Jul 2026 22:27:31 -0400 Subject: [PATCH 7/8] Add log-ratio histograms to Section 8 maps, declutter correlation scatter Each Section 8 metric-family figure gains a third row: a histogram of the log-ratio map's pixel distribution, colour-matched to the log-ratio panel above it, making the A-vs-B divergence easier to read at a glance. The adjacent per-pixel correlation scatter plots drop their overall/tail regression fit lines (kept the 1:1 reference) in favour of coloring each point by its own log-ratio value on the same bwr scale, linking the scatter back to the map instead of adding two more lines to interpret. Co-Authored-By: Claude Sonnet 5 --- CLAUDE.md | 4 +- analysis/image_filters.py | 144 ++++++++++++++++++++------------------ core/models.py | 1 - report/report_builder.py | 66 ++++++++--------- 4 files changed, 112 insertions(+), 103 deletions(-) diff --git a/CLAUDE.md b/CLAUDE.md index 3ebda62..81743bf 100644 --- a/CLAUDE.md +++ b/CLAUDE.md @@ -60,8 +60,9 @@ synthetic/ | `_annular_rdf(log_data, xc, yc, radius)` | `gui/halo_dialog.py` | 1-px annular mean/std in log10 space; mirrors `HaloAnalyzer._annular_stats` | | `_power_ratio_db(freq_a, rp_a, freq_b, rp_b)` | `report_builder.py` | 10·log10 dB ratio between two radial power curves; returns `None` on missing data or misaligned frequency bins | | `_log_ratio_map(a, b)` | `analysis/image_filters.py` | Per-pixel `log10(\|A\|/\|B\|)` map with percentile-based epsilon floor and defensive shape crop — the Section 8 replacement for plain `A − B` diff | +| `_log_ratio_color_range(diff)` | `analysis/image_filters.py` | Symmetric `(vmin, vmax)` for the `bwr` log-ratio colormap, shared by the log-ratio map panel, its histogram, and the correlation scatter dot coloring | | `_plot_mask_illustration(base, mask_neb, mask_bg)` | `analysis/image_filters.py` | Translucent steelblue/tomato mask overlay on a grayscale base image | -| `_plot_metric_correlation(map_a, map_b, mask_neb, mask_bg, ...)` | `analysis/image_filters.py` | 1×2 masked-region scatter (A vs B) with a 1:1 line plus overall + tail-restricted dual linear fits | +| `_plot_metric_correlation(map_a, map_b, log_ratio, mask_neb, mask_bg, ...)` | `analysis/image_filters.py` | 1×2 masked-region scatter (A vs B) with a 1:1 line; each point colored by its pixel's log-ratio value using the same `bwr` scale as the adjacent map figure | | `_family_figs_with_corr(rows, map_key_fn)` | `report_builder.py` | Emits a Section 8 family's map figure immediately followed by its `corr_*` correlation scatter, one scale at a time, in numeric order (`_SPATIAL_CORR_ROWS`) — the pattern to follow when adding any new per-scale Section 8 figure pair | --- @@ -374,7 +375,6 @@ pytest tests/ --cov=analysis,core,synthetic,report --cov-report=html | Adding a float64 cast in analysis code | Don't. All image data is float32 after `AstroImage.load()`. The only float64 exception is `astroalign` in `gui/analysis_thread.py`. Redundant float64 casts waste memory and defeat the float32 performance gains. | | Mixed float32/float64 arithmetic silently widens to float64 | NumPy upcasts when operands differ (e.g. `float32_array - float64_scalar`). If photutils ever returns a float64 background model, `background_subtracted()` will silently return float64. Guard by adding `.astype(np.float32)` at the end of `background_subtracted()` in `astro_image.py` if this is observed. | | `_section_snr` crashed when SNR metric is unchecked | `_plot_snr_pair` (`report_builder.py`) built a `panels` list filtered to non-`None` entries but never checked whether it was empty before calling `plt.subplots(1, len(panels), ...)` — 0 columns raised `ValueError: Number of columns must be a positive integer, not 0`. Hit whenever SNR is unchecked while another metric (e.g. Power Spectrum) is run. Fixed with an early `if not panels: return None` guard, matching `_plot_radial_overlay`/`_plot_radial_ratio_db`; both call sites already pipe the result through `_img_tag`, which turns `None` into `""`. | -| Whole-population regression slope hides tail-specific divergence | A single OLS fit over an entire A-vs-B scatter is dominated by leverage (distance from the mean), not point count — it can show near-zero correlation even when a sparse high-value tail diverges sharply (or vice versa), because it blends "the bulk agrees" and "the tail diverges" into one ambiguous number. When the signal of interest lives specifically in the tail (Section 8g: "which filter shows more detail in the brightest/most-structured pixels"), fit two lines — an overall fit and a second fit restricted to the top N% by combined magnitude (`_plot_metric_correlation`, `SECTION8_SCATTER_TAIL_PERCENTILE`) — and report slope + R² for each. High tail R² with slope far from 1 is a real, systematic effect; low tail R² means the apparent divergence is mostly noise. | | New Section 8 panel key doesn't need Report Inspector code changes | `gui/report_inspector.py` is fully generic — driven entirely by a companion `_inspector.npz` (raw float32/uint8 arrays) plus an embedded `catalog_json` built in `report_builder.py::_write_inspector_file`. `_panel_display_name`/`_panel_concept` dynamically parse any `panels` dict key prefix, so a new `SpatialDetailAnalyzer` panel family auto-appears in the inspector with zero inspector-side changes. A genuinely new *visual type* is a different story: the inspector only knows how to `imshow` 2D/RGB arrays (side-by-side or slider-reveal), so scatter-style plots (Section 8's `corr_*` correlation figures, interleaved into 8b–8f right after each map figure via `_family_figs_with_corr`) must stay static-HTML-only unless new inspector canvas code is written. | | Renumbering a Section 8 subsection misses caption cross-references | Section 8's sub-heading letters (8a–8g) are referenced by literal string in caption/info-box text scattered throughout `_section_spatial` — not just in the `

` tags (e.g. "see 8g for…", "(8b–8f, 8g)"). After adding, removing, or renumbering a subsection, `grep` the function for every old *and* new heading letter — HTML renders a stale cross-reference without error, it just silently misdirects the reader to the wrong subsection. | | Stale ROI crashes Section 8 with "index -1 is out of bounds for axis 0 with size 0" | `MainWindow._roi` (`gui/main_window.py`) is never cleared when a new image is loaded into either panel. If the user draws an ROI on one image pair, then loads a smaller replacement pair without clearing it, the stale coordinates go out of bounds for the new image. NumPy doesn't raise on an out-of-range slice — `norm_a[ry0:ry1, rx0:rx1]` silently returns a zero-size array — so the crash surfaces much later and far from the real cause: `SpatialDetailAnalyzer._plot_mask_illustration → _stretch_for_display → np.percentile(empty_array, ...)`. The same unguarded `bgsub[y0:y1, x0:x1]` pattern exists in `power_spectrum.py::_extract_roi` and `edge_analyzer.py::analyze`, so a stale ROI can corrupt those sections too (with different, equally misleading errors) if they happen to run. Fixed at the single real boundary — `MainWindow._on_run()`, which is the only path that constructs `AnalysisThread` — by validating `self._roi` against every loaded image's `data.shape` right before `settings["roi"]` is set; an out-of-bounds ROI is cleared (falls back to auto-detect/full-image) with a `QMessageBox` explaining why, rather than patching each analyzer's slice individually. | diff --git a/analysis/image_filters.py b/analysis/image_filters.py index ce81827..2d52d87 100644 --- a/analysis/image_filters.py +++ b/analysis/image_filters.py @@ -12,7 +12,6 @@ import matplotlib.pyplot as plt from matplotlib.patches import Patch from scipy.ndimage import generic_filter, gaussian_filter, gaussian_laplace, gaussian_gradient_magnitude, map_coordinates, zoom, maximum_filter, minimum_filter, median_filter -from scipy.stats import linregress import pywt from core.astro_image import AstroImage @@ -21,8 +20,7 @@ WEBER_KERNEL_SIZES, XS_LINE_ALPHA, SECTION8_BORDER_CROP_FRACTION, SECTION8_ANALYSIS_CMAP, XS_SNR_REGION_WIDTH, SECTION8_DIFF_DIST_MAX_SAMPLES, - SECTION8_LOGRATIO_EPS_PERCENTILE, SECTION8_SCATTER_MAX_SAMPLES, - SECTION8_SCATTER_TAIL_PERCENTILE) + SECTION8_LOGRATIO_EPS_PERCENTILE, SECTION8_SCATTER_MAX_SAMPLES) MAX_DIM_FOR_STD = 2048 # downsample to this before generic_filter (performance) _DISPLAY_SMOOTH_SIGMA = 1.0 # applied to maps before plotting; does NOT affect metrics @@ -412,7 +410,7 @@ def _std_analysis(self, norm_a, norm_b, diff, mask_neb_shared, mask_bg_shared, diff_dist_rng) if not single: corr_fig = self._plot_metric_correlation( - std_a, std_b, mask_neb_shared, mask_bg_shared, + std_a, std_b, diff, mask_neb_shared, mask_bg_shared, label_a, label_b, f"Local σ (kernel {ks}px)", diff_dist_rng) if corr_fig is not None: figures[f"corr_std_{ks}px"] = corr_fig @@ -562,6 +560,13 @@ def _log_ratio_map(a: np.ndarray, b: np.ndarray) -> np.ndarray: ratio = np.clip(abs_a, eps, None) / np.clip(abs_b, eps, None) return np.log10(ratio).astype(np.float32) + @staticmethod + def _log_ratio_color_range(diff: np.ndarray) -> tuple[float, float]: + """Symmetric (vmin, vmax) for the bwr log-ratio colormap, shared by the + log-ratio map panel, its histogram, and the correlation scatter dots.""" + d_max = float(np.percentile(np.abs(diff), 99.5)) or 1.0 + return -d_max, d_max + @staticmethod def _diff_distribution(diff_map: np.ndarray, mask_neb_shared: np.ndarray | None, @@ -644,7 +649,7 @@ def _log_analysis(self, norm_a, norm_b, sigmas, diff, mask_neb_shared, mask_bg_shared, diff_dist_rng) if not single: corr_fig = self._plot_metric_correlation( - log_a, log_b, mask_neb_shared, mask_bg_shared, + log_a, log_b, diff, mask_neb_shared, mask_bg_shared, label_a, label_b, f"|LoG| (σ={sigma}px)", diff_dist_rng) if corr_fig is not None: figures[f"corr_log_{sigma}"] = corr_fig @@ -747,7 +752,7 @@ def _gradient_analysis(self, norm_a, norm_b, sigmas, diff, mask_neb_shared, mask_bg_shared, diff_dist_rng) if not single: corr_fig = self._plot_metric_correlation( - gm_a, gm_b, mask_neb_shared, mask_bg_shared, + gm_a, gm_b, diff, mask_neb_shared, mask_bg_shared, label_a, label_b, f"Gradient |G| (σ={sigma}px)", diff_dist_rng) if corr_fig is not None: figures[f"corr_gradient_{sigma}"] = corr_fig @@ -882,7 +887,7 @@ def _wavelet_analysis(self, norm_a, norm_b, wavelet, levels, # sign-discarding log-ratio map, shows whether band-pass detail # flips sign between the two filters at a given pixel. corr_fig = self._plot_metric_correlation( - rec_a, rec_b, mask_neb_shared, mask_bg_shared, + rec_a, rec_b, diff, mask_neb_shared, mask_bg_shared, label_a, label_b, f"Wavelet level {display_level}", diff_dist_rng) if corr_fig is not None: figures[f"corr_wavelet_{display_level}"] = corr_fig @@ -1019,7 +1024,7 @@ def _weber_analysis(self, norm_a, norm_b, kernel_sizes, diff, mask_neb_shared, mask_bg_shared, diff_dist_rng) if not single: corr_fig = self._plot_metric_correlation( - wc_a, wc_b, mask_neb_shared, mask_bg_shared, + wc_a, wc_b, diff, mask_neb_shared, mask_bg_shared, label_a, label_b, f"Weber contrast (kernel {ks}px)", diff_dist_rng) if corr_fig is not None: figures[f"corr_weber_{ks}px"] = corr_fig @@ -1143,21 +1148,32 @@ def _plot_side_by_side(self, arr_a: np.ndarray, arr_b: np.ndarray, diff = self._log_ratio_map(arr_a, arr_b) # Symmetric about zero so the "bwr" midpoint (white) always means no difference. - d_max = float(np.percentile(np.abs(diff), 99.5)) or 1.0 - dvmin, dvmax = -d_max, d_max + dvmin, dvmax = self._log_ratio_color_range(diff) - # 2×2 grid: A|B on top, log-ratio diff | cross-section on bottom. A/B/diff - # share the source array's pixel aspect ratio (aspect="equal" imshow); the - # cross-section panel is a line plot with no such constraint but occupies - # an equal-size grid cell so the other three panels stay geometrically - # identical whether or not a crosshair (and thus xs_data) is set. + # Dark-mode-aware reference-line color (project convention). + is_dark = matplotlib.rcParams.get("figure.facecolor", "white") not in ("white", "#ffffff", 1.0) + orig_color = "white" if is_dark else "black" + + # 3-row grid: A|B on top, log-ratio diff | cross-section in the middle, + # a log-ratio pixel-value histogram spanning both columns on the bottom. + # A/B/diff share the source array's pixel aspect ratio (aspect="equal" + # imshow); the cross-section panel is a line plot with no such constraint + # but occupies an equal-size grid cell so the other three panels stay + # geometrically identical whether or not a crosshair (and thus xs_data) + # is set. h, w = arr_a.shape[:2] aspect_ratio = h / max(w, 1) panel_w = 5.0 # half the old single-column width — 2 columns now share it panel_h = panel_w * aspect_ratio - fig, axes = plt.subplots(2, 2, figsize=(panel_w * 2, panel_h * 2 + 1.5), - constrained_layout=True) - (ax_a, ax_b), (ax_diff, ax_xs) = axes + hist_h = 2.2 + fig = plt.figure(figsize=(panel_w * 2, panel_h * 2 + hist_h + 1.5), + constrained_layout=True) + gs = fig.add_gridspec(3, 2, height_ratios=[panel_h, panel_h, hist_h]) + ax_a = fig.add_subplot(gs[0, 0]) + ax_b = fig.add_subplot(gs[0, 1]) + ax_diff = fig.add_subplot(gs[1, 0]) + ax_xs = fig.add_subplot(gs[1, 1]) + ax_hist = fig.add_subplot(gs[2, :]) for ax, arr, title in zip([ax_a, ax_b], [arr_a, arr_b], [title_a, title_b]): im = ax.imshow(arr, origin="upper", cmap=cmap, @@ -1183,6 +1199,26 @@ def _plot_side_by_side(self, arr_a: np.ndarray, arr_b: np.ndarray, else: ax_xs.axis("off") + # Histogram of the log-ratio map's pixel distribution, colored to match + # the diff panel above: full data range on the x-axis (no pixels hidden), + # but each bin's fill color is clipped to [dvmin, dvmax] so extreme-tail + # bins saturate to the same end colors imshow already uses for its own + # outliers. + counts, bin_edges, patches = ax_hist.hist(diff.ravel(), bins=60) + hist_norm = mcolors.Normalize(vmin=dvmin, vmax=dvmax, clip=True) + hist_cmap = plt.get_cmap("bwr") + bin_centers = 0.5 * (bin_edges[:-1] + bin_edges[1:]) + for patch, center in zip(patches, bin_centers): + patch.set_facecolor(hist_cmap(hist_norm(center))) + ax_hist.axvline(0.0, color=orig_color, linestyle="--", linewidth=1.0, label="A = B") + ax_hist.set_yscale("log") + ax_hist.set_xlabel("Log ratio, log10(|A|/|B|)", fontsize=8) + ax_hist.set_ylabel("Pixel count", fontsize=8) + ax_hist.tick_params(labelsize=7) + ax_hist.legend(fontsize=6.5, loc="upper right") + ax_hist.grid(True, alpha=0.3) + ax_hist.set_title("Log-ratio pixel distribution", fontsize=9) + return fig def _plot_single(self, arr_a: np.ndarray, title_a: str, @@ -1358,31 +1394,31 @@ def _draw_cross_section(ax, pos: np.ndarray, prof_a: np.ndarray, prof_b: np.ndar @staticmethod def _plot_metric_correlation(map_a: np.ndarray, map_b: np.ndarray, + log_ratio: np.ndarray, mask_neb_shared: np.ndarray, mask_bg_shared: np.ndarray, label_a: str, label_b: str, metric_title: str, rng: np.random.Generator, max_samples: int = SECTION8_SCATTER_MAX_SAMPLES) -> plt.Figure | None: """1x2 correlation scatter (Nebula | Background): metric value in A (y) - vs. metric value in B (x), with a dashed 1:1 reference line plus two - linear fits — an overall OLS fit over all points, and a second fit - restricted to the top (100 - SECTION8_SCATTER_TAIL_PERCENTILE)% of - points by combined magnitude ("tail fit"). The tail fit isolates the - high-detail regime's slope from the bulk's, since a whole-population - slope is dominated by leverage rather than point count and conflates - the two regimes into one ambiguous number. + vs. metric value in B (x), with a dashed 1:1 reference line. Each point + is colored by that pixel's log-ratio value (log_ratio — the same array + driving the adjacent log-ratio map panel), using the same bwr colormap + and range, so the scatter visually links back to the map. Axis limits reflect the FULL pooled population range per subplot (not percentile-clipped like the Section 8a violin plots) so upper-tail divergence from the 1:1 line — the signal this plot exists to surface — stays visible. Point clouds are randomly subsampled (up to max_samples) - for render cost only; the axis range and both fits are always computed - from the full, unsampled population. Returns None if both subplots have - too few points. + for render cost only; the axis range is always computed from the full, + unsampled population. Returns None if both subplots have too few points. """ - h = min(map_a.shape[0], map_b.shape[0], mask_neb_shared.shape[0], mask_bg_shared.shape[0]) - w = min(map_a.shape[1], map_b.shape[1], mask_neb_shared.shape[1], mask_bg_shared.shape[1]) + h = min(map_a.shape[0], map_b.shape[0], log_ratio.shape[0], + mask_neb_shared.shape[0], mask_bg_shared.shape[0]) + w = min(map_a.shape[1], map_b.shape[1], log_ratio.shape[1], + mask_neb_shared.shape[1], mask_bg_shared.shape[1]) map_a = map_a[:h, :w] map_b = map_b[:h, :w] + log_ratio = log_ratio[:h, :w] mask_neb_shared = mask_neb_shared[:h, :w] mask_bg_shared = mask_bg_shared[:h, :w] @@ -1390,12 +1426,17 @@ def _plot_metric_correlation(map_a: np.ndarray, map_b: np.ndarray, is_dark = matplotlib.rcParams.get("figure.facecolor", "white") not in ("white", "#ffffff", 1.0) orig_color = "white" if is_dark else "black" + # Shared across both panels so a given color always means the same + # log-ratio value, matching the adjacent map figure's own scale. + dvmin, dvmax = SpatialDetailAnalyzer._log_ratio_color_range(log_ratio) + fig, axes = plt.subplots(1, 2, figsize=(9, 4.5)) any_data = False for ax, region_name, mask in zip(axes, ("Nebula", "Background"), (mask_neb_shared, mask_bg_shared)): a_vals = map_a[mask] b_vals = map_b[mask] + c_vals = log_ratio[mask] n = a_vals.size if n < 3: ax.set_visible(False) @@ -1410,46 +1451,15 @@ def _plot_metric_correlation(map_a: np.ndarray, map_b: np.ndarray, if n > max_samples: idx = rng.choice(n, max_samples, replace=False) - a_plot, b_plot = a_vals[idx], b_vals[idx] + a_plot, b_plot, c_plot = a_vals[idx], b_vals[idx], c_vals[idx] else: - a_plot, b_plot = a_vals, b_vals - - # Linear fits computed on the FULL population (not the render subsample) - # for statistical accuracy. Two fits, not one: an overall OLS fit over - # all points is dominated by leverage, not point count, so it conflates - # "the bulk agrees near 1:1" with "the tail diverges" into one ambiguous - # number. A second fit restricted to the top (100 - - # SECTION8_SCATTER_TAIL_PERCENTILE)% of points by combined magnitude - # (A + B, a proxy for "how far out along the diagonal") isolates the - # high-detail regime specifically — that slope is the more direct - # answer to "does one filter show disproportionately more detail/ - # contrast in the brightest/most-structured regions". - overall_fit = linregress(b_vals, a_vals) - m_all, c_all, r_all = overall_fit.slope, overall_fit.intercept, overall_fit.rvalue - - combined = a_vals + b_vals - tail_thresh = np.percentile(combined, SECTION8_SCATTER_TAIL_PERCENTILE) - tail_mask = combined >= tail_thresh - tail_fit = None - if np.count_nonzero(tail_mask) >= 3: - tf = linregress(b_vals[tail_mask], a_vals[tail_mask]) - tail_fit = (tf.slope, tf.intercept, tf.rvalue) - - ax.scatter(b_plot, a_plot, alpha=0.55, s=8, color="steelblue", - zorder=3, edgecolors="none", rasterized=True) + a_plot, b_plot, c_plot = a_vals, b_vals, c_vals + + sc = ax.scatter(b_plot, a_plot, c=c_plot, cmap="bwr", vmin=dvmin, vmax=dvmax, + alpha=0.55, s=8, zorder=3, edgecolors="none", rasterized=True) ax.plot([lo, hi], [lo, hi], color=orig_color, linestyle="--", linewidth=1.2, zorder=4, label="Slope = 1 (A = B)") - fit_x = np.array([lo, hi]) - ax.plot(fit_x, m_all * fit_x + c_all, color="darkorange", linestyle="-", - linewidth=1.6, zorder=5, - label=f"Fit: y={m_all:.2f}x{c_all:+.2f} (R²={r_all ** 2:.2f})") - if tail_fit is not None: - m_t, c_t, r_t = tail_fit - tail_pct = 100 - SECTION8_SCATTER_TAIL_PERCENTILE - ax.plot(fit_x, m_t * fit_x + c_t, color="crimson", linestyle=":", - linewidth=1.8, zorder=6, - label=f"Tail fit (top {tail_pct:.0f}%): " - f"y={m_t:.2f}x{c_t:+.2f} (R²={r_t ** 2:.2f})") + fig.colorbar(sc, ax=ax, fraction=0.046, pad=0.04, label="log10(|A|/|B|)") ax.set_xlim(lo, hi) ax.set_ylim(lo, hi) ax.set_aspect("equal") diff --git a/core/models.py b/core/models.py index 31fd51d..f413f67 100644 --- a/core/models.py +++ b/core/models.py @@ -43,7 +43,6 @@ SECTION8_DIFF_DIST_MAX_SAMPLES = 1000000 # per masked population, per scale — caps violin/KDE cost on full-res diff maps SECTION8_LOGRATIO_EPS_PERCENTILE = 1.0 # percentile of pooled positive |A|,|B| values used as the epsilon floor in log10(|A|/|B|) SECTION8_SCATTER_MAX_SAMPLES = 50000 # per masked population, per scale — caps render cost of Section 8g correlation scatter plots -SECTION8_SCATTER_TAIL_PERCENTILE = 80 # percentile of (A+B) combined magnitude above which points are fit separately as the "tail" regression in Section 8g PSF_SPATIAL_MAP_SIZE = 150 # px; long-axis resolution of FWHM / eccentricity spatial maps PSF_SPATIAL_MAP_SMOOTH_SIGMA = 5.0 # Gaussian smoothing sigma (px) applied to spatial maps before display diff --git a/report/report_builder.py b/report/report_builder.py index f40a1ed..b6d9c91 100644 --- a/report/report_builder.py +++ b/report/report_builder.py @@ -4159,21 +4159,14 @@ def _family_figs_with_corr(rows, map_key_fn) -> str: 'A = B). Axis limits span the full data range for that panel — ' 'unlike the 8a violin plots, they are not clipped to a percentile — so ' 'the behaviour of the upper tail is always visible. ' - '

Two fitted lines are drawn on every panel: a solid ' - 'orange overall fit through every ' - 'point, and a dotted crimson tail fit ' - 'restricted to the brightest/most-structured points (the top fraction by ' - 'combined A+B magnitude — see the legend for the exact percentage). A ' - 'whole-population slope is dominated by leverage (distance from the mean), not ' - 'point count, so it blends "the bulk agrees near 1:1" and "the tail diverges" ' - 'into one ambiguous number. The tail fit\'s slope is the more direct ' - 'answer to which filter shows disproportionately more detail/contrast in the ' - 'brightest, most-structured pixels — a tail slope far from 1 with a ' - 'high R² is strong, consistent evidence of a real difference; a tail slope near ' - '1 (even if the point cloud looks scattered) means the apparent divergence is ' - 'mostly noise, not a systematic effect.

' - 'Point clouds are randomly subsampled for rendering; the axis range, both fits, ' - 'and each panel\'s point count (n) are always computed from the full, unsampled ' + 'Each point is colored by that pixel\'s log-ratio value ' + '(log10(|A|/|B|)), using the same bwr colour scale as the ' + 'adjacent log-ratio map panel and its histogram — red points are pixels where A ' + 'is stronger, blue where B is stronger, white where the two agree — so the ' + 'scatter\'s spatial pattern (visible in the map) links directly to its position ' + 'in this plot. ' + 'Point clouds are randomly subsampled for rendering; the axis range and each ' + 'panel\'s point count (n) are always computed from the full, unsampled ' 'population.', title="Per-pixel correlation methodology") if _has_any_corr else "" @@ -4182,9 +4175,9 @@ def _family_figs_with_corr(rows, map_key_fn) -> str: has_crosshair = sm.get("crosshair") is not None xs_note = _info_box( 'ℹ When a cross-section line is set in the viewer, its profile is embedded ' - 'as the bottom-right panel of each map-pair figure below (8b–8f): Image A ' - '(top-left), Image B (top-right), log-ratio map (bottom-left), cross-section ' - 'profile (bottom-right, steelblue = A, tomato = B).', + 'as the middle-right panel of each map-pair figure below (8b–8f): Image A ' + '(top-left), Image B (top-right), log-ratio map (middle-left), cross-section ' + 'profile (middle-right, steelblue = A, tomato = B).', title="Cross-section profiles", open=True, ) if has_crosshair else "" @@ -4252,15 +4245,17 @@ def _family_figs_with_corr(rows, map_key_fn) -> str: {smooth_note} {_info_box('All maps below are computed on mean-signal-normalised data ' '(each image divided by its own mean signal), making them dimensionless and comparable ' - 'across different filter bandwidths. Each figure is a 2×2 grid: Image A ' - '(top-left) and Image B (top-right) share a colour scale; the bottom-left panel shows ' + 'across different filter bandwidths. Each figure has Image A (top-left) and Image B ' + '(top-right) sharing a colour scale; the middle-left panel shows ' 'the per-pixel log-ratio log10(A / B): red (positive) means A is stronger at ' 'that pixel, blue (negative) means B is stronger, white means the two are equal. ±0.3 ≈ ' 'a 2× difference, ±1.0 ≈ a 10× difference. For the Original and Wavelet maps ' '— the only two families whose values can be negative — the ratio uses |A| and |B|, so ' 'it compares the magnitude of structure rather than signed brightness. The ' - 'bottom-right panel shows the cross-section profile along the selected line when one is ' - 'set, otherwise left blank.', + 'middle-right panel shows the cross-section profile along the selected line when one is ' + 'set, otherwise left blank. A third row spans both columns with a histogram of the ' + 'log-ratio map\'s pixel values, coloured bin-by-bin with the same bwr scale and range as ' + 'the log-ratio panel above it.', title="Spatial detail maps overview")} {dist_html} {nc_methodology_box} @@ -4277,7 +4272,7 @@ def _family_figs_with_corr(rows, map_key_fn) -> str: 'window. Brighter regions contain more local variation — typically nebula filaments, ' 'star halos, or noise. A filter with higher std values in targeted emission regions ' 'preserves more structure; higher std in blank sky regions indicates more photon noise. ' - 'When a cross-section line is set, its profile is embedded in the bottom-right panel of ' + 'When a cross-section line is set, its profile is embedded in the middle-right panel of ' 'each map figure below, showing how local detail amplitude varies along the selected line.', title="Local standard deviation")} @@ -4290,8 +4285,9 @@ def _family_figs_with_corr(rows, map_key_fn) -> str:
{_family_figs_with_corr(_std_rows, lambda k: k)}

Local σ maps at each kernel size (shared colour scale): Image A (top-left), -Image B (top-right), log-ratio map (bottom-left) highlighting where one filter preserves more -local variation, and — when a cross-section line is set — its profile (bottom-right). Each map +Image B (top-right), log-ratio map (middle-left) highlighting where one filter preserves more +local variation, and — when a cross-section line is set — its profile (middle-right), plus a +bottom-row histogram of the log-ratio map's pixel values (same colour scale). Each map is immediately followed by its per-pixel correlation scatter (see methodology above).

Noise-normalised (× noise floor) — shared colour scale is a fair A/B comparison.

{figs_for("nrm_std_")} @@ -4313,8 +4309,9 @@ def _family_figs_with_corr(rows, map_key_fn) -> str: {_family_figs_with_corr(_log_rows, lambda k: "log_sigma" + k.split("_", 1)[1])}

|LoG| maps at σ = 1.5, 3, and 6 px (shared colour scale per figure): -Image A (top-left), Image B (top-right), log-ratio map (bottom-left), and — when a -cross-section line is set — its profile (bottom-right). A filter preserving more fine +Image A (top-left), Image B (top-right), log-ratio map (middle-left), and — when a +cross-section line is set — its profile (middle-right), plus a bottom-row histogram of +the log-ratio map's pixel values (same colour scale). A filter preserving more fine detail shows brighter, more defined boundaries at small σ. Each map is immediately followed by its per-pixel correlation scatter (see methodology above).

Noise-normalised (× noise floor) — shared colour scale is a fair A/B comparison.

@@ -4337,10 +4334,11 @@ def _family_figs_with_corr(rows, map_key_fn) -> str: {_family_figs_with_corr(_wavelet_rows, lambda k: "wavelet_level" + k.split("_", 1)[1])}

Reconstructed detail images at levels 2 and 3 (shared colour scale, -diverging colourmap): Image A (top-left), Image B (top-right), log-ratio panel (bottom-left) +diverging colourmap): Image A (top-left), Image B (top-right), log-ratio panel (middle-left) showing where fine structure differs between the two filters (sign discarded — |A|/|B| — since wavelet reconstructions can be negative), and — when a cross-section line is set — -its profile (bottom-right). Each map is immediately followed by its per-pixel correlation +its profile (middle-right), plus a bottom-row histogram of the log-ratio map's pixel values +(same colour scale). Each map is immediately followed by its per-pixel correlation scatter (see methodology above).

Noise-normalised (× noise floor) — shared colour scale is a fair A/B comparison.

{figs_for("nrm_wavelet_")} @@ -4378,9 +4376,10 @@ def _family_figs_with_corr(rows, map_key_fn) -> str: {_family_figs_with_corr(_weber_rows, lambda k: k)}

Per-pixel Weber fraction contrast maps (c = ΔL / L, square-root colour -scale, viridis): Image A (top-left), Image B (top-right), log-ratio panel (bottom-left) +scale, viridis): Image A (top-left), Image B (top-right), log-ratio panel (middle-left) showing where one image achieves greater relative contrast, and — when a cross-section -line is set — its profile (bottom-right). Brighter regions have higher Weber contrast — +line is set — its profile (middle-right), plus a bottom-row histogram of the log-ratio +map's pixel values (same colour scale). Brighter regions have higher Weber contrast — the local intensity range is large relative to the local median luminance. High values over dark-sky regions are expected; use a nebula ROI for meaningful filter comparison. Each map is immediately followed by its per-pixel correlation scatter (see methodology above).

@@ -4407,8 +4406,9 @@ def _family_figs_with_corr(rows, map_key_fn) -> str: {_family_figs_with_corr(_gradient_rows, lambda k: k)}

Gradient magnitude maps at σ = 1.5, 3, and 6 px (shared colour scale per -figure): Image A (top-left), Image B (top-right), log-ratio map (bottom-left), and — when a -cross-section line is set — its profile (bottom-right). A filter preserving sharper +figure): Image A (top-left), Image B (top-right), log-ratio map (middle-left), and — when a +cross-section line is set — its profile (middle-right), plus a bottom-row histogram of the +log-ratio map's pixel values (same colour scale). A filter preserving sharper boundaries shows brighter, more defined gradient response. Each map is immediately followed by its per-pixel correlation scatter (see methodology above).

Noise-normalised (× noise floor) — shared colour scale is a fair A/B comparison.

From 22d9bd4bca83e9593ccfe5365495bd5433be1c67 Mon Sep 17 00:00:00 2001 From: Brent <52629076+brentmantooth@users.noreply.github.com> Date: Wed, 15 Jul 2026 22:27:48 -0400 Subject: [PATCH 8/8] setup --- .claude/settings.json | 5 ++++- AstroImageLab.py | 4 ++-- 2 files changed, 6 insertions(+), 3 deletions(-) diff --git a/.claude/settings.json b/.claude/settings.json index b15e19f..ff6952b 100644 --- a/.claude/settings.json +++ b/.claude/settings.json @@ -4,7 +4,10 @@ "PowerShell(cmd /c \"type C:\\\\Users\\\\bmant\\\\anaconda3\\\\envs\\\\astrolab\\\\qt6.conf\")", "PowerShell(cmd /c \"C:\\\\Users\\\\bmant\\\\anaconda3\\\\Scripts\\\\conda.exe run -n astrolab python -c \"\"import struct, sys; data = open\\(r'C:\\\\Users\\\\bmant\\\\anaconda3\\\\envs\\\\astrolab\\\\Lib\\\\site-packages\\\\PyQt6\\\\Qt6\\\\bin\\\\Qt6Core.dll', 'rb'\\).read\\(\\); print\\('DLL size:', len\\(data\\)\\)\"\" 2>&1\")", "PowerShell(cmd /c \"C:\\\\Users\\\\bmant\\\\anaconda3\\\\Scripts\\\\conda.exe run -n astrolab pip install pyqt6==6.11.0 pyqt6-qt6==6.11.0 pyqt6-sip==13.11.1 --force-reinstall --ignore-installed 2>&1\")", - "WebSearch" + "WebSearch", + "PowerShell(conda activate astrolab)", + "PowerShell(python -m pytest tests/ -m \"not slow\" -q 2>&1)", + "PowerShell(Get-ChildItem -Path \"$env:USERPROFILE\" -Filter \"condabin\" -Directory -ErrorAction SilentlyContinue -Recurse -Depth 2 | Select-Object -First 5 FullName; Get-ChildItem -Path \"C:\\\\ProgramData\" -Filter \"condabin\" -Directory -ErrorAction SilentlyContinue | Select-Object FullName)" ] } } diff --git a/AstroImageLab.py b/AstroImageLab.py index a5117b8..85ad46e 100644 --- a/AstroImageLab.py +++ b/AstroImageLab.py @@ -4,8 +4,8 @@ # # PR → merge to main (CI runs tests + build to verify everything works) # Tag the merge commit on main → triggers the release workflow -# git tag v0.0.7 -# git push origin v0.0.7 +# git tag v0.0.8 +# git push origin v0.0.8 import sys import os