Home > wafo > wavemodels > ltwcpdf.m

ltwcpdf

PURPOSE

Long Term Wave Climate PDF of significant wave height and wave period

SYNOPSIS

p = ltwcpdf(t,h,Pt1,Ph1, def,condon)

DESCRIPTION

``` LTWCPDF Long Term Wave Climate PDF of significant wave height and wave period

CALL: f = ltwcpdf(T,Hs,Pt,Ph,def)

f      = density
T      = peak period, Tp  or mean zerodowncrossing period, Tz (sec)
Hs     = significant wave height (meter)
Pt     = [A1 A2 A3 B1 B2 B3 B4] parameters for the conditional
distribution of T given Hs
(default [1.59 0.42 2 0.005 0.09 0.13 1.34])
Ph     = [C1 C2 C3 Hc M1 S1], parameters for the marginal
distribution for Hs (default [2.82 1.547 0 3.27  0.836 0.376])
def    = defines the parametrization used (default 0)

The probability distribution is given as the product:

f(T,Hs)=f(T|Hs)*f(Hs)
where

f(T|Hs) = wlognpdf(T,my(Hs),S(Hz)) with mean and variance
my(Hs)  = A1 + A2*log(Hs+A3)          if mod(def,2)==0
A1 + A2*Hs^A3               if mod(def,2)==1
S(Hs)   = B1+B2*exp(-B3*Hs.^B4)

For def <=1:
f(Hs)    = wlognpdf(Hs,M1,S1)       for Hs <= Hc
= wweibpdf(Hs-C3,C1,C2)    for Hs >  Hc
For def >1:
f(Hs)    = wlognpdf(Hs,M1,S1)       for Hs <= Hc
= wggampdf(Hs,C1,C2,C3)    for Hs >  Hc  (Suggested by Ochi)

The default values for T and Hs are suitable for the Northern North
Sea for peak period, Tp and significant waveheight, Hs.
With a suitable change of parameters for Pt this model fit mean wave
period, Tz, also reasonably well.

The size of f is the common size of T and Hs

NOTE: - by specifying nan's  in the vectors Pt or Ph default values
will be used.
- if length(Pt) or length(Ph) is shorter than the parameters
needed then the default values are used for the parameters
not specified.
- For tables of fitted parameter values to specific sites see
the end of this file (type ltwcpdf.m)

Example: % Set  C1 = 2.73,  Hc = 3.95  and the rest to their default values
Ph = [2.73, NaN,NaN, 3.95]; x = linspace(0, 15);
[T,H] = meshgrid(x);
f = ltwcpdf(T,H,[],Ph);
contour(x,x,f)

CROSS-REFERENCE INFORMATION

This function calls:
 wggampdf Generalized Gamma probability density function wlogncdf Lognormal cumulative distribution function wlognpdf Lognormal probability density function wweibpdf Weibull probability density function comnsize Check if all input arguments are either scalar or of common size. error Display message and abort function.
This function is called by:

SOURCE CODE

```001 function p = ltwcpdf(t,h,Pt1,Ph1, def,condon)
002 %LTWCPDF Long Term Wave Climate PDF of significant wave height and wave period
003 %
004 % CALL: f = ltwcpdf(T,Hs,Pt,Ph,def)
005 %
006 %    f      = density
007 %    T      = peak period, Tp  or mean zerodowncrossing period, Tz (sec)
008 %    Hs     = significant wave height (meter)
009 %    Pt     = [A1 A2 A3 B1 B2 B3 B4] parameters for the conditional
010 %             distribution of T given Hs
011 %             (default [1.59 0.42 2 0.005 0.09 0.13 1.34])
012 %    Ph     = [C1 C2 C3 Hc M1 S1], parameters for the marginal
013 %             distribution for Hs (default [2.82 1.547 0 3.27  0.836 0.376])
014 %    def    = defines the parametrization used (default 0)
015 %
016 %   The probability distribution is given as the product:
017 %
018 %    f(T,Hs)=f(T|Hs)*f(Hs)
019 %   where
020 %
021 %       f(T|Hs) = wlognpdf(T,my(Hs),S(Hz)) with mean and variance
022 %        my(Hs)  = A1 + A2*log(Hs+A3)          if mod(def,2)==0
023 %                  A1 + A2*Hs^A3               if mod(def,2)==1
024 %        S(Hs)   = B1+B2*exp(-B3*Hs.^B4)
025 %
026 %    For def <=1:
027 %       f(Hs)    = wlognpdf(Hs,M1,S1)       for Hs <= Hc
028 %                = wweibpdf(Hs-C3,C1,C2)    for Hs >  Hc
029 %    For def >1:
030 %       f(Hs)    = wlognpdf(Hs,M1,S1)       for Hs <= Hc
031 %                = wggampdf(Hs,C1,C2,C3)    for Hs >  Hc  (Suggested by Ochi)
032 %
033 %    The default values for T and Hs are suitable for the Northern North
034 %    Sea for peak period, Tp and significant waveheight, Hs.
035 %    With a suitable change of parameters for Pt this model fit mean wave
036 %    period, Tz, also reasonably well.
037 %
038 %   The size of f is the common size of T and Hs
039 %
040 %   NOTE: - by specifying nan's  in the vectors Pt or Ph default values
041 %           will be used.
042 %         - if length(Pt) or length(Ph) is shorter than the parameters
043 %           needed then the default values are used for the parameters
044 %           not specified.
045 %         - For tables of fitted parameter values to specific sites see
046 %           the end of this file (type ltwcpdf.m)
047 %
048 % Example: % Set  C1 = 2.73,  Hc = 3.95  and the rest to their default values
049 %   Ph = [2.73, NaN,NaN, 3.95]; x = linspace(0, 15);
050 %   [T,H] = meshgrid(x);
051 %   f = ltwcpdf(T,H,[],Ph);
052 %   contour(x,x,f)
053 %
054 % See also  wweibpdf, wlognpdf, wggampdf
055
056 %   References:
057 %   Haver, S (1980)
058 %  'Analysis of uncertainties related to the stochastic modelling of
059 %    Ocean waves'
060 %   Ph.D. thesis, Norwegian Institute of Technology, NTH, Trondheim, Norway
061 %
062 %   Haver,S and Nyhus, K. A.  (1986)
063 %   'A wave climate description for long term response calculation.'
064 %   In Proc. of OMAE'86, tokyo, Japan, pp. 27-34
065 %
066 %  Sagli, Gro (2000)
067 %  "Model uncertainty and simplified estimates of long term extremes of
068 %  hull girder loads in ships"
069 %  Ph.D. thesis, Norwegian University of Science and Technology, NTNU,
070 %  Trondheim, Norway, pp 45--47
071 %
072 %   Michel K. Ochi (1998),
073 %   "OCEAN WAVES, The stochastic approach",
074 %   OCEAN TECHNOLOGY series 6, Cambridge, pp 127.
075 % (Generalized Gamma distribution for wave height)
076 %
077 %  Bitner-Gregersen, E.M. and Hagen, ุ. (2000)
078 % "Aspects of joint distribution for metocean Phenoma at the Norwegian
079 % Continental shelf."
080 % In Proc. of OMAE'2000,
081
082 % tested on: matlab 5.2
083 % history:
084 % revised pab 04.11.2000
085 %  -updated calls to new wstats functions
086 % by  Per A. Brodtkorb 13.05.2000
087
088
089
090
091
092 error(nargchk(2,6,nargin))
093
094 if nargin < 5|isempty(def),  def=0; end
095 if nargin < 6|isempty(condon),  condon=0; end
096
097 Ph=[2.82 1.547 0 3.27  0.836 0.376]; % default values
098 if nargin < 4|isempty(Ph),
099
100 else
101   nh=length(Ph1);
102   ind=find(~isnan(Ph1(1:min(nh,6))));
103   if any(ind) % replace default values with those from input data
104     Ph(ind)=Ph1(ind);
105   end
106 end
107
108 Pt=[1.59 0.42 2 0.005 0.09 0.13 1.34]; % default values
109 if nargin < 3|isempty(Pt),
110 else
111   nt=length(Pt1);
112   ind=find(~isnan(Pt1(1:min(nt,7))));
113   if any(ind) % replace default values with those from input data
114     Pt(ind) = Pt1(ind);
115   end
116 end
117
118 [errorcode,t,h] = comnsize(t,h);
119 if errorcode > 0
120   error('Requires non-scalar arguments to match in size.');
121 end
122
123
124 % Model parameters for the marginal distribution for Hz:
125 %~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
126 % scale, shape and location parameter for weibull distribution
127 % or shape,shape and scale parameter for the generalized Gamma distribution
128 C1=Ph(1);C2=Ph(2);C3=Ph(3);
129 % Split value, mean and standard deviation for Lognormal distribution
130 Hc=Ph(4);myHz=Ph(5);sHz=Ph(6);
131
132 % Model parameters for the conditional distribution of Tp or Tz given Hz
133 %~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
134 A1=Pt(1);A2=Pt(2);A3=Pt(3);
135 B1=Pt(4);B2=Pt(5);B3=Pt(6);B4=Pt(7);
136 % Mean and Variance of the lognormal distribution
137 if mod(abs(def),2)==0
138   myTp = A1+A2*log(h+A3);
139 else
140   myTp = A1+A2*h.^A3;
141 end
142 sTp =(B1+B2*exp(-B3*h.^B4));
143
144 p=zeros(size(h));
145 k0=find(t>0 & h>0);
146 if any(k0)
147   if condon==5,
148      p(k0) = wlogncdf(t(k0),myTp(k0),sTp(k0));
149    elseif condon ==10
150      p(k0) = 1;
151    else
152     p(k0) = wlognpdf(t(k0),myTp(k0),sTp(k0));
153   end
154 end
155
156
157
158 if condon==0 | condon==5 | condon ==10
159   k=find(0<h & h<=Hc );
160   if any(k)
161     p(k)=p(k).*wlognpdf(h(k),myHz,sHz);
162   end
163   k1=find(h>Hc);
164   if any(k1)
165     if def<=1,
166       % NB! weibpdf must be modified to correspond to
167       % pdf=x^(b-1)/a^b*exp(-(x/a)^b)
168       p(k1)=p(k1).*wweibpdf(h(k1)-C3,C1,C2);
169     else
170       p(k1)=p(k1).*wggampdf(h(k1),C1,C2,C3);
171     end
172   end
173 end
174
175 switch condon,
176   case 0, % regular pdf is returned
177   case 1, %pdf conditioned on Hz ie. p(Tp|Hz)
178   case 3, % secret option  used by ltwcstat: returns Tp*p(Tp|Hz)
179     p = t.*p
180   case 4, % secret option  used by ltwcstat: returns Tp.^2*p(Tp|Hz)
181     p = t.^2.*p
182   case 5, % p(Hz)*P(Tp|Hz) is returned special case used by ltwccdf
183   case 10 % p(Hz) is returned
184
185   otherwise error('unknown option')
186 end
187
188
189
190
191 % The following is taken from Sagli (2000):
192 % Parameters for the long term description of the wave climate (Tp|Hs):
193 %~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
194 %
195 %   Area            Pt = [  A1     A2    A3     B1    B2    B3   B4]  def
196 %   จจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจ
197 %   Northern North Sea:  [1.59,  0.42, 2.00, 0.005, 0.09, 0.13, 1.34]  0
198 %   Aasgard           :  [1.72,  0.34, 0.46, 0.005, 0.10, 0.29, 1.00]  1
199 %   Sleipner          :  [0.23,  1.69, 0.15, 0.005, 0.12, 0.40, 1.00]  1
200 %   Troms I           :  [1.35,  0.61, 0.34, 0.005, 0.12, 0.36, 1.00]  1
201 %   Statfjord/Brent   :  [0.40,  1.59, 0.15, 0.005, 0.12, 0.34, 1.00]  1
202 %   Ekofisk           :  [1.00,  0.90, 0.25, 0.005, 0.10, 0.44, 1.00]  1
203 %   Ekofisk Zone II   :  [0.03,  1.81, 0.15, 0.005, 0.16, 0.58, 1.00]  1
204 %   Ekofisk Zone III  :  [0.03,  1.81, 0.15, 0.005, 0.16, 0.58, 1.00]  1
205 %   Ekofisk Zone IV   :  [1.41,  0.38, 0.54, 0.010, 0.08, 0.41, 1.00]  1
206 %   Ekofisk Zone VI   :  [0.00,  1.72, 0.15, 0.010, 0.14, 1.00, 1.00]  1
207 %
208 %
209 % Parameters for the long term description of the wave climate (Hs):
210 %~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
211 %
212 %   Area            Ph = [  C1     C2    C3    Hc    M1    S1  ]  def
213 %   จจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจ
214 %   Northern North Sea:  [2.82,  1.55, 0.00, 3.27, 0.84, 0.61^2  ]  0
215 %   Aasgard           :  [2.73,  1.45, 0.00, 3.95, 0.84, 0.59^2  ]  1
216 %   Sleipner          :  [2.38,  1.44, 0.00, 3.00, 0.69, 0.60^2  ]  1
217 %   Troms I           :  [2.26,  1.36, 0.00, 4.25, 0.72, 0.56^2  ]  1
218 %   Statfjord/Brent   :  [2.82,  1.55, 0.00, 3.27, 0.84, 0.61^2  ]  1
219 %   Ekofisk           :  [2.08,  1.37, 0.00, 2.65, 0.52, 0.67^2  ]  1
220 %   Ekofisk Zone II   :  [2.12,  1.37, 0.00, 2.70, 0.54, 0.67^2  ]  1
221 %   Ekofisk Zone III  :  [2.05,  1.55, 0.00, 2.15, 0.50, 0.65^2  ]  1
222 %   Ekofisk Zone IV   :  [1.56,  1.46, 0.00, 2.00, 0.25, 0.62^2  ]  1
223 %   Ekofisk Zone VI   :  [1.18,  1.42, 0.00, 1.25, 0.07, 0.70^2  ]  1
224 %
225
226
227 % The following is taken from Bitner-Gregersen and Hagen (2000):
228 % Description of sites:
229 % V๘ring plateau:  (position 67, 27'N 5, 58'E water depth 1460m) NORWAVE
230 %                   buoy every 3rd hour for a period of 2.3 years
231 %
232
233 % Parameters for the long term description of the wave climate (Tp|Hs):
234 %~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
235 %
236 %   Area           Pt = [  A1     A2      A3      B1     B2       B3     B4]  def
237 %   จจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจ
238 %   V๘ring M         :  [1.718,  0.194, 0.642, -288.28, 288.52,0.279e-5 , 1]  1
239 %   V๘ring M/H       :  [1.670,  0.220, 0.617, - 93.13,  93.37,0.330e-4 , 1]  1
240 %   V๘ring H         :  [0.474,  1.438, 0.186,   0.047,   0.46,   0.476 , 1]  1
241 %   Haltenbanken     :  [1.921,  0.171, 0.665,   0.064,  0.327,   0.227 , 1]  1
242
243 % Parameters for the long term description of the wave climate (Tz|Hs)
244 % (mean zero-down crossing period given Hz)                           :
245 %~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
246 %
247 %   Area            Pt = [  A1     A2      A3     B1    B2     B3   B4]  def
248 %   จจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจ
249 %   Statfjord&Gullfaks:  [1.790,  0.110, 0.759, 0.106, 0.361,0.969 , 1]  1
250 %   Statfjord Platform:  [1.771,  0.080, 0.902, 0.001, 0.239,0.169 , 1]  1
251 %   Ekofisk           :  [0.611,  0.902, 0.289, 0.071, 0.211,0.606 , 1]  1
252
253
254 % Parameters for the long term description of the wave climate (Hs):
255 %~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
256 %
257 %   Area            Ph = [  C1     C2    C3    Hc  ]  def
258 %   จจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจ
259 %   V๘ring M          :  [2.371,  1.384, 0.755, 0  ]  1
260 %   V๘ring M/H        :  [2.443,  1.405, 0.820, 0  ]  1
261 %   V๘ring H          :  [2.304,  1.383, 0.899, 0  ]  1
262 %   Haltenbanken      :  [2.154,  1.273, 0.763, 0  ]  1
263 %   Statfjord&Gullfaks:  [2.264,  1.398, 0.969, 0  ]  1
264 %   Statfjord Platform:  [2.502,  1.527, 0.657, 0  ]  1
265 %   Ekofisk           :  [1.597,  1.184, 0.852, 0  ]  1
266
267
268 % Fitted parameters by splitting the data into wind Sea and swell:
269
270 % Parameters for the long term description of the wave climate (Tp|Hs):
271 %~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
272 %
273 %   Area            Pt = [  A1     A2    A3     B1    B2    B3   B4]  def
274 %   จจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจ
275 %   V๘ring M         :  [-0.707, 2.421, 0.126, 0.076, 0.063, 0.162, 1]  1    (Wind Sea)
276 %   V๘ring M         :  [ 1.501, 0.882, 0,135, 0.126, 0.000, 0.000, 1]  1    (Swell)
277 %   Haltenbanken     :  [-4383.9, 4385.6, 0.948e-4, 0.045, 0.369,0.780,1]  1 (Wind Sea)
278 %   Haltenbanken     :  [ 2.257, 0.088, 0.746, 0.010, 0.229, 0.167, 1]  1    (Swell)
279
280 % Parameters for the long term description of the wave climate (Hs):
281 %~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
282 %
283 %   Area            Ph = [  C1     C2    C3    Hc  ]  def
284 %   จจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจจ
285 %   V๘ring M          :  [2.170,  1.382, 0.203, 0  ]  1  (Wind sea)
286 %   V๘ring M          :  [1.606,  1.252, 0.229, 0  ]  1  (Swell)
287 %   Haltenbanken      :  [1.473,  1.032, 0.241, 0  ]  1  (Wind sea)
288 %   Haltenbanken      :  [1.661,  1.117, 0.240, 0  ]  1  (Swell)
289
290
291
292
293```

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