astrololo/services/logic/app/engine/houses.py

"""Osie i domy — czysta matematyka sferyczna (LOG-05).

Bezstanowe funkcje: z lokalnego czasu gwiazdowego (RAMC), nachylenia ekliptyki
(ε) i szerokości geograficznej (φ) wyliczają Ascendent i MC, a stąd cusps domów
dla prostych systemów (Whole Sign, Equal, Porphyry). Niezależne od silnika —
silnik dostarcza tylko RAMC i ε.
"""
from __future__ import annotations

import math

from app.engine.formats import SIGN_ABBR, norm360, sign_index  # noqa: F401

WHOLE_SIGN = "whole_sign"
EQUAL = "equal"
PORPHYRY = "porphyry"
SYSTEMS = (WHOLE_SIGN, EQUAL, PORPHYRY)


def mean_obliquity(tt_jd: float) -> float:
    """Średnie nachylenie ekliptyki [°] dla daty (Julian TT). Wystarcza do domów."""
    t = (tt_jd - 2451545.0) / 36525.0
    arcsec = 84381.448 - 46.8150 * t - 0.00059 * t * t + 0.001813 * t ** 3
    return arcsec / 3600.0


def compute_mc(ramc_deg: float, eps_deg: float) -> float:
    r, e = math.radians(ramc_deg), math.radians(eps_deg)
    mc = math.atan2(math.sin(r), math.cos(r) * math.cos(e))
    return norm360(math.degrees(mc))


def compute_asc(ramc_deg: float, eps_deg: float, lat_deg: float) -> float:
    """Ascendent — punkt ekliptyki wschodzący na horyzoncie.

    KOREKTA GAŁĘZI (błąd wykryty przez porównanie z wyrocznią, tests/oracle):
    ekliptyka przecina horyzont w DWÓCH punktach — wschodzącym (Asc) i zachodzącym
    (Dsc). `atan2` wybiera jeden z nich, ale powyżej koła podbiegunowego potrafi
    wskazać ten NIEWŁAŚCIWY: dla szerokości 67°+ i szerokiego zakresu RAMC
    zwracaliśmy Descendent, czyli Ascendent przesunięty o 180°. Skutek nie był
    subtelny — planety lądowały w PRZECIWNYCH domach dla całej Skandynawii
    północnej (Tromsø, Rovaniemi, Murmańsk).

    Rozstrzyga położenie względem MC: punkt wschodzący leży zawsze w półkolu
    (0°, 180°) na wschód od MC. Reguła zweryfikowana na 46 800 przypadkach wobec
    Swiss Ephemeris — zero rozbieżności.
    """
    r, e, phi = math.radians(ramc_deg), math.radians(eps_deg), math.radians(lat_deg)
    asc = norm360(math.degrees(math.atan2(
        math.cos(r),
        -(math.sin(r) * math.cos(e) + math.tan(phi) * math.sin(e)),
    )))
    mc = compute_mc(ramc_deg, eps_deg)
    return norm360(asc + 180.0) if (asc - mc) % 360.0 > 180.0 else asc


def _trisect(a: float, b: float) -> tuple[float, float]:
    """Dwa punkty dzielące łuk a→b (w kierunku zodiaku) na trzy równe części."""
    span = (b - a) % 360.0
    return norm360(a + span / 3.0), norm360(a + 2.0 * span / 3.0)


# Poniżej tej odległości od granicy znaku traktujemy Ascendent jak leżący DOKŁADNIE
# na niej. 1e-9° to 3,6 mikrosekundy łuku — o rzędy wielkości poniżej jakiejkolwiek
# realnej dokładności danych urodzeniowych, więc nie zmienia to żadnego horoskopu.
_SIGN_SNAP_DEG = 1e-9


def _snap_to_sign_boundary(lon: float) -> float:
    """Przyciąga długość do granicy znaku, gdy jest od niej o włos.

    Whole sign jest NIECIĄGŁY na granicach znaków: różnica 10⁻¹¹° w Ascendencie
    przerzuca cały dom I o 30°. Bez tego przyciągania ten sam horoskop policzony
    na innej maszynie mógłby dać inny wynik (wykryte przez porównanie z wyrocznią:
    nasz Asc = 359,999999999976, swissepha = 1e-10 — ta sama wartość po dwóch
    stronach granicy Ryby/Baran). Determinizm jest tu ważniejszy niż dosłowność
    zmiennoprzecinkowa."""
    nearest = round(lon / 30.0) * 30.0
    return norm360(nearest) if abs(lon - nearest) < _SIGN_SNAP_DEG else lon


def cusps(asc: float, mc: float, system: str) -> list[float]:
    """Zwraca 12 cusps (długości) domów 1..12."""
    if system == WHOLE_SIGN:
        start = sign_index(_snap_to_sign_boundary(asc)) * 30.0
        return [norm360(start + 30.0 * i) for i in range(12)]
    if system == EQUAL:
        return [norm360(asc + 30.0 * i) for i in range(12)]
    if system == PORPHYRY:
        dsc, ic = norm360(asc + 180.0), norm360(mc + 180.0)
        c = [0.0] * 12
        c[0], c[3], c[6], c[9] = asc, ic, dsc, mc
        c[1], c[2] = _trisect(asc, ic)      # domy 2,3
        c[4], c[5] = _trisect(ic, dsc)      # domy 5,6
        c[7], c[8] = _trisect(dsc, mc)      # domy 8,9
        c[10], c[11] = _trisect(mc, asc)    # domy 11,12
        return c
    raise ValueError(f"nieznany system domów: {system}")


def polar_circle(eps_deg: float) -> float:
    """Szerokość koła podbiegunowego [°] dla danego nachylenia ekliptyki.

    NIE jest to stała 66,56°: ε zmienia się z datą (ok. 23,71° w 370 p.n.e.,
    23,44° dziś), więc granica przesuwa się o ~0,3° w zakresie dat programu.
    Powyżej niej stopnie ekliptyki bywają okołobiegunowe — nie wschodzą ani nie
    zachodzą — przez co systemy oparte na łuku dobowym (Placidus, Koch) tracą
    definicję."""
    return 90.0 - abs(eps_deg)


def cusps_for(ramc: float, eps: float, lat: float, system: str) -> list[float]:
    """Kanoniczne wejście: (RAMC, ε, φ) → 12 cusps.

    Systemy proste (whole sign / equal / porphyry) potrzebują tylko Asc i MC,
    ale systemy egzotyczne dzielą inne koła wielkie i wymagają pełnego zestawu
    (RAMC, ε, φ). Ta funkcja jest wspólnym punktem wejścia dla obu rodzajów —
    i to ją porównuje z wyrocznią framework testowy (tests/oracle).
    """
    asc = compute_asc(ramc, eps, lat)
    mc = compute_mc(ramc, eps)
    return cusps(asc, mc, system)


def assign_house(lon: float, cusp_list: list[float]) -> int:
    """Numer domu (1..12), w którym leży dana długość ekliptyczna."""
    lon = norm360(lon)
    for i in range(12):
        start = cusp_list[i]
        end = cusp_list[(i + 1) % 12]
        span = (end - start) % 360.0
        offset = (lon - start) % 360.0
        if offset < span:
            return i + 1
    return 12