基于新循环算法的低旁瓣非周期单模序列设计

Design of Low Sidelobe Aperiodic Unimodular Sequences Based on New Cycling Algorithm

摘要: 针对雷达波形设计中旁瓣过高的问题，提出了一种基于快速傅里叶变换的新循环算法的单模非周期序列优化方法。首先，构建了以最小化积分旁瓣水平(integrated-sidelobe-level, ISL)为优化准则，以单模序列为约束条件的优化问题。由于该优化问题的非凸性，提出一种基于新循环算法的迭代优化算法对其进行求解。此外，还在理论上改进其算法的收敛性，确保它在迭代过程中能稳定地收敛到最优解或局部最优解。仿真实验结果表明，该方法有效降低非周期单模序列自相关旁瓣，提高单模非周期序列的优点因子。

Abstract: Because of the issue of high sidelobes in radar waveform design, we proposed a new optimization method for unimodular aperiodic sequences based on fast Fourier transform cyclic algorithm. Firstly, an optimization problem is constructed with the minimization of integrated sidelobe level (ISL) as the optimization criterion and single-mode sequences as constraint conditions. Due to the non-convex nature of the optimization problem, a new iterative optimization algorithm based on a new cycling algorithm is proposed to solve it. Furthermore, improvements are being made in the theoretical convergence of its algorithm, ensuring that it can stably converge to the optimal or local optimal solution during the iterative process. Simulation experiment results show that the method effectively reduces the autocorrelation sidelobes of single-mode non-periodic sequences, and improves the merit factor of single-mode non-periodic sequences.

1. 引言

雷达系统中存在干扰、强杂波、多路径等多种因素，信号在进行多目标检测时，会出现微弱目标被掩盖在较强的目标或者强杂波旁瓣中的情况 [1] 。所以雷达波形设计一般以降低旁瓣电平为主 [2] 。

在雷达波形设计中，典型的波形有线性调频波形、相位编码波形等，相位编码波形又分为多相编码序列和二相编码序列。多相编码序列也叫做单模序列 [3] 。常见的单模序列有Chu序列、Px序列、Frank序列、golomb序列等 [4] 。二相编码序列有m序列、Gold序列、Kasami序列等。对于非周期情况，两类序列都存在自相关旁瓣过高的缺点，但是单模序列Chu序列要比m序列的自相关旁瓣低，并且单模序列的相位自由度比m序列要高，设计性强。在进行雷达探测时，不易被侦察和识别，隐蔽性好。

在不同的主动传感应用中，需要设计具有良好自相关特性的序列。例如，在雷达和声纳系统中使用低自相关峰值旁瓣水平(peak-sidelobe-level, PSL)序列进行距离压缩，可以缓解强回波掩蔽弱目标的旁瓣问题 [5] [6] 。使用低自相关ISL序列可以减少接近感兴趣目标的分布式杂波返回的有害影响 [7] 。根据使用场景的不同，在许多领域都需要低周期或低非周期自相关旁瓣的序列。所以我们设计的目标是使序列的ISL或PSL尽可能低。

文献 [8] 通过代数构造，构造出巴克序列，它具有良好的非周期自相关性，但是长度不超过13，不适用于实际应用 [8] 。文献 [9] 中基于勒让德多项式和旋转构建了勒让德序列 [9] ，并且长度满足特定需求。然而，这些序列的长度并非任意可选，因为它们必须是质数。还有一些其他的构造序列，如鲁丁–夏皮罗序列 [10] 、雅克布序列 [11] 、Gold序列 [12] 。这些通过构造出来的低旁瓣序列，都存在着长度受限，多样性受限的问题。于是一部分学者通过算法来设计足够长度、多样性的良好非周期自相关序列。如文献 [13] 中通过穷举搜索最优的序列 [13] ，对于短长度的序列来说，搜索算法比较简便，但是搜索较长的序列，需要花费较多时间。文献 [14] 提出新循环算法即CAN算法，将现有的良好的单模序列进行初始化，经过算法迭代后，能优化该序列。因为使用快速傅里叶变换的方法来优化序列，所以计算效率高，能设计较长长度的序列，复杂度低，实现方便。但CAN算法收敛性不能保证，文献 [15] 中使用Golomb序列也存在局限性，长度只能为奇数。不能有效解决序列存在旁瓣过高，会造成微弱目标丢失的问题。

