Jyotisha/scripts/gulika.py

#!/usr/bin/env python3
"""Swiss-Ephemeris Gulika calculator using the Prasna Marga Ghatika table."""
from __future__ import annotations

from datetime import datetime
from typing import Any

import swisseph as swe

try:
    from saham_daynight import determine_daytime
except ImportError:
    from scripts.saham_daynight import determine_daytime
try:
    from ayanamsa_utils import DEFAULT_AYANAMSA_NAME, ayanamsa_display_name, normalize_ayanamsa_name, temporary_ayanamsa
except ImportError:
    from scripts.ayanamsa_utils import DEFAULT_AYANAMSA_NAME, ayanamsa_display_name, normalize_ayanamsa_name, temporary_ayanamsa


# Monday=0, matching datetime.weekday(). Values are the end of Saturn's share
# measured in Ghatika from the relevant sunrise/sunset (30 Ghatika per period).
GHATIKA_END = {
    0: {"day": 22, "night": 6},
    1: {"day": 18, "night": 2},
    2: {"day": 14, "night": 26},
    3: {"day": 10, "night": 22},
    4: {"day": 6, "night": 18},
    5: {"day": 2, "night": 14},
    6: {"day": 26, "night": 10},
}

# Saturn's zero-based share in eight equal day/night parts. Python weekday:
# Monday=0. This is the variant used by PyJHora's public black-box oracle.
SATURN_PART_START = {
    0: {"day": 5, "night": 1},
    1: {"day": 4, "night": 0},
    2: {"day": 3, "night": 6},
    3: {"day": 2, "night": 5},
    4: {"day": 1, "night": 4},
    5: {"day": 0, "night": 3},
    6: {"day": 6, "night": 2},
}

# Planet indices follow the local seven-graha ordering: Sun through Saturn.
# Each row assigns the eight daylight/nighttime shares for a weekday where
# Python's datetime.weekday() is Monday=0.
DAY_RULERS = (
    (1, 2, 3, 4, 5, 6, -1, 0),
    (2, 3, 4, 5, 6, -1, 0, 1),
    (3, 4, 5, 6, -1, 0, 1, 2),
    (4, 5, 6, -1, 0, 1, 2, 3),
    (5, 6, -1, 0, 1, 2, 3, 4),
    (6, -1, 0, 1, 2, 3, 4, 5),
    (0, 1, 2, 3, 4, 5, 6, -1),
)
NIGHT_RULERS = (
    (5, 6, -1, 0, 1, 2, 3, 4),
    (6, -1, 0, 1, 2, 3, 4, 5),
    (-1, 0, 1, 2, 3, 4, 5, 6),
    (0, 1, 2, 3, 4, 5, 6, -1),
    (1, 2, 3, 4, 5, 6, -1, 0),
    (2, 3, 4, 5, 6, -1, 0, 1),
    (4, 5, 6, -1, 0, 1, 2, 3),
)
TEMPORAL_UPAGRAHAS = {
    "Kaala": (0, "middle"),
    "Mrityu": (2, "middle"),
    "Artha_Praharaka": (3, "middle"),
    "Yama_Ghantaka": (4, "middle"),
    "Gulika": (6, "begin"),
    "Maandi": (6, "middle"),
}


def _sidereal_ascendant(jd_ut: float, lat: float, lon: float, ayanamsa: str = DEFAULT_AYANAMSA_NAME) -> float:
    with temporary_ayanamsa(ayanamsa, swe):
        _, ascmc = swe.houses_ex(jd_ut, lat, lon, b"P", swe.FLG_SWIEPH | swe.FLG_SIDEREAL)
        return float(ascmc[0]) % 360


def calculate_gulika(
    moment: datetime,
    *,
    lat: float,
    lon: float,
    tz: float,
    method: str = "saturn_part_start",
    ayanamsa: str = DEFAULT_AYANAMSA_NAME,
) -> dict[str, Any]:
    """Return Gulika from local moment/location using Swiss sunrise and sunset."""
    ayanamsa_name = normalize_ayanamsa_name(ayanamsa)
    daynight = determine_daytime(moment, lat=lat, lon=lon, tz=tz)
    is_day = bool(daynight["is_daytime"])
    period = "day" if is_day else "night"
    start_jd = daynight["sunrise_jd_ut"] if is_day else daynight["sunset_jd_ut"]
    end_jd = daynight["sunset_jd_ut"] if is_day else daynight["sunrise_jd_ut"] + 1.0
    if end_jd <= start_jd:
        end_jd += 1.0
    start_year, start_month, start_day, _ = swe.revjul(start_jd + float(tz) / 24.0)
    period_weekday = datetime(int(start_year), int(start_month), int(start_day)).weekday()
    if method == "saturn_part_start":
        part_index = SATURN_PART_START[period_weekday][period]
        segment_fraction = part_index / 8.0
        ghatika_end = None
    elif method == "legacy_ghatika_end":
        part_index = None
        ghatika_end = GHATIKA_END[period_weekday][period]
        segment_fraction = ghatika_end / 30.0
    else:
        raise ValueError("method must be saturn_part_start or legacy_ghatika_end")
    segment_jd = start_jd + (end_jd - start_jd) * segment_fraction
    longitude = _sidereal_ascendant(segment_jd, float(lat), float(lon), ayanamsa_name)
    return {
        "scope": "gulika_prasna_marga",
        "status": "partial",
        "longitude": round(longitude, 6),
        "sign_idx": int(longitude / 30) % 12,
        "degree_in_sign": round(longitude % 30, 6),
        "period": period,
        "weekday": period_weekday,
        "moment_weekday": moment.weekday(),
        "method": method,
        "part_index": part_index,
        "segment_fraction": segment_fraction,
        "ghatika_end": ghatika_end,
        "segment_jd_ut": segment_jd,
        "daynight_evidence": daynight,
        "ayanamsa": ayanamsa_name,
        "ayanamsa_display": ayanamsa_display_name(ayanamsa_name),
        "rule_source": "references/prashna-complete-guide.md#3.5",
        "boundary": "PyJHora-aligned Saturn-part-start variant. Numeric parity is evidence only and does not enable Prashna verdict layers.",
    }


