线性调频(linear frequency modulated，LFM)信号由于易实现大时、宽带、宽积、低截获等优点，广泛应用于雷达、声呐与通信工程等领域。在电子战中，对LFM信号的参数估计是一项重要课题，因此成为研究人员的研究热点。目前常用的LFM信号参数估计方法主要包括时频分析法^[1-2]、CFCR域分析方法^[3-5]以及分数阶变换法^[6-10]。时频变换法受到窗函数的限制，很难同时满足时域和频域分辨率需求；分数阶变换法为了满足精度，需要进行大量的匹配搜索计算，因而，计算复杂度比较高，难于进行工程实现；CFCR域分析方法由于信号采样点数的限制，部分情况难以满足估计精度。本文提出一种LVD与Zoom-FRFT联合的多LFM信号参数估计方法。该方法通过LVD提供的粗估计信息减少了FRFT计算复杂度，利用FRFT提高了LVD的估计精度，并且可在较低的信噪比下对LFM信号进行有效的参数估计。

1 LVD和Zoom-FRFT原理 1.1 LVD原理

LVD是一种新型的时频分析方法^[5]，使LFM信号在中心频率-调频率域具有良好的能量聚集性，采用尺度变换原理，使多LFM信号的交叉项问题得到解决，能够有效地展现出多LFM信号的分量个数，每个分量都可以通过峰值检测过程被找到，其频率参数可以直接从坐标值中读取。

为了更好地解释LVD的工作原理，我们将连续无限时间多分量LFM信号设为

$s(t) = \sum\limits_{i = 1}^K {{A_i}\exp \left[ {{\rm{j}}2{\rm{{\text{π}} }}\left( {{f_i}t + \frac{1}{2}{\gamma _i}{t^2}} \right)} \right]} $

式中：t∈(−T_ob/2,T_ob/2),T_ob为观测时长；A_i为幅度；f_i为信号分量的中心频率；γ_i为信号分量的调频斜率；K为信号分量的数目。

其参数对称自相关函数定义为

$\begin{split} R_{\rm{s}}^c(t,\tau ) = s\left( {t + \dfrac{{\tau + a}}{2}} \right){s^ * }\left( {t - \dfrac{{\tau + a}}{2}} \right){\rm{ = }} \\ \displaystyle\sum\limits_{i = 1}^K {R_{{{\rm{s}}_i}}^c(t,\tau )} + \displaystyle\sum\limits_{i = 1}^{K - 1} {\displaystyle\sum\limits_{j = i + 1}^K {\left( {R_{{{\rm{s}}_i}{{\rm{s}}_j}}^c(t,\tau ) + R_{{{\rm{s}}_j}{{\rm{s}}_i}}^c(t,\tau )} \right)} } {\rm{ = }} \\ \displaystyle\sum\limits_{i = 1}^K {A_i^2\exp ({\rm{j2{\text{π}} }}{f_i}(\tau + a) + } {\rm{2{\text{π}} }}{\gamma _i}(\tau + a)t) + {\rm{ }} \\ \displaystyle\sum\limits_{i = 1}^{K - 1} {\displaystyle\sum\limits_{j = i + 1}^K {\left( {R_{{{\rm{s}}_i}{{\rm{s}}_j}}^c(t,\tau ) + R_{{{\rm{s}}_j}{{\rm{s}}_i}}^c(t,\tau )} \right)} } \end{split}$

(1)

式中：τ为时延变量；a为一个与尺度变换相关的时延常量； $R_{{{\rm{s}}_i}}^c(t,\tau )$ 和 $R_{{{\rm{s}}_i}{{\rm{s}}_j}}^c$ 分别为信号项分量与交叉项分量。

由式(1)可以很容易地看出，时间变量t和时延变量τ在信号项的相位项上相互耦合。二维线性相位分量中时间变量t和时延变量τ的耦合是线性频率偏移导致线性调频信号在二维参数平面上模糊表示的根本机制。换言之，如果能够解除时间变量t和时延变量τ的耦合，LFM信号的每个分量将在中心频率−调频斜率域聚焦为峰值。

