Identify zooplankton school from glider echosounder data with external location and environmental data incorporated

Identify zooplankton school from glider echosounder data with external location and environmental data incorporated#

This notebook accompanies the following paper:

Echopype: Interoperable and scalable echosounder data processing with Echopype
Authors: Wu-Jung Lee, Landung Setiawan, Caesar Tuguinay, Emilio Mayorga, and Valentina Staneva
ICES Journal of Marine Science 2024: https://doi.org/10.1093/icesjms/fsae133

Introduction#

Description#

This notebook uses echosounder data from a glider to illustrate how Echopype can use external location (latitude and longitude) and environmental (temperature, salinity and pressure) data in the processing, if the data do not already exist in the echosounder raw files or if users would like to udpdate them. The data are provided by Delphine Mossman from the Department of Marine and Coastal Sciences at Rutgers University. The data were a 15-min section of echosounder data collected by an Acoustic Zooplankton and Fisher Profiler (AZFP) mounted on a Slocum glider deployed in the Southern Ocean off Antarctica (Ross Sea) in 2018.

In this notebook we attempt to correct for the glider orientation and identify zooplankton swarms in the echogram following steps detailed in the 2021 M.S. thesis from Ailey Sheehan, another member of the same Rutgers team. The computations are available in the associated GitHub repository.

Outline#

Open and convert raw AZFP data
Update EchoData Platform group with glider data and metadata
Calibrate backscatter data with glider mean environmental data
Remove background noise
Align echosounder data according to depth
Compute and apply shoal mask on corrected Sv

Running the notebook#

This notebook can be run with a conda environment created using the conda environment file. It uses data files found in the notebooks/example_data/glider_AZFP directory in this repository.

Note#

We encourage importing Echopype as ep for consistency.

from pathlib import Path

import geopandas as gpd
import pandas as pd
import xarray as xr
from scipy.spatial.transform import Rotation as R
import numpy as np

import matplotlib.pyplot as plt
import cartopy.crs as ccrs

import echopype as ep
import echopype.colormap

import warnings
warnings.simplefilter("ignore", category=DeprecationWarning)

Path setup#

# Set path to example data
data_path = Path('./example_data/glider_AZFP')

Data processing workflow#

Update EchoData `Platform` group with glider data and metadata#

Here we update the EchoData object with metadata entries for this deployment, and the GPS locations from the glider netCDF data.

Update EchoData `Top-level` and `Platform` group attributes#

While not required for data processing, it’s always a good idea to include as much metadata as possible. We’ll take advantage of the metadata found in the glider netCDF file to populate some of these metadata.

# Manually populate additional metadata about the dataset and the platform
# -- SONAR-netCDF4 Top-level Group attributes
ed["Top-level"] = ed["Top-level"].assign_attrs(
    title="2018 Ross Sea Slocum glider AZFP echosounder data from Rutgers University",
    summary=ds_glider.attrs.get("summary", ""),
)

# -- SONAR-netCDF4 Platform Group attributes
ed["Platform"] = ed["Platform"].assign_attrs(
    platform_type=ds_glider.attrs.get("platform_type", ""),
    platform_name="Rutgers r32 Slocum Webb G2 glider",
    platform_code_ICES=ds_glider.attrs.get("wmo_id", ""),
)

Here are the updated attributes in the “Top-level” group:

# Show top level attributes
ed['Top-level'].attrs

{'conventions': 'CF-1.7, SONAR-netCDF4-1.0, ACDD-1.3',
 'keywords': 'AZFP',
 'sonar_convention_authority': 'ICES',
 'sonar_convention_name': 'SONAR-netCDF4',
 'sonar_convention_version': '1.0',
 'summary': 'This project integrated an Acoustic Zooplankton and Fish Profiler (AZFP) multi-frequency echo sounder into a Slocum Webb G2 glider. The AZFP is complemented with existing glider sensors including a CTD, a WET Labs BB2FL ECO puck configured for simultaneous chlorophyll fluorescence (phytoplankton biomass) and optical backscatter measurements, and an Aanderaa Optode for measuring dissolved oxygen. This glider deployment is located in the polynya of Terra Nova Bay (western Ross Sea, Antarctica), and is focused on investigating relationships between phytoplankton-zooplankton-fish distributions and the physical drivers of zooplankton and silverfish species and size distributions.',
 'title': '2018 Ross Sea Slocum glider AZFP echosounder data from Rutgers University',
 'date_created': '2018-01-11T07:44:40Z'}

The external glider file name, glider_nc_fpath.name, is used in update_platform only to store it as provenance information recording the origin of the data.

Plot glider location#

Extract and join together the latitude and longitude variables from the Platform group in the ed EchoData object. Convert to a Pandas DataFrame first, then to a GeoPandas GeoDataFrame for convenient viewing and manipulation.

# Create gps dataframe
gps_df = ed['Platform'].latitude.to_dataframe().join(ed['Platform'].longitude.to_dataframe())
gps_df.head(3)

	latitude	longitude
