Matlab script Decomp

18 April 2005

 

Separate spectrum of irregular wave in incident and reflected part.

 

1.       Calculate distance of wave gauges needed

 

To be able to separate the waves into incident and reflected wave, 2 or 3 wave gauges are needed (3 is preferable). The signals from the wave gauges lead to a system of equations, which will enable the calculation of incident and reflected wave. If the distance is a multiple of the wavelength, the system of equations will be singular.

If three wave gauges are used, distances of 0.3 m between gauge 1 and gauge 2, and 0.4 m between gauge 2 and gauge 3 will do in most cases, where gauges 1 is supposed to be nearest to the wave board. This will work if the peek period is between 1 to 3 seconds.

If two wave gauges are used, the program Distance of WL | Delft Hydraulics can be used:

 

DelftAUKE, program Distance (DOS), 2 wave gauges.

 

Input to Distance:

 

Tp   peak period

h    water depth

 

Output:

 

fmin - fmax                  frequency range of valid solution

LminLmax                wave lengths according to frequencies at given water depth

dxlow - dxhigh                  range of distances between which a valid solution is possible.

 

The minimum distance of the wave gauges specified in the manual is 0.2 m. If this distance is not in the range found by the program Distance, specify 0.2 m as input to the program Distance. In that case Distance calculates the frequncy range of valid solution if the distace is 0.2 m.

 

 

2.       Measurements -> Dasylab file

 

3.       Input parameters to Decomp (input from .txt-file):

 

About measurements:

 

parameter

description

dt (s)

time step used in data file

Tp (s)

approximated peak period (only used to calculate the frequency resolution in the variance spectra to be calculated)

fres

desired frequency resolution of spectrum/peak frequency

fmin - fmax (Hz)

frequency range where the spectrum of incident and reflected wave must be calculated (fmax=-1 to calculate this interval automatically from the distances in case of two gauges, and takes the default 0 – 2.5 Hz in case of 3 gauges).

Some times, fmax=-1 will select a part that doesn’t cover the whole Jonswap spectrum. Decomp will show a message box on the screen in that case.

thresh

threshold to determine the part of the variance spectrum where the decomposition is valid, E>=threshEmax.

Most of the times, the frequency range of an ideal Jonswap spectrum with the same threshold will extend this part, caused by noise in the calculated spectrum. The program will print both ranges in the output file.

 

 

col_t

column number in the data file where the time has been registrated(0 if not registrated);

col_t is used by the input function only; it will not be used by this program, but it is needed to locate the colum numbers of the wave gauges correctly.

 

     

about channels to be analysed:

 

column number

position of gauge (m)

scale from registered data value to meters

 

about the water depth as a function of the position:

 

x (m), h (m)            position and water depth in at least P+1 points, where P is the number of wave gauges

 

 

 

4.       Program Decomp

 

Method:                  based on method of Zelt en Skjelbreia (1992)

Matlab version:        Henk Jan Bakkenes (2002)

 

Suppose N points, measured at time t=n∆t, n=1,2,…N. In this case a wave can be described by N spectral terms.

The data registered by wave gauge p can be represented by a sum of incident and reflected wave:

 

                                             (1)

 

where ai,j and ar,j are complex (i.e. containing the phase shift of the incident and reflected wave respectively). The mean value of η is zero, so j=0 disappears.

Furthermore, ωj=2π/T, whereT=N∆t.

 

Standard Fourier analysis of η(xp,t ) gives:

 

                                                                   (2)

 

Combination of (1) and (2) -> system of P linear equations and two unknown complex values (ai,j and ar,j ) for each spectral component j  (P : number of gauges).

A least-squares solution can be calculated using a standard Matlab operation.

Equations must be well conditioned: it is not allowed to put wave gauges in points of equal phase.

To calculate the distances needed, a DOS program Distance is available (from the Delft Hydraulics package DELFT-AUKE).

 

Steps of Decomp:

 

(1)    Ask which output the user desires (see list on p.5).

 

(2)    Read data file, parameter file, and file containing profile of water depth.

 

(3)    If asked by user: calculate frequency range of valid decomposition according to the distances of the wave gauges. Otherwise, use the values supplied by the user.

 

(4)    Ask part of input data to be analysed (begin and end point).

 