本文针对此问题，选取长度不受限制的单模非周期序列Chu序列作为优化对象，提出了基于快速傅里叶变换的新循环算法的单模非周期序列的方法，并改进了算法的收敛性，确保它在迭代过程中能稳定的收敛到最优解或局部最优解。该方法能有效降低单模非周期序列的旁瓣，提高优点因子以及弱目标的分辨能力。

2. 单模非周期序列

相位编码信号如下：

$s\left(t\right)={\displaystyle \underset{n=1}{\overset{N}{\sum }}x\left(n\right)\sqrt{\frac{1}{t_p}}rect\left(\frac{t-\left(n-1\right)t_p}{t_p}\right),n=1,\cdots ,N.}$ (1)

$\left\{x\left(n\right)\right\}$ 是相位编码序列， $t_p$ 是脉宽。实际传输的波形是由 $s\left(t\right){\text{e}}^{j2\pi f_{c}t}$ 的同相分量和正交分量组成的，其中 $f_c$ 为载波频率。假设信号解调已经在接收端完成，所以在分析中，我们就忽略 ${\text{e}}^{j2\pi f_{c}t}$ 。

为了避免非线性的副作用和最大化传输在系统中可用的功率，发射的相位编码序列为单模序列或者低峰均功率比的序列。

$x\left(n\right)=\exp \left(j\varphi \left(n\right)\right),n=1,\cdots ,N$ (2)

式(2)为相位编码序列表达式，对其施加单位模量约束，即 $\left|x\left(n\right)\right|=1$ 。满足这个条件的序列叫做单模序列。那么降低单模非周期序列的旁瓣的问题，实际上就是以最小化ISL为目标函数，单位模量为约束条件的优化问题。

$\varphi \left(n\right)$ 是相位，二相编码序列的相位只有0或 $\pi$ 两个取值，可以用二进制相位序列 $\left\{\varphi \left(n\right)=0,\pi \right\}$ 表示，也可以用二进制序列 $\left\{x\left(n\right)=+1,-1\right\}$ 来代替。如巴克码、m序列、Gold序列等；巴克码副瓣均匀，具有良好的自相关性，但是长度限制在13；Gold序列长度限制在 $2^{n-1}$ ，n为整数；m序列是周期自相关性比较好，但是缺乏多样性并且非周期自相关性不太好；多相编码序列的相位是 $\left[0,2\pi \right]$ ，相位范围更大，设计更具有自由度。如Frank序列、Px序列、Golomb序列、Chu序列等。

其中Chu序列表达式为

