Newer
Older
SmTIAS-LightSim / tests / evaluation / test_oetf.py
"""evaluation/oetf.py の仕様ベーステスト.

仕様出典:
- docs/04_SPEC/SPEC_03_基準設計の物理構成.md(「手法訂正」「luma vs luminance」節)
- output/sweep_reflectance/apply_oetf.py(検証済み参照実装)
- docs/01_GUIDE/GUIDE_08_テスト方針.md(テスト方針)
- IEC 61966-2-1(sRGB 規格)

テスト対象:
  srgb_oetf(c): IEC 61966-2-1 OETF の順方向変換
  solve_exposure(lin_roi, target_mean): OETF 後 ROI mean を target に収束させる露出ソルバ
  apply_srgb_oetf(lin, scale, quantize_uint8): スケール適用 + OETF + 量子化

テスト方針:
  - 既知の入力・出力で数値精度を検証する(実装詳細でなく振る舞いを検証)
  - ファイル I/O はモックせず実 ndarray を使う (GUIDE_08)
  - 乱数を使う場合は seed を必ず固定する (GUIDE_08)
"""

from __future__ import annotations

import numpy
import numpy.testing
import pytest

from sm_tias_lightsim.evaluation.oetf import apply_srgb_oetf, solve_exposure, srgb_oetf


# ---------------------------------------------------------------------------
# 1. srgb_oetf — IEC 61966-2-1 の既知点
# ---------------------------------------------------------------------------


def test_srgb_oetf_zero_input_returns_zero():
    """入力 0.0 → 出力 0.0 (OETF(0) = 0 は IEC 61966-2-1 の自明な既知点).

    仕様: IEC 61966-2-1: 線形部 c <= 0.0031308 → 12.92 * c
    c=0 では 12.92 * 0 = 0.0
    """
    result = srgb_oetf(numpy.array([0.0]))
    assert result[0] == pytest.approx(0.0, abs=1e-12)


def test_srgb_oetf_one_input_returns_one():
    """入力 1.0 → 出力 1.0 (OETF(1) = 1 は IEC 61966-2-1 の既知点).

    仕様: IEC 61966-2-1: 非線形部 c > 0.0031308 → 1.055 * c^(1/2.4) - 0.055
    c=1 では 1.055 * 1^(1/2.4) - 0.055 = 1.055 - 0.055 = 1.0
    """
    result = srgb_oetf(numpy.array([1.0]))
    assert result[0] == pytest.approx(1.0, abs=1e-12)


def test_srgb_oetf_midpoint_exact():
    """中間点 0.5 の出力が標準式(非線形部)と一致すること.

    仕様: c=0.5 > 0.0031308 → 1.055 * 0.5^(1/2.4) - 0.055 (IEC 61966-2-1 非線形部)
    """
    c = 0.5
    expected = 1.055 * (c ** (1.0 / 2.4)) - 0.055
    result = srgb_oetf(numpy.array([c]))
    assert result[0] == pytest.approx(expected, rel=1e-9)


def test_srgb_oetf_linear_segment_exact():
    """線形部(c ≤ 0.0031308)の出力が 12.92 * c であること.

    仕様: IEC 61966-2-1: c <= 0.0031308 → 12.92 * c(線形トイ部)
    """
    c_vals = numpy.array([0.0, 0.001, 0.0031308])
    result = srgb_oetf(c_vals)
    expected = 12.92 * c_vals
    numpy.testing.assert_allclose(result, expected, rtol=1e-9)


def test_srgb_oetf_nonlinear_segment_at_boundary():
    """非線形部(c > 0.0031308)の出力が 1.055 * c^(1/2.4) - 0.055 であること.

    仕様: IEC 61966-2-1: c > 0.0031308 → 1.055 * c^(1/2.4) - 0.055
    """
    c_vals = numpy.array([0.003131, 0.01, 0.1, 0.5, 0.9])
    result = srgb_oetf(c_vals)
    expected = 1.055 * numpy.power(c_vals, 1.0 / 2.4) - 0.055
    numpy.testing.assert_allclose(result, expected, rtol=1e-9)


