"""Shared constants, paths and small helpers for the anomalous_small_bodies study."""
import os
import numpy as np

STUDY = os.path.dirname(os.path.dirname(os.path.abspath(__file__)))
DATA = os.path.join(STUDY, "data")
FIGS = os.path.join(STUDY, "figures")
RAW = os.environ.get("ASB_RAW", "fermi_raw/anomalous_small_bodies")
for _d in (DATA, FIGS, RAW):
    os.makedirs(_d, exist_ok=True)

# ---------------------------------------------------------------- constants
AU_M = 1.495978707e11            # m
DAY_S = 86400.0
AUD2_TO_MS2 = AU_M / DAY_S**2    # 1 au/d^2 in m/s^2  (= 20.04)
GAUSS_K = 0.01720209895
MU_SUN = GAUSS_K**2              # au^3/d^2
C_AUD = 299792.458 * DAY_S / (AU_M / 1e3)   # speed of light, au/d (173.14)
SOLAR_CONST = 1361.0             # W/m^2 at 1 au
C_MS = 299792458.0
P_SRP = SOLAR_CONST / C_MS       # N/m^2 at 1 au, perfect absorber (4.54e-6)
G_SUN_1AU = MU_SUN * AUD2_TO_MS2 # solar gravity at 1 au in m/s^2 (5.93e-3)

# A1 [au/d^2] produced by solar radiation pressure on a body with
# area-to-mass ratio AMR [m^2/kg] and radiation-pressure coefficient C_R:
#     A1 = C_R * P_SRP * AMR / AUD2_TO_MS2
A1_PER_AMR = P_SRP / AUD2_TO_MS2   # au/d^2 per (m^2/kg), C_R = 1  (2.27e-7)

# size bins: D = 1329 km * 10^(-H/5) / sqrt(p_V)
PV_NOM, PV_LO, PV_HI = 0.14, 0.04, 0.45


def diam_m(H, pv=PV_NOM):
    return 1.329e6 * 10 ** (-np.asarray(H, float) / 5.0) / np.sqrt(pv)


def H_of_diam(D_m, pv=PV_NOM):
    return -5.0 * np.log10(D_m * np.sqrt(pv) / 1.329e6)


H_10M = float(H_of_diam(10.0))     # 27.75
H_100M = float(H_of_diam(100.0))   # 22.75
H_1M = float(H_of_diam(1.0))       # 32.75


def size_bin(H):
    H = np.asarray(H, float)
    return np.where(H > H_10M, "1-10 m", np.where(H > H_100M, "10-100 m", ">100 m"))


def amr_of_A1(A1_aud2, cr=1.0):
    """Equivalent area-to-mass ratio [m^2/kg] of a radial 1/r^2 acceleration."""
    return np.asarray(A1_aud2, float) / (A1_PER_AMR * cr)


def amr_sphere(D_m, rho):
    """Area-to-mass of a homogeneous sphere, m^2/kg."""
    return 1.5 / (rho * np.asarray(D_m, float))


def a2_yark_max(D_m, rho, bond=0.0):
    """Upper bound on the Yarkovsky transverse A2 [au/d^2] (1/r^2 law) for a
    sphere: A2 = (4(1-A)/9) * Phi(1au) * f(Theta) * cos(gamma), with
    Phi = 3 S/(2 D rho c) and max f(Theta) = sqrt(2)/(2+sqrt(2)) = 0.414
    (Farnocchia et al. 2013, Icarus 224, 1; Vokrouhlicky et al. 2000)."""
    phi = 1.5 * P_SRP / (rho * np.asarray(D_m, float))        # m/s^2
    fmax = np.sqrt(2.0) / (2.0 + np.sqrt(2.0))
    return (4.0 * (1.0 - bond) / 9.0) * phi * fmax / AUD2_TO_MS2


# ------------------------------------------------- packed designations (MPC)
_B62 = "0123456789ABCDEFGHIJKLMNOPQRSTUVWXYZabcdefghijklmnopqrstuvwxyz"
_HALF = "ABCDEFGHJKLMNOPQRSTUVWXY"      # half-month letters (no I, no Z)
_SECOND = "ABCDEFGHJKLMNOPQRSTUVWXYZ"   # second letters (no I)


def pack_prov(des):
    """'2009 BD' -> 'K09B00D';  '2016 GW221' -> 'K16GM1W'.  Returns None for
    designations that are not standard provisional ones."""
    try:
        yr, rest = des.strip().split()
        year = int(yr)
        h, s = rest[0], rest[1]
        cyc = int(rest[2:]) if len(rest) > 2 else 0
    except Exception:
        return None
    if h not in _HALF or s not in _SECOND or year < 1800 or year > 2099:
        return None
    if cyc <= 619:
        cent = {18: "I", 19: "J", 20: "K"}[year // 100]
        return "%s%02d%s%s%d%s" % (cent, year % 100, h, _B62[cyc // 10], cyc % 10, s)
    # extended packed form for cycle counts >= 620
    n = cyc * 25 + _SECOND.index(s) - 15500
    out = ""
    for _ in range(4):
        out = _B62[n % 62] + out
        n //= 62
    return "_" + _B62[year % 100] + h + out


def mjd_from_ymd(y, m, d):
    """Julian-calendar-free MJD for Gregorian dates (vector-safe, float day)."""
    y = np.asarray(y); m = np.asarray(m); d = np.asarray(d, float)
    a = (14 - m) // 12
    yy = y + 4800 - a
    mm = m + 12 * a - 3
    jdn = (153 * mm + 2) // 5 + 365 * yy + yy // 4 - yy // 100 + yy // 400 - 32045
    return jdn + d - 2400001.0   # integer-day part d -> MJD at 0h + fraction