time3
2018-01-11 07:26:30.393129984	-75.018896	165.499855
2018-01-11 07:48:05.826169856	-75.017433	165.500057
2018-01-11 07:48:06.832060160	-75.017432	165.500058

# Create geodataframe from gps_df
gps_gdf = gpd.GeoDataFrame(
    gps_df,
    geometry=gpd.points_from_xy(gps_df['longitude'], gps_df['latitude']), 
)

try:
    gps_gdf = gps_gdf.set_crs("EPSG:4326")
except Exception:
    gps_gdf = gps_gdf.set_crs("+proj=longlat +datum=WGS84 +no_defs +type=crs")

# Plot lat lon bounding box points
print("Minimum Longitude:", gps_df.longitude.min())
print("Maximum Longitude:", gps_df.longitude.max())
print("Minimum Latitude:", gps_df.latitude.min())
print("Maximum Latitude:", gps_df.latitude.max())

Minimum Longitude: 165.49765868921702
Maximum Longitude: 165.50005833333333
Minimum Latitude: -75.01889590423247
Maximum Latitude: -75.01737333333334

Using matplotlib and cartopy, plot the individual GPS points with a South Polar Stereo map projection.

# Plot on Polar Map
_, ax = plt.subplots(figsize=(5, 5), subplot_kw={"projection": ccrs.SouthPolarStereo()})
ax.set_extent([-180, 180, -90, -65], crs=ccrs.PlateCarree())
ax.coastlines(resolution='110m')
ax.stock_img()
ax.gridlines()
gps_gdf.plot(ax=ax, markersize=10, color='red', transform=ccrs.PlateCarree(), aspect=None)

<GeoAxes: >

_images/0b3106e51c23354434f6ea6c1f81834545682bab19265cd1fd88aac802c9e688.png

Remove background noise#

Following De Robertis and Higginbottom (2007), we can remove background noise by estimating it from the mean calibrated power of a collection of pings.

# Remove background noise
ds_Sv = ep.clean.remove_background_noise(
    ds_Sv,
    ping_num=5,
    range_sample_num=5,
    # Threshold value from https://github.com/a-sheehan/Echopype-Processing-Pipeline-for-AZFP-and-Glider-Data/blob/main/Sheehan_azfp_forloop.ipynb
    SNR_threshold="2.0dB",
)

# Plot corrected Sv
g = ds_Sv["Sv_corrected"].plot(
    col="channel", yincrease=False, y="range_sample", vmin=-100, vmax=-60, cmap="ep.ek500"
)

# Iterate through each facet
for ax in g.axes.flatten():
    labels = [label.get_text() for label in ax.get_xticklabels()]
    # Set the current ticks
    ax.set_xticks(ax.get_xticks())
    # Apply rotation
    ax.set_xticklabels(labels, rotation=45)
    # Set x label
    ax.set_xlabel("Ping Time (HH:MM)")

plt.show()

/tmp/ipykernel_250/3887182110.py:7: FutureWarning: self.axes is deprecated since 2022.11 in order to align with matplotlibs plt.subplots, use self.axs instead.
  for ax in g.axes.flatten():

_images/1faebaabbf553a6650e67caf3037ca78388b90bf215c9bac39bc09e6ca8e7a02.png

Since we’re most interested in 125 kHz for capturing zooplankton schools, let’s look at the Sv histograms (with and without background noise removed) for the 125 kHz channels:

# Extract and flatten data from 125 kHz channels
data_Sv = ds_Sv["Sv"].sel(channel="59006-125-2").values.ravel() # Flatten
data_Sv_corrected = ds_Sv["Sv_corrected"].sel(channel="59006-125-2").values.ravel() # Flatten

# Create a figure and axis
fig, ax = plt.subplots(figsize=(8, 5))

# Plot histograms
ax.hist(data_Sv, bins=100, color='blue', alpha=0.2, label='125 kHz Sv')
ax.hist(data_Sv_corrected, bins=100, color='red', alpha=0.2, label='Corrected 125 kHz Sv')

# Add labels, title, and legend
ax.set_xlabel("Volume backscattering strength (Sv re 1 m-1)[dB]")
ax.set_ylabel("Sv Value Frequency")
ax.set_title("Histogram Comparison between original Sv and background noise removed Sv")
ax.legend()

plt.show()

_images/c9a860d7bcc1aa53e8eac1fa4feb83e623d7228f04500acd7b6ea656699bd66e.png

From both the histogram and the echogram images above, we can see that a significant amount of Sv values, between -120 and -90 dB, surrounding the orange and yellow regions in the 125 kHz channel have been removed. Being that the majority of values removed were so weak, we suspect that it was mostly background noise that was removed.

Compute and apply shoal mask on corrected Sv#

