Skip to content

Commit faaf1f7

Browse files
author
Abhinab Bhattacharjee
committed
Minor comments and refactoring
1 parent 826d6dc commit faaf1f7

24 files changed

+181
-67
lines changed

‎.gitignore

+1
Original file line numberDiff line numberDiff line change
@@ -4,3 +4,4 @@
44
*.fig
55
*.pdf
66
*.mat
7+
*.html

‎docs/CONTRIBUTING.md

+17
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,17 @@
1+
# Contributing to DATools
2+
3+
General Guidelines for adding new filters/algorithms to DATools
4+
5+
## Submitting new pull requests
6+
7+
To request a pull request, please follow the steps below
8+
9+
* Create your own branch if you are starting new or modify the existing one
10+
* Commit your changes to local branch and push to Github
11+
* Create a draft pull request
12+
* Do not merge or modify the master branch
13+
14+
## Style Conventions
15+
16+
In order to prevent confusions and maintain a consistent coding style, the following styles should be strictly followed. The conventions match most of the MATLAB's own conventions. Failing to follow the conventions can delay in processing and approving new pull requests.
17+

‎src/+datools/+error/Gaussian.m

+1-1
Original file line numberDiff line numberDiff line change
@@ -22,7 +22,7 @@
2222
end
2323

2424
function xp = adderr(obj, ~, x)
25-
A = obj.CovarianceSqrt * randn(size(x, 1), size(x, 2));
25+
A = obj.CovarianceSqrt * randn(size(x), 'like', x);
2626
xp = x + A + obj.Bias;
2727
end
2828

‎src/+datools/+error/Laplace.m

+1-1
Original file line numberDiff line numberDiff line change
@@ -20,7 +20,7 @@
2020
end
2121

2222
function xp = adderr(obj, ~, x)
23-
X = obj.CovarianceSqrt * randn(size(x, 1), size(x, 2));
23+
X = obj.CovarianceSqrt * randn(size(x), 'like', x);
2424
Z = exprnd(1, 1, size(x, 2));
2525
Y = sqrt(Z) .* X;
2626
xp = x + Y + obj.Bias;

‎src/+datools/+observation/+operator/Indexed.m renamed to ‎src/+datools/+observation/+operator/IndexedObservation.m

+2-2
Original file line numberDiff line numberDiff line change
@@ -1,4 +1,4 @@
1-
classdef Indexed < datools.observation.operator.Observation
1+
classdef IndexedObservation < datools.observation.operator.Observation
22
%INDEXED Defines the indexed observation operator
33

44
properties
@@ -7,7 +7,7 @@
77

88
methods
99

10-
function obj = Indexed(nvars, varargin)
10+
function obj = IndexedObservation(nvars, varargin)
1111
p = inputParser;
1212
p.KeepUnmatched = true;
1313
addParameter(p, 'Indices', 1);

‎src/+datools/+observation/+operator/Linear.m renamed to ‎src/+datools/+observation/+operator/LinearObservation.m

+2-2
Original file line numberDiff line numberDiff line change
@@ -1,4 +1,4 @@
1-
classdef Linear < datools.observation.operator.Observation
1+
classdef LinearObservation < datools.observation.operator.Observation
22
%LINEAR Defines the linear observation operator
33
% H is already linearized
44

@@ -8,7 +8,7 @@
88

99
methods
1010

11-
function obj = Linear(nvars, varargin)
11+
function obj = LinearObservation(nvars, varargin)
1212
p = inputParser;
1313
p.KeepUnmatched = true;
1414
addParameter(p, 'H', speye(nvars));

‎src/+datools/+observation/+operator/Nonlinear.m renamed to ‎src/+datools/+observation/+operator/NonlinearObservation.m

+2-2
Original file line numberDiff line numberDiff line change
@@ -1,4 +1,4 @@
1-
classdef Nonlinear < datools.observation.operator.Observation
1+
classdef NonlinearObservation < datools.observation.operator.Observation
22
%NONLINEAR Defines the non linear observation operator
33

