Home > wafo > wavemodels > mk87pdf.m

mk87pdf

PURPOSE ^

Myrhaug and Kjeldsen (1987) joint (Scf,Hd) PDF.

SYNOPSIS ^

p = mk87pdf(Hd,Scf,Hs,Tz, condon,normalizedInput)

DESCRIPTION ^

 MK87PDF Myrhaug and Kjeldsen (1987) joint (Scf,Hd) PDF. 
  
  CALL:  f = mk87pdf(Hd,Scf,Hs,Tz) 
  
     f  = density  
    Hd  = zero down crossing wave height 
    Scf = crest front steepness 
    Hs  = significant wave height. 
    Tz  = average zero down crossing period.   
  
  MK87PDF returns the joint PDF of (Scf, Hd) given Hs and Tz, 
  i.e., crest front steepness (2*pi*Ac/(g*Td*Tcf)) and wave height, given 
  the seastate. The root mean square values of Hd and Scf (Hrms,Erms) are 
  related to the significant waveheight and the average zero down 
  crossing period by: 
              Hrms = 0.715*Hs; 
              Erms = 0.0202+0.826*Hs/(Tz^2); 
  
    The size of f is the common size of  the input arguments 
  
  Example: 
   Hs = 7;Tz=10;   
   h  = linspace(0,3*Hs)';  
   s  = linspace(0,4*Hs/Tz^2)'; 
   [S,H] = meshgrid(s,h); 
   contour(s,h,mk87pdf(H,S,Hs,Tz)) 
  
  See also  mk87pdf2, wweibpdf, lognpdf

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SUBFUNCTIONS ^

SOURCE CODE ^

001 function p = mk87pdf(Hd,Scf,Hs,Tz, condon,normalizedInput) 
002 %MK87PDF Myrhaug and Kjeldsen (1987) joint (Scf,Hd) PDF. 
003 % 
004 % CALL:  f = mk87pdf(Hd,Scf,Hs,Tz) 
005 % 
006 %    f  = density  
007 %   Hd  = zero down crossing wave height 
008 %   Scf = crest front steepness 
009 %   Hs  = significant wave height. 
010 %   Tz  = average zero down crossing period.   
011 % 
012 % MK87PDF returns the joint PDF of (Scf, Hd) given Hs and Tz, 
013 % i.e., crest front steepness (2*pi*Ac/(g*Td*Tcf)) and wave height, given 
014 % the seastate. The root mean square values of Hd and Scf (Hrms,Erms) are 
015 % related to the significant waveheight and the average zero down 
016 % crossing period by: 
017 %             Hrms = 0.715*Hs; 
018 %             Erms = 0.0202+0.826*Hs/(Tz^2); 
019 % 
020 %   The size of f is the common size of  the input arguments 
021 % 
022 % Example: 
023 %  Hs = 7;Tz=10;   
024 %  h  = linspace(0,3*Hs)';  
025 %  s  = linspace(0,4*Hs/Tz^2)'; 
026 %  [S,H] = meshgrid(s,h); 
027 %  contour(s,h,mk87pdf(H,S,Hs,Tz)) 
028 % 
029 % See also  mk87pdf2, wweibpdf, lognpdf 
030  
031 %   References: 
032 %   Myrhaug, D. and Kjelsen S.P. (1987)  
033 %  'Prediction of occurences of steep and high waves in deep water'. 
034 %   Journal of waterway, Port, Coastal and Ocean Engineers,  
035 %   Vol. 113, pp 122--138 
036    
037 % 
038 %   Myrhaug & Dahle (1984)  
039 %   'Parametric modelling of joint probability density  
040 %   distributions for steepness and asymmetry in deep water waves' 
041 % 
042  
043 % tested on: matlab 5.1 
044 % history: 
045 % revised pab 09.08.2003   
046 % Changed input + updated help header 
047 % revised pab 23.01.2001 
048 % - no longer dependent on stats toolbox only wstats 
049 % revised pab 10.02.2000 
050 % - fixed normalization 
051 % by  Per A. Brodtkorb 1998 
052  
053  
054 error(nargchk(3,6,nargin)) 
055    
056 if nargin < 6|isempty(normalizedInput),  normalizedInput=0; end 
057 if nargin < 5|isempty(condon),  condon=0; end % regular pdf is returned 
058 if nargin < 4|isempty(Tz),  Tz=8;end 
059 if nargin < 3|isempty(Hs),  Hs=6;end 
060  
061 [errorcode, Scf,Hd,Tz,Hs] = comnsize(Scf,Hd,Tz,Hs); 
062 if errorcode > 0 
063     error('Requires non-scalar arguments to match in size.'); 
064 end 
065 if (normalizedInput>0) 
066   Hrms = 1; 
067   Erms = 1; 
068 else 
069   Hrms = 0.715*Hs; 
070   Erms = 0.0202 + 0.826*Hs./(Tz.^2);  
071 end 
072  
073 s = Scf./Erms; 
074 h = Hd./Hrms; 
075  
076 sig = (-0.21*atan(2*(h-1.4))+0.325); 
077 sig = max(sig,eps); % pab fix: make sure it is positive 
078 p   = zeros(size(h)); 
079  
080 % NB! weibpdf must be modified to correspond to 
081 % pdf=x^(b-1)/a^b*exp(-(x/a)^b) or else insert 
082 % weibpdf=2.39.*h.^1.39/(1.05^2.39).*exp(-(h./1.05).^2.39); 
083 switch condon, 
084  case 0, % regular pdf is returned  
085   p = wweibpdf(h,1.05,2.39).*wlognpdf(s,my(h),sig); 
086  case 1, %pdf conditioned on x1 ie. p(x2|x1)  
087   p = wlognpdf(s,my(h),sig); 
088  case 3, % secret option  used by mk87stat: returns x2*p(x2|x1)  
089   p = s.*wlognpdf(s,my(h),sig); 
090  case 4, % secret option  used by mk87stat: returns x2.^2*p(x2|x1)  
091   p = s.^2.*wlognpdf(s,my(h),sig); 
092  case 5, % p(h)*P(V|h) is returned special case used by mk87cdf 
093   p = wweibpdf(h,1.05,2.39).*wlogncdf(s,my(h),sig); 
094  case 6, % P(V|h) is returned special case used by mk87cdf 
095   p = wlogncdf(s,my(h),sig); 
096  case 7, % p(h)*(1-P(V|h)) is returned special case used by mk87cdf 
097   p = wweibpdf(h,1.05,2.39).*(1-wlogncdf(s,my(h),sig)); 
098  otherwise error('unknown option') 
099 end 
100 if condon~=6 
101   p = p./Hrms./Erms; 
102 end 
103 return 
104  
105 function y = my(h) 
106 y   = zeros(size(h)); 
107 ind = (h <= 1.7); 
108 h1  = h(ind); 
109 y(ind) = 0.024-1.065.*h1+0.585.*h1.^2; 
110  
111 %ind=find(h > 1.7); 
112 %h2=h(~ind); 
113 y(~ind) = 0.32*atan(3.14*(h(~ind)-1.7)) - 0.096; 
114 return 
115

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