"""Unit evidence for phase-continuous Galileo E1 CBOC replicas."""

from __future__ import annotations

import numpy as np

from galileo_rx.signals.boc import SquareSubcarrierNco, boc_subcarrier
from galileo_rx.signals.e1_modulation import (
    CBOC_ALPHA,
    CBOC_BETA,
    E1Modulation,
    e1_subcarrier_replica,
    e1_transmitted_subcarrier_replica,
)


def _phases(sample_count: int, sampling_frequency_hz: float) -> np.ndarray:
    return np.arange(sample_count, dtype=np.float64) * 1.023e6 / sampling_frequency_hz


def test_boc11_subcarrier_sequence() -> None:
    wave, ending = boc_subcarrier(8, 8.0, 1.0)
    np.testing.assert_array_equal(wave, [1.0, 1.0, 1.0, 1.0, -1.0, -1.0, -1.0, -1.0])
    assert ending == 0.0


def test_boc61_subcarrier_sequence() -> None:
    wave, ending = boc_subcarrier(24, 24.0, 6.0)
    np.testing.assert_array_equal(wave, np.tile([1.0, 1.0, -1.0, -1.0], 6))
    assert ending == 0.0


def test_subcarrier_phase_continuity_across_arbitrary_chunks() -> None:
    fs = 31.25e6  # deliberately nonintegral samples per BOC(6,1) sign interval
    total = 10_003
    reference = SquareSubcarrierNco(fs, 6.138e6, 0.137).samples(total)
    chunked_nco = SquareSubcarrierNco(fs, 6.138e6, 0.137)
    chunks = [1, 7, 31, 2, 997, 4093, 3, 4869]
    chunked = np.concatenate([chunked_nco.samples(count) for count in chunks])
    np.testing.assert_array_equal(chunked, reference)


def test_cboc_power_normalization_and_component_sign() -> None:
    phase = _phases(4092 * 32, 32.736e6)
    sc11 = e1_subcarrier_replica("E1B", E1Modulation.BOC11, phase)
    sc61 = np.where((np.floor(2.0 * np.mod(6.0 * phase, 1.0) + 1e-12).astype(int) & 1) == 0, 1.0, -1.0)
    e1b = e1_subcarrier_replica("E1B", E1Modulation.CBOC, phase, 6.0 * phase)
    e1c = e1_subcarrier_replica("E1C", E1Modulation.CBOC, phase, 6.0 * phase)
    assert abs(float(np.mean(e1b**2)) - 1.0) < 1e-12
    assert abs(float(np.mean(e1c**2)) - 1.0) < 1e-12
    np.testing.assert_allclose((e1b + e1c) / 2.0, CBOC_ALPHA * sc11, atol=1e-15)
    np.testing.assert_allclose((e1b - e1c) / 2.0, CBOC_BETA * sc61, atol=1e-15)
    transmitted_e1b = e1_transmitted_subcarrier_replica(
        "E1B", E1Modulation.CBOC, phase, 6.0 * phase
    )
    transmitted_e1c = e1_transmitted_subcarrier_replica(
        "E1C", E1Modulation.CBOC, phase, 6.0 * phase
    )
    # Equation (12): low-rate terms oppose, high-rate terms have equal sign.
    np.testing.assert_allclose((transmitted_e1b - transmitted_e1c) / 2.0, CBOC_ALPHA * sc11, atol=1e-15)
    np.testing.assert_allclose((transmitted_e1b + transmitted_e1c) / 2.0, CBOC_BETA * sc61, atol=1e-15)


def test_cboc_spectrum_contains_the_boc61_lobe() -> None:
    fs = 32.736e6
    phase = _phases(4092 * 32, fs)
    boc11 = e1_subcarrier_replica("E1B", E1Modulation.BOC11, phase)
    cboc = e1_subcarrier_replica("E1B", E1Modulation.CBOC, phase, 6.0 * phase)
    frequencies = np.fft.rfftfreq(phase.size, 1.0 / fs)
    boc_psd = np.abs(np.fft.rfft(boc11)) ** 2
    cboc_psd = np.abs(np.fft.rfft(cboc)) ** 2
    high_lobe = (frequencies > 5.8e6) & (frequencies < 6.5e6)
    fundamental = (frequencies > 0.8e6) & (frequencies < 1.25e6)
    assert cboc_psd[high_lobe].sum() > 1_000.0 * boc_psd[high_lobe].sum()
    assert cboc_psd[fundamental].sum() > cboc_psd[high_lobe].sum()


def test_cboc_autocorrelation_side_structure_and_epl_symmetry() -> None:
    fs = 32.736e6
    samples_per_chip = 32
    phase = _phases(4092 * samples_per_chip, fs)
    cboc = e1_subcarrier_replica("E1B", E1Modulation.CBOC, phase, 6.0 * phase)

    def circular_correlation(offset_samples: int) -> float:
        return float(np.mean(cboc * np.roll(cboc, offset_samples)))

    assert circular_correlation(0) > 0.999999
    assert circular_correlation(samples_per_chip // 2) < -0.80
    for offset in (2, 4, 6, 8, 10, 12):
        assert abs(circular_correlation(offset) - circular_correlation(-offset)) < 2e-4
    # E and L at +/-0.25 chip see equal magnitude at the correct prompt.
    assert abs(abs(circular_correlation(8)) - abs(circular_correlation(-8))) < 2e-4


def test_cboc_sample_rate_independence() -> None:
    powers = []
    quarter_chip_correlations = []
    for fs in (32.736e6, 30.0e6, 27.5e6):
        count = int(round(4092 * fs / 1.023e6))
        phase = _phases(count, fs)
        cboc = e1_subcarrier_replica("E1C", E1Modulation.CBOC, phase, 6.0 * phase)
        powers.append(float(np.mean(cboc**2)))
        quarter_late = e1_subcarrier_replica(
            "E1C", E1Modulation.CBOC, phase + 0.25, 6.0 * (phase + 0.25)
        )
        quarter_chip_correlations.append(float(np.mean(cboc * quarter_late)))
    assert max(abs(power - 1.0) for power in powers) < 2e-3
    assert np.ptp(quarter_chip_correlations) < 0.04