基于压缩采样阵列的空频二维谱估计方法

摘要

本发明公开了一种基于压缩采样阵列的空频二维谱估计方法，主要解决在大采样率和多通道时，由于数据量非常惊人，使得数据存储变得异常困难的问题。其实现步骤是：1)将用Nl个天线的接收机采集到的信号作为输入信号矩阵X；2)根据输入信号矩阵X生成观测矩阵Y；3)根据观测矩阵Y构造频域基矩阵F、设置时域压缩矩阵Φb和联合稀疏矩阵Zf，并建立信号模型YT＝ΦbFZf；4)求解上述信号模型得到联合稀疏矩阵Zf，根据联合稀疏矩阵得到非零支撑集合Γ；5)根据非零支撑集合Γ解得空频二维谱Zp。本发明降低了系统的通道数和时域采样速率，提高了瞬时信号接收带宽和测向精度，可用于雷达、侦察一体化系统。

说明书

技术领域

本发明属于信号处理技术领域，特别涉及一种空频二维谱估计方法，可用于雷达、 侦察一体化系统。

背景技术

在雷达、侦察一体化系统中，天线阵列接收到的信号，分别进入各个通道的射频 前端，而射频前端主要完成低噪放大、可变增益控制、混频、中频放大等功能，由于 系统有测向功能，射频前端需要有较好的相位一致性。

中频信号从射频前端进入数字信号处理分系统，首先，经过模数转换A/D变换， 再送给现场可编程门阵列FPGA，在现场可编程门阵列FPGA中进行高速缓存和预处 理；由于系统带宽要达到1GHz，因此先对中频信号进行信道化，对信道化后的信号进 行自相关，求得信号的幅度；再与自适应门限比较进行信号检测，进而将非弱信号段 提取出来，送给数字信号处理器DSP进行波达方向DOA估计；然后，由数字信号处 理器DSP将计算的测向结果和用于波束形成的权值传输给现场可编程门阵列FPGA， 在现场可编程门阵列FPGA中进行数字波束合成；最后，进行传统的信号处理，包括 信号的分选、信号脉内特性分析等。

在传统的阵列信号处理系统框架中，每个通道都包含了低噪声放大、增益控制、 下变频、带通滤波、模数转换ADC等电路，为了提高波达方向DOA估计的精度，通 道的数量会不断的增大，随着通道的增加，模拟部分的电路设备量会急剧增大。系统 通常具有很高的采样率和很多的通道数，这就需要使用大量的高速模数转换ADC和更 大规模的现场可编程门阵列FPGA，从而导致系统的功耗、体积急剧增大，并且在采 样阶段出现瓶颈，尤其是在大采样率和多通道时，由于数据量非常惊人，使得数据存 储也变得异常困难。

发明内容

本发明的目的在于针对上述已有技术的不足，提出基于基于压缩采样阵列的宽带 空频二维谱估计方法，以降低采样速率和通道的数量，减小数字信号的数据量。

本发明的技术关键是同时对空域-频域进行压缩，其实现步骤包括如下：

1、一种基于压缩采样阵列的空频二维谱估计方法，包括以下步骤：

1)采用具有N_l个天线的接收机采集到的N_l×N_t维信号作为输入信号矩阵X，定 义每个天线为一个阵元，记为i；设阵元i到阵元1的间距为d_i，并假设有R个信号同 时入射到该天线接收机形成的随机线性阵列上，其中，i＝1,2,…,N_l，N_t是时间总长度， 1≤R≤N_l；

2)将输入信号矩阵X依次经过空域压缩、模拟信息转换和模数转换，生成观测矩 阵Y；

3)构造频域基矩阵F：

3a)将载频搜索范围分为N_f份，每一份用表示为：

$F_{f_q}={\left(\begin{array}{cccc}1& e^{j{\omega }_q}& ...& e^{j{\omega }_q(N_t-1)}\end{array}\right)}^T/\sqrt{N_t}$

其中，是信号频率f_q的一个傅里叶基向量，q＝1,2,…,N_f；ω_q是归一化信号角频率， ω_q＝2π(f_q/F_s)，F_s是奈奎斯特采样频率，且N_t＝TF_s，T是时间观察窗，j表示虚数， [·]^Τ表示向量的转置；

3b)在频域设置一个大小为N_t×N_f维的傅里叶基矩阵F：

$F=[F_{f_1}...F_{f_q}...F_{f_{N_f}}],N_f\ge N_t$

其中，$f_q=f_1,f_2,...,f_{N_f},N_f\ge 1;$

4)信号模型的建立及求解：

