求根MUSIC算法源程序
合集下载
相关主题
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
f=a2(I(1:source_number));
%计算信号到达方向角
source_doa=[asin(angle(f(1))/pi)*180/pi asin(angle(f(2))/pi)*180/pi];
source_doa=sort(source_doa);
disp('source_doa');
s=10.^((snr/2)/10)*exp(j*w*[0:N_x-1]);%仿真信号
%x=awgn(s,snr);
x=A*s+(1/sqrt(2))*(randn(sensor_number,N_x)+j*randn(sensor_number,N_x));%加了高斯白噪声后的阵列接收信号
disp(source_doa);
a(i+sensor_number) = sum( diag(Gn,i) );
end
%使用ROOTS函数求出多项式的根
a1=roots(a);
%除掉大于1的增根（方程式在此条件下变形的）
a2=a1(abs(a1)<1);
%挑选出最接近单位圆的N个根
[lamda,I]=sort(abs(abs(a2)-1));
%ROOT_MUSIC ALOGRITHM
%DOA ESTIMATION BY ROOT_MUSIC
clear all;
close all;
clc;
source_number=2;%信元数
sensor_number=8;%阵元数
N_x=1024; %信号长度
snapshot_number=N_x;%快拍数
w=[pi/4 pi/6].';%信号频率
l=((2*pi*3e8)/w(1)+(2*pi*3e8)/w(2))/2;%信号波长
d=0.5*l;%阵元间距
snr=0;%信噪比
theta1=45;theta2=60;%两个信号的入射角度
A=[exp(-j*(0:sensor_number-1)*d*2*pi*sin(theta1*pi/180)/l);exp(-j*(0:sensor_number-1)*d*2*pi*sin(theta2*pi/180)/l)].';%阵列流型
R=(x*x')/N_x;
%对协方差矩阵进行特征分解
[V,D]=eig(R);
Un=V(:,1:sensor_number-source_number);
Gn=Un*Un';
%找出多项式的系数，并按阶数从高至低排列
a = zeros(2*sensor_number-1,1)';
for i=-(sensor_number-1):(sensor_number-1)
%计算信号到达方向角
source_doa=[asin(angle(f(1))/pi)*180/pi asin(angle(f(2))/pi)*180/pi];
source_doa=sort(source_doa);
disp('source_doa');
s=10.^((snr/2)/10)*exp(j*w*[0:N_x-1]);%仿真信号
%x=awgn(s,snr);
x=A*s+(1/sqrt(2))*(randn(sensor_number,N_x)+j*randn(sensor_number,N_x));%加了高斯白噪声后的阵列接收信号
disp(source_doa);
a(i+sensor_number) = sum( diag(Gn,i) );
end
%使用ROOTS函数求出多项式的根
a1=roots(a);
%除掉大于1的增根（方程式在此条件下变形的）
a2=a1(abs(a1)<1);
%挑选出最接近单位圆的N个根
[lamda,I]=sort(abs(abs(a2)-1));
%ROOT_MUSIC ALOGRITHM
%DOA ESTIMATION BY ROOT_MUSIC
clear all;
close all;
clc;
source_number=2;%信元数
sensor_number=8;%阵元数
N_x=1024; %信号长度
snapshot_number=N_x;%快拍数
w=[pi/4 pi/6].';%信号频率
l=((2*pi*3e8)/w(1)+(2*pi*3e8)/w(2))/2;%信号波长
d=0.5*l;%阵元间距
snr=0;%信噪比
theta1=45;theta2=60;%两个信号的入射角度
A=[exp(-j*(0:sensor_number-1)*d*2*pi*sin(theta1*pi/180)/l);exp(-j*(0:sensor_number-1)*d*2*pi*sin(theta2*pi/180)/l)].';%阵列流型
R=(x*x')/N_x;
%对协方差矩阵进行特征分解
[V,D]=eig(R);
Un=V(:,1:sensor_number-source_number);
Gn=Un*Un';
%找出多项式的系数，并按阶数从高至低排列
a = zeros(2*sensor_number-1,1)';
for i=-(sensor_number-1):(sensor_number-1)