(5)    Calculate phase shift at positions of the wave gauges, and if needed, calculate the change of the amplitude caused by changes in water depth.

 

(6)    For each wave gauge, calculate FFT A , variance spectrum E, and averaged spectrum Es, averaged over frequency bands of length specified by input.

 

Points (7) through (9) describe parameter calculations per wave gauge. The parameters are chosen according to a program Reflecmf from the DELFT-AUKE package (manual July 2000).

 

(7)    From averaged variance spectrum: peak period Tp , period TpD in middle of connected part of the spectrum where Es >= 0.8 Esmax , and the frequency range of a connected part of the spectrum where Es >= thresh Emax . This frequency range will be used to calculate the parameters in points (8) and (9).

 

(8)    From "raw" variance spectrum in frequency range Es >= threshEmax : spectral moments (m0, m2, and m4).

 

(9)    From spectral moments: mean period between successive zero crossings and between successive maxima, and correlation between successive wave heights.

 

Continuation of program Decomp:

 

As a valid frequency range, the frequencies must be valid for the used distances (point (3)), and the spectrum must be above the threshold (point (7)), otherwise the equations will be of the type 0.x = 0.

Most of the times, an ideal Jonswap spectrum will have a wider frequency range above the threshold. This range will be given in the output as well. Compare the range used by Decomp with the frequecy range of an ideal Jonswap spectrum to evaluate the results.

In the figures, more frequencies than the above range are shown. This will help the user to inspect if the valid frequency range contains a good representation of the signal. Vertical lines with arrows are used to show which part of the spectrum is valid.

If the coefficient matrix of the equations is nearly singular, the program will stop, and show an error message. Sometimes the conditions in the equations are reasonable, but the spectrum will show an unexpected peak. If this peak is in the valid part of the spectrum, the distances may be chosen incorrectly. Check if the values are in accordance with this program.

More information about distances: DELFT-AUKE Process Document (2000).

 

 

(10)   Calculate amplitude spectra of incident and reflected wave from FFT-values A for each spectral component.

 

(11)   Variance spectrum and spectrum averaged over frequency bands of incident and reflected wave in specified frequency range.

 

(12)   Spectral parameters from "raw" variance spectra of incident and reflected wave: moment m0, significant wave height, and mean factor of reflection.

 

(13)   Reconstruction of the incident and reflected wave as a time series by inverse FFT. Only the "valid" part of the spectrum will be used for this computation.

 

 


Parameters calculated in steps (8) and (9):

 

parameter

description

spectral moments

T0,1 = m0/m1

average period

average period between two zero crossings

average period between two successive local maximum values

narrowness parameter

broadness parameter

κ

correlation value, defined in DELFT-AUKE manual (2000), p.4-6

γ

coefficient of linear correlation, defined in DELFT-AUKE manual (2000), p.4-7 as ρHH

Hm0 = 4√m0

significant wave height

 

                          

Usage of the program:

 

N.B.:     This program uses dialog boxes. In Matlab 7, dialog boxes will not be deleted immediately after clicking OK, but the program will continue; please don’t click OK again.

 

Input files:

·         data file (Dasylab file),

 

·         parameter file: .txt file

Contents of parameter file:

 

time step used in data file

0.02

Tp (s),  f-resol/Tp  fmin (Hz),  fmax (Hz),  thresh

1.27,    0.025,      0.25,       1.65,       0.005

column with time (0 if none)

1

columno gauge, pos (m),    scale fact to m

2,             0,         0.025

3,             0.25,       0.025

 

Line number 1, 3, 5, and 7 are text lines; the other lines are data lines; if more values are stored in one line, the values must be separated by commas. The values are examples, and have to be replaced by the user. If fmax = -1, fmin and fmax will be calculated by the program (2 gauges), or replaced by 0.02 and 2.5 respectively (3 or more gauges).

 

·         profile of depth

Contents of file containing profile of depth:

 

pos (m),       depth (m)        

0.0            0.6

0.125          0.6

0.25           0.6

0.5                                                0.6

 

Line 1 contains text; the other lines contain data, separated by blanks or tabs. The values are examples.

If P gauges are used, at least P+1 data lines are needed (this restriction has not yet been checked).

 