def test_srgb_oetf_continuity_at_threshold():
    """境界点 0.0031308 の前後で出力が連続(差が小さい)こと.

    仕様: OETF は [0,1] で連続な単調関数(IEC 61966-2-1)
    線形部と非線形部の接続点での不連続(ジャンプ)がないことを確認する.
    """
    threshold = 0.0031308
    eps = 1e-7
    lo = srgb_oetf(numpy.array([threshold - eps]))[0]
    hi = srgb_oetf(numpy.array([threshold + eps]))[0]
    # 連続性:差が ε の 100 倍程度以下であること
    assert abs(hi - lo) < 1e-4


def test_srgb_oetf_monotonically_increasing():
    """OETF が [0, 1] で単調増加であること.

    仕様: OETF は単調増加 → solve_exposure の二分法が成立する根拠(SPEC_03)
    """
    c_vals = numpy.linspace(0.0, 1.0, 1000)
    result = srgb_oetf(c_vals)
    diffs = numpy.diff(result)
    assert numpy.all(diffs >= 0.0)


def test_srgb_oetf_clips_negative_input():
    """入力が負のとき 0.0 にクリップされること.

    仕様: 入力は [0, 1] にクリップしてから変換する(oetf.py docstring)
    """
    result = srgb_oetf(numpy.array([-1.0, -0.5, -0.001]))
    numpy.testing.assert_allclose(result, 0.0, atol=1e-12)


def test_srgb_oetf_clips_above_one_input():
    """入力が 1.0 を超えるとき 1.0 にクリップされること.

    仕様: 入力は [0, 1] にクリップしてから変換する(oetf.py docstring)
    """
    result = srgb_oetf(numpy.array([1.001, 2.0, 10.0]))
    numpy.testing.assert_allclose(result, 1.0, atol=1e-12)


def test_srgb_oetf_output_range_is_zero_to_one():
    """任意の入力に対して出力が [0, 1] の範囲に収まること.

    仕様: sRGB OETF: 線形 [0,1] → 符号化 [0,1] へ変換する(oetf.py docstring)
    """
    rng = numpy.random.RandomState(42)
    # [-1, 2] の広い範囲の入力
    c_vals = rng.uniform(-1.0, 2.0, 1000)
    result = srgb_oetf(c_vals)
    assert numpy.all(result >= 0.0)
    assert numpy.all(result <= 1.0)


def test_srgb_oetf_output_dtype_is_float():
    """出力 dtype が float 型(float32 または float64)であること.

    仕様: Returns: sRGB 符号化値(float64, 値域 [0, 1])(oetf.py docstring)
    """
    result = srgb_oetf(numpy.array([0.5], dtype=numpy.float64))
    assert numpy.issubdtype(result.dtype, numpy.floating)


def test_srgb_oetf_matches_reference_implementation():
    """参照実装(apply_oetf.py の srgb_oetf)と数値が一致すること.

    仕様: output/sweep_reflectance/apply_oetf.py は検証済み参照実装(SPEC_03)
    参照実装の式:
      lo (c <= 0.0031308): 12.92 * c
      hi (c >  0.0031308): 1.055 * c^(1/2.4) - 0.055
    """
    # 参照実装と同じ計算
    def ref_srgb_oetf(c: numpy.ndarray) -> numpy.ndarray:
        c = numpy.clip(c, 0.0, 1.0)
        lo = c <= 0.0031308
        out = numpy.empty_like(c)
        out[lo] = 12.92 * c[lo]
        out[~lo] = 1.055 * numpy.power(c[~lo], 1.0 / 2.4) - 0.055
        return out

    rng = numpy.random.RandomState(0)
    c_vals = rng.uniform(0.0, 1.0, 500)
    expected = ref_srgb_oetf(c_vals)
    result = srgb_oetf(c_vals)
    numpy.testing.assert_allclose(result, expected, rtol=1e-12)


# ---------------------------------------------------------------------------
# 2. solve_exposure — 露出ソルバの収束と仕様
# ---------------------------------------------------------------------------


def test_solve_exposure_roi_mean_converges_to_target():
    """solve_exposure 後に OETF 済み ROI mean が target_mean に収束すること.

    仕様: OETF 後の ROI mean が target_mean に収束する(許容誤差内)(SPEC_03)
    OETF は単調増加 → mean(scale) も単調 → 二分法で一意解(SPEC_03)
    """
    rng = numpy.random.RandomState(1)
    # 均一に近い合成線形輝度(値域 [0.1, 0.3])
    lin_roi = rng.uniform(0.1, 0.3, (100, 100)).astype(numpy.float64)
    target_mean = 146.7146

    scale, clipped = solve_exposure(lin_roi, target_mean)

    # 得られた scale で実際にエンコードして mean を確認する
    encoded_mean = (srgb_oetf(lin_roi * scale) * 255.0).mean()
    assert encoded_mean == pytest.approx(target_mean, abs=0.01)


