Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Dcp uvfits rdate fix #1436

Closed
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Jump to
Jump to file
Failed to load files.
Diff view
Diff view
40 changes: 24 additions & 16 deletions pyuvdata/uvdata/uvfits.py
Original file line number Diff line number Diff line change
Expand Up @@ -11,6 +11,7 @@
from astropy import constants as const
from astropy.io import fits
from astropy.time import Time
from astropy.coordinates import EarthLocation
from docstring_parser import DocstringStyle

from .. import utils as uvutils
Expand Down Expand Up @@ -59,9 +60,11 @@ def _get_parameter_data(
# angles in uvfits files are stored in degrees, so convert to radians
self.lst_array = np.deg2rad(vis_hdu.data.par("lst"))
if run_check_acceptability:
(latitude, longitude, altitude) = (
self.telescope_location_lat_lon_alt_degrees
)
(
latitude,
longitude,
altitude,
) = self.telescope_location_lat_lon_alt_degrees
uvutils.check_lsts_against_times(
jd_array=self.time_array,
lst_array=self.lst_array,
Expand Down Expand Up @@ -1158,8 +1161,18 @@ def write_uvfits(
hdu.header["PSCAL" + str(i + 1) + " "] = pscal_dict[key]
hdu.header["PZERO" + str(i + 1) + " "] = pzero_dict[key]

# ISO string of first time in self.time_array
hdu.header["DATE-OBS"] = Time(self.time_array[0], scale="utc", format="jd").isot
# Create an astropy.time.Time object from first timestamp.
# This is used to generate DATE-OBS and RDATE keywords
eloc = EarthLocation(
self.telescope_location[0],
self.telescope_location[1],
self.telescope_location[2],
unit="m",
)
obs_date0 = Time(self.time_array[0], format="jd", scale="utc", location=eloc)

# Per AIPS memo 117, DATE-OBS is the YYYY-MM-DD string
hdu.header["DATE-OBS"] = obs_date0.strftime("%Y-%m-%d")

hdu.header["CTYPE2 "] = "COMPLEX "
hdu.header["CRVAL2 "] = 1.0
Expand Down Expand Up @@ -1338,26 +1351,21 @@ def write_uvfits(
else:
ant_hdu.header["FREQ"] = self.freq_array[0, 0]

# RDATE: Reference date when obs started, YYYY-MM-DD
if (self.rdate is None) or (self.rdate == ""):
rdate_obj = Time(np.floor(self.time_array[0]), format="jd", scale="utc")
else:
try:
rdate_obj = Time(self.rdate, scale="utc")
except ValueError:
rdate_obj = Time(np.floor(self.time_array[0]), format="jd", scale="utc")

if self.rdate is None:
ant_hdu.header["RDATE"] = rdate_obj.strftime("%Y-%m-%d")
ant_hdu.header["RDATE"] = obs_date0.strftime("%Y-%m-%d")
else:
ant_hdu.header["RDATE"] = self.rdate

# GSTIA0: Greenwich sidereal time in degrees at zero hours on RDATE
if self.gst0 is None:
ant_hdu.header["GSTIA0"] = rdate_obj.sidereal_time("apparent", "tio").deg
ant_hdu.header["GSTIA0"] = obs_date0.sidereal_time("apparent", "tio").deg
else:
ant_hdu.header["GSTIA0"] = self.gst0

# UT1UTC: the difference between UT1 and UTC in seconds on RDATE
if self.dut1 is None:
ant_hdu.header["UT1UTC"] = float(rdate_obj.delta_ut1_utc)
ant_hdu.header["UT1UTC"] = float(obs_date0.delta_ut1_utc)
else:
ant_hdu.header["UT1UTC"] = self.dut1

Expand Down