高动态环境下载波频率的精确估计算法

金 磊

中国西南电子技术研究所，成都 610036

当飞机、火箭和导弹等载体作大机动运动时，接收信号的载波频率会产生极大的多普勒频移和多普勒变化率[1]，而当载体处于障碍遮挡、多径严重和信号干扰等环境时，信号功率由于受到各种信道损耗而严重衰减[2]，因此，高动态、弱信号等复杂因素在很大程度上制约着卫星导航、航天测控及雷达探测等通信产品的应用范围。随着众多军民用通信产品对其性能需求的不断提升，如何在高动态、弱信号等复杂环境下实现载波频率的精确快速估计已成为研究焦点，目前国内外广大科研学者和机构对其进行了持久深入的研究。

国外学者从上世纪已开始研究并给出多种解决办法，Hurd W等采用最大似然估计(MLE)算法得到多普勒频移估计值，适用于高动态环境，具有较高精度[3]；Kumar R等提出了自适应最小均方误差(ALS)估计算法，测量精度高，能够适应快速变化的频率[4]；Jack M等采用信号包络的平方偏差实时估计多普勒频移，此算法的复杂度较高但可适用于低信噪比环境[5]。

国内学者开展研究较晚但也取得了大量研究成果，杨昂研究了高动态环境下各种经典的多普勒频移估计算法，并在此基础上提出了改进的多普勒频移估计算法[6]；黄富彪等结合类Rife频率修正算法提出了适用于数字方法实现的卫星信号载波频率高精度估计算法[7]；易辉等通过非线性变换将BPSK信号转化点频信号并采用相位差校正法提高了估计精度[8]；郑纪民等对非线性变换后的信号进行正交下变频和2级DFT，精确估计了MPSK信号的载波频率[9]。

上述各研究均只对多普勒频移进行了高精度快速估计，并未考虑复杂环境下的多普勒变化率，且难以在动态范围、信号强弱、搜索时间、测频精度和资源消耗等因素之间达到均衡，为了克服上述现存缺陷，本文提出了一种改进的载波频率估计算法。该算法在快速傅里叶变换前对降采样数据进行载波频率的双重补偿，实现多普勒频移、多普勒变化率的双重估计；在不同的测频状态，根据信号调制类型进行模式识别与控制完成载波恢复，可提高载波频率的搜索范围和测量精度；在粗测频、精测频状态对第一级降采样数据频率预补偿后进行第二级降采样处理，可适当缩短测频状态的搜索时间。

1 算法原理

改进的载波频率估计算法的具体原理结构如图1所示，其主要由主控单元(1)、第一级降采样(2)、频率补偿(4)、第二级降采样(5)、模式识别与控制(7)、快速傅里叶变换(8)和峰值搜索频率解算(9)等单元组成。

图1 改进的载波频率估计算法原理结构图

下面分别对该算法原理结构各功能模块的工作原理进行详细说明[10-11]：

1.1 主控单元

主控单元控制整个测频状态，根据输入信号的测频范围和调制类型确定各测频状态中降采样频率、存储容量、频率补偿、模式识别、峰值搜索和频率解算等状态参数，控制其他功能模块的数据传输，通过资源复用实现通用化。

1.2 降采样单元

两级降采样单元均采用ID积分滤波器结构，主要包括DDS、积分清零和数据截位等3个部分，根据主控单元的测频状态确定积分频率控制字，由DDS产生清零脉冲，对输入数据进行积分清零操作以完成ID积分滤波，最后将积分滤波数据有效截位输出至存储单元。

1.3 频率补偿单元

频率补偿单元根据每个频率细分槽的多普勒频移fdopl、多普勒变化率frate得到多普勒频移控制字Kdopl和多普勒变化率控制字Krate，分别通过地址映射、查表生成2路本地载波，经过复乘运算后得到本地复合载波，然后将采样数据和本地复合载波进行复乘运算，对采样数据完成多普勒频移、多普勒变化率的双重补偿，其具体结构如图2所示。

图2 频率补偿单元结构图

图2中多普勒频移控制字Kdopl、多普勒变化率控制字Krate可由下式计算得到，

