陣列信號處理--基礎概念

假設一平面波傳播方向為\(\mathbf{a}\)，頻率為\(w\),那么經過空間中的麥克風陣列系統處理后，可以得到該麥克風陣列對平面波的響應為

\[y(t,\mathbf{k})=\mathbf{H}^T\mathbf{v}_{k}(\mathbf{k})e^{jwt} \]

其中\(\mathbf{H}\)表示的是濾波器沖激響應的傅里葉變換，\(\mathbf{k}\)表示波數，其幅度最大值為\(\lvert\mathbf{k}\rvert=\frac{2\pi}{\lambda}\)，\(\mathbf{v}_{k}\)表示陣列波形矢量，可以表示成

\[\mathbf{v}_{k}(\mathbf{k})=\left[\begin{matrix} e^{-\mathbf{k^Tp_{0}}}\\e^{-\mathbf{k^Tp_{1}}}\\\vdots\\e^{-\mathbf{k^Tp_{N-1}}} \end{matrix} \right] \]

然后對\(y(t,\mathbf{k})\)進行傅里葉變換，可以得到

\[\mathbf{Y}(w,\mathbf{k})=\mathbf{H}^{T}(w)\mathbf{v}_{k}(\mathbf{k}) \]

將后面這一項定義為\(\Upsilon(w,\mathbf{k})\),即陣列的頻率-波數響應函數，描述了一個陣列對於波數 \(\mathbf{k}\)，時頻頻率\(w\)的輸入平面波的復增益。

\[\Upsilon(w,\mathbf{k})=\mathbf{H}^T(w)\mathbf{v}_{k}(\mathbf{k}) \]

波束方向圖是用入射方向表示的頻率-波數響應函數，可以將其表示成

\[B(w;\theta,\phi)=\Upsilon(w,\mathbf{k})\rvert_{\mathbf{k}=\frac{2\pi}{\lambda}\mathbf{a}(\theta,\phi)} \]

定義一個復數加權矢量為

\[\mathbf{w}^{H}[w_{0}^{*},\cdots,w_{N-1}^{*}] \]

那么\(y(t,\mathbf{k})\)可以表示成

\[y(t,\mathbf{k})=\mathbf{w}^{H}\mathbf{v}_{k}(\mathbf{k})e^{jwt} \]

此時

\[\Upsilon(w,\mathbf{k})=\mathbf{w}^T(w)\mathbf{v}_{k}(\mathbf{k}) \]

對於均勻線陣而言

，\(\mathbf{v}_{k}(\mathbf{k})\)可以表示成

\[\mathbf{v}_{k}(k)=[e^{j\frac{N-1}{2}\mathbf{k}_{z}d},e^{j(\frac{N-1}{2}-1)\mathbf{k}_{z}d},\cdots,e^{-j\frac{N-1}{2}\mathbf{k}_{z}d}]^T \]

其中\(\mathbf{k}_z=-\frac{2\pi}{\lambda}cos\theta\),那么此時的\(\Upsilon(w,k_z)\)可以表示成

\[\Upsilon(w,k_z)=\mathbf{w}^H\mathbf{v}_{k}(k)\\=\sum_{n=1}^{N-1}w_{n}^{*}e^{-j(n-\frac{N-1}2)k_zd} \]

令\(\psi=-k_zd=\frac{2\pi}{\lambda}dcos\theta=\frac{2\pi}{\lambda}du_z\),進而可以將\(\Upsilon(w,k_z)\)表示成

\[\Upsilon(w,k_z)=e^{-j\frac{N-1}{2}\psi}\sum_{n=1}^{N-1}w_{n}^{*}e^{jn\psi} \]

波數方向圖可以表示為

\[B(\psi)=e^{-j\frac{N-1}{2}\psi}\sum_{n=1}^{N-1}w_{n}^{*}e^{jn\psi} \]

對於均勻加權線陣而言，此時\(w_n=\frac{1}{N}\),那么\(\Upsilon(w,k_z)\)和\(B(\psi)\)可以表示成

\[\Upsilon(\psi)=\frac{1}{N}\frac{sin(N\frac{\psi}{2})}{sin(\psi/2)}\\B(u)=\frac{1}{N}\frac{sin(\frac{\pi Nd}{\lambda}u)}{sin(\frac{\pi d}{\lambda}u)} -1\le u\le1 \]

