Download

Oceanlyz Matlab Toolbox

Oceanlyz User Manual

Download
oceanlyz_user_manual_ver1.4.pdf
Adobe Acrobat Document 537.0 KB


Ocean Wave Data Analysis:

Introduction to Time Series Analysis, Signal Processing, and Wave Prediction

Download
Contents.pdf
Adobe Acrobat Document 238.4 KB


License Agreement and Disclaimer

 

 

 

Copyright (c) 2018 Arash Karimpour

 

 

 

Permission is hereby granted, free of charge, to any person obtaining a copy

 

of this software and associated documentation files (the "Software"), to deal

 

in the Software without restriction, including without limitation the rights

 

to use, copy, modify, merge, publish, distribute, sublicense, and/or sell

 

copies of the Software, and to permit persons to whom the Software is

 

furnished to do so, subject to the following conditions:

 

 

 

The above copyright notice and this permission notice shall be included in all

 

copies or substantial portions of the Software.

 

 

 

THE SOFTWARE 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 THE

 

AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER

 

LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,

 

OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE

 

SOFTWARE.

 

 

 

 

 

Citation

 

 

 

Cite this toolbox as:

 

 

 

Karimpour, A., & Chen, Q. (2017). Wind Wave Analysis in Depth Limited Water Using OCEANLYZ, a MATLAB toolbox. Computers & Geosciences.

 

 

Recommended Book

 

For more details on coastal and ocean wave data analysis refer to:

 

Karimpour A., (2018), Ocean Wave Data Analysis: Introduction to Time Series Analysis, Signal Processing, and Wave Prediction, KDP.

 

Book link:

www.amazon.com/dp/0692109978

 

 

 

 

Introduction

 

 

 

OCEANLYZ, Ocean Wave Analyzing Toolbox, is a toolbox for analyzing the wave time series data collected by sensors in open body of water such as ocean, sea, and lake or in a laboratory.

 

 

 

This toolbox contains functions that each one is suitable for a particular purpose. Both spectral and zero-crossing methods are offered for wave analysis. This toolbox can calculate wave properties such as zero-moment wave height, significant wave height, mean wave height, peak wave period and mean period. This toolbox can correct and account for the pressure attention (pressure loss) in the water column for data collected by a pressure sensor. This toolbox can separate wind sea and swell energies and reports their properties.

 

 

 

 

 

Download

 

 

 

To download Oceanlyz, visit www.arashkarimpour.com

 

 

 

 

 

Installation

 

 

 

To use this toolbox, download and unzip it in any location you choose such as “C:\”.

 

 

 

 

 

Operating System

 

 

 

This code can be run on Windows, Mac and Linux. However, make sure any given path is compatible with a running operating system. In particular, “\” is used in Windows path, while “/” is used in Mac or Linux path. For example, if a path is “C:\” on Windows machine, it would be “C:/” on Mac or Linux.

 

 

 

 

 

Required Programing Language

 

 

 

This toolbox can be run by using MATLAB (www.mathworks.com) or GNU Octave (www.gnu.org/software/octave).

 

 

 

 

 

Required Package for MATLAB

 

 

 

MATLAB users need to install MATLAB Signal Processing Toolbox for running the Oceanlyz spectral analysis. It gives Oceanlyz access to MATLAB Welch's power spectral density calculation. However, MATLAB Signal Processing Toolbox it is not required for zero-crossing analysis.

 

 

 

 

 

Required Package for GNU Octave

 

 

 

GNU Octave users need to install/load GNU Octave Signal Package for running the Oceanlyz spectral analysis. It gives Oceanlyz access to GNU Octave Welch's power spectral density calculation. However, GNU Octave Signal Package it is not required for zero-crossing analysis.

 

 

 

GNU Octave Signal Package can be loaded inside GNU Octave by using a following command in a command window (This should be done every time GNU Octave is opened):

 

 

 

>> pkg load signal

 

 

 

If GNU Octave Signal Package is not already installed, it should be first installed from Octave Forge (octave.sourceforge.io), and then get loaded by using the following commands in a command window:

 

 

 

>> pkg install -forge signal

 

>> pkg load signal

 

 

 

 

 

Quick Start

 

 

 

·       Open MATLAB or GNU Octave

 

·       Change a current folder (current directory) to a folder that contains OCEANLYZ toolbox, for example “C:\oceanlyz”, in MATLAB or GNU Octave.

 

·       Open a file named “oceanlyzinput.m” in MATLAB or GNU Octave editor and modify it based on the properties of the collected dataset and required analysis.

 

·       Run a file named “RunOceanlyz.m” in MATLAB or GNU Octave to start calculations.

 

 

 

For more information, refer to “Tutorial” section.

 

 

 

 


 

 

Tutorial

 

 

 

In order to use Oceanlyz toolbox, first, one of the MATLAB or GNU Octave programming language should be installed with required packages (Refer to required packages for MATLAB/GNU Octave).

 

 

 

Then, the Oceanlyz toolbox should be downloaded and unzip in a folder of your choice. For example, if the Oceanlyz is unzipped on C: derive, then Oceanlyz folder would be “C:\oceanlyz” on a windows machine.

 

 

 

Before using the Oceanlyz, these steps should be followed first:

 

 

 

