The DOLfYN API

This is the DOLfYN API. It is a high-level object-oriented library composed of a set of data-object classes (types) that contain data from a particular measurement instrument, and a collection of functions that manipulate those data objects to accomplish data processing and data analysis tasks.

Contents

• The DOLfYN API

  • DOLfYN data objects

  • DOLfYN tools and functions

DOLfYN data objects

dolfyn.Velocity	This is the base class for velocity data objects.
dolfyn.ADPdata	The acoustic Doppler profiler (ADP) data type.
dolfyn.ADVdata	The acoustic Doppler velocimeter (ADV) data type.

class dolfyn.Velocity(*args, **kwargs)

This is the base class for velocity data objects.

All ADP and ADV data objects inherit from this base class.

See also

dict

Notes

DOLfYN Velocity objects are based on Python dicts, but have fancy interactive printing properties and indexing properties.

First, the interactive printing:

>>> import dolfyn as dlfn
>>> dat = dlfn.read_example('BenchFile01.ad2cp')

In an interactive interpreter, view the contents of the data object by:

>>> dat
<ADP data object>
. 9.11 minutes (started: Feb 24, 2017 10:01)
. BEAM-frame
. (38 bins, 1094 pings @ 2Hz)
*------------
| mpltime                  : <time_array; (1094,); float64>
| range                    : <array; (38,); float64>
| range_b5                 : <array; (38,); float64>
| vel                      : <array; (4, 38, 1094); float32>
| vel_b5                   : <array; (1, 38, 1094); float32>
+ alt                      : + DATA GROUP
+ altraw                   : + DATA GROUP
+ config                   : + DATA GROUP
+ env                      : + DATA GROUP
+ orient                   : + DATA GROUP
+ props                    : + DATA GROUP
+ signal                   : + DATA GROUP
+ sys                      : + DATA GROUP

You can view the contents of a ‘DATA GROUP’ by:

>>> dat['env']
<class 'dolfyn.data.base.TimeData'>: Data Object with Keys:
*------------
| c_sound                  : <array; (1094,); float32>
| press                    : <array; (1094,); float32>
| temp                     : <array; (1094,); float32>

Or you can also use attribute-style syntax:

>>> dat.signal
<class 'dolfyn.data.base.TimeData'>: Data Object with Keys:
*------------
| amp                      : <array; (4, 38, 1094); float16>
| amp_b5                   : <array; (1, 38, 1094); float16>
| corr                     : <array; (4, 38, 1094); uint8>
| corr_b5                  : <array; (1, 38, 1094); uint8>

You can directly access an item in a subgroup by:

>>> dat['env.c_sound']
array([1520.9   , 1520.8501, 1520.8501, ..., 1522.3   , 1522.3   ,
       1522.3   ], dtype=float32)

# And you can test for the presence of a variable by:

>>> 'signal.amp' in dat
True

property U

Horizontal velocity as a complex quantity.

property U_angle

Angle of horizontal velocity vector (radians clockwise from east/X/streamwise).

property U_mag

Horizontal velocity magnitude.

append(self, other)

Join two data objects together.

For example, two data objects d1 and d2 (which must contain the same variables, with the same array dimensions) can be joined together by:

>>> dat = d1.append(d2)

copy(self)

Create a copy of the data object.

get(self, key, default=None, /)

Return the value for key if key is in the dictionary, else default.

iter_data(self, include_hidden=False)

Generate the keys for all data items in this data object, including walking through sub-data objects.

Parameters

include_hiddenbool (Default: False)

Whether entries starting with ‘_’ should be included in the iteration.

iter_subgroups(self, include_hidden=False)

Generate the keys for all sub-groups in this data object, including walking through sub-groups.

Parameters

include_hiddenbool (Default: False)

Whether entries starting with ‘_’ should be included in the iteration.

keys()

property n_time

The number of timesteps in the data object.

pop(self, indx, d=<object object at 0x105342060>)

If key is not found, d is returned if given, otherwise KeyError is raised

rotate2(self, out_frame, inplace=False)

Rotate the data object into a new coordinate system.

Parameters

out_framestring {‘beam’, ‘inst’, ‘earth’, ‘principal’}

The coordinate system to rotate the data into.

inplacebool

Operate on self (True), or return a copy that has been rotated (False, default).

Returns

objoutVelocity

The rotated data object. This is self if inplace is True.

See also

dolfyn.rotate2()

set_declination(self, declination)

