1.原始文档
《用于正弦波频率估计的修正I-Rife算法》,王哲文,2024
DOI: 10. 16337/j. 1004‑9037. 2024. 02. 019
1.1 这篇论文所属的自科基金U21A20500:近5年所承担的重要科研项目表-智能感知系统与安全教育部重点实验室(湖北大学)
2.问题速记
2.1 摘要部分
- 量化频点的因素
- 为啥接近量化频谱点处会有更大的误差?他不应该更小吗?
- 文章的中心是Rif算法的性能和误差分析以及针对性改进
- 0.5 pt???how?
- 论文复杂度与I-Rif相当
- I-Rif算法
- 移频因子
- Rif算法区别
- 修正的基本逻辑
- 克拉美-罗下界
- 均方根误差
- Cramer-Rao lower Bound
- 关键词:
- 频率估计Freq estimation
- 细化谱 spectrum subdivision
- CRLB
- 中国分类TN957.51
- 文献标志码:A
- Sinusoid Wave
2.2 引言部分
- 调频连续波Frequency modulation continuous waveFMCW)雷达
- 混频
- 中频信号
- 距离速度
3.频率估计的途径
- Abatzoglou 最大似然估计最大似然估计(Maximum likelihood ML)算法
- FFT 频谱泄露、栅栏效应
- 信噪比较低时误差很大
- 基于FFT的插值算法
- 三谱线插值法
- 迭代差值法
- 抛物线内插
- 双线性对称频率内插
- Rif算法★
- 异频相位处理
- 全相位时移相位差频率估计
- 能量重心矫正法
- 相位差分算法
- 频谱细化方法的Zoom‑FFT算法★
- 线性调频变换CZT算法★
ML慢而精确。RIf快,Zoom-FFT, CZT性能很好但无法满足实时解析场景。
3.1 Rif算法
- Rife1970年提出,利用幅度最高谱线与相邻谱线,差值,确定真实频率。
- 当信噪比过低或其真实频率接近量化频率时,易出现插值方向的错误从而引起较大的
误差。 - Quinn提出用最大谱线,次大谱线的之比的实部来替代幅值之比
- 邓帧淼M-Rif,引入移频原理,将估计的频率移到两个量化频率中心,再差值,移频量固定1/3
- 有时需要二次移频
- 王宏位I-Rif,频谱细化,移频
- 孙宏君P-Rif,相角
- A-P-Rif幅值相角判据的改进频移门限值来确定是否使用相角判据,得到了精度更高的估计性能
- 叶茂Z-Rif
- 在A-P-Rif上,幅值相角规律
- NianA-Rif
- 基于估计
3.2 Rif算法的实际性能与误差
4.Python重现Rif算法和上述现象
4.1 准备信号源
我的基本信号处理对象是振动信号,这里生成了一个冲击性振动信号,很简单:
叠加了幅度为信号-20db的高斯噪声。生成代码参见附录B 冲击信号模拟源
4.2 使用不同的采样率去采集信号源,生成多组采样数据
附录A 一些可用的三方库
1. specutils
这是一个用于天文观测雷达信号的解析库,作为一个完整领域的信号分析库,可以用来对信号处理的可能问题和算法有一个大致的第一印象。
download from: specutils · PyPId
doc:Specutils Documentation — specutils v1.15.1.dev0+gce3f43c.d20240501
附录B 源码
1. 冲击信号模拟源
#!/usr/bin/env python3
# -*- coding: utf-8 -*-
# 获取当前脚本文件所在目录的父目录,并构建相对路径
import os
import sys
current_dir = os.path.dirname(os.path.abspath(__file__))
project_path = os.path.join(current_dir, '..')
sys.path.append(project_path)
sys.path.append(current_dir)
import numpy as np
import matplotlib.pyplot as plt
from matplotlib import rcParams
from scipy.signal import get_window# 设置中文字体
rcParams['font.sans-serif'] = ['SimHei'] # 指定中文字体
rcParams['axes.unicode_minus'] = False # 解决负号问题# 参数设置
fs = 10000 # 采样频率 (Hz)
f_signal = 181.5 # 信号频率 (Hz)
rms_signal = 5e-3 # 信号 RMS (m/s)
duration = 1 # 信号持续时间 (秒)
pulse_duration = 1 / (100 * f_signal) # 脉冲持续时间 (秒)# 生成时间序列
t = np.arange(0, duration, 1/fs)# 生成冲击性振动信号(脉冲信号)
pulse_start = int(fs * (duration / 2 - pulse_duration / 2))
pulse_end = int(fs * (duration / 2 + pulse_duration / 2))
signal = np.zeros_like(t)
signal[pulse_start:pulse_end] = rms_signal * np.sqrt(2)# 生成高斯噪声
rms_noise = rms_signal * 0.1
noise = rms_noise * np.random.randn(len(t))# 将信号和噪声相加
signal_with_noise = signal + noise# 绘图
plt.figure(figsize=(12, 6))
plt.subplot(1, 2, 1)
plt.plot(t, signal, label='冲击性振动信号')
plt.title('冲击性振动信号')
plt.xlabel('时间 (秒)')
plt.ylabel('幅值')
plt.legend()plt.subplot(1, 2, 2)
plt.plot(t, signal_with_noise, label='含噪声信号', color='red')
plt.title('含噪声的冲击性振动信号')
plt.xlabel('时间 (秒)')
plt.ylabel('幅值')
plt.legend()plt.tight_layout()
plt.show()