/
correlation_3rd_order_comparison_v2.m
124 lines (115 loc) · 3.4 KB
/
correlation_3rd_order_comparison_v2.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
close all
clear all
clc
load FTLE_wrf_short
load error_comparison
s1 = s1(:,:,1);
l1 = l1(:,:,1);
a1 = a1(:,:,1);
a2 = a2(:,:,1);
T
i=0;
s1 = s1((i+1):end-i,(i+1):end-i,1);
l1 = l1((i+1):end-i,(i+1):end-i,1);
a1 = a1((i+1):end-i,(i+1):end-i,1);
a2 = a2((i+1):end-i,(i+1):end-i,1);
sigma = sigma((i+1):end-i,(i+1):end-i,:);
sigma(:,:,1)=-s1;
n_t = length(T);
for i =1:n_t
i
ftle_t = squeeze(sigma(:,:,i));
sig_true = reshape(ftle_t,[],1);
sig_approx = reshape(-s1-T(i)*(-s1.^2+0.5*l1+T(i)*(4/3*s1.^3-s1.*l1+0.25*a1)),[],1);
sig_approx2 = reshape(-s1-T(i)*(-s1.^2+0.5*l1+T(i)*(4/3*s1.^3-s1.*l1+0.25*a2)),[],1);
ind = ~isnan(sig_true) & ~isnan(sig_approx) & ~isnan(sig_approx2);
sig_true = sig_true(ind);
sig_approx = sig_approx(ind);
rmsea1(i) = sqrt(mean((sig_approx-sig_true).^2));
sa_bar = mean(sig_approx);
st_bar = mean(sig_true);
n=length(sig_true);
numerator = sum(sig_approx.*sig_true)-(n*sa_bar*st_bar);
den1 = sqrt(sum(sig_approx.^2)-n*sa_bar.^2);
den2 = sqrt(sum(sig_true.^2)-n*st_bar.^2);
denominator = den1*den2;
corra1(i) = numerator./denominator;
sig_approx = sig_approx2(ind);
rmsea2(i) = sqrt(mean((sig_approx-sig_true).^2));
sa_bar = mean(sig_approx);
st_bar = mean(sig_true);
n=length(sig_true);
numerator = sum(sig_approx.*sig_true)-(n*sa_bar*st_bar);
den1 = sqrt(sum(sig_approx.^2)-n*sa_bar.^2);
den2 = sqrt(sum(sig_true.^2)-n*st_bar.^2);
denominator = den1*den2;
corra2(i) = numerator./denominator;
sig_approx = reshape(-s1-T(i)*(-s1.^2+0.5*l1),[],1);
sig_approx = sig_approx(ind);
rmse2(i) = sqrt(mean((sig_approx-sig_true).^2));
sa_bar = mean(sig_approx);
st_bar = mean(sig_true);
numerator = sum(sig_approx.*sig_true)-(n*sa_bar*st_bar);
den1 = sqrt(sum(sig_approx.^2)-n*sa_bar.^2);
den2 = sqrt(sum(sig_true.^2)-n*st_bar.^2);
denominator = den1*den2;
corr2(i) = numerator./denominator;
sig_approx = reshape(-s1,[],1);
sig_approx = sig_approx(ind);
rmse1(i) = sqrt(mean((sig_approx-sig_true).^2));
sa_bar = mean(sig_approx);
st_bar = mean(sig_true);
numerator = sum(sig_approx.*sig_true)-(n*sa_bar*st_bar);
den1 = sqrt(sum(sig_approx.^2)-n*sa_bar.^2);
den2 = sqrt(sum(sig_true.^2)-n*st_bar.^2);
denominator = den1*den2;
corr1(i) = numerator./denominator;
end
T=T/3600;
figure
subplot(121)
hold on
plot(abs(T),corr1,'r.-')
plot(abs(T),corr2,'b.-')
plot(abs(T),corra1,'k.-')
legend('-s1','-s1-O(T)','-s1-O(T^2)','Location','northwest')
ylabel('RMSE s^{-1}')
xlabel('|T| s')
title('lambda2 = X0^T*Q*X0 + X0^T*B*X1 - X0^T*S*X1')
subplot(122)
hold on
plot(abs(T),corr1,'r.-')
plot(abs(T),corr2,'b.-')
plot(abs(T),corra2,'k.-')
legend('-s1','-s1-O(T)','-s1-O(T^2)','Location','northwest')
ylabel('RMSE s^{-1}')
xlabel('|T| s')
title('lambda2 = X0^T*Q*X0-d.^2/m')
%
figure
subplot(121)
hold on
plot(abs(T),corr1,'r.-')
plot(abs(T),corr2,'b.-')
plot(abs(T),corra1,'k.-')
legend('-s1-T*corr','-s1','Location','southeast')
ylabel('RMSE s^{-1}')
xlabel('|T| s')
set(gca, 'YScale', 'log')
set(gca, 'XScale', 'log')
axis('equal')
subplot(122)
hold on
plot(abs(T),corr1,'r.-')
plot(abs(T),corr2,'b.-')
plot(abs(T),corra2,'k.-')
legend('-s1-T*cor','-s1','Location','southeast')
ylabel('RMSE s^{-1}')
xlabel('|T| s')
set(gca, 'YScale', 'log')
set(gca, 'XScale', 'log')
axis('equal')
%}
corr3 = corra1;
time = T;
save wrf_plot_corr_data corr1 corr2 corr3 time