Merge pull request #5 from lofar-astron/lbasparse

Add LBA sparse support
master
maaijke 2020-12-03 16:42:46 +01:00 committed by GitHub
commit 9672b52bb9
No known key found for this signature in database
GPG Key ID: 4AEE18F83AFDEB23
5 changed files with 102 additions and 80 deletions

View File

@ -61,7 +61,8 @@
"import warnings\n",
"from mpl_toolkits.axes_grid1 import make_axes_locatable\n",
"import matplotlib.axes as maxes\n",
"from matplotlib.ticker import FormatStrFormatter"
"from matplotlib.ticker import FormatStrFormatter\n",
"import re"
]
},
{
@ -118,8 +119,7 @@
"outputs": [],
"source": [
"data_dir = \"./test\"\n",
"caltable_dir = \"./test/CalTables\" # Root directory under which station information is stored in subdirectories DE602C/etc/, RS106/etc/, ...\n",
"station_name = 'DE603'"
"caltable_dir = \"./test/CalTables\" # Root directory under which station information is stored in subdirectories DE602C/etc/, RS106/etc/, ..."
]
},
{
@ -139,7 +139,7 @@
"metadata": {},
"outputs": [],
"source": [
"station_type = get_station_type(station_name)"
"os.makedirs('results', exist_ok=True)"
]
},
{
@ -147,22 +147,13 @@
"execution_count": 11,
"metadata": {},
"outputs": [],
"source": [
"os.makedirs('results', exist_ok=True)"
]
},
{
"cell_type": "code",
"execution_count": 12,
"metadata": {},
"outputs": [],
"source": [
"start1 = time.time()"
]
},
{
"cell_type": "code",
"execution_count": 13,
"execution_count": 12,
"metadata": {},
"outputs": [
{
@ -170,7 +161,7 @@
"output_type": "stream",
"text": [
"Searching for available files in ./test\n",
"0: ./test/20200525_084652_mode_5_xst_sb100.dat\n",
"0: ./test/20170621_072634_mode_sparseeven_xst_sb350.dat\n",
"1: ./test/20170720_095816_mode_3_xst_sb297.dat\n"
]
}
@ -185,27 +176,58 @@
},
{
"cell_type": "code",
"execution_count": 58,
"execution_count": 33,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"File selected: ./test/20200525_084652_mode_5_xst_sb100.dat\n"
"File selected: ./test/20170720_095816_mode_3_xst_sb297.dat\n"
]
}
],
"source": [
"# Select a file\n",
"xst_filename = files[0]\n",
"xst_filename = \"./test/20170720_095816_mode_3_xst_sb297.dat\" # files[0]\n",
"\n",
"print(\"File selected:\", xst_filename)"
]
},
{
"cell_type": "code",
"execution_count": 15,
"execution_count": 34,
"metadata": {},
"outputs": [],
"source": [
"station_name = \"DE603\"\n",
"\n",
"rcu_mode = re.search(\"mode_([^_]*)\", xst_filename).groups(0)[0]\n",
"subband = int(re.search(\"sb([0-9]*)\", xst_filename).groups(0)[0])\n",
"obsdatestr, obstimestr, *_ = os.path.basename(xst_filename).rstrip(\".dat\").split(\"_\")"
]
},
{
"cell_type": "code",
"execution_count": 35,
"metadata": {},
"outputs": [],
"source": [
"rcu_mode = rcu_mode.replace(\"sparseeven\", \"sparse_even\").replace(\"sparseodd\", \"sparse_odd\")"
]
},
{
"cell_type": "code",
"execution_count": 36,
"metadata": {},
"outputs": [],
"source": [
"station_type = get_station_type(station_name)"
]
},
{
"cell_type": "code",
"execution_count": 37,
"metadata": {},
"outputs": [
{
@ -225,10 +247,6 @@
}
],
"source": [
"# Distill metadata from filename\n",
"cubename = os.path.basename(xst_filename)\n",
"obsdatestr, obstimestr, _, rcu_mode, _, subbandname = cubename.rstrip(\".dat\").split(\"_\")\n",
"subband = int(subbandname[2:])\n",
"station_name = get_full_station_name(station_name, rcu_mode)\n",
"\n",
"# Get the data\n",
@ -246,8 +264,7 @@
"obstime = datetime.datetime.strptime(obsdatestr + \":\" + obstimestr, '%Y%m%d:%H%M%S')\n",
"\n",
"# Confirm the data has been read correctly\n",
"print(f\"\"\"Filename: {cubename}\n",
"Station: {station_name}\n",
"print(f\"\"\"Station: {station_name}\n",
"Station type: {station_type}\n",
"Subband: {subband}\n",
"Timestep: {timestep}\n",
@ -260,27 +277,7 @@
},
{
"cell_type": "code",
"execution_count": 16,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"'20170720_095816_DE603LBA_SB297'"
]
},
"execution_count": 16,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"fname"
]
},
{
"cell_type": "code",
"execution_count": 17,
"execution_count": 40,
"metadata": {},
"outputs": [],
"source": [
@ -290,7 +287,7 @@
},
{
"cell_type": "code",
"execution_count": 18,
"execution_count": 41,
"metadata": {},
"outputs": [],
"source": [
@ -299,7 +296,7 @@
},
{
"cell_type": "code",
"execution_count": 19,
"execution_count": 42,
"metadata": {},
"outputs": [],
"source": [
@ -315,7 +312,7 @@
},
{
"cell_type": "code",
"execution_count": 20,
"execution_count": 43,
"metadata": {},
"outputs": [
{
@ -324,7 +321,7 @@
"((29, 96, 96), (29, 96, 96), (96, 96), dtype('complex128'))"
]
},
"execution_count": 20,
"execution_count": 43,
"metadata": {},
"output_type": "execute_result"
}
@ -408,7 +405,7 @@
},
{
"cell_type": "code",
"execution_count": 23,
"execution_count": 46,
"metadata": {},
"outputs": [],
"source": [
@ -418,7 +415,7 @@
},
{
"cell_type": "code",
"execution_count": 24,
"execution_count": 47,
"metadata": {},
"outputs": [],
"source": [
@ -452,7 +449,7 @@
},
{
"cell_type": "code",
"execution_count": 27,
"execution_count": 49,
"metadata": {},
"outputs": [],
"source": [
@ -487,7 +484,7 @@
},
{
"cell_type": "code",
"execution_count": 29,
"execution_count": 51,
"metadata": {},
"outputs": [
{

View File

@ -101,6 +101,9 @@ def write_hdf5(filename: str, xst_data: np.ndarray, visibilities: np.ndarray, sk
"""
short_station_name = station_name[:5]
if subtracted is None:
subtracted = []
with h5py.File(filename, 'a') as h5file:
new_obsname = get_new_obsname(h5file)
obs_group = h5file.create_group(new_obsname)

View File

@ -100,16 +100,18 @@ def freq_from_sb(sb: int, rcu_mode: Union[str, int] = 1):
58007812.5
"""
clock = 200e6
if int(rcu_mode) == 6:
clock = 160e6
freq_offset = 0
if int(rcu_mode) == 5:
freq_offset = 100e6
elif int(rcu_mode) == 6:
freq_offset = 160e6
elif int(rcu_mode) == 7:
freq_offset = 200e6
if 'sparse' not in str(rcu_mode):
if int(rcu_mode) == 6:
clock = 160e6
if int(rcu_mode) == 5:
freq_offset = 100e6
elif int(rcu_mode) == 6:
freq_offset = 160e6
elif int(rcu_mode) == 7:
freq_offset = 200e6
sb_bandwidth = 0.5 * clock / 512.
freq = (sb * sb_bandwidth) + freq_offset
@ -150,6 +152,10 @@ def find_caltable(field_name: str, rcu_mode: Union[str, int], caltable_dir='calt
filename += "-HBA-170_230.dat"
elif str(rcu_mode) == '7':
filename += "-HBA-210_250.dat"
elif str(rcu_mode) == 'sparse_even':
filename += "-LBA_SPARSE_EVEN-10_90.dat"
elif str(rcu_mode) == 'sparse_odd':
filename += "-LBA_SPARSE_ODD-10_90.dat"
else:
raise RuntimeError("Unexpected mode: " + str(rcu_mode) + " for field_name " + str(field_name))
@ -317,15 +323,24 @@ def get_station_pqr(station_name: str, rcu_mode: Union[str, int], db):
full_station_name = get_full_station_name(station_name, rcu_mode)
station_type = get_station_type(full_station_name)
if 'LBA' in station_name or str(rcu_mode) in ('1', '2', '3', '4', 'inner', 'outer'):
# Get the PQR positions for an individual station
station_pqr = db.antenna_pqr(full_station_name)
# Exception: for Dutch stations (sparse not yet accommodated)
if (station_type == 'core' or station_type == 'remote') and int(rcu_mode) in (3, 4):
station_pqr = station_pqr[0:48, :]
elif (station_type == 'core' or station_type == 'remote') and int(rcu_mode) in (1, 2):
station_pqr = station_pqr[48:, :]
if 'LBA' in station_name or str(rcu_mode) in ('1', '2', '3', '4', 'inner', 'outer', 'sparse_even', 'sparse_odd', 'sparse'):
if (station_type == 'core' or station_type == 'remote'):
if str(rcu_mode) in ('3', '4', 'inner'):
station_pqr = db.antenna_pqr(full_station_name)[0:48, :]
elif str(rcu_mode) in ('1', '2', 'outer'):
station_pqr = db.antenna_pqr(full_station_name)[48:, :]
elif rcu_mode in ('sparse_even', 'sparse'):
all_pqr = db.antenna_pqr(full_station_name)
# Indices 0, 49, 2, 51, 4, 53, ...
station_pqr = np.ravel(np.column_stack((all_pqr[:48:2], all_pqr[49::2]))).reshape(48, 3)
elif rcu_mode == 'sparse_odd':
all_pqr = db.antenna_pqr(full_station_name)
# Indices 1, 48, 3, 50, 5, 52, ...
station_pqr = np.ravel(np.column_stack((all_pqr[1:48:2], all_pqr[48::2]))).reshape(48, 3)
else:
raise RuntimeError("Cannot select subset of LBA antennas for mode " + rcu_mode)
else:
station_pqr = db.antenna_pqr(full_station_name)
elif 'HBA' in station_name or str(rcu_mode) in ('5', '6', '7', '8'):
selected_dipole_config = {
'intl': GENERIC_INT_201512, 'remote': GENERIC_REMOTE_201512, 'core': GENERIC_CORE_201512
@ -549,14 +564,13 @@ def get_full_station_name(station_name: str, rcu_mode: Union[str, int]) -> str:
return station_name
if str(rcu_mode) in ('1', '2', 'outer'):
if len(station_name) == 5:
station_name += "LBA"
station_name += "LBA"
elif str(rcu_mode) in ('3', '4', 'inner'):
if len(station_name) == 5:
station_name += "LBA"
station_name += "LBA"
elif 'sparse' in str(rcu_mode):
station_name += "LBA"
elif str(rcu_mode) in ('5', '6', '7'):
if len(station_name) == 5:
station_name += "HBA"
station_name += "HBA"
else:
raise Exception("Unexpected rcu_mode: ", rcu_mode)
@ -645,6 +659,14 @@ def make_xst_plots(xst_data: np.ndarray,
>>> type(leafletmap)
<class 'folium.folium.Map'>
>>> xst_data = read_acm_cube("test/20170621_072634_sb350_xst.dat", "remote")[0]
>>> obstime = datetime.datetime(2017, 6, 21, 7, 26, 34)
>>> sky_fig, ground_fig, leafletmap = make_xst_plots(xst_data, "RS509", obstime, 350, \
'sparse_even', \
caltable_dir="test/CalTables", \
hdf5_filename="test/test.h5")
Maximum at 2m east, -2m north of station center (lat/long 53.40884, 6.78531)
"""
if extent is None:
extent = [-150, 150, -150, 150]
@ -853,7 +875,7 @@ def reimage_sky(h5: h5py.File, obsnum: str, db: lofarantpos.db.LofarAntennaDatab
sky_data = h5[obsnum]["sky_img"]
freq = h5[obsnum].attrs['frequency']
marked_bodies_lmn = dict(zip(h5[obsnum].attrs["source_names"], h5[obsnum].attrs["source_lmn"]))
visibilities = h5[obsnum]['calibrated_data'].value
visibilities = h5[obsnum]['calibrated_data'][:]
visibilities_xx = visibilities[0::2, 0::2]
visibilities_yy = visibilities[1::2, 1::2]
# Stokes I
@ -902,7 +924,7 @@ def reimage_nearfield(h5: h5py.File, obsnum: str, db: lofarantpos.db.LofarAntenn
rcu_mode = h5[obsnum].attrs['rcu_mode']
freq = h5[obsnum].attrs['frequency']
marked_bodies_lmn = dict(zip(h5[obsnum].attrs["source_names"], h5[obsnum].attrs["source_lmn"]))
visibilities = h5[obsnum]['calibrated_data'].value
visibilities = h5[obsnum]['calibrated_data'][:]
visibilities_xx = visibilities[0::2, 0::2]
visibilities_yy = visibilities[1::2, 1::2]
# Stokes I

Binary file not shown.