Newer
Older
SmTIAS-LightSim / scripts / build_reflectance_panel.py
"""反射率スイープの比較パネル生成.

scripts/sweep_reflectance.py の出力から以下を生成:
    - 1×5 輝度マップパネル (shared vmin/vmax + 独立スケール 2 段比較)
    - 既存 4 指標の収束カーブ
    - 新 5 指標 (shape_metrics) の変化カーブ
    - 数値表

実行方法:
    uv run python scripts/build_reflectance_panel.py

出力:
    - output/_slide_panels/reflectance_sweep.png         (1×5 輝度マップ 2 段)
    - output/_slide_panels/reflectance_convergence.png   (既存 4 指標カーブ)
    - output/_slide_panels/reflectance_shape_metrics.png (新指標カーブ)
"""

from __future__ import annotations

import csv
import sys
from pathlib import Path

import matplotlib
import numpy as np

matplotlib.use("Agg")
import matplotlib.font_manager as fm
import matplotlib.pyplot as plt

from sm_tias_lightsim.evaluation.roi_presets import SIM_ROI, extract_roi_xywh

ROOT = Path(__file__).resolve().parents[1]
SWEEP_BASE = ROOT / "output" / "reflectance_sweep"
SUMMARY_CSV = SWEEP_BASE / "summary.csv"
OUT_LUM_PANEL = ROOT / "output" / "_slide_panels" / "reflectance_sweep.png"
OUT_CONV_PANEL = ROOT / "output" / "_slide_panels" / "reflectance_convergence.png"
OUT_SHAPE_PANEL = ROOT / "output" / "_slide_panels" / "reflectance_shape_metrics.png"

TARGET = {"cov": 0.0695, "gradient": 8.35, "cp": 1.091, "rmm": 0.860}


def _setup_japanese_font() -> None:
    for p in ("C:/Windows/Fonts/meiryo.ttc", "C:/Windows/Fonts/YuGothM.ttc"):
        if Path(p).exists():
            fm.fontManager.addfont(p)
            plt.rcParams["font.family"] = fm.FontProperties(fname=p).get_name()
            return


