-- orbital_elements regression tests
--
-- Tests orbital_elements type I/O, accessors, MPC parser,
-- small_body_heliocentric(), and small_body_observe().
CREATE EXTENSION IF NOT EXISTS pg_orrery;

\set boulder '''40.015N 105.270W 1655m'''::observer

-- ============================================================
-- Test 1: Type I/O round-trip
-- Construct orbital_elements from text literal, verify output matches.
-- Angles in degrees, stored as radians internally.
-- ============================================================
SELECT 'io_roundtrip' AS test,
       '(2460600.000000,2.5478000000,0.0789126000,10.586640,73.429370,80.268600,2460319.000000,3.33,0.12)'::orbital_elements AS oe;

-- ============================================================
-- Test 2: Accessor functions (direct field access)
-- ============================================================
SELECT 'accessor_epoch' AS test,
       round(oe_epoch('(2460600.0,2.5478,0.0789126,10.58664,73.42937,80.2686,2460319.0,3.33,0.12)'::orbital_elements)::numeric, 1) AS epoch_jd;

SELECT 'accessor_q' AS test,
       round(oe_perihelion('(2460600.0,2.5478,0.0789126,10.58664,73.42937,80.2686,2460319.0,3.33,0.12)'::orbital_elements)::numeric, 4) AS q_au;

SELECT 'accessor_e' AS test,
       round(oe_eccentricity('(2460600.0,2.5478,0.0789126,10.58664,73.42937,80.2686,2460319.0,3.33,0.12)'::orbital_elements)::numeric, 7) AS ecc;

SELECT 'accessor_inc' AS test,
       round(oe_inclination('(2460600.0,2.5478,0.0789126,10.58664,73.42937,80.2686,2460319.0,3.33,0.12)'::orbital_elements)::numeric, 5) AS inc_deg;

SELECT 'accessor_omega' AS test,
       round(oe_arg_perihelion('(2460600.0,2.5478,0.0789126,10.58664,73.42937,80.2686,2460319.0,3.33,0.12)'::orbital_elements)::numeric, 5) AS omega_deg;

SELECT 'accessor_Omega' AS test,
       round(oe_raan('(2460600.0,2.5478,0.0789126,10.58664,73.42937,80.2686,2460319.0,3.33,0.12)'::orbital_elements)::numeric, 4) AS Omega_deg;

SELECT 'accessor_tp' AS test,
       round(oe_tp('(2460600.0,2.5478,0.0789126,10.58664,73.42937,80.2686,2460319.0,3.33,0.12)'::orbital_elements)::numeric, 1) AS tp_jd;

SELECT 'accessor_h' AS test,
       round(oe_h_mag('(2460600.0,2.5478,0.0789126,10.58664,73.42937,80.2686,2460319.0,3.33,0.12)'::orbital_elements)::numeric, 2) AS h_mag;

SELECT 'accessor_g' AS test,
       round(oe_g_slope('(2460600.0,2.5478,0.0789126,10.58664,73.42937,80.2686,2460319.0,3.33,0.12)'::orbital_elements)::numeric, 2) AS g_slope;

-- ============================================================
-- Test 3: Computed accessors
-- ============================================================

-- Semi-major axis: a = q / (1 - e) = 2.5478 / (1 - 0.0789126) ~ 2.766 AU
SELECT 'sma_elliptic' AS test,
       round(oe_semi_major_axis('(2460600.0,2.5478,0.0789126,10.58664,73.42937,80.2686,2460319.0,3.33,0.12)'::orbital_elements)::numeric, 3) AS a_au;

-- Period: a^1.5 years.  a ~ 2.766, period ~ 4.599 years
SELECT 'period_elliptic' AS test,
       round(oe_period_years('(2460600.0,2.5478,0.0789126,10.58664,73.42937,80.2686,2460319.0,3.33,0.12)'::orbital_elements)::numeric, 2) AS period_yr;

-- Hyperbolic orbit (e=1.5): semi-major axis and period should be NULL
SELECT 'sma_hyperbolic' AS test,
       oe_semi_major_axis('(2460600.0,1.0,1.5,10.0,73.0,80.0,2460319.0,NaN,NaN)'::orbital_elements) IS NULL AS is_null;

SELECT 'period_hyperbolic' AS test,
       oe_period_years('(2460600.0,1.0,1.5,10.0,73.0,80.0,2460319.0,NaN,NaN)'::orbital_elements) IS NULL AS is_null;

