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
Lognormal cumulative distribution function | |
Lognormal probability density function | |
Weibull probability density function | |
Check if all input arguments are either scalar or of common size. | |
Display message and abort function. |
Myrhaug and Kjeldsen (1987) joint (Scf,Hd) CDF. | |
Myrhaug and Kjeldsen (1987) joint (Scf,Hd) PDF. |
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
Comments or corrections to the WAFO group