4a)在时域范围内设置一个服从高斯随机分布的M_t×N_t维时域压缩矩阵Φ_b，其 中，Φ_b表示以M_t/N_t的奈奎斯特采样率的模拟信息转换器AIC采样，M_t是输入信号 矩阵X经空域压缩后的列数，且M_t＜N_t；

4b)设Z_f是大小为N_f×M_l的联合稀疏矩阵：$Z_f=[z\left(f_1\right)...z\left(f_l\right)...z\left(f_{M_l}\right)]$其中，z(f_l)对应一个信号的频域表示，f_l是经过以M_t/N_t为奈奎斯特的AIC模拟信息 转换器采样后的频率，且

4c)在频域范围内，依据宽带空频二维压缩采样矩阵得到信号模型：

Y^T＝Φ_bFZ_f；

4d)将联合稀疏矩阵Z_f的恢复问题转化为求解如下最小范数矩阵方程：

$\underset{Z_f}{\min }\left(\begin{array}{ccc}{\left|\right|Z_f\left|\right|}_{\mathrm{2,1}}& s.t.& Y^T\end{array}\right)={\Phi }_{b}F{Z}_f$

其中||·||_2,1表示矩阵的范数；

4e)求解上述最小范数矩阵方程，得到联合稀疏矩阵Z_f；

4f)对联合稀疏矩阵Z_f进行压缩重构，得到Z_f的非零支撑集合Γ，其中，f_s∈Γ， f_s＝f₁,f₂…f_Γ，f_s是经过以M_t/N_t为奈奎斯特的AIC模拟信息转换器采样后的非零频 率；

；

5)根据非零支撑集合Γ解得空频二维谱Z_p：

5a)在空域设置一个服从高斯随机分布的M_l×N_l维空域压缩矩阵Φ_a，且 M_l＜N_l；

5b)对于联合稀疏矩阵Z_f的某个幅度非零的信号频率f_s，定义该频率信号的波 达方向DOA谱为z_θ(f_s)，构造角度域基矩阵Α(f_s)，则该信号的频域值z(f_s)的空域稀 疏表示为：

z(f_s)＝Φ_aΑ(f_s)z_θ(f_s)；

5c)将空频二维谱估计问题转化为求解如下最小范数矩阵方程：

$\underset{z_{\theta }\left(f_s\right)}{\min }\left(\begin{array}{ccc}{\left|\right|z_{\theta }\left(f_s\right)\left|\right|}_1& s.t.& z\left(f_s\right)={\Phi }_a\end{array}A\left(f_s\right)z_{\theta }\left(f_s\right)\right)$

其中，||·||₁表示矩阵的范数；

5d)求解上述最小范数矩阵方程，得到波达方向DOA谱z_θ(f_s)；

5e)对所有非零支撑集合Γ求波达方向DOA谱，得到空频二维谱Z_p:

Z_p＝[z_θ(f₁)…z_θ(f_s)…z_θ(f_Γ)]。

本发明与现有技术相比具有以下优点：

1)本发明由于在时域引入时域压缩矩阵Φ_b，与传统的模拟数字转换相比，降低 了系统采样率，同时使系统获得了非常大的瞬时信号接收带宽；

2)本发明由于在空域引入空域压缩矩阵Φ_a，与传统的滤波法相比，降低了系统通 道数，进而降低了硬件设计的复杂度与器件的要求，从而降低了系统功耗，同时减小 了系统体积；

3)本发明采用随机天线阵列模型，同时在空域引入空域压缩矩阵Φ_a，因此使系 统有更好的测向精度，从而为射频前端提供了相对较好的相位一致性；

4)本发明与其他的压缩感知的方法相比，不仅利用了数据在频域上的稀疏性同时 也利用了信号在空域上的稀疏性，使数据的得到最大程度的压缩，在相同总压缩率下， 本发明的系统重构性优于传统方法。

附图说明

图1是本发明的实现流程图；

图2是本发明的实现场景图；

图3是本发明中生成观测矩阵的示意图；

图4是本的发明的仿真结果图。

具体实施方式

以下参照附图，对本发明的技术方案和效果作进一步的详细说明。

参照图1，本发明的实现步骤如下：

步骤1：采用具有N_l个天线的接收机采集到的信号作为输入信号矩阵X。

如图2所示，本发明的实现场景包括具有N_l个天线的接收机，R个入射信号。

