分别通过LMS,RLS,SMI三种自适应滤波方法实现雷达信号的波束形成matlab仿真,包含仿真操作录像,代码注释

clc; clear; close all; warning off; addpath(genpath(pwd)); rng('default') f1=0; % 第一个信号的频率，初始化为0 f2=0; % 第二个信号的频率，初始化为0 f3=2; % 第三个信号的频率，初始化为2 f4=130; % 第四个信号的频率，初始化为130 global f0;% 设定f0的值为100e6 f0=100e6; d1=60*pi/180;% 第一个信号的入射角，单位是弧度 d2=-60*pi/180;% 第二个信号的入射角，单位是弧度 d3=-40*pi/180;% 第三个信号的入射角，单位是弧度 d4=20*pi/180;% 第四个信号的入射角，单位是弧度 SP=2000;%采样点数 N=16;%阵元数 A=[conj(exp(j*(0:N-1)*pi*sin(d1)));conj(exp(j*(0:N-1)*pi*sin(d2)));conj(exp(j*(0:N-1)*pi*sin(d3)));conj(exp(j*(0:N-1)*pi*sin(d4)))]';%A为阵列 Ts=1/(300e3);%采样率300kHz，即每个采样点的时间间隔为1/(300e3)秒 y1=50*Ts; % 第一个信号的延迟，单位是秒 y2=100*Ts; % 第2个信号的延迟，单位是秒 y3=120*Ts; % 第3个信号的延迟，单位是秒 y4=200*Ts; % 第4个信号的延迟，单位是秒 t=linspace(Ts,SP*Ts,SP);% 生成一个从Ts到SP*Ts的等差数列，用于模拟时间轴 Ps4=0;% 初始化S4的功率为0 K=200;% 设定一个参数K为200 % 循环从1到SP for n=1:SP S1(1,n)=(t(1,n)-y1)*exp(-j*2*pi*f1*t(1,n));%目标反射 S2(1,n)=(t(1,n)-y2)*exp(-j*2*pi*f2*t(1,n));%地反射 S3(1,n)=(t(1,n)-y3)*exp(-j*2*pi*f3*t(1,n));%地物反射 S4(1,n)=(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=[S1;S2;S3;S4];% 将S1，S2，S3，S4拼接起来形成S矩阵 M=[1e5 1e4 1e4 1];%幅度，设定四个信号的幅度分别为1e5，1e4，1e4，1 M=diag(M);% 将M转换为对角矩阵 SNR=2;% 信噪比设定为2 Pn=Ps/(10^(N*SNR/10));%噪声功率 for i=1:200% 循环从1到200 %噪声功率，根据公式计算噪声功率 noise=wgn(N,SP,10*log10(Pn),'complex');% 生成复高斯白噪声 X=A*M*S+noise;% 根据阵列响应、信号幅度和信号源以及噪声生成接收信号 I=eye(N); p=(1/0.004)*I;% 生成一个N阶单位矩阵I w1=zeros(N,1); % 定义一个矩阵p，初始值为I/0.004 e=zeros(1,SP);% 定义一个N*1的矩阵w1，初始值为0 e(1,1)=S4(1,1); for n=2:SP% 从2开始循环，直到SP kg=p*X(:,n)/(1+X(:,n)'*p*X(:,n));% 计算卡尔曼增益 e(n)=conj(S4(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;%求均值% 对Wc的每一列求和，然后除以200，得到平均值 Y=Wc'*X;% 计算Y cc=0;% 初始化cc为0 for d=-pi/2:0.01:pi/2% 从-pi/2到pi/2，步长为0.01，进行循环 cc=cc+1;% 计算A1的每一列 A1(:,cc)=[exp(j*(0:N-1)*pi*sin(d))].'; % 计算y_d的每一列 y_d(1,cc)=Wc'*A1(:,cc); end d=-pi/2:0.01:pi/2;% 将d从弧度转换为角度 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)');