1.       Setting values representing the properties of the dataset needs to be analyzed

 

2.       Defining what type of analysis is required to be done

 

3.       Preparing dataset in a format that the Oceanlyz can use

 

 

 

These steps are explained in detail in following sections.

 

 

 

 

 

Set up Input, Output, and Calculation parameters

 

 

 

Values that define input, output, and calculation parameters are defined in a file named “oceanlyzinput.m”. To set up these parameters, first, open file “oceanlyzinput.m” by following these steps:

 

 

 

·       Open MATLAB or GNU Octave

 

·       Change a current folder (current directory) to a folder that contains OCEANLYZ toolbox in MATLAB or GNU Octave. For example, open “C:\oceanlyz” directory.

 

·       Open a file named “oceanlyzinput.m” in MATLAB or GNU Octave editor and modify it based on the properties of the collected dataset and required analysis.

 

 

 

Parameters that can be defined in “oceanlyzinput.m” are described in following.

 

 

 

 

 

Input, Output, and Calculation Parameters

 

 

 

InputFileName

 

It defines a name of the input data file

 

 

 

Notes:

 

Enter name of the file contains water depth or surface elevation (Eta) time series

 

Data should be in a single column vector file without any text, each burst of data should follow the previous burst

 

Most file formats such as '.mat', '.txt', '.csv' can be imported

 

Example:

 

InputFileName='waterdepthsample.mat';

 

 

 

 

 

InputFileFolder

 

It defines a location of the input data file

 

 

 

Notes:

 

Enter a location of the file contains water depth or surface elevation (Eta) time series

 

For Linux or Mac, use '/' instead of '\'

 

Example:

 

InputFileFolder='C:\oceanlyz\Sample';

 

 

 

SaveOutput

 

It defines if to save output(s) and a location the output(s) to be saved

 

 

 

SaveOutput='off';     :Does not save output(s)

 

SaveOutput='on';     :Save output(s)

 

 

 

OutputFileFolder

 

Location of the output file to be saved

 

 

 

OutputFileFolder='C:\oceanlyz';

 

 

 

Notes:

 

For Linux or Mac, use '/' instead of '\'

 

Output is a structure array

 

Name of output file is 'wave.mat'

 

After wave.mat is loaded in Matlab, field(s) in structure array can be called by using '.'

 

Example: Output peak wave period is : "wave.Tp"

 

 

 

AnalysisMethod

 

It defines analysis method

 

 

 

AnalysisMethod='spectral';     :Use spectral analysis method / Fast Fourier Transform

 

AnalysisMethod='zerocross';     :Use zero-crossing method

 

 

 

WaveParameterCalc

 

It defines if to calculate wave parameters or not

 

 

 

WaveParameterCalc='off';     :Does not calculate/report wave properties

 

WaveParameterCalc='on';     :Calculates/reports wave properties

 

 

 

Notes:

 

Use WaveParameterCalc='off' to only report corrected water level data measured by a pressure sensor

 

If WaveParameterCalc='off', then pressureattenuation should be pressureattenuation='on' or 'all'

 

If WaveParameterCalc='off', it reports corrected water level data measured by a pressure sensor

 

by accounting for pressure attenuation (pressure loss) in depth

 

If WaveParameterCalc='on' and pressureattenuation='on' or 'all', it reports both waves and corrected water level data measured by a pressure sensor

 

 

 

 

 

burst

 

it defines number of burst(s) in the input file

 

 

 

burst=5;     :Number of burst(s) in the input file

 

 

 

Notes:

 

burst=total number of data points/(duration*fs)

 

 

 

Example:

 

For 12 bursts of data, which each burst has a duration of 30 minutes, and collected at sampling frequency of 10 Hz:

 

duration=(30*60)

 

total number of data points=number of burst*duration of each burst*sampling frequency=12*(30*60)*10

 

 

 

duration

 

it defines duration time that data collected in each burst in (second)

 

 

 

duration=1024;     :Duration time that data collected in each burst in (second)

 

 

 

fs

 

It defines sampling frequency that data are collected at in (Hz)

 

 

 

fs=10;     :Sampling frequency that data are collected at in (Hz)

 

 

 

heightfrombed

 

It defines pressure sensor height from a bed in (m)

 

 

 

heightfrombed=0.05;     :Pressure sensor height from a bed in (m)

 

 

 

Notes:

 

Leave heightfrombed=0.0; if data is not measured by a pressure sensor or if a sensor sits on the seabed

 

 

 

nfft

 

It defines NFFT for Fast Fourier Transform

 

 

 

nfft=2^10;     :NFFT for Fast Fourier Transform

 

 

 

Notes:

 

Results will be reported for frequency range of 0 <= f <= (fs/2) with (nfft/2+1) data points

 

 

 

Example:

 

If fs=4 Hz and nfft=1024, then output frequency has a range of 0 <= f <= 2 with 513 data points

 

 

 

seaswellCalc

 

it defines if to separate wind sea and swell waves or not

 

 

 

seaswellCalc='off';     :Does not separate wind sea and swell waves

 

seaswellCalc='on';     :Separates wind sea and swell waves

 

 

 

fminswell

 

It defines minimum frequency that swell can have

 

 

 

