Coverage for PyFHD/pyfhd_tools/unit_conv.py: 100%
44 statements
« prev ^ index » next coverage.py v7.9.1, created at 2025-07-01 10:58 +0800
1from astropy.coordinates import SkyCoord, EarthLocation, AltAz, ICRS
2from astropy.wcs import WCS
3from astropy.time import Time
4from astropy import units as u
5import numpy as np
6from numpy.typing import NDArray
9def altaz_to_radec(
10 alt: float,
11 az: float,
12 lat: float,
13 lon: float,
14 height: float,
15 time: float,
16 time_format="jd",
17) -> tuple[float, float]:
18 """
19 Turn AltAz coordinates into the equatorial/celestial coordinates RA and DEC.
20 The exact location and time must given in order for the coordinates to be calculated.
22 Parameters
23 ----------
24 alt : float
25 The altitude in degrees
26 az : float
27 The azimuth in degrees
28 lat : float
29 The latitude of the location (default is MWA's latitude)
30 lon : float
31 The longitude of the location (default is MWA's longitude)
32 height : float
33 The altitude of the location in metres above sea level (default is MWA's altitude above sea level)
34 time : float
35 The time from the UVFITS file
36 time_format : str
37 The time format given, as per AstroPy's Time Object, by default jd (Julian)
39 Returns
40 -------
41 ra : float
42 Right Ascension from the given location and time with altitude and azimuth
43 dec : float
44 Declination from the given location and time with altitude and azimuth
45 """
47 loc = EarthLocation.from_geodetic(lon * u.deg, lat * u.deg, height=height * u.meter)
48 altaz = AltAz(
49 alt=alt * u.deg,
50 az=az * u.deg,
51 location=loc,
52 obstime=Time(time, format=time_format),
53 )
54 return altaz.transform_to(ICRS()).ra.deg, altaz.transform_to(ICRS()).dec.deg
57def radec_to_altaz(
58 ra: float,
59 dec: float,
60 lat: float,
61 lon: float,
62 height: float,
63 time: float,
64 time_format="jd",
65) -> tuple[float, float]:
66 """
67 Turn Celestial/Equatorial coordinates into AltAz at the given location and time.
68 Time Format by default is Julian, but you can use any of the formats provided by the AstroPy Time classes.
70 Parameters
71 ----------
72 ra : float
73 Right Ascension in degrees
74 dec : float
75 Declination in degrees
76 lat : float
77 The latitude in degrees (default is MWA's latitude)
78 lon : float
79 The longtitude in degrees (default is MWA's longitude)
80 height : float
81 The height of the location in metres above sea level (default is MWA's altitude above sea level)
82 time : float
83 The time
84 time_format : str, optional
85 The format of the time, by default 'jd' (Julian)
87 Returns
88 -------
89 alt : float
90 Altitude coordinates for the given location, time and celestial coordinates
91 az : float
92 Azimuth coordinates for the given location, time and celestial coordinates
93 """
95 # Create the Earth Location
96 loc = EarthLocation.from_geodetic(lon * u.deg, lat * u.deg, height=height * u.meter)
97 loc_time = Time(time, format=time_format)
98 radec = SkyCoord(ra=ra * u.deg, dec=dec * u.deg)
99 altaz = radec.transform_to(AltAz(location=loc, obstime=loc_time))
100 return altaz.alt.deg, altaz.az.deg
103def radec_to_pixel(ra: float, dec: float, astr: dict) -> tuple[float, float]:
104 """
105 Turn Celestial Coordinates into Pixel coordinates (X & Y). The astr dictionary should contain
106 cdelt, ctype, crpix and crval as per the WCS standard when naxis is 2.
108 Parameters
109 ----------
110 ra : float
111 The Right Ascension in degrees
112 dec : float
113 The Declination in degrees
114 astr : dict
115 The astrometry dictionary used in PyFHD.
117 Returns
118 -------
119 x : float
120 The pixel x coordinate for the given celestial coordinate.
121 y : float
122 The pixel y coordinate for the given celestial coordinate.
123 """
125 # Create WCS object with astr
126 wcs_astr = WCS(naxis=2)
127 wcs_astr.wcs.cdelt = astr["cdelt"]
128 # AstroPy ctype requires a list of python string objects, NumPy string objects will not work
129 wcs_astr.wcs.ctype = [str(projection) for projection in astr["ctype"]]
130 wcs_astr.wcs.crpix = astr["crpix"]
131 wcs_astr.wcs.crval = astr["crval"]
132 wcs_astr.wcs.equinox = astr["equinox"]
133 wcs_astr.wcs.mjdobs = astr["mjdobs"]
134 wcs_astr.wcs.dateobs = astr["dateobs"]
135 wcs_astr.wcs.radesys = astr["radecsys"]
136 wcs_astr.wcs.set_pv(
137 [
138 (1, 1, astr["pv1"][1]),
139 (1, 2, astr["pv1"][2]),
140 (1, 3, astr["pv1"][3]),
141 (1, 4, astr["pv1"][4]),
142 (2, 1, astr["pv2"][0]),
143 (2, 2, astr["pv2"][1]),
144 ]
145 )
146 # Now use world_to_pixel function for WCS objects
147 x, y = wcs_astr.world_to_pixel(
148 SkyCoord(ra=ra * u.deg, dec=dec * u.deg, frame=astr["radecsys"].lower())
149 )
150 # AstroPy returns values as an array, cast to float
151 if np.size(x) == 1:
152 return float(x), float(y)
153 else:
154 return x, y
157def pixel_to_radec(
158 x: float | NDArray[np.float64], y: float | NDArray[np.float64], astr: dict
159) -> tuple[float | NDArray[np.float64], float | NDArray[np.float64]]:
160 """
161 Turn Pixel coordinates (X & Y) into Celestial Coordinates based off a WCS.
162 The astr dictionary should contain cdelt, ctype, crpix and crval as per the
163 WCS standard when naxis is 2.
165 Parameters
166 ----------
167 x : float | np.ndarray
168 x coordinate(s)
169 y : float | np.ndarray
170 y coordinate(s)
171 astr : dict
172 The astrometry dictionary from an obs dictionary
174 Returns
175 -------
176 ra : float | NDArray[np.float64]
177 Right Ascension from the given pixel coordinates and WCS
178 dec : float | NDArray[np.float64]
179 Declination from the given pixel coordinates and WCS
180 """
181 wcs_astr = WCS(naxis=2)
182 wcs_astr.wcs.cdelt = astr["cdelt"]
183 # AstroPy ctype requires a list of python string objects, NumPy string objects will not work
184 wcs_astr.wcs.ctype = [str(projection) for projection in astr["ctype"]]
185 wcs_astr.wcs.crpix = astr["crpix"]
186 wcs_astr.wcs.crval = astr["crval"]
187 wcs_astr.wcs.equinox = astr["equinox"]
188 wcs_astr.wcs.mjdobs = astr["mjdobs"]
189 wcs_astr.wcs.dateobs = astr["dateobs"]
190 wcs_astr.wcs.radesys = astr["radecsys"]
191 wcs_astr.wcs.set_pv(
192 [
193 (1, 1, astr["pv1"][1]),
194 (1, 2, astr["pv1"][2]),
195 (1, 3, astr["pv1"][3]),
196 (1, 4, astr["pv1"][4]),
197 (2, 1, astr["pv2"][0]),
198 (2, 2, astr["pv2"][1]),
199 ]
200 )
201 radec = wcs_astr.pixel_to_world(x, y)
203 return radec.ra, radec.dec