为了解除类似的耦合，借助Keystone变换的思想，设Γ为尺度变换操作，对时间变量t进行尺度变换，令

$ t=t_{n} /[h(\tau+a)] $

式中：t_n为尺度时间；h为尺度因子，则对参数对称自相关函数的尺度变换结果如下：

$\begin{split} \varGamma (R_{\rm{s}}^c(t,\tau ){\rm{ ) = }}&\displaystyle\sum\limits_{i = 1}^K {A_i^2\exp \Bigg({\rm{j}}2{\rm{{\text{π}} }}{f_i}(\tau + a) + } 2{\rm{{\text{π}} }}{\gamma _i}\frac{{{t_n}}}{h}\Bigg) + \\& \displaystyle\sum\limits_{i = 1}^{K - 1} {\displaystyle\sum\limits_{j = i + 1}^K {\left( {R_{{{\rm{s}}_i}{{\rm{s}}_j}}^c({t_n},\tau ) + R_{{{\rm{s}}_j}{{\rm{s}}_i}}^c({t_n},\tau )} \right)} } \end{split}$

(2)

由式(2)可知，因为信号项引入了新的时间变量t_n，消除了时间变量t和时延变量τ的耦合。再分别沿尺度时间t_n和时延变量τ的方向做FFT运算可得：

$\begin{array}{c} L(\gamma ,f) = \displaystyle\sum\limits_{i = 1}^K {A_i^2\exp ({\rm{j}}2{\rm{{\text{π}} }}af)} \delta (f - {f_i}) \cdot \\ {\rm{ }}\delta (\gamma - {\gamma _i}/h) + {L_{{\rm{cross}}}}(\gamma ,f) \end{array} $

各LFM信号分量在CFCR域聚焦为峰值，根据峰值下标可获得对各分量的参数估计，图1为包含3个分量的LFM信号的LVD谱。

LVD的核心步骤是尺度变换。为了快速地完成离散的尺度变换过程，需要借助尺度傅里叶−快速傅里叶反变换(scaling Fourier transform-inverse fast Fourier transform，SFT-IFFT)来进行。设多分量LFM信号采样间隔为T_s，采样点数为N，则式(1)中信号项分量 $R_{\rm{s}}^c(t,\tau )$ 的离散形式为

$\begin{array}{c} R_{\rm{s}}^c(n,m) = \displaystyle\sum\limits_{i = 1}^K {A_i^2\exp ({\rm{j2{\text{π}} }}{f_i}(2m{T_{\rm{s}}} + a) + } \\ {\rm{ 2{\text{π}} }}{k_i}(2m{T_{\rm{s}}} + a)n{T_{\rm{s}}}) \end{array} $

尺度傅里叶变换(SFT)可以表示为

$S(\eta ,f) = SF(s(t)) = \int {s(t)\exp ( - {\rm{j2{\text{π}} }}\eta ft){\rm{d}}t} $

SFT能够通过chirp变换快速实现^[4]，如图2所示。

基于SFT-IFFT的尺度变换离散流程如图3所示。

