Fix merge conflicts

Signed-off-by: Cees Bassa <bassa@astron.nl>

bstplot.py

@@ -8,7 +8,7 @@ import matplotlib.pyplot as plt
 import matplotlib.dates as mdates
 from matplotlib.ticker import FormatStrFormatter, MultipleLocator
 
-from lofty.io import read_bst
+from lofty.io import read_bst, read_minio_bst
 from lofty.plotting import plot_bst_timeseries, plot_bst_dynamic_spectrum
 
 import astropy.units as u
@@ -17,6 +17,7 @@ from astropy.time import Time
 if __name__ == "__main__":
     # Read command line arguments
     parser = argparse.ArgumentParser(description="Plot BST data")
+    parser.add_argument("-m", "--minio", help="Minio HDF5", action="store_true")
     parser.add_argument("-d", "--dynspec", help="Plot dynamic spectrum", action="store_true")
     parser.add_argument("-t", "--timeseries", help="Plot timeseries", action="store_true")
     parser.add_argument("-s", "--subband", help="Subband to plot timeseries [int, default: 300]", default=300, type=int)
@@ -25,9 +26,10 @@ if __name__ == "__main__":
     args = parser.parse_args()
 
     # Read data
+    if args.minio:
+        metadata, data = read_minio_bst(args.filenames)
+    else:
         metadata, data = read_bst(args.filenames)
-    metadata["frequency_bands"] = [["LBA_10_90", "LBA_10_90"]]
-    metadata["subbands"] = range(12, 500)
 
     print(data.shape)
     print(metadata)

lofty/io.py

@@ -197,6 +197,98 @@ def read_xst(xst_filenames, fill_empty=False):
     return metadata, data
 
 
+def read_single_minio_xst(xst_filename, fill_empty=False):
+    # Open HDF5 file
+    h5 = h5py.File(xst_filename, "r")
+
+    # Store main header meta data
+    metadata = {}
+    for key in ["antenna_names", "antenna_status", "antenna_type", "antenna_usage_mask", "station_name"]:
+        metadata[key] = h5["/"].attrs[key]
+
+    # Antennafield
+    metadata["antennafield"] = h5["/"].attrs["antennafield_device"].upper().replace("STAT/AFH/", "").replace("STAT/AFL/", "")
+    
+    # Number of subintegrations
+    nsub = len(h5.keys())
+    metadata["nsub"] = nsub
+
+    # Loop over keys
+    data_dict = []
+    for key in h5.keys():
+        # Keys in group
+        attribute_keys = list(h5[key].attrs.keys())
+
+        # Keys to store in dict
+        keys_to_store = []
+        
+        # Create dict
+        d = {"data": np.array(h5[key]["data"]).squeeze()}
+        for key_to_store in keys_to_store:
+            if key_to_store in attribute_keys:
+                d[key_to_store] = h5[key].attrs[key_to_store]
+            else:
+                print(f"Failed to find {key_to_store} in {key}")
+                d[key_to_store] = None
+
+        # Add fake data
+        if "LBA" in metadata["antennafield"]:
+            d["frequency_band"] = [["LBA_10_90"] * 2] * 96
+        else:
+            d["frequency_band"] = [["HBA_110_190"] * 2] * 48
+        d["subbands"] = h5[key]["subband"].squeeze()
+        d["integration_interval"] = 1.0
+
+        # Timestamp
+        d["timestamp"] = key.replace("XST_", "")
+        
+        # Store dict
+        data_dict.append(d)
+
+    # Close file
+    h5.close()
+
+    # Store arrays
+    metadata["timestamps"] = np.array([d["timestamp"] for d in data_dict])
+    metadata["durations"] = np.array([d["integration_interval"] for d in data_dict])
+    metadata["subbands"] = np.array([d["subbands"] for d in data_dict])
+    
+    # Assume first entry holds for all (ASSUMPTION!!)
+    metadata["frequency_bands"] = data_dict[0]["frequency_band"]
+    if "HBA" in metadata["antennafield"]:
+        metadata["tile_beam_pointing_direction"] = data_dict[0]["tile_beam_pointing_direction"][0]
+    else:
+        metadata["tile_beam_pointing_direction"] = ""
+
+    # Data shape
+    nrcu, _ = data_dict[0]["data"].shape
+    nant = nrcu // 2
+    metadata["nrcu"] = nrcu
+    metadata["nant"] = nant
+    
+    # Store as nsub x nrcu x nrcu
+    data = np.zeros((nsub, nrcu, nrcu), dtype=data_dict[0]["data"].dtype)
+    for isub in range(nsub):
+        data[isub] = data_dict[isub]["data"]
+
+    # Fill empty parts
+    if fill_empty:
+        data += np.conj(np.transpose(data, (0, 2, 1))) * (data==0)
+            
+    return metadata, data
+    
+def read_minio_xst(xst_filenames, fill_empty=False):
+    first = True
+    for xst_filename in xst_filenames:
+        # Get data
+        if first:
+            metadata, data = read_single_minio_xst(xst_filename, fill_empty)
+            first = False
+
+    metadata["station_name"] = metadata["station_name"].upper()
+    return metadata, data
+
+
 def read_single_sst(sst_filename):
     # Open HDF5 file
     h5 = h5py.File(sst_filename, "r")
