Home > wafo > spec > wallop.m

wallop

PURPOSE ^

Calculates (and plots) a Wallop spectral density.

SYNOPSIS ^

S1=wallop(w1,sdata,plotflag)

DESCRIPTION ^

 WALLOP  Calculates (and plots) a Wallop spectral density.
  
  CALL:  S = wallop(w,data,plotflag); 
         S = wallop(wc,data,plotflag);
 
         S    = a struct containing the spectral density, see datastructures.
         w    = angular frequency (default linspace(0,3,257))
         wc   = angular cutoff frequency (default 33/Tp)
         data = [Hm0 Tp M]
                Hm0 = significant wave height (default 7 (m))
                Tp  = peak period (default 11 (sec))
                M   = shape factor, i.e. slope for the high frequency
                      part (default depending on Hm0 and Tp, see below)
     plotflag = 0, do not plot the spectrum (default).
                1, plot the spectrum.
 
   The WALLOP spectrum parameterization used is 
 
      S(w)=Bw*Hm0^2/wp*(wp./w).^M.*exp(-M/4*(wp./w).^4);
  where
      Bw = normalization factor
      M  = abs((log(2*pi^2)+2*log(Hm0/4)-2*log(Lp))/log(2));
      Lp = wave length corresponding to the peak frequency, wp.
 
   If M=5 it becomes the same as the JONSWAP spectrum with 
   peak enhancement factor gamma=1 or the Pierson-Moskowitz spectrum. 
 
  Example: 
    S = wallop(1.1,[6.5 10]), wspecplot(S)
   
  See also  jonswap, torsethaugen, simpson

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

001 function S1=wallop(w1,sdata,plotflag)
002 %WALLOP  Calculates (and plots) a Wallop spectral density.
003 % 
004 % CALL:  S = wallop(w,data,plotflag); 
005 %        S = wallop(wc,data,plotflag);
006 %
007 %        S    = a struct containing the spectral density, see datastructures.
008 %        w    = angular frequency (default linspace(0,3,257))
009 %        wc   = angular cutoff frequency (default 33/Tp)
010 %        data = [Hm0 Tp M]
011 %               Hm0 = significant wave height (default 7 (m))
012 %               Tp  = peak period (default 11 (sec))
013 %               M   = shape factor, i.e. slope for the high frequency
014 %                     part (default depending on Hm0 and Tp, see below)
015 %    plotflag = 0, do not plot the spectrum (default).
016 %               1, plot the spectrum.
017 %
018 %  The WALLOP spectrum parameterization used is 
019 %
020 %     S(w)=Bw*Hm0^2/wp*(wp./w).^M.*exp(-M/4*(wp./w).^4);
021 % where
022 %     Bw = normalization factor
023 %     M  = abs((log(2*pi^2)+2*log(Hm0/4)-2*log(Lp))/log(2));
024 %     Lp = wave length corresponding to the peak frequency, wp.
025 %
026 %  If M=5 it becomes the same as the JONSWAP spectrum with 
027 %  peak enhancement factor gamma=1 or the Pierson-Moskowitz spectrum. 
028 %
029 % Example: 
030 %   S = wallop(1.1,[6.5 10]), wspecplot(S)
031 %  
032 % See also  jonswap, torsethaugen, simpson
033 
034 % References:
035 % Huang, N.E., Long, S.R., Tung, C.C, Yuen, Y. and Bilven, L.F. (1981)
036 % "A unified two parameter wave spectral model for a generous sea state"
037 % J. Fluid Mechanics, Vol.112, pp 203-224
038 
039 % Tested on: matlab 6.0, 5.3
040 % History:
041 % revised jr 03.04.2001
042 %  - added wc to input 
043 %  - updated information
044 % revised pab 18.02.2000
045 %  - normalization so that int S(w) dw = m0
046 % revised pab 24.01.2000
047 %  - updated note to 'Wallop Hm0='....
048 % by pab 01.12.99
049 
050 monitor=0;
051 
052 if nargin<3|isempty(plotflag)
053   plotflag=0;
054 end
055 
056 % Old call  
057 %if nargin<1|isempty(w1)
058 %  w=linspace(0,3,257).';
059 %end
060 
061 M=[];
062 if nargin<2|isempty(sdata)
063   sdata=[7 11];
064 else
065   switch length(sdata)
066     case 1, sdata=[sdata 11];
067     case 3, M=sdata(3); sdata=sdata(1:2);
068   end
069 end %
070 
071 if nargin<1|isempty(w1), wc = 33/sdata(2)
072 elseif length(w1)==1,    wc = w1; 
073 else w = w1 ; end
074 nw = 257;
075 if isempty(w), w = linspace(0,wc,nw).'; end
076  
077 n=length(w);
078 
079 S1=createspec;
080 S1.S=zeros(n,1);
081 S1.w=w;
082 S1.norm=0; % The spectrum is not normalized
083 
084 Hm0=sdata(1);
085 Tp=sdata(2);
086 S1.note=['Wallop, Hm0 = ' num2str(Hm0)  ', Tp = ' num2str(Tp)];
087 wp=2*pi/Tp;
088 
089 if monitor
090   disp(['Hm0, Tp      = ' num2str([Hm0 Tp])])
091 end
092 
093 if isempty(M),
094   kp=w2k(wp,0,inf); % wavenumber at peak frequency
095   Lp=2*pi/kp; % wave length at the peak frequency
096   M=abs((log(2*pi^2)+2*log(Hm0/4)-2*log(Lp))/log(2));
097 end
098 
099 %Bw=0.06238*M^((M-1)/4)/(4^((M-5)/4)*gamma((M-1)/4))*(1+0.7458*(M+2)^(-1.057));
100 
101 Bw = M^((M-1)/4)/(4^((M-5)/4)*gamma((M-1)/4))/16;
102 
103 % for w>0 % avoid division by zero
104 k=find(w>0);
105 S1.S(k)=Bw*Hm0^2/wp*(wp./w(k)).^M.*exp(-M/4*(wp./w(k)).^4);
106 
107 if plotflag
108   wspecplot(S1,plotflag)
109 end
110 
111

Mathematical Statistics
Centre for Mathematical Sciences
Lund University with Lund Institute of Technology

Comments or corrections to the WAFO group


Generated on Thu 06-Oct-2005 02:21:16 for WAFO by m2html © 2003