$x\left(n\right)=\{\begin{array}{l}\exp \left[jQ\pi \frac{n\left(n-1\right)}{N}\right],N为奇数,\\ \exp \left[jQ\pi \frac{n^2}{N}\right],N为偶数,\end{array}$ (3)

其中Q是与N互素的任意整数。Chu序列相比于Golomb序列和Frank序列，对于长度N的选择范围更广，所以在接下来的优化中，我们选择Chu序列作为优化对象。

3. 基于快速傅里叶变换的优化旁瓣算法

3.1. 优化目标函数

设 ${\left\{x\left(n\right)\right\}}_{n=1}^N$ 是长度为N的非周期离散序列，其自相关函数如下：

$r\left(k\right)={\displaystyle \underset{n=k+1}{\overset{N}{\sum }}x\left(n\right)x^{*}\left(n-k\right)}=r^{*}\left(-k\right),k=0,\cdots ,N-1.$ (4)

式(1)中， $r\left(0\right)$ 是零时延的自相关函数，等于信号能量N，我们也称作主瓣。除去在零时延处的自相关函数外，其余的自相关函数，即 ${\left\{r\left(k\right)\right\}}_{k=1}^{N-1}$ 称为自相关旁瓣。但是旁瓣峰值过高，在检测目标时会造成多目标检测的丢失，以及弱目标的检测性能差等问题。所以我们要尽可能的降低旁瓣，也就是使 ${\left\{r\left(k\right)\right\}}_{k=1}^{N-1}$ 尽可能小。

在本文中我们专注于最小化ISL，其定义如下：

$\text{ISL}={\displaystyle \underset{\begin{array}{c}k=-\left(N-1\right)\\ k\ne 0\end{array}}{\overset{N-1}{\sum }}{\left|r\left(k\right)\right|}^2}=2{\displaystyle \underset{k=1}{\overset{N-1}{\sum }}{\left|r\left(k\right)\right|}^2}.$ (5)

优点因子(merit factor, MF)与ISL成反比，定义

$\text{MF}=\frac{{\left|r\left(0\right)\right|}^2}{{\displaystyle \underset{\begin{array}{c}k=-\left(N-1\right)\\ k\ne 0\end{array}}{\overset{N-1}{\sum }}{\left|r\left(k\right)\right|}^2}}=\frac{N^2}{\text{ISL}}.$ (6)

根据式(5)，通过最小化ISL来设计具有良好自相关特性的序列。根据式(6)，最小化ISL就相当于最大化MF。

为了更好的研究非周期序列的自相关性，我们用归一化非周期自相关函数来表示

$\text{NAAF}\left(\text{dB}\right)=20{\log }_{10}\left|\frac{r\left(k\right)}{r\left(0\right)}\right|\text{,}k=1,\cdots ,N-1.$ (7)

3.2. 优化旁瓣算法

CAN算法是专门针对单模非周期序列的ISL问题所设计的算法 [14] 。其核心思想是基于快速傅里叶变换。

首先将时域上的自相关函数转化成频域上的表达，令

$\begin{array}{c}\Phi \left(\omega \right)\triangleq {\left|{\displaystyle \underset{n}{\overset{N}{\sum }}x\left(n\right){\text{e}}^{-j\omega n}}\right|}^2\\ ={\displaystyle \underset{k=-\left(N-1\right)}{\overset{N-1}{\sum }}r\left(k\right){\text{e}}^{-j\omega k}},\omega \in \left[0,2\pi \right]\end{array}$ (8)

那么ISL可以写成 [15] ：

$\text{ISL}=\frac{1}{2N}{\displaystyle \underset{p=1}{\overset{2N}{\sum }}\left[\Phi \left({\omega }_p\right)-N\right]}$ (9)

其中

${\omega }_p=\frac{\pi }{N}p,p=1,\cdots ,2N.$ (10)

则ISL的优化问题可以表述为：

$\begin{array}{l}\underset{{\left\{x\left(n\right)\right\}}_{n=1}^N}{\min }{\displaystyle \underset{p=1}{\overset{2N}{\sum }}{\left({\left|{\displaystyle \underset{n=1}{\overset{N}{\sum }}x}\left(n\right){\text{e}}^{-j{\omega }_{p}n}\right|}^2-N\right)}^2}\\ \text{s}\text{.t}.\left|x\left(n\right)\right|=1,n=1,\cdots ,N.\end{array}$ (11)

式(11)又近似等价于以下二次优化问题：

$\begin{array}{l}\underset{{\left\{x\left(n\right)\right\}}_{n=1}^N;{\left\{{\psi }_p\right\}}_{p=1}^{2N}}{\min }{\displaystyle \underset{p=1}{\overset{2N}{\sum }}{\left|{\displaystyle \underset{n=1}{\overset{N}{\sum }}x}\left(n\right){\text{e}}^{-j{\omega }_{p}n}-\sqrt{N}{\text{e}}^{j{\psi }_p}\right|}^2}\\ \text{s}\text{.t}.{\psi }_p\in \left[0,2\pi \right],p=1,\cdots ,2N,\\ \left|x\left(n\right)\right|=1,n=1,\cdots ,N.\end{array}$ (12)

为了简化目标函数，将式(12)表示成矩阵形式：

$\begin{array}{l}\underset{{\left\{x\left(n\right)\right\}}_{n=1}^N;{\left\{{\psi }_p\right\}}_{p=1}^{2N}}{\min }\Vert A^{H}z-v\Vert \\ \text{s}\text{.t}\text{.}{\psi }_p\in \left[0,2\pi \right],p=1,\cdots ,2N,\\ \left|x\left(n\right)\right|=1,n=1,\cdots ,N.\end{array}$ (13)

所以最小化ISL问题就转化成当目标函数 $\Vert A^{H}z-v\Vert$ 最小时，求 ${\left\{x\left(n\right)\right\}}_{n=1}^N$ 的值的二次函数优化问题。此时的序列 ${\left\{x\left(n\right)\right\}}_{n=1}^N$ 就是优化后的序列。

其中

${\alpha }_p^H=\left[{\text{e}}^{-j{\omega }_p}\cdots {\text{e}}^{-j2N{\omega }_p}\right],$ (14)

$A^H=\frac{1}{\sqrt{2N}}\left[\begin{array}{l}{\alpha }_1^H\\ ⋮\\ {\alpha }_{2N}^H\end{array}\right],$ (15)

$z={\left[x\left(1\right)\cdots x\left(N\right)0\cdots 0\right]}_{2N\times 1}^{\text{T}},$ (16)

$v=\frac{1}{\sqrt{2}}\left[{\text{e}}^{j{\psi }_1}\cdots {\text{e}}^{j{\psi }_{2N}}\right].$ (17)

$A^H$ 是 $2N\times 2N$ 的离散傅里叶变换半酉矩阵。令 $f=A^{H}z$ ，则 $f$ 是 $z$ 的快速傅里叶变换矩阵。 $\psi$ 为 $f$ 的相位角，所以有 ${\psi }_p=\arg \left(f_p\right),p=1,\cdots ,2N$ ，计算出 ${\left\{{\psi }_p\right\}}_{p=1}^{2N}$ 的最小值。令 $g=Av$ ，则 $g$ 是 $v$ 的逆快速傅里叶变换。因为 ${\Vert A^{H}z-v\Vert }^2={\Vert z-Av\Vert }^2={\Vert z-g\Vert }^2$ ，所以有 $x\left(n\right)={\text{e}}^{j\arg \left(g_n\right)},n=1,\cdots ,N$ ，可以通过该式计算序列的最小值 ${\left\{x\left(n\right)\right\}}_{n=1}^N$ 。

3.3. 改进后的优化旁瓣算法

CAN算法是通过交替优化 ${\left\{x\left(n\right)\right\}}_{n=1}^N$ 和 $v$ 来求解式(13)。文献 [15] 中结束迭代的条件是满足 $\Vert X^{i+1}-X^i\Vert <\xi$ ， $X^i$ 是第i次迭代后生成的序列， $\xi$ 是一个事先设置好的阈值。如10⁻³。然而，优化变量 $X^i$ 的收敛性在理论上不能保证。我们的目标是最小化ISL，目标函数是 $\Vert A^{H}z-v\Vert$ 。如果使用 $\Vert X^{i+1}-X^i\Vert <\xi$ 作为终止标准，则目标函数可能在满足所述终止标准之前就已经停止改变。所以我们将目标函数作为迭代结束的条件。接下来是证明其收敛性。

令 $g=Av$ ， $z^i={\left[x^i\left(1\right)\cdots x^i\left(N\right)0\cdots 0\right]}_{2N\times 1}^{\text{T}}$ ， $D_i$ 表示在i次迭代时 ${\Vert A^H{z}^i-v^i\Vert }^2$ 的值，我们有

$\begin{array}{c}{D}_i={\Vert A^H{z}^i-v^i\Vert }^2\ge {\Vert A^H{z}^i-v^{i+1}\Vert }^2\\ \ge {\Vert A^H{z}^{i+1}-v^{i+1}\Vert }^2=D_{i+1}.\end{array}$ (18)

计算 $D_i$ 就是计算 $\Vert z^i-g^i\Vert$ 。由于 $g^i$ 在第i次迭代中已经储存。 $D_i$ 是一个二范数，计算复杂度与 $\Vert X^{i+1}-X^i\Vert$ 一样。以 $\left|D_i-D_{i-1}\right|<\xi$ 作为终止准则，符合优化问题的目标(即降低目标函数的值)，理论上保证了 $D_i$ 的单调减小并且计算复杂度保持与原CAN算法不变。改进后的CAN算法步骤如下：

4. 数值结果与分析

4.1. 初始序列的选择

从图1可以看出，Chu序列的自相关函数旁瓣总体比m序列的自相关旁瓣要低很多，但是未被优化的Chu序列两端的旁瓣峰值距离主瓣的旁瓣也很高，很容易导致在目标检测时，相邻的目标不能被分辨以及小目标被掩盖在大目标下的问题 [3] 。所以我们选择对Chu序列进行优化，来降低其自相关旁瓣。

Figure 1. The normalized aperiodic auto-correlation function

图1. 归一化非周期自相关函数

4.2. 对Chu序列用改进后的CAN算法进行优化

图2是长度为127的Chu序列和经过CAN算法优化后的Chu序列的NAAF，根据前文的公式(4)，当 $k=0$ 时， $r\left(0\right)=N$ ， $\text{NAAF}=20{\log }_{10}1=0$ 。从图中可以看出k在零处，对应为主瓣，相应的剩下的部分为旁瓣。当k接近0和 $N-1$ 时，Chu序列的自相关旁瓣比起其他旁瓣都要高，称作峰值旁瓣。经过CAN算法优化后的Chu序列在接近0和 $N-1$ 时，其峰值旁瓣有了很明显的降低。Chu序列的峰值旁瓣为−27.8 dB，经过优化后的峰值旁瓣为−45.6 dB。从整体来看，经过CAN算法优化后的Chu序列的自相关旁瓣比Chu的自相关旁瓣都很大程度的降低了。有些地方的旁瓣都达到了−80 dB以下。

Figure 2. The NAAF of the Chu and CAN(C) sequences

图2. 优化序列后的归一化非周期自相关函数

验证了改进后的CAN算法能够优化已有良好非周期单模序列的旁瓣。接下来用MF数值来进一步验证改进后的CAN算法能降低其旁瓣。

4.3. 优点因子比较

从图3可以很清晰的看出经过CAN算法优化Chu序列的优点因子比m序列、Chu序列、随机相位序列都要高。所以很容易观察到，CAN(C)序列比m序列或随机相位序列提供了更低的相关性旁瓣。随着长度的增加，优点因子也在不断的增加。在长度为 $2^5-1$ 的Chu序列和CAN(C)序列的优点因子非常接近。通过表1，当长度为 $2^{13}-1$ 时，CAN(C)的优点因子为1431.5，而未经优化的Chu序列的优点因子为142.15。在优点因子方面，CAN算法可以将优点因子提升了将近十倍。

Figure 3. The merit factors of four types of sequences

图3. 四种序列的优点因子

Table 1. Merit factor

表1. 优点因子

4.4. 模糊函数对比

前文所提出的CAN算法能够降低Chu序列的自相关旁瓣，解决被旁瓣掩盖的弱目标问题。但是只考虑了固定目标，在实践中更多的是移动目标，所以在进行波形设计时，还要考虑多普勒频移问题。模糊函数是时延和多普勒频移的相关函数，其零多普勒切割就是序列的自相关函数。所以模糊函数可以看作是序列自相关函数的二维拓展。

(a) (b)

Figure 4. Ambiguity function. (a) The AF of a length-127 Chu sequence; (b) The AF of a length-127 CAN(C) sequence

图4. 模糊函数。(a) 长度127的Chu序列模糊函数；(b)长度127的CAN(C)序列模糊函数

(a) (b)

Figure 5. 2D-Ambiguity function. (a) The AF of a length-127 Chu sequence; (b) The AF of a length-127 CAN(C) sequence

图5. 2D模糊函数。(a) Chu序列的2D模糊函数；(b) CAN(C)序列的2D模糊函数

图4(a)是长度为127的Chu序列的模糊函数图，图4(b)是相同长度下用Chu序列进行优化后的序列模糊函数图。可以看出呈斜刃状，周围的相关旁瓣都有所降低。图5是2D的模糊函数图，用颜色的深浅来表示能量的高低。从图5中可以看出相比于未进行优化前的Chu序列，无论是主瓣还是旁瓣，它们的能量都较之前更为分散了，由密集的斜刃峰值变成了稀疏的斜刃峰值，也在一定程度上降低了旁瓣，验证了CAN算法的有效性。

5. 结论

提出了基于快速傅里叶变换的算法优化方法，并改进了算法的收敛性。改进后能优化现有的单模序列Chu序列，通过数值实验，发现经过优化后的Chu序列能够获得更高的优点因子，达到了优化的目的。最后的模糊函数图也进一步表明，由该算法优化后的序列，能获得较低的旁瓣，在一体化波形设计时，可以考虑使用其优化后的序列。