所有天线采集到的信号是一个N_l×N_t的输入信号矩阵X，定义每个天线为一个阵 元，记为i，设阵元i到阵元1的间距为d_i，并假设有R个信号同时入射到该N_l个天线 构成的阵列，入射信号的波达方向DOA角度为θ_k，其中，k＝0,1,2,…,R-1，1≤R≤N_l， i＝1,2,…,N_l，N_t是时间总长度；

输入信号矩阵X表示如下：

其中，x(i,t)表示第i个阵元在t时刻接收到的复基带射频信号，

$x(i,t)=\underset{k=0}{\overset{R-1}{\Sigma }}{\beta }_k{e}^{j2\pi (f_k/F_s)t}{e}^{-j2\pi (d_i/\lambda )\sin {\theta }_k}+n(i,t),$

其中，k表示同时入射到具有N_l个天线的接收机的第k个信号，k＝0,1,…,R-1，β_k是 复振幅，θ_k是波达方向DOA，λ是入射信号波长，f_k表示同时入射到输入信号矩阵X 的第k个信号的频率，f_k＝f₀,f₁,…,f_R-1，n(i,t)是加性噪声，t＝1,2,…,N_t。

步骤2：根据输入信号矩阵X生成观测矩阵Y。

参照图3，本步骤的具体实现如下：

2a)用一个服从高斯随机分布的M_l×N_l维的空域压缩矩阵Φ_a对输入信号矩阵X 进行空域压缩，使输入信号矩阵X的行数由N_l行减少到M_l行，即N_l＞M_l≥1，以减少 系统的通道数，得到M_l×N_t维空域压缩信号；

2b)对M_l×N_t维空域压缩信号通过以M_t/N_t为奈奎斯特采样率的模拟信息转换 器AIC进行转换，使输入信号矩阵X的列数由N_t列减少到M_t列，即N_t＞M_t≥1，以 降低系统的时域采样速率，使系统获得非常大的瞬时接收带宽，得到M_l×M_t维空频域 压缩模拟采样信号；

2c)对M_l×M_t维空频域压缩采样信号通过模数转换器A/D进行模数转换，生成 观测矩阵Y：

其中，y(m,n)表示接收到入射信号的第m个阵元在n时刻接收到的数字信号， m＝1,2,…,M_l，n＝1,2,…,M_t，M_l是输入信号矩阵X经过空域压缩后的行数，M_t是输 入信号矩阵X经过模拟信息转换后的列数。

步骤3：构造频域基矩阵F。

3a)将载频搜索范围分为N_f份，每一份用F_fq表示为：

$F_{f_q}={\left(\begin{array}{cccc}1& e^{j{\omega }_q}& ...& e^{j{\omega }_q(N_t-1)}\end{array}\right)}^T/\sqrt{N_t}$

其中，是信号频率f_q的一个傅里叶基向量，q＝1,2,…,N_f；ω_q是归一化信号角频率， ω_q＝2π(f_q/F_s)，F_s是奈奎斯特采样频率，且N_t＝TF_s，T是时间观察窗，j表示虚数， [·]^Τ表示向量的转置，N_f是大于等于1的整数；

3b)在频域设置一个大小为N_t×N_f维的傅里叶基矩阵F：

$F=[F_{f_1}...F_{f_q}...F_{f_{N_f}}],N_f\ge N_t$

其中，$f_q=f_1,f_2,...,f_{N_f}.$

步骤4：信号模型的建立及求解。

4a)在频域范围内设置一个服从高斯随机分布的M_t×N_t维时域压缩矩阵Φ_b，其 中，Φ_b表示以M_t/N_t的奈奎斯特采样率的模拟信息转换器AIC采样，M_t是输入信号 矩阵X经过模拟信息转换器AIC转换后的列数，且N_t＞M_t≥1；

4b)设Z_f是大小为N_f×M_l的联合稀疏矩阵：$Z_f=[z\left(f_1\right)...z\left(f_l\right)...z\left(f_{M_l}\right)]$其中，z(f_l)对应一个信号的频域表示，f_l是经过以M_t/N_t为奈奎斯特的AIC模拟信息 转换器采样后的频率，

4c)依据宽带空频二维压缩采样矩阵得到信号模型：Y^T＝Φ_bFZ_f，求得联合稀 疏矩阵Z_f关于观测矩阵Y的表达式，即将观测矩阵Y由时域转换到频域，其中，联合 稀疏矩阵Z_f是信号的频域表示；

4d)将联合稀疏矩阵Z_f的求解问题转化为求解如下最小范数矩阵方程，以 得到观测矩阵Y的频域表示联合稀疏矩阵Z_f：

$\underset{Z_f}{\min }\left(\begin{array}{ccc}{\left|\right|Z_f\left|\right|}_{\mathrm{2,1}}& s.t.& Y^T\end{array}\right)={\Phi }_{b}F{Z}_f$

