您所在位置：主题：【求助】雷达跟踪滤波器，里面有源代码，可否帮忙看下丰富下。
vito86
[专家分：0] 发布于 2009-05-13 17:58:00
大家好，这个是根据文献《基于MATLAB的雷达跟踪滤波器设计》的源代码，是做αβγ滤波器的程序及仿真的，但是现在被定，研究内容太少了，现在才50多行，可否帮忙丰富一下，让这个程序达到100+行就行，我不知道再做点什么，或者丰富什么，只求程序能运行及仿真。谢谢大家了！
另外朋友给我了一段 信杂比的程序代码，但是我不知道怎么能把这两个程序联系在一起，请达人帮忙！不知道怎么把他那段程序改改，添加到我程序里。
原程序代码：
function abc
clear all
eps = 0.0000001;
npts = 5000;
del = 1./ 5000.;
t = 0.:del: 1.;
% 产生输入序列
inp = 1+ t.^3 + 0.5 * t.^2 + cos(2 * pi * 5 * t);
x0 = [2, 0.1, 0.01]';
% 在几秒钟内更新间隔
T = del;
smcof = .8; nvar=0.05;
[e, y] = abr_filter(x0, smcof, inp, npts, T, nvar);
figure(1)
NN = 3499.;
n = 1: NN;
subplot(2,1,1)
plot(n * T , inp(1: NN),'k')
ylabel ('目标真实位置')
axis tight
subplot(2, 1, 2)
plot(n * T, y(1: NN),'k')
xlabel ('time(s)'), ylabel ('预测目标位置')
axis tight
figure(2)
plot([1: 500]* T, e(1: 500)','k')
axis([0 inf -0.5 0.5])
xlabel ('time(s)'), ylabel ('预留误差')
grid on
function [e, y] = abr_filter(x0, smcof, inp, npts, T , nvar)
rn=1.;
% 状态向量赋初值
X = x0;
theta = smcof;
% 计算αβγ的值
alpha = 1. - (theta^3);
beta= 1.5 * (1. + theta) * ((1. - theta)^2)/T;
gamma= ((1. - theta)^3)/(T^2);
% 装入状态转移矩阵PHI
PHI = [1. T (T^2)/2.;0. 1. T;0. 0. 1.];
while rn < npts;
% 用状态转移矩阵预测下一状态
XN = PHI * X;
error = (inp(rn) + normrnd(0,nvar)) - XN(1);
e(rn) = error;
tmp1 = alpha * error;
tmp2 = beta * error;
tmp3 = gamma * error;
% 预测下一状态
X(1) = XN(1) + tmp1;
X(2) = XN(2) + tmp2;
X(3) = XN(3) + tmp3;
y(rn) = X(1);
rn = rn +1.;
end
附加 信杂比 代码，帮忙看下能否加我程序里，怎么改？
clc;
clear all;
close all;
load 'Input.mat'
N=500;
x= 1:1:N;
Y=zeros(1,N);
Yd=zeros(1,N);
Yp=zeros(1,N);
M_IN=M(row,:);
Y(1:24)=M_IN(1:24);
Yd(1:24)=1;
gama=0.6;
alpha=0:0.01:0.9;
beta=0:0.001:0.9;
SCR_input=zeros(length(beta),length(alpha));
SCR_output=zeros(length(beta),length(alpha));
SCR_improve=zeros(length(beta),length(alpha));
% SCR_max=zeros(1,length(alpha))
% for i=1:length(alpha)
% for i=83
for k=83
for j=1:length(beta)
for m=24:N
Yp(m)=Y(m-1)+alpha(k)*Yd(m-1);
Y(m)=Yp(m)+(M_IN(m)-Yp(m))*beta(j);
Yd(m)=Yd(m-1)+(M_IN(m)-Yp(m)-Yd(m-1))*gama;
end
Output=M_IN-Y;
Input=M_IN;
SCR_input(j,k)=sum(Input(col-1:col+1).^2)/mean(Input.^2);
SCR_output(j,k)=sum(Output(col-1:col+1).^2)/mean(Output.^2);
SCR_improve(j,k)=10*log(SCR_output(j,k)/SCR_input(j,k));
end
end
[c,b]=max(SCR_improve);
[e,f]=max(max(SCR_improve));
% SCR_max(i)=e;
% end
figure;
plot(beta,SCR_improve(:,83))
% mesh(beta,alpha,SCR_improve);
xlabel('beta');
% ylabel('alpha');
title('the improvement of SCR (dB)');
axis([0 1 -40 20])
[c,b]=max(SCR_improve);
[e,f]=max(max(SCR_improve));
row=b(f);
col=f;
beta(row)
alpha(col)
e
另外朋友给我了一段 信杂比的程序代码，但是我不知道怎么能把这两个程序联系在一起，请达人帮忙！不知道怎么把他那段程序改改，添加到我程序里。
原程序代码：
function abc
clear all
eps = 0.0000001;
npts = 5000;
del = 1./ 5000.;
t = 0.:del: 1.;
% 产生输入序列
inp = 1+ t.^3 + 0.5 * t.^2 + cos(2 * pi * 5 * t);
x0 = [2, 0.1, 0.01]';
% 在几秒钟内更新间隔
T = del;
smcof = .8; nvar=0.05;
[e, y] = abr_filter(x0, smcof, inp, npts, T, nvar);
figure(1)
NN = 3499.;
n = 1: NN;
subplot(2,1,1)
plot(n * T , inp(1: NN),'k')
ylabel ('目标真实位置')
axis tight
subplot(2, 1, 2)
plot(n * T, y(1: NN),'k')
xlabel ('time(s)'), ylabel ('预测目标位置')
axis tight
figure(2)
plot([1: 500]* T, e(1: 500)','k')
axis([0 inf -0.5 0.5])
xlabel ('time(s)'), ylabel ('预留误差')
grid on
function [e, y] = abr_filter(x0, smcof, inp, npts, T , nvar)
rn=1.;
% 状态向量赋初值
X = x0;
theta = smcof;
% 计算αβγ的值
alpha = 1. - (theta^3);
beta= 1.5 * (1. + theta) * ((1. - theta)^2)/T;
gamma= ((1. - theta)^3)/(T^2);
% 装入状态转移矩阵PHI
PHI = [1. T (T^2)/2.;0. 1. T;0. 0. 1.];
while rn < npts;
% 用状态转移矩阵预测下一状态
XN = PHI * X;
error = (inp(rn) + normrnd(0,nvar)) - XN(1);
e(rn) = error;
tmp1 = alpha * error;
tmp2 = beta * error;
tmp3 = gamma * error;
% 预测下一状态
X(1) = XN(1) + tmp1;
X(2) = XN(2) + tmp2;
X(3) = XN(3) + tmp3;
y(rn) = X(1);
rn = rn +1.;
end
附加 信杂比 代码，帮忙看下能否加我程序里，怎么改？
clc;
clear all;
close all;
load 'Input.mat'
N=500;
x= 1:1:N;
Y=zeros(1,N);
Yd=zeros(1,N);
Yp=zeros(1,N);
M_IN=M(row,:);
Y(1:24)=M_IN(1:24);
Yd(1:24)=1;
gama=0.6;
alpha=0:0.01:0.9;
beta=0:0.001:0.9;
SCR_input=zeros(length(beta),length(alpha));
SCR_output=zeros(length(beta),length(alpha));
SCR_improve=zeros(length(beta),length(alpha));
% SCR_max=zeros(1,length(alpha))
% for i=1:length(alpha)
% for i=83
for k=83
for j=1:length(beta)
for m=24:N
Yp(m)=Y(m-1)+alpha(k)*Yd(m-1);
Y(m)=Yp(m)+(M_IN(m)-Yp(m))*beta(j);
Yd(m)=Yd(m-1)+(M_IN(m)-Yp(m)-Yd(m-1))*gama;
end
Output=M_IN-Y;
Input=M_IN;
SCR_input(j,k)=sum(Input(col-1:col+1).^2)/mean(Input.^2);
SCR_output(j,k)=sum(Output(col-1:col+1).^2)/mean(Output.^2);
SCR_improve(j,k)=10*log(SCR_output(j,k)/SCR_input(j,k));
end
end
[c,b]=max(SCR_improve);
[e,f]=max(max(SCR_improve));
% SCR_max(i)=e;
% end
figure;
plot(beta,SCR_improve(:,83))
% mesh(beta,alpha,SCR_improve);
xlabel('beta');
% ylabel('alpha');
title('the improvement of SCR (dB)');
axis([0 1 -40 20])
[c,b]=max(SCR_improve);
[e,f]=max(max(SCR_improve));
row=b(f);
col=f;
beta(row)
alpha(col)
e