BMR00CDF Brodtkorb et.al. (2000) joint (Scf,Hd) CDF from North Sea. CALL: F = bmr00cdf(Hd,Scf,Hs,Tz,tail) F = CDF Hd = zero down crossing wave height. Scf = crest front steepness. Hs = significant wave height. Tz = average zero down crossing period. tail = 1 if upper tail is calculated 0 if lower tail is calulated (default) BMR00CDF returns the joint CDF 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 = Hs/sqrt(2); Erms = 5/4*Hs/(Tz^2); The size of F is the common size of the input arguments Examples: Hs = 5.5; Tz = 8.5; sc = 0.25; hc = 4; p = bmr00cdf(hc,sc,Hs,Tz) % Prob(Hd<=hc,Scf<=sc|Hs,Tz) = 0.68 % Conditional probability of steep and high waves given seastates % i.e., Prob(Hd>hc,Scf>sc|Hs,Tz) upperTail = 1; Hs = linspace(2.5,18.5,17); Tz = linspace(4.5,19.5,16); [T,H] = meshgrid(Tz,Hs); p = bmr00cdf(hc,sc,H,T,upperTail); v = 10.^(-6:-1); contourf(Tz,Hs,log10(p),log10(v)) xlabel('Tz') ylabel('Hs') fcolorbar(log10(v)) See also mk87cdf, gaussq
Brodtkorb et.al (2000) joint (Scf,Hd) PDF from North Sea. | |
Numerically evaluates a integral using a Gauss quadrature. | |
Check if all input arguments are either scalar or of common size. | |
Display message and abort function. |
001 function [p, eps2, a, b] = bmr00cdf(Hd,Scf,Hs,Tz,tail) 002 %BMR00CDF Brodtkorb et.al. (2000) joint (Scf,Hd) CDF from North Sea. 003 % 004 % CALL: F = bmr00cdf(Hd,Scf,Hs,Tz,tail) 005 % 006 % F = CDF 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 % tail = 1 if upper tail is calculated 012 % 0 if lower tail is calulated (default) 013 % 014 % BMR00CDF returns the joint CDF of (Scf, Hd) given Hs and Tz, 015 % i.e., crest front steepness (2*pi*Ac/(g*Td*Tcf)) and wave height, given 016 % the seastate. The root mean square values of Hd and Scf (Hrms,Erms) are 017 % related to the significant waveheight and the average zero down 018 % crossing period by: 019 % Hrms = Hs/sqrt(2); 020 % Erms = 5/4*Hs/(Tz^2); 021 % 022 % The size of F is the common size of the input arguments 023 % 024 % Examples: 025 % Hs = 5.5; 026 % Tz = 8.5; 027 % sc = 0.25; 028 % hc = 4; 029 % p = bmr00cdf(hc,sc,Hs,Tz) % Prob(Hd<=hc,Scf<=sc|Hs,Tz) = 0.68 030 % 031 % % Conditional probability of steep and high waves given seastates 032 % % i.e., Prob(Hd>hc,Scf>sc|Hs,Tz) 033 % upperTail = 1; 034 % Hs = linspace(2.5,18.5,17); 035 % Tz = linspace(4.5,19.5,16); 036 % [T,H] = meshgrid(Tz,Hs); 037 % p = bmr00cdf(hc,sc,H,T,upperTail); 038 % v = 10.^(-6:-1); 039 % contourf(Tz,Hs,log10(p),log10(v)) 040 % xlabel('Tz') 041 % ylabel('Hs') 042 % fcolorbar(log10(v)) 043 % 044 % See also mk87cdf, gaussq 045 046 % Reference 047 % Brodtkorb, P.A. and Myrhaug, D. and Rue, H. (2000) 048 % "Joint Distributions of Wave Height and Wave Steepness Parameters", 049 % In Proc. 27'th Int. Conf. on Coastal Eng., ICCE, Sydney, Australia }, 050 % vol. 1, pp. 545--558, Paper No. 162 051 052 053 %tested on matlab 5.2 054 %history: 055 % by Per A. brodtkorb July 2004 056 057 058 error(nargchk(3,5,nargin)) 059 060 if (nargin < 5|isempty(tail)),tail = 0; end 061 if (nargin < 4|isempty(Tz)),Tz = 8; end 062 if (nargin < 3|isempty(Hs)), Hs = 6; end 063 064 multipleSeaStates = any(prod(size(Hs))>1|prod(size(Tz))>1); 065 if multipleSeaStates 066 [errorcode, Scf,Hd,Tz,Hs] = comnsize(Scf,Hd,Tz,Hs); 067 else 068 [errorcode, Scf,Hd] = comnsize(Scf,Hd); 069 end 070 if errorcode > 0 071 error('Requires non-scalar arguments to match in size.'); 072 end 073 074 Hrms = Hs/sqrt(2); 075 Erms = 5/4*Hs./(Tz.^2); 076 %Erms = (0.0202 + 0.826*Hs./(Tz.^2)); 077 078 s = Scf./Erms; 079 hMax = 20; 080 h = min(Hd./Hrms,hMax); 081 082 eps2 = 1e-6; 083 084 085 p = zeros(size(Hd)); 086 087 %k0 = find((Hs<=(Tz-4)*13/6+4)); 088 upLimit = 6.5/1.4; 089 loLimit = 2.5/1.26;; 090 k0 = find((loLimit*sqrt(Hs)<Tz)); 091 if any(k0) 092 if multipleSeaStates 093 h = h(k0); 094 s = s(k0); 095 Hs = Hs(k0); 096 Tz = Tz(k0); 097 else 098 k0 = 1:prod(size(Hd)); 099 end 100 hlim = h; 101 102 h0 = 3.79937912772197; 103 normalizedInput = 1; 104 lowerTail = 0; 105 106 107 if 0 108 % This is a trick to get the html documentation correct. 109 k = bmr00pdf(1,1,2,3); 110 end 111 112 113 if (tail == lowerTail) 114 k = find(h>h0); 115 hlim(k) = h0; 116 p(k0) = gaussq('bmr00pdf',0,hlim,eps2/2,[],s,Hs,Tz,5,normalizedInput)... 117 + gaussq('bmr00pdf',hlim,h,eps2/2,[],s,Hs,Tz,5,normalizedInput); 118 else 119 k = find(h<h0); 120 hlim(k) = h0; 121 p(k0) = gaussq('bmr00pdf',h,hlim,eps2/2,[],s,Hs,Tz,7,normalizedInput)... 122 + gaussq('bmr00pdf',hlim,hMax,eps2/2,[],s,Hs,Tz,7,normalizedInput); 123 end 124 end 125 return 126 127
Comments or corrections to the WAFO group