fminswell=0.1;     :Minimum frequency that swell can have (it is used for Tpswell calculation) in (Hz)

 

 

 

fmaxswell

 

It defines maximum frequency that swell can have

 

 

 

fmaxswell=0.25;     :Maximum frequency that swell can have (It is about 0.2 in Gulf of Mexico) in (Hz)

 

 

 

pressureattenuation

 

It defines if to apply pressure attenuation factor or not

 

 

 

pressureattenuation='off';     :No pressure attenuation is applied

 

pressureattenuation='on';     :Pressure attenuation is applied without correction after fmaxpcorr

 

pressureattenuation='all';     :Pressure attenuation is applied with constant correction after fmaxpcorr

 

 

 

Notes:

 

Pressure attenuation factor is used to account for pressure attenuation (pressure loss) in depth

 

For pressureattenuation='on' or 'all', input data should be water depth

 

 

 

autofmaxpcorr

 

It defines if to calculate fmaxpcorr and ftailcorrection based on water depth or not

 

 

 

autofmaxpcorr='off':     :Does not calculate fmaxpcorr and ftailcorrection based on water depth

 

autofmaxpcorr='on':     :Calculate fmaxpcorr and ftailcorrection based on water depth

 

 

 

Notes:

 

Code calculate a maximum frequency that a pressure attenuation factor should be applied up to that

 

 

 

fminpcorr

 

It defines minimum frequency that automated calculated fmaxpcorr can have if autofmaxpcorr='on' in (Hz)

 

 

 

fminpcorr=0.15;     :Minimum frequency that automated calculated fmaxpcorr can have if autofmaxpcorr='on' in (Hz)

 

 

 

Notes:

 

If autofmaxpcorr='on', then fmaxpcorr will be checked to be larger or equal to fminpcorr

 

 

 

fmaxpcorr

 

It defines maximum frequency for applying pressure attenuation factor in (Hz)

 

 

 

fmaxpcorr=0.55;     :Maximum frequency for applying pressure attenuation factor in (Hz)

 

 

 

 

 

Notes:

 

Pressure attenuation factor is not applied on frequency larger than fmaxpcorr

 

 

 

mincutoff

 

It defines if to cut off the spectrum below fmin, i.e. where f<fmin, or not

 

 

 

mincutoff='off';     : Does not cut off spectrum below fmin

 

mincutoff='on';     : Cuts off spectrum below fmin

 

 

 

fmin

 

It defines minimum frequency to cut off the lower part of spectrum in (Hz)

 

 

 

fmin=0.04;     :Minimum frequency to cut off the lower part of spectrum in (Hz)

 

 

 

Notes:

 

If mincutoff='on', then results with frequency f<fmin will be removed from analysis

 

It is a simple high pass filter

 

 

 

maxcutoff

 

It defines if to cut off the spectrum beyond fmax, i.e. where f>fmax, or not

 

 

 

maxcutoff='off';     : Does not cut off spectrum beyond fmax

 

maxcutoff='on';     : Cut off spectrum beyond fmax

 

 

 

fmax

 

It defines maximum frequency to cut off the upper part of spectrum in (Hz)

 

 

 

fmax=1;     :Maximum frequency to cut off the upper part of spectrum in (Hz)

 

 

 

Notes:

 

If maxcutoff='on', then results with frequency f>fmax will be removed from analysis

 

It is a simple low pass filter

 

 

 

tailcorrection

 

It defines if to apply diagnostic tail correction or not

 

 

 

tailcorrection='off';     :Does not apply diagnostic tail

 

tailcorrection='jonswap';     :Applies JONSWAP Spectrum tail

 

tailcorrection='tma';     :Applies TMA Spectrum tail

 

 

 

Notes

 

For tailcorrection='tma', input data should be water depth

 

 

 

ftailcorrection

 

It defines frequency that diagnostic tail applies after that in (Hz)

 

 

 

ftailcorrection=0.9;     :Frequency that diagnostic tail applies after that in (Hz)

 

 

 

Notes:

 

ftailcorrection is typically set at 2.5fm where fm=1/Tm01

 

 

 

tailpower

 

It defines power that a diagnostic tail will be applied based on that

 

 

 

tailpower=-5;     :Power that a diagnostic tail will be applied based on that

 

 

 

Notes:

 

Diagnostic tail will be proportional with (f^tailpower)

 

tailpower=-3 for shallow water, tailpower=-5 for deep water

 

 

 

dispout

 

It defines if to plot spectrum or not

 

 

 

dispout='off';     :Does not plot

 

dispout='on';     :Plot

 

 

 

 

 

Required Parameters for Spectral Analysis

 

 

 

All parameters mentioned in a previous section might be required for the spectral analysis (depending on which module is on or off). In other words, if AnalysisMethod='spectral'; then all mentioned parameters above might be required. If a parameter is not required, it is ignored by Oceanlyz if defined.

 

 

 

Required Parameters for Zero-Crossing Method

 

 

 

Not all parameters mentioned in previous section are required for the zero-crossing method. If AnalysisMethod= 'zerocross'; then only following parameters are required. All other parameters, if defined, are ignored by Oceanlyz.

 

 

 

InputFileName='waterdepthsample.mat';

 

InputFileFolder=pwd;

 