Set the declination of the data object.

Parameters

declinationfloat

The value of the magnetic declination in degrees (positive values specify that Magnetic North is clockwise from True North)

Notes

This method modifies the data object in the following ways:

• If the data-object is in the earth reference frame at the time of setting declination, it will be rotated into the “True-East, True-North, Up” (hereafter, ETU) coordinate system

• dat['orient']['orientmat'] is modified to be an ETU to instrument (XYZ) rotation matrix (rather than the magnetic-ENU to XYZ rotation matrix). Therefore, all rotations to/from the ‘earth’ frame will now be to/from this ETU coordinate system.

• The value of the specified declination will be stored in dat.props['declination']

• dat['orient']['heading'] is adjusted for declination (i.e., it is relative to True North).

• If dat['props']['principal_heading'] is set, it is adjusted to account for the orientation of the new ‘True’ earth coordinate system (i.e., calling set_declination on a data object in the principal coordinate system, then calling dat.rotate2(‘earth’) will yield a data object in the new ‘True’ earth coordinate system)

property shape

The shape of ‘scalar’ data in this data object.

property subset

Subset is an indexer for creating subsets of the data object using Python slice syntax.

For example, you can do:

dat2 = dat.subset[10:500]

Which returns a new data object that contains the 10th through 500th entry for each variable in the data object. This also operates recursively through sub-objects.

to_hdf5(self, buf, chunks=True, compression='gzip')

Write the data in this object to an hdf5 file.

property u

The first velocity component.

This is simply a shortcut to self[‘vel’][0]. Therefore, depending on the coordinate system of the data object (self[‘props’][‘coord_sys’]), it is:

• beam: beam1

• inst: x

• earth: east

• principal: streamwise

property v

The second velocity component.

This is simply a shortcut to self[‘vel’][1]. Therefore, depending on the coordinate system of the data object (self[‘props’][‘coord_sys’]), it is:

• beam: beam2

• inst: y

• earth: north

• principal: cross-stream

property w

The third velocity component.

This is simply a shortcut to self[‘vel’][2]. Therefore, depending on the coordinate system of the data object (self[‘props’][‘coord_sys’]), it is:

• beam: beam3

• inst: z

• earth: up

• principal: up

class dolfyn.ADPdata(*args, **kwargs)

The acoustic Doppler profiler (ADP) data type.

See also

dolfyn.Velocity

property S2

The horizontal shear-squared.

See also

dvdz(), dudz()

Notes

This is actually (dudz)^2 + (dvdz)^2. So, if those variables are not actually vertical derivatives of the horizontal velocity, then this is not the ‘horizontal shear-squared’.

property dudz

The shear in the first velocity component.

Notes

The derivative direction is along the profiler’s ‘z’ coordinate (‘dz’ is actually diff(self[‘range’])), not necessarily the ‘true vertical’ direction.

property dvdz

The shear in the second velocity component.

Notes

The derivative direction is along the profiler’s ‘z’ coordinate (‘dz’ is actually diff(self[‘range’])), not necessarily the ‘true vertical’ direction.

property dwdz

The shear in the third velocity component.

Notes

The derivative direction is along the profiler’s ‘z’ coordinate (‘dz’ is actually diff(self[‘range’])), not necessarily the ‘true vertical’ direction.

class dolfyn.ADVdata(*args, **kwargs)

The acoustic Doppler velocimeter (ADV) data type.

See also

dolfyn.Velocity

DOLfYN tools and functions

dolfyn.read	Read a binary Nortek (e.g., .VEC, .wpr, .ad2cp, etc.) or RDI (.000, .PD0, etc.) data file.
dolfyn.rotate2	Rotate a data object to a new coordinate system.
dolfyn.adv.motion.CorrectMotion	This object performs motion correction on an IMU-ADV data object.
dolfyn.VelBinner	This is the base binning (averaging) tool.
dolfyn.adv.turbulence.TurbBinner	Computes various averages and turbulence statistics from cleaned ADV data.
dolfyn.calc_principal_heading	Compute the principal angle of the horizontal velocity.

dolfyn.read(fname, userdata=True, nens=None)

Read a binary Nortek (e.g., .VEC, .wpr, .ad2cp, etc.) or RDI (.000, .PD0, etc.) data file.

Parameters

filenamestring

Filename of Nortek file to read.

userdataTrue, False, or string of userdata.json filename

(default True) Whether to read the ‘<base-filename>.userdata.json’ file.