图3中，n=[−N/2],[−N/2+1],…,[(N−1)/2]，m=[−N/2], [−N/2+1],…,[(N−1)/2]，[·]表示取整运算，η=(2mT_s+a)h/ N， $R_{\rm{s}}^c(n',m)$ 是 $R_{\rm{s}}^c(n,m)$ 经过尺度变换操作后的结果。

1.2 Zoom-FRFT原理

Zoom-FRFT^[11]是一种FRFT的高分辨计算方法，它解决了FRFT的固定分辨率问题，可以提供分数阶域上的任意段局部谱的细节。

假定要计算FRFT在局部谱区间[u₁,u₂]上的M点等间隔取样值(u₁、u₂和M的取值任意)，将分数阶傅里叶域变量离散化为

$ u = {u_0} + m\Delta I $

(3)

式中：−M/2≤m≤M/2；u₀=(u₂−u₁)/2为区间中点；ΔI=(u₂−u₁)/(M−1)为采样间隔。然后将式(3)代入FRFT离散化公式中，得到

$\begin{array}{c} {X_p}({u_0} + m\Delta I) = \\ \dfrac{{{A_\alpha }}}{{2\Delta x}}\exp ({\rm{j{\text{π}} }}\gamma {({u_0} + m\Delta I)^2})\displaystyle\sum\limits_{n = - N}^N {\exp \left( {{\rm{j{\text{π}} }}\beta \dfrac{{\Delta I}}{{2\Delta x}}( - 2mn)} \right)} \cdot \\ \exp \left( { - {\rm{j2{\text{π}} }}\beta {u_0}\left( {\dfrac{n}{{2\Delta x}}} \right)} \right)\exp \left( {{\rm{j{\text{π}} }}\gamma {{\left( {\dfrac{n}{{2\Delta x}}} \right)}^2}} \right)x\left( {\dfrac{n}{{2\Delta x}}} \right)\\[-16pt] \end{array} $

(4)

将恒等式−2mn=(m−n)²−m²−n²代入式(3)，替换得到

$\begin{array}{c} {X_{\rm{p}}}({u_0} + m\Delta I) = \\ \dfrac{{{A_\alpha }}}{{2\Delta x}}\exp \left( {{\rm{j{\text{π}} }}\gamma {{({u_0} + m\Delta I)}^2}} \right)\exp \left( { - {\rm{j{\text{π}} }}\beta \dfrac{1}{{2\Delta I\Delta x}}{{(m\Delta I)}^2}} \right) \cdot \\ \displaystyle\sum\limits_{n = - N}^N {\exp \left( {{\rm{j2{\text{π}} }}\beta \Delta I\Delta x{{\left( {\dfrac{{m - n}}{{2\Delta x}}} \right)}^2}} \right)} \cdot \\ \exp \left( {{\rm{j{\text{π}} }}\left( { - 2\beta {u_0}\left( {\dfrac{n}{{2\Delta x}}} \right) + } \right.} \right.\left. {\left. {(\gamma - 2\beta \Delta I\Delta x){{\left( {\dfrac{n}{{2\Delta x}}} \right)}^2}} \right)} \right)x\left( {\dfrac{n}{{2\Delta x}}} \right) \end{array} $

(5)

式(4)的计算需要直接确定u₀和ΔI这2个参数的取值，为此引入了2个相对因子，一个是“平移因子”λ=u₀/x，−0.5≤λ≤0.5，它表示局部谱中心u₀在整个谱范围[−Δx/2，Δx/2]中的相对位置；另一个是“变焦因子”P=1/(2ΔIΔx)，它表示局部谱的分辨率ΔI相对于标准分辨率(1/2Δx)的放大倍数，一般取P为大于1的自然数。将关系式u₀=λΔx和ΔI=1/(2PΔx)代入到式(4)得到

$\begin{array}{c} {X_p}\left( {\lambda \Delta x + \dfrac{m}{{2P\Delta x}}} \right) = \dfrac{{{A_\alpha }}}{{2\Delta x}}\exp \left( {{\rm{j{\text{π}} }}\gamma \Delta {x^2}{\lambda ^2}} \right) \cdot \\ \exp \left( {{\rm{j2{\text{π}} }}\gamma \lambda \Delta x\left( {\dfrac{m}{{2P\Delta x}}} \right)} \right)\exp \left( {{\rm{j{\text{π}} }}(\gamma - P\beta ){{\left( {\dfrac{m}{{2P\Delta x}}} \right)}^{\rm{2}}}} \right) \cdot \\ \displaystyle\sum\limits_{n = - N}^N {\exp \left( {{\rm{j{\text{π}} }}\dfrac{{\rm{1}}}{P}\beta {{\left( {\dfrac{{m - n}}{{2\Delta x}}} \right)}^2}} \right)} \exp \left( {{\rm{j{\text{π}} }}\left( {\gamma - \dfrac{{\rm{1}}}{P}\beta } \right){{\left( {\dfrac{n}{{2\Delta x}}} \right)}^2}} \right) \cdot \\ {\rm{ }}\exp \left( {{\rm{ - j2{\text{π}} }}\beta \lambda \Delta x\left( {\dfrac{n}{{2\Delta x}}} \right)} \right)x\left( {\dfrac{n}{{2\Delta x}}} \right) \\ \end{array} $

记 $\bar{x}[n]=x\left(\dfrac{n}{2 \Delta x}\right),\; \bar{X}_{\alpha, \lambda, P}[m]=X_{p}\left(\lambda \Delta x+\dfrac{m}{2 P \Delta x}\right)$ ，式(5)可简化为

$\begin{array}{c} {{\bar X}_{\alpha ,\lambda ,P}}\left[ m \right] = \dfrac{{{A_\alpha }\exp ({\rm{j{\text{π}} }}\gamma N{\lambda ^2})}}{{2\sqrt N }} \cdot \\ \exp \left( {{\rm{j{\text{π}} }}\left( {\dfrac{{\gamma \lambda }}{P}} \right)m} \right)\exp {\rm{j{\text{π}} }}\left( {j{\text{π}} \left( {\dfrac{{\gamma - P\beta }}{{4PN}}} \right){m^{\rm{2}}}} \right) \cdot \\ \displaystyle\sum\limits_{n = - N}^N {\exp \left( {{\rm{j{\text{π}} }}\left( {\dfrac{\beta }{{4PN}}} \right){{(m - n)}^2}} \right)} \cdot \\ \exp \left( {{\rm{j{\text{π}} }}\left( {\dfrac{{P\gamma - \beta }}{{4PN}}} \right){n^2}} \right)\exp \left( { - {\rm{j{\text{π}} }}(\beta \lambda )n} \right)\bar x\left[ n \right] \\ \end{array} $

图4为高分辨局部谱对比曲线。由2条谱线对比可见，Zoom-FRFT弥补了离散计算过程中损失的信息，能够更好地将细节信息展现出来。

2 FRFT对参数估计的影响因素

基于FRFT的信号二维搜索检测法的基本思路可叙述为：令变换阶次p遍历，对观测信号连续进行FRFT，形成信号能量在参数空间(p,u)平面上的二维分布X_p(u)。在此平面上搜寻最大峰值点的坐标 $\left(\hat{p}_{0}, \hat{u}_{0}\right)$ ，由此可实现信号的检测与参数估计。该过程可描述为

$\{ {\hat p_0},{\hat u_0}\} = \arg \mathop {\max }\limits_{p,u} |{X_p}{\rm{(}}u{\rm{)|}}$

(6)

$\left\{ {\begin{array}{*{20}{c}} {\hat \gamma = - \cot ({{\hat p}_0}{\text{π}} /2)} \\ {{{\hat f}_0} = {{\hat u}_0}\csc ({{\hat p}_0}{\text{π}} /2)} \end{array}} \right.$

(7)

由式(6)、(7)可知，变换阶次 ${\hat p_0}$ 决定了调频斜率估计值 $\hat \gamma $ 的准确性，变换阶次 ${\hat p_0}$ 与FRFT域坐标 ${\hat u_0}$ 共同决定了中心频率估计值 ${\hat f_0}$ 的准确性。因为离散过程所导致的峰值偏差会对参数估计产生不可避免的误差，所以需要同时提高变换阶次p的搜索精度与FRFT域的分辨率。

3 LVD与Zoom-FRFT联合的参数估计方法

针对FRFT中搜索步长所导致的计算量激增问题，本文借助LVD能够有效地展现LFM信号中多个分量的优点。首先采用LVD对多分量信号的参数进行粗略估计，然后采用FRFT与Zoom-FRFT对各分量进行精估计，在保证估计精度的前提下，提高估计效率。

LVD与Zoom-FRFT联合的多分量LFM信号参数估计步骤如下：

1)对多分量LFM信号进行LVD，在中心频率−调频斜率域中提取谱峰，确定分量个数。

2)对提取到的分量的粗估计调频斜率 ${\hat \gamma _c}$ 进行标记，代入式(7)，计算各分量粗估计变换阶数 ${\hat p_c}$ 。

