frequency_differencing#
- echopype.mask.frequency_differencing(source_Sv: Union[xarray.core.dataset.Dataset, str, pathlib.Path], storage_options: Optional[dict] = {}, freqAB: Optional[List[float]] = None, chanAB: Optional[List[str]] = None, operator: str = '>', diff: Optional[Union[float, int]] = None) xarray.core.dataarray.DataArray #
Create a mask based on the differences of Sv values using a pair of frequencies. This method is often referred to as the “frequency-differencing” or “dB-differencing” method.
- Parameters
- source_Sv: xr.Dataset or str or pathlib.Path
If a Dataset this value contains the Sv data to create a mask for, else it specifies the path to a zarr or netcdf file containing a Dataset. This input must correspond to a Dataset that has the coordinate
channel
and variablesfrequency_nominal
andSv
.- storage_options: dict, optional
Any additional parameters for the storage backend, corresponding to the path provided for
source_Sv
- freqAB: list of float, optional
The pair of nominal frequencies to be used for frequency-differencing, where the first element corresponds to
freqA
and the second element corresponds tofreqB
. Only one offreqAB
andchanAB
should be provided, and not both.- chanAB: list of strings, optional
The pair of channels that will be used to select the nominal frequencies to be used for frequency-differencing, where the first element corresponds to
freqA
and the second element corresponds tofreqB
. Only one offreqAB
andchanAB
should be provided, and not both.- operator: {“>”, “<”, “<=”, “>=”, “==”}
The operator for the frequency-differencing
- diff: float or int
The threshold of Sv difference between frequencies
- Returns
- xr.DataArray
A DataArray containing the mask for the Sv data. Regions satisfying the thresholding criteria are filled with
True
, else the regions are filled withFalse
.
- Raises
- ValueError
If neither
freqAB
orchanAB
are given- ValueError
If both
freqAB
andchanAB
are given- TypeError
If any input is not of the correct type
- ValueError
If either
freqAB
orchanAB
are provided and the list does not contain 2 distinct elements- ValueError
If
freqAB
contains values that are not contained infrequency_nominal
- ValueError
If
chanAB
contains values that not contained inchannel
- ValueError
If
operator
is not one of the following:">", "<", "<=", ">=", "=="
- ValueError
If the path provided for
source_Sv
is not a valid path- ValueError
If
freqAB
orchanAB
is provided and the Dataset produced bysource_Sv
does not contain the coordinatechannel
and variablefrequency_nominal
Notes
This function computes the frequency differencing as follows:
Sv_freqA - Sv_freqB operator diff
. Thus, ifoperator = "<"
anddiff = "5"
the following would be calculated:Sv_freqA - Sv_freqB < 5
.Examples
Compute frequency-differencing mask using a mock Dataset and channel selection:
>>> n = 5 # set the number of ping times and range samples ... >>> # create mock Sv data >>> Sv_da = xr.DataArray(data=np.stack([np.arange(n**2).reshape(n,n), np.identity(n)]), ... coords={"channel": ['chan1', 'chan2'], ... "ping_time": np.arange(n), "range_sample":np.arange(n)}) ... >>> # obtain mock frequency_nominal data >>> freq_nom = xr.DataArray(data=np.array([1.0, 2.0]), ... coords={"channel": ['chan1', 'chan2']}) ... >>> # construct mock Sv Dataset >>> Sv_ds = xr.Dataset(data_vars={"Sv": Sv_da, "frequency_nominal": freq_nom}) ... >>> # compute frequency-differencing mask using channel names >>> echopype.mask.frequency_differencing(source_Sv=mock_Sv_ds, storage_options={}, freqAB=None, ... chanAB = ['chan1', 'chan2'], ... operator = ">=", diff=10.0) <xarray.DataArray 'mask' (ping_time: 5, range_sample: 5)> array([[False, False, False, False, False], [False, False, False, False, False], [ True, True, True, True, True], [ True, True, True, True, True], [ True, True, True, True, True]]) Coordinates: * ping_time (ping_time) int64 0 1 2 3 4 * range_sample (range_sample) int64 0 1 2 3 4