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

Cabled EK60 and EK80 updates #5

Open
wants to merge 2 commits into
base: main
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 1 addition & 1 deletion VERSION.txt
Original file line number Diff line number Diff line change
@@ -1 +1 @@
0.2.0
0.2.1
5 changes: 2 additions & 3 deletions environment.yml
Original file line number Diff line number Diff line change
@@ -1,14 +1,13 @@
name: echogram
channels:
- defaults
- conda-forge
dependencies:
- python>=3.9
- python>=3.10
- argcomplete
- bottleneck
- cmocean
- conda-build
- echopype>=0.8.1
- echopype>=0.8.4
- h5netcdf
- matplotlib
- numpy
Expand Down
85 changes: 55 additions & 30 deletions ooi_zpls_echograms/zpls_echogram.py
Original file line number Diff line number Diff line change
Expand Up @@ -382,7 +382,8 @@ def generate_echogram(data, site, long_name, deployed_depth, output_directory, f
# populate the subplots
im = []
for index in range(len(frequency_list)):
im.append(data.isel(frequency_nominal=index).Sv.plot(x='ping_time', y='echo_range', vmin=v_min, vmax=v_max,
# TODO: It feels like I need to dropna here as well...
im.append(data.isel(frequency_nominal=index).Sv.sortby('ping_time').plot(x='ping_time', y='echo_range', vmin=v_min, vmax=v_max,
ax=ax[index], cmap=my_cmap, add_colorbar=False))
ax_config(ax[index], frequency_list[index])

Expand Down Expand Up @@ -468,7 +469,7 @@ def ek60_file_list(data_directory, dates):
file_list = []
for i in range(delta.days + 1):
day = sdate + timedelta(days=i)
ek60_files = glob.glob(os.path.join(data_directory, day.strftime('%m'), day.strftime('%d')) + '/*.raw')
ek60_files = glob.glob(os.path.join(data_directory, day.strftime('%Y'), day.strftime('%m'), day.strftime('%d')) + '/*.raw')
file_list.append(ek60_files)

return file_list
Expand Down Expand Up @@ -516,7 +517,7 @@ def process_sonar_data(site, data_directory, output_directory, dates, zpls_model
file_list.sort()

# convert and process the raw files using echopype
desc = 'Converting and processing %d raw %s data files' % (len(file_list), zpls_model)
desc = 'Converting and processing %d raw %s data files [%s]' % (len(file_list), zpls_model, dates)
echo = [_process_file(file, site, output_directory, zpls_model, xml_file, tilt_correction)
for file in tqdm(file_list, desc=desc)]

Expand Down Expand Up @@ -597,7 +598,17 @@ def _process_file(file, site, output_directory, zpls_model, xml_file, tilt_corre
ds.to_netcdf(Path(output_directory))

# process the data, calculating the volume acoustic backscatter strength and the vertical range
ds_sv = ep.calibrate.compute_Sv(ds, env_params=env_params) # calculate Sv
'''
The EK80 echosounder can be configured to transmit either broadband (waveform_mode="BB") or narrowband (waveform_mode="CW") signals.
When transmitting in broadband mode, the returned echoes are encoded as complex samples (encode_mode="complex").
When transmitting in narrowband mode, the returned echoes can be encoded either as complex samples (encode_mode="complex")
or as power/angle combinations (encode_mode="power") in a format similar to those recorded by EK60 echosounders.
'''
if zpls_model == 'EK80':
ds_sv = ep.calibrate.compute_Sv(ds, env_params=env_params, waveform_mode='CW',
encode_mode='power') # calculate Sv
else:
ds_sv = ep.calibrate.compute_Sv(ds, env_params=env_params) # calculate Sv

# calculate the depth from the range and convert the channels dimension to frequency
ds_sv = ep.consolidate.add_depth(ds_sv, depth_offset=depth_offset, tilt=tilt_correction, downward=downward)
Expand All @@ -623,7 +634,7 @@ def _process_file(file, site, output_directory, zpls_model, xml_file, tilt_corre
return data


def zpls_echogram(site, data_directory, output_directory, dates, zpls_model, xml_file, **kwargs):
def zpls_echogram(site, data_directory, output_directory, dates, zpls_model, xml_file, make_echogram, **kwargs):
"""
Main processing function to convert and process data from either the ASL
AZFP or the Kongsberg Simrad EK60. Uses echopype to convert the raw data
Expand All @@ -645,6 +656,7 @@ def zpls_echogram(site, data_directory, output_directory, dates, zpls_model, xml
:param xml_file: If the model is an AZFP, the xml_file (with instrument
calibration coefficients needed for the conversion) must also be
specified (usually just one per deployment)
:param make_echogram: Boolean whether to process the echogram PNG output
:kwargs tilt_correction: Tilt of the sonar transducers (typically 15
degrees for the uncabled sensors to avoid interference from the
riser elements)
Expand All @@ -671,8 +683,8 @@ def zpls_echogram(site, data_directory, output_directory, dates, zpls_model, xml
'configuration and calibration parameters.')

# convert and process the data
if zpls_model not in ['AZFP', 'EK60']:
raise ValueError('The ZPLS model must be set as either AZFP or EK60 (case sensitive)')
if zpls_model not in ['AZFP', 'EK60', 'EK80']:
raise ValueError('The ZPLS model must be set as either AZFP, EK60, or EK80 (case sensitive)')
else:
data = process_sonar_data(site, data_directory, output_directory, dates, zpls_model, xml_file, tilt_correction)

Expand All @@ -684,6 +696,10 @@ def zpls_echogram(site, data_directory, output_directory, dates, zpls_model, xml
# save the full resolution processed data to daily NetCDF files
file_name = set_file_name(site, dates)

# clean the interwoven NaN values out due to the different sampling rates on the channels
if zpls_model == 'EK80':
data['Sv'] = data["Sv"].dropna(dim="ping_time", how="all")

# reset data types (helps to control size of NetCDF files)
data['range_sample'] = data['range_sample'].astype(np.int32)
data['echo_range'] = data['echo_range'].astype(np.float32)
Expand All @@ -694,6 +710,11 @@ def zpls_echogram(site, data_directory, output_directory, dates, zpls_model, xml
# split the data into daily records
days, datasets = zip(*data.groupby("ping_time.day"))

# TODO: maybe a EK80 thing...or from processing one day of data, doesn't seem correct to do this
# fixes dataset count to match day/file _Full count
if zpls_model == 'EK80':
datasets = datasets[:-1]

# create a list of file names based on the day of the record
start = datetime.strptime(dates[0], '%Y%m%d')
stop = datetime.strptime(dates[1], '%Y%m%d')
Expand All @@ -718,38 +739,42 @@ def zpls_echogram(site, data_directory, output_directory, dates, zpls_model, xml
xr.save_mfdataset(datasets, nc_files, mode='w', format='NETCDF4', engine='h5netcdf')

# clean up any NaN's in the range values (some EK60 files seem to have this problem)
data = data.dropna('range_sample', subset=['echo_range'])
if zpls_model != 'EK80':
data = data.dropna('range_sample', subset=['echo_range'])

# if a global mooring, create hourly averaged data records, otherwise create 15-minute records
if 'HYPM' in site:
# resample the data into a 60 minute, median averaged record, filling gaps less than 180 minutes
avg = data.resample(ping_time='60Min', skipna=True).median(dim='ping_time', keep_attrs=True)
avg = avg.interpolate_na(dim='ping_time', max_gap='180Min')
if zpls_model == 'EK80':
avg = data.resample(ping_time='15Min', skipna=True).median(dim='ping_time', keep_attrs=True)
# avg = avg.interpolate_na(dim='ping_time', max_gap='45Min')
else:
# resample the data into a 15 minute, median averaged record, filling gaps less than 45 minutes
avg = data.resample(ping_time='15Min', skipna=True).median(dim='ping_time', keep_attrs=True)
avg = avg.interpolate_na(dim='ping_time', max_gap='45Min')

# generate the echogram
long_name = site_config[site]['long_name']
generate_echogram(avg, site, long_name, deployed_depth, output_directory, file_name, dates,
vertical_range=vertical_range, colorbar_range=colorbar_range)

# add the OOI logo as a watermark
echogram = os.path.join(output_directory, file_name + '.png')
echo_image = Image.open(echogram)
# noinspection PyTypeChecker
ooi_image = Image.open(files('ooi_zpls_echograms').joinpath('ooi-logo.png'))
width, height = echo_image.size
transparent = Image.new('RGBA', (width, height), (0, 0, 0, 0))
transparent.paste(echo_image, (0, 0))
if max(vertical_range) > 99:
transparent.paste(ooi_image, (96, 15), mask=ooi_image)
else:
transparent.paste(ooi_image, (80, 15), mask=ooi_image)

# re-save the echogram with the added logo
transparent.save(echogram)
if make_echogram:
long_name = site_config[site]['long_name']
generate_echogram(avg, site, long_name, deployed_depth, output_directory, file_name, dates, vertical_range=vertical_range, colorbar_range=colorbar_range)

# add the OOI logo as a watermark
echogram = os.path.join(output_directory, file_name + '.png')
echo_image = Image.open(echogram)
# noinspection PyTypeChecker
ooi_image = Image.open(files('ooi_zpls_echograms').joinpath('ooi-logo.png'))
width, height = echo_image.size
transparent = Image.new('RGBA', (width, height), (0, 0, 0, 0))
transparent.paste(echo_image, (0, 0))
if max(vertical_range) > 99:
transparent.paste(ooi_image, (96, 15), mask=ooi_image)
else:
transparent.paste(ooi_image, (80, 15), mask=ooi_image)

# re-save the echogram with the added logo
transparent.save(echogram)

# save the averaged data
avg['ping_time'] = avg['ping_time'].values.astype(np.float64) / 10.0 ** 9
Expand Down Expand Up @@ -777,7 +802,7 @@ def main(argv=None):
help=('Date range to plot as either YYYYMM or YYYYMMDD. Specifying an end date is optional, '
'it will be assumed to be 1 month or 1 day depending on input.'))
parser.add_argument('-zm', '--zpls_model', dest='zpls_model', type=str, required=True,
help='Specifies the ZPLS instrument model, either AZFP or EK60.')
help='Specifies the ZPLS instrument model, either AZFP, EK60, EK80.')
parser.add_argument('-xf', '--xml_file', dest='xml_file', type=str, required=False,
help='The path to .XML file used to process the AZFP data in the .01A files')
parser.add_argument('-tc', '--tilt_correction', dest='tilt_correction', type=int, required=False,
Expand Down Expand Up @@ -826,8 +851,8 @@ def main(argv=None):
os.makedirs(output_directory, exist_ok=True)

# convert and process the data
if zpls_model not in ['AZFP', 'EK60']:
raise ValueError('The ZPLS model must be set as either AZFP or EK60 (case sensitive)')
if zpls_model not in ['AZFP', 'EK60', 'EK80']:
raise ValueError('The ZPLS model must be set as either AZFP, EK60, or EK80 (case sensitive)')
else:
zpls_echogram(site, data_directory, output_directory, dates, zpls_model, xml_file,
deployed_depth=deployed_depth, tilt_correction=tilt_correction,
Expand Down