其中，fs是降采样数据的采样频率。

1.4 模式识别与控制单元

模式识别与控制单元根据主控单元的测频状态和输入信号的调制类型确定控制状态，在初测频状态，将输入信号统一划分为单频模式进行处理；在粗测频、精测频状态，根据输入信号的调制类型完成信号的载波恢复，其中将BPSK信号划分为二倍频模式，将QPSK，SQPSK和UQPSK信号等划分为四倍频模式。

二倍频模式为了恢复BPSK信号被抑制的载波，在进行快速傅里叶变换之前需对BPSK信号进行平方运算完成载波恢复，而平方运算带来了处理信号的信噪比恶化[12]，其恶化损失可由下式计算得到

其中，SNRi为输入信噪比。

同理，四倍频模式为了恢复QPSK，SQPSK以及UQPSK信号等被抑制的载波，需对QPSK，SQPSK和UQPSK信号等进行四次方运算完成载波恢复，而四次方运算带来了处理信号更严重的的信噪比恶化，通过仿真可得到其恶化损失。

综上分析，二倍频模式和四倍频模式带来的处理信号信噪比恶化损失曲线如图3所示

图3 处理信号信噪比恶化损失曲线图

1.5 峰值搜索与频率解算单元

峰值搜索与频率解算单元将快速傅里叶变换输出的处理结果经过非相干积分累加后进行峰值判决，单频模式采用频谱重心法，其余模式均采用频谱比较法，得到每个频率槽的峰值，然后对当前测频状态所有频率槽进行峰值比较，找到最大峰值以及其所对应的频率槽，完成当前测频状态采样数据的频率搜索，最后计算并输出所得的多普勒频移、多普勒变化率。

2 流程控制

主控单元将整个测频过程分为初测频、粗测频、精测频等3个状态，根据输入信号的测频范围、调制类型等确定各测频状态的状态参数和控制各功能模块间的数据传输，具体流程控制如图4所示，图中详细说明了各测频状态的工作原理。

图4 改进的载波频率估计算法流程控制图

2.1 初测频状态

初测频状态根据测频范围对输入信号完成第一级降采样处理，对第一级降采样数据进行多普勒变化率的频率补偿后以单频模式进行测频，得到粗略的多普勒频移测量值。

设最大多普勒频移为fdopl_max，最大多普勒变化率为frate_max，系统工作时钟为fclk，则第一级降采样频率为

fs≥4fdopl_max

为了检测输入信号被抑制的载波，第一级降采样数据量为Nfull，则第一级降采样时间为

ts=Nfull/fs

快速傅里叶变换测频点数为COH，则多普勒频移、多普勒变化率的测量精度为

初测频状态对第一级降采样数据以单频模式进行测频，测频精度较差，则频率补偿的多普勒变化率步长可为frate_max≤4frate_res1、其载波频率搜索轮次分别为

slotfull≥2frate_max/frate_step1

则初测频状态的搜索时间为

tfull=ts+Nfull/fclk×slotfull

2.2 粗测频状态

粗测频状态采用初测频状态得到的多普勒频移测量值对第一级降采样数据完成频率预补偿处理后进行第二级降采样处理，对第二级降采样数据完成多普勒频移、多普勒变化率的双重频率补偿后进行输入信号的模式识别，得到较精确的多普勒频移测量值和粗略的多普勒变化率测量值。

粗测频状态第二级降采样的ID积分滤波点数为ID2coarse，则第二级降采样数据量为，

Ncoarse=Nfull/ID2coarse

则其第二级采样时间为

tid 2=Nfull/fclk

粗测频状态的多普勒频移、多普勒变化率的测量精度为

式中，K为模式识别的控制指数，二倍频模式时为2，四倍频模式时为4。

粗测频状态根据输入信号的调制类型确定模式识别的控制状态，其测频精度较好，则频率补偿的多普勒变化率步长可为frate_step2≤2frate_res2，其载波频率搜索轮次为

slotcoarse≥4frate_step1/frate_step2

则粗测频状态的搜索时间为

