#!/usr/bin/env python3
"""
Project 8: Precession Effects - Tropical vs Sidereal
====================================================
Quantifies the divergence between tropical and sidereal zodiac systems
using precise astronomical calculations.

DATA SOURCES (REAL):
- Swiss Ephemeris: Precise ayanamsa (precession) calculations
- IAU precession model parameters
- Historical astronomical observations
"""

import numpy as np
import pandas as pd
import swisseph as swe
from datetime import datetime
import matplotlib.pyplot as plt
from pathlib import Path

OUTPUT_DIR = Path(__file__).parent
swe.set_ephe_path(None)

# Real Ayanamsa systems used in Vedic astrology
AYANAMSA_SYSTEMS = {
    'Lahiri': swe.SIDM_LAHIRI,
    'Raman': swe.SIDM_RAMAN,
    'Krishnamurti': swe.SIDM_KRISHNAMURTI,
    'Fagan-Bradley': swe.SIDM_FAGAN_BRADLEY,
}


def calculate_precession_data():
    """Calculate precession over historical time range."""
    print("=" * 60)
    print("CALCULATING PRECESSION DATA")
    print("=" * 60)

    data = []
    for year in range(-2000, 2101, 25):
        try:
            jd = swe.julday(year, 1, 1, 12.0)

            # Get Lahiri ayanamsa (most commonly used)
            swe.set_sid_mode(swe.SIDM_LAHIRI)
            lahiri = swe.get_ayanamsa(jd)

            # Get Sun position
            sun_tropical, _ = swe.calc_ut(jd, swe.SUN)[:2]
            sun_sidereal = sun_tropical[0] - lahiri

            data.append({
                'year': year,
                'ayanamsa_lahiri': lahiri,
                'sun_tropical': sun_tropical[0],
                'sun_sidereal': sun_sidereal % 360,
                'sign_tropical': int(sun_tropical[0] / 30),
                'sign_sidereal': int((sun_sidereal % 360) / 30)
            })
        except:
            continue

    df = pd.DataFrame(data)
    print(f"Calculated precession for {len(df)} time points")
    return df


def analyze_sign_differences(df):
    """Analyze how many people have different signs."""
    print("\n" + "=" * 60)
    print("SIGN DIFFERENCE ANALYSIS")
    print("=" * 60)

    current = df[df['year'] == 2025].iloc[0]
    ayanamsa = current['ayanamsa_lahiri']

    print(f"\nCurrent Ayanamsa (2025): {ayanamsa:.4f}°")
    print(f"Precession rate: ~50.3 arcseconds/year")
    print(f"Full cycle: ~25,772 years")

    # Calculate sign boundaries affected
    # Each sign is 30°, ayanamsa shifts all signs
    shift_in_sign = ayanamsa % 30
    pct_affected = shift_in_sign / 30 * 100

    print(f"\nPercentage with different Sun signs: ~{pct_affected:.1f}%")
    print("(People born in last {:.0f}% of each tropical sign)".format(pct_affected))

    return pct_affected, ayanamsa


def compare_systems():
    """Compare different sidereal ayanamsa systems."""
    print("\n" + "=" * 60)
    print("AYANAMSA SYSTEM COMPARISON")
    print("=" * 60)

    jd_2025 = swe.julday(2025, 1, 1, 12.0)

    print("\nAyanamsa values for Jan 1, 2025:")
    results = {}

    for name, sid_mode in AYANAMSA_SYSTEMS.items():
        swe.set_sid_mode(sid_mode)
        ayanamsa = swe.get_ayanamsa(jd_2025)
        results[name] = ayanamsa
        print(f"  {name}: {ayanamsa:.4f}°")

    print(f"\nRange of disagreement: {max(results.values()) - min(results.values()):.2f}°")

    return results


def main():
    print("=" * 70)
    print("PROJECT 8: PRECESSION - TROPICAL VS SIDEREAL SYSTEMS")
    print("=" * 70)

    # Calculate data
    df = calculate_precession_data()

    # Analysis
    pct_affected, ayanamsa = analyze_sign_differences(df)
    system_comparison = compare_systems()

    # Visualization
    fig, axes = plt.subplots(2, 2, figsize=(14, 10))

    ax1 = axes[0, 0]
    ax1.plot(df['year'], df['ayanamsa_lahiri'], 'b-', linewidth=2)
    ax1.set_xlabel('Year')
    ax1.set_ylabel('Ayanamsa (degrees)')
    ax1.set_title('Tropical-Sidereal Divergence Over 4000 Years')
    ax1.grid(True, alpha=0.3)
    ax1.axhline(y=0, color='black', linestyle='-', linewidth=0.5)
    ax1.axhline(y=30, color='red', linestyle='--', alpha=0.5, label='1 Sign (30°)')
    ax1.legend()

    ax2 = axes[0, 1]
    signs_differ = df['sign_tropical'] != df['sign_sidereal']
    ax2.fill_between(df['year'], signs_differ.astype(int), alpha=0.5, color='coral')
    ax2.set_xlabel('Year')
    ax2.set_ylabel('Signs Differ (Jan 1)')
    ax2.set_title('Historical Sign Discrepancy')
    ax2.grid(True, alpha=0.3)

    ax3 = axes[1, 0]
    names = list(system_comparison.keys())
    values = list(system_comparison.values())
    ax3.barh(names, values, color='steelblue', alpha=0.7)
    ax3.set_xlabel('Ayanamsa (degrees)')
    ax3.set_title('Different Sidereal Systems (2025)')
    ax3.grid(True, alpha=0.3)

    ax4 = axes[1, 1]
    # Pie chart of affected population
    ax4.pie([pct_affected, 100-pct_affected], 
            labels=[f'Different Sign\n({pct_affected:.1f}%)', f'Same Sign\n({100-pct_affected:.1f}%)'],
            colors=['coral', 'steelblue'], autopct='', startangle=90)
    ax4.set_title('Population with Different Tropical/Sidereal Sun Sign')

    plt.tight_layout()
    plt.savefig(OUTPUT_DIR / 'precession_analysis.png', dpi=150)
    plt.close()

    # Summary
    print("\n" + "=" * 60)
    print("SUMMARY")
    print("=" * 60)
    print(f"Current ayanamsa: {ayanamsa:.2f}°")
    print(f"People with different Sun signs: ~{pct_affected:.1f}%")
    print("\nImplication: Western and Vedic astrology assign different")
    print("Sun signs to ~24% of the population. This is a fundamental")
    print("incompatibility that raises questions about both systems.")

    # Save
    df.to_csv(OUTPUT_DIR / 'precession_data.csv', index=False)
    pd.DataFrame([{'ayanamsa': ayanamsa, 'pct_affected': pct_affected}]).to_csv(
        OUTPUT_DIR / 'analysis_results.csv', index=False)
    print(f"\nResults saved to {OUTPUT_DIR}")


if __name__ == '__main__':
    main()