def _temporal_upagraha(
    moment: datetime,
    *,
    lat: float,
    lon: float,
    tz: float,
    planet_index: int,
    segment_part: str,
    ayanamsa: str = DEFAULT_AYANAMSA_NAME,
) -> dict[str, Any]:
    ayanamsa_name = normalize_ayanamsa_name(ayanamsa)
    daynight = determine_daytime(moment, lat=lat, lon=lon, tz=tz)
    is_day = bool(daynight["is_daytime"])
    period = "day" if is_day else "night"
    start_jd = daynight["sunrise_jd_ut"] if is_day else daynight["sunset_jd_ut"]
    end_jd = daynight["sunset_jd_ut"] if is_day else daynight["sunrise_jd_ut"] + 1.0
    if end_jd <= start_jd:
        end_jd += 1.0
    start_year, start_month, start_day, _ = swe.revjul(start_jd + float(tz) / 24.0)
    weekday = datetime(int(start_year), int(start_month), int(start_day)).weekday()
    ruler_table = DAY_RULERS if is_day else NIGHT_RULERS
    part_index = ruler_table[weekday].index(planet_index)
    fraction = (part_index + (0.5 if segment_part == "middle" else 0.0)) / 8.0
    segment_jd = start_jd + (end_jd - start_jd) * fraction
    longitude = _sidereal_ascendant(segment_jd, float(lat), float(lon), ayanamsa_name)
    return {
        "longitude": round(longitude, 6),
        "sign_idx": int(longitude / 30) % 12,
        "degree_in_sign": round(longitude % 30, 6),
        "period": period,
        "weekday": weekday,
        "planet_index": planet_index,
        "part_index": part_index,
        "segment_part": segment_part,
        "segment_fraction": fraction,
        "segment_jd_ut": segment_jd,
    }


def _solar_upagrahas(sun_longitude: float) -> dict[str, dict[str, Any]]:
    sun = float(sun_longitude) % 360.0
    longitudes = {
        "Dhuma": (sun + 133.0 + 20.0 / 60.0) % 360.0,
        "Vyatipata": (360.0 - ((sun + 133.0 + 20.0 / 60.0) % 360.0)) % 360.0,
    }
    longitudes["Parivesha"] = (longitudes["Vyatipata"] + 180.0) % 360.0
    longitudes["Indrachapa"] = (360.0 - longitudes["Parivesha"]) % 360.0
    longitudes["Upaketu"] = (sun - 30.0) % 360.0
    return {
        name: {
            "longitude": round(longitude, 6),
            "sign_idx": int(longitude / 30) % 12,
            "degree_in_sign": round(longitude % 30, 6),
            "source": "solar_longitude_formula",
        }
        for name, longitude in longitudes.items()
    }


def calculate_upagrahas(
    moment: datetime,
    *,
    lat: float,
    lon: float,
    tz: float,
    sun_longitude: float,
    ayanamsa: str = DEFAULT_AYANAMSA_NAME,
) -> dict[str, Any]:
    """Calculate the D1 Upagraha packet using local Swiss Ephemeris inputs."""
    ayanamsa_name = normalize_ayanamsa_name(ayanamsa)
    temporal = {
        name: _temporal_upagraha(
            moment,
            lat=lat,
            lon=lon,
            tz=tz,
            planet_index=planet_index,
            segment_part=segment_part,
            ayanamsa=ayanamsa_name,
        )
        for name, (planet_index, segment_part) in TEMPORAL_UPAGRAHAS.items()
    }
    return {
        "status": "computed_pending_validation",
        "raw": {**temporal, **_solar_upagrahas(sun_longitude)},
        "ayanamsa": ayanamsa_name,
        "ayanamsa_display": ayanamsa_display_name(ayanamsa_name),
        "scope": "d1_only",
        "boundary": (
            "Local D1 formulas are retained with segment evidence. Public PyJHora observations are "
            "used only for calibration; D-N Upagraha mapping and external numeric closure remain pending."
        ),
    }