With the calibrated backscatter volume data in an easily manipulable array format, here we use echopype’s shoal detection functionality to extract the zooplankton schools seen in the echogram above. There are multiple options to perform shoal detection, and below we use the one used by Echoview, with the implementation adapted from another open-source fisheries acoustics package, Echopy.

from echopype.mask import detect_shoal

shoal_mask = detect_shoal(
    ds_Sv,
    method="echoview",
    params={
        "var_name": "Sv_corrected",
        "channel": "59006-125-2",
        "idim": np.arange(len(ds_Sv["range_sample"])+1),
        "jdim": np.arange(len(ds_Sv["ping_time"])+1),
        "thr": -90,
        "mincan": (5, 5),
        "maxlink": (30, 2),
        "minsho": (30, 2),
    }
)

The mask is binary, where 1 is within a zooplankton school, and 0 otherwise:

# Plot the school mask
shoal_mask.plot(yincrease=False, y="range_sample", cmap="gray")

<matplotlib.collections.QuadMesh at 0x7f46e6d69010>

_images/0a4a5f6a3d573f075c1bcce055c0e7138df3306a44b2ac0e9f415e3b2a5b429e.png

# Apply shoal mask on the corrected Sv
ds_shoal_masked_Sv = ep.mask.apply_mask(ds_Sv, shoal_mask, var_name="Sv_corrected")

Plot the 125 kHz channel of the shoal masked Sv dataset using depth as a 2D coordinate grid:

# Create a single subplot
fig, ax = plt.subplots(1,1,figsize=(12,6))

# Plot shoal masked 38 kHz Sv channel as a 2D mesh with depth and ping time coords
pcolormesh = plt.pcolormesh(
    ds_shoal_masked_Sv["ping_time"].broadcast_like(ds_shoal_masked_Sv["range_sample"]).values.T,
    ds_shoal_masked_Sv["depth"].sel(channel="59006-125-2").values, 
    ds_shoal_masked_Sv["Sv_corrected"].sel(channel="59006-125-2").values,
    shading='auto',
    vmin=-100,
    vmax=-60,
    cmap="ep.ek500",
)
plt.gca().invert_yaxis()

# Add title and axis labels
plt.title("125 kHz Sv (with background noise removed and shoal mask applied)")
plt.xlabel("Ping Time (DD HH:MM)")
plt.ylabel("Depth (m)")

# Add colorbar
plt.colorbar(pcolormesh, label="Volume backscattering strength (Sv re 1 m-1)[dB]")

# Show plot
plt.show()

/tmp/ipykernel_250/2654804491.py:5: UserWarning: The input coordinates to pcolormesh are interpreted as cell centers, but are not monotonically increasing or decreasing. This may lead to incorrectly calculated cell edges, in which case, please supply explicit cell edges to pcolormesh.
  pcolormesh = plt.pcolormesh(

_images/2c1bc04c221167fc96ac4f8f18e92f12e394bd9febb30149be853c14318262e7.png

Package versions#

import datetime
print(f"echopype: {ep.__version__}, xarray: {xr.__version__}, geopandas: {gpd.__version__}")

print(f"\n{datetime.datetime.utcnow()} +00:00")

echopype: 0.11.1, xarray: 2026.2.0, geopandas: 1.1.2

2026-04-18 19:42:11.009823 +00:00

Identify zooplankton school from glider echosounder data with external location and environmental data incorporated

Contents

Identify zooplankton school from glider echosounder data with external location and environmental data incorporated#

Introduction#

Description#

Outline#

Running the notebook#

Note#

Path setup#

Data processing workflow#

Open and convert raw AZFP data#

Update EchoData `Platform` group with glider data and metadata#

Load external glider netCDF file#

Update EchoData `Top-level` and `Platform` group attributes#

Update latitude and longitude in EchoData `Platform` group#

Plot glider location#

Calibrate backscatter data with glider mean environmental data#

Remove background noise#

Align echosounder data according to depth#

What are `range_sample` and `echo_range`?#

Adding depth to Sv#

Compute and apply shoal mask on corrected Sv#

Package versions#

Identify zooplankton school from glider echosounder data with external location and environmental data incorporated

Contents

Identify zooplankton school from glider echosounder data with external location and environmental data incorporated#

Introduction#

Description#

Outline#

Running the notebook#

Note#

Path setup#

Data processing workflow#

Open and convert raw AZFP data#

Update EchoData Platform group with glider data and metadata#

Load external glider netCDF file#

Update EchoData Top-level and Platform group attributes#

Update latitude and longitude in EchoData Platform group#

Plot glider location#

Calibrate backscatter data with glider mean environmental data#

Remove background noise#

Align echosounder data according to depth#

What are range_sample and echo_range?#

Adding depth to Sv#

Compute and apply shoal mask on corrected Sv#

Package versions#

Update EchoData `Platform` group with glider data and metadata#

Update EchoData `Top-level` and `Platform` group attributes#

Update latitude and longitude in EchoData `Platform` group#

What are `range_sample` and `echo_range`?#