CN105842713A - 一种基于sft的ins辅助bds信号快速捕获方法 - Google Patents

Abstract

本发明是一种基于稀疏傅里叶变换的INS辅助BDS信号快速捕获方法，属于卫星导航领域。本方法包括：通过天线接收北斗卫星信号，经过射频模块处理后得到中频信号；利用INS的惯性信息和卫星星历计算多普勒频率，并引入到本地载波NCO中，调整本地载波NCO的搜索频率范围，输出本地载波正余弦两路信号；将中频信号与本地载波混频，得到的同相I支路和正交Q支路信号构成复数序列，对复数序列和本地伪码求快速傅里叶变换，将结果进行复数相乘；将相乘结果进行稀疏傅里叶变换得到的信号峰值，判断是否成功捕获。本发明可以大幅缩小BDS信号的频率预估范围，简化IFFT运算的同时，保证了捕获性能，大幅降低运算量，减少硬件消耗。

Description

一种基于SFT的INS辅助BDS信号快速捕获方法

技术领域

本发明涉及一种基于稀疏傅里叶变换(SFT)的惯性导航***(INS)辅助BDS信号快速捕获方法，尤其适用于高动态环境下的捕获，属于卫星导航领域。

背景技术

目前，卫星导航技术已成为世界各国在现代化战争中的重要武器组成部分。北斗卫星导航***(BDS)是中国自行研制、自主控制的卫星定位导航***，已经实现亚太地区的定位、导航、授时、短报文等服务；BDS将在我国经济***会的发展中发挥重大作用。

卫星导航***导航技术应用十分广泛，具有高精度、低成本、定位误差不随时间积累等诸多优点，但是在高动态环境下存在先天的局限性——不易捕获和跟踪卫星载波信号，甚至出现所谓的“周跳”现象。然而，在航空、航天、精确武器打击等应用领域，所研究的飞行器大都属于高动态，如何满足高动态环境的需要，对导弹等进行精确的导航，成为全球卫星导航***(GNSS)在高动态领域研究的重点。

随着导航技术的发展和日益成熟，GNSS/INS组合导航可以很好的改善***的动态性能，基于INS辅助的GNSS捕获技术可以快速地完成卫星信号的捕获。尤其，在卫星接收机冷启动以及高动态环境下卫星失锁后的快速重捕获方面，有非常好的优势。利用INS的惯性信息得到高速运动载体相对卫星的多普勒频移，可以大幅缩小捕获过程中频率搜索范围，从而缩短卫星信号的捕获时间。

卫星捕获方法的选择直接影响着捕获的时间和性能，目前常用的捕获方法有：滑动相关法、匹配滤波法、并行频率法、并行码相位法等。前两种为时域法，捕获时间较长；后两种为频域法，捕获时间短，但复杂度高运算量大。

发明内容

针对现有捕获方法存在的问题，本发明提出了一种基于稀疏傅里叶变换的INS辅助BDS信号快捕获方法。本发明方法在运用并行码相位法捕获时利用SFT简化离散傅里叶逆变换(IFFT)运算，将运算的复杂度大大降低，结合INS惯性信息，可以实现高动态环境下BDS信号的快速捕获。

本发明的一种基于稀疏傅里叶变换的INS辅助BDS信号快速捕获方法，具体实现步骤为：

步骤一：通过天线接收北斗卫星信号，经过射频模块处理后得到中频信号S_IF(t)。

$S_{IF}\left(t\right)=\sqrt{2P_s}C(t-\mathrm{\τ})D(t-\mathrm{\τ})cos\[{\mathrm{\ω}}_{IF}t+\mathrm{\φ}\left(t\right)\]+n\left(t\right)$

其中，t表示时间，P_s表示信号功率；C(t)表示测距码，其值为土1；τ表示在传输过程中的时间延时；D(t)表示导航电文，比特率是50bps；ω_IF表示中频载波频率；φ(t)表示初始载波相位；n(t)表示白噪声。