3)以 ${\hat p_c}$ 为中心，限定搜索范围r，设置搜索步长δ，在限定区间内[ ${\hat p_c}$ −r/2, ${\hat p_c}$ +r/2]，对信号进行FRFT，提取峰值对应的变换阶次 ${\hat p_0}$ 与FRFT域坐标 ${\hat u_c}$ 。

4)在分数阶次 ${\hat p_0}$ 下，选定FRFT域坐标 ${\hat u_c}$ 为中心的区间，进行Zoom-FRFT，提取峰值的FRFT域坐标 ${\hat u_0}$ 。

5)将变换阶次 ${\hat p_0}$ 与FRFT域坐标 ${\hat u_0}$ 代入式(7)，计算各分量的调频斜率估计值 $\hat \gamma $ 、中心频率估计值 ${\hat f_0}$ 。

4 仿真实验 4.1 估计精度分析 4.1.1 调频率

以2个分量的情况进行分析，采样率f_s=128 Hz，采样时长T_p=2 s。分量1的参数为A₁=1，f_c=4 Hz，K₁=16 Hz/s；分量2的参数为A₁=1，f_c=8 Hz，K₁=8 Hz/s。加入高斯白噪声，信噪比20 dB，取搜索精度范围[10⁻²,10⁻⁷]，递减10倍，分别进行100次蒙特卡洛实验，统计调频率相对误差，如图5所示。