N.B.:     In the file containing the depth profile, spaces or tabs ar allowed as a separation. In the parameter file, commas are allowed only, since the data are read in a different way. The program will send a message if a comma is missing.

 

Usage:

 

·         Open Matlab and go to the folder where Decomp and subfolder .\Functions are stored.

Type Decomp.

 

Remark:      the user must make the parameter file and the file containing the profile of depth in advance. The contents are described above.

 

·         Desired output:

 

The program shows a list of possible outputs:

 

plot of raw variance spectrum of original time signal

plot of averaged variance spectrum of original time signal

plot of averaged variance spectrum of incident and reflected wave

plot of reflection factor

file of statistical results (spectral moments, refl. coeff)

file of variance spectra of inc. and refl. waves

file of time signal of inc. and refl. wave

 

For each item of the list, fill in 1 for yes or 0 for no.

 

·         The input files are asked from dialog boxes.

 

·         The time interval to be used in the calculations is asked by the program. A plot of the measurements will be put on the screen, which will help the user to select suitable data.

The begin time and end time in seconds must be given in a dialog box.

Close the figure to start the calculation.

 

·         If the frequency range selected in the input file covers less than 90% of the frequency range of an ideal Jonswap spectrum with peak period TpD (mean value of all wave gauges), a message box will be shown on the screen.

The frequency range of an ideal Jonswap spectrum will be printed to the output file of statistical results. The algorithme to calculated this ideal Jonswap spectrum has been copied from Ap van Dongeren, WL | Delft Hydraulics (see function FN_jonswap).

 

·         The desired figures will be produced, and dialog boxes will ask the names of the desired output files. If the spectra of incident and reflected waves show peeks, the equations probably have a singularity at those frequencies.

 

 

The incident and reflected wave as time series are calculated using the inverse FFT of the spectra of incident and reflected wave. Only frequencies for which the spectrum has been detected as "valid" will be used for the reconstruction. If the output file of time series of incident and reflected wave is asked by the user (see "Desired Output"), these series can be used as input for the DELFT-AUKE program Waves, which determines special wave parameters. A special Matlab script Asc2wav has to be run to make the files readable by Waves (and other DELFT-AUKE programs). Information from the programmer of Stevin III, room 0.07, tel. 85974.

How to handle if a file cannot be read

 

The next errors may occur:

 

  • If the number of header lines specified in case of a non-Dasylab file is less than the actual number, the read function used in the scripts gives a message. However, in some cases Matlab gives the next error:

 

Number of colums in line -- of ASCII file ... must be the same as the previous line.

 

In such a case, remove all spaces, signs, or numbers from the first position of each header line.

 

If the number of header lines specified is more than the actual number, no message will occur, but less data lines will be read.

 

  • If a Dasylab file gives an error like above, check if the right file has been chosen. If yes, copy this file, delete the first line, and read the copy specifying the remaining number of header lines.

 

  • The maximum number of samples to be read is exceeded. This number is set to 500000/nchan,

where nchan is the number of channels used. If more data are needed, change the value of maxnum in \Functions\FN_indata7, line 86, or read the file in parts, using different sample numbers. Sample numbers to be read can be specified via the dialog box:

 

   read all samples?

           

            Specify No in this case.   

 

 

In a Dasylab file, it is preferable to store the time in decimal format. The time format with a ':' in Dasylab is not read as a time in Matlab.

 

 

 

 

References

 

Bakkenes, Henk Jan (2002), Observation and separation of bound and free low-frequency waves in the nearshore zone, MSc-Thesis, Delft University of Technology, Faculty of Civil Engineering and Geosciences, Fluidmechanics section, June 2002.

 

Battjes,J.A., Windgolven. Collegehandleiding  b78, uitgave 1984, 5e herdruk 1992 (Dutch version).

 

Goda, Y. (1985). Random seas and design of maritime structures. Univ. Tokyo Press, Tokyo, Japan.

 

WL Delft Hydraulics (2000), Delft-Auke Process Description, July 2000.

 

Zelt, J.A. en J.E. Skjelbreia (1992), Estimating Incident and Reflected Wave Fields Using an Arbitrary Number of Wave Gauges, Proc. 23rd Int. Conf. Coastal Eng., ASCE, pp. 777-789.