How do I calculate the static Kalman gain (= steady-state Kalman gain)?

151 次查看（过去 30 天）
Joonas
Joonas 2015-11-3
回答： Bill Tubbs 2022-9-18
Dear Matlab users,
how do I calculate the static Kalman gain in advance? Is there any functions that do this for me if I provide the system matrices and covariance data Q and R of the process and measurement noises respectively? Please inform me ways to do this.
Thank you!

请先登录，再回答此问题。

回答（4 个）

Nitin Khola
Nitin Khola 2015-11-5
编辑：Nitin Khola 2015-11-5
A popular FileExchange submission "Kalman Filter Package" might have some functionality coded for you. The following is the link to the submission: http://www.mathworks.com/matlabcentral/fileexchange/38302-kalman-filter-package
I encourage you to use the trial version of our product, "Control System Toolbox" which offers a wide array of functionalities related to analyzing, designing, and tuning linear control systems. Below is the link to the homepage of this product: http://www.mathworks.com/products/control/
Roger Labbe
Roger Labbe 2015-11-9
There are two accepted ways to do this.
First, write a numerical simulation and just see what value it converges to.
Second, any of several textbooks will give you the procedure to compute it directly. I am not going to try to type several pages of instruction here, not the least because I am not solid on the theory itself. Crassidis' "Optimal Estimation of Dynamic Systems" gives the clearest presentation that I know of in section 5.3.4. Dan Simon's book talks about computing this for different filters (Crassidis is for a linear, discrete filter), and sites the research that discusses how to determine if your filter will converge or not.
Me, since I'm not sending rockets into space or doing other life-critical work, I just do the numerical simulation!
Bill Tubbs
Bill Tubbs 2022-9-18
编辑：Bill Tubbs 2022-9-18
For the record, I found a way to do it using the idare function (see this material online: hw5sol.pdf from S. Boyd).
This is for the "prediction form" of the Kalman Filter where the correction gain is used as follows:
The steady-state correction gain is found by solving this equation for the steady-state error covariance:
% Example system matrices
A = [0.7 1;
0 1];
B = [1; 0];
C = [0.3 0];
D = 0;
% Error covariance matrices
Q = diag([0.0001 0.01]);
R = 0.01;
% Solve discrete-time algebraic Riccati equation
[P,K] = idare(A',C',Q,R,[],[]);
K = K';
% Check results satisfy eqn.s
assert(all(A*(P - P*C'*(C*P*C' + R)^-1*C*P)*A' + Q - P < 1e-14, [1 2]))
assert(all(K - A*P*C'*(C*P*C' + R)^-1 < 1e-14, [1 2]))
K
K =
1.7674
0.7514
You can also check this by using the built-in Kalman filter object:
% Construct model
n = 2; ny = 1; Ts = 1;
N = zeros(n, ny);
G = eye(n); % apply process noises to all states
H = zeros(ny, n); % no direct transmission of noises
Gmodel = ss(A, [B G], C, [D H], Ts);
% Use MATLAB's Kalman filter function
[~, K, P] = kalman(Gmodel, Q, R, N, 'delayed');
K
K =
1.7674
0.7514
If anyone knows a better way to do this I would be interested.
Bill Tubbs
Bill Tubbs 2022-9-18
In addition to my answer above, I just discovered a special function for Kalman estimator design called dlqe in the Control System Toolbox:
G = eye(2);
[Kf,P,~] = dlqe(A,G,C,Q,R);
Kf
Kf =
1.4515
0.7514
I can't find the documentation page for this function online but the help text says:
[M,P,Z,E] = dlqe(A,G,C,Q,R) returns the gain matrix M such
that the discrete, stationary Kalman filter with observation
and time update equations
x[n|n] = x[n|n-1] + M(y[n] - Cx[n|n-1] - Du[n])
x[n+1|n] = Ax[n|n] + Bu[n]
produces an optimal state estimate x[n|n] of x[n] given y[n] and
the past measurements.
I think Kf here is the gain for the "filtering" version of the Kalman filter:
So, then you can calculate the prediction form gain simply as follows:
K = A * Kf
K =
1.7674
0.7514
I think this function may be obsolete though:
>> which dlqe
/Applications/MATLAB_R2021b.app/toolbox/control/ctrlobsolete/dlqe.m

类别

在 Help Center 和 File Exchange 中查找有关 State-Space Control Design and Estimation 的更多信息

标签

产品

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by