nensNone (default: read entire file), int, or

2-element tuple (start, stop)

Number of pings to read from the file

Returns

dat<~dolfyn.data.velocity.Velocity>

A DOLfYN velocity data object.

dolfyn.rotate2(obj, out_frame='earth', inplace=False)

Rotate a data object to a new coordinate system.

Parameters

objVelocity

The dolfyn velocity-data (ADV or ADP) object to rotate.

out_framestring {‘beam’, ‘inst’, ‘earth’, ‘principal’}

The coordinate system to rotate the data into.

inplacebool

Operate on the input data object (True), or return a copy that has been rotated (False, default).

Returns

objoutVelocity

The rotated data object. Note that when inplace=True, the input object is modified in-place and returned (i.e., objout is obj).

Notes

This function rotates all variables in obj.props['rotate_vars'].

class dolfyn.adv.motion.CorrectMotion(accel_filtfreq=None, vel_filtfreq=None, separate_probes=False)

This object performs motion correction on an IMU-ADV data object. The IMU and ADV data should be tightly synchronized and contained in a single data object.

Parameters

accel_filtfreqfloat

the frequency at which to high-pass filter the acceleration signal to remove low-frequency drift.

vel_filtfreqfloat (optional)

a second frequency to high-pass filter the integrated acceleration. (default: 1/3 of accel_filtfreq)

separate_probesbool (optional: False)

a flag to perform motion-correction at the probe tips, and perform motion correction in beam-coordinates, then transform back into XYZ/earth coordinates. This correction seems to be lower than the noise levels of the ADV, so the defualt is to not use it (False).

Notes

Acceleration signals from inertial sensors are notorious for having a small bias that can drift slowly in time. When integrating these signals to estimate velocity the bias is amplified and leads to large errors in the estimated velocity. There are two methods for removing these errors,

1. high-pass filter the acceleration signal prior and/or after integrating. This implicitly assumes that the low-frequency translational velocity is zero.

2. provide a slowly-varying reference position (often from a GPS) to an IMU that can use the signal (usually using Kalman filters) to debias the acceleration signal.

Because method (1) removes real low-frequency acceleration, method (2) is more accurate. However, providing reference position estimates to undersea instruments is practically challenging and expensive. Therefore, lacking the ability to use method (2), this function utilizes method (1).

For deployments in which the ADV is mounted on a mooring, or other semi-fixed structure, the assumption of zero low-frequency translational velocity is a reasonable one. However, for deployments on ships, gliders, or other moving objects it is not. The measured velocity, after motion-correction, will still hold some of this contamination and will be a sum of the ADV motion and the measured velocity on long time scales. If low-frequency motion is known separate from the ADV (e.g. from a bottom-tracking ADP, or from a ship’s GPS), it may be possible to remove that signal from the ADV signal in post-processing. The accuracy of this approach has not, to my knowledge, been tested yet.

Examples

>>> from dolfyn.adv import api as avm
>>> dat = avm.read_nortek('my_data_file.vec')
>>> mc = avm.CorrectMotion(0.1)
>>> corrected_data = mc(dat)

__call__(self, advo, to_earth=True)

Perform motion correction on an IMU-equipped ADV object.

Parameters

advoADVdata

The adv object on which to perform motion correction. It must contain the following data attributes:

• vel : The velocity array.

• accel : The translational acceleration array.

• angrt : The rotation-rate array.

• orientmat : The orientation matrix.

• props : a dictionary that has ‘inst2head_vec’ and ‘coord_sys’.

to_earthbool (optional, default: True)

A boolean that specifies whether the data should be rotated into the earth frame.

Notes

After calling this function, advo will have velrot and velacc data attributes. The velocity vector attribute vel will be motion corrected according to:

u_corr = u_raw + velacc + velrot

Therefore, to recover the ‘raw’ velocity, subtract velacc and velrot from vel.

This method does not return a data object, it operates on (motion corrects) the input advo.

class dolfyn.VelBinner(n_bin, fs, n_fft=None, n_fft_coh=None)

This is the base binning (averaging) tool.

All DOLfYN binning tools derive from this base class.

Examples

The VelBinner class is used to compute averages and turbulence statistics from ‘raw’ (unaveraged) ADV or ADP measurements, for example:

# First read or load some data.
rawdat = dlfn.read_example('BenchFile01.ad2cp')

# Now initialize the averaging tool:
binner = dlfn.VelBinner(n_bin=600, fs=rawdat['props']['fs'])

