% CHAPTER3 Demonstrates distributions of wave characteristics Chapter3 contains the commands used in Chapter3 in the tutorial. Some of the commands are edited for fast computation. Each set of commands is followed by a 'pause' command.
Cavanie et al. (1976) approximation of the density (Tc,Ac) | |
Estimate one-sided spectral density from data. | |
Estimate one-sided spectral density, version 2. | |
Extracts waveheights and steepnesses from data. | |
Extracts troughs and crests from data. | |
Extracts turning points from data, | |
Estimate transformation, g, from data. | |
Extracts sequence of wavelengths from data. | |
Extracts exact level v crossings | |
Computes and plots the empirical CDF | |
Transforms xx using the inverse of transformation g. | |
Calculate transformation, g, proposed by Winterstein | |
Calculates (and plots) a JONSWAP spectral density | |
Kernel Density Estimator. | |
Binned Kernel Density Estimator. | |
Create or alter KDE OPTIONS structure. | |
Longuet-Higgins (1983) approximation of the density (Tc,Ac) | |
Make a directional spectrum | |
Plot contents of pdf structures | |
reconstruct the spurious/missing points of timeseries | |
Create or alter RIND OPTIONS structure. | |
Rotate spectrum anti-clockwise around the origin. | |
Numerical integration with the Simpson method | |
Evaluates some spectral bandwidth and irregularity factors | |
Evaluates spectral characteristics and their covariance | |
Calculates joint density of Maximum, minimum and period. | |
Calculates spectral moments from spectrum | |
Estimates the moments of 2'nd order non-linear waves | |
Evaluates densities of wave period Tcc, wave lenght Lcc. | |
Joint density of amplitude and period/wave-length characteristics | |
Evaluates densities for crest-,trough-period, length. | |
Directional spreading functions | |
Plots specified waves from a timeseries | |
Calculates a double peaked (swell + wind) spectrum | |
Calculates min2Max and Max2min cycles from a sequence of turning points | |
Transforms process X and up to four derivatives | |
Calculates transformed Rayleigh approximation for amplitudes | |
Prints a caption "made by WAFO" in current figure. | |
Plots a histogram | |
Plots data on a Normal distribution paper |
0001 %% CHAPTER3 Demonstrates distributions of wave characteristics 0002 % 0003 % Chapter3 contains the commands used in Chapter3 in the tutorial. 0004 % 0005 % Some of the commands are edited for fast computation. 0006 % Each set of commands is followed by a 'pause' command. 0007 % 0008 0009 % Tested on Matlab 5.3 0010 % History 0011 % Revised pab sept2005 0012 % Added sections -> easier to evaluate using cellmode evaluation. 0013 % Revised by pab Feb 2005 0014 % -updated calls to kdetools+spec2XXpdf programs 0015 % Created by GL July 12, 2000 0016 % from commands used in Chapter 3, written by IR 0017 % 0018 0019 0020 %% Section 3.2 Estimation of wave characteristics from data 0021 %% Example 1 0022 pstate = 'off'; 0023 0024 xx = load('sea.dat'); 0025 xx(:,2) = detrend(xx(:,2)); 0026 rate = 8; 0027 Tcrcr = dat2wa(xx,0,'c2c','tw',rate); 0028 Tc = dat2wa(xx,0,'u2d','tw',rate); 0029 disp('Block = 1'), pause(pstate) 0030 0031 %% Histogram of crestperiod compared to the kernel density estimate 0032 clf 0033 mean(Tc) 0034 max(Tc) 0035 t = linspace(0.01,8,200); 0036 L2 = 0; 0037 kopt = kdeoptset('L2',L2); 0038 ftc1 = kde(Tc,kopt,t); 0039 pdfplot(ftc1), hold on 0040 whisto(Tc,[],[],1) 0041 axis([0 8 0 0.5]), hold off 0042 wafostamp([],'(ER)') 0043 disp('Block = 2'), pause(pstate) 0044 0045 clf 0046 ftc2 = kdebin(Tc,kopt); 0047 disp('Block = 3'), pause(pstate) 0048 0049 %% Extreme waves - model check: the highest and steepest wave 0050 clf 0051 method = 0; 0052 rate = 8; 0053 [S, H, Ac, At, Tcf, Tcb, z_ind, yn] = ... 0054 dat2steep(xx,rate,method); 0055 disp('Block = 4'), pause(pstate) 0056 0057 clf 0058 [Smax indS]=max(S) 0059 [Amax indA]=max(Ac) 0060 spwaveplot(yn,[indA indS],'k.') 0061 wafostamp([],'(ER)') 0062 disp('Block = 5'), pause(pstate) 0063 0064 %% Does the highest wave contradict a transformed Gaussian model? 0065 clf 0066 inds1 = (5965:5974)'; 0067 Nsim = 10; 0068 [y1, grec1, g2, test, tobs, mu1o, mu1oStd] = ... 0069 reconstruct(xx,inds1,Nsim); 0070 spwaveplot(y1,indA-10) 0071 hold on 0072 plot(xx(inds1,1),xx(inds1,2),'+') 0073 lamb = 2.; 0074 muLstd = tranproc(mu1o-lamb*mu1oStd,fliplr(grec1)); 0075 muUstd = tranproc(mu1o+lamb*mu1oStd,fliplr(grec1)); 0076 plot (y1(inds1,1), [muLstd muUstd],'b-') 0077 wafostamp([],'(ER)') 0078 disp('Block = 6'), 0079 pause(pstate) 0080 0081 %% 0082 clf 0083 plot(xx(inds1,1),xx(inds1,2),'+'), hold on 0084 mu = tranproc(mu1o,fliplr(grec1)); 0085 plot(y1(inds1,1), mu) 0086 disp('Block = 7'), pause(pstate) 0087 0088 %% Crest height 0089 clf 0090 L2 = 0.6; 0091 wnormplot(Ac.^L2) 0092 0093 fac = kde(Ac,{'L2',L2},linspace(0.01,3,200)); 0094 pdfplot(fac) 0095 wafostamp([],'(ER)') 0096 simpson(fac.x{1},fac.f) 0097 disp('Block = 8'), pause(pstate) 0098 0099 %% Empirical wave amplitude CDF 0100 clf 0101 Fac = flipud(cumtrapz(fac.x{1},flipud(fac.f))); 0102 Fac = [fac.x{1} 1-Fac]; 0103 Femp = empdistr(Ac,Fac); 0104 axis([0 2 0 1]) 0105 wafostamp([],'(ER)') 0106 disp('Block = 9'), pause(pstate) 0107 0108 %% 0109 facr = trraylpdf(fac.x{1},'Ac',grec1); 0110 Facr = cumtrapz(facr.x{1},facr.f); 0111 hold on 0112 plot(facr.x{1},Facr,'.') 0113 axis([1.25 2.25 0.95 1]) 0114 wafostamp([],'(ER)') 0115 disp('Block = 10'), pause(pstate) 0116 0117 %% Joint pdf of crest period and crest amplitude 0118 clf 0119 kopt2 = kdeoptset('L2',0.5,'inc',256); 0120 Tc = Tcf+Tcb; 0121 fTcAc = kdebin([Tc Ac],kopt2); 0122 fTcAc.labx={'Tc [s]' 'Ac [m]'} 0123 pdfplot(fTcAc) 0124 hold on 0125 plot(Tc,Ac,'k.') 0126 hold off 0127 wafostamp([],'(ER)') 0128 disp('Block = 11'), pause(pstate) 0129 0130 %% Example 4: Simple wave characteristics obtained from Jonswap spectrum 0131 clf 0132 S = jonswap([],[5 10]); 0133 [m, mt]= spec2mom(S,4,[],0); 0134 disp('Block = 12'), pause(pstate) 0135 0136 clf 0137 spec2bw(S) 0138 [ch Sa2] = spec2char(S,[1 3]) 0139 disp('Block = 13'), pause(pstate) 0140 0141 %% Section 3.3.2 Explicit form approximations of wave characteristic densities 0142 %% Longuett-Higgins model for Tc and Ac 0143 clf 0144 t = linspace(0,15,100); 0145 h = linspace(0,6,100); 0146 flh = lh83pdf(t,h,[m(1),m(2),m(3)]); 0147 disp('Block = 14'), pause(pstate) 0148 0149 %% Transformed Longuett-Higgins model for Tc and Ac 0150 clf 0151 [sk, ku ]=spec2skew(S); 0152 sa = sqrt(m(1)); 0153 gh = hermitetr([],[sa sk ku 0]); 0154 flhg = lh83pdf(t,h,[m(1),m(2),m(3)],gh); 0155 disp('Block = 15'), pause(pstate) 0156 0157 %% Cavanie model for Tc and Ac 0158 clf 0159 t = linspace(0,10,100); 0160 h = linspace(0,7,100); 0161 fcav = cav76pdf(t,h,[m(1) m(2) m(3) m(5)],[]); 0162 disp('Block = 16'), pause(pstate) 0163 0164 %% Example 5 Transformed Rayleigh approximation of crest height 0165 clf 0166 xx = load('sea.dat'); 0167 x = xx; 0168 x(:,2) = detrend(x(:,2)); 0169 SS = dat2spec2(x); 0170 [sk, ku, me, si ] = spec2skew(SS); 0171 gh = hermitetr([],[si sk ku me]); 0172 Hs = 4*si; 0173 r = (0:0.05:1.1*Hs)'; 0174 fac_h = trraylpdf(r,'Ac',gh); 0175 fat_h = trraylpdf(r,'At',gh); 0176 h = (0:0.05:1.7*Hs)'; 0177 facat_h = trraylpdf(h,'AcAt',gh); 0178 pdfplot(fac_h) 0179 hold on 0180 pdfplot(fat_h) 0181 hold off 0182 wafostamp([],'(ER)') 0183 disp('Block = 17'), pause(pstate) 0184 0185 %% 0186 clf 0187 TC = dat2tc(xx, me); 0188 tc = tp2mm(TC); 0189 Ac = tc(:,2); 0190 At = -tc(:,1); 0191 AcAt = Ac+At; 0192 disp('Block = 18'), pause(pstate) 0193 0194 %% 0195 clf 0196 Fac_h = [fac_h.x{1} cumtrapz(fac_h.x{1},fac_h.f)]; 0197 subplot(3,1,1) 0198 Fac = empdistr(Ac,Fac_h); 0199 hold on 0200 plot(r,1-exp(-8*r.^2/Hs^2),'.') 0201 axis([1. 2. 0.9 1]) 0202 Fat_h = [fat_h.x{1} cumtrapz(fat_h.x{1},fat_h.f)]; 0203 subplot(3,1,2) 0204 Fat = empdistr(At,Fat_h); 0205 hold on 0206 plot(r,1-exp(-8*r.^2/Hs^2),'.') 0207 axis([1. 2. 0.9 1]) 0208 Facat_h = [facat_h.x{1} cumtrapz(facat_h.x{1},facat_h.f)]; 0209 subplot(3,1,3) 0210 Facat = empdistr(AcAt,Facat_h); 0211 hold on 0212 plot(r,1-exp(-2*r.^2/Hs^2),'.') 0213 axis([1.5 3.5 0.9 1]) 0214 wafostamp([],'(ER)') 0215 disp('Block = 19'), pause(pstate) 0216 0217 %% Section 3.4 Exact wave distributions in transformed Gaussian Sea 0218 %% Section 3.4.1 Density of crest period, crest length or encountered crest period 0219 clf 0220 S1 = torsethaugen([],[6 8],1); 0221 D1 = spreading(101,'cos',pi/2,[15],[],0); 0222 D12 = spreading(101,'cos',0,[15],S1.w,1); 0223 SD1 = mkdspec(S1,D1); 0224 SD12 = mkdspec(S1,D12); 0225 disp('Block = 20'), pause(pstate) 0226 0227 %% Crest period 0228 clf 0229 f_tc = spec2tpdf(S1,[],'Tc',[0 11 56],[],4); 0230 pdfplot(f_tc) 0231 wafostamp([],'(ER)') 0232 simpson(f_tc.x{1},f_tc.f) 0233 disp('Block = 21'), pause(pstate) 0234 0235 %% Crest length 0236 clf 0237 disp('NIT=5 may take time, running with NIT=2 in the following') 0238 %f_Lc = spec2tpdf2(S1,[],'Lc',[0 200 81],[],opt1); 0239 f_Lc = spec2tpdf(S1,[],'Lc',[0 200 81],[],2); 0240 % f_Lc = spec2tpdf(S1,[],'Lc',[0 200 81],[],5); 0241 pdfplot(f_Lc,'-.') 0242 wafostamp([],'(ER)') 0243 disp('Block = 22'), pause(pstate) 0244 0245 % f_Lc_1 = spec2tpdf(S1,[],'Lc',[0 200 81],1.5,5); 0246 f_Lc_1 = spec2tpdf(S1,[],'Lc',[0 200 81],1.5,2); 0247 %f_Lc_1 = spec2tpdf(S1,[],'Lc',[0 200 81],1.5,opt1); 0248 0249 hold on 0250 pdfplot(f_Lc_1) 0251 wafostamp([],'(ER)') 0252 disp('Block = 23'), pause(pstate) 0253 %% 0254 clf 0255 simpson(f_Lc.x{1},f_Lc.f) 0256 simpson(f_Lc_1.x{1},f_Lc_1.f) 0257 0258 disp('Block = 24'), pause(pstate) 0259 %% 0260 clf 0261 % f_Lc_d1 = spec2tpdf(spec2spec(SD1,'rotdir',pi/2),[],... 0262 % 'Lc',[0 300 121],[],5); 0263 f_Lc_d1 = spec2tpdf(rotspec(SD1,pi/2),[],... 0264 'Lc',[0 300 121],[],2); 0265 % f_Lc_d1 = spec2tpdf2(spec2spec(SD1,'rotdir',pi/2),[],... 0266 % 'Lc',[0 300 121],[],opt1); 0267 pdfplot(f_Lc_d1,'-.') 0268 hold on 0269 % f_Lc_d12 = spec2tpdf(SD12,[],'Lc',[0 200 81],[],5); 0270 f_Lc_d12 = spec2tpdf(SD12,[],'Lc',[0 200 81],[],2); 0271 % f_Lc_d12 = spec2tpdf2(SD12,[],'Lc',[0 200 81],[],opt1); 0272 pdfplot(f_Lc_d12) 0273 hold off 0274 wafostamp([],'(ER)') 0275 disp('Block = 25'), pause(pstate) 0276 0277 %% 0278 disp('Run the following example only if you want a check on computing time') 0279 disp('Edit the command file and remove %') 0280 % clf 0281 % f_Lc_d1_5 = spec2tpdf(spec2spec(SD1,'rotdir',pi/2),[],... 0282 % 'Lc',[0 300 121],[],5); 0283 % f_Lc_d1_3 = spec2tpdf(spec2spec(SD1,'rotdir',pi/2),[],... 0284 % 'Lc',[0 300 121],[],3); 0285 % f_Lc_d1_2 = spec2tpdf(spec2spec(SD1,'rotdir',pi/2),[],... 0286 % 'Lc',[0 300 121],[],2); 0287 % f_Lc_d1_0 = spec2tpdf(spec2spec(SD1,'rotdir',pi/2),[],... 0288 % 'Lc',[0 300 121],[],0); 0289 % f_Lc_d1_n4 = spec2tpdf2(spec2spec(SD1,'rotdir',pi/2),[],... 0290 % 'Lc',[0 400 161],-4); 0291 % pdfplot(f_Lc_d1_5) 0292 % hold on 0293 % pdfplot(f_Lc_d1_2) 0294 % pdfplot(f_Lc_d1_0) 0295 % pdfplot(f_Lc_d1_n4) 0296 % simpson(f_Lc_d1_n4.x{1},f_Lc_d1_n4.f) 0297 disp('Block = 26'), pause(pstate) 0298 0299 %% Section 3.4.2 Density of wave period, wave length or encountered wave period 0300 %% Example 7: Crest period and high crest waves 0301 clf 0302 xx = load('sea.dat'); 0303 x = xx; 0304 x(:,2) = detrend(x(:,2)); 0305 SS = dat2spec2(x); 0306 si = sqrt(spec2mom(SS,1)); 0307 SS.tr = dat2tr(x); 0308 Hs = 4*si 0309 method = 0; 0310 rate = 2; 0311 [S, H, Ac, At, Tcf, Tcb, z_ind, yn] = dat2steep(x,rate,method); 0312 t = linspace(0.01,8,200); 0313 ftc1 = kde(Tc,{'L2',0},t); 0314 pdfplot(ftc1) 0315 hold on 0316 % f_t = spec2tpdf(SS,[],'Tc',[0 8 81],0,4); 0317 f_t = spec2tpdf(SS,[],'Tc',[0 8 81],0,2); 0318 simpson(f_t.x{1},f_t.f) 0319 pdfplot(f_t,'-.') 0320 hold off 0321 wafostamp([],'(ER)') 0322 disp('Block = 27'), pause(pstate) 0323 0324 %% 0325 clf 0326 % f_t2 = spec2tpdf(SS,[],'Tc',[0 8 81],[Hs/2],4); 0327 f_t2 = spec2tpdf(SS,[],'Tc',[0 8 81],[Hs/2],2); 0328 Pemp = sum(Ac>Hs/2)/sum(Ac>0) 0329 simpson(f_t2.x{1},f_t2.f) 0330 index = find(Ac>Hs/2); 0331 ftc1 = kde(Tc(index),{'L2',0},t); 0332 ftc1.f = Pemp*ftc1.f; 0333 pdfplot(ftc1) 0334 hold on 0335 pdfplot(f_t2,'-.') 0336 hold off 0337 wafostamp([],'(ER)') 0338 disp('Block = 28'), pause(pstate) 0339 0340 0341 % clf 0342 % f_tcc2 = spec2tccpdf(SS,[],'t>',[0 12 61],[Hs/2],[0],-1); 0343 % simpson(f_tcc2.x{1},f_tcc2.f) 0344 % f_tcc3 = spec2tccpdf(SS,[],'t>',[0 12 61],[Hs/2],[0],3,5); 0345 % f_tcc3 = spec2tccpdf(SS,[],'t>',[0 12 61],[Hs/2],[0],1,5); 0346 % simpson(f_tcc3.x{1},f_tcc3.f) 0347 % pdfplot(f_tcc2,'-.') 0348 % hold on 0349 % pdfplot(f_tcc3) 0350 % hold off 0351 disp('Block = 29'), pause(pstate) 0352 0353 %% 0354 clf 0355 [TC tc_ind v_ind] = dat2tc(yn,[],'dw'); 0356 N = length(tc_ind); 0357 t_ind = tc_ind(1:2:N); 0358 c_ind = tc_ind(2:2:N); 0359 Pemp = sum(yn(t_ind,2)<-Hs/2 & yn(c_ind,2)>Hs/2)/length(t_ind) 0360 ind = find(yn(t_ind,2)<-Hs/2 & yn(c_ind,2)>Hs/2); 0361 spwaveplot(yn,ind(2:4)) 0362 wafostamp([],'(ER)') 0363 disp('Block = 30'), pause(pstate) 0364 0365 %% 0366 clf 0367 Tcc = yn(v_ind(1+2*ind),1)-yn(v_ind(1+2*(ind-1)),1); 0368 t = linspace(0.01,14,200); 0369 L2 = 0; 0370 ftcc1 = kde(Tcc,{'kernel' 'epan','L2',L2},t); 0371 ftcc1.f = Pemp*ftcc1.f; 0372 pdfplot(ftcc1,'-.') 0373 wafostamp([],'(ER)') 0374 disp('Block = 31'), pause(pstate) 0375 0376 disp('The rest of this chapter deals with joint densities.') 0377 disp('Some calculations may take some time.') 0378 disp('You could experiment with other NIT.') 0379 %return 0380 0381 %% Example 8: Wave period for high crest waves 0382 clf 0383 f_tcc22_1 = spec2tccpdf(SS,[],'t>',[0 12 61],[Hs/2],[Hs/2],-1); 0384 simpson(f_tcc22_1.x{1},f_tcc22_1.f) 0385 hold on 0386 pdfplot(f_tcc22_1) 0387 hold off 0388 wafostamp([],'(ER)') 0389 disp('Block = 32'), pause(pstate) 0390 0391 %% Section 3.4.3 Joint density of crest period and crest height 0392 %% Example 9. Some preliminary analysis of the data 0393 clf 0394 yy = load('gfaksr89.dat'); 0395 SS = dat2spec(yy); 0396 si = sqrt(spec2mom(SS,1)); 0397 SS.tr = dat2tr(yy); 0398 Hs = 4*si 0399 v = gaus2dat([0 0],SS.tr); 0400 v = v(2) 0401 disp('Block = 33'), pause(pstate) 0402 0403 %% 0404 clf 0405 [TC, tc_ind, v_ind] = dat2tc(yy,v,'dw'); 0406 N = length(tc_ind); 0407 t_ind = tc_ind(1:2:N); 0408 c_ind = tc_ind(2:2:N); 0409 v_ind_d = v_ind(1:2:N+1); 0410 v_ind_u = v_ind(2:2:N+1); 0411 T_d = ecross(yy(:,1),yy(:,2),v_ind_d,v); 0412 T_u = ecross(yy(:,1),yy(:,2),v_ind_u,v); 0413 0414 % Old call 0415 %T_d = yy(v_ind_d,1)- yy(v_ind_d,2)* ... 0416 % (yy(2,1)-yy(1,1))./(yy(v_ind_d+1,2)-yy(v_ind_d,2)); 0417 %T_u = yy(v_ind_u,1)- yy(v_ind_u,2)* ... 0418 % (yy(2,1)-yy(1,1))./(yy(v_ind_u+1,2)-yy(v_ind_u,2)); 0419 Tc = T_d(2:end)-T_u(1:end); 0420 Tt = T_u(1:end)-T_d(1:end-1); 0421 Tcf = yy(c_ind,1)-T_u; 0422 Ac = yy(c_ind,2)-v; 0423 At = v-yy(t_ind,2); 0424 disp('Block = 34'), pause(pstate) 0425 0426 %% 0427 clf 0428 t = linspace(0.01,15,200); 0429 kopt3 = kdeoptset('hs',0.25,'L2',0); 0430 ftc1 = kde(Tc,kopt3,t); 0431 ftt1 = kde(Tt,kopt3,t); 0432 pdfplot(ftt1,'k') 0433 hold on 0434 pdfplot(ftc1,'k-.') 0435 f_tc4 = spec2tpdf(SS,[],'Tc',[0 12 81],0,4,5); 0436 f_tc2 = spec2tpdf(SS,[],'Tc',[0 12 81],0,2,5); 0437 f_tc = spec2tpdf(SS,[],'Tc',[0 12 81],0,-1); 0438 pdfplot(f_tc,'b') 0439 hold off 0440 wafostamp([],'(ER)') 0441 disp('Block = 35'), pause(pstate) 0442 0443 %% Example 10: Joint characteristics of a half wave: 0444 %% position and height of a crest for a wave with given period 0445 clf 0446 ind = find(4.4<Tc & Tc<4.6); 0447 f_AcTcf = kde([Tcf(ind) Ac(ind)],{'L2',[1 .5]}); 0448 plot(Tcf(ind), Ac(ind),'.'); 0449 hold on 0450 pdfplot(f_AcTcf) 0451 wafostamp([],'(ER)') 0452 disp('Block = 36'), pause(pstate) 0453 0454 %% 0455 clf 0456 %opt1 = rindoptset('speed',5,'method',3); 0457 %opt2 = rindoptset('speed',5,'nit',2,'method',0); 0458 opt1 = rindoptset('speed',9,'method',3); 0459 opt2 = rindoptset('speed',7,'nit',2,'method',0); 0460 0461 0462 f_tcfac1 = spec2thpdf(SS,[],'TcfAc',[4.5 4.5 46],[0:0.25:8],opt1); 0463 f_tcfac2=spec2thpdf(SS,[],'TcfAc',[4.5 4.5 46],[0:0.25:8],opt2); 0464 0465 plot(Tcf(ind), Ac(ind),'.'); 0466 hold on 0467 pdfplot(f_tcfac1,'-.') 0468 pdfplot(f_tcfac2) 0469 0470 simpson(f_tcfac1.x{1},simpson(f_tcfac1.x{2},f_tcfac1.f,1)) 0471 simpson(f_tcfac2.x{1},simpson(f_tcfac2.x{2},f_tcfac2.f,1)) 0472 f_tcf4=spec2tpdf(SS,[],'Tc',[4.5 4.5 46],[0:0.25:8],6); 0473 f_tcf4.f(46) 0474 wafostamp([],'(ER)') 0475 disp('Block = 37'), pause(pstate) 0476 0477 0478 0479 %% 0480 clf 0481 f_tcac_s = spec2thpdf(SS,[],'TcAc',[0 12 81],[Hs/2:0.1:2*Hs],opt1); 0482 disp('Block = 38'), pause(pstate) 0483 0484 clf 0485 mom = spec2mom(SS,4,[],0); 0486 t = f_tcac_s.x{1}; 0487 h = f_tcac_s.x{2}; 0488 flh_g = lh83pdf(t',h',[mom(1),mom(2),mom(3)],SS.tr); 0489 clf 0490 ind=find(Ac>Hs/2); 0491 plot(Tc(ind), Ac(ind),'.'); 0492 hold on 0493 pdfplot(flh_g,'k-.') 0494 pdfplot(f_tcac_s) 0495 wafostamp([],'(ER)') 0496 disp('Block = 39'), pause(pstate) 0497 0498 %% 0499 clf 0500 % f_tcac = spec2thpdf(SS,[],'TcAc',[0 12 81],[0:0.2:8],opt1); 0501 % pdfplot(f_tcac) 0502 disp('Block = 40'), pause(pstate) 0503 0504 %% Section 3.4.4 Joint density of crest and trough height 0505 %% Section 3.4.5 Min-to-max distributions Markov method 0506 %% Example 11. (min-max problems with Gullfaks data) 0507 %% Joint density of maximum and the following minimum 0508 clf 0509 tp = dat2tp(yy); 0510 Mm = fliplr(tp2mm(tp)); 0511 fmm = kde(Mm); 0512 f_mM = spec2mmtpdf(SS,[],'mm',[],[-7 7 51],opt2); 0513 clf 0514 pdfplot(f_mM,'-.') 0515 hold on 0516 pdfplot(fmm,'k-') 0517 hold off 0518 wafostamp([],'(ER)') 0519 disp('Block = 41'), pause(pstate) 0520 0521 %% The joint density of still water separated maxima and minima. 0522 clf 0523 ind = find(Mm(:,1)>v & Mm(:,2)<v); 0524 Mmv = abs(Mm(ind,:)-v); 0525 fmmv = kde(Mmv); 0526 f_vmm = spec2mmtpdf(SS,[],'vmm',[],[-7 7 51],opt2); 0527 clf 0528 pdfplot(fmmv,'k-') 0529 hold on 0530 pdfplot(f_vmm,'-.') 0531 hold off 0532 wafostamp([],'(ER)') 0533 disp('Block = 42'), pause(pstate) 0534 0535 0536 %% 0537 clf 0538 facat = kde([Ac At]); 0539 f_acat = spec2mmtpdf(SS,[],'AcAt',[],[-7 7 51],opt2); 0540 clf 0541 pdfplot(f_acat,'-.') 0542 hold on 0543 pdfplot(facat,'k-') 0544 hold off 0545 wafostamp([],'(ER)') 0546 disp('Block = 43'), pause(pstate) 0547 0548
Comments or corrections to the WAFO group