其中||·||_2,1表示矩阵的范数；

4e)求解上述最小范数矩阵方程，得到联合稀疏矩阵Z_f；

4f)对联合稀疏矩阵Z_f进行压缩重构，得到Z_f的非零支撑集合Γ，其中，f_s∈Γ， f_s＝f₁,f₂…f_Γ，f_s是经过以M_t/N_t为奈奎斯特的AIC模拟信息转换器采样后的非零频 率。

步骤5：根据非零支撑集合Γ解得空频二维谱Z_p。

5a)在空域设置一个服从高斯随机分布的M_l×N_l维空域压缩矩阵Φ_a，且 M_l＜N_l，其中，空域压缩矩阵Φ_a的作用是将输入信号矩阵X的行数由N_l行减少到M_l行，以减少系统的通道数；

5b)构造角度域基矩阵Α(f_s)：

5b1)将波达方向DOA的搜索范围分为N_θ份，每一份记为θ_p，用α(f_s,θ_p)表 示频率为f_s、波达方向DOA为θ_p的信号阵列导向矢量：

$\alpha (f_s,{\theta }_p)={\left(\begin{array}{cccc}1& e^{-j2\pi f_s{d}_2\sin {\theta }_p/c}& ...& e^{-j2\pi f_s{d}_{N_l}\sin {\theta }_p/c}\end{array}\right)}^T,$

其中，p＝1,2,…,N_θ，c为光速，c＝3.0×10⁸m/s；

5b2)用所有波达方向DOA的信号阵列导向矢量构成角度域基矩阵A(f_s)：

$A\left(f_s\right)={[\alpha (f_s,{\theta }_1)\alpha (f_s,{\theta }_2)...\alpha (f_s,{\theta }_{N_{\theta }})]}_{N_t\times N_{\theta }}$

其中，N_θ＞＞N_l；

5c)对于联合稀疏矩阵Z_f的某个幅度非零的信号频率f_s，定义该频率信号的波 达方向DOA谱为z_θ(f_s)，则该信号的频域值z(f_s)的空域稀疏表示为： z(f_s)＝Φ_aΑ(f_s)z_θ(f_s)；

5d)将波达方向DOA谱z_θ(f_s)的求解问题转化为求解如下最小范数矩阵方 程：

$\underset{z_{\theta }\left(f_s\right)}{\min }\left(\begin{array}{ccc}{\left|\right|z_{\theta }\left(f_s\right)\left|\right|}_1& s.t.& z\left(f_s\right)={\Phi }_a\end{array}A\left(f_s\right)z_{\theta }\left(f_s\right)\right)$

其中，||·||₁表示矩阵的范数；

5e)求解上述最小范数矩阵方程，得到波达方向DOA谱z_θ(f_s)；

5f)对所有非零支撑集合Γ求波达方向DOA谱，得到空频二维谱Z_p:

Z_p＝[z_θ(f₁)…z_θ(f_s)…z_θ(f_Γ)]。

本发明的效果可通过以下仿真说明：

1.仿真条件：

采用具有10个天线的接收机形成的随机线性阵列，每个天线为一个阵元，记为i， 阵元i到阵元1的间距为d_i，观测空域角度范围为[-80°,80°]，入射信号波长为λ，频 域范围为[-0.5F_s,0.5F_s]，入射信号的频率f和波达方向DOA角度分别为(0，0°)、 (-0.15，5°)和(-0.25，20°)，信号在接收信道中加入了均值为0的高斯白噪声， 其中，i＝1,2,...,10，F_s是奈奎斯特采样率。

2.仿真内容与结果：

采用本发明对同时入射到上述随机线性阵列的射频信号的频率f和波达方向 DOA角度θ进行估计，结果如图4所示，其中图4中的纵坐标表示入射信号的频率f， 横坐标表示入射信号的波达方向DOA角度，白色区域表示该区域没有信号入射，且只 有3个区域有值，而其他区域均没有入射信号，这三个有值区域的坐标分别为(0.02, 0°)、(-0.13,4°)和(-0.24,19°)。

从图4中可以看出：1)对于同时入射到随机线性阵列的射频信号，在经过空域压 缩和时域压缩矩阵后，使系统的通道数减少和时域采样速率降低的同时，还可以估算 出入射信号的频率f和波达方向DOA角度；2)本发明估算出的入射信号的波达方向 DOA角度的误差为仅为1°，提高了测向精度，能为射频前端提供相对较好的相位一致 性。