步骤二：利用INS的惯性信息和卫星星历计算多普勒频率，将得到的多普勒频率引入到本地载波数字控制振荡器(NCO)中，调整本地载波NCO的搜索频率范围为[-f_max,f_max]，f_max表示最大搜索频率。

输出本地载波的正余弦两路信号分别表示为:

$I_{ca}\left(t\right)=\sqrt{2}cos\[({\mathrm{\ω}}_{IF}+\mathrm{\Δ}\widehat{\mathrm{\ω}})t+{\widehat{\mathrm{\φ}}}_0\];Q_{ca}\left(t\right)=\sqrt{2}sin\[({\mathrm{\ω}}_{IF}+\mathrm{\Δ}\widehat{\mathrm{\ω}})t+{\widehat{\mathrm{\φ}}}_0\];$

其中，I_ca(t)为正弦信号，Q_ca(t)为余弦信号，为本地复制码初始载波相位，为本地搜索频移，为本地载波频率。

步骤三：将中频信号与本地载波混频，得到的同相I支路和正交Q支路信号构成复数序列信号，对其求快速傅里叶变换(FFT)；对本地伪码求FFT，并将两个FFT结果进行复数相乘。

将本地载波正弦信号I_ca(t)、余弦信号Q_ca(t)分别与中频信号S_IF(t)相乘，得到信号I(t)和Q(t)，再经过积分器低通滤波后，I(t)和Q(t)简化为：

$I\left(t\right)=\sqrt{P_s}C(t-\mathrm{\τ})D(t-\mathrm{\τ})cos\[\mathrm{\φ}\left(t\right)-\mathrm{\Δ}\widehat{\mathrm{\ω}}t-{\widehat{\mathrm{\φ}}}_0\],$

$Q\left(t\right)=\sqrt{P_s}C(t-\mathrm{\τ})D(t-\mathrm{\τ})sin\[\mathrm{\φ}\left(t\right)-\mathrm{\Δ}\widehat{\mathrm{\ω}}t-{\widehat{\mathrm{\φ}}}_0\];$

得到复数序列信号H(t)＝I(t)+jQ(t)，对H(t)进行FFT得到频域值H(k)。

设本地伪码为 为本地伪码相位，经FFT获得的频域值C(k)。

将两个FFT结果进行复数相乘，得到信号Z(k)＝H(k)C(k)。

步骤四：将相乘结果Z(k)进行稀疏傅里叶变换(SFT)。

对Z(k)进行频域重排，得到结果Z_p(k)。然后通过SFT估计得到Z(k)的信号峰值。

步骤五：根据步骤四得到的信号峰值判断是否成功捕获。

将步骤四中得到的信号峰值，与已设定的捕获阈值进行比较，大于阈值则信号捕获成功，转入信号跟踪环路，小于阈值则调整载波频率和码频率重新捕获，直到信号捕获成功为止。

本发明的优点和积极效果在于：本发明利用INS惯性信息和历书信息获得载体高速运动产生的多普勒频移和卫星运动产生的多普勒频移，可以大幅缩小BDS信号的频率预估范围。根据捕获相关峰的稀疏特性，在运用并行码相位法时利用SFT简化IFFT运算，可以大幅降低运算量，减少硬件消耗，同时又保证了捕获性能，从而实现了高动态环境下BDS信号的快速捕获。

附图说明

图1是本发明基于SFT的INS辅助BDS信号快速捕获结构示意图；

图2是本发明基于SFT的INS辅助BDS信号快速捕获方法的流程示意图；

图3是导航卫星信号软件接收机结构图；

图4是SFT运算过程示意图。

具体实施方式

下面将结合附图及具体实施方式对本发明加以详细说明，需要指出的是，所描述的实施例仅旨在便于对本发明的理解，而不起任何限定作用。

本发明基于稀疏傅里叶变换的INS辅助BDS信号快速捕获方法的一个实施例的实现结构，如图1所示，实现流程如图2所示。本发明方法可基于卫星信号软件接收机很方便的实现，图3是软件接收机的结构图。如图3所示，除了射频前端与数字采样模块外，在通用的基础硬件平台上，将接收机的基带信号处理(捕获、跟踪)与导航信息处理(解码、导航解算)部分都用软件进行，软件接收机具有低价格、小型化、方便灵活、便于扩展等优点。

