一种基于CZT的频率估计方法与流程

一种基于czt的频率估计方法
技术领域
1.本发明涉及频率估计领域，更具体地涉及一种czt的频率估计方法。


背景技术：

2.信号的频率估计是工程应用中常常存在的问题，许多场景需要对信号的频率进行精确估计。例如，在频率调制连续波(fmcw)雷达测距系统中，其可以根据发送波和反射波的差频信号来得到相关的距离信息，即通过估计差频信号的频率来得到估计距离，因此对于差频信号的频率估计准确程度直接影响到了估计距离的准确程度，差频信号频率的测量精度直接决定了测距的精度。
3.为了提升频率估计的精度，许多算法在fft变换(傅里叶变换)的基础上进行了改进，例如选带傅里叶变换(zoom-fft)、rife算法、线性调频z变换(czt)等。zoom-fft算法(赵强,侯孝民,廉昕.基于zoom-fft的改进rife正弦波频率估计算法[j].数据采集与处理,2017,32(04):731-736.doi:10.16337/j.1004-9037.2017.04.009)对于信号进行移频，再进行降采样，通过降采样提升了频谱分辨率，从而提升了频率估计的精度；rife算法利用频谱幅值最大的谱线及其相邻谱线中幅值较大谱线来计算实际频率与最大谱线对应的频率的偏差值，从而提升频率估计的精度。但是，rife算法存在估计精度受噪声影响大、偏差值较小时估计误差较大的问题。孙宏军等(孙宏军,王小威.基于幅值-相角判据的修正rife正弦波频率估计算法[j].天津大学学报(自然科学与工程技术版),2018,51(08):810-816)通过设定频移门限值，对频移因子满足门限值的信号使用rife算法，否则便使用相角来估计频移因子，提升了算法的性能。quinn(quinn b g.estimation of frequency，amplitude and phase from the dft of a time series[j].ieee trans on signal processing，1997，45(3)：814-817)提出利用fft变换主瓣内次大值与最大值之比的实部代替幅值之比的频率估计方法，减少了由于频移因子方向判断错误引起的误差。czt算法(欧阳鑫信,贺青.基于czt的多子带信号时差估计方法[j].现代雷达,2019,41(02):31-34.doi:10.16592/j.cnki.1004-7859.2019.02.007)对于部分频段进行细化，相较于fft算法具有更高的频率估计精度，但是需要大量的计算量。刘上(刘上,朱国富,王玲,等.lfmcw雷达高精度测距相位差改进算法[j].雷达科学与技术,2021,19(03):287-291+297)提出了基于czt的相位差算法，对于信号进行czt变换，求出信号相位，利用相位差来求解信号的频率，解决了相位差算法受噪声影响大缺陷，提升了频率估计的精度。
[0004]
综上所述，现有的频率估计算法均存在受噪声影响大、计算量大等缺陷，无法满足目前工业领域对于频率估计的精度高、运算快的要求。


技术实现要素：