-- ============================================================
-- Test 4: MPC parser -- (1) Ceres
--
-- MPCORB.DAT line for Ceres (epoch 2024 Oct 17.0 = K24AM):
-- Packed epoch K24AM -> K=2000, 24, A=Oct, M=22 -> 2024-10-22
-- (Actually: K=2000, year=24, month=A=10, day=M=22)
-- JD for 2024-10-22 = 2460605.5
--
-- a = 2.7660961 AU, e = 0.0789126
-- q = a*(1-e) = 2.7660961*(1-0.0789126) = 2.5478...
-- M = 60.07966 deg
-- ============================================================
SELECT 'mpc_ceres_sma' AS test,
       round(oe_semi_major_axis(oe_from_mpc(
           '00001    3.33  0.12 K24AM  60.07966  73.42937   80.26860  10.58664  0.0789126  0.21406048   2.7660961  0 MPO838504  8738 115 1801-2024 0.65 M-v 30k MPCLINUX   0000      (1) Ceres              20240825'
       ))::numeric, 4) AS a_au;

SELECT 'mpc_ceres_q' AS test,
       round(oe_perihelion(oe_from_mpc(
           '00001    3.33  0.12 K24AM  60.07966  73.42937   80.26860  10.58664  0.0789126  0.21406048   2.7660961  0 MPO838504  8738 115 1801-2024 0.65 M-v 30k MPCLINUX   0000      (1) Ceres              20240825'
       ))::numeric, 4) AS q_au;

SELECT 'mpc_ceres_ecc' AS test,
       round(oe_eccentricity(oe_from_mpc(
           '00001    3.33  0.12 K24AM  60.07966  73.42937   80.26860  10.58664  0.0789126  0.21406048   2.7660961  0 MPO838504  8738 115 1801-2024 0.65 M-v 30k MPCLINUX   0000      (1) Ceres              20240825'
       ))::numeric, 7) AS ecc;

SELECT 'mpc_ceres_inc' AS test,
       round(oe_inclination(oe_from_mpc(
           '00001    3.33  0.12 K24AM  60.07966  73.42937   80.26860  10.58664  0.0789126  0.21406048   2.7660961  0 MPO838504  8738 115 1801-2024 0.65 M-v 30k MPCLINUX   0000      (1) Ceres              20240825'
       ))::numeric, 5) AS inc_deg;

SELECT 'mpc_ceres_h' AS test,
       round(oe_h_mag(oe_from_mpc(
           '00001    3.33  0.12 K24AM  60.07966  73.42937   80.26860  10.58664  0.0789126  0.21406048   2.7660961  0 MPO838504  8738 115 1801-2024 0.65 M-v 30k MPCLINUX   0000      (1) Ceres              20240825'
       ))::numeric, 2) AS h_mag;

-- ============================================================
-- Test 5: small_body_heliocentric -- compare to kepler_propagate
-- Both should produce the same heliocentric position for Ceres.
-- ============================================================
SELECT 'helio_vs_kepler' AS test,
       round(helio_x(small_body_heliocentric(
           oe_from_mpc('00001    3.33  0.12 K24AM  60.07966  73.42937   80.26860  10.58664  0.0789126  0.21406048   2.7660961  0 MPO838504  8738 115 1801-2024 0.65 M-v 30k MPCLINUX   0000      (1) Ceres              20240825'),
           '2025-01-01 00:00:00+00'
       ))::numeric, 6) AS x_au,
       round(helio_y(small_body_heliocentric(
           oe_from_mpc('00001    3.33  0.12 K24AM  60.07966  73.42937   80.26860  10.58664  0.0789126  0.21406048   2.7660961  0 MPO838504  8738 115 1801-2024 0.65 M-v 30k MPCLINUX   0000      (1) Ceres              20240825'),
           '2025-01-01 00:00:00+00'
       ))::numeric, 6) AS y_au,
       round(helio_z(small_body_heliocentric(
           oe_from_mpc('00001    3.33  0.12 K24AM  60.07966  73.42937   80.26860  10.58664  0.0789126  0.21406048   2.7660961  0 MPO838504  8738 115 1801-2024 0.65 M-v 30k MPCLINUX   0000      (1) Ceres              20240825'),
           '2025-01-01 00:00:00+00'
       ))::numeric, 6) AS z_au;

-- Verify heliocentric distance is reasonable (Ceres orbits ~2.5-3.0 AU)
SELECT 'helio_distance' AS test,
       round(helio_distance(small_body_heliocentric(
           oe_from_mpc('00001    3.33  0.12 K24AM  60.07966  73.42937   80.26860  10.58664  0.0789126  0.21406048   2.7660961  0 MPO838504  8738 115 1801-2024 0.65 M-v 30k MPCLINUX   0000      (1) Ceres              20240825'),
           '2025-01-01 00:00:00+00'
       ))::numeric, 2) AS dist_au,
       helio_distance(small_body_heliocentric(
           oe_from_mpc('00001    3.33  0.12 K24AM  60.07966  73.42937   80.26860  10.58664  0.0789126  0.21406048   2.7660961  0 MPO838504  8738 115 1801-2024 0.65 M-v 30k MPCLINUX   0000      (1) Ceres              20240825'),
           '2025-01-01 00:00:00+00'
       )) BETWEEN 2.5 AND 3.0 AS in_range;