SaveOutput='on';

 

OutputFileFolder=pwd;

 

AnalysisMethod='spectral';

 

WaveParameterCalc='on';

 

burst=1;

 

duration=1024;

 

fs=1;

 

heightfrombed=0.0;

 

pressureattenuation='off';

 

dispout='off';

 

 

 

 

 

 

 

Default Values for Input, Output, and Calculation Parameters

 

 

 

Default values are set as follow:

 

 

 

InputFileName='waterdepthsample.mat';

 

InputFileFolder=pwd;

 

SaveOutput='on';

 

OutputFileFolder=pwd;

 

AnalysisMethod='spectral';

 

WaveParameterCalc='on';

 

burst=1;

 

duration=1024;

 

fs=1;

 

heightfrombed=0.0;

 

nfft=2^10;

 

seaswellCalc='off';

 

fminswell=0.1;

 

fmaxswell=0.25;

 

pressureattenuation='off';

 

autofmaxpcorr='off';

 

fminpcorr=0.15;

 

fmaxpcorr=0.55;

 

mincutoff='off';

 

fmin=0.05;

 

maxcutoff='off';

 

fmax=fs/2;

 

tailcorrection='off';

 

ftailcorrection=0.9;

 

tailpower=-5;

 

dispout='off';

 

 

 

 

 

Run Oceanlyz

 

 

 

To run Oceanlyz follow these steps:

 

 

 

 

 

·       Open MATLAB or GNU Octave

 

·       Change a current folder (current directory) to a folder that contains OCEANLYZ toolbox inside MATLAB or GNU Octave. For example, open “C:\oceanlyz” directory.

 

·       Run a file named “RunOceanlyz.m” in MATLAB or GNU Octave to start calculations.

 

 

 

 

 

Outputs

 

 

 

Output(s) of the wave properties are reported based on the selected parameters as a structure array named "wave". Field(s) in the structure array "wave" can be called by using ".". For example, an output a peak period is "wave.Tp", an output for zero-moment wave height is "wave.Hm0", and an output for a water surface elevation power spectral density is "wave.Syy".

 

 

 

In general, output(s) for each time step (each burst) is reported in one row. For example, if an input file contains 5 bursts of data, then outputs has 5 rows, each row contains output for one burst. For this example, wave.Tp(1,1) or wave.Syy(1,:) are outputs for the first burst. Similarly, wave.Tp(5,1) or wave.Syy(5,:) are outputs for the fifth burst.

 

 

 

If SaveOutput='on', then the output(s) is saved in a file named "wave.mat" as a structure array in a folder defined by OutputFileFolder.

 

 

 

 

 

Notes

 

 

 

Note1:

 

If data are collected in continuous mode you can choose burst and duration as follow:

 

 

 

Duration is equal to a period of time that you want data averaged over that. For example, if you need wave properties reported every 15 min, then the duration would be 15*60 second

 

 

 

Burst is equal to the total length of the time series divided by the duration. Burst should be a rounded number. So, if the total length of the time series divided by the duration leads to a decimal number, then data should be shortened to avoid that.

 

 

 

Note2:

 

In a calculation, NFFT value that is set in “oceanlyzinput.m” file will be used. However, a user can set NFFT to be calculated automatically. This should be done inside each function. In that case, NFFT will be set equal to the smallest power of two that is larger than or equal to the absolute value of the total number of data points in each burst. This should be done manually inside each function.

 

 

 

Note3:

 

Welch spectrum is used to calculate a power spectral density. In all spectral calculation, a default window function with a default overlap window between segments is used. If any other values are required, it should be changed manually inside each function.

 

 

 

Note4:

 

If autofmaxpcorr='on', then the package calculates fmaxpcorr and ftailcorrection based on water depth and a sensor height from seabed (refer to Applying Pressure Response Factor section). A maximum value for calculated fmaxpcorr and ftailcorrection will be limited to the ones user set in “oceanlyzinput.m” file.  

 

 


 

 

Preparing Input Data File

 

 

 

This section explains about a format of input data file that is being used in a calculation. An input data file contains a time series of water surface elevation or water depth. An input data file should be prepared as a single file, with only one column. An input data file should not contain any text. Most of file formats such as ".mat", ".txt", ".xlsx", and ".csv" can be imported by the program. The schematics of input data file for data recorded in a burst mode and a continuous mode are demonstrated in following sections.

 

 

 

 

 

Data Recorded in Burst Mode

 

 

 

If data are recorded in a burst mode, then an instrument records data for a specific duration, such as 20 minutes, and then it goes to sleep to save the battery for a predefined duration, such as 60 minutes, before it wakes up again and repeats the recording/sleeping procedure again. Assuming there are M bursts of data, and each burst contains N data points, then a schematic of input data file is:

 

 

 

 

 

Burst 1, data point 1

Burst 1, data point 2

Burst 1, data point 3

.

.

.

Burst 1, data point N-2

Burst 1, data point N-1

Burst 1, data point N

Burst 2, data point 1

Burst 2, data point 2

Burst 2, data point 3

.

.

.

Burst 2, data point N-2

Burst 2, data point N-1

Burst 2, data point N

.

.

.

 

 

 

Burst M, data point 1

Burst M, data point 2

Burst M, data point 3

.