tcoarse=tid2+Ncoarse/fclk×slotcoarse

2.3 精测频状态

精测频状态采用粗测频状态得到的载波频率测量值对第一级降采样数据完成较精确的频率预补偿处理，并重复粗测频状态的操作流程，从而得到更精确的多普勒频移和多普勒变化率测量值。

精测频状态第二级降采样的ID积分滤波点数为ID2fine(ID2fine>ID2coarse)，其测频状态的操作流程同上所述，可得精测频状态的第二级降采样数据量为Nfine、多普勒频移测量精度fdopl_res3、多普勒变化率测量精度frate_res3及多普勒变化率步长frate_step3，频率预补偿的多普勒变化率为概略测量值，需进一步提高其测量精度，则其载波频率搜索轮次为

slotfine≥2frate_step2/frate_step3

则精测频状态的搜索时间为

tfine=tid 2+Nfine/fclk×slotfine

3 试验与分析

为检验该算法的测频性能，采用某FPGA+DSP平台进行测频试验，其参数设置如下：系统工作时钟fclk为120MHz，数据调制速率Rb为10kbps，多普勒频移范围fdopl_max为±250kHz，多普勒变化范围frate_max为±25kHz/s。

3.1 测频精度

若接收信号为BPSK信号时，在初测频状态，第一级降频率采样为fs=1000kHz，为了检测接收信号中已完全淹没在背景噪声的载波频率，需延长积分时间以增加积分增益，则采样时间ts为0.256s，FFT测频点数COH选取2048，其载波频率搜索轮次slotfull为30；在粗测频状态，第二级降采样的ID积分滤波点数ID2coarse为4，其载波频率搜索轮次slotcoarse为100；在精测频状态，第二级降采样的ID积分滤波点数ID2fine为16，其载波频率搜索轮次slotfine为40；模式识别与控制单元在初测频状态采用单频模式，在粗测频、精测频状态采用二倍频模式。则根据上述理论分析可知整个测频过程的搜索时间、多普勒频移精度、多普勒变化率精度分别约为0.386s、±10.0Hz和±30Hz/s。

若接收信号为QPSK信号时，初测频、粗测频状态流程同上，在精测频状态，第二级降采样的ID积分滤波点数ID2fine为8，其载波频率搜索轮次slotfine为60；模式识别与控制单元在初测频状态采用单频模式，在粗测频、精测频状态采用四倍频模式。则根据上述理论分析可知整个测频过程的搜索时间、多普勒频移精度及多普勒变化率精度分别约为0.397s、±10.0Hz和±30Hz/s。

3.2 检测概率

根据参数设置，在不同频率动态范围、不同信噪比条件下对该算法的检测概率Pd进行试验分析，以多普勒频移误差Δfdopl≤20Hz、多普勒变化率误差Δfrate≤50Hz/s为成功检测条件，可得试验结果如图5～6所示。

图5 BPSK信号测频试验检测概率图

图6 QPSK信号测频试验检测概率图

图5为BPSK信号测频试验的检测概率统计，当输入信噪比Eb/N0≥-3.0dB时，该算法在不同频率动态范围的检测概率Pd≥90%；图6为QPSK信号测频试验的检测概率统计，当输入信噪比Eb/N0≥4.0dB时，该算法在不同频率动态范围的检测概率Pd≥90%；当频率动态逐渐增大时，相干积分的扇贝效应和FFT运算的栅栏效应让输出增益产生衰减，相应地降低检测概率，但该算法在高动态环境下仍具有明显的性能优势。

4 结论

研究了高动态环境下载波频率的精确估计问题，在总结国内外研究现状的基础上，提出了一种改进的载波频率精确估计算法，详细说明该算法的原理结构和流程控制，实现了多普勒频移、多普勒变化率的同时精确估计，可适应BPSK、QPSK、SQPSK和UQPSK等多种不同调制类型和调制速率的输入信号，并完成该算法在高动态低信噪比条件下的试验实验。通过试验分析可以看出，该算法具有更大的搜索范围、更高的估计精度及更好的测量性能，在复杂环境下有更明显的应用优势。