-- ============================================================
-- Test 6: small_body_observe -- topocentric observation
-- Verify we get reasonable az/el/range for Ceres from Boulder.
-- Range should be on the order of 1.5-4.5 AU in km.
-- ============================================================
SELECT 'observe_range' AS test,
       topo_range(small_body_observe(
           oe_from_mpc('00001    3.33  0.12 K24AM  60.07966  73.42937   80.26860  10.58664  0.0789126  0.21406048   2.7660961  0 MPO838504  8738 115 1801-2024 0.65 M-v 30k MPCLINUX   0000      (1) Ceres              20240825'),
           :boulder, '2025-01-01 00:00:00+00'
       )) BETWEEN 1.5 * 149597870.7 AND 4.5 * 149597870.7 AS range_reasonable;

-- Elevation should be a finite number
SELECT 'observe_elevation' AS test,
       topo_elevation(small_body_observe(
           oe_from_mpc('00001    3.33  0.12 K24AM  60.07966  73.42937   80.26860  10.58664  0.0789126  0.21406048   2.7660961  0 MPO838504  8738 115 1801-2024 0.65 M-v 30k MPCLINUX   0000      (1) Ceres              20240825'),
           :boulder, '2025-01-01 00:00:00+00'
       )) BETWEEN -90 AND 90 AS el_reasonable;

-- ============================================================
-- Test 7: Parabolic/hyperbolic elements via text I/O
-- Verify type handles e >= 1 without error.
-- ============================================================
SELECT 'parabolic_io' AS test,
       oe_eccentricity('(2460600.0,1.0,1.0,45.0,90.0,180.0,2460500.0,12.0,0.15)'::orbital_elements) AS ecc;

SELECT 'hyperbolic_io' AS test,
       oe_eccentricity('(2460600.0,0.5,2.5,30.0,60.0,120.0,2460400.0,8.0,0.04)'::orbital_elements) AS ecc;

-- ============================================================
-- Test 8: Error paths
-- ============================================================

-- Invalid text input (wrong number of fields)
SELECT 'bad_text' AS test, '(1.0,2.0,3.0)'::orbital_elements;

-- Negative perihelion distance
SELECT 'negative_q' AS test, '(2460600.0,-1.0,0.5,10.0,73.0,80.0,2460319.0,3.33,0.12)'::orbital_elements;

-- Negative eccentricity
SELECT 'negative_e' AS test, '(2460600.0,1.0,-0.1,10.0,73.0,80.0,2460319.0,3.33,0.12)'::orbital_elements;

-- MPC line too short
SELECT 'mpc_short' AS test, oe_from_mpc('too short');

-- ============================================================
-- Test 9: MPC packed date decoding
-- Verify K24AM -> 2024-10-22 -> JD 2460605.5
-- ============================================================
SELECT 'mpc_epoch' AS test,
       round(oe_epoch(oe_from_mpc(
           '00001    3.33  0.12 K24AM  60.07966  73.42937   80.26860  10.58664  0.0789126  0.21406048   2.7660961  0 MPO838504  8738 115 1801-2024 0.65 M-v 30k MPCLINUX   0000      (1) Ceres              20240825'
       ))::numeric, 1) AS epoch_jd;

-- ============================================================
-- Test 10: Consistency -- small_body_observe vs comet_observe
-- Both pipelines should produce the same topocentric result
-- when fed the same elements and Earth position.
-- ============================================================
WITH ceres AS (
    SELECT oe_from_mpc(
        '00001    3.33  0.12 K24AM  60.07966  73.42937   80.26860  10.58664  0.0789126  0.21406048   2.7660961  0 MPO838504  8738 115 1801-2024 0.65 M-v 30k MPCLINUX   0000      (1) Ceres              20240825'
    ) AS oe
), earth AS (
    SELECT planet_heliocentric(3, '2025-06-15 12:00:00+00') AS pos
)
SELECT 'pipeline_match' AS test,
       round(abs(
           topo_elevation(small_body_observe(c.oe, :boulder, '2025-06-15 12:00:00+00'))
           - topo_elevation(comet_observe(
               oe_perihelion(c.oe), oe_eccentricity(c.oe),
               oe_inclination(c.oe), oe_arg_perihelion(c.oe), oe_raan(c.oe),
               oe_tp(c.oe),
               helio_x(e.pos), helio_y(e.pos), helio_z(e.pos),
               :boulder, '2025-06-15 12:00:00+00'
           ))
       )::numeric, 6) AS el_diff_deg
FROM ceres c, earth e;
