Lekshmi B S
Member level 1
i want to do channel estimation using kalman filter. this is the matlab code i formulated using the equations. the noise added is gaussian. but the problem is that i am getting a constant error instead of error getting reduced with the number of iterations. can anyone please help me to find the correction needed in this...
its urgent..if anyone can find the mistakes, help...
Thank you..
N=20;
F=0.4;
H=0.07;
Q=0.4;
R=0.8;
x=zeros(1,N);
z=zeros(1,N);
residual=zeros(1,N);
K=zeros(1,N);
apriori_state=zeros(1,N);
apriori_error=ones(1,N);
aposteriori_state=zeros(1,N);
aposteriori_error=ones(1,N);
x_0=1;
aposteriori_state_0=1.5;
aposteriori_error_0=1;
w(1)=(sqrt(1/(((2*pi)^N)*det(Q))))*(exp((-x_0'*inv(Q)*x_0)/2));
v(1)=(sqrt(1/(((2*pi)^N)*det(R))))*(exp((-x_0'*inv(R)*x_0)/2));
x(1)=F*x_0+w(1);
z(1)=H*x(1)+v(1);
apriori_state(1)=F*aposteriori_state_0;
apriori_error(1)=F*F'*aposteriori_error_0+Q;
residual(1)=z(1)-H*apriori_state(1);
K(1)=H*apriori_error(1)/(H*H'*apriori_error(1)+R);
aposteriori_state(1)=apriori_state(1)+K(1)*residual(1);
aposteriori_error(1)=apriori_error(1)*(1-(H*K(1)));
for i=2:N
w(i)=(sqrt(1/(((2*pi)^N)*det(Q))))*(exp((-x(i-1)'*inv(Q)*x(i-1))/2));
v(i)=(sqrt(1/(((2*pi)^N)*det(R))))*(exp((-x(i-1)'*inv(R)*x(i-1))/2));
x(i)=F*x(i-1)+w(i);
z(i)=H*x(i)+v(i);
apriori_state(i)=F*aposteriori_state(i-1);
apriori_error(i)=F*F'*aposteriori_error(i-1)+Q;
residual(i)=z(i)-H*apriori_state(i);
K(i)=H*apriori_error(i)/(H*H'*apriori_error(i)+R);
aposteriori_state(i)=apriori_state(i)+K(i)*residual(i);
aposteriori_error(i)=apriori_error(i)*(1-(H*K(i)));
end
error1=x-apriori_state;
error2=x-aposteriori_state;
figure(1);
t=1:N;
plot(t,K);
xlabel('Observation number');
ylabel('Kalman gain');
title('Kalman gain');
figure(2);
t=1:N;
subplot(2,1,1);
plot(t,error1);
xlabel('Observation number');
ylabel('Error1');
title('Error between actual value and predicted value');
subplot(2,1,2);
plot(t,error2);
xlabel('Observation number');
ylabel('Error2');
title('Error between actual value and updated value');
figure(3);
subplot(2,1,1);
plot(error1,K);
xlabel('Error1');
ylabel('Kalman gain');
title('Error1 Vs Kalman gain');
subplot(2,1,2);
plot(error2,K);
xlabel('Error2');
ylabel('Kalman gain');
title('Error2 Vs Kalman gain');
its urgent..if anyone can find the mistakes, help...
Thank you..
N=20;
F=0.4;
H=0.07;
Q=0.4;
R=0.8;
x=zeros(1,N);
z=zeros(1,N);
residual=zeros(1,N);
K=zeros(1,N);
apriori_state=zeros(1,N);
apriori_error=ones(1,N);
aposteriori_state=zeros(1,N);
aposteriori_error=ones(1,N);
x_0=1;
aposteriori_state_0=1.5;
aposteriori_error_0=1;
w(1)=(sqrt(1/(((2*pi)^N)*det(Q))))*(exp((-x_0'*inv(Q)*x_0)/2));
v(1)=(sqrt(1/(((2*pi)^N)*det(R))))*(exp((-x_0'*inv(R)*x_0)/2));
x(1)=F*x_0+w(1);
z(1)=H*x(1)+v(1);
apriori_state(1)=F*aposteriori_state_0;
apriori_error(1)=F*F'*aposteriori_error_0+Q;
residual(1)=z(1)-H*apriori_state(1);
K(1)=H*apriori_error(1)/(H*H'*apriori_error(1)+R);
aposteriori_state(1)=apriori_state(1)+K(1)*residual(1);
aposteriori_error(1)=apriori_error(1)*(1-(H*K(1)));
for i=2:N
w(i)=(sqrt(1/(((2*pi)^N)*det(Q))))*(exp((-x(i-1)'*inv(Q)*x(i-1))/2));
v(i)=(sqrt(1/(((2*pi)^N)*det(R))))*(exp((-x(i-1)'*inv(R)*x(i-1))/2));
x(i)=F*x(i-1)+w(i);
z(i)=H*x(i)+v(i);
apriori_state(i)=F*aposteriori_state(i-1);
apriori_error(i)=F*F'*aposteriori_error(i-1)+Q;
residual(i)=z(i)-H*apriori_state(i);
K(i)=H*apriori_error(i)/(H*H'*apriori_error(i)+R);
aposteriori_state(i)=apriori_state(i)+K(i)*residual(i);
aposteriori_error(i)=apriori_error(i)*(1-(H*K(i)));
end
error1=x-apriori_state;
error2=x-aposteriori_state;
figure(1);
t=1:N;
plot(t,K);
xlabel('Observation number');
ylabel('Kalman gain');
title('Kalman gain');
figure(2);
t=1:N;
subplot(2,1,1);
plot(t,error1);
xlabel('Observation number');
ylabel('Error1');
title('Error between actual value and predicted value');
subplot(2,1,2);
plot(t,error2);
xlabel('Observation number');
ylabel('Error2');
title('Error between actual value and updated value');
figure(3);
subplot(2,1,1);
plot(error1,K);
xlabel('Error1');
ylabel('Kalman gain');
title('Error1 Vs Kalman gain');
subplot(2,1,2);
plot(error2,K);
xlabel('Error2');
ylabel('Kalman gain');
title('Error2 Vs Kalman gain');