import os
os.environ['HDF5_USE_FILE_LOCKING'] = 'FALSE'
import h5py
import hdf5plugin
import matplotlib
import matplotlib.pyplot as plt
from matplotlib.ticker import FormatStrFormatter
import matplotlib.patches as mpatches
import matplotlib.gridspec as gridspec
import numpy as np
from pyuvdata import UVCal, UVData
import sys
import glob
import uvtools as uvt
from astropy.time import Time
from astropy.coordinates import EarthLocation, AltAz, Angle
from astropy.coordinates import SkyCoord as sc
import pandas
import warnings
import copy
from hera_notebook_templates import utils
import hera_qm
from hera_mc import cm_hookup
import importlib
from scipy import stats
from IPython.display import display, HTML
#warnings.filterwarnings('ignore')
%matplotlib inline
%config InlineBackend.figure_format = 'retina'
display(HTML("<style>.container { width:100% !important; }</style>"))
#get data location
data_path = os.environ['DATA_PATH']
print(f'DATA_PATH = "{data_path}"')
statuses = os.environ['APRIORI_STATUSES']
print(f'APRIORI_STATUSES = {statuses}')
JD = os.environ['JULIANDATE']
print(f'JULIANDATE = {JD}')
utc = Time(JD, format='jd').datetime
print(f'Date = {utc.month}-{utc.day}-{utc.year}')
DATA_PATH = "/mnt/sn1/data2/2460902" APRIORI_STATUSES = dish_maintenance,dish_ok,RF_maintenance,RF_ok,digital_ok,digital_maintenance,calibration_maintenance,calibration_triage,calibration_ok JULIANDATE = 2460902 Date = 8-14-2025
# Load in data
HHfiles, difffiles, HHautos, diffautos, uvdx, uvdy = utils.load_data(data_path,JD)
uvd = UVData()
unread = True
readInd=0
while unread and readInd<len(HHautos):
try:
uvd.read(HHautos[readInd])
unread = False
except:
readInd += 1
continue
use_ants = utils.get_use_ants(uvd,statuses,JD)
print(f'This day contains {len(use_ants)} antennas of the given status category.')
uvd.read(HHautos[::10], skip_bad_files=True, antenna_nums = use_ants)
lsts = uvd.lst_array
uvdx.select(antenna_nums=use_ants)
uvdy.select(antenna_nums=use_ants)
377 sum files found between JDs 2460902.21076 and 2460902.44867 377 diff files found between JDs 2460902.21076 and 2460902.44867 377 sum auto files found between JDs 2460902.21076 and 2460902.44867 377 diff auto files found between JDs 2460902.21076 and 2460902.44867
This day contains 298 antennas of the given status category.
Sky Coverage Map¶
Map of the sky (made using the Haslam 408MHz map). The RA/DEC range covered by this night of observation is shaded based on a 12 degree FWHM of the beam. Horizontal dashed lines represent the stripe that HERA can observe, while the shaded region is what was observed on this night. Vertical lines represent the beginning and ending LSTs of this observation. Selected sources are labelled, sources included are those in the GLEAM 4Jy catalog with a flux >10.9 Jy. Note that the map is clipped at the northern horizon.
sources = utils.gather_source_list()
utils.plot_sky_map(uvd,dec_pad=55,ra_pad=55,clip=False,sources=sources)
LST Coverage¶
Shows the LSTs (in hours) and JDs for which data is collected. Green represents data, red means no data.
utils.plot_lst_coverage(uvd)
Autocorrelations for a single file¶
This plot shows autocorrelations for one timestamp of each antenna that is active and each polarization. For each node, antennas are ordered by SNAP number, and within that by SNAP input number. The antenna number label color corresponds to the a priori status of that antenna.
### plot autos
utils.plot_autos(uvdx, uvdy)
Waterfalls of Autocorrelation Amplitudes for each Antenna and Each polarization¶
These plots show autocorrelation waterfalls of each antenna that is active and whose status qualifies for this notebook. For each node, antennas are ordered by SNAP number, and within that by SNAP input number. The antenna number label color corresponds to the a priori status of that antenna.
utils.plot_wfs(uvd, pol = 0)
utils.plot_wfs(uvd, pol = 1)
Correlation Metrics¶
The first plot shows the correlation metric (described below) for a set of baseline types, as calculated at several times throughout the night. It is expected that longer baselines (darker color) will exhibit lower values than the short baselines.
The matrices show the phase correlation between antennas. Using the even and odd visibilities, each pixel is calculated as (even/abs(even)) * (conj(odd)/abs(odd)), and then averaged across time and frequency. If the phases are noise-like, this value will average down to zero. If the antennas are well correlated, the phases should not be noise-like, and this value should average to 1. The lines denoting node boundaries are intended to help confirm that inter-node correlations are functioning - if they aren't, this plot will appear block-diagonal.
This metric has shown to be LST locked - when comparing to other nights, be sure to compare for the same LST. It is expected that some LSTs will look much better or worse than others.
Note: Within each node, the order of antennas is determined by snap, and within that by snap input number.
badAnts = []
badAnts = utils.plotNodeAveragedSummary(uvd,HHfiles,JD,use_ants,mat_pols=['xx','yy','xy','yx'])
WARNING: unable to read /mnt/sn1/data2/2460902/zen.2460902.41757.sum.uvh5
Visibility amplitude spectra for a set of redundant baselines, labeled by inter vs. intranode baselines. The red and blue should exhibit the same bandpass shape - if the red are consistently different from the blue, this indicates an issue with internode correlations.
Note: antennas that were identified as bad by the correlation matrix have been removed from this plot.
utils.plotVisibilitySpectra(HHfiles[len(HHfiles)//2+1], JD, use_ants, badAnts=[])
<Figure size 640x480 with 0 Axes>
Even and Odd File Checks¶
A waterfall showing the ratio between the even and odd visibilities. The purpose of this is to highlight xengine failures, which will cause this value to fall to zero or go to infinity. If things are working properly, this value should be stable at 1. The boundaries between different x-engines are shown by the vertical white lines.
if len(HHautos) == len(diffautos):
uvd_diff = UVData()
uvd_diff.read(diffautos[::10], skip_bad_files=True, antenna_nums=use_ants)
rat = utils.plotEvenOddWaterfalls(uvd,uvd_diff)
else:
uvd_diff = UVData()
use_diffs = [f for f in diffautos if '%s/zen.%s.%s.sum.autos.uvh5' % (data_path,f.split('.')[1],f.split('.')[2]) in HHautos[::10]]
uvd_diff.read(use_diffs, skip_bad_files=True, antenna_nums = use_ants)
uvd_sum = uvd.select(times=np.unique(uvd_diff.time_array),inplace=False)
rat = utils.plotEvenOddWaterfalls(uvd_sum,uvd_diff)
Crossed Antenna Check¶
These are differences between different panels of the correlation matrices shown above (see panel titles for specifics). Antennas showing as consistently blue are ones which are correlating stronger in the cross pols than in the auto pols, indicating that the antenna polarizations are likely crossed.
crossedAnts = utils.plotNodeAveragedSummary(uvd,HHfiles,JD,use_ants,mat_pols=['xx','yy','xy','yx'],plotRatios=True,
plotSummary=False)
Antenna Positions¶
Antennas outlined in black here have been identified by the correlation matrix as bad antennas. Antennas with a colorful outline correspond to their status as identified by ant_metrics (see above plot). Faded antennas are those not meeting the apriori status requirement for this notebook run. Gold stars are node box locations.
uvd1 = UVData()
uvd1.read(HHfiles[readInd], skip_bad_files=True)
utils.plot_antenna_positions(uvd1, badAnts=badAnts,use_ants=use_ants)
Observer Inspection Plots¶
Antennas of status digital_OK or better that are flagged as bad by any of the above metrics are plotted here so observers can inspect their failures in more detail. Additionally, a 'good' template has been used to identify outliers. The upper line plots are averages over the whole observation, and the lower line plots are slices of a single time in the middle of the observation. These plots are recommended diagnostics for demoting antennas to lower statuses or reporting issues. If the plots below look OK, check other plots in notebook to hunt why the antenna was flagged. NOTE: The colorbar/power scales in these plots are NOT locked between antennas OR polarizations so that the detail will be visible on all plots. Be sure to check for reasonable power levels, as this may be the reason the antenna was flagged for inspection.
d, tempAnts = utils.flag_by_template(uvd,HHautos,JD,use_ants=use_ants,pols=['XX','YY'],plotMap=False)
--------------------------------------------------------------------------- ValueError Traceback (most recent call last) Cell In[14], line 1 ----> 1 d, tempAnts = utils.flag_by_template(uvd,HHautos,JD,use_ants=use_ants,pols=['XX','YY'],plotMap=False) File ~/mambaforge/envs/RTP/lib/python3.12/site-packages/hera_notebook_templates/utils.py:112, in flag_by_template(uvd, HHfiles, jd, use_ants, pols, polDirs, temp_norm, plotMap) 110 for i,lst in enumerate(use_lsts): 111 hdat = UVData() --> 112 hdat.read(use_files[i],antenna_nums=use_ants) 113 for p,pol in enumerate(pols): 114 ant_dfs[pol][lst] = {} File ~/mambaforge/envs/RTP/lib/python3.12/site-packages/pyuvdata/uvdata/uvdata.py:11125, in UVData.read(self, filename, axis, file_type, read_data, skip_bad_files, background_lsts, astrometry_library, ignore_name, use_future_array_shapes, antenna_nums, antenna_names, ant_str, bls, catalog_names, frequencies, freq_chans, times, time_range, lsts, lst_range, polarizations, blt_inds, phase_center_ids, keep_all_metadata, run_check, check_extra, run_check_acceptability, strict_uvw_antpos_check, check_autos, fix_autos, projected, correct_lat_lon, calc_lst, fix_old_proj, fix_use_ant_pos, params_file, obs_file, flags_file, layout_file, settings_file, data_column, pol_order, ignore_single_chan, raise_error, read_weights, allow_flex_pol, multidim_index, remove_flex_pol, blt_order, blts_are_rectangular, time_axis_faster_than_bls, data_array_dtype, use_aoflagger_flags, remove_dig_gains, remove_coarse_band, correct_cable_len, correct_van_vleck, cheby_approx, flag_small_auto_ants, propagate_coarse_flags, flag_init, edge_width, start_flag, end_flag, flag_dc_offset, remove_flagged_ants, phase_to_pointing_center, nsample_array_dtype, corrchunk, receivers, sidebands, mir_select_where, apply_tsys, apply_flags, apply_dedoppler, pseudo_cont, rechunk, compass_soln, swarm_only, codes_check, recompute_nbls) 11106 self.read_ms( 11107 filename, 11108 data_column=data_column, (...) 11121 astrometry_library=astrometry_library, 11122 ) 11124 elif file_type == "uvh5": > 11125 self.read_uvh5( 11126 filename, 11127 antenna_nums=antenna_nums, 11128 antenna_names=antenna_names, 11129 ant_str=ant_str, 11130 bls=bls, 11131 frequencies=frequencies, 11132 freq_chans=freq_chans, 11133 times=times, 11134 time_range=time_range, 11135 lsts=lsts, 11136 lst_range=lst_range, 11137 polarizations=polarizations, 11138 blt_inds=blt_inds, 11139 phase_center_ids=phase_center_ids, 11140 catalog_names=catalog_names, 11141 read_data=read_data, 11142 data_array_dtype=data_array_dtype, 11143 keep_all_metadata=keep_all_metadata, 11144 multidim_index=multidim_index, 11145 remove_flex_pol=remove_flex_pol, 11146 background_lsts=background_lsts, 11147 run_check=run_check, 11148 check_extra=check_extra, 11149 run_check_acceptability=run_check_acceptability, 11150 strict_uvw_antpos_check=strict_uvw_antpos_check, 11151 fix_old_proj=fix_old_proj, 11152 fix_use_ant_pos=fix_use_ant_pos, 11153 check_autos=check_autos, 11154 fix_autos=fix_autos, 11155 time_axis_faster_than_bls=time_axis_faster_than_bls, 11156 blts_are_rectangular=blts_are_rectangular, 11157 recompute_nbls=recompute_nbls, 11158 astrometry_library=astrometry_library, 11159 ) 11160 select = False 11162 if select: File ~/mambaforge/envs/RTP/lib/python3.12/site-packages/pyuvdata/uvdata/uvdata.py:9816, in UVData.read_uvh5(self, filename, **kwargs) 9809 raise ValueError( 9810 "Reading multiple files from class specific " 9811 "read functions is no longer supported. " 9812 "Use the generic `uvdata.read` function instead." 9813 ) 9815 uvh5_obj = uvh5.UVH5() -> 9816 uvh5_obj.read_uvh5(filename, **kwargs) 9817 self._convert_from_filetype(uvh5_obj) 9818 del uvh5_obj File ~/mambaforge/envs/RTP/lib/python3.12/site-packages/pyuvdata/uvdata/uvh5.py:1074, in UVH5.read_uvh5(self, filename, antenna_nums, antenna_names, ant_str, bls, frequencies, freq_chans, times, time_range, lsts, lst_range, polarizations, blt_inds, phase_center_ids, catalog_names, keep_all_metadata, read_data, data_array_dtype, multidim_index, remove_flex_pol, background_lsts, run_check, check_extra, run_check_acceptability, strict_uvw_antpos_check, fix_old_proj, fix_use_ant_pos, check_autos, fix_autos, use_future_array_shapes, blt_order, blts_are_rectangular, time_axis_faster_than_bls, recompute_nbls, astrometry_library) 1063 self._read_header( 1064 meta, 1065 run_check=run_check, (...) 1069 astrometry_library=astrometry_library, 1070 ) 1072 if read_data: 1073 # Now read in the data -> 1074 self._get_data( 1075 meta.datagrp, 1076 antenna_nums=antenna_nums, 1077 antenna_names=antenna_names, 1078 ant_str=ant_str, 1079 bls=bls, 1080 frequencies=frequencies, 1081 freq_chans=freq_chans, 1082 times=times, 1083 time_range=time_range, 1084 lsts=lsts, 1085 lst_range=lst_range, 1086 polarizations=polarizations, 1087 blt_inds=blt_inds, 1088 phase_center_ids=phase_center_ids, 1089 catalog_names=catalog_names, 1090 data_array_dtype=data_array_dtype, 1091 keep_all_metadata=keep_all_metadata, 1092 multidim_index=multidim_index, 1093 ) 1094 if close_meta: 1095 meta.close() File ~/mambaforge/envs/RTP/lib/python3.12/site-packages/pyuvdata/uvdata/uvh5.py:749, in UVH5._get_data(self, dgrp, antenna_nums, antenna_names, ant_str, bls, frequencies, freq_chans, times, time_range, lsts, lst_range, polarizations, blt_inds, phase_center_ids, catalog_names, data_array_dtype, keep_all_metadata, multidim_index) 744 raise ImportError( 745 "hdf5plugin is not installed but is required to read this dataset" 746 ) from hdf5plugin_error 748 # figure out what data to read in --> 749 blt_inds, freq_inds, pol_inds, history_update_string = self._select_preprocess( 750 antenna_nums=antenna_nums, 751 antenna_names=antenna_names, 752 ant_str=ant_str, 753 bls=bls, 754 frequencies=frequencies, 755 freq_chans=freq_chans, 756 times=times, 757 time_range=time_range, 758 lsts=lsts, 759 lst_range=lst_range, 760 polarizations=polarizations, 761 blt_inds=blt_inds, 762 phase_center_ids=phase_center_ids, 763 catalog_names=catalog_names, 764 ) 766 # figure out which axis is the most selective 767 if blt_inds is not None: File ~/mambaforge/envs/RTP/lib/python3.12/site-packages/pyuvdata/uvdata/uvdata.py:6719, in UVData._select_preprocess(self, antenna_nums, antenna_names, ant_str, bls, frequencies, freq_chans, times, time_range, lsts, lst_range, polarizations, blt_inds, phase_center_ids, catalog_names) 6710 if bls is not None: 6711 bls, polarizations = utils.bls._extract_bls_pol( 6712 bls=bls, 6713 polarizations=polarizations, (...) 6717 nants_telescope=self.telescope.Nants, 6718 ) -> 6719 blt_inds, blt_selections = utils.bltaxis._select_blt_preprocess( 6720 select_antenna_nums=antenna_nums, 6721 select_antenna_names=antenna_names, 6722 bls=bls, 6723 times=times, 6724 time_range=time_range, 6725 lsts=lsts, 6726 lst_range=lst_range, 6727 blt_inds=blt_inds, 6728 phase_center_ids=phase_center_ids, 6729 antenna_names=self.telescope.antenna_names, 6730 antenna_numbers=self.telescope.antenna_numbers, 6731 ant_1_array=self.ant_1_array, 6732 ant_2_array=self.ant_2_array, 6733 baseline_array=self.baseline_array, 6734 time_array=self.time_array, 6735 time_tols=self._time_array.tols, 6736 lst_array=self.lst_array, 6737 lst_tols=self._lst_array.tols, 6738 phase_center_id_array=self.phase_center_id_array, 6739 ) 6740 selections.extend(blt_selections) 6742 freq_inds, freq_selections = utils.frequency._select_freq_helper( 6743 frequencies=frequencies, 6744 freq_chans=freq_chans, (...) 6749 obj_spw_id_array=self.flex_spw_id_array, 6750 ) File ~/mambaforge/envs/RTP/lib/python3.12/site-packages/pyuvdata/utils/bltaxis.py:386, in _select_blt_preprocess(select_antenna_nums, select_antenna_names, bls, times, time_range, lsts, lst_range, blt_inds, phase_center_ids, antenna_names, antenna_numbers, ant_1_array, ant_2_array, baseline_array, time_array, time_tols, lst_array, lst_tols, phase_center_id_array) 381 ant_check = np.logical_or( 382 np.isin(select_antenna_nums, ant_1_array), 383 np.isin(select_antenna_nums, ant_2_array), 384 ) 385 if not np.all(ant_check): --> 386 raise ValueError( 387 f"Antenna number {select_antenna_nums[~ant_check]} is not present " 388 "in the ant_1_array or ant_2_array" 389 ) 390 ant_blt_inds = np.where( 391 np.logical_and( 392 np.isin(ant_1_array, select_antenna_nums), 393 np.isin(ant_2_array, select_antenna_nums), 394 ) 395 )[0] 396 ant_blt_inds = np.asarray(ant_blt_inds, dtype=np.int64) ValueError: Antenna number [157 158 176 177 178 329 333 156 155 136 179 299 302 311 314 343 346 135 300 312 347 342 301 313] is not present in the ant_1_array or ant_2_array
inspectAnts = utils.plot_inspect_ants(uvd,JD,badAnts=badAnts,use_ants=use_ants,
tempAnts=tempAnts,crossedAnts=crossedAnts)
--------------------------------------------------------------------------- NameError Traceback (most recent call last) Cell In[15], line 2 1 inspectAnts = utils.plot_inspect_ants(uvd,JD,badAnts=badAnts,use_ants=use_ants, ----> 2 tempAnts=tempAnts,crossedAnts=crossedAnts) NameError: name 'tempAnts' is not defined
Mean-Subtracted Waterfalls¶
Here the mean value in each frequency bin has been subtracted out. This effectively subtracts out the bandpass shape, making time variations more visible.
utils.plot_wfs(uvd,0,mean_sub=True,jd=JD)
utils.plot_wfs(uvd,1,mean_sub=True,jd=JD)