可以得出\(\Upsilon(\psi)\)是一個周期函數，當N是奇數的時候，其周期為\(2\pi\)，當N是偶數的時候，其周期為\(4\pi\)，對於任意的N值，其周期為\(2\pi\)，下圖給出了當\(N=11\)的時候，\(\Upsilon(\psi)\)的變化趨勢。

下面分析波束方向圖的幾個重要參數

• 3dB帶寬(HPBW)
• 到第一零點的距離
• 到第一個旁瓣的距離
• 第一旁瓣的高度
• 其余零點的位置旁瓣衰減的速率
• 柵瓣

對於波束方向圖

\[B(u)=\frac{1}{N}\frac{sin(\frac{\pi Nd}{\lambda}u)}{sin(\frac{\pi d}{\lambda}u)} -1\le u\le1 \]

（1）3dB帶寬：定義為使\(\lvert B(\psi)\rvert^2=0.5\)的點。

（2）第一零點：方向圖的零點出現在\(B(\psi)\)的分子為0而分母不為0的時候，則零點的出現滿足下面的兩個條件

\[\frac{\pi Nd}{\lambda}u=m\pi,m=1,2,\cdots.\\ \Rightarrow u=m\frac{\lambda}{Nd}\\u \neq \frac{\lambda}{d} \]

那么第一個零點出現的位置就是\(u=\frac{\lambda}{Nd}\),那么定義\(\Delta u_2=2\frac{\lambda}{Nd}\)為零點-零點波束寬度(BWNN)，該值的一半(0.5BWNN)是到第一零點的距離，這個值衡量了陣列分辨兩個不同平面波的能力，即瑞利限。如果第二個波束方向圖的峰值在第一個波束方向圖的第一零點之外，那么就稱這個平面波是可以被分辨的。

(3)旁瓣位置：旁瓣極大值出現的位置在當分子出現極大值的時候，即

\[sin(\frac{\pi Nd}{\lambda}u)=1 \]

可以得到解

\[u=\pm \frac{2m+1}{N}\frac{\lambda}{2d} \]

(4)柵瓣：柵瓣就是和主瓣波束一樣高的波瓣。柵瓣出現分子和分母均為0的時候，柵瓣出現的間隔為

\[u=m.\frac{\lambda}{d} \]

如果陣列的間距大於\(\lambda\)，那么柵瓣的峰值出現在信號傳播區域內，即\(\lvert u \rvert \le1\)的區域。此時就會出現峰值響應模糊的問題。另外標准線陣指\(d=\lambda/2\)的均勻線陣。

上面三幅圖分別對應着\(d=\lambda/4;d=\lambda/2;d=\lambda\),可以看到當\(d=\lambda\)的時候波束方向圖出現了明顯的柵瓣。

上面三幅圖對應的代碼依次如下

clc
% 陣元個數 
N = 11;
psi = (-5 : 1/ 400 : 5) * pi;
beam = zeros(1, length(psi(:)));
y = sin(.5 * psi);
i = find( abs(y) > 1e-12);
j = 1 : length(psi(:));
j(i) = [];
beam(i) = sin((N/2)*psi(i))./(N*y(i));
beam(j) = sign(cos(psi(j) * ((N+1)/2)));

plot(psi / pi, beam);
grid on;
xlabel('{\it \psi}/\pi', 'Fontsize', 9)
ylabel('頻率波數響應函數','Fontsize', 9);
axis([-5 5 -0.4 1])