44
properties
@@ -9,7 +9,7 @@
99

1010
methods
1111

12-
function obj = Nonlinear(nvars, varargin)
12+
function obj = NonlinearObservation(nvars, varargin)
1313
p = inputParser;
1414
p.KeepUnmatched = true;
1515
addParameter(p, 'F', @(~, x) x);

‎src/+datools/+observation/StateObservation.m

+6-6
Original file line numberDiff line numberDiff line change
@@ -3,10 +3,10 @@
33
% Observed state are noisy snapshots of reality (truth)
44

55
properties
6-
Y % Current Observation of States (maybe sparse)
7-
R % Observation Error Covariance
8-
Noise % Noise
9-
NumObs % Number of Observations(Make sure > 0)
6+
Y % Current Observation of States (maybe sparse)
7+
Covariance % Observation Error Covariance
8+
Noise % Noise
9+
NumObs % Number of Observations(Make sure > 0)
1010
end
1111

1212
methods
@@ -30,11 +30,11 @@
3030
unMatched = p.Unmatched;
3131

3232
p = inputParser;
33-
addParameter(p, 'R', speye(obj.NumObs));
33+
addParameter(p, 'Covariance', speye(obj.NumObs));
3434

3535
parse(p, unMatched);
3636
s = p.Results;
37-
obj.R = s.R;
37+
obj.Covariance = s.Covariance;
3838
end
3939

4040
function updateY(obj, y)

‎src/+datools/+statistical/+ensemble/DEnKF.m

+4-2
Original file line numberDiff line numberDiff line change
@@ -2,14 +2,16 @@
22

33
methods
44

5-
function analysis(obj, R, y)
5+
function analysis(obj)
66

77
inflation = obj.Inflation;
88

99
tc = obj.Model.TimeSpan(1);
1010

1111
xf = obj.Ensemble;
1212
ensN = obj.NumEnsemble;
13+
14+
R = obj.Observation.Covariance;
1315

1416
xfm = mean(xf, 2);
1517

@@ -43,7 +45,7 @@ function analysis(obj, R, y)
4345

4446
Aa = Af - 0.5 * (PfHt * (dS \ HAf));
4547

46-
d = y - Hxfm;
48+
d = obj.Observation.Y - Hxfm;
4749
xam = xfm + (PfHt * (dS \ d));
4850

4951
obj.Ensemble = repmat(xam, 1, ensN) + Aa;

‎src/+datools/+statistical/+ensemble/ETKF.m

+7-2
Original file line numberDiff line numberDiff line change
@@ -3,15 +3,20 @@
33
methods
44

55
function analysis(obj)
6-
6+
%ANALYSIS Method to overload the analysis function
7+
%
8+
% ANALYSIS(OBJ) assimilates the current observation with the
9+
% background/prior information to get a better estimate
10+
% (analysis/posterior)
11+
712
inflation = obj.Inflation;
813

914
tc = obj.Model.TimeSpan(1);
1015

1116
xf = obj.Ensemble;
1217
ensN = obj.NumEnsemble;
1318

14-
R = obj.Observation.R;
19+
R = obj.Observation.Covariance;
1520

1621
xfm = mean(xf, 2);
1722
Af = xf - repmat(xfm, 1, ensN);

‎src/+datools/+statistical/+ensemble/ETPF.m

+11-4
Original file line numberDiff line numberDiff line change
@@ -2,17 +2,24 @@
22

33
methods
44

5-
function analysis(obj, R, y)
6-
5+
function analysis(obj)
6+
%ANALYSIS Method to overload the analysis function
7+
%
8+
% ANALYSIS(OBJ) assimilates the current observation with the
9+
% background/prior information to get a better estimate
10+
% (analysis/posterior)
11+
712
tau = obj.Rejuvenation;
813

914
tc = obj.Model.TimeSpan(1);
15+
16+
R = obj.Observation.Covariance;
1017

