OHHSSPDF Joint (Scf,Hd) PDF linear waves in space with Ochi-Hubble spectra. CALL: f = ohhsspdf(Hd,Scf,Hm0,Tp) f = pdf evaluated at (Scf,Hd). Hd = zero down crossing wave height Scf = crest front steepness Hm0 = significant wave height [m]. def = defines the parametrization of the spectral density (default 1) 1 : The most probable spectrum (default) 2,3,...11 : gives 95% Confidence spectra OHHSSPDF approximates the joint distribution of (Scf, Hd), i.e., crest front steepness (Ac/Lcf) and wave height in space, for a Gaussian process with a bimodal Ochi-Hubble spectral density. The empirical parameters of the model is fitted by least squares to simulated (Scf,Hd) data for 24 classes of Hm0. Between 50000 and 300000 zero-downcrossing waves were simulated for each class of Hm0. OHHSSPDF is restricted to the following range for Hm0: 0.5 < Hm0 [m] < 12 The size of f is the common size of the input arguments, Hd, Scf and Hm0. Example: Hm0 = 6;Tp = 8;def= 2; h = linspace(0,4*Hm0/sqrt(2))'; v = linspace(0,4*2*Hm0/Tp)'; [V,H] = meshgrid(v,h); f = ohhsspdf(H,V,Hm0,def); w = linspace(0,10,2*1024+1).'; S = ohspec2(w,[Hm0 def]); Sk = spec2spec(specinterp(S,.55),'k1d'); dk = 1; x = spec2sdat(Sk,80000,dk); rate = 8; [vi,hi] = dat2steep(x,rate,1); fk = kdebin([vi,hi],'epan',[],[],.5,128); plot(vi,hi,'.'), hold on contour(v,h,f,fk.cl,'g'), pdfplot(fk,'r'), hold off See also ohhpdf, thvpdf
PDF class constructor | |
Wave height, Hd, distribution parameters for Ochi-Hubble spectra. | |
Calculates the difference between the maximum and minimum values. | |
Calculates a smoothing spline. | |
Generalized Gamma probability density function | |
Weibull cumulative distribution function | |
Weibull probability density function | |
Check if all input arguments are either scalar or of common size. | |
Display message and abort function. | |
1-D interpolation (table lookup) | |
2-D interpolation (table lookup). | |
X and Y arrays for 3-D plots. | |
Write formatted data to string. | |
Remove singleton dimensions. | |
Set unique. | |
Display warning message; disable or enable warning messages. |
Joint (Scf,Hd) CDF for linear waves in space with Ochi-Hubble spectra. | |
Joint (Scf,Hd) PDF linear waves in space with Ochi-Hubble spectra. |
001 function [f,Hrms,Vrms,fA,fB] = ohhsspdf(Hd,Scf,Hm0,def,normalizedInput,condon) 002 %OHHSSPDF Joint (Scf,Hd) PDF linear waves in space with Ochi-Hubble spectra. 003 % 004 % CALL: f = ohhsspdf(Hd,Scf,Hm0,Tp) 005 % 006 % f = pdf evaluated at (Scf,Hd). 007 % Hd = zero down crossing wave height 008 % Scf = crest front steepness 009 % Hm0 = significant wave height [m]. 010 % def = defines the parametrization of the spectral density (default 1) 011 % 1 : The most probable spectrum (default) 012 % 2,3,...11 : gives 95% Confidence spectra 013 % 014 % OHHSSPDF approximates the joint distribution of (Scf, Hd), i.e., crest 015 % front steepness (Ac/Lcf) and wave height in space, for a Gaussian 016 % process with a bimodal Ochi-Hubble spectral density. The empirical 017 % parameters of the model is fitted by least squares to simulated 018 % (Scf,Hd) data for 24 classes of Hm0. 019 % Between 50000 and 300000 zero-downcrossing waves were simulated for 020 % each class of Hm0. 021 % OHHSSPDF is restricted to the following range for Hm0: 022 % 0.5 < Hm0 [m] < 12 023 % The size of f is the common size of the input arguments, Hd, Scf and 024 % Hm0. 025 % 026 % Example: 027 % Hm0 = 6;Tp = 8;def= 2; 028 % h = linspace(0,4*Hm0/sqrt(2))'; 029 % v = linspace(0,4*2*Hm0/Tp)'; 030 % [V,H] = meshgrid(v,h); 031 % f = ohhsspdf(H,V,Hm0,def); 032 % w = linspace(0,10,2*1024+1).'; 033 % S = ohspec2(w,[Hm0 def]); 034 % Sk = spec2spec(specinterp(S,.55),'k1d'); 035 % dk = 1; 036 % x = spec2sdat(Sk,80000,dk); rate = 8; 037 % [vi,hi] = dat2steep(x,rate,1); 038 % fk = kdebin([vi,hi],'epan',[],[],.5,128); 039 % plot(vi,hi,'.'), hold on 040 % contour(v,h,f,fk.cl,'g'), 041 % pdfplot(fk,'r'), hold off 042 % 043 % See also ohhpdf, thvpdf 044 045 046 % Reference 047 % P. A. Brodtkorb (2004), 048 % The Probability of Occurrence of Dangerous Wave Situations at Sea. 049 % Dr.Ing thesis, Norwegian University of Science and Technolgy, NTNU, 050 % Trondheim, Norway. 051 052 % History 053 % revised pab jan2004 054 % By pab 20.12.2000 055 056 057 058 error(nargchk(3,6,nargin)) 059 if (nargin < 6|isempty(condon)), condon = 0;end 060 if (nargin < 5|isempty(normalizedInput)), normalizedInput = 0;end 061 if (nargin < 4|isempty(def)), def=1;end 062 063 multipleSeaStates = any(prod(size(Hm0))>1); 064 if multipleSeaStates 065 [errorcode, Scf,Hd,Hm0] = comnsize(Scf,Hd,Hm0); 066 else 067 [errorcode, Scf,Hd] = comnsize(Scf,Hd); 068 end 069 if errorcode > 0 070 error('Requires non-scalar arguments to match in size.'); 071 end 072 073 if any(Hm0>12| Hm0<=0.5) 074 disp('Warning: Hm0 is outside the valid range') 075 disp('The validity of the Hd distribution is questionable') 076 end 077 078 if def>11|def<1 079 Warning('DEF is outside the valid range') 080 def = mod(def-1,11)+1; 081 end 082 083 global OHHSSPAR 084 if isempty(OHHSSPAR) 085 OHHSSPAR = load('ohhsspar.mat'); 086 end 087 088 089 %w = linspace(0,100,16*1024+1).'; %original spacing 090 %ch = spec2char(ohspec2(w,[Hm0,def]),{'Tm02','eps2'}); 091 %Tm02 = ch(1); 092 %Vrms = 2*Hm0/Tm02; 093 Hm00 = OHHSSPAR.Hm0; 094 095 if normalizedInput, 096 Hrms = 1; 097 Vrms = 1; 098 else 099 method = 'cubic'; %'spline' 100 Tm020 = OHHSSPAR.Tm02; 101 Hrms = Hm0/sqrt(2); 102 Tm02 = interp1(Hm00,Tm020(:,def),Hm0,method); 103 Vrms = 2*Hm0./Tm02; % Srms 104 end 105 106 107 h = Hd./Hrms; 108 v = Scf./Vrms; 109 cSize = size(h); % common size 110 111 % Logarithm of Weibull distribution parameters as a function of Hm0 and h2 112 A11 = squeeze(OHHSSPAR.A11s(:,def,:)); 113 B11 = squeeze(OHHSSPAR.B11s(:,def,:)); 114 115 h2 = OHHSSPAR.h2; %Hd/Hrms 116 [E2 H2] = meshgrid(Hm00,h2); 117 method = '*cubic'; %'spline' 118 if multipleSeaStates 119 h = h(:); 120 v = v(:); 121 Hm0 = Hm0(:); 122 A1 = zeros(length(h),1); 123 B1 = A1; 124 [Hm0u,ix,jx] = unique(Hm0); 125 numSeaStates = length(ix); 126 Nh2 = length(h2); 127 Hm0i = zeros(Nh2,1); 128 for iz=1:numSeaStates 129 k = find(jx==iz); 130 Hm0i(:) = Hm0u(iz); 131 A1(k) = exp(smooth(h2,interp2(E2,H2,log(A11.'),Hm0i,h2,method),... 132 1,h(k),1)); 133 B1(k) = (smooth(h2,interp2(E2,H2,B11.',Hm0i,h2,method),... 134 1,h(k),1)); 135 end 136 else 137 Hm0i = Hm0(ones(size(h2))); 138 A1 = exp(smooth(h2,interp2(E2,H2,log(A11.'), ... 139 Hm0i,h2,method),1,h,1)); 140 B1 = (smooth(h2,interp2(E2,H2,B11.', ... 141 Hm0i,h2,method),1,h,1)); 142 end 143 %fh = ohhpdf(h(:)/Hrms,Hm0,def,'time',1); 144 % Fh = fh.f; 145 % Waveheight distribution in time 146 % Generalized gamma distribution parameters as a function of Hm0 147 [A0 B0 C0] = ohhgparfun(Hm0,def,'space'); 148 149 150 switch condon, 151 case 0, % regular pdf is returned 152 f = wggampdf(h,A0,B0,C0).*wweibpdf(v,A1,B1); 153 case 1, %pdf conditioned on x1 ie. p(x2|x1) 154 f = wweibpdf(v,A1,B1); 155 case 3, % secret option used by XXstat: returns x2*p(x2|x1) 156 f = v.*wweibpdf(v,A1,B1); 157 case 4, % secret option used by XXstat: returns x2.^2*p(x2|x1) 158 f = v.^2.*wweibpdf(v,A1,B1); 159 case 5, % p(h)*P(V|h) is returned special case used by thscdf 160 f = wggampdf(h,A0,B0,C0).*wweibcdf(v,A1,B1); 161 case 6, % P(V|h) is returned special case used by thscdf 162 f = wweibcdf(v,A1,B1); 163 case 7,% p(h)*(1-P(V|h)) is returned special case used by thscdf 164 f = wggampdf(h,A0,B0,C0).*(1-wweibcdf(v,A1,B1)); 165 otherwise error('unknown option') 166 end 167 if multipleSeaStates 168 f = reshape(f,cSize); 169 end 170 171 if condon~=6 172 f = f./Hrms./Vrms; 173 end 174 f(find(isnan(f)|isinf(f) ))=0; 175 if any(size(f)~=cSize) 176 disp('Wrong size') 177 end 178 if nargout>3, 179 fA = createpdf(2); 180 fA.f = A11; 181 fA.x = {Hm00, h2}; 182 fA.labx = {'Hm0', 'h'}; 183 fA.note = sprintf('%s \n',... 184 'The conditonal Weibull distribution Parameter',... 185 'A(h)/Hrms for wave heigth as a function of',... 186 'h=Hd/Hrms and Hm0 for the ohspec2 spectrum, ',... 187 sprintf('given def = %d.',def)); 188 189 ra = range(A11(:)); 190 st = round(min(A11(:))*100)/100; 191 en = max(A11(:)); 192 fA.cl = st:ra/20:en; 193 end 194 if nargout>4, 195 fB = createpdf(2); 196 fB.f = B11; 197 fB.x = {Hm00,h2}; 198 fB.labx = {'Hm0','h'}; 199 fB.note = sprintf('%s \n',... 200 'The conditonal Weibull distribution Parameter',... 201 'B(h)/Hrms for wave heigth as a function of',... 202 'h=Hd/Hrms and Hm0 for the ohspec2 spectrum, ',... 203 sprintf('given def = %d.',def)); 204 ra = range(B11(:)); 205 st = round(min(B11(:))*100)/100; 206 en = max(B11(:)); 207 fB.cl = st:ra/20:en; 208 end 209 return 210 211
Comments or corrections to the WAFO group