def test_solve_exposure_returns_positive_scale():
    """solve_exposure が正の scale を返すこと.

    仕様: scale は線形ゲイン(正の実数でなければ意味がない)
    """
    rng = numpy.random.RandomState(2)
    lin_roi = rng.uniform(0.05, 0.2, (50, 50)).astype(numpy.float64)
    scale, _ = solve_exposure(lin_roi, 120.0)
    assert scale > 0.0


def test_solve_exposure_returns_tuple():
    """solve_exposure が (scale, clipped_fraction) のタプルを返すこと.

    仕様: Returns: (scale, clipped_fraction) のタプル(oetf.py docstring)
    """
    lin_roi = numpy.full((10, 10), 0.1, dtype=numpy.float64)
    result = solve_exposure(lin_roi, 100.0)
    assert isinstance(result, tuple)
    assert len(result) == 2


def test_solve_exposure_clipped_fraction_is_zero_for_small_roi():
    """白飛びしない(clip が 0 付近の)線形入力で clipped_fraction ≈ 0 になること.

    仕様: clipped_fraction: OETF 入力が 1.0 を超えるピクセルの割合(oetf.py docstring)
    入力が小さい場合、スケールを掛けても 1.0 を超えないためクリップ率は 0.
    """
    # 非常に小さい linear 値(スケール後も確実に 1.0 未満)
    lin_roi = numpy.full((20, 20), 1e-5, dtype=numpy.float64)
    _, clipped = solve_exposure(lin_roi, 10.0)
    assert clipped == pytest.approx(0.0, abs=1e-6)


def test_solve_exposure_clipped_fraction_is_positive_for_large_roi():
    """大きな linear 値でターゲット mean が高い場合、clipped_fraction > 0 になること.

    仕様: clipped_fraction: OETF 入力が 1.0 を超えるピクセルの割合(oetf.py docstring)
    """
    # 均一な輝度でターゲット mean をほぼ 255 に設定する
    # → scale が大きくなり、一部ピクセルが 1.0 を超える
    lin_roi = numpy.full((50, 50), 0.8, dtype=numpy.float64)
    _, clipped = solve_exposure(lin_roi, 254.0)
    # scale が大きいので clip が発生するはず
    assert clipped >= 0.0  # 0 または正


def test_solve_exposure_clipped_fraction_range():
    """clipped_fraction が [0, 1] の範囲内に収まること.

    仕様: OETF 入力が 1.0 を超えるピクセルの割合 [0, 1](oetf.py docstring)
    """
    rng = numpy.random.RandomState(3)
    lin_roi = rng.uniform(0.01, 1.0, (30, 30)).astype(numpy.float64)
    _, clipped = solve_exposure(lin_roi, 180.0)
    assert 0.0 <= clipped <= 1.0


def test_solve_exposure_different_targets_give_different_scales():
    """異なる target_mean に対して異なる scale が返ること(単調性).

    仕様: OETF は単調増加 → mean(scale) も単調 → target が大きいほど scale が大きい
    """
    lin_roi = numpy.full((20, 20), 0.1, dtype=numpy.float64)
    scale_low, _ = solve_exposure(lin_roi, 100.0)
    scale_high, _ = solve_exposure(lin_roi, 180.0)
    assert scale_high > scale_low


def test_solve_exposure_uniform_roi_converges():
    """均一な ROI(全画素同値)でも収束すること.

    仕様: 均一入力は最もシンプルなケース;収束が成功すること
    """
    lin_roi = numpy.full((64, 64), 0.05, dtype=numpy.float64)
    target = 146.7146
    scale, clipped = solve_exposure(lin_roi, target)
    encoded_mean = (srgb_oetf(lin_roi * scale) * 255.0).mean()
    assert encoded_mean == pytest.approx(target, abs=0.1)


# ---------------------------------------------------------------------------
# 3. apply_srgb_oetf — スケール適用 + OETF + 量子化
# ---------------------------------------------------------------------------


