Multibeam fish detection filter

The following algorithm attempts to identify fish in a multibeam echogram as a Boolean result. It generates a bitmap to filter out background signal in multibeam and imaging sonars. It is inspired by the filter described in Balk et al. 2009.

The Window size (pings) setting governs the background noise estimation. In this example we set it to 17.

Use a Mask operator in conjunction with this filter and the original data. Then use the Multibeam Target Detection operator to detect fish targets.

List of operands

Operand 1: a Sound Metrics DIDSON or ARIS multibeam variable, or other multibeam variables.

Code

"""
Echoview Code Operator source file
========================================================================

Created by Echoview(R) 10.0.218 on Monday, 29 April 2019

See the Echoview help file for Code operator documentation and
examples.

NumPy User Guide: https://docs.scipy.org/doc/numpy/user/
NumPy Reference: https://docs.scipy.org/doc/numpy/Reference/
SciPy Reference Guide: https://docs.scipy.org/doc/scipy/Reference/

Echoview(R) is a registered trademark of Echoview Software Pty Ltd.
"""

# Authorship information
__author__ = "Echoview Software Pty Ltd. 2019."
__disclaimer__ = (
    "This example code is provided AS IS, without warranty of any "
    "kind, express or implied, including but not limited to the "
    "warranties of merchantability, fitness for a particular purpose "
    "and noninfringement. In no event shall Echoview Software Pty Ltd "
    "be liable for any claim, damages or other liability, arising "
    "from, out of or in connection with the use of this example code."
)
__version__ = "1.0"

# System Imports
from typing import List

# Libraries
from echoview import OperatorBase, MeasurementType, OperandInput, Error
import numpy as np
from scipy import signal, ndimage


class Operator(OperatorBase):
    """
    Fish detection filter (background suppression) for DIDSON data
    ====================================================================

    Generates a bitmap representing potential fish detections.

    Operands
    ---------------------

    * Operand 1 - DIDSON Multibeam Sv
    """

    def result_type(self, input_types):
        return MeasurementType.MULTIBEAM_BOOLEAN

    def __init__(self):
        self.foreground_filter = np.full((7, 3, 3), 1.)  # boxcar (range, beams, pings)
        self.foreground_filter /= np.sum(self.foreground_filter)

        self.threshold = 2.2  # Filter threshold in dB

    def eval(self, inputs: List[OperandInput]):
        input = inputs[0]
        data = np.concatenate([ping.data[..., np.newaxis]
                               for ping in input.window_measurements], axis=2)

        # Generate a mask that reflects all valid values (i.e. not
        # nodata), this is then used to convert output to nodata
        # Note: That as signal.convolve doesn't support nodata values we
        # mask them to -999db (0 linearly)
        mask = np.isfinite(input.measurement.data)
        if input.measurement.type.is_db:
            data[~mask] = -999  # 0 linearly
            data = to_linear(data)
        else:
            data[~mask] = 0

        # smooth data using the foreground_filter, this cuts out
        # the highest frequency noise and is considered the foreground
        foreground_linear = signal.convolve(data, self.foreground_filter, mode='same')
        foreground = to_db(foreground_linear[..., input.window_index])

        # calculate the background signal, this is taken to be mean over
        # the ping/time axis
        background = to_db(np.mean(data, axis=(2)))

        # calculate the filtered signal, defined as the difference
        # between the foreground and background signals plus some offset
        filtered_data = foreground - background > self.threshold

        # no data values go to False since Echoview doesn't support
        # nodata for boolean data
        filtered_data[~mask] = False

        return filtered_data


# helper functions

THRESHOLD_LOWER_DB = -999
THRESHOLD_UPPER_DB = 999
THRESHOLD_UPPER_LINEAR = 10 ** (THRESHOLD_UPPER_DB / 10)


def to_linear(data: np.ma.MaskedArray, out=None):
    """Converts the specified input to linear values"""

    # Equivalent to ``10 ** (data / 10)`` but more efficient due to
    # reusing same output array

    if out is not None and out.shape != data.shape:
        raise ValueError(
            'out must be a NumPy array with the same shape as the input array')

    below_threshold = data <= THRESHOLD_LOWER_DB
    above_upper_threshold = data >= THRESHOLD_UPPER_DB
    out = np.divide(data, 10, out=out)
    out = np.power(10, out, out=out)

    out[below_threshold] = 0.
    out[above_upper_threshold] = THRESHOLD_UPPER_LINEAR

    return out


def to_db(data: np.ma.MaskedArray):
    """Converts the specified linear input to dB values"""

    thresholded_values = data <= 0
    out = 10 * np.ma.log10(data)
    out[thresholded_values] = THRESHOLD_LOWER_DB

    return out

Multibeam fish detection filter

List of operands

Code

See also