"""
    NUVMAP Projection Engine — Julia Port

Mirrors `python/nuvmap/projection_engine.py` (SilverSight Q16_16 version).

Cross-validated against:
  * Python Q16_16 (SilverSight)
  * Python float (Research Stack)
  * R Q16_16 (SilverSight)

All three agree within Q16_16 quantization tolerance.

References:
  * `formal/BraidStateN.lean` — Rossby/Kelvin boundary correspondence
  * `julia/CoreFormalism/Q16_16.jl` — Q16_16 fixed-point library
"""
module NUVMAP

using SHA
using Dates
using ..Q16_16

export NUVMAPCell, NUVMAPSurface, NUVMAPProjectionEngine,
       build_nuvmap_from_eigenmass

# ── Schema ──────────────────────────────────────────────────────────

const SCHEMA_VERSION = 1

# ── Data structures ─────────────────────────────────────────────────

"""
    NUVMAPCell

A single NUVMAP address cell. All numeric fields are Q16_16 raw values.
Mirrors `python/nuvmap/projection_engine.NUVMAPCell`.
"""
struct NUVMAPCell
    u_i::Int32
    v_i::Int32
    k_i::Int32
    E_i::Q16_16Value
    R_i::Q16_16Value
    chi_i::Q16_16Value
    S_i::Q16_16Value
    L_i::Q16_16Value
    q_i::Int32
    admissible::Bool
    equation_id::Int32
    fingerprint::String
end

"""
    NUVMAPSurface

The full NUVMAP address surface. Mirrors `python/nuvmap/projection_engine.NUVMAPSurface`.
"""
struct NUVMAPSurface
    cells::Vector{NUVMAPCell}
    total_qubits::Int32
    bekenstein_bound::Q16_16Value
    area_utilization::Q16_16Value
    root_fingerprint::String
    timestamp::String
end

# ── Engine ──────────────────────────────────────────────────────────

"""
    NUVMAPProjectionEngine

Projects eigenmass data into a NUVMAP address surface using Q16_16.
Mirrors `python/nuvmap/projection_engine.NUVMAPProjectionEngine`.
"""
struct NUVMAPProjectionEngine
    total_qubit_budget::Int32
    chi_max_q16::Q16_16Value
    R_max_q16::Q16_16Value
    landauer_threshold_q16::Q16_16Value
end

function NUVMAPProjectionEngine(;
    total_qubit_budget::Integer = 0,
    chi_max_q16::Q16_16Value = Q16_16.half,
    R_max_q16::Q16_16Value = Q16_16.half,
    landauer_threshold_q16::Union{Nothing, Q16_16Value} = nothing
)
    lt = landauer_threshold_q16 === nothing ? div_op(Q16_16.one, of_raw_int(10)) : landauer_threshold_q16
    NUVMAPProjectionEngine(
        Int32(total_qubit_budget), chi_max_q16, R_max_q16, lt
    )
end

