#!/usr/bin/env python
"""Chase the on-plane IRAS gate-passers lacking a WISE-HII distance: cross-match
to Urquhart+2018 ATLASGAL distance catalogue (J/MNRAS/473/1059/table5, kinematic
dist + VLSR + Lbol). Report any genuine <3 kpc & >1e4 Lsun not already in sample."""
import numpy as np, csv
from astropy.table import Table
from astropy.coordinates import SkyCoord
import astropy.units as u
from astroquery.vizier import Vizier

miss = Table.read('onplane_iras_notRMS.fits')
mc = SkyCoord(miss['_ra'], miss['_dec'], unit=u.deg)

# which lacked a WISE-HII distance (reproduce the 120" match w/ valid Dist)
Vizier.ROW_LIMIT = -1
wise = Vizier(columns=['GLON', 'GLAT', 'Dist'], row_limit=-1).query_constraints(
    catalog='J/ApJS/212/1/wisecat')[0]
wc = SkyCoord(wise['GLON'], wise['GLAT'], unit=u.deg, frame='galactic').icrs
iw, dw, _ = mc.match_to_catalog_sky(wc)
has_wise = np.array([(dw[i].arcsec < 120) and np.isfinite(np.nan_to_num(float(wise['Dist'][iw[i]])))
                     and float(wise['Dist'][iw[i]]) > 0 for i in range(len(miss))])
nod = miss[~has_wise]
print(f'on-plane not-in-RMS lacking WISE distance: {len(nod)}')

# Urquhart 2018 ATLASGAL table8 (GLON, GLAT, Dist)
t5 = Vizier(columns=['*'], row_limit=-1).get_catalogs('J/MNRAS/473/1059/table8')[0]
gl = [c for c in t5.colnames if c.lower() in ('glon', 'glong', '_glon')][0]
gb = [c for c in t5.colnames if c.lower() in ('glat', '_glat')][0]
print('Urquhart table8 cols:', t5.colnames)
uc = SkyCoord(t5[gl], t5[gb], unit=u.deg, frame='galactic').icrs
nc = SkyCoord(nod['_ra'], nod['_dec'], unit=u.deg)
iu, du, _ = nc.match_to_catalog_sky(uc)

# sample positions (for "already in sample" check)
samp = []
for o in csv.DictReader(open('master_distances.csv')):
    samp.append((float(o['GLON']), float(o['GLAT'])))
for o in csv.DictReader(open('candidate_distances_ext.csv')):
    samp.append((float(o['GLON']), float(o['GLAT'])))
sl, sb = map(np.array, zip(*samp))
scs = SkyCoord(sl * u.deg, sb * u.deg, frame='galactic').icrs
isx, dsx, _ = nc.match_to_catalog_sky(scs)

dcol = 'Dist'
nclose = 0; nmatch = 0; nofar = 0
print(f"\n{'IRAS':14s}{'l':>7s}{'b':>6s}{'L@3kpc':>9s} {'AGAL_d':>7s}  verdict")
for j, row in enumerate(nod):
    insamp = dsx[j].arcsec < 180
    d = None
    if du[j].arcsec < 60:
        try:
            d = float(t5[dcol][iu[j]])
        except (ValueError, TypeError):
            d = None
    if d is not None and np.isfinite(d) and d > 0:
        nmatch += 1
        L = row['L_at3kpc'] * (d / 3.0) ** 2
        if d < 3 and L >= 1e4 and not insamp:
            nclose += 1
            v = f'*** <3kpc&>1e4 NOT in sample: d={d:.1f} L={L:.1e}'
        elif d >= 3:
            nofar += 1; v = f'AGAL d={d:.1f} far (L={L:.1e})'
        elif insamp:
            v = f'AGAL d={d:.1f} but in sample'
        else:
            v = f'AGAL d={d:.1f} -> L={L:.1e} (<1e4)'
        ds = f'{d:.1f}'
    else:
        v = 'still no distance'; ds = '--'
    if '***' in v or 'far' not in v:  # print interesting only (skip the many far)
        print(f"  IRAS{row['IRAS']:12s}{row['GLON']:7.2f}{row['GLAT']:+6.2f}{row['L_at3kpc']:9.1e} {ds:>7s}  {v}")
print(f'\nAGAL distance recovered: {nmatch}/{len(nod)} | far: {nofar} | '
      f'still nodist: {len(nod)-nmatch} | *** new <3kpc&>1e4: {nclose}')