// Copyright (C) 2026 Kiyotsugu Arai
// SPDX-License-Identifier: LGPL-3.0-or-later

/**
 * @file tonemap.hpp
 * @brief トーンマッピング (DRW-65) + サンプリング (DRW-66) + デノイザー (DRW-67) + AOV (DRW-68)
 */

#ifndef CALX_IMAGE_RENDER_TONEMAP_HPP
#define CALX_IMAGE_RENDER_TONEMAP_HPP

#include <image/Image.hpp>
#include <cmath>
#include <algorithm>
#include <vector>
#include <string>

namespace calx::image {

// ========================================================================
// DRW-65: トーンマッピング
// ========================================================================

/// Reinhard トーンマッピング (グローバル)
template<typename T>
Image<T> tonemapReinhard(const Image<T>& hdr, T exposure = T(1)) {
    Image<T> result(hdr.height(), hdr.width(), hdr.numChannels());
    for (std::size_t y = 0; y < hdr.height(); ++y)
        for (std::size_t x = 0; x < hdr.width(); ++x)
            for (int c = 0; c < hdr.numChannels(); ++c) {
                T v = hdr.at(y, x, c) * exposure;
                result.at(y, x, c) = v / (T(1) + v);
            }
    return result;
}

/// Filmic トーンマッピング (Uncharted 2)
template<typename T>
Image<T> tonemapFilmic(const Image<T>& hdr, T exposure = T(1)) {
    auto filmic = [](T x) -> T {
        T A=T(0.15), B=T(0.50), C=T(0.10), D=T(0.20), E=T(0.02), F=T(0.30);
        return ((x*(A*x+C*B)+D*E) / (x*(A*x+B)+D*F)) - E/F;
    };
    T W = T(11.2);
    T whiteScale = T(1) / filmic(W);
    Image<T> result(hdr.height(), hdr.width(), hdr.numChannels());
    for (std::size_t y = 0; y < hdr.height(); ++y)
        for (std::size_t x = 0; x < hdr.width(); ++x)
            for (int c = 0; c < hdr.numChannels(); ++c) {
                T v = hdr.at(y, x, c) * exposure;
                result.at(y, x, c) = std::clamp(filmic(v) * whiteScale, T(0), T(1));
            }
    return result;
}

/// ACES トーンマッピング (簡易フィット)
template<typename T>
Image<T> tonemapACES(const Image<T>& hdr, T exposure = T(1)) {
    Image<T> result(hdr.height(), hdr.width(), hdr.numChannels());
    for (std::size_t y = 0; y < hdr.height(); ++y)
        for (std::size_t x = 0; x < hdr.width(); ++x)
            for (int c = 0; c < hdr.numChannels(); ++c) {
                T v = hdr.at(y, x, c) * exposure;
                T a = T(2.51), b = T(0.03), cc = T(2.43), d = T(0.59), e = T(0.14);
                result.at(y, x, c) = std::clamp((v*(a*v+b)) / (v*(cc*v+d)+e), T(0), T(1));
            }
    return result;
}

// ========================================================================
// DRW-66: 適応的サンプリング
// ========================================================================

/// ピクセルごとの分散を計算し、追加サンプルが必要な領域を判定
template<typename T>
Image<T> computeVarianceMap(const std::vector<Image<T>>& samples) {
    if (samples.empty()) return {};
    int h = static_cast<int>(samples[0].height());
    int w = static_cast<int>(samples[0].width());
    int ch = samples[0].numChannels();
    int n = static_cast<int>(samples.size());

    Image<T> variance(h, w, 1, T(0));
    for (int y = 0; y < h; ++y)
        for (int x = 0; x < w; ++x) {
            T var = T(0);
            for (int c = 0; c < ch; ++c) {
                T mean = T(0);
                for (int s = 0; s < n; ++s) mean += samples[s].at(y, x, c);
                mean /= T(n);
                T v = T(0);
                for (int s = 0; s < n; ++s) {
                    T d = samples[s].at(y, x, c) - mean;
                    v += d * d;
                }
                var += v / T(n);
            }
            variance(y, x) = var / T(ch);
        }
    return variance;
}

// ========================================================================
// DRW-67: デノイザー (統計的)
// ========================================================================

/// 簡易デノイザー: 分散に基づくバイラテラルフィルタ
template<typename T>
Image<T> denoise(const Image<T>& noisy, const Image<T>& varianceMap,
                  int radius = 3, T strengthFactor = T(1)) {
    int h = static_cast<int>(noisy.height());
    int w = static_cast<int>(noisy.width());
    int ch = noisy.numChannels();
    Image<T> result(h, w, ch, T(0));

    for (int y = 0; y < h; ++y)
        for (int x = 0; x < w; ++x) {
            T localVar = varianceMap(y, x) * strengthFactor;
            T sigmaR2 = std::max(localVar, T(1e-6));

            for (int c = 0; c < ch; ++c) {
                T sumW = T(0), sumV = T(0);
                T center = noisy.at(y, x, c);
                for (int dy = -radius; dy <= radius; ++dy)
                    for (int dx = -radius; dx <= radius; ++dx) {
                        int ny = std::clamp(y+dy, 0, h-1);
                        int nx = std::clamp(x+dx, 0, w-1);
                        T v = noisy.at(ny, nx, c);
                        T spatialW = std::exp(T(-(dx*dx + dy*dy)) / T(2*radius*radius));
                        T rangeW = std::exp(-(v-center)*(v-center) / (T(2)*sigmaR2));
                        T w = spatialW * rangeW;
                        sumV += v * w;
                        sumW += w;
                    }
                result.at(y, x, c) = sumV / (sumW + T(1e-10));
            }
        }
    return result;
}

// ========================================================================
// DRW-68: AOV / レンダーパス出力
// ========================================================================

/// レンダーパスの種類
enum class RenderPass { Beauty, Diffuse, Specular, Normal, Depth, ObjectId, Albedo, Emission };

/// AOV バッファ (Arbitrary Output Variable)
template<typename T>
struct AOVBuffer {
    Image<T> beauty;    // 最終画像
    Image<T> diffuse;   // 拡散成分
    Image<T> specular;  // スペキュラ成分
    Image<T> normal;    // ワールド法線 (RGB = XYZ)
    Image<T> depth;     // 深度 (1ch)
    Image<T> objectId;  // オブジェクト ID (1ch)
    Image<T> albedo;    // アルベド (テクスチャ色)
    Image<T> emission;  // 発光成分

    void init(std::size_t h, std::size_t w) {
        beauty = Image<T>(h, w, 3, T(0));
        diffuse = Image<T>(h, w, 3, T(0));
        specular = Image<T>(h, w, 3, T(0));
        normal = Image<T>(h, w, 3, T(0));
        depth = Image<T>(h, w, 1, T(0));
        objectId = Image<T>(h, w, 1, T(0));
        albedo = Image<T>(h, w, 3, T(0));
        emission = Image<T>(h, w, 3, T(0));
    }
};

} // namespace calx::image

#endif