clc;
N = 10;                                % Elements in array
d = 0.5;                               % spacing wrt wavelength
beamwidth = 2/(N*d);                   % null-to-null (LL BW is half this)
D=d*[-(N-1)/2:1:(N-1)/2];              % element locations
u   = [-0.45:0.01:0.45]; 
AS  = exp(-j*2*pi*D'*0);               % BP points to 0
Au  = exp(-j*2*pi*D'*u);
B   = real(AS'*Au)/N;                          % BP

h=plot(u,B,'-');
set(h,'LineWidth',1.5)
hold on
% HPBW
I=find(B>=0.707);
plot([u(min(I)-1) u(max(I)+1)],[B(min(I)-1) B(max(I)+1)],'-')
plot(u(min(I)-1)*[1 1],B(min(I)-1)*[1 1]+0.03*[-1 1],'-')
plot(u(max(I)+1)*[1 1],B(max(I)+1)*[1 1]+0.03*[-1 1],'-')
plot([0.03 0.13],[0.72 0.845],'-')
h=text(0.155,0.84,'D');
set(h,'FontName','Symbol')
text(0.17,0.84,'{\itu}   = HPBW')
h=text(0.185,0.82,'1');
set(h,'Fontsize',10)
text(0.12, 0.7,'0.707')
% BW-NN
plot([-0.2 0.2],[-0.3 -0.3],'-')
plot(-0.2*[1 1],-0.3*[1 1]+0.03*[-1 1],'-')
plot(0.2*[1 1],-0.3*[1 1]+0.03*[-1 1],'-')
h=text(-0.05,-0.35,'D');
set(h,'FontName','Symbol')
text(-0.035,-0.35,'{\itu   = BWNN}')
% h=text(-0.02,-0.37,'2');
% set(h,'Fontsize',10)
% h=text(0.04,-0.37,'\itNN');
% set(h,'Fontsize',10)
% axes
plot([-1 1 ],[0 0],'-')
plot([0 0],[-0 1.1],'-')
% tick marks
plot(-0.4*[1 1],0.03*[-1 1])
plot(-0.2*[1 1],0.03*[-1 1])
plot(0.2*[1 1],0.03*[-1 1])
plot(0.4*[1 1],0.03*[-1 1])
text(-0.005,-0.05,'0')
%tick labels
xx=0.2;
yy=-0.1;
h=text(xx-0.01,yy+0.02,'l');
set(h,'FontName','Symbol')
plot([xx-0.02 xx+0.01],[yy-0.01 yy-0.01],'-')
h=text(xx-0.02,yy-0.05,'\itNd');

xx=-0.18;
yy=-0.1;
h=text(xx-0.01,yy+0.02,'l');
set(h,'FontName','Symbol')
plot([xx-0.02 xx+0.01],[yy-0.01 yy-0.01],'-')
h=text(xx-0.02,yy-0.05,'\itNd');
text(xx-0.04, yy-0.005,'-')

xx=0.41;
yy=-0.1;
h=text(xx-0.01,yy+0.02,'l');
plot([xx-0.02 xx+0.01],[yy-0.01 yy-0.01],'-')
h=text(xx-0.02,yy-0.05,'\itNd');
text(xx-0.035, yy-0.005,'2')

xx=-0.4;
yy=-0.1;
h=text(xx-0.01,yy+0.02,'l');
plot([xx-0.02 xx+0.01],[yy-0.01 yy-0.01],'-')
h=text(xx-0.02,yy-0.05,'\itNd');
text(xx-0.05, yy-0.005,'-2')

text(0.05,1.1,'\itB(u)')
hold off
axis([-0.42 0.42 -0.4 1.2])

%set(gca,'Box','off')
set(gca,'Visible','off')

clc;
N = 10;
n = (-(N-1)/2:(N-1)/2).';
psi = pi * (-3 : 0.001 : 3);
w = 1/N * [ones(N,1)];
d = [1/4 1/2 1];
figure
for k = 1 : length(d)
    vv = exp(1j * 2 * d(k) * n * psi);
    B(k,:) = 20 * log10(abs(w'*vv));
    subplot(3, 1, k);
    plot(psi / pi, B(k, :));
    hold on
    plot([-1 -1],[-25 0],'--');
    plot([1 1],[-25 0],'--');
    axis([-3 3 -25 5])
    set(gca,'YTick',[-25 -20 -15 -10 -5 0])
    grid on
    xlabel('\itu')
    if k == 1
        plot([-1 1],[2.5 2.5])
        plot([-1 -0.9],[2.5 0.5])
        plot([-1 -0.9],[2.5 4.5])
        plot([0.9 1],[0.5 2.5])
        plot([0.9 1],[4.5 2.5])
        text(-0.4,7,'Visible region')
    end
    hold off
end