下面说明本发明方法的各个步骤的具体实现。

步骤一：天线接收北斗卫星信号，经过射频模块处理后得到中频信号。

通过天线接收北斗卫星导航信号，经过射频前端的滤波器、放大器、下变频器、数模转换器等处理，输出中频信号S_IF(t)表示为：

$S_{IF}\left(t\right)=\sqrt{2P_s}C(t-\mathrm{\τ})D(t-\mathrm{\τ})cos\[{\mathrm{\ω}}_{IF}t+\mathrm{\φ}\left(t\right)\]+n\left(t\right)---\left(1\right)$

式中，P_s表示信号功率；C(t)表示测距码，其值为土1；τ表示在传输过程中的时间延时；D(t)表示导航电文，比特率是50bps；ω_IF表示中频载波频率，由射频前端和多普勒频移决定；φ(t)表示初始载波相位；n(t)表示白噪声。

步骤二：利用INS的惯性和星历、历书信息计算多普勒频率，并根据其调整搜索环路的本地载波数字控制振荡器(NCO)的载波频率。

卫星与载体的相对运动是北斗卫星信号产生多普勒频移的主要因素，利用卫星星历、历书和INS的惯性信息可以计算出其产生的多普勒频移为：

其中，f_rec表示载体运动产生的多普勒频率；f_s表示卫星运动产生的多普勒频率；λ表示载波的波长，表示载体的速度，表示卫星的速度，表示卫星到用户连线方向的单位矢量。可由卫星历书估算得到，由INS提供。

将INS和卫星历书估算得到的多普勒频率引入到载波NCO中，缩小频率方向的搜索范围。设搜索频率范围为[-f_max,f_max]，f_max表示最大搜索频率。

捕获环路使用INS数据和卫星星历、历书数据提供的载体相对于卫星运动产生的多普勒频率来不断的调整复制码和复制的载波频率。

NCO输出本地载波的正余弦两路信号分别表示为:

$I_{ca}\left(t\right)=\sqrt{2}cos\[({\mathrm{\ω}}_{IF}+\mathrm{\Δ}\widehat{\mathrm{\ω}})t+{\widehat{\mathrm{\φ}}}_0\]---\left(3\right)$

$Q_{ca}\left(t\right)=\sqrt{2}\sin \[({\mathrm{\ω}}_{IF}+\mathrm{\Δ}\widehat{\mathrm{\ω}})t+{\widehat{\mathrm{\φ}}}_0\]---\left(4\right)$

式中，I_ca(t)表示输出的正弦信号，Q_ca(t)表示输出的余弦信号，表示本地复制码初始载波相位，表示本地载波频率，表示本地搜索频移。

步骤三：将步骤一得到的中频信号与步骤二得到的本地载波混频，得到的同相I支路和正交Q支路信号构成复数序列，对其求FFT；对本地伪码求FFT，并将两个FFT结果进行复数相乘。

将本地载波的正余弦两路信号分别送给混频器和输入的中频信号相乘，得到：

$I\left(t\right)=\sqrt{2P_s}C(t-\mathrm{\τ})D(t-\mathrm{\τ})cos\[{\mathrm{\ω}}_{IF}t+\mathrm{\φ}\left(t\right)\]\×\sqrt{2}cos\[({\mathrm{\ω}}_{IF}+\mathrm{\Δ}\widehat{\mathrm{\ω}})t+{\widehat{\mathrm{\φ}}}_0\]---\left(5\right)$

$Q\left(t\right)=\sqrt{2P_s}C(t-\mathrm{\τ})D(t-\mathrm{\τ})cos\[{\mathrm{\ω}}_{IF}t+\mathrm{\φ}\left(t\right)\]\×\sqrt{2}sin\[({\mathrm{\ω}}_{IF}+\mathrm{\Δ}\widehat{\mathrm{\ω}})t+{\widehat{\mathrm{\φ}}}_0\]---\left(6\right)$

