阵列信号处理作业

航海学院 阵列信号处理 阵列信号处理课程2011年作业第1题假定半波长间隔均匀分布线列阵的阵元数N16，若入射平面波为62.5Hz的正弦信号，信号持续时间为0.4s，系统采样频率为1kHz，阵列加权方式为均匀加权。分别给出1. 当平面波信号分别从0,10, 20, 30, 40, 50, 60, 70, 80, 90, 100度方向入射时，指向90度的波束形成器的输出序列。2. 当平面波信号分别从0:1:180度方向入射时，指向90度的波束形成器的输出序列经过平方求和后的分贝数输出。（把所有181个输出绘制在同一幅图中）1）仿真图图一：所求角度入射信号输出序列三维表示注：1. 为信号入射角度，取值从0度到100度，每10为一个间隔；2. T为整个阵元采样时间，对于不同的入射角度，t的取值范围不同；3. 输出信号幅度表示所有阵元的求和输出幅度，为有噪声情况。结论：从图一可以看出：从90入射的信号输出序列没有得到衰减，而其它角度入射的都得到了衰减；从100入射的信号和从80入射的信号输出序列关于90方向是对称的；整个阵列对噪声有很好的抑制作用。图二：入射信号0到50的输出序列图三：入射信号60到100的输出序列结论： 从图二和图三可以看出：图一的所有结论；90方向入射信号没有时延，0方向入射信号时延最长；在不同角度，信号衰减倍数不同。图四：所求角度输出序列分贝数结论：从图四可以看出：从90方向入射的信号输出序列分贝数最大，高出旁瓣近13dB；整个图形关于90方向对称。2）仿真程序clear all; close all; clc%以最先有信号的阵元为参考，信号采用正弦，考虑有加性白噪声的情况N = 16; fc = 62.5; fs = 1000; T = 0.4; %阵元数目、信号频率、采样频率和信号%持续时间% 第1题 （a）for theta = 0:10:100; %入射信号角度 t_delta = abs(cos(theta/180*pi)/(2*fc); %相邻阵元的时延 t_noise = 0:1/fs:T+(N-1)*t_delta; %整个信号+噪声采样时间 X = zeros(N+1,length(t_noise); for n = 1:N nn = (N+1)*(theta90) + (-1)(theta90)*n; %判断角度是否大于90 temp = floor(n-1)*t_delta*fs); if (n-1)*t_delta*fs-temp 90) + (-1)(theta90)*n; %判断角度是否大于90 temp = floor(n-1)*t_delta*fs); if (n-1)*t_delta*fs-temp = 10(-12) %判断信号起点是否位于采样点 t_signal = 0:1/fs:T; %信号采样点 X(nn,temp+1:temp+length(t_signal) = sin(2*pi*fc*t_signal); %存储信号 else t_signal = (temp+1-(n-1)*t_delta*fs)/fs:1/fs:T; %信号采样点 X(nn,temp+2:temp+1+length(t_signal) = sin(2*pi*fc*t_signal); %存储信号 end noise = 0.01*randn(size(t_noise); %生成噪声，如果不需要，将方差设为0 X(nn,:) = X(nn,:) + noise; %信号+噪声 X(N+1,:) = X(N+1,:) + X(nn,:); %所有阵元采样求和 end X_out = X(N+1,:)/16; X_sum = 0; %采用均匀加权 X_dB(theta+1) = 10*log10(sum(X_out.2); %求分贝（没有进行归一）endfigure; plot(0:180,X_dB); grid on;xlabel(it theta); ylabel(输出信号求和：dB);注：1. 所给程序为有噪声的情况。限于篇幅，画图部分程序稍有改变，但不影响最终结果。2. 程序包含第一问的三维输出序列（图一）和第二问的经过平方求和取分贝的图形（图四），在程序中有说明。3. 对于图二和图三的程序，可以很容易的根据图一的程序得出。但由于篇幅所限，没有给出具体程序；第2题考虑一个10元标准线列阵（阵元间隔为1000Hz对应波长的一半），该阵列接收从远场入射的两个平面波连续信号，它们的入射角和频率如下表所示。设这两个信号的信噪比均为5dB，请针对表中所示的入射角和信号频率，分别给出常规波束形成法和MUSIC方法的空间方位谱。计算时假定观测时间为0.1秒，采样频率为8192Hz，阵元噪声为空间白噪声，且各阵元上噪声功率相等。说明1：对于每一种入射角和信号频率情况，把两种方法给出的空间方位谱放在同一图中。说明2：信噪比统一定义为信号功率和噪声功率之比。信号1频率（Hz）960999960999入射角（度）85859090信号2频率（Hz）1040100110401001入射角（度）959595951）仿真图图五：入射角度分别为85和95频率分别为960Hz和1040Hz信号的常规波束搜索和MUSIC算法仿真结果图六：入射角度分别为85和95频率分别为999Hz和1001Hz信号的常规波束搜索和MUSIC算法仿真结果图七：入射角度分别为90和95频率分别为960Hz和1040Hz信号的常规波束搜索和MUSIC算法仿真结果图八：入射角度分别为90和95频率分别为999Hz和1001Hz信号的常规波束搜索和MUSIC算法仿真结果结论：从上图五、图六、图七和图八可以看出：对于所给频率和入射角度的两个信号，MUSIC算法在图五、图六和图七信号源位置准确的显示了两个分离的谱峰。在图八中也出现了两个谱峰，但谱峰之间的凹陷非常小，以至于信号可能无法被分辨；常规波束搜索算法在图五和图六也出现了两个谱峰，但在图五谱峰之间的凹陷也非常小，以至于信号可能无法被分辨。从图五、图六、图七和图八对比中我们可以得出：在相同的信噪比情况下，对于常规波束算法，如果两个信号频率间距不变，但入射角度间距缩小，或入射角度间距不变，频率间距扩大，或入射角度间距缩小，频率间距扩大，可能导致信号从可被分辨变为无法分辨；对于MUSIC算法，如果信号入射角度间距不变，但频率间距缩小，可能导致信号从可被分辨变为无法分辨。在信号可被分辨的情况下，MUSIC算法表现出了比常规波束更好的性能。2）仿真程序clear all; close all; clc%以第一个阵元为坐标原点，信号采用复指数信号，分别使用常规波束形成算法和MUSIC（多重%信号分类）算法N = 10; fs = 8192; T = 0.1; %阵元数目、采样频率和信号持续时间SNR = 5; f0 = 1000; Nfft = 2048; %信噪比、阵列间隔对应频率和FFT点数f1 = 960 999 960 999; 85 85 90 90;f2 = 1040 1001 1040 1001; 95 95 95 95; %所给信号条件t_sample = 0:1/fs:T; phi = linspace(-1,1,Nfft); A = sqrt(2*10(SNR/10);noise = randn(N,length(t_sample) + 1j*randn(N,length(t_sample); %噪声S_1 = zeros(N,length(t_sample); S_2 = zeros(N,length(t_sample);for i = 1:4 for n = 0:N-1 %生成所给频率和入射角的信号 S_1(n+1,:) = A*exp(1j*(2*pi*f1(1,i)*t_sample + pi*f1(1,i)/f0*cos(f1(2,i)/180*pi)*n); S_2(n+1,:) = A*exp(1j*(2*pi*f2(1,i)*t_sample + pi*f2(1,i)/f0*cos(f2(2,i)/180*pi)*n); end signal = S_1 + S_2 + noise; P_b = 0; for n = 1:length(t_sample) %采用常规波束搜索 P_b = P_b + abs(fftshift(fft(signal(:,n).,Nfft).2; %使用FFT%来提高速度 end P_b = P_b/length(t_sample); figure; plot(acos(phi)/pi*180,10*log10(P_b),-); hold on; xlabel(it theta); ylabel(it SNR_0(dB); grid on; C_x = 1/length(t_sample)*(signal*signal); %采样谱矩阵 v d = eig(C_x); P_music = 0; %得到特征值和特征向量 for n = 1:N-2 %采用MUSIC方法 P_music = P_music + abs(fftshift(fft(v(:,n).,Nfft).2; %使用%FFT来提高速度 end P_music = 1./P_music; plot(acos(phi)/pi*180,10*log10(P_music),-); hold off; legend(波束搜索,MUSIC);end注：1. 程序包含常规波束搜索算法和MUSIC（多重信号分类）算法（在程序相应部分有说明）。2. 两种方法都采用FFT算法来提高运算速度，简化程序编程（在程序相应部分有说明）。3. 程序仿真结果对应图五、图六、图七和图八。4. 限于篇幅，画图部分程序稍有改变，但不影响最终结果。2011年5月29日完成9