.

.

Burst M, data point N-2

Burst M, data point N-1

Burst M, data point N

 

 

 

 

 

If data recorded with a sampling frequency of fs, then an input data file has a size of:

 

 

 

(fs x duration of each burst in seconds x total number of recorded bursts,1)

 

 

 

 

 

In this mode, the total number of bursts is M. The total number of samples in each burst is:

 

 

 

N= fs x duration of each burst in seconds

 

 

 

 

 

Data Recorded in Continuous Mode:

 

 

 

If data are recorded in a burst mode, and the total number of recorded data points is equal to N, then a schematic of input data file is:

 

 

 

 

 

data point 1

data point 2

data point 3

.

.

.

 

 

 

data point N-2

data point N-1

data point N

 

 

 

 

 

If data recorded with a sampling frequency of fs for a total duration equal to hr hours, then an input data file has a size of:

 

 

 

(fs x hr x 3600,1)

 

 

 

 

 

To define a number of bursts when data are collected in a continuous mode, a length (duration) of each burst should be defined. The total number of bursts is equal to the total duration that data are collected divided by a duration one burst. For example, consider data that were recorded continuously for 24 hours. Now, if 30 minutes of data is selected to represent one burst, then the total number of bursts is equal to the 24 hours divided by 0.5 hour (30 minutes), i.e. (24 / 0.5)=48. In this example, the total number of samples in each burst is (fs x 30 x 60).

 

 

 

 

 

 


 

 

Sample Files

 

 

 

With this package, two sample files are provided. The first file is named “input.mat”, and contains one burst of data recorded by a pressure sensor. It can be used as an input for any of function. The second file is named “waterdepthsample.mat”, and contains five bursts of data recorded by a pressure sensor. It can be used as an input file while using “RunOceanlyz.m” file.

 

 

 

Measurement properties for a file named “input.mat” is:

 

 

 

Properties

Value

File name

input.mat

Number of recorded burst (burst)

1

Sampling frequency (fs)

10 (Hz)

Recording duration (duration). This is duration of one burst.

1024 (second)

Pressure sensor height from bed (heightfrombed)

0.05 (m)

Mean water depth (h)

1.07 (m)

 

 

 

 

 

Measurement properties for a file named “waterdepthsample.mat” is:

 

 

 

Properties

Value

File name

waterdepthsample.mat

Number of recorded burst (burst)

5

Sampling frequency (fs)

10 (Hz)

Recording duration (duration). This is duration of one burst.

1024 (second)

Pressure sensor height from bed (heightfrombed)

0.05 (m)

 

 


 

 

Applying Pressure Response Factor

 

 

 

Dynamic pressure from water wave will attenuate by moving from a water surface toward a seabed. Because of that, data collected by a pressure sensor have lower magnitudes compared to actual pressure values. To fix this issue and account for the pressure loss in the water column, pressure data collected by a pressure sensor should be divided by a pressure response factor. To do that, first, pressure data should be converted to water depth as:

 

 

 

 

 

 

 

Then, data are corrected as:

 

 

 

 

 

 

 

 

 

 

 

where,  is water depth read by a pressure sensor,  is pressure,  is a pressure response factor,  is a wavenumber, is a mean water depth, and  is a distance of the sensor location from a seabed.

 

 

 

One method to use these equations to correct the data is to use them on every single wave. To do that, a zero-crossing method is used to calculate  for every single wave. Then,  can be applied to each wave height or on water level time series corresponded to each wave.

 

 

 

Another method to use these equations to correct the data is to use Fast Fourier transform and convert a time series to the frequency domain, and then applying  on each frequency. In that case, calculated  in the frequency domain decreases from 1 to 0 as frequency, , increases (Figure 1). It means that, a value of  will be 1 for , and increases toward infinity as a frequency increases. Because of that, if  applies to an entire frequency domain, it results in will have very large values in the high frequency region. To avoid that, lower and upper limits for frequencies should be considered, where  is not applied to raw data outside of that range.

 

 

 

 

 

Figure 1: Schematic trend of  versus

 

 

 

A lower-frequency limit for applying  is the smallest frequency that is available in a time series such as  to . But an upper-frequency limit is depending on the deployment situations. As it was pointed out, wave pressure attenuates from a water surface toward a seabed (Figure 2). As a wave gets smaller (its frequency increases) it damps in water depth sooner. Because of that, there is a limit that a wave smaller than that (wave with larger frequency than that limit) cannot be sensed by a sensor, as their effect is damped completely before reaching a sensor depth. In that case,  should only apply to a range of frequencies that their effects reach a sensor location. So, for any case, depending on deployment situations, a range of frequencies that are large enough to reach a sensor should be calculated, and  only applies to that range.

 

 

 

 

 

Figure2: Schematic of sensor deployment setup

 

 

 

To calculate a maximum frequency that  should be applied, deep-water condition can be considered as  or . Considering a sensor setup as described in Figure 2, it can be written:

 

 

 

  or    so   

 

 

 

 

 

 

 

Then, a maximum frequency that  should apply to data up to that frequency is:

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

For the case that a pressure sensor sits on a seabed (i.e.   and ), the  and  are:

 

 

 

 

 

 

 

 

 

 

 