def main() -> int:
    _setup_japanese_font()

    if not SUMMARY_CSV.exists():
        print(f"エラー: {SUMMARY_CSV} が見つかりません", file=sys.stderr)
        return 1

    with SUMMARY_CSV.open(encoding="utf-8", newline="") as f:
        rows = list(csv.DictReader(f))

    reflectances = [float(r["reflectance"]) for r in rows]
    covs = [float(r["cov"]) for r in rows]
    grads = [float(r["gradient_magnitude"]) for r in rows]
    cps = [float(r["center_periphery_ratio"]) for r in rows]
    rmms = [float(r["radial_min_max_ratio"]) for r in rows]
    means = [float(r["mean"]) for r in rows]
    # 新指標
    centroid_y = [float(r["intensity_centroid_y_norm"]) for r in rows]
    ang_ratios = [float(r["angular_ratio"]) for r in rows]
    eccs = [float(r["pca_eccentricity"]) for r in rows]

    # npy 読み込み
    lums = []
    for r in rows:
        npy = SWEEP_BASE / r["job"] / "luminance.npy"
        if not npy.exists():
            lums.append(None)
            continue
        arr = extract_roi_xywh(np.load(npy), SIM_ROI)
        lums.append(arr.mean(axis=-1) if arr.ndim == 3 else arr)

    valid = [x for x in lums if x is not None]
    stacked = np.stack(valid)
    vmin_shared = float(np.percentile(stacked, 1))
    vmax_shared = float(np.percentile(stacked, 99))

    # ============ 輝度マップ 2 段比較 (shared / independent) ============
    n = len(rows)
    fig, axes = plt.subplots(2, n, figsize=(4.5 * n, 10))
    for j, (r, lum, cov, grad) in enumerate(zip(reflectances, lums, covs, grads, strict=False)):
        # 上段: shared
        ax_t = axes[0, j]
        im_t = ax_t.imshow(lum, cmap="hot", vmin=vmin_shared, vmax=vmax_shared)
        ax_t.set_title(f"R = {r}\nCoV={cov:.4f}  勾配={grad:.2f}%", fontsize=11)
        ax_t.set_xticks([])
        ax_t.set_yticks([])
        plt.colorbar(im_t, ax=ax_t, fraction=0.046, pad=0.04)
        # 下段: independent
        ax_b = axes[1, j]
        vmin_i = float(np.percentile(lum, 1))
        vmax_i = float(np.percentile(lum, 99))
        im_b = ax_b.imshow(lum, cmap="hot", vmin=vmin_i, vmax=vmax_i)
        ax_b.set_title(f"独立スケール [{vmin_i:.0f}, {vmax_i:.0f}]", fontsize=10)
        ax_b.set_xticks([])
        ax_b.set_yticks([])
        plt.colorbar(im_b, ax=ax_b, fraction=0.046, pad=0.04)

    axes[0, 0].set_ylabel("① 共有スケール", fontsize=12, fontweight="bold")
    axes[1, 0].set_ylabel("② 独立スケール\n(形比較用)", fontsize=12, fontweight="bold")
    fig.suptitle("反射率スイープ (d=16) — 共有 vs 独立スケール 2 段比較",
                 fontsize=14, fontweight="bold")
    plt.tight_layout()
    OUT_LUM_PANEL.parent.mkdir(parents=True, exist_ok=True)
    plt.savefig(OUT_LUM_PANEL, dpi=110, bbox_inches="tight")
    print(f"輝度マップパネル保存: {OUT_LUM_PANEL}")

    # ============ 既存 4 指標の収束カーブ ============
    fig, axes = plt.subplots(2, 2, figsize=(13, 9))
    metrics_existing = [
        ("CoV", covs, TARGET["cov"], "CoV (低いほど均一)"),
        ("勾配 [%]", grads, TARGET["gradient"], "勾配 [%]"),
        ("C/P 比", cps, TARGET["cp"], "C/P 比"),
        ("動径 min/max", rmms, TARGET["rmm"], "動径 min/max"),
    ]
    for ax, (label, vals, target, ylabel) in zip(axes.flat, metrics_existing, strict=False):
        ax.plot(
            reflectances, vals, marker="o", linewidth=2, markersize=10, color="C0", label="シミュ"
        )
        ax.axhline(
            target, color="red", linestyle="--", linewidth=1.5, label=f"実機 TARGET ({target})"
        )
        for r, v in zip(reflectances, vals, strict=False):
            ax.annotate(f"{v:.4f}", (r, v), textcoords="offset points", xytext=(0, 10),
                        ha="center", fontsize=9)
        ax.set_xlabel("反射率 R")
        ax.set_ylabel(ylabel)
        ax.set_title(label)
        ax.grid(alpha=0.3)
        ax.legend(loc="best", fontsize=9)

    fig.suptitle("反射率スイープ — 既存 4 指標の変化 (赤破線=実機 TARGET)",
                 fontsize=14, fontweight="bold")
    plt.tight_layout()
    plt.savefig(OUT_CONV_PANEL, dpi=110, bbox_inches="tight")
    print(f"既存指標カーブ保存: {OUT_CONV_PANEL}")

    # ============ 新指標カーブ ============
    fig, axes = plt.subplots(1, 3, figsize=(16, 5))
    new_metrics = [
        ("輝度重心 Y (norm)", centroid_y, 0.0, "centroid_y_norm"),
        ("angular_ratio (0=同心円, 1=角度依存)", ang_ratios, 0.0, "angular_ratio"),
        ("PCA eccentricity", eccs, 0.0, "eccentricity"),
    ]
    for ax, (label, vals, ideal, ylabel) in zip(axes, new_metrics, strict=False):
        ax.plot(
            reflectances, vals, marker="o", linewidth=2, markersize=10, color="C2", label="シミュ"
        )
        ax.axhline(ideal, color="gray", linestyle=":", linewidth=1.5, label=f"理想 ({ideal})")
        for r, v in zip(reflectances, vals, strict=False):
            ax.annotate(f"{v:.4f}", (r, v), textcoords="offset points", xytext=(0, 10),
                        ha="center", fontsize=9)
        ax.set_xlabel("反射率 R")
        ax.set_ylabel(ylabel)
        ax.set_title(label)
        ax.grid(alpha=0.3)
        ax.legend(loc="best", fontsize=9)

    fig.suptitle("反射率スイープ — 新指標 (shape_metrics) の変化",
                 fontsize=14, fontweight="bold")
    plt.tight_layout()
    plt.savefig(OUT_SHAPE_PANEL, dpi=110, bbox_inches="tight")
    print(f"新指標カーブ保存: {OUT_SHAPE_PANEL}")

    # ============ サマリ ============
    print()
    print("=== サマリ ===")
    h = ("R", "mean", "CoV", "勾配%", "C/P", "rmin/max", "centroid_y", "ang_ratio", "ecc")
    print(
        f"  {h[0]:<6}{h[1]:<10}{h[2]:<10}{h[3]:<10}{h[4]:<10}"
        f"{h[5]:<10}{h[6]:<12}{h[7]:<12}{h[8]:<10}"
    )
    for r, mean, cov, grad, cp, rmm, cy, ar, ec in zip(
            reflectances, means, covs, grads, cps, rmms, centroid_y, ang_ratios, eccs,
            strict=False):
        print(f"  {r:<6}{mean:<10.2f}{cov:<10.4f}{grad:<10.2f}{cp:<10.4f}{rmm:<10.4f}"
              f"{cy:<12.4f}{ar:<12.4f}{ec:<10.4f}")
    print(f"  TARGET: -          {TARGET['cov']:<10}{TARGET['gradient']:<10}"
          f"{TARGET['cp']:<10}{TARGET['rmm']:<10}")
    return 0


if __name__ == "__main__":
    sys.exit(main())