"""
    project(engine, eigenmass_data; eigenvalue_q16=nothing)

Project eigenmass data into a NUVMAP surface.

`eigenmass_data`: Vector of Dicts with keys:
  `:equation_id` (Int), `:amvr_q16`, `:avmr_q16`, `:chiral_residual_q16`
  (all Q16_16 raw), `:chiral_state` (String).

All numeric values MUST already be Q16_16 raw integers.
Convert at the outermost call boundary via `Q16_16.of_float`.
"""
function project(engine::NUVMAPProjectionEngine,
                 eigenmass_data::Vector{Dict{Symbol, Any}};
                 eigenvalue_q16::Union{Nothing, Q16_16Value} = nothing)
    isempty(eigenmass_data) && return NUVMAPSurface(
        NUVMAPCell[], Int32(0), Q16_16.zero, Q16_16.zero, "", ""
    )

    n = length(eigenmass_data)
    cells = NUVMAPCell[]

    # max eigenmass
    max_eigenmass = Q16_16.zero
    for d in eigenmass_data
        raw_q = q_div(q_add(d[:amvr_q16], d[:avmr_q16]), of_raw_int(Q16_16.to_int(Q16_16.one) * 2))
        if raw_q > max_eigenmass
            max_eigenmass = raw_q
        end
    end
    if max_eigenmass == Q16_16.zero
        max_eigenmass = Q16_16.one
    end

    for (i, d) in enumerate(eigenmass_data)
        amvr_q = d[:amvr_q16]
        avmr_q = d[:avmr_q16]
        cr_q   = d[:chiral_residual_q16]
        eq_id  = Int32(d[:equation_id])
        cs     = String(d[:chiral_state])

        raw_eigenmass = q_div(q_add(amvr_q, avmr_q), of_raw_int(Q16_16.to_int(Q16_16.one) * 2))
        E_norm = q_div(raw_eigenmass, max_eigenmass)

        one_minus = q_sub(Q16_16.one, E_norm)
        R_i = if one_minus > q16_pct1; one_minus else q16_pct1 end

        if cs == "chiral_scarred"
            R_i = q_mul(R_i, q16_150pct)
        end

        S_i = if cs == "achiral_stable"
            Q16_16.one
        elseif cs in ("left_handed_mass_bias", "right_handed_vector_bias")
            q16_pct70
        else
            q16_pct30
        end

        L_i = if E_norm > engine.landauer_threshold_q16
            Q16_16.one
        else
            q_div(E_norm, engine.landauer_threshold_q16)
        end

        lam_q = eigenvalue_q16 === nothing ? Q16_16.one : eigenvalue_q16
        v_abs = E_norm

        E_i = q_div(q_mul(q_mul(q_mul(lam_q, v_abs), S_i), L_i),
                    q_add(R_i, Q16_16.epsilon))

        chi_i = cr_q

        is_R_ok  = R_i ≤ engine.R_max_q16
        is_chi_ok = chi_i ≤ engine.chi_max_q16
        admissible = is_R_ok && is_chi_ok

        fp_payload = string(eq_id, '\0', Q16_16.to_int(amvr_q), '\0',
                           Q16_16.to_int(avmr_q), '\0', Q16_16.to_int(cr_q), '\0', cs)
        fp = bytes2hex(SHA.sha256(fp_payload))

        push!(cells, NUVMAPCell(
            Int32(i-1), Int32(i-1), Int32(i-1),
            E_i, R_i, chi_i, S_i, L_i,
            Int32(0), admissible, eq_id, fp
        ))
    end

    # Qubit allocation proportional to E_i / (R_i + eps)
    total_weight = Q16_16.zero
    for c in cells
        total_weight = q_add(total_weight, q_div(c.E_i, q_add(c.R_i, Q16_16.epsilon)))
    end
    if total_weight == Q16_16.zero
        total_weight = Q16_16.one
    end

    budget = if engine.total_qubit_budget > 0
        Int64(engine.total_qubit_budget)
    else
        b = Int64(0)
        for c in cells
            if c.admissible
                b += div(Int64(Q16_16.to_int(c.E_i)) * 100, Q16_16.q16_scale)
            end
        end
        b
    end
    n_admissible = count(c -> c.admissible, cells)
    budget = max(budget, Int64(n_admissible))

    for (i, c) in enumerate(cells)
        if c.admissible
            weight = q_div(c.E_i, q_add(c.R_i, Q16_16.epsilon))
            total_weight_int = Q16_16.to_int(total_weight)
            weight_int = Q16_16.to_int(weight)
            raw_q = div(budget * Int64(weight_int), Int64(total_weight_int))
            target_q = max(Int32(1), Int32(raw_q))
            cells[i] = NUVMAPCell(
                c.u_i, c.v_i, c.k_i, c.E_i, c.R_i, c.chi_i,
                c.S_i, c.L_i, target_q, c.admissible,
                c.equation_id, c.fingerprint
            )
        else
            cells[i] = NUVMAPCell(
                c.u_i, c.v_i, c.k_i, c.E_i, c.R_i, c.chi_i,
                c.S_i, c.L_i, Int32(0), c.admissible,
                c.equation_id, c.fingerprint
            )
        end
    end

    total_qubits = Int32(sum(c.q_i for c in cells))

    # Bekenstein-like bound
    bekenstein = if !isempty(cells)
        sum_E = Q16_16.zero
        for c in cells
            sum_E = q_add(sum_E, c.E_i)
        end
        q_div(sum_E, of_raw_int(length(cells) * Q16_16.to_int(Q16_16.one)))
    else
        Q16_16.zero
    end

    area_util = if bekenstein > Q16_16.zero
        q_div(of_raw_int(Q16_16.to_int(total_qubits) * Q16_16.to_int(Q16_16.one)),
              q_add(bekenstein, Q16_16.epsilon))
    else
        Q16_16.zero
    end

    root_fp = _compute_surface_root(cells)
    timestamp = string(Dates.now(Dates.UTC))

    NUVMAPSurface(cells, total_qubits, bekenstein, area_util, root_fp, timestamp)
