阵列通道不一致现象仿真-MATLAB源码

function MUSIC_Loss_Analy(X_AP_Ident,X_AP_Inconf,N,N_Targ,Array_Sel,theta_limit,phi_limit,theta_step,phi_step,lambda,Param_ULA,Param_UPA,Param_UCA) % Create Time :2022 03 18 % Finish Time :2022 03 18 % Author :Xu Y. B.（CSDN ID:在路上，正出发） % MATLAB Edition:R2021a % 函数功能：分析加入幅相误差后 MUSIC 算法的性能损失情况 % 参数声明： % 输入： % X_AP_Ident : 各通道幅相完全一致的阵列信号 % X_AP_Inconf: 各通道存在幅相误差的阵列信号 % N : 快拍点数 % N_Targ : 目标数目 % Array_Sel : 阵列类型 % 均匀线阵 'ULA' % 均匀面阵 'UPA' % 均匀圆阵 'UCA' % theta_limit: 俯仰角搜索范围 单位：° % phi_limit ：方位角搜索范围 单位：° % theta_step ：俯仰角搜索步长 单位：° % phi_step ：方位角搜索步长 单位：° % lambda ：电磁波波长 单位：m % Param_ULA ：均匀线阵的参数集合（数组） % 数组包含 个元素，分别的定义如下： % 索引号 1：阵元个数 % 索引号 2：阵元间隔 单位：m % Param_UPA ：均匀面阵的参数集合（数组） % 数组包含 个元素，分别的定义如下： % 索引号 1：x轴阵元个数 % 索引号 2：y轴阵元个数 % 索引号 2：阵元间隔 单位：m % Param_UCA ：均匀圆阵的参数集合（数组） % 数组包含 个元素，分别的定义如下： % 索引号 1：阵元个数 % 索引号 2：圆阵半径 单位：m % 输出： % %% 实现：↓ % -1- 对幅相无差别通道阵列数据 求噪声子空间 % -1.1- 求协方差矩阵 Rx_AP_Ident = X_AP_Ident*X_AP_Ident'/N; % -1.2- 特征值分解 [Rx_AP_Id_v,Rx_AP_Id_eig] = eig(Rx_AP_Ident); Rx_AP_Id_eig = diag(Rx_AP_Id_eig); % 取特征值 转为列向量 % -1.3- 特征值降序排列 [~,Rx_eig_index_Id] = sort(abs(Rx_AP_Id_eig),'descend'); % -1.4- 构建噪声子空间 Un_Id = Rx_AP_Id_v(:,Rx_eig_index_Id(N_Targ+1:end)); % -2- 对幅相有差别通道进行 MUSIC 算法处理 % -2.1- 求协方差矩阵 Rx_AP_Inconf = X_AP_Inconf*X_AP_Inconf'/N; % -2.2- 特征值分解 [Rx_AP_In_v,Rx_AP_In_eig] = eig(Rx_AP_Inconf); Rx_AP_In_eig = diag(Rx_AP_In_eig); % 取特征值 转为列向量 % -2.3- 特征值降序排列 [~,Rx_eig_index_In] = sort(abs(Rx_AP_In_eig),'descend'); % -2.4- 构建噪声子空间 Un_In = Rx_AP_In_v(:,Rx_eig_index_In(N_Targ+1:end)); % -3- MUSIC搜索算法 % -3.1- 搜索配置 theta_search = theta_limit(1):theta_step:theta_limit(2); phi_search = phi_limit(1):phi_step:phi_limit(2); % -3.2- 空间谱函数初始化 if(strcmp(Array_Sel,'ULA')) P_MUSIC = zeros(2,length(theta_search)); else P_MUSIC = zeros(length(phi_search),length(theta_search),2); end % -3.2- 搜索 if(strcmp(Array_Sel,'ULA')) for i = 1:2 if(i==1) Un = Un_Id; else Un = Un_In; end for j = 1:length(theta_search) a = exp(-1j*2*pi/lambda*Param_ULA(2)*((0:Param_ULA(1)-1).')*sin(theta_search(j)/180*pi)); P_MUSIC(i,j) = 1/(a'*(Un*Un')*a); end end elseif(strcmp(Array_Sel,'UPA')) for k = 1:2 for i = 1:length(phi_search) if(k == 1) Un = Un_Id; else Un = Un_In; end for j = 1:length(theta_search) u = sin(phi_search(i)/180*pi).*sin(theta_search(j)/180*pi); % 矢量 u v = cos(phi_search(i)/180*pi).*sin(theta_search(j)/180*pi); % 矢量 v ax = exp(-1j*2*pi*Param_UPA(3)*((0:Param_UPA(1)-1).')*v/lambda); ay = exp(-1j*2*pi*Param_UPA(3)*((0:Param_UPA(2)-1).')*u/lambda); a = Kronecker_Mult_Cal(ay,ax); P_MUSIC(i,j,k) = 1/(a'*(Un*Un')*a); end end end else for k = 1:2 for i = 1:length(phi_search) if(k == 1) Un = Un_Id; else Un = Un_In; end for j = 1:length(theta_search) a = exp(1j*2*pi*Param_UCA(2)/lambda*... (cos(phi_search(i)/180*pi-(((0:Param_UCA(1)-1).')*2*pi/Param_UCA(1))))*... sin(theta_search(j)/180*pi)); P_MUSIC(i,j,k) = 1/(a'*(Un*Un')*a); end end end end % 作图分析 if(strcmp(Array_Sel,'ULA')) figure; subplot(211) plot(theta_search,abs(P_MUSIC(1,:)),'r') axis tight xlabel('俯仰角/°') ylabel('谱函数数值') title([Array_Sel,' 无通道误差阵列信号处理——DOA估计图']) subplot(212) plot(theta_search,abs(P_MUSIC(2,:)),'b.-') axis tight xlabel('俯仰角/°') ylabel('谱函数数值') title([Array_Sel,' 有通道误差阵列信号处理——DOA估计图']) else figure; colormap('jet') subplot(121) imagesc(theta_search,phi_search,abs(P_MUSIC(:,:,1))) colorbar axis tight axis xy xlabel('方位角/°') ylabel('俯仰角/°') title([Array_Sel,' 无通道误差阵列信号处理——二维DOA估计图']) subplot(122) imagesc(theta_search,phi_search,abs(P_MUSIC(:,:,2))) colorbar axis tight axis xy xlabel('方位角/°') ylabel('俯仰角/°') title([Array_Sel,' 有通道误差阵列信号处理——二维DOA估计图']) end end