# This computes the basic averages
avg = binner(rawdat)

calc_acov(self, indat, n_bin=None)

Calculate the auto-covariance of the raw-signal indat.

As opposed to calc_xcov, which returns the full cross-covariance between two arrays, this function only returns a quarter of the full auto-covariance. It computes the auto-covariance over half of the range, then averages the two sides (to return a ‘quartered’ covariance).

This has the advantage that the 0 index is actually zero-lag.

calc_omega(self, fs=None, coh=False)

Calculate the radial-frequency vector for the psd’s.

Parameters

fsfloat (optional)

The sample rate (Hz).

cohbool

Calculate the frequency vector for coherence/cross-spectra (default: False) i.e. use self.n_fft_coh instead of self.n_fft.

calc_stress(self, veldat)

Calculate the stresses (cross-covariances of u,v,w).

Parameters

veldata velocity data array. The last dimension is assumed

to be time.

Returns

outAn array of stress values.

calc_tke(self, veldat, noise=[0, 0, 0])

Calculate the tke (variances of u,v,w).

Parameters

veldata velocity data array. The last dimension is assumed

to be time.

noisea three-element vector of the noise levels of the

velocity data for ach component of velocity.

Returns

outAn array of tke values.

calc_vel_cpsd(self, veldat, rotate_u=False, fs=None, window='hann', n_bin=None, n_fft=None, n_pad=None, step=None)

Calculate the cross-spectra of velocity components.

Parameters

veldatnp.ndarray

The raw velocity data.

rotate_ubool (optional)

If True, each ‘bin’ of horizontal velocity is rotated into its principal axis prior to calculating the psd. (default: False).

fsfloat (optional)

The sample rate (default: from the binner).

windowstring or array

Specify the window function.

n_binint (optional)

The bin-size (default: from the binner).

n_fftint (optional)

The fft size (default: from the binner).

n_padint (optional)

The number of values to pad with zero (default: 0)

stepint (optional)

Controls amount of overlap in fft (default: the step size is chosen to maximize data use, minimize nens, and have a minimum of 50% overlap.).

Returns

CSpecnp.ndarray (3, M, N_FFT)

The first-dimension of the cross-spectrum is the three different cross-spectra: ‘uv’, ‘uw’, ‘vw’ (in that order).

calc_vel_psd(self, veldat, rotate_u=False, fs=None, window='hann', noise=[0, 0, 0], n_bin=None, n_fft=None, n_pad=None, step=None)

Calculate the psd of velocity.

Parameters

veldatnp.ndarray

The raw velocity data.

rotate_ubool (optional)

If True, each ‘bin’ of horizontal velocity is rotated into its principal axis prior to calculating the psd. (default: False).

fsfloat (optional)

The sample rate (default: from the binner).

windowstring or array

Specify the window function.

noiselist(3 floats) (optional)

Noise level of each component’s velocity measurement (default 0).

n_binint (optional)

The bin-size (default: from the binner).

n_fftint (optional)

The fft size (default: from the binner).

n_padint (optional)

The number of values to pad with zero (default: 0)

stepint (optional)

Controls amount of overlap in fft (default: the step size is chosen to maximize data use, minimize nens, and have a minimum of 50% overlap.).

Returns

Specnp.ndarray (3, M, N_FFT)

The first-dimension of the spectrum is the three different spectra: ‘uu’, ‘vv’, ‘ww’.

calc_xcov(self, indt1, indt2, npt=None, n_bin1=None, n_bin2=None, normed=False)

Calculate the cross-covariance between arrays indt1 and indt2 for each bin.

cohere(self, dat1, dat2, window='hann', debias=True, noise=0, 0, n_fft=None, n_bin1=None, n_bin2=None)

Calculate coherence between dat1 and dat2.

cpsd(self, dat1, dat2, fs=None, window='hann', n_fft=None, n_bin1=None, n_bin2=None)

Calculate the ‘cross power spectral density’ of dat.

Parameters

dat1np.ndarray

The first raw-data array of which to calculate the cpsd.

dat2np.ndarray

The second raw-data array of which to calculate the cpsd.

windowstring

String indicating the window function to use (default: ‘hanning’).

Returns

outnp.ndarray

The cross-spectral density of dat1 and dat2

demean(self, dat, n_pad=0, n_bin=None)

Reshape the array dat and remove the mean from each ensemble.