信号I(t)和Q(t)经过积分器(LPF)实现低通滤波，忽略高频成分，简化为：

$I\left(t\right)=\sqrt{P_s}C(t-\mathrm{\τ})D(t-\mathrm{\τ})cos\[\mathrm{\φ}\left(t\right)-\mathrm{\Δ}\widehat{\mathrm{\ω}}t-{\widehat{\mathrm{\φ}}}_0\]---\left(7\right)$

$Q\left(t\right)=\sqrt{P_s}C(t-\mathrm{\τ})D(t-\mathrm{\τ})sin\[\mathrm{\φ}\left(t\right)-\mathrm{\Δ}\widehat{\mathrm{\ω}}t-{\widehat{\mathrm{\φ}}}_0\]---\left(8\right)$

同相I支路信号I(t)和正交Q支路信号Q(t)构成复数序列H(t)＝I(t)+jQ(t)，j为虚部单位。对H(t)利用FFT计算得到信号H(t)的频域值H(k)，k表示频点。

假设本地伪码发生器产生的本地伪码表示为 表示某一本地伪码相位，则本地伪码和接收信号的伪码相位相差利用FFT计算信号的频域值C(k)。

将两个FFT结果进行复数相乘，得到相乘结果Z(k)，如下所示：

Z(k)＝H(k)C(k) (9)

步骤四：将相乘结果Z(k)进行稀疏傅里叶变换。

根据捕获相关峰的稀疏特性，在计算相关时利用SFT简化IFFT运算。

首先对Z(k)进行频域重排得到Z_p(k)。然后通过平滑窗函数和随机排列处理，在控制参量c_t与m的辅助下运用SFT算法估计Z(k)的相关峰。对Z(k)进行SFT运算的过程如图4所示。

SFT主要包括内循环和外循环两个环节。

外循环的主要步骤为：

(1.1)选择参数ε、δ、B、d′、L_s，满足下式：