def test_apply_srgb_oetf_returns_uint8_range_when_quantized():
    """quantize_uint8=True のとき出力が [0, 255] の整数値(float64)であること.

    仕様: quantize_uint8=True のとき np.clip(np.round(x * 255), 0, 255)(oetf.py docstring)
    """
    lin = numpy.linspace(0.0, 1.0, 100)
    result = apply_srgb_oetf(lin, scale=1.0, quantize_uint8=True)
    assert numpy.all(result >= 0.0)
    assert numpy.all(result <= 255.0)
    # 整数値かどうか確認
    numpy.testing.assert_allclose(result, numpy.round(result), atol=1e-10)


def test_apply_srgb_oetf_zero_input_returns_zero_with_quantize():
    """入力 0 → 出力 0(uint8 量子化時).

    仕様: OETF(0) = 0 → round(0 * 255) = 0
    """
    result = apply_srgb_oetf(numpy.array([0.0]), scale=1.0, quantize_uint8=True)
    assert result[0] == pytest.approx(0.0, abs=1e-10)


def test_apply_srgb_oetf_one_input_returns_255_with_quantize():
    """入力 1.0 → 出力 255(uint8 量子化時).

    仕様: OETF(1) = 1 → round(1 * 255) = 255
    """
    result = apply_srgb_oetf(numpy.array([1.0]), scale=1.0, quantize_uint8=True)
    assert result[0] == pytest.approx(255.0, abs=1e-10)


def test_apply_srgb_oetf_no_quantize_returns_continuous():
    """quantize_uint8=False のとき連続値(量子化なし)が返ること.

    仕様: quantize_uint8=True のとき uint8 量子化(oetf.py docstring)
    → False のときは量子化しないので連続値になる
    """
    lin = numpy.array([0.1, 0.2, 0.3, 0.5, 0.8])
    result = apply_srgb_oetf(lin, scale=1.0, quantize_uint8=False)
    # 連続値のため整数ちょうどにはならないケースがほとんど
    # OETF 結果 * 255 の少数部が 0 でないものが存在すること
    fractional_parts = numpy.abs(result - numpy.round(result))
    assert numpy.any(fractional_parts > 1e-3)


def test_apply_srgb_oetf_scale_amplifies_input():
    """scale > 1 のとき出力が scale=1 より大きくなること.

    仕様: apply_srgb_oetf は lin * scale に OETF を掛ける(OETF は単調増加)
    scale > 1 → 入力が大きい → OETF 出力が大きい
    """
    lin = numpy.array([0.1, 0.2, 0.3])
    result_s1 = apply_srgb_oetf(lin, scale=1.0, quantize_uint8=False)
    result_s2 = apply_srgb_oetf(lin, scale=2.0, quantize_uint8=False)
    assert numpy.all(result_s2 >= result_s1)


def test_apply_srgb_oetf_matches_manual_computation():
    """apply_srgb_oetf の結果が手動計算(OETF → *255 → round → clip)と一致すること.

    仕様: 内部演算は float64 で行い,quantize_uint8=True のとき
          np.clip(np.round(x * 255), 0, 255)(oetf.py docstring)
    """
    lin = numpy.array([0.0, 0.01, 0.1, 0.5, 0.9, 1.0], dtype=numpy.float64)
    scale = 1.5
    result = apply_srgb_oetf(lin, scale=scale, quantize_uint8=True)

    # 手動計算
    expected = numpy.clip(numpy.round(srgb_oetf(lin * scale) * 255.0), 0.0, 255.0)
    numpy.testing.assert_allclose(result, expected, atol=1e-10)


def test_apply_srgb_oetf_2d_input_shape_preserved():
    """2D 入力に対して形状が保持されること.

    仕様: lin: linear 輝度配列(任意次元 float64)(oetf.py docstring)
    """
    lin = numpy.full((32, 32), 0.2, dtype=numpy.float64)
    result = apply_srgb_oetf(lin, scale=1.0, quantize_uint8=True)
    assert result.shape == (32, 32)


def test_apply_srgb_oetf_default_quantize_is_true():
    """quantize_uint8 のデフォルトが True であること(引数省略時).

    仕様: quantize_uint8: bool = True(oetf.py apply_srgb_oetf デフォルト引数)
    デフォルトで量子化されるため、出力は整数値になる
    """
    lin = numpy.array([0.5], dtype=numpy.float64)
    result = apply_srgb_oetf(lin, scale=1.0)
    fractional_part = abs(result[0] - round(result[0]))
    assert fractional_part < 1e-10


