"""Только программные нарушения IQ; эталоны сюда не передаются приёмнику."""

from __future__ import annotations

import json
import math
from pathlib import Path
import statistics
from typing import cast

import numpy as np
from scipy.optimize import brentq
from scipy.signal import lfilter

from protocol.video003_modem import SAMPLE_RATE_HZ, samples_per_symbol
from protocol.video005_iq_waveform import NOMINAL_CARRIER_HZ, SYMBOL_RATE


ROOT = Path(__file__).resolve().parents[1]
LAB066_RESULT = ROOT / "data/processed/lab066/rf_result.json"


def impair(
    iq: np.ndarray,
    *,
    timing_offset_symbols: float = 0.0,
    cfo_hz: float = 0.0,
    phase_degrees: float = 0.0,
    ebn0_db: float | None = None,
    seed: int = 20260917,
    symbol_rate: int = SYMBOL_RATE,
) -> np.ndarray:
    """Сдвиг времени, CFO, фаза и AWGN с фиксированной затравкой."""
    if not all(math.isfinite(value) for value in (timing_offset_symbols, cfo_hz, phase_degrees)):
        raise ValueError("параметры канала должны быть конечными")
    if ebn0_db is not None and not math.isfinite(ebn0_db):
        raise ValueError("Eb/N0 должно быть конечным")
    transmitted = np.asarray(iq, dtype=np.complex64)
    received = np.empty_like(transmitted)
    delay = timing_offset_symbols * samples_per_symbol(symbol_rate)
    # Ограничиваем рабочие массивы: полная видеозапись занимает миллионы IQ-отсчётов.
    random_generator = np.random.default_rng(seed)
    sigma = (math.sqrt(1 / (2 * (2 / 3) * 10 ** (ebn0_db / 10)) / 2)
             if ebn0_db is not None else 0.0)
    for start in range(0, transmitted.size, 1_000_000):
        stop = min(start + 1_000_000, transmitted.size)
        indexes = np.arange(start, stop, dtype=np.float64)
        if delay:
            positions = indexes - delay
            lower = np.floor(positions).astype(np.int64)
            valid = (lower >= 0) & (lower + 1 < transmitted.size)
            part = np.zeros(stop - start, dtype=np.complex64)
            fraction = positions[valid] - lower[valid]
            part[valid] = (transmitted[lower[valid]] * (1 - fraction)
                           + transmitted[lower[valid] + 1] * fraction)
        else:
            part = transmitted[start:stop].copy()
        if cfo_hz or phase_degrees:
            phase = math.radians(phase_degrees) + 2 * math.pi * cfo_hz * indexes / SAMPLE_RATE_HZ
            part *= np.exp(1j * phase).astype(np.complex64)
        if ebn0_db is not None:
            part += (sigma * (random_generator.standard_normal(stop - start)
                              + 1j * random_generator.standard_normal(stop - start))
                     ).astype(np.complex64)
        received[start:stop] = part
    return received


def lab066_like(iq: np.ndarray, *, seed: int = 20260917) -> tuple[np.ndarray, dict]:
    """Диагностический суррогат по двум измеренным сводным характеристикам.

    Частотная рябь и добавочный остаток не воспроизводят точную OTA-передачу.
    """
    report = json.loads(LAB066_RESULT.read_text(encoding="utf-8"))
    rows = report["physical_diagnostics"]["rows"]
    ripple_db = statistics.median(row["band_gain_nonuniformity_p95_p05_db"] for row in rows)
    useful_to_residual_db = statistics.mean(row["matched_residual_db"] for row in rows)
    frequencies = np.linspace(NOMINAL_CARRIER_HZ - 337_500,
                              NOMINAL_CARRIER_HZ + 337_500, 2048)
    delay_samples = 8

    def modeled_ripple(alpha: float) -> float:
        response = np.abs(1 + alpha * np.exp(-2j * np.pi * frequencies
                                             * delay_samples / SAMPLE_RATE_HZ))
        return 20 * math.log10(float(np.percentile(response, 95)
                                     / np.percentile(response, 5)))

    alpha = cast(float, brentq(lambda value: modeled_ripple(value) - ripple_db,
                              0.0, 0.8, full_output=False))
    coefficients = np.zeros(delay_samples + 1)
    coefficients[0], coefficients[-1] = 1.0, alpha
    filtered = cast(np.ndarray, lfilter(coefficients, [1.0],
                                       np.asarray(iq, dtype=np.complex64)))
    residual_ratio = 10 ** (useful_to_residual_db / 10)
    sigma = math.sqrt(1 / (2 * residual_ratio))
    random_generator = np.random.default_rng(seed)
    result = filtered + sigma * (
        random_generator.standard_normal(filtered.size)
        + 1j * random_generator.standard_normal(filtered.size)
    )
    return result.astype(np.complex64), {
        "label": "LAB066-LIKE DIAGNOSTIC CHANNEL MODEL",
        "measured_band_p95_p05_db": ripple_db,
        "modeled_band_p95_p05_db": modeled_ripple(alpha),
        "measured_matched_residual_db": useful_to_residual_db,
        "filter_delay_samples": delay_samples,
        "filter_second_tap": alpha,
        "limitation": "FIR ripple and white residual surrogate do not reproduce OTA propagation",
    }