[0005]
本发明的目的在于提供一种基于czt的频率估计方法，其计算量低、抗噪性能好。
[0006]
本发明提供一种基于czt的频率估计方法，包括以下步骤：
[0007]
s1：对待估计信号进行采样，并对采样后的信号进行fft变换，获得最大谱线频率
索引；
[0008]
s2：根据fft的最大谱线频率索引确定频率细化区间，对采样后的信号进行czt变换；
[0009]
s3：对czt变换后的频谱进行搜索，获得最大谱线的幅值及其左侧和右侧谱线的幅值；
[0010]
s4：根据最大谱线右侧谱线的幅值与最大谱线的幅值的比值以及最大谱线左侧谱线的幅值与最大谱线的幅值的比值确定待估计信号的瞬时频率与待估计信号在czt变换后的频率的误差；
[0011]
s5：根据所述误差确定待估计信号的瞬时频率。
[0012]
进一步地，所述待估计信号为fmwc雷达的差频信号s(t)，fmwc雷达的差频信号s(t)满足如下关系式：
[0013]
s(t)＝a1e
jψ(t)
[0014]
其中，a1为待估计信号的振幅，j为虚数单位，ψ(t)为待估计信号的相位。
[0015]
进一步地，采样后的信号满足如下关系式：
[0016][0017]
其中，fc为待估计信号的瞬时频率，fs为采样频率，n为经过采样的信号的点数，为采样后的信号的相位。
[0018]
进一步地，步骤s2中，频率细化区间为其中，fs为采样频率，q用于控制频率细化区间的大小，m0为fft的最大谱线频率索引。
[0019]
进一步地，czt变换后的频谱中的最大谱线的幅值为：
[0020][0021]
其中，m1为czt的最大谱线的频率索引，m为频率细化倍数，δ为待估计信号的瞬时频率与czt变换后的频率的误差。
[0022]
进一步地，czt变换后的频谱中的最大谱线左侧谱线的幅值为：
[0023][0024]
进一步地，czt变换后的频谱中的最大谱线右侧谱线的幅值为：
[0025][0026]
进一步地，待估计信号的瞬时频率与待估计信号在czt变换后的频率的误差δ满足如下关系式：
[0027]
[0028]
其中，a1为最大谱线右侧谱线的幅值与最大谱线的幅值的比值，a2为最大谱线左侧谱线的幅值与最大谱线的幅值的比值。
[0029]
进一步地，所述待估计信号的瞬时频率满足如下关系式：
[0030]
fc＝(m
0-q)δf1+(m1+δ)δf2[0031]
其中，m0为fft的最大谱线频率索引，q用于控制频率细化区间的大小，δf1为fft的频率分辨率，m1为czt的最大谱线的频率索引，δ为待估计信号的瞬时频率与待估计信号在czt变换后的频率的误差，δf2为czt的频率分辨率。
[0032]
本发明的基于czt的频率估计方法，先使用fft对待估计信号进行粗估计，再使用czt变换对粗估计频率附近的频段进行细化，利用细化后频谱的最大谱线及其左右谱线与误差之间的函数关系，求解出误差的表达式，从而提升频率估计的精度。本发明能够以较低的细化倍数和较小的计算量，获得较高精度的频率估计结果，从而准确测量目标与雷达之间的距离。
附图说明
[0033]
图1为根据本发明实施例的基于czt的频率估计方法的流程图。
具体实施方式
[0034]
下面结合附图，给出本发明的较佳实施例，并予以详细描述。
[0035]
在fmcw雷达发射系统中，fmcw雷达的发射信号为：
[0036][0037]
其中a
tx
是发射信号的振幅，f0是初始频率，是线性调频的斜率，b是雷达带宽，tc代表扫频周期，初始相位为
[0038]
雷达的接收信号为：
[0039][0040]
其中，a
rx
是接收信号的振幅，φ
rx
是接收信号的相位，是接收信号相对于发射信号的时延，r表示目标与雷达之间的距离，c为电磁波的传播速度。
[0041]
将发射信号与接收信号进行混频操作，将混频过后的信号通过低通滤波器，从而获得同向分量：
[0042][0043]
将发射信号偏移90
°
后与接收引号进行混频，再将信号通过低通滤波器后，获得正交分量：
[0044][0045]
将同向分量与正交分量结合获得复信号，即为差频信号s(t)。由于τ的值较小，可
以忽略τ2的项，因此fmwc雷达的差频信号s(t)为：
[0046][0046][0047]
信号的瞬时频率由相位求得：
[0048][0049]
其中，是差频信号的相位。
[0050]
根据可以推导出雷达物体之间的距离与瞬时频率成正比：
[0051][0052]
因此，为了得到准确度高的距离，必须得到准确度高的差频信号的瞬时频率。
[0053]
如图1所示，本发明提供一种基于czt的频率估计方法，包括以下步骤：
[0054]
s1：对待估计信号进行采样，并对采样后的信号进行fft变换，获得最大谱线频率索引；
[0055]
在本实施例中，待估计信号为fmwc雷达的差频信号s(t)。
[0056]
如上文所述，fmwc雷达的差频信号s(t)满足如下关系式：
[0057]
s(t)＝a1e
jψ(t)
[0058]
其中，a1为待估计信号的振幅，j为虚数单位，ψ(t)为待估计信号的相位。
[0059]
以fs的频率对所述待估计信号进行采样，得到采样后的信号：
[0060][0061]
其中，fc为待估计信号的瞬时频率，fs为采样频率，为采样后的信号的相位，n是经过采样的信号的点数。
[0062]
对采样后的信号进行fft运算，得到其频谱：
[0063][0064]
对频谱的幅值进行搜索，幅值最大处即为频率的粗估计值
[0065][0066]
其中，m0为最大谱线频率索引，即幅值最大处对应的频率索引值，其中，m0为最大谱线频率索引，即幅值最大处对应的频率索引值，即为频率分辨率。
[0067]
s2：根据fft的最大谱线频率索引确定频率细化区间，对采样后的信号进行czt变
换；
[0068]
频率细化区间(f1，f2)以为中心，其两端为：
[0069][0070]
其中，fs为采样频率，q用于控制频率细化区间的大小，可根据需要进行设置，例如可设置为q＝1，2，3或其他值等，m0为fft的最大谱线频率索引。
[0071]
信号的czt变换表示为：
[0072][0073]
其中zk＝aw-k
，是起始采样角度，是两相邻采样点之间的角度，m是频率细化倍数。
[0074]
在本实施例中，a0＝1、w0＝1，信号s(n)的czt变换为：
[0075][0076][0077][0078][0079]
s3：对czt变换后的频谱进行搜索，获得最大谱线的幅值及其左右两根谱线的幅值；
[0080]
对czt变换后的频谱进行搜索后，求得其最大谱线的频率索引m1，从而估计出频率：
[0081][0082]
其中，即为czt变换的频率分辨率。
[0083]
信号的瞬时频率与czt细化后频谱的采样点之间往往还存在一定的误差δ，瞬时频率fc表示为：
[0084]
fc＝(m
0-q)δf1+(m1+δ)δf2ꢀꢀ
(15)
[0085]
其中，m0为fft的最大谱线的频率索引，q用于控制频率细化区间的大小，δf1为fft的频率分辨率，m1为czt的最大谱线的频率索引，δ为待估计信号的瞬时频率与待估计信号在czt变换后的频率的误差，δf2为czt的频率分辨率，δ∈[-0.5，0.5]。
[0086]
czt变换后的频谱中的最大谱线的幅值|s
czt
(m1)|为：
[0087][0088]
其中，m1为czt的最大谱线的频率索引，m为频率细化倍数，δ为待估计信号的瞬时
频率与czt变换后的频率的误差。
[0089]
czt变换后的频谱中的最大谱线左侧谱线的幅值|s
czt
(m
1-1)|为：
[0090][0091]
czt变换后的频谱中的最大谱线右侧谱线的幅值|s
czt
(m1+1)|为：
[0092][0093]
s4：根据最大谱线右侧谱线的幅值与最大谱线的幅值的比值以及最大谱线左侧谱线的幅值与最大谱线的幅值的比值确定待估计信号的瞬时频率与待估计信号在czt变换后的频率的误差δ；
[0094]
|s
czt
(m1+1)|和|s
czt
(m1)|的比值记为a1：
[0095][0096]
n为经过采样的信号点数，一般为2的整数次幂，例如512、1024、2048等，对于较大的mn(mn≥256)时，a1可以表示为：
[0097][0098]
整理后可得：
[0099][0100]
|s
czt
(m
1-1)|和|s
czt
(m1)|的比值记为a2：
[0101][0102]
对于较大的mn(mn≥256)时，a2可以表示为：
[0103][0104]
整理后可得：
[0105][0106]
将式(21)和式(24)相减，可得：
[0107][0107][0108]
利用三角函数的和差化积公式对于式(25)进行化简可得：
[0109][0110]
对于式(26)进行求解，可得待估计信号的瞬时频率与待估计信号在czt变换后的频率的误差δ的表达式：
[0111][0112]
s5：根据待估计信号的瞬时频率与czt变换后的频率的误差确定待估计信号的瞬时频率。
[0113]
将误差δ代入式(15)中即可获得待估计信号的瞬时频率fc。
[0114]
然后再将瞬时频率代入式(7)中即可得到fmwc雷达的距离。
[0115]
本发明实施例提供的基于czt的频率估计方法，先使用fft对待估计信号进行粗估计，再使用czt变换对粗估计频率附近的频段进行细化，利用细化后频谱的最大谱线及其左右谱线与误差δ之间的函数关系，求解出误差δ的表达式，从而提升频率估计的精度。本发明能够以较低的细化倍数和较小的计算量，获得较高精度的频率估计结果，从而准确测量目标与雷达之间的距离。
[0116]
以上所述的，仅为本发明的较佳实施例，并非用以限定本发明的范围，本发明的上述实施例还可以做出各种变化。即凡是依据本发明申请的权利要求书及说明书内容所作的简单、等效变化与修饰，皆落入本发明专利的权利要求保护范围。本发明未详尽描述的均为常规技术内容。