LMS、RLS和SMI这三个程序该如何改
现有三个程序，都是基于均匀线阵的，大家讨论一下该如何改成非均匀的%LMS
clear all;
close all;
clc;
f1=0;
f2=0;
f3=2;
f4=130;
global f0;
f0=100e6;
d1=60*pi/180;
d2=-60*pi/180;
d3=-40*pi/180;
d4=20*pi/180;
SP=1000;%采样点数
N=16;%阵元数
A=';%A为阵列
Ts=1/(300e3);%采样率300kHz
y1=50*Ts;
y2=100*Ts;
y3=120*Ts;
y4=200*Ts;
t=linspace(Ts,SP*Ts,SP);
Ps4=0;
K=200;
%w=zeros(N,K);
for n=1:SP
S1(1,n)=s(t(1,n)-y1)*exp(-j*2*pi*f1*t(1,n));%直达波
S2(1,n)=s(t(1,n)-y2)*exp(-j*2*pi*f2*t(1,n));%地反射
S3(1,n)=s(t(1,n)-y3)*exp(-j*2*pi*f3*t(1,n));%地物反射
S4(1,n)=s(t(1,n)-y4)*exp(-j*2*pi*f4*t(1,n));%目标反射
Ps4=Ps4+(abs(S4(1,n)))^2;
end
Ps=N^2*Ps4/SP;
S=;
M=;%幅度
M=diag(M);
SNR=2;
Pn=Ps/(10^(N*SNR/10));%噪声功率
Wc=zeros(N,200);
e=zeros(SP,1);%%%%
y=zeros(SP,1);%%%%
S1=1e5*S1;
for i=1:200
noise=wgn(N,SP,10*log10(Pn),'complex');
X=A*M*S+noise;
%d=zeros(1,SP);
w1=zeros(N,1);
y(1)=w1'*X(:,1);
e(1)=S1(:,1)-y(1);
for n=2:SP
w1=w1+0.00000000000005*X(:,n-1)*conj(e(n-1));
y(n)=w1'*X(:,n);
e(n)=S1(1,n)-y(n);
end
Wc(:,i)=w1;
% figure;%%%
% plot(1:SP,abs(e));%%%
% title('error');%%%
end
Wc=sum(Wc,2)/200;%求均值
Y=Wc'*X;
cc=0;
for d=-pi/2:0.01:pi/2
cc=cc+1;
A1(:,cc)=.';
y_d(1,cc)=Wc'*A1(:,cc);
end
d=-pi/2:0.01:pi/2;
d=d.*(180/pi);
subplot(121);
plot(d,10*log10(abs(y_d).^2));grid on;
title('阵列输出功率增益方向图');
xlabel('到达角(deg)');
ylabel('阵列输出功率增益(dB)');
subplot(122);
plot(d,abs(y_d));grid on;
title('阵列输出幅度增益方向图');
xlabel('到达角(deg)');
ylabel('阵列输出幅度增益(dB)');
[ 本帖最后由 eight 于 2007-11-15 13:03 编辑 ]
%RLS,forgeting facter==1
%RLS,forgeting facter==1clear all;
close all;
clc;
f1=0;
f2=0;
f3=2;
f4=130;
global f0;
f0=100e6;
d1=60*pi/180;
d2=-60*pi/180;
d3=-40*pi/180;
d4=20*pi/180;
SP=2000;%采样点数
N=16;%阵元数
A=';%A为阵列
Ts=1/(300e3);%采样率300kHz
y1=500*Ts;
y2=100*Ts;
y3=120*Ts;
y4=200*Ts;
t=linspace(Ts,SP*Ts,SP);
Ps1=0;
K=200;
%w=zeros(N,K);
for n=1:SP
S1(1,n)=s(t(1,n)-y1)*exp(-j*2*pi*f1*t(1,n));%目标反射
S2(1,n)=s(t(1,n)-y2)*exp(-j*2*pi*f2*t(1,n));%地反射
S3(1,n)=s(t(1,n)-y3)*exp(-j*2*pi*f3*t(1,n));%地物反射
S4(1,n)=s(t(1,n)-y4)*exp(-j*2*pi*f4*t(1,n));% 直达波
Ps1=Ps1+(abs(S1(1,n)))^2;
end
Ps=N^2*Ps1/SP;
S=;
M=;%幅度
M=diag(M);
SNR=2;
Pn=Ps/(10^(N*SNR/10));%噪声功率
S1=1e5*S1;
for i=1:200
noise=wgn(N,SP,10*log10(Pn),'complex');
X=A*M*S+noise;
I=eye(N);
p=(1/0.004)*I;
w1=zeros(N,1);
e=zeros(1,SP);
e(1,1)=S1(1,1);
for n=2:SP
kg=p*X(:,n)/(1+X(:,n)'*p*X(:,n));
e(n)=conj(S1(1,n))-X(:,n)'*w1;
Wc(:,i)=w1+kg*(e(n));
p=p-kg*X(:,n)'*p;
w1=Wc(:,i);
end
end
Wc=sum(Wc,2)/200;%求均值
Y=Wc'*X;
cc=0;
for d=-pi/2:0.01:pi/2
cc=cc+1;
A1(:,cc)=.';
y_d(1,cc)=Wc'*A1(:,cc);
end
d=-pi/2:0.01:pi/2;
d=d*180/pi;
subplot(121);
plot(d,10*log10(abs(y_d).^2));grid on;
title('阵列输出功率增益方向图');
xlabel('到达角(deg)');
ylabel('阵列输出功率增益(dB)');
subplot(122);
plot(d,abs(y_d));grid on;
title('阵列输出幅度增益方向图');
xlabel('到达角(deg)');
ylabel('阵列输出幅度增益(dB)');
%SMI算法
%SMI算法clear all;
close all;
clc;
f1=0;
f2=0;
f3=2;
f4=130;
global f0;
f0=100e6;
d1=60*pi/180;
d2=-60*pi/180;
d3=-40*pi/180;
d4=20*pi/180;
SP=2000;%采样点数
N=16;%阵元数
A=';%A为阵列
Ts=1/(300e3);%采样率300kHz
y1=500*Ts;
y2=100*Ts;
y3=120*Ts;
y4=200*Ts;
t=linspace(Ts,SP*Ts,SP);
Ps1=0;
for n=1:SP
S1(1,n)=s(t(1,n)-y1)*exp(-j*2*pi*f1*t(1,n));%目标直达
S2(1,n)=s(t(1,n)-y2)*exp(-j*2*pi*f2*t(1,n));%地反射
S3(1,n)=s(t(1,n)-y3)*exp(-j*2*pi*f3*t(1,n));%地物反射
S4(1,n)=s(t(1,n)-y4)*exp(-j*2*pi*f4*t(1,n));%反射波
Ps1=Ps1+(abs(S1(1,n)))^2;
end
Ps=N^2*Ps1/SP;
S=;
M=;%幅度
M=diag(M);
SNR=2;
S1=1e5*S1;
for c=1:200
Pn=Ps/(10^(N*SNR/10));%噪声功率
W=wgn(N,SP,10*log10(Pn),'complex');
X=A*M*S+W;
R=X*X'/SP;
Rd=X*S1'/SP;
Wc(:,c)=inv(R)*Rd;%权向量
end
Wc=sum(Wc,2)/200;%求均值
Y=Wc'*X;
cc=0;
for d=-pi/2:0.01:pi/2
cc=cc+1;
A1(:,cc)=.';
y_d(1,cc)=Wc'*A1(:,cc);
end
d=-pi/2:0.01:pi/2;
d=d.*(180/pi);
subplot(121);
plot(d,10*log10(abs(y_d).^2));grid on;
title('阵列输出功率增益方向图');
xlabel('到达角(deg)');
ylabel('阵列输出功率增益(dB)');
subplot(122);
plot(d,abs(y_d));grid on;
title('阵列输出幅度增益方向图');
xlabel('到达角(deg)');
ylabel('阵列输出幅度增益(dB)');
function output=s(t)
function output=s(t)%%%%%%%%%%%%%%%%%LFM%%%%%%%%%%%%%%%%%%%%%
% global b;
% global T;
% F=1;
% global f0;
% ift>0&t<T
% output=exp(j*(2*pi*f0*t+b*t^2));%信号源
% else
% output=0;
% end
%%%%%%%%%%%%%%%%%FM调频广播%%%%%%%%%%%%%
F=1;
Mf=3;
W=2*pi*15e3;%带宽180kHz
global f0;
ift>0
output=F*exp(j*(2*pi*f0*t+Mf*sin(W*t)));
else
output=0;
end
% %%%%%%%%%%%%%%%电视伴音信号%%%%%%%%%%%
%F=1;
% Mf=50/15;
% global f0;
% W=2*pi*15e3;
% if t>0
% output=F*exp(j*(2*pi*f0*t+Mf*sin(W*t)));
% else
% output=0;
% end
% %%%%%%%%%%%%%%电视图象信号%%%%%%%%%%%
% F=1;
% global f0;
% W1=2*pi*3e6;W2=2*pi*0.5e3;
% if t>=0
% output=F*(1+sin(W1*t)+sin(W2*t))*exp(j*(2*pi*f0*t));
% else
% output=0;
% end
谢谢你的程序
谢谢楼主的程序，我毕业设计要写智能天线的波束形成算法仿真，就有你这三个程序的比较，哈哈，真是太感谢你了但是，你的程序是不是有点bug，就是那个输入信号，就是我们的有用信号，怎么运行时老报错！ 建议楼主多看看本版的置顶帖
页:
[1]