gusucode.com > 均匀线阵和均匀圆阵的DOA估计matlab源码程序 > code/ULCAMUS.m
%-------JIN Wei %----ULA and UCA estimation based on the MUSIC method %Remark 1: The UCA estimation has a worse performance than ULA %estimation under the MUSIC ,espatially in the low snr,though it has the same resolution in %360 degree. The reason is that the vector of the UCA convariance matrix is %not orthogonal to each other. %Remark 2: The program is based on the assumption of narrowband signal %model. %Remark 3:If we increase the element number or snr, we will find the the %merit of UCA array. clc; clear; clear all; msens = 16; nsamp = 256; fc = 10000;%the carrier frequency f = [123,183,263,333];%the signal frequency ts = 1/2/fc;%the sampling time inteval v=3*10^8;%the speed of electromagnetic wave lambda=v/fc;% msource=4; bearing = [20 40 65 85]*pi/180; %--------------ULA,uniform linear array DOA estimation------------- dspace = 0.5; for i=1:length(bearing) smat_ULA(i,:)=hilbert(cos(2*pi*f(i)*(0:nsamp-1)*ts).*cos(2*pi*fc*(0:nsamp-1)*ts));%signal matrix of the sources end %-----------------generate data-------------------------------- amat_ULA = exp(sqrt(-1) * [0:msens-1]' * (2*pi*dspace*sin(bearing)));%array manifold matrix smat_ULA = amat_ULA*smat_ULA;%signal matrix of the array elements theta=0:pi/100:pi/2; for run=1:30 %-----------signal + noise------------------------------------- snr=20; for i=1:msens, noise= hilbert(randn(1,nsamp)); confi=(norm(smat_ULA(i,:))^2)/(10^(snr/10))/(norm(noise)^2); smat_ULA(i,:)=smat_ULA(i,:)+noise.*sqrt(confi); end R=smat_ULA*smat_ULA'/nsamp; [V,eigvalue]=eig(R); [Y,i]=sort(diag(eigvalue)); amat1 = exp(sqrt(-1) * [0:msens-1]' * (2*pi*dspace*sin(theta))); P=zeros(size(theta')); for k=1:msens-msource, P=P+abs(amat1'*V(:,i(k))).^2; end P=1./P; % P=P/max(P); P=10*log10(P); figure(1) plot(theta*180/pi,P,'b') hold on end %--------------UCA,uniform circular array DOA estimation------------- R = lambda/4/sin(pi/msens); beta=2*pi*R/lambda; for i=1:length(bearing) smat_UCA(i,:)=hilbert(cos(2*pi*f(i)*(0:nsamp-1)*ts).*cos(2*pi*fc*(0:nsamp-1)*ts));%signal matrix of the sources end %-----------------generate data-------------------------------- for m=1:length(bearing) amat_UCA(:,m) = exp(sqrt(-1)*(beta*cos(bearing(m)-[0:msens-1]'*2*pi/msens)));%array manifold matrix end a=size(amat_UCA) b=size(smat_UCA) smat_UCA = amat_UCA*smat_UCA;%signal matrix of the array elements theta=0:pi/100:pi/2; for run=1:30 %-----------signal + noise------------------------------------- snr=20; for i=1:msens, noise= hilbert(randn(1,nsamp)); confi=(norm(smat_UCA(i,:))^2)/(10^(snr/10))/(norm(noise)^2); smat_UCA(i,:)=smat_UCA(i,:)+noise.*sqrt(confi); end R=smat_UCA*smat_UCA'/nsamp; [V,eigvalue]=eig(R); [Y,i]=sort(diag(eigvalue)); for m=1:length(theta) amat1(:,m) = exp(sqrt(-1)*(beta*cos(theta(m)-[0:msens-1]'*2*pi/msens))); end P=zeros(size(theta')); for k=1:msens-msource, P=P+abs(amat1'*V(:,i(k))).^2; end P=1./P; % P=P/max(P); P=10*log10(P); figure(2) plot(theta*180/pi,P,'r') hold on end