@@ -292,6 +384,95 @@ def read_sst(sst_filenames):
 
     return metadata, data
 
+def read_single_minio_sst(sst_filename):
+    # Open HDF5 file
+    h5 = h5py.File(sst_filename, "r")
+
+    # Store main header meta data
+    metadata = {}
+    for key in ["antenna_names", "antenna_status", "antenna_type", "antenna_usage_mask", "station_name"]:
+        metadata[key] = h5["/"].attrs[key]
+
+    # Antennafield
+    metadata["antennafield"] = h5["/"].attrs["antennafield_device"].upper().replace("STAT/AFH/", "").replace("STAT/AFL/", "")
+    
+    # Number of subintegrations
+    nsub = len(h5.keys())
+    metadata["nsub"] = nsub
+
+    # Loop over keys
+    data_dict = []
+    for key in h5.keys():
+        # Keys in group
+        attribute_keys = list(h5[key].attrs.keys())
+
+        # Keys to store in dict
+        keys_to_store = []
+        
+        # Create dict
+        d = {"data": np.array(h5[key])}
+        for key_to_store in keys_to_store:
+            if key_to_store in attribute_keys:
+                d[key_to_store] = h5[key].attrs[key_to_store]
+            else:
+                print(f"Failed to find {key_to_store} in {key}")
+                d[key_to_store] = None
+
+        # Add fake data
+        if "LBA" in metadata["antennafield"]:
+            d["frequency_band"] = [["LBA_10_90"] * 2] * 96
+        else:
+            d["frequency_band"] = [["HBA_110_190"] * 2] * 48
+        d["subbands"] = list(range(0, 512))
+        d["integration_interval"] = 1.0
+
+        # Timestamp
+        d["timestamp"] = key.replace("SST_", "")
+        
+        # Store dict
+        data_dict.append(d)
+
+    # Close file
+    h5.close()
+
+    # Store arrays
+    metadata["timestamps"] = np.array([d["timestamp"] for d in data_dict])
+    metadata["durations"] = np.array([d["integration_interval"] for d in data_dict])
+    
+    # Assume first entry holds for all (ASSUMPTION!!)
+    metadata["subbands"] = data_dict[0]["subbands"]
+    metadata["frequency_bands"] = data_dict[0]["frequency_band"]
+    if "HBA" in metadata["antennafield"]:
+        metadata["tile_beam_pointing_direction"] = data_dict[0]["tile_beam_pointing_direction"][0]
+    else:
+        metadata["tile_beam_pointing_direction"] = ""
+
+    # Data shape
+    nrcu, nchan = data_dict[0]["data"].shape
+    nant = nrcu // 2
+    metadata["nrcu"] = nrcu
+    metadata["nant"] = nant
+    metadata["nchan"] = nchan
+    
+    # Store as nsub x npol x nant x nchan
+    data = np.zeros((nsub, 2, nant, nchan))
+    for isub in range(nsub):
+        data[isub, 0] = data_dict[isub]["data"][0::2, :]
+        data[isub, 1] = data_dict[isub]["data"][1::2, :]
+    
+    return metadata, data
+    
+def read_minio_sst(sst_filenames):
+    first = True
+    for sst_filename in sst_filenames:
+        # Get data
+        if first:
+            metadata, data = read_single_minio_sst(sst_filename)
+            first = False
+
+    return metadata, data
+
+
 def read_single_bst(bst_filename):
     # Open HDF5 file
     h5 = h5py.File(bst_filename, "r")
