How would I compute the offsets between tle and ephemeris RA and Dec in the cross track and along track directions?

Question

I have computed the offsets in ra and Dec between tle predicted values and those from an ephemeris file for a particular gps satellite , however to fully scrutinise these offsets , I’d like to see how they look when transformed in the cross track and along track directions. Although I’m not sure how to go about this mathematically. An explanation on how to do this and the transformations themselves would be much appreciated.

I’m not really sure where to start. I have those code snippet but I would need more explanation/context to understand whether I’m going the right way to get what I desire.

from astropy.time import Time

# Initial start time and final end time
initial_start_time = Time('2023-12-10T03:31:51.488')
final_end_time = Time('2023-12-10T05:59:19.634')

# Number of sets
num_sets = 73
start_times = np.linspace(initial_start_time.jd, final_end_time.jd, num_sets + 1)[:-1]
mid_times = start_times + (final_end_time.jd - initial_start_time.jd) / num_sets / 2
end_times = np.linspace(initial_start_time.jd, final_end_time.jd, num_sets + 1)[1:]

start_times = Time(start_times, format='jd')
mid_times = Time(mid_times, format='jd')
end_times = Time(end_times, format='jd')

loader = Loader('/var/tmp/')
timescale = loader.timescale()
de421 = loader('de421.bsp')

sat = Satrec()
line1 = '1 28190U 04009A  23343.04034598 -.00000057 00000-0 00000-0 0 9994'
line2 = '2 28190 55.6161 312.6614 0092109 141.9301 239.7489 2.00564352144478'
satellite = EarthSatellite(line1,line2)
earth = de421['earth']
sun = de421['sun']
site_topos = wgs84.latlon(28.7606, -17.8795, elevation_m=2350)

unrotated_2 = []
ROT = []
for (start_time,midtime,end_time) in zip(start_times,mid_times,end_times):
  ra, dec, _ = (satellite - site_topos).at(ts.from_astropy(Time([start_time, end_time, midtime]))).radec()
  dra = (http://ra.to(u.deg)[1] - http://ra.to(u.deg)[0]).value * np.cos(http://dec.to(u.rad)[2].value)
  ddec = (http://dec.to(u.deg)[1] - http://dec.to(u.deg)[0]).value
  aperture_length = np.sqrt(dra**2 + ddec**2) * 3600 / 0.98

  angle = -np.degrees(np.arctan2(ddec, dra))
  c, s = np.cos(np.radians(angle)), np.sin(np.radians(angle))
  rot_matrix = np.array([[c, s], [-s, c]])
  unrotated = rot_matrix @ np.asarray([dra,ddec])