detrend(self, dat, n_pad=0, n_bin=None)

Reshape the array dat and remove the best-fit trend line.

… Need to fix this to deal with NaNs…

do_avg(self, rawdat, outdat=None, names=None, n_time=None)

Parameters

rawdatraw_data_object

The raw data structure to be binned

outdatavg_data_object

The bin’d (output) data object to which averaged data is added.

nameslist of strings

The names of variables to be averaged. If names is None, all data in rawdat will be binned.

do_var(self, rawdat, outdat=None, names=None, suffix='_var')

Calculate the variance of data attributes.

Parameters

rawdatraw_data_object

The raw data structure to be binned.

outdatavg_data_object

The bin’d (output) data object to which variance data is added.

nameslist of strings

The names of variables of which to calculate variance. If names is None, all data in rawdat will be binned.

mean(self, dat, axis=- 1, n_bin=None, mask_thresh=None)

Average an array object.

Parameters

n_binint (default is self.n_bin)

mask_threshfloat (between 0 and 1)

if the input data is a masked array, and mask_thresh is not None mask the averaged values where the fraction of bad points is greater than mask_thresh

mean_angle(self, dat, axis=- 1, units='radians', n_bin=None, mask_thresh=None)

Average an angle array.

Parameters

units{‘radians’ | ‘degrees’}

n_binint (default is self.n_bin)

mask_threshfloat (between 0 and 1)

if the input data is a masked array, and mask_thresh is not None mask the averaged values where the fraction of bad points is greater than mask_thresh

phase_angle(self, dat1, dat2, window='hann', n_fft=None, n_bin1=None, n_bin2=None)

Calculate the phase difference between two signals as a function of frequency (complimentary to coherence).

Parameters

dat1np.ndarray

The first raw-data array of which to calculate the cpsd.

dat2np.ndarray

The second raw-data array of which to calculate the cpsd.

windowstring

String indicating the window function to use (default: ‘hanning’).

Returns

outnp.ndarray

The phase difference between signal dat1 and dat2.

psd(self, dat, fs=None, window='hann', noise=0, n_bin=None, n_fft=None, step=None, n_pad=None)

Calculate ‘power spectral density’ of dat.

Parameters

datdata_object

The raw-data array of which to calculate the psd.

windowstring

String indicating the window function to use (default: ‘hanning’).

noisefloat

The white-noise level of the measurement (in the same units as dat).

reshape(self, arr, n_pad=0, n_bin=None)

Reshape the array arr to shape (…,n,n_bin+n_pad).

Parameters

arrnp.ndarray

n_padint

Is used to add n_pad/2 points from the end of the previous ensemble to the top of the current, and n_pad/2 points from the top of the next ensemble to the bottom of the current. Zeros are padded in the upper-left and lower-right corners of the matrix (beginning/end of timeseries). In this case, the array shape will be (…,`n`,`n_pad`+`n_bin`)

n_binfloat, int (optional)

Override this binner’s n_bin.

Notes

n_bin can be non-integer, in which case the output array size will be n_pad`+`n_bin, and the decimal will cause skipping of some data points in arr. In particular, every mod(n_bin,1) bins will have a skipped point. For example: - for n_bin=2048.2 every 1/5 bins will have a skipped point. - for n_bin=4096.9 every 9/10 bins will have a skipped point.

class dolfyn.adv.turbulence.TurbBinner(n_bin, fs, n_fft=None, n_fft_coh=None)

Computes various averages and turbulence statistics from cleaned ADV data.

Parameters

n_binint

The length of bin s, in number of points, for this averaging operator.

n_fftint (optional, default: n_fft = n_bin)

The length of the FFT for computing spectra (must be < n_bin)

__call__(self, advr, out_type=None, omega_range_epsilon=[6.28, 12.57], Itke_thresh=0, window='hann')

Compute a suite of turbulence statistics for the input data advr, and return a binned data object.

Parameters

advrADVdata

The raw adv data-object to bin, average and compute turbulence statistics of.

omega_range_epsiloniterable(2)

The frequency range (low, high) over which to estimate the dissipation rate epsilon [rad/s].

Itke_threshThe threshold for velocity magnitude for

computing the turbulence intensity. Values of Itke where U_mag < Itke_thresh are set to NaN. (default: 0).

window1, None, ‘hann’

The window to use for psds.

Returns

advbbase.ADVbinned

Returns an ‘binned’ (i.e. ‘averaged’) data object. All fields of the input data object are averaged in n_bin chunks. This object also computes the following items over those chunks:

