From a5debc7329f4c63222b4bdab6a542b561c90fe34 Mon Sep 17 00:00:00 2001 From: "Adam Ginsburg (keflavich)" Date: Mon, 20 Dec 2021 14:32:21 -0500 Subject: [PATCH] add HISTORY parsing for CASA cubes --- spectral_cube/base_class.py | 4 ++++ spectral_cube/io/casa_image.py | 44 ++++++++++++++++++++++++++++------ spectral_cube/spectral_cube.py | 6 ++++- 3 files changed, 46 insertions(+), 8 deletions(-) diff --git a/spectral_cube/base_class.py b/spectral_cube/base_class.py index 566b22604..f8d0b0a34 100644 --- a/spectral_cube/base_class.py +++ b/spectral_cube/base_class.py @@ -107,6 +107,10 @@ def header(self): header.insert(3+ind, Card(keyword='NAXIS{0:1d}'.format(ind+1), value=sh)) + if self._history is not None: + for row in self._history: + header.add_history(row) + return header def check_jybeam_smoothing(self, raise_error_jybm=True): diff --git a/spectral_cube/io/casa_image.py b/spectral_cube/io/casa_image.py index 2e0fce9b5..40f6008ea 100644 --- a/spectral_cube/io/casa_image.py +++ b/spectral_cube/io/casa_image.py @@ -36,6 +36,30 @@ def is_casa_image(origin, filepath, fileobj, *args, **kwargs): return filepath is not None and filepath.lower().endswith('.image') +def make_fits_history(logtable): + """ + Convert a CASA logtable to a FITS history + + https://github.com/casacore/casacore/blob/dbf28794ef446bbf4e6150653dbe404379a3c429/fits/FITS/FITSHistoryUtil.cc + """ + from astropy.time import Time + + history_rows = [] + for row in logtable: + head_row = f"{Time(row['TIME']*u.s, format='mjd').isot} {row['PRIORITY']} SRCCODE::{row['LOCATION']}" + msg = row['MESSAGE'] + if len(msg) <= 72: + datarows = [msg] + else: + datarows = [] + datarows.append(msg[:72]) + for ii in range(72, len(msg), 71): + datarows.append('>' + msg[ii:ii+71]) + history_rows.append(head_row) + history_rows.extend(datarows) + return history_rows + + def load_casa_image(filename, skipdata=False, memmap=True, skipvalid=False, skipcs=False, target_cls=None, use_dask=None, target_chunksize=None, **kwargs): @@ -75,6 +99,9 @@ def load_casa_image(filename, skipdata=False, memmap=True, desc = getdesc(filename) + logtable = Table.read(filename + "/logtable") + history = make_fits_history(logtable) + casa_cs = desc['_keywords_']['coords'] if 'units' in desc['_keywords_']: @@ -168,12 +195,14 @@ def load_casa_image(filename, skipdata=False, memmap=True, mask = None if 'beam' in locals(): - cube = DaskSpectralCube(data, wcs_slice, mask, meta=meta, beam=beam) + cube = DaskSpectralCube(data, wcs_slice, mask, meta=meta, beam=beam, history=history) elif 'beams' in locals(): - cube = DaskVaryingResolutionSpectralCube(data, wcs_slice, mask, meta=meta, - beams=beams['I']) + cube = DaskVaryingResolutionSpectralCube(data, wcs_slice, mask, + meta=meta, + beams=beams['I'], + history=history) else: - cube = DaskSpectralCube(data, wcs_slice, mask, meta=meta) + cube = DaskSpectralCube(data, wcs_slice, mask, meta=meta, history=history) # with #592, this is no longer true # we've already loaded the cube into memory because of CASA # limitations, so there's no reason to disallow operations @@ -196,16 +225,17 @@ def load_casa_image(filename, skipdata=False, memmap=True, if 'beam' in locals(): data[component] = DaskSpectralCube(data_, wcs_slice, mask[component], - meta=meta, beam=beam) + meta=meta, beam=beam, history=history) elif 'beams' in locals(): data[component] = DaskVaryingResolutionSpectralCube(data_, wcs_slice, mask[component], meta=meta, - beams=beams[component]) + beams=beams[component], + history=history) else: data[component] = DaskSpectralCube(data_, wcs_slice, mask[component], - meta=meta) + meta=meta, history=history) data[component].allow_huge_operations = True diff --git a/spectral_cube/spectral_cube.py b/spectral_cube/spectral_cube.py index ccca59f3b..10bf5994f 100644 --- a/spectral_cube/spectral_cube.py +++ b/spectral_cube/spectral_cube.py @@ -171,11 +171,15 @@ class BaseSpectralCube(BaseNDClass, MaskableArrayMixinClass, HeaderMixinClass): def __init__(self, data, wcs, mask=None, meta=None, fill_value=np.nan, - header=None, allow_huge_operations=False, wcs_tolerance=0.0): + header=None, allow_huge_operations=False, wcs_tolerance=0.0, + history=None): # Deal with metadata first because it can affect data reading self._meta = meta or {} + # add the "history", which needs to be read specially for some cubes + self._history = history + # must extract unit from data before stripping it if 'BUNIT' in self._meta: self._unit = cube_utils.convert_bunit(self._meta["BUNIT"])