|
| 1 | +clear; |
| 2 | +clc; |
| 3 | +addpath('../utilities'); |
| 4 | + |
| 5 | +%% get measurement data |
| 6 | +calmstd = load('../data/hovsore_wind/data/offshore_detrend.mat').calmstd; |
| 7 | +x1 = calmstd.u; |
| 8 | +x2 = calmstd.u_stdv; |
| 9 | +save('res/data_measurement.mat', 'x1','x2'); |
| 10 | + |
| 11 | + |
| 12 | +%-------------------------------------------------------------------------- |
| 13 | +%% marginal distrbution fitting |
| 14 | +% marginal distribution of x1 |
| 15 | +[x_pd1,pd1,x_cd1,cd1,x_icd1,icd1] = epf(x1, 100); % histogram method |
| 16 | + |
| 17 | + % fit x1 to lognormal distribution |
| 18 | + mu_i = mean(x1); |
| 19 | + sigma_i = std(x1); |
| 20 | + ln_mu_i = log(mu_i^2/sqrt(sigma_i^2+mu_i^2)); |
| 21 | + ln_sigma_i = sqrt(log(sigma_i^2/mu_i^2+1)); |
| 22 | + pd1_ln = pdf('Lognormal',x_pd1,ln_mu_i,ln_sigma_i); |
| 23 | + % fit x1 to Weibull distribution |
| 24 | + pd_f1 = fitdist(x1, 'Weibull'); %fitted pd1 |
| 25 | + pd1_wb = pdf('Weibull',x_pd1,pd_f1.A,pd_f1.B); |
| 26 | + % [M,V] = wblstat(pd_f1.A,pd_f1.B); |
| 27 | + % fit x2 to Normal distribution |
| 28 | + pd1_normal = pdf('Normal',x_pd1, mu_i,sigma_i); |
| 29 | + |
| 30 | + |
| 31 | +% marginal distribution of x2 |
| 32 | +[x_pd2,pd2,x_cd2,cd2,x_icd2,icd2] = epf(x2, 100); % histogram method |
| 33 | + |
| 34 | + % fit x2 to lognormal distribution |
| 35 | + mu_j = mean(x2); |
| 36 | + sigma_j = std(x2); |
| 37 | + ln_mu_j = log(mu_j^2/sqrt(sigma_j^2+mu_j^2)); |
| 38 | + ln_sigma_j = sqrt(log(sigma_j^2/mu_j^2+1)); |
| 39 | + pd2_ln = pdf('Lognormal',x_pd2,ln_mu_j,ln_sigma_j); |
| 40 | + |
| 41 | + % fit x2 to Weibull distribution |
| 42 | + pd_f2 = fitdist(x2, 'Weibull'); |
| 43 | + pd2_wb = pdf('Weibull',x_pd2,pd_f2.A,pd_f2.B); |
| 44 | + % [M,V] = wblstat(pd_f2.A,pd_f2.B); |
| 45 | + |
| 46 | + % fit x2 Normal distribution |
| 47 | + pd2_normal = pdf('Normal',x_pd2, mu_j,sigma_j); |
| 48 | + |
| 49 | + save('res/dis_measurement.mat', 'x_pd1', 'pd1', 'x_cd1', 'cd1', 'x_icd1', 'icd1', ... |
| 50 | + 'x_pd2', 'pd2', 'x_cd2', 'cd2', 'x_icd2', 'icd2'); |
| 51 | + |
| 52 | +%-------------------------------------------------------------------------- |
| 53 | +%% plot marginal distributions |
| 54 | +xPmin = [0 0]; % lower bound for figure plot |
| 55 | +xPmax = [40 4]; % upper bound for figure plot |
| 56 | +figure; |
| 57 | + % plot(x_pd1,pd1,'DisplayName', 'x_1'); |
| 58 | + semilogy(x_pd1, pd1,'DisplayName', 'x_1'); |
| 59 | + hold on; |
| 60 | + plot(x_pd1,pd1_ln,'-.','DisplayName', 'Lognormal','LineWidth',1.5) |
| 61 | + plot(x_pd1,pd1_wb,':', 'DisplayName', 'Weibull','LineWidth',1.5) |
| 62 | + plot(x_pd1,pd1_normal,'k--', 'DisplayName', 'Gaussian','LineWidth',1.5); |
| 63 | + ylim([1e-6 0.12]) |
| 64 | + xlabel('{\it u} (m/s)') |
| 65 | + ylabel('Probability density') |
| 66 | + legend() |
| 67 | + xlim([xPmin(1), xPmax(1)]) |
| 68 | + % yticks([1e-6 1e-5 1e-4 1e-3 1e-2 1e-1]) |
| 69 | + % yticklabels({'10^{-6}','10^{-5}','10^{-4}','10^{-3}','10^{-2}','10^{-1}'}) |
| 70 | + % ylim([1e-7, 10]) |
| 71 | + set(gca,'FontSize',14,'FontName','Times New Roman') |
| 72 | + hold off |
| 73 | + |
| 74 | +figure; |
| 75 | + % plot(x_pd2,pd2,'DisplayName', 'x_2'); |
| 76 | + semilogy(x_pd2, pd2,'DisplayName', 'x_2'); |
| 77 | + hold on; |
| 78 | + plot(x_pd2,pd2_ln,'-.','DisplayName', 'Lognormal','LineWidth',1.5); |
| 79 | + plot(x_pd2,pd2_wb,':', 'DisplayName', 'Weibull','LineWidth',1.5) |
| 80 | + plot(x_pd2,pd2_normal,'k--', 'DisplayName', 'Gaussian','LineWidth',1.5); |
| 81 | + ylim([1e-6 1.4]) |
| 82 | + xlabel('\sigma_u (m/2)') |
| 83 | + ylabel('Probability density') |
| 84 | + legend() |
| 85 | + xlim([xPmin(2), xPmax(2)]) |
| 86 | + % yticks([1e-6 1e-5 1e-4 1e-3 1e-2 1e-1]) |
| 87 | + % yticklabels({'10^{-6}','10^{-5}','10^{-4}','10^{-3}','10^{-2}','10^{-1}'}) |
| 88 | + % ylim([1e-7, 10]) |
| 89 | + set(gca,'FontSize',14,'FontName','Times New Roman') |
| 90 | + hold off |
0 commit comments