88
99from astropy .time import Time
1010from astropy .coordinates import AltAz , SkyCoord
11+ from astropy .coordinates .erfa_astrom import erfa_astrom , ErfaAstromInterpolator
1112import astropy .units as u
1213
1314from fact .io import read_h5py , to_h5py
3132from ..preprocessing import calc_true_disp
3233from ..logging import setup_logging
3334
35+ # use interpolation for all coordinate transforms
36+ # 100x speed increase with no precision lost (~uas)
37+ erfa_astrom .set (ErfaAstromInterpolator (5 * u .min ))
38+
3439
3540dl3_columns = [
3641 'run_id' ,
4651 'theta_deg_off_5' ,
4752 'pointing_position_az' ,
4853 'pointing_position_zd' ,
54+ 'az_prediction' ,
55+ 'zd_prediction' ,
4956]
5057dl3_columns_sim_read = [
5158 'corsika_run_header_run_number' ,
@@ -92,10 +99,14 @@ def to_altaz(obstime, source):
9299
93100
94101def concat_results_altaz (results ):
95- obstime = np .concatenate ([s .obstime for s in results ])
102+ jd1 = np .concatenate ([s .obstime .jd1 for s in results ])
103+ jd2 = np .concatenate ([s .obstime .jd2 for s in results ])
104+ obstime = Time (jd1 , jd2 , format = 'jd' , copy = False )
105+
106+ alt = u .Quantity (np .concatenate ([s .alt .deg for s in results ]), u .deg , copy = False )
107+ az = u .Quantity (np .concatenate ([s .az .deg for s in results ]), u .deg , copy = False )
96108 return SkyCoord (
97- alt = np .concatenate ([s .alt .deg for s in results ]) * u .deg ,
98- az = np .concatenate ([s .az .deg for s in results ]) * u .deg ,
109+ alt = alt , az = az ,
99110 frame = AltAz (location = LOCATION , obstime = obstime )
100111 )
101112
@@ -109,6 +120,12 @@ def calc_source_features_common(
109120 pointing_position_az ,
110121):
111122 result = {}
123+
124+ k_zd , k_az = 'zd_prediction' , 'az_prediction'
125+ result [k_zd ], result [k_az ] = camera_to_horizontal (
126+ prediction_x , prediction_y ,
127+ pointing_position_zd , pointing_position_az ,
128+ )
112129 result ['theta_deg' ] = calc_theta_camera (
113130 prediction_x ,
114131 prediction_y ,
0 commit comments