end

function _compute_surface_root(cells::Vector{NUVMAPCell})
    payload = join(map(c -> string(c.u_i, ':', Q16_16.to_int(c.E_i), ':',
                                   Q16_16.to_int(c.chi_i), ':', c.q_i),
                       sort(cells, by = c -> c.u_i)), "|")
    bytes2hex(SHA.sha256(payload))
end

"""
    quantum_storage_admissible(engine, surface, node_i, tau_q16)

Implements the Lean-safe gate from BraidStateN.lean.
"""
function quantum_storage_admissible(engine::NUVMAPProjectionEngine,
                                     surface::NUVMAPSurface,
                                     node_i::Integer, tau_q16::Q16_16Value)
    if node_i < 1 || node_i > length(surface.cells)
        return false
    end
    c = surface.cells[node_i]
    lhs = q_mul(c.E_i, q_add(c.R_i, Q16_16.epsilon))
    rhs = q_mul(tau_q16, q_add(c.R_i, Q16_16.epsilon))
    return lhs ≤ rhs && c.chi_i ≤ engine.chi_max_q16 && c.admissible
end

"""
    get_density_map(surface)

Return the eigenmass density distribution.
"""
function get_density_map(surface::NUVMAPSurface)
    isempty(surface.cells) && return Dict(
        :E_i => Int32[], :q_i => Int32[], :chi_i => Q16_16Value[],
        :R_i => Q16_16Value[], :equation_ids => Int32[]
    )
    return Dict(
        :E_i => [c.E_i for c in surface.cells],
        :q_i => [c.q_i for c in surface.cells],
        :chi_i => [c.chi_i for c in surface.cells],
        :R_i => [c.R_i for c in surface.cells],
        :equation_ids => [c.equation_id for c in surface.cells],
    )
end

"""
    summary(surface)
"""
function summary(surface::NUVMAPSurface)
    admiss = filter(c -> c.admissible, surface.cells)
    return Dict{String, Any}(
        "num_cells" => Int32(length(surface.cells)),
        "num_admissible" => Int32(length(admiss)),
        "num_rejected" => Int32(length(surface.cells) - length(admiss)),
        "total_qubits" => surface.total_qubits,
        "avg_qubits_per_cell" => div(surface.total_qubits, max(Int32(1), Int32(length(admiss)))),
        "bekenstein_bound_q16" => Q16_16.to_int(surface.bekenstein_bound),
        "area_utilization_q16" => Q16_16.to_int(surface.area_utilization),
        "max_eigenmass_q16" => isempty(surface.cells) ? 0 : Q16_16.to_int(maximum(c -> c.E_i, surface.cells)),
        "max_chiral_q16" => isempty(surface.cells) ? 0 : Q16_16.to_int(maximum(c -> c.chi_i, surface.cells)),
        "surface_root" => surface.root_fingerprint,
    )
end

# ── Convenience ─────────────────────────────────────────────────────

function build_nuvmap_from_eigenmass(eigenmass_data::Vector{Dict{Symbol, Any}};
                                     qubit_budget::Integer = 0)
    engine = NUVMAPProjectionEngine(total_qubit_budget = qubit_budget)
    return project(engine, eigenmass_data)
end

# ── Internal aliases ────────────────────────────────────────────────

const q16_pct1   = Q16_16.q16_pct1
const q16_pct70  = Q16_16.q16_pct70
const q16_pct30  = Q16_16.q16_pct30
const q16_150pct = Q16_16.q16_150pct
const of_raw_int = Q16_16.of_raw_int
const q_mul = Q16_16.q_mul
const q_div = Q16_16.q_div
const q_add = Q16_16.q_add
const q_sub = Q16_16.q_sub

end # module NUVMAP