由5可知，当搜索精度达到10⁻⁴时，搜索精度对调频率产生的影响趋于平稳。在该条件下，搜索精度对中心频率估计的影响同样趋于平稳，需要对FRFT域进行细化。

4.1.2 中心频率

采用本文估计步骤，信号分量条件与4.1.1节相同，加入高斯白噪声，信噪比20 dB，固定搜索精度为10⁻⁴，Zoom-FRFT细化倍数为[10:80]，递进10，分别进行100次蒙特卡洛实验，统计中心频率相对误差，如图6所示，Zoom-FRFT对FRFT域的细化使相对误差得到改善。

4.2 抗噪性能分析

采样率f_s=256 Hz、采样时长T_p=1 s，信号分量条件与4.1.1节相同，加入高斯白噪声，取信噪比范围[−15 dB,20 dB]，递进1 dB，分别进行300次蒙特卡洛实验，统计分量估计参数的均方误差。

采用参考文献[3]中算法与本文算法进行比较，由图7可见，在本文实验条件下，所提算法的信噪比门限为−11 dB，当信噪比大于−11 dB时，对参数的估计性能远高于文献[3]中算法。

4.3 计算量分析

本文方法的计算量主要源于LVD、FRFT和Zoom-FRFT，对信号做LVD的计算次数为6N²logN，N为输入点数。对信号做FRFT的计算量为L(2N+NlogN)，L为变换阶次离散点数，N为输入点数。对信号做N点输入M点输出的Zoom-FRFT的计算次数为3(2N+M)log(4N+8M)+8M，文中取M=N，计算次数为9Nlog(12N)+8N。

根据估计精度分析，搜索精度采用10⁻⁴时，对误差的影响趋于平稳，在该条件下，采用FRFT的计算次数为20 000(2N+NlogN)。本文利用LVD粗估计完成对FRFT的搜索范围限定，FRFT的计算量缩减到200(2N+NlogN)，再进行Zoom-FRFT，对FRFT域进行细化，计算量为9Nlog(12N)+8N。设分量个数K，本文算法总计算量为6N²logN+K(200×(2N+NlogN)+ 9Nlog(12N)+8N)，当K≤3时，信号采样点数小于3 770，本文算法计算量小于FRFT的计算量。

5 结论

通过在不同信噪比下进行统计实验与计算量分析，验证了本文算法的有效性并得出结论：

1)在信噪比满足的情况下，可以快速精准地对多LFM信号参数进行估计。

2)通过LVD算法的结合，减小了总运算次数，提高FRFT算法的运算效率。