Figures 3 and 4 show the maximum frequency that  should not apply beyond that frequency.

 

 

 

For more details on this topic refer to Karimpour and Chen (2017) and Karimpour (2018).

 

 

 

 

 

 

 

Figure 3: Maximum frequency that  should not applied beyond that frequency.

 

 

 

 

 

Figure 4: Maximum frequency that  should not apply beyond that frequency, for a sensor sitting on a seabed.

 

 


 

 

Applying Diagnostic Tail

 

 

 

A diagnostic tail is used to correct a high frequency section of the wave spectrum. It is not recommended to use a diagnostic tail for measured data, unless data are recorded with a low sampling frequency which leads to no data in higher frequency, or in case that the noise in higher frequency contaminated the data. In these cases, a higher section of the spectrum can be replaced by a diagnostic tail as (Siadatmousavi et al. 2012):

 

 

 

 

 

 

 

where,  is a water surface elevation power spectral density,  is a frequency,  is a frequency that a tail applied after that, and  is a power coefficient. The  typically set at , where  is a mean frequency (Ardhuin et al. 2010). A value of  depends on deployment conditions, however, typically it is -5 for deep and -3 for shallow water (e.g. Kaihatu et al. 2007, Siadatmousavi et al. 2012).

 

 

 

For more details on this topic refer to Karimpour and Chen (2017) and Karimpour (2018).

 

 

 

 


 

 

Functions

 

 

 

This code contains 5 main functions to analyze wave data (water level data) time series. Water level could be in a form of the total depth or it can be in the form of surface elevation. If a pressure sensor is used for data collection, then it is required to correct the data for pressure attenuation in the water column. In such cases, only water depth data should be used as an input. In any other cases, either of water surface elevation or water depth time series can be used for data analysis.

 

 

 

The provided “RunOceanlyz.m” file analyzes the the entire time series by calling an appropriate function(s). Note that any of these functions may be called by a user in his/her own code as well.

 

 

 

The main functions that are used in this toolbox are:

 

 

 

Function

Description

WaveSpectraFun

Calculates wave properties from spectral analysis

WaveZerocrossingFun

Calculates wave properties from zero-crossing

PcorFFTFun

Corrects pressure attenuation effect using Fast Fourier transform (FFT)

PcorZerocrossingFun

Corrects pressure attenuation effect using zero-crossing

SeaSwellFun

Separates wind sea and swell

 

 

 

 

 

Calculating Wave Parameters Using Spectral Analysis Method (WaveSpectraFun)

 

 

 

Function “WaveSpectraFun.m” calculates wave properties using wave surface elevation power spectral density.

 

 

 

Inputs for this function are:

 

 

 

Value

Description

input

wave data in (m). Input file should be in column format (n*1)

fs

sampling frequency that data collected at in (Hz)

duration

duration time that data collected in input in each burst in second

nfft

NFFT for Fast Fourier Transform

heightfrombed

sensor height from bed

fmin

minimum frequency for cut off the lower part of spectra

fmax

maximum frequency for cut off the upper part of spectra

ftailcorrection

frequency that diagnostic tail apply after that (typically set at 2.5fm, fm=1/Tm01)

tailpower

power that diagnostic tail apply based on that (-3 for shallow to -5 for deep water)

mincutoff

defines if to cut off the spectra below fmin or not (0: cutoff off, 1: cutoff on)

maxcutoff

defines if to cut off the spectra beyond fmax or not (0: cutoff off, 1: cutoff on)

tailcorrection

defines if to apply diagnostic tail correction (0: not apply, 1: JONSWAP tail, 2: TMA tail)

dispout

defines to display outputs or not ('off': not display, 'on': display)

 

 

 

 

 

Outputs of this function are:

 

 

 

Value

Description

Hm0

Zero-Moment Wave Height (m)

Tm01

Wave Period from m01 (second), Mean Wave Period

Tm02

Wave Period from m02 (second), Mean Zero Crossing Period

Tp

Peak Wave Period (second)

fp

Peak Wave Frequency (Hz)

f

Frequency (Hz)

Syy

Wave Surface Elevation Power Spectral Density (m^2s)

 

 

 

 

 

In case that data are measured by a pressure sensor, data should be corrected for pressure attenuation. In that case, function “PcorFFTFun.m” or “PcorZerocrossingFun.m” should be called first to correct pressure (water depth) data before calculating wave properties. Results from either of those functions can be imported to this function for spectral analysis of wave parameters. If provided “RunOceanlyz.m” file is used, it will do this procedure if proper input parameters are selected in “oceanlyzinput.m” file. Please read the note in “PcorFFTFun.m”.

 

 

 

This function can be used as a standalone command in Matlab/GNU Octave command line or it can be embedded in Matalb/GNU Octave script file (.m file) as:

 

 

 

[Hm0,Tm01,Tm02,Tp,fp,f,Syy]=WaveSpectraFun(input,fs,duration,nfft,h,heightfrombed,fmin,fmax,ftailcorrection,tailpower,mincutoff,maxcutoff,tailcorrection,dispout);

 

 

 

Example for using a provided sample input file:

 

 

 

[Hm0,Tm01,Tm02,Tp,fp,f,Syy]= WaveSpectraFun(input,10,1024,2^10, 0.05,0.04,2,0.9,-4,’on’, ’on’, ’off’, ’on’);

 

 

 

 