$\left\{\begin{array}{c}\mathrm{\&epsiv;}<1\\ \mathrm{\&delta;}=\frac{1}{n^{O\left(1\right)}}\\ B=O\left(\sqrt{\frac{nK}{\mathrm{\&epsiv;}}\log \frac{n}{\mathrm{\&delta;}}}\right)\\ d^{\&prime;}=O\left(\frac{1}{\mathrm{\&epsiv;}}\right)\\ L_s=O\left(\log n\right)\end{array}\right.---\left(10\right)$

其中，ε规定了阻带起始频率与采样频率的比值，δ规定了阻带和通带共同的起伏特性，B为频域重排分组数，d′表示稀疏性因子，L_s表示内循环执行次数，n表示信号长度，K表示稀疏度。O(.)表示空间复杂度。

(1.2)对于参数r∈[1,2,...L_s]，以满足supp(G)＜σ′＜N/K的参数σ′执行定位内循环，得到向量I_r。其中，σ′是一个在区间[0,n-1]上随机选取的奇数，N表示采样点数，G表示平滑窗函数，supp(·)表示函数的支集。

(1.3)令I＝I₁UI₂U...UI_r，对于i∈I，令s_i＝|{r|i∈I_r}|。

(1.4)计算

(1.5)对于r∈[1,2,...L_s]，在I″上计算信号估计谱的幅值得到估计结果

其中，x表示对Z(k)频域重排后得到的Z_p(k)，为x的估计谱，为估计的x在第i个频点处幅值，为x进行降采样后的序列的估计谱。h_σ′(i)为映射函数，参见公式(12)。

表示在参数r下估计谱在第i个频点处幅值。

(1.6)对于集合I″中的所有元素，计算中位值

(1.7)对于上述中位值构成的集合求出其中所有元素中的最大值

(1.8)设集合中不小于共有m个，其中c_t是算法选择的一个比例系数，其数值根据实际需要确定。

(1.9)令对于r∈[1,2,...L_s]，在I′上运行估计循环，得到估计结果

(1.10)对于集合I′中的所有元素，计算中位值

内循环的主要步骤为：

(2.1)随机选取参数σ′∈[n]与τ_s∈[n]，使得σ′是奇数。n表示信号长度，[n]表示一个数集[0,n-1]。

(2.2)记x经过窗函数滤波和随机排列的结果为y，即满足条件：

$\left\{\begin{array}{c}{y}_i=G_i{x}_{{\mathrm{\&sigma;}}^{\&prime;}i+{\mathrm{\&tau;}}_s}\\ \sup p\left(y\right)\&SubsetEqual;\sup p\left(G\right)=\&lsqb;\mathrm{\&omega;}\&rsqb;\end{array}\right.---\left(11\right)$

其中，表示随机排列方式，ω表示滤波器长度，[ω]表示一个数集[0,ω-1]，i∈[B]，B为频域重排分组数，[B]表示一个数集[0,B-1]。表示x在第σ′i+τ_s个频点处信号，y_i表示y的第i个频点处信号，G_i为第i个频点对应的平滑窗函数。

(2.3)对于i∈[B]，定义序列对序列z_i进行傅里叶变换得到

(2.4)设映射函数h_σ′(i)和o_σ′(i)满足下面条件：

$\left\{\begin{array}{c}{h}_{{\mathrm{\&sigma;}}^{\&prime;}}\left(i\right)=round\&lsqb;\frac{{\mathrm{\&sigma;}}^{\&prime;}iB}{n}\&rsqb;\\ o_{{\mathrm{\&sigma;}}^{\&prime;}}\left(i\right)={\mathrm{\&sigma;}}^{\&prime;}i-h_{{\mathrm{\&sigma;}}^{\&prime;}}\left(i\right)\frac{n}{B}\end{array}\right.---\left(12\right)$

其中round(·)表示取整操作。

(2.5)定位内循环：将中幅值最大的d′K个谱峰表示为集合J，对于J中的每个元素，令集合H满足条件：H＝{i∈[n]|h_σ′(i)∈J}。

(2.6)估计内循环：对于i∈H，对输入信号谱进行估计，按下式计算

${\tilde{x}}_i^{\&prime;}={\hat{z}}_{h_{{\mathrm{\&sigma;}}^{\&prime;}}\left(i\right)}{e}^{-j\frac{2\mathrm{\&pi;}}{n}\mathrm{\&tau;}i}/{\hat{G}}_{o_{{\mathrm{\&sigma;}}^{\&prime;}}\left(i\right)}---\left(13\right)$

其中，为x在映射o_σ′(i)下平滑窗函数的谱。

就是最后得到的信号峰值。

外循环是为了检验内循环的结果。由于重排具备随机性，且整个过程受能量泄漏与降采样的影响，为了使计算的峰值位置与峰值幅值结果具有较高的置信度，增加了外循环过程。

步骤五：根据步骤四中结果判断是否成功捕获。

将步骤四中得到的信号峰值，与已设定的捕获阈值进行比较，大于阈值则信号捕获成功，得到码相位估计和粗略的载波频率，进行码环跟踪。根据码相位估计调整本地伪码发生器，同输入的中频信号调制将B1码剥离；然后根据粗略的载波频率调整本地载波NCO，在粗略的载波频率附近窄带频率范围利用Zoom-FFT频率细化方法进行频率细化处理，得到精细的载波频率，转入信号跟踪环路。若小于阈值，则转步骤二继续执行，调整本地载波频率和码频率重新进行信号捕获，直到信号捕获成功为止。

本发明在高动态环境下，利用INS惯性信息可以计算出载体高速运动产生的多普勒频移，利用历书信息可以得到卫星运动产生的多普勒频移，利用以上信息可以大幅缩小BDS信号的频率预估范围，运用并行码相位捕获可以实现2046个码相位的并行搜索，因此本发明方法可以明显减少捕获时间。经过计算，本发明提出的方法相较于传统的基于FFT的方法可以把运算复杂度降低1个数量级以上。综述，本发明即实现了高动态环境下的快速捕获有可以大幅降低运算复杂度，具有很好的工程应用价值。

以上所述仅为本发明的具体实施例，并不用于限定本发明的保护范围，凡在本发明的精神和原则之内，所做的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (3)

1.一种基于稀疏傅里叶变换的INS辅助BDS信号快速捕获方法，其特征在于，包括如下步骤：

步骤一，天线接收北斗卫星信号，经过射频模块处理后得到中频信号S_IF(t)，如下：

$S_{IF}\left(t\right)=\sqrt{2P_s}C(t-\mathrm{\&tau;})D(t-\mathrm{\&tau;})cos\&lsqb;{\mathrm{\&omega;}}_{IF}t+\mathrm{\&phi;}\left(t\right)\&rsqb;+n\left(t\right)$

其中，t表示时间，P_s表示信号功率；C(t)表示测距码，其值为土1；τ表示在传输过程中的时间延时；D(t)表示导航电文，比特率是50bps；ω_IF表示中频载波频率；φ(t)表示初始载波相位；n(t)表示白噪声；

步骤二：利用INS的惯性信息和卫星历书计算多普勒频率，将得到的多普勒频率引入到本地载波NCO中，调整本地载波NCO的搜索频率范围；NCO表示数字控制振荡器；

输出本地载波正余弦两路信号分别表示为：

$I_{ca}\left(t\right)=\sqrt{2}cos\&lsqb;({\mathrm{\&omega;}}_{IF}+\mathrm{\&Delta;}\widehat{\mathrm{\&omega;}})t+{\widehat{\mathrm{\&phi;}}}_0\&rsqb;,Q_{ca}\left(t\right)=\sqrt{2}sin\&lsqb;({\mathrm{\&omega;}}_{IF}+\mathrm{\&Delta;}\widehat{\mathrm{\&omega;}})t+{\widehat{\mathrm{\&phi;}}}_0\&rsqb;;$

其中，I_ca(t)为正弦信号，Q_ca(t)为余弦信号，为本地复制码初始载波相位，为本地搜索频移，为本地载波频率；

步骤三：将中频信号与本地载波混频，得到的同相I支路和正交Q支路信号构成复数序列，对其求FFT；对本地伪码求FFT，将两个FFT结果进行复数相乘；FFT表示快速傅里叶变换；

将本地载波正弦信号、余弦信号分别与中频信号相乘，得到同相I支路信号I(t)和正交Q支路信号Q(t)，经过积分器低通滤波后，简化为：

$I\left(t\right)=\sqrt{P_s}C(t-\mathrm{\&tau;})D(t-\mathrm{\&tau;})cos\&lsqb;\mathrm{\&phi;}\left(t\right)-\mathrm{\&Delta;}\widehat{\mathrm{\&omega;}}t-{\widehat{\mathrm{\&phi;}}}_0\&rsqb;,$

$Q\left(t\right)=\sqrt{P_s}C(t-\mathrm{\&tau;})D(t-\mathrm{\&tau;})sin\&lsqb;\mathrm{\&phi;}\left(t\right)-\mathrm{\&Delta;}\widehat{\mathrm{\&omega;}}t-{\widehat{\mathrm{\&phi;}}}_0\&rsqb;;$

得到的复数序列H(t)＝I(t)+jQ(t)，对H(t)进行FFT得到频域值H(k)；

设本地伪码为 为本地伪码相位，经FFT获得的频域值C(k)；

将两个FFT结果进行复数相乘，得到信号Z(k)＝H(k)C(k)；

步骤四：将信号Z(k)进行频域重排得到信号Z_p(k)，再进行SFT，获得Z(k)的信号峰值；SFT表示稀疏傅里叶变换；

步骤五：将步骤四中得到的信号峰值，与已设定的捕获阈值进行比较，大于阈值则信号捕获成功，转入信号跟踪环路，小于阈值则转步骤二执行，调整本地载波频率和码频率重新捕获，直到信号捕获成功为止。

2.根据权利要求1所述的基于稀疏傅里叶变换的INS辅助BDS信号快速捕获方法，其特征在于，所述的步骤二中，利用卫星历书和INS的惯性信息计算得到的多普勒频移为：

其中，f_rec表示载体运动产生的多普勒频率；f_s表示卫星运动产生的多普勒频率；λ表示载波的波长，表示载体的速度，表示卫星的速度，表示卫星到用户连线方向的单位矢量。

3.根据权利要求1所述的基于稀疏傅里叶变换的INS辅助BDS信号快速捕获方法，其特征在于，所述的步骤四中，进行SFT，包括内循环和外循环两个环节；

实现外循环的步骤包括(1.1)～(1.10)；

(1.1)选择参数ε、δ、B、d′、L_s，满足下式：

$\left\{\begin{array}{c}\mathrm{\&epsiv;}<1\\ \mathrm{\&delta;}=\frac{1}{n^{O\left(1\right)}}\\ B=O\left(\sqrt{\frac{nK}{\mathrm{\&epsiv;}}log\frac{n}{\mathrm{\&delta;}}}\right)\\ d^{\&prime;}=O\left(\frac{1}{\mathrm{\&epsiv;}}\right)\\ L_s=O\left(\log n\right)\end{array}\right.$

其中，ε规定了阻带起始频率与采样频率的比值，δ规定了阻带和通带共同的起伏特性，B为频域重排分组数，d′表示稀疏性因子，L_s表示内循环执行次数，n表示信号长度，K表示稀疏度，O(.)表示空间复杂度；

(1.2)对于参数r∈[1,2,...L_s]，以满足supp(G)＜σ′＜N/K的参数σ′执行定位内循环，得到向量I_r；其中，σ′是一个在区间[0,n-1]上随机选取的奇数，N表示采样点数，G表示平滑窗函数，supp(·)表示函数的支集；

(1.3)令集合I＝I₁UI₂U...UI_r，对于参数i∈I，令参数s_i＝|{r|i∈I_r}|；

(1.4)计算

(1.5)对于r∈[1,2,...L_s]，在I″上计算信号x的估计谱的幅值，得到估计结果其中，x＝Z_p(k)， 表示在参数r下估计谱在第i个频点处幅值；

(1.6)对于集合I″中的所有元素，计算中位值

(1.7)对于中位值构成的集合求出其中的最大值

(1.8)设集合中不小于共有m个，其中c_t为设定的比例系数；

(1.9)令对于r∈[1,2,...L_s]，在I′上运行估计循环，得到估计结果

(1.10)对于集合I′中的所有元素，计算中位值

实现内循环的步骤包括(2.1)～(2.6)；

(2.1)随机选取参数σ′∈[n]与τ_s∈[n]，且σ′是奇数；

(2.2)记x进行随机排列和频域滤波的结果为y，满足条件：

$\left\{\begin{array}{c}{y}_i=G_i{x}_{{\mathrm{\&sigma;}}^{\&prime;}i+{\mathrm{\&tau;}}_s}\\ \sup p\left(y\right)\&SubsetEqual;\sup p\left(G\right)=\&lsqb;\mathrm{\&omega;}\&rsqb;\end{array}\right.$

其中，表示随机排列方式，ω表示滤波器长度，i∈[B]，y_i表示y的第i个频点处信号，表示x在第σ′i+τ_s个频点处信号，G_i为第i个频点对应的平滑窗函数；

(2.3)对于i∈[B]，定义序列对序列z_i进行傅里叶变换得到信号

(2.4)设映射函数h_σ′(i)和o_σ′(i)满足下面条件：

其中round(·)表示取整操作；

(2.5)定位内循环：将中幅值最大的d′K个谱峰表示为集合J，对于J中的每个元素，令集合H满足条件：H＝{i∈[n]|h_σ′(i)∈J}；

(2.6)估计内循环：对于i∈H，按下式计算信号峰值如下：

${{\tilde{x}}_i}^{\&prime;}={\hat{z}}_{h_{{\mathrm{\&sigma;}}^{\&prime;}}\left(i\right)}{e}^{-j\frac{2\mathrm{\&pi;}}{n}\mathrm{\&tau;}i}/{\hat{G}}_{o_{{\mathrm{\&sigma;}}^{\&prime;}}\left(i\right)}$

其中，τ表示在传输过程中的时间延时，为x在映射o_σ′(i)下平滑窗函数的谱。