@@ -362,6 +543,7 @@ def read_single_bst(bst_filename):
 
     return metadata, data
 
+    
 def read_bst(bst_filenames):
     first = True
     for bst_filename in bst_filenames:
@@ -385,6 +567,92 @@ def read_bst(bst_filenames):
 
     return metadata, data
 
+def read_single_minio_bst(bst_filename):
+    # Open HDF5 file
+    h5 = h5py.File(bst_filename, "r")
+
+    # Store main header meta data
+    metadata = {}
+    for key in ["antenna_names", "antenna_status", "antenna_type", "antenna_usage_mask", "station_name"]:
+        metadata[key] = h5["/"].attrs[key]
+
+    # Antennafield
+    metadata["antennafield"] = h5["/"].attrs["antennafield_device"].upper().replace("STAT/AFH/", "").replace("STAT/AFL/", "")
+    
+    # Number of subintegrations
+    nsub = len(h5.keys())
+    metadata["nsub"] = nsub
+
+    # Loop over keys
+    data_dict = []
+    for key in h5.keys():
+        # Keys in group
+        attribute_keys = list(h5[key].attrs.keys())
+
+        # Keys to store in dict
+        keys_to_store = []
+        
+        # Create dict
+        d = {"data": np.array(h5[key])}
+        for key_to_store in keys_to_store:
+            if key_to_store in attribute_keys:
+                d[key_to_store] = h5[key].attrs[key_to_store]
+            else:
+                print(f"Failed to find {key_to_store} in {key}")
+                d[key_to_store] = None
+
+        # Add fake data
+        if "LBA" in metadata["antennafield"]:
+            d["frequency_band"] = [["LBA_10_90"] * 2] * 96
+        else:
+            d["frequency_band"] = [["HBA_110_190"] * 2] * 48
+        d["subbands"] = list(range(12, 500))
+        d["integration_interval"] = 1.0
+
+        # Timestamp
+        d["timestamp"] = key.replace("BST_", "")
+        
+        # Store dict
+        data_dict.append(d)
+
+    # Close file
+    h5.close()
+
+    # Store arrays
+    metadata["timestamps"] = np.array([d["timestamp"] for d in data_dict])
+    metadata["durations"] = np.array([d["integration_interval"] for d in data_dict])
+    
+    # Assume first entry holds for all (ASSUMPTION!!)
+    metadata["subbands"] = data_dict[0]["subbands"]
+    metadata["frequency_bands"] = data_dict[0]["frequency_band"]
+    if "HBA" in metadata["antennafield"]:
+        metadata["tile_beam_pointing_direction"] = data_dict[0]["tile_beam_pointing_direction"][0]
+    else:
+        metadata["tile_beam_pointing_direction"] = ""
+    
+    # Data shape
+    nchan, npol = data_dict[0]["data"].shape
+    metadata["nchan"] = nchan
+    
+    # Store as nsub x npol x nchan
+    data = np.zeros((nsub, 2, nchan))
+    for isub in range(nsub):
+        data[isub, 0] = data_dict[isub]["data"][:, 0::2].T
+        data[isub, 1] = data_dict[isub]["data"][:, 1::2].T
+    
+    return metadata, data
+    
+def read_minio_bst(bst_filenames):
+    first = True
+    for bst_filename in bst_filenames:
+        # Get data
+        if first:
+            metadata, data = read_single_minio_bst(bst_filename)
+            first = False
+
+    return metadata, data
+
+
 def read_acm_cube(xst_filename, fill=True):
     # Open file
     h5 = h5py.File(xst_filename, "r")

lofty/plotting.py

@@ -152,7 +152,7 @@ def plot_sst_timeseries(metadata, data, subband):
     fig.autofmt_xdate(rotation=0, ha="center")
 
     # Select subband
-    c = metadata["subbands"] == subband
+    c = np.array(metadata["subbands"]) == subband
     px = decibel(data[:, 0, :, c].squeeze())
     py = decibel(data[:, 1, :, c].squeeze())
 
