确知波束形成matlab仿真

阵列信号处理中的导向矢量

假设一均匀线性阵列，有N个阵元组成，满足：远场、窄带假设。

图1. 均匀线性阵模型

假设信源发射信号，来波方向为$\theta$，第一个阵元接收到的信号为$x(t)$，则第$n$个阵元接收的信号可以表示为$X=x(t)e^{j2\pi \ast (n-1)\frac{dsin\theta }{\lambda }}$，可以得到阵列对信号的响应，即导向矢量（steering vector）：
$$a(\theta )=[1,e^{j2\pi \ast \frac{dsin\theta }{\lambda }},e^{j2\pi \ast 2\ast \frac{dsin\theta }{\lambda }},...,e^{j2\pi \ast (N-1)\ast \frac{dsin\theta }{\lambda }}{]}^T$$
那么阵列输出为：$y=a(\theta )x(t)$
matlab中可以一行代码来实现：

% Input : 
%   * N: 阵元数
%   * theta：来波方向
%   * dLembda：阵元间距与波长的比值，通常为1/2
% Output : 
%    * aTheta：导向矢量，列向量
aTheta = exp(1i*2*pi*dLembda*[0:N-1]'*sin(theta/180*pi));
1
2
3
4
5
6
7

确知波束形成

确知波束形成就是通过对各阵元接收到的信号进行加权求和的方式，使的波束指向一个确定的方向，经过波束形成之后，阵列接收的信号就由一个矩阵变为了一个向量，此时只是一个信号。假设来波方向为$\theta$，确知波束形成的权重，就是导向矢量的埃尔米特转置，即共轭转置：
$$W=a(\theta {)}^H=[1,e^{-j2\pi \ast \frac{dsin\theta }{\lambda }},e^{-j2\pi \ast 2\ast \frac{dsin\theta }{\lambda }},...,e^{-j2\pi \ast (N-1)\ast \frac{dsin\theta }{\lambda }}]$$
波束形成之后，输出信号为：
$$y=WX=Wa(\theta )x(t)$$
静态波束形成其实就是空间滤波，它对指定方向上的信号进行加强和选择性通过，而对其他方向的来波信号进行抑制，例如有一个信源发出信号为$x_1(t)$的波达方向为${\theta }_1=10^{\circ }$,而波束指向${\theta }_2=20^{\circ }$,此时
导向矢量（steering vector）：
$$a({\theta }_1)=[1,e^{j2\pi \ast \frac{dsin{\theta }_1}{\lambda }},e^{j2\pi \ast 2\ast \frac{dsin{\theta }_1}{\lambda }},...,e^{j2\pi \ast (N-1)\ast \frac{dsin{\theta }_1}{\lambda }}{]}^T$$
那么阵列输出为：$y_{=}a({\theta }_1)x_1(t)$。
指向$20^{\circ }$确知波束对阵元的加权的权重为：
$$W=a({\theta }_2{)}^H=[1,e^{-j2\pi \ast \frac{dsin{\theta }_2}{\lambda }},e^{-j2\pi \ast 2\ast \frac{dsin{\theta }_2}{\lambda }},...,e^{-j2\pi \ast (N-1)\ast \frac{dsin{\theta }_2}{\lambda }}]$$
波束形成之后，输出信号为：
$$y=WX=Wa({\theta }_1)x_1(t)=(a({\theta }_2){)}^{T}a({\theta }_1)x_1(t)=(1\ast 1+e^{j\ast 2\pi \frac{d}{\lambda }(sin{\theta }_1-sin{\theta }_2)}+e^{j\ast 2\pi \ast 2\ast \frac{d}{\lambda }(sin{\theta }_1-sin{\theta }_2)}+...+e^{j\ast 2\pi \ast (N-1)\ast \frac{d}{\lambda }(sin{\theta }_1-sin{\theta }_2)})x_1(t)$$
我们对其中的权重取绝对值进行分析，因为${\theta }_1\ne {\theta }_2$ ，所以每一项进行加权后的绝对值都小于1，而只有第一项等于1：
$$|e^{j\ast 2\pi \frac{d}{\lambda }\ast n\ast (sin{\theta }_1-sin{\theta }_2)}|<1,n=2,3,...,N-1$$
但是当波束指向方向与来波方向重合的时候，即${\theta }_1={\theta }_2=20^{\circ }$的时候
$$y=WX=Wa({\theta }_1)x_1(t)=(a({\theta }_2){)}^{T}a({\theta }_1)x_1(t)=(1+1+...+1)x_1(t)=N{x}_1(t)$$
其中$N$为阵元个数，即信号$x_1(t)$被增强$N$倍。

阵列方向图

我们可以通过仿真，绘制方向图，来直观展示阵列的空间滤波响应，方向图定义为阵列输入为平面波的时候，阵列输出（通常只是考虑输出幅度或者功率与平面波入射角的关系）。一维线阵在阵列基线上的幅度方向图为：
$$G(\theta )=|W^{H}a(\theta )|$$
功率方向图为：
$$G(\theta )=|W^{H}a(\theta ){|}^2$$
不管是幅度方向图还是功率方向图，都和阵列的权矢量有关，在工程中，分贝形式的方向图更为普遍，幅度方向图和功率方向图的分贝形式分别为：
$$G(\theta )=10lg|W^{H}a(\theta )|$$
功率方向图为：
$$G(\theta )=20lg||W^{H}a(\theta ){|}^2$$

阵列方向图仿真

下面采用matlab仿真，绘制方向图，编写的方向图自定义函数

function Gtheta = DBF(elementNum,observaTheta,pointTheta,dLmbda,plotEnableHigh)
% 该函数用于产生方向图
%  Author:huasir 2023.9.24 @Beijing
% Input : 
%   * elementNum: 阵元数
%   * observaTheta：观测角度范围，例如[-90:0.1:90]
%   * theta：来波方向
%   * dLembda：阵元间距与波长的比值，通常为1/2
% Output : 
%    * GTheta：方向图
N = elementNum; %阵元数
theta = observaTheta; %观测角度范围
theta0 = pointTheta; %波形形成的方向
theta0 = theta0*pi/180;
theta = theta*pi/180;
d_lembda = dLmbda; %阵元间距比波长
a = exp(j*2*pi*d_lembda*(0:N-1)'*sin(theta)); %导向矢量
w = exp(j*2*pi*d_lembda*(0:N-1)'*sin(theta0)); %加权权向量
y = abs(w'*a);
%win = taylorwin(N); % 泰勒窗
%win = chebwin(N,30); %切比雪夫窗，抑制旁瓣
% ww = w.*win; %加窗操作
% y1 = abs(ww'*a); %加窗后的方向图
if plotEnableHigh == 1
    figure;
    %plot(theta*180/pi,20*log10(y/max(y)),'--',theta*180/pi,20*log10(y1/max(y1)),'-');
    plot(theta*180/pi,20*log10(y/max(y)),'linewidth',1);
    xlabel('方位角/°');ylabel('归一化功率方向图/dB');
    axis tight;% axistight 使得图形框图靠近数据
    grid on; %添加栅格线
    ylim([-40, 0]); % 为了限制y值范围，使得图像显示的更加合理
    title(sprintf('阵元数：%d，波束方向：%.0f°',N,pointTheta))
end
Gtheta = y; %返回方向图
end
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35

在主函数中进行调用，主要参数如下：

主函数中进行调用：

close all; clear all; clc;
Gtheta = DBF(8,(-90:0.1:90),10,1/2,1);
1
2

绘制的方向图如下：

图1. 功率方向图