TESSPoint/tesswcs discrepancy

main.py

#!/usr/bin/env -S uv run --script
# /// script
# dependencies = [
#   "astropy==7.1.0",
#   "astroquery==0.4.10",
#   "tess-point==0.9.2",
#   "tesswcs==1.6.3",
#   "pandas==2.2.3"
# ]
# ///

import sys
import astropy.units as u
from tesswcs import latest_sector
from tesswcs.locate import get_pixel_locations
from astroquery.mast import Catalogs
from astropy.coordinates import SkyCoord
from tess_stars2px import tess_stars2px_function_entry as tesspoint

tics = [
    43361132,
    152184857,
    423357559,
    140686801,
    309702890,
    152831893,
    139257315,
    290714823,
    209453258,
    408014108,
    197569385,
    439403832,
    32070860,
    402129745,
    156084952,
    31309834,
    237335281,
    304391892,
    139225273,
    123266870,
    197599172,
    197868378,
    238230649,
    393488969,
    38460940, 
    7316458
]

def log(m):
    print(m, file=sys.stderr)

ras, decs = [], []
log("Retrieving coordinates from MAST")
for tic in tics:
    data = Catalogs.query_object(f"TIC {tic}", catalog="TIC", radius=3 * u.arcsec) # type:ignore
    ras.append(data["ra"][0])
    decs.append(data["dec"][0])

tic2coords = {tic: (ras[i], decs[i]) for i, tic in enumerate(tics)}

log("Locating targets with TESSPoint")
outID, outEclipLong, outEclipLat, outSec, outCam, outCcd, \
            outColPix, outRowPix, scinfo = tesspoint(tics, ras, decs)

wrong_sectors = []
wrong_ccds = []
log("Checking results against tesswcs")
for i, (tic, sec, cam, ccd) in enumerate(zip(outID, outSec, outCam, outCcd)):
    if sec > latest_sector:
        continue

    ra, dec = tic2coords[tic]
    tw_res = get_pixel_locations(SkyCoord(ra, dec, unit="deg"), sector=sec).to_pandas()
    if not len(tw_res) > 0:
        wrong_sectors.append((tic, sec, cam, ccd))
    elif tw_res["Camera"][0].item() != cam or tw_res["CCD"][0].item() != ccd:
        actual = tw_res.iloc[0]
        _, _, acam, accd, _, _ = actual
        wrong_ccds.append((tic, sec, cam, ccd, int(acam), int(accd)))

    if i % 5 == 0 or i == len(outID) - 1:
        print(f"\rProgress: {i}/{len(outID)}", end="", flush=True, file=sys.stderr)

def print_row(items):
    print("|" + "|".join([str(item) for item in items]) + "|")

def print_table_header(headers):
    print_row(headers)
    print_row(["-" * len(i) for i in headers])

print("\n## Incorrect Results")
headers = ["TIC", "RA", "Dec", "Sector", "Camera", "CCD"]
print_table_header(headers)
for tic, sec, cam, ccd in wrong_sectors:
    ra, dec = tic2coords[tic]
    print_row([tic, f"{ra:.5f}", f"{dec:.5f}", sec, cam, ccd])

print("\n## Wrong CCD")
headers = ["TIC", "RA", "Dec", "Sector", "Pred. Cam", "Pred. CCD", "Actual Cam", "Actual CCD"]
print_table_header(headers)
for tic, sec, cam, ccd, acam, accd in wrong_ccds:
    ra, dec = tic2coords[tic]
    print_row([tic, f"{ra:.5f}", f"{dec:.5f}", sec, cam, ccd, acam, accd])

results.md

Incorrect Results

TIC	RA	Dec	Sector	Camera	CCD
43361132	9.45108	-27.50636	2	1	2
152184857	308.90991	-30.39838	1	1	4
423357559	14.54043	-5.47595	3	1	3
140686801	72.98291	-77.63662	3	4	2
309702890	79.70681	-65.02435	3	4	1
152831893	345.08576	-39.27944	2	2	4
139257315	347.33058	-44.82874	2	2	4
290714823	84.26332	-53.08898	1	4	4
209453258	46.05911	-34.41177	3	2	2
408014108	3.68331	-11.82457	3	1	4
197569385	325.85889	-37.13034	1	1	1
439403832	18.24109	-15.53774	3	1	1
32070860	57.96999	-67.25703	2	3	3
402129745	11.53561	-30.80440	3	2	4
156084952	11.41810	-50.59947	2	2	1
31309834	88.34894	-66.07577	3	4	1
237335281	347.30704	-62.52071	1	2	1
304391892	352.83435	-25.48495	2	1	1
139225273	346.39947	-42.68354	2	2	4
123266870	28.35180	-34.43606	3	2	1
197599172	327.33022	-40.08492	1	1	1
197868378	55.79128	-58.92492	3	3	3
197868378	55.79128	-58.92492	65	4	4
38460940	61.13322	-62.98722	12	4	4

Wrong CCD

TIC	RA	Dec	Sector	Pred. Cam	Pred. CCD	Actual Cam	Actual CCD
238230649	59.46272	-72.40573	1	4	1	3	4
238230649	59.46272	-72.40573	28	4	2	3	3
393488969	35.70007	-44.55174	3	2	1	3	2
38460940	61.13322	-62.98722	2	4	1	3	4
7316458	34.68630	-44.99545	3	3	2	2	1