# ---------------------------------------------------------------------------
# 4. 参照実装との数値一致(R0.85 sRGB 行 apply_oetf.py 記載値)
# ---------------------------------------------------------------------------


def test_srgb_oetf_at_known_luminance_value():
    """SPEC_03 の R0.85 sRGB mean≈146.71 の元となる OETF 計算が正しく機能すること.

    仕様: OETF 適用後の比較(SPEC_03「手法訂正」節)
      R0.85 (d12,spp256) sRGB: mean≈146.71, clip 0%
    テスト: 正規化輝度を OETF 変換して *255 した値の範囲が物理的に妥当であること

    参照実装(apply_oetf.py)の encode_with_exposure は
    srgb_oetf(lin * scale) * 255.0 で計算する
    """
    # sRGB OETF の典型的な動作点の検証
    # 例: 線形 0.3 の OETF 値が約 0.597 であることを確認(標準式の計算値)
    c = 0.3
    expected = 1.055 * (c ** (1.0 / 2.4)) - 0.055  # ≈ 0.5969
    result = float(srgb_oetf(numpy.array([c]))[0])
    assert result == pytest.approx(expected, rel=1e-9)
    # この値に 255 を掛けると ≈ 152.2 程度(正常な sRGB 空間の値)
    assert 140.0 < result * 255.0 < 165.0


def test_solve_exposure_with_r085_like_luminance():
    """R0.85 データ相当の合成輝度でソルバが target_mean=146.7146 に収束すること.

    仕様: SPEC_03「R0.85 (d12,spp256) sRGB mean≈146.71, clip 0%」
    R0.85 の線形空間 mean ≈ 0.1~0.3 程度を仮定した合成データで検証する.
    (実際のデータを使わず合成で近似する.実データはレンダが必要なため slow テスト)
    """
    rng = numpy.random.RandomState(42)
    # R0.85 データの linear mean が約 0.15 程度であることを想定
    lin_roi = rng.uniform(0.05, 0.25, (100, 100)).astype(numpy.float64)
    target_mean = 146.7146

    scale, clipped = solve_exposure(lin_roi, target_mean)
    encoded_mean = (srgb_oetf(lin_roi * scale) * 255.0).mean()

    # target_mean への収束精度(200 回対数二分なので高精度)
    assert encoded_mean == pytest.approx(target_mean, abs=0.01)
    # clip 率は 0% のはず(R0.85 データの既知値と整合)
    assert clipped >= 0.0
    assert clipped <= 1.0


# ---------------------------------------------------------------------------
# 5. 無彩色シーン: luma vs luminance の等価性(SPEC_03「luma vs luminance」節)
# ---------------------------------------------------------------------------


def test_grayscale_oetf_channel_vs_luminance_equivalence():
    """無彩色(R=G=B)入力で「チャネル別 OETF → Rec.709 加重和」と「luminance に直接 OETF」が厳密一致すること.

    仕様: SPEC_03「luma vs luminance と演算順序」節
      無彩色画素(R=G=B=x)では:
        Σwᵢ·f(x) = f(x·Σwᵢ) = f(x)(Rec.709 重みの和は 1.0)
      →「チャネル別 OETF→加重和」=「加重和→OETF」= 「luminance に直接 OETF」
    """
    rng = numpy.random.RandomState(99)
    # 無彩色: R=G=B(ランダム値)
    x = rng.uniform(0.0, 1.0, (50, 50)).astype(numpy.float64)
    rgb = numpy.stack([x, x, x], axis=-1)  # shape (50, 50, 3)

    # Rec.709 係数
    w = numpy.array([0.2126, 0.7152, 0.0722])
    assert abs(w.sum() - 1.0) < 1e-10  # 重みの和は 1

    # 方法1: チャネル別 OETF → Rec.709 加重和(実測の演算順)
    r_prime = srgb_oetf(rgb[..., 0])
    g_prime = srgb_oetf(rgb[..., 1])
    b_prime = srgb_oetf(rgb[..., 2])
    luma_from_channels = w[0] * r_prime + w[1] * g_prime + w[2] * b_prime

    # 方法2: luminance(= x)に直接 OETF(シミュの実装)
    luminance_direct = srgb_oetf(x)

    # 無彩色条件(R=G=B=x)では厳密一致する
    numpy.testing.assert_allclose(luma_from_channels, luminance_direct, rtol=1e-12)