Calculating Wave Parameters Using Zero-Crossing Method (WaveZerocrossingFun)

 

 

 

Function “WaveZerocrossingFun.m” calculates wave properties using an upward zero-crossing method.

 

 

 

Inputs for this function are: 

 

 

 

Value

Description

input

wave data in (m). Input file should be in column format (n*1).

fs

sampling frequency that data collected at in (Hz)

duration

duration time that data collected in input in each burst in second

dispout

defines to display outputs or not ('off': not display, 'on': display)

 

 

 

 

 

Outputs of this function are: 

 

 

 

Value

Description

Hs

Significant Wave Height (m)

Hz

Zero Crossing Mean Wave Height (m)

Tz

Zero Crossing Mean Wave Period (second)

Ts

Significant Wave Period (second)

H

Wave Height Data Series (m)

T

Wave Period Data Series (second)

 

 

 

 

 

In case that data are measured by a pressure sensor, data should be corrected for pressure attenuation. In that case, function “PcorFFTFun.m” or “PcorZerocrossingFun.m” should be called first to correct pressure (water depth) data before calculating wave properties. Results from either of those functions can be imported to this function for zero-crossing calculation of wave parameters. If provided “RunOceanlyz.m” file is used, it will do this procedure if proper input parameters are selected in “oceanlyzinput.m” file. Please read the note in “PcorFFTFun.m”.

 

 

 

This function can be used as a standalone command in Matlab/GNU Octave command line or it can be embedded in Matalb/GNU Octave script file (.m file) as:

 

 

 

 [Hs,Hz,Tz,Ts,H,T]=WaveZerocrossingFun(input,fs,duration,dispout);

 

 

 

Example for using a provided sample input file:

 

 

 

[Hs,Hz,Tz,Ts,H,T]=WaveZerocrossingFun(input,10,1024,’on’);

 

 

 

 

 

Correcting Wave Pressure Data Using Fast Fourier Transform (PcorFFTFun)

 

 

 

Function “PcorFFTFun.m” corrects water level data measured by a pressure sensors using a pressure response factor, or a pressure attenuation coefficient, and the Fast Fourier Transform.

 

 

 

Inputs for this function are:

 

 

 

Value

Description

input

wave data in (m). Input file should be in column format (n*1).

fs

sampling frequency that data collected at in (Hz)

duration

duration time that data collected in input in each burst in second

nfft

NFFT for Fast Fourier Transform

h

mean water depth in (m)

heightfrombed

sensor height from bed

fminpcorr

minimum frequency for applying pressure attenuation factor

fmaxpcorr

maximum frequency for applying pressure attenuation factor

tailcorrection

define if to apply diagnostic tail correction (0: not apply, 1: JONSWAP tail, 2: TMA tail)

pressureattenuation

define if to apply pressure attenuation factor or not (0:off, 1:on, without correction after fmaxpcorr, 2:on, with constant correction after fmaxpcorr)

autofmaxpcorr

define if to calculate fmaxpcorr and ftailcorrection based on water depth or not (0:off, 1:on)

dispout

defines to display outputs or not ('off': not display, 'on': display)

 

 

 

 

 

Output of this function is: 

 

 

 

Value

Description

Eta

Corrected Water Surface Level Time Series (m)

 

 

 

 

 

This function can be used as a standalone command in Matlab/GNU Octave command line or it can be embedded in Matalb/GNU Octave script file (.m file) as:

 

 

 

[Eta,ftailcorrection]=PcorFFTFun(input,fs,duration,nfft,h,heightfrombed,fminpcorr,fmaxpcorr,ftailcorrection,pressureattenuation,autofmaxpcorr,dispout);

 

 

 

Example for using a provided sample input file:

 

 

 

[Eta]=PcorFFTFun(input,10,1024,2^10,1.07,0.05,0.04,0.8,’all’,’off’,’on’);

 

 

 

 

 

Note:

 

In case of using a pressure attenuation correction function, choosing a proper upper limit for applying a correction is essential. If the upper limit is chosen unreasonably high, it will lead to a wrong over-estimation of the results. The under-estimation also can happen, if the upper limit is chosen unreasonably low. The waves that are large enough that can be sensed by a sensor should be included in attenuation correction. Note that pressure from small waves with high frequency can be damped in the water column and therefore may not reach a sensor depth. In that case, these small waves should not be included in attenuation correction. Wave height, wave frequency, water depth, height of a sensor above a seabed all play roles on if a wave effect reaches down to a sensor depth or not. Because of that, a deployment situation and wave properties should be used to define the highest frequency that a sensor senses its effect. Pressure correction should not be applied beyond that point. (refer to Applying Pressure Response Factor section)

 

 

 

 

Correcting Wave Pressure Data Using Linear Wave Theory (PcorZerocrossingFun)

 

 

 

Function “PcorZerocrossingFun.m” corrects water level data measured by a pressure sensors using a pressure response factor, or a pressure attenuation coefficient, and a zero-crossing method.

 

 

 

Inputs for this function are: 

 

 

 

Value

Description

input

wave data in (m). Input file should be in column format (n*1).

fs

sampling frequency that data collected at in (Hz)