@@ -236,6 +236,7 @@ def plot_sst_bandpass(metadata, data):
 
     px = decibel(data[0, 0, :, :].squeeze())
     py = decibel(data[0, 1, :, :].squeeze())
+
     nant = metadata["nant"]
     mask = metadata["antenna_usage_mask"][:nant] # Array is longer by default
     for iant in range(nant):
@@ -400,6 +401,7 @@ def sky_plot(img, title, subtitle, source_names, source_lmn, fig, vmin=None, vma
     ax.add_artist(circle1)
 
     # Plot sources
+    if source_lmn is not None:
         for name, lmn in zip(source_names, source_lmn):
             if lmn[2] >= -1:
                 ax.plot(lmn[0], lmn[1], marker="x", color="k")

sstplot.py

@@ -8,12 +8,13 @@ import matplotlib.pyplot as plt
 import matplotlib.dates as mdates
 from matplotlib.ticker import FormatStrFormatter, MultipleLocator
 
-from lofty.io import read_sst
+from lofty.io import read_sst, read_minio_sst
 from lofty.plotting import plot_sst_timeseries, plot_sst_bandpass, plot_sst_dynamic_spectrum
 
 if __name__ == "__main__":
     # Read command line arguments
     parser = argparse.ArgumentParser(description="Plot SST data")
+    parser.add_argument("-m", "--minio", help="Minio HDF5", action="store_true")
     parser.add_argument("-d", "--dynspec", help="Plot dynamic spectrum", action="store_true")
     parser.add_argument("-t", "--timeseries", help="Plot timeseries", action="store_true")
     parser.add_argument("-b", "--bandpass", help="Plot bandpass", action="store_true")    
@@ -27,6 +28,9 @@ if __name__ == "__main__":
     print(args.pol)
 
     # Read data
+    if args.minio:
+        metadata, data = read_minio_sst(args.filenames)
+    else:
         metadata, data = read_sst(args.filenames)        
 
     print(metadata)

xstimager.py

@@ -9,8 +9,7 @@ import matplotlib.pyplot as plt
 
 from lofarantpos.db import LofarAntennaDatabase
 
-from lofty.io import read_npz
-from lofty.io import read_xst
+from lofty.io import read_xst, read_minio_xst
 from lofty.imaging import sky_imager
 from lofty.sources import lmn_sources
 from lofty.plotting import sky_plot
@@ -19,12 +18,15 @@ from lofty.utils import get_station_type, get_full_station_name, freq_from_sb, g
 if __name__ == "__main__":
     # Read command line arguments
     parser = argparse.ArgumentParser(description="Create XST images")
+    parser.add_argument("-m", "--minio", help="Minio HDF5", action="store_true")
     parser.add_argument("filenames", help="XST filenames", nargs="*", metavar="FILE")
     args = parser.parse_args()
 
     # Read data
+    if args.minio:
+        metadata, data = read_minio_xst(sorted(args.filenames), fill_empty=True)
+    else:
         metadata, data = read_xst(sorted(args.filenames), fill_empty=True)
-    #metadata, data = read_npz(sorted(args.filenames), fill_empty=True)
 
     # Get information
     station_type = get_station_type(metadata["station_name"].upper())

xstplot.py

@@ -8,13 +8,13 @@ import matplotlib.pyplot as plt
 import matplotlib.dates as mdates
 from matplotlib.ticker import FormatStrFormatter, MultipleLocator
 
-from lofty.io import read_xst
-from lofty.io import read_npz
+from lofty.io import read_xst, read_minio_xst
 from lofty.plotting import plot_xst_matrix
 
 if __name__ == "__main__":
     # Read command line arguments
     parser = argparse.ArgumentParser(description="Plot XST data")
+    parser.add_argument("-m", "--minio", help="Minio HDF5", action="store_true")
     parser.add_argument("filenames", help="XST filenames", nargs="*", metavar="FILE")
     args = parser.parse_args()
 
@@ -26,5 +26,8 @@ if __name__ == "__main__":
     for key in metadata.keys():
         print(key, metadata[key])
 
+    print(metadata["subbands"])
+    print(metadata["nsub"])
+    
     for isub in range(metadata["nsub"]):
         plot_xst_matrix(metadata, data, isub)