• tke_vec : The energy in each component (components are also accessible as upup_, vpvp_, wpwp_)

• stress : The Reynolds stresses (each component is accessible as upwp_, vpwp_, upvp_)

• sigma_Uh : The standard deviation of the horizontal velocity.

• Spec : A data group containing the spectra of the velocity in radial frequency units. The data group contains:

  • vel : the velocity spectra array (m^2/s/radian))

  • omega : the radial frequncy (radian / s)

calc_Lint(self, corr_vel, U_mag, fs=None)

Calculate integral length scales.

Parameters

corr_velnumpy.ndarray

The auto-covariance array (i.e. computed using calc_acov).

U_magnumpy.ndarray (…, n_time)

The velocity magnitude for this bin.

fsfloat

The raw sample rate.

Returns

Lintnumpy.ndarray (…, n_time)

The integral length scale (Tint*U_mag).

Notes

The integral time scale (Tint) is the lag-time at which the auto-covariance falls to 1/e.

calc_epsilon_LT83(self, spec, omega, U_mag, omega_range=[6.28, 12.57])

Calculate the dissipation rate from the spectrum.

Parameters

specnumpy.ndarray (…,n_time,n_f)

The spectrum array [m^2/s/radian]

omeganumpy.ndarray (n_f)

The frequency array [rad/s]

U_magnumpy.ndarray (…,n_time)

The velocity magnitude [m/s]

omega_rangeiterable(2)

The range over which to integrate/average the spectrum.

Returns

epsilonnp.ndarray (…,n_time)

The dissipation rate.

Notes

This uses the standard formula for dissipation:

\[S(k) = \alpha \epsilon^{2/3} k^{-5/3}\]

where \(\alpha = 0.5\) (1.5 for all three velocity components), k is wavenumber and S(k) is the turbulent kinetic energy spectrum.

With \(k \rightarrow \omega / U\) then–to preserve variance– \(S(k) = U S(\omega)\) and so this becomes:

\[S(\omega) = \alpha \epsilon^{2/3} \omega^{-5/3} U^{2/3}\]

LT83 : Lumley and Terray “Kinematics of turbulence convected by a random wave field” JPO, 1983, 13, 2000-2007.

calc_epsilon_SF(self, veldat, umag, fs=None, freq_rng=[0.5, 5.0])

Calculate epsilon using the “structure function” (SF) method.

Parameters

veldatnumpy.ndarray (…, n_time, n_bin)

The raw velocity signal (last dimension time) upon which to perform the SF technique.

umagnumpy.ndarray (…, n_time)

The bin-averaged horizontal velocity magnitude.

fsfloat

The sample rate of veldat [hz].

freq_rngiterable(2)

The frequency range over which to compute the SF [hz].

Returns

epsilonnumpy.ndarray (…, n_time)

The dissipation rate.

calc_epsilon_TE01(self, advbin, advraw, omega_range=[6.28, 12.57])

Calculate the dissipation according to TE01.

Parameters

advbinADVbinned

The binned adv object. The following spectra and basic turbulence statistics must already be computed.

advrawADVdata

The raw adv object.

Notes

TE01Trowbridge, J and Elgar, S, “Turbulence measurements in

the Surf Zone” JPO, 2001, 31, 2403-2417.

up_angle(self, Uh_complex)

Calculate the angle of the turbulence fluctuations.

Parameters

Uh_complexnumpy.ndarray (…, n_time * n_bin)

The complex, raw horizontal velocity (non-binned).

Returns

thetanumpy.ndarray (…, n_time)

The angle of the turbulence [rad]

dolfyn.calc_principal_heading(vel, tidal_mode=True)

Compute the principal angle of the horizontal velocity.

Parameters

velnp.ndarray (2,…,Nt), or (3,…,Nt)

The 2D or 3D velocity array (3rd-dim is ignored in this calculation)

tidal_modebool (default: True)

Returns

p_headingfloat or ndarray

The principal heading(s) in degrees clockwise from North.

Notes

The tidal mode follows these steps:

1. rotates vectors with negative v by 180 degrees

2. then doubles those angles to make a complete circle again

3. computes a mean direction from this, and halves that angle again.

4. The returned angle is forced to be between 0 and 180. So, you may need to add 180 to this if you want your positive direction to be in the western-half of the plane.

Otherwise, this function simply compute the average direction using a vector method.