duration

duration time that data collected in input in each burst in second

h

mean water level (m)

heightfrombed

sensor height from bed

dispout

defines to display outputs or not ('off': not display, 'on': display)

 

 

 

 

 

Output of this function is: 

 

 

 

Value

Description

Eta

Corrected Water Surface Level Time Series (m)

 

 

 

This function can be used as a standalone command in Matlab/GNU Octave command line or it can be embedded in Matalb/GNU Octave script file (.m file) as:

 

 

 

 [Eta]=PcorZerocrossingFun(input,fs,duration,h,heightfrombed,dispout);

 

 

 

Example for using a provided sample input file:

 

 

 

[Eta]=PcorZerocrossingFun(input,10,1024,1.07,0.05,’on’);

 

 

 

 

Separating Sea and Swell (SeaSwellFun)

 

 

 

Function “SeaSwellFun.m” separates sea wave from swell wave using a wave surface elevation power spectral density.

 

 

 

Inputs for this function are:  

 

 

 

Value

Description

input

wave data in (m). Input file should be in column format (n*1).

fs

sampling frequency that data collected at in (Hz)

duration

duration time that data collected in input in each burst in second

nfft

NFFT for Fast Fourier Transform

h

mean water depth in (m)

fmin

minimum frequency for cut off the lower part of spectra

fmax

maximum frequency for cut off the upper part of spectra

ftailcorrection

frequency that diagnostic tail apply after that (typically set at 2.5fm, fm=1/Tm01)

tailpower

power that diagnostic tail apply based on that (-3 for shallow to -5 for deep water)

fminswell

minimum frequency that is used for Tpswell calculation

fmaxswell

maximum frequency that swell can have, It is about 0.2 in Gulf of Mexico

mincutoff

define if to cut off the spectra below fmin (0: cutoff off, 1: cutoff on)

maxcutoff

define if to cut off the spectra beyond fmax (0: cutoff off, 1: cutoff on)

tailcorrection

define if to apply diagnostic tail correction (0: not apply, 1: JONSWAP tail, 2: TMA tail)

dispout

defines to display outputs or not ('off': not display, 'on': display)

 

 

 

 

 

Outputs of this function are:  

 

 

 

Value

Description

Hm0

Zero-Moment Wave Height (m)

Hm0sea

Sea Zero-Moment Wave Height (m)

Hm0swell

Swell Zero-Moment Wave Height (m)

Tp

Peak Wave Period (second)

Tpsea

Peak Sea Period (second)

Tpswell

Peak Swell Period (second)

fp

Peak Wave Frequency (Hz)

f

Frequency (Hz)

fseparation

Sea and Swell Separation Frequency (Hz)

Syy

Wave Surface Elevation Power Spectral Density (m^2s)

 

 

 

 

 

In case that data are measured by a pressure sensor, data should be corrected for pressure attenuation. In that case, function “PcorFFTFun.m” or “PcorZerocrossingFun.m” should be called first to correct pressure (water depth) data before calculating wave properties. Results from either of those functions can be imported to this function for spectral analysis of wave parameters. If provided “RunOceanlyz.m” file is used, it will do this procedure if proper input parameters are selected in “oceanlyzinput.m” file. Please read the note in “PcorFFTFun.m”.

 

 

 

This function can be used as a standalone command in Matlab/GNU Octave command line or it can be embedded in Matalb/GNU Octave script file (.m file) as:

 

 

 

[Hm0,Hm0sea,Hm0swell,Tp,Tpsea,Tpswell,fp,fseparation,f,Syy]=SeaSwellFun(input,fs,duration,nfft,h,fmin,fmax,ftailcorrection,tailpower,fminswell,fmaxswell,mincutoff,maxcutoff,tailcorrection,dispout);

 

 

 

Example for using a provided sample input file:

 

 

 

[Hm0,Hm0sea,Hm0swell,Tp,Tpswell,fp,fseparation,f,Syy]= SeaSwellFun(input,10,1024,2^10,1.07,0.05,2,0.9,-4,0.1,0.25,’on’, ’on’, ’off’, ’on’);

 

 

 

 

 

 

References

 

 

 

Ardhuin, F., Rogers, E., Babanin, A. V., Filipot, J. F., Magne, R., Roland, A., ... & Collard, F. (2010). Semiempirical dissipation source functions for ocean waves. Part I: Definition, calibration, and validation. Journal of Physical Oceanography, 40(9), 1917-1941.

 

 

 

Kaihatu, J. M., Veeramony, J., Edwards, K. L. & Kirby, J. T. 2007 Asymptotic behaviour of frequency and wave number spectra of nearshore shoaling and breaking waves. J. Geophys. Res. 112, C06016

 

 

 

Karimpour, A., & Chen, Q. (2017). Wind Wave Analysis in Depth Limited Water Using OCEANLYZ, a MATLAB toolbox. Computers & Geosciences, 106,181-189.

 

 

 

Karimpour A., (2018), Ocean Wave Data Analysis: Introduction to Time Series Analysis, Signal Processing, and Wave Prediction, KDP.

 

 

 

Siadatmousavi, S. M., Jose, F., & Stone, G. W. (2011). On the importance of high frequency tail in third generation wave models. Coastal Engineering.