1118
dR = decomposition(R, 'chol');
1219

1320
xf = obj.Ensemble;
1421
ensN = obj.NumEnsemble;
15-
22+
1623
AeqT = [kron(speye(ensN), ones(1, ensN)); kron(ones(ensN, 1), speye(ensN)).'];
1724
lbT = zeros((ensN)*ensN, 1);
1825
optsT = optimoptions('linprog', 'Display', 'off');
@@ -26,7 +33,7 @@ function analysis(obj, R, y)
2633
xdist(i, :) = vecnorm(xtemp).^2;
2734
end
2835

29-
t0 = Hxf - y;
36+
t0 = Hxf - obj.Observation.Y;
3037

3138
% more efficient way of calculating weights
3239
as = (-0.5 * sum(t0.*(dR \ t0), 1)).';

‎src/+datools/+statistical/+ensemble/ETPF2.m

+9-3
Original file line numberDiff line numberDiff line change
@@ -2,14 +2,20 @@
22

33
methods
44

5-
function analysis(obj, R, y)
6-
5+
function analysis(obj)
6+
%ANALYSIS Method to overload the analysis function
7+
%
8+
% ANALYSIS(OBJ) assimilates the current observation with the
9+
% background/prior information to get a better estimate
10+
% (analysis/posterior).
711
tau = obj.Rejuvenation;
812

913
tc = obj.Model.TimeSpan(1);
1014

1115
xf = obj.Ensemble;
1216
ensN = obj.NumEnsemble;
17+
18+
R = obj.Observation.Covariance;
1319

1420
dR = decomposition(R, 'chol');
1521

@@ -26,7 +32,7 @@ function analysis(obj, R, y)
2632
xdist(i, :) = vecnorm(xtemp).^2;
2733
end
2834

29-
t0 = Hxf - y;
35+
t0 = Hxf - obj.Observation.Y;
3036

3137
% more efficient way of calculating weights
3238
as = (-0.5 * sum(t0.*(dR \ t0), 1)).';

‎src/+datools/+statistical/+ensemble/EnGMF.m

+12-5
Original file line numberDiff line numberDiff line change
@@ -2,14 +2,21 @@
22

33
methods
44

5-
function analysis(obj, R, y)
6-
5+
function analysis(obj)
6+
%ANALYSIS Method to overload the analysis function
7+
%
8+
% ANALYSIS(OBJ) assimilates the current observation with the
9+
% background/prior information to get a better estimate
10+
% (analysis/posterior)
11+
712
tc = obj.Model.TimeSpan(1);
813

914
xf = obj.Ensemble;
1015
ensN = obj.NumEnsemble;
1116
n = size(xf, 1);
12-
17+
18+
R = obj.Observation.Covariance;
19+
1320
beta = ((4/(n + 2))^(2/(n + 4)))*((ensN)^(-2/(n + 4)));
1421

1522
xfm = mean(xf, 2);
@@ -23,11 +30,11 @@ function analysis(obj, R, y)
2330

2431
dS = decomposition((HAf * HAf.') + R, 'chol');
2532

26-
xtilde = xf - (Af*HAf.')*(dS\(Hxf - y));
33+
xtilde = xf - (Af*HAf.')*(dS\(Hxf - obj.Observation.Y));
2734
Aa = xtilde - mean(xtilde, 2);
2835
Aa = sqrt(beta)*Aa/sqrt(ensN - 1);
2936

30-
t0 = (Hxf - y);
37+
t0 = (Hxf - obj.Observation.Y);
3138

3239
as = (-0.5 * sum(t0.*(dS \ t0), 1)).';
3340
m = max(as);

‎src/+datools/+statistical/+ensemble/EnKF.m

+10-3
Original file line numberDiff line numberDiff line change
@@ -4,13 +4,20 @@
44
methods
55

66
function analysis(obj)
7-
7+
%ANALYSIS Method to overload the analysis function
8+
%
9+
% ANALYSIS(OBJ) assimilates the current observation with the
10+
% background/prior information to get a better estimate
11+
% (analysis/posterior)
12+
813
inflation = obj.Inflation;
914

1015
tc = obj.Model.TimeSpan(1);
1116

1217
xf = obj.Ensemble;
1318
ensN = obj.NumEnsemble;
19+
20+
R = obj.Observation.Covariance;
1421

1522
xfm = mean(xf, 2);
1623

@@ -38,12 +45,12 @@ function analysis(obj)
3845
PfHt = rhoHt .* ((1 / (ensN - 1)) * (Af * (HAf.')));
3946
HPfHt = HrhoHt .* ((1 / (ensN - 1)) * (HAf * (HAf.')));
4047

41-
S = HPfHt + obj.Observation.R;
48+
S = HPfHt + R;
4249
dS = decomposition(S, 'chol');
4350
d = obj.Observation.Y - Hxfm;
4451

4552
xam = xfm + PfHt * (dS \ d);
46-
Aa = Af + PfHt * (dS \ (sqrtm(obj.Observation.R) *...
53+
Aa = Af + PfHt * (dS \ (sqrtm(R) *...
4754
randn(size(HAf)) - HAf));
4855

4956
obj.Ensemble = repmat(xam, 1, ensN) + Aa;

‎src/+datools/+statistical/+ensemble/FETKF.m

+10-5
Original file line numberDiff line numberDiff line change
@@ -35,16 +35,21 @@
3535

3636
end
3737

38-
function analysis(obj, R, y)
39-
40-
38+
function analysis(obj)
39+
%ANALYSIS Method to overload the analysis function
40+
%
41+
% ANALYSIS(OBJ) assimilates the current observation with the
42+
% background/prior information to get a better estimate
43+
% (analysis/posterior)
44+
4145
inflation = obj.Inflation;
4246

4347
tc = obj.Model.TimeSpan(1);
4448

4549
xf = obj.Ensemble;
4650
ensN = obj.NumEnsemble;
47-
51+
52+
R = obj.Observation.Covariance;
4853

4954
xfm = mean(xf, 2);
5055
n = numel(xfm);
@@ -101,7 +106,7 @@ function analysis(obj, R, y)
101106

102107
TT = (eye(ensN+ensNW) - ZZ.' / S * ZZ);
103108

104-
d = y - Hxfm;
109+
d = obj.Observation.Y - Hxfm;
105110

106111
AAa = [Afs, ws] * real(sqrtm(TT));
107112
Aa = sqrt(ensN-1) * AAa(:, 1:ensN) / sqrt(1-gamma);

‎src/+datools/+statistical/+ensemble/FETPF.m

+10-3
Original file line numberDiff line numberDiff line change
@@ -42,8 +42,13 @@
4242

4343
methods
4444

45-
function analysis(obj, R, y)
46-
45+
function analysis(obj)
46+
%ANALYSIS Method to overload the analysis function
47+
%
48+
% ANALYSIS(OBJ) assimilates the current observation with the
49+
% background/prior information to get a better estimate
50+
% (analysis/posterior)
51+
4752
inflation = obj.Inflation;
4853
tau = obj.Rejuvenation;
4954

@@ -53,6 +58,8 @@ function analysis(obj, R, y)
5358

5459
xf = obj.Ensemble;
5560
ensN = obj.NumEnsemble;
61+
62+
R = obj.Observation.Covariance;
5663

5764
dR = decomposition(R, 'chol');
5865

@@ -138,7 +145,7 @@ function analysis(obj, R, y)
138145
xdist(i, :) = vecnorm(X_temp, 2).^2;
139146
end
140147

141-
t0 = Hchif - y;
148+
t0 = Hchif - obj.Observation.Y;
142149

143150
% more efficient way of calculating weights
144151
as = (-0.5 * sum(t0.*(dR \ t0), 1)).';

0 commit comments

Comments
 (0)