Węzły księżycowe, mean Lilith i wykrywanie stacji (LOG-02/03)

Punkty wirtualne (LOG-02): - engine/points.py: mean Node (Ω) i mean Lilith (apogeum) wzorami Meeusa; prędkości numerycznie. SN = NN + 180° (ta sama prędkość), zawsze Rx. - DEFAULT_OBJECTS + North Node / South Node / Lilith — automatycznie dostają domy, aspekty i A/S. Parzystość silnika B: swe.MEAN_NODE / swe.MEAN_APOG (uwaga: stała pyswisseph to MEAN_APOG, nie MEAN_APOGEE). - significators: tokeny [NN / [SN / [Lilith (zgodne z SIGNIFICATORS KEY). Stacje (LOG-03): - engine/stations.py: skan prędkości (krok 4 dni, okno ±800 dni — pokrywa najdłuższe przerwy Marsa/Wenus) + bisekcja; klasyfikacja SD/SR; poprzednia/ następna stacja (dni, data, stopień w znaku) + flaga station_soon (<7 dni). - /chart/positions: opt-in stations:true; UI: checkbox + tabela stacji. Walidacja: - mean NN vs astro-seek (Gem 8°09'24"): Δ=0,3'; vs swisseph: Δ=17"; mean Lilith vs swisseph: Δ=1,5'. NN dom 12 / SN dom 6 zgodnie z astro-seek. - Stacje Marsa 1984 trafiają w historię: SR 5.04.1984, SD 19.06.1984; samospójność |speed|<0,01°/d w znalezionych momentach; flaga "blisko" działa (Merkury +5,3d, Jowisz -0,6d). - E2E na realnej bazie: [SN 134 rekordy, trafienie w 6. domu. 54 testy przechodzą. Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
2026-07-17 07:02:09 +00:00 · 2026-07-08 10:51:29 +02:00
parent 93932246f3
commit a8c3072e62
13 changed files with 330 additions and 9 deletions
@@ -13,10 +13,12 @@ from typing import Any

 from app.engine import formats

-# kanoniczny zestaw i kolejność obiektów (LOG-02: światła + 7 klasycznych + 3 nowożytne)
+# kanoniczny zestaw i kolejność obiektów (LOG-02: światła + 7 klasycznych +
+# 3 nowożytne + punkty wirtualne: węzły mean i mean Lilith)
 DEFAULT_OBJECTS = [
    "Sun", "Moon", "Mercury", "Venus", "Mars",
    "Jupiter", "Saturn", "Uranus", "Neptune", "Pluto",
+    "North Node", "South Node", "Lilith",
 ]


@@ -0,0 +1,43 @@
+"""Punkty wirtualne liczone analitycznie (LOG-02): mean Node i mean Lilith.
+
+Wzory Meeusa (Astronomical Algorithms) w stuleciach juliańskich od J2000 (TT):
+- Ω — średni węzeł wstępujący orbity Księżyca (mean ascending node). Porusza się
+  zawsze wstecz (~−0,053°/dobę) — stąd węzły są wiecznie Rx.
+- średnie perygeum orbity Księżyca; mean Lilith (Black Moon) = średnie APOGEUM
+  = perygeum + 180° (~+0,111°/dobę).
+
+Wersje TRUE (oskulacyjne) — osobny, późniejszy krok (notatki: mean to
+historyczny standard i domyślne ustawienie programów).
+"""
+from __future__ import annotations
+
+from app.engine.formats import norm360
+
+_DAYS_PER_CENTURY = 36525.0
+
+
+def _t(tt_jd: float) -> float:
+    return (tt_jd - 2451545.0) / _DAYS_PER_CENTURY
+
+
+def mean_lunar_node(tt_jd: float) -> float:
+    """Długość ekliptyczna średniego Węzła Północnego (Ω) [°]."""
+    t = _t(tt_jd)
+    omega = (125.0445479 - 1934.1362891 * t + 0.0020754 * t * t
+             + t ** 3 / 467441.0 - t ** 4 / 60616000.0)
+    return norm360(omega)
+
+
+def mean_lilith(tt_jd: float) -> float:
+    """Długość ekliptyczna mean Lilith (średnie apogeum Księżyca) [°]."""
+    t = _t(tt_jd)
+    perigee = (83.3532465 + 4069.0137287 * t - 0.0103200 * t * t
+               - t ** 3 / 80053.0 + t ** 4 / 18999000.0)
+    return norm360(perigee + 180.0)
+
+
+def point_speed(fn, tt_jd: float, dt_days: float = 0.1) -> float:
+    """Prędkość [°/dobę] punktu analitycznego — różnica po małym kroku."""
+    a = fn(tt_jd)
+    b = fn(tt_jd + dt_days)
+    return (((b - a + 180.0) % 360.0) - 180.0) / dt_days
@@ -60,6 +60,27 @@ class SkyfieldEngine(EphemerisEngine):
        lat, lon, _dist = astrometric.ecliptic_latlon(epoch="date")
        return lon.degrees, lat.degrees

+    def _virtual_point(self, name: str, tt_jd: float) -> ObjectPosition:
+        """Punkty analityczne (LOG-02): mean Node (NN/SN) i mean Lilith.
+
+        Liczone wzorami Meeusa, nie z jądra JPL. SN = NN + 180° (ta sama prędkość).
+        Punkty leżą na ekliptyce (latitude = 0).
+        """
+        from app.engine.points import mean_lilith, mean_lunar_node, point_speed
+
+        if name in ("North Node", "South Node"):
+            lon = mean_lunar_node(tt_jd)
+            if name == "South Node":
+                lon = norm360(lon + 180.0)
+            speed = point_speed(mean_lunar_node, tt_jd)
+        else:  # Lilith
+            lon = mean_lilith(tt_jd)
+            speed = point_speed(mean_lilith, tt_jd)
+        return ObjectPosition(
+            name=name, longitude=float(lon), latitude=0.0,
+            speed=float(speed), retrograde=bool(speed < 0),
+        )
+
    def positions(
        self, moment: ChartMoment, objects: list[str] | None = None
    ) -> list[ObjectPosition]:
@@ -70,6 +91,9 @@ class SkyfieldEngine(EphemerisEngine):

        out: list[ObjectPosition] = []
        for name in names:
+            if name not in _TARGETS:                 # punkt wirtualny (NN/SN/Lilith)
+                out.append(self._virtual_point(name, t.tt))
+                continue
            target = self.eph[_TARGETS[name]]
            lon, lat = self._ecliptic_lon_lat(target, t)
            lon2, _ = self._ecliptic_lon_lat(target, t2)
@@ -0,0 +1,100 @@
+"""Wykrywanie stacji planet (LOG-03): poprzednia/następna stacja, SD/SR, flaga <7 dni.
+
+Stacja ścisła = moment, w którym prędkość zodiakalna przechodzi przez zero.
+Metoda: próbki prędkości co 1 dzień w oknie ± SEARCH_DAYS → zmiana znaku →
+bisekcja do dokładności ~1 minuty. Klasyfikacja: prędkość przed<0 i po>0 → SD
+(stationary direct), odwrotnie → SR (stationary retrograde).
+
+Pomijamy Słońce/Księżyc (nigdy Rx) i punkty mean (NN/SN/Lilith — ruch jednostajny).
+"""
+from __future__ import annotations
+
+from datetime import timedelta
+
+from app.engine.formats import in_sign
+from app.engine.models import ChartMoment
+
+# obiekty bez stacji
+NO_STATIONS = {"Sun", "Moon", "North Node", "South Node", "Lilith"}
+
+# Okno musi pokryć najdłuższą przerwę między stacjami (Mars/Wenus: ~700 dni),
+# a krok skanu musi być krótszy niż najkrótsza retrogradacja (Merkury ~21 dni).
+SEARCH_DAYS = 800          # okno poszukiwań w każdą stronę
+SCAN_STEP_DAYS = 4.0       # krok zgrubnego skanu (potem bisekcja)
+STATION_SOON_DAYS = 7.0    # próg flagi "stacja blisko" (konfigurowalny, notes2)
+
+
+def _speed_fn(engine, name: str):
+    """Zwraca funkcję: dni_od_momentu_bazowego -> prędkość [°/dobę]."""
+    def speed(base_moment: ChartMoment, offset_days: float) -> float:
+        m = ChartMoment(
+            when_utc=base_moment.when_utc + timedelta(days=offset_days),
+            lat=base_moment.lat, lon=base_moment.lon,
+        )
+        return engine.positions(m, [name])[0].speed
+    return speed
+
+
+def _bisect_zero(speed, moment: ChartMoment, lo: float, hi: float, iters: int = 20) -> float:
+    """Bisekcja miejsca zerowego prędkości między dniami lo i hi."""
+    s_lo = speed(moment, lo)
+    for _ in range(iters):
+        mid = (lo + hi) / 2.0
+        s_mid = speed(moment, mid)
+        if (s_lo < 0) == (s_mid < 0):
+            lo, s_lo = mid, s_mid
+        else:
+            hi = mid
+    return (lo + hi) / 2.0
+
+
+def _station_info(engine, moment: ChartMoment, name: str, day: float, speed) -> dict:
+    """Opis stacji w danym dniu (offset od momentu bazowego)."""
+    before = speed(moment, day - 0.5)
+    kind = "SD" if before < 0 else "SR"
+    when = moment.when_utc + timedelta(days=day)
+    m = ChartMoment(when_utc=when, lat=moment.lat, lon=moment.lon)
+    lon = engine.positions(m, [name])[0].longitude
+    return {
+        "type": kind,
+        "date": when.strftime("%Y-%m-%d %H:%M"),
+        "days": round(day, 1),               # ujemne = w przeszłości
+        "degree": in_sign(lon),
+    }
+
+
+def find_stations(engine, moment: ChartMoment, name: str, step_days: float = SCAN_STEP_DAYS) -> dict | None:
+    """Poprzednia i następna stacja obiektu względem momentu horoskopu."""
+    if name in NO_STATIONS:
+        return None
+    speed = _speed_fn(engine, name)
+
+    prev_day = next_day = None
+    # w przeszłość
+    s_right = speed(moment, 0.0)
+    d = 0.0
+    while d > -SEARCH_DAYS:
+        s_left = speed(moment, d - step_days)
+        if (s_left < 0) != (s_right < 0):
+            prev_day = _bisect_zero(speed, moment, d - step_days, d)
+            break
+        d, s_right = d - step_days, s_left
+    # w przyszłość
+    s_left = speed(moment, 0.0)
+    d = 0.0
+    while d < SEARCH_DAYS:
+        s_right = speed(moment, d + step_days)
+        if (s_left < 0) != (s_right < 0):
+            next_day = _bisect_zero(speed, moment, d, d + step_days)
+            break
+        d, s_left = d + step_days, s_right
+
+    result: dict = {}
+    if prev_day is not None:
+        result["prev"] = _station_info(engine, moment, name, prev_day, speed)
+    if next_day is not None:
+        result["next"] = _station_info(engine, moment, name, next_day, speed)
+    result["station_soon"] = any(
+        abs(x["days"]) < STATION_SOON_DAYS for x in result.values() if isinstance(x, dict)
+    )
+    return result or None