CN110471018B - 一种频谱校正方法 - Google Patents
一种频谱校正方法 Download PDFInfo
- Publication number
- CN110471018B CN110471018B CN201910888489.3A CN201910888489A CN110471018B CN 110471018 B CN110471018 B CN 110471018B CN 201910888489 A CN201910888489 A CN 201910888489A CN 110471018 B CN110471018 B CN 110471018B
- Authority
- CN
- China
- Prior art keywords
- value
- frequency
- amplitude
- corrected
- phase
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Expired - Fee Related
Links
- 238000000034 method Methods 0.000 title claims abstract description 29
- 238000001228 spectrum Methods 0.000 title description 23
- 230000003595 spectral effect Effects 0.000 claims abstract description 32
- 238000005070 sampling Methods 0.000 claims description 14
- 238000004422 calculation algorithm Methods 0.000 claims description 3
- 238000006243 chemical reaction Methods 0.000 abstract description 3
- 238000010586 diagram Methods 0.000 description 9
- 238000010183 spectrum analysis Methods 0.000 description 7
- 238000004364 calculation method Methods 0.000 description 3
- 230000005540 biological transmission Effects 0.000 description 2
- 239000006185 dispersion Substances 0.000 description 2
- 230000000694 effects Effects 0.000 description 2
- 238000000819 phase cycle Methods 0.000 description 2
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000015572 biosynthetic process Effects 0.000 description 1
- 238000004883 computer application Methods 0.000 description 1
- 238000003745 diagnosis Methods 0.000 description 1
- 238000005516 engineering process Methods 0.000 description 1
- 238000000605 extraction Methods 0.000 description 1
- 238000007670 refining Methods 0.000 description 1
- 239000004576 sand Substances 0.000 description 1
- 230000001629 suppression Effects 0.000 description 1
- 238000003786 synthesis reaction Methods 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01R—MEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
- G01R23/00—Arrangements for measuring frequencies; Arrangements for analysing frequency spectra
- G01R23/16—Spectrum analysis; Fourier analysis
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01R—MEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
- G01R35/00—Testing or calibrating of apparatus covered by the other groups of this subclass
- G01R35/005—Calibrating; Standards or reference devices, e.g. voltage or resistance standards, "golden" references
Landscapes
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Mathematical Physics (AREA)
- Radar Systems Or Details Thereof (AREA)
- Complex Calculations (AREA)
Abstract
本发明利用变窗长方法,求解得到系列幅值、频率和相位。根据主瓣中心与最近谱线的距离随窗长变化的原理,选取不同窗长的窗函数进行截断,获得多组FFT变换结果,则最近谱线相对于主瓣中心的位置将会存在规律性的接近和远离。理论上,只要遍历各种窗函数长度,最近谱线与主瓣中心的距离就可以无限接近。此时,最大幅值对应的谱线位置即与主瓣中心最近,即最接近幅值真值,并且主瓣中心与最近谱线的距离为零均值的分布。因此,可以将系列频率值的均值作为校正的频率值,将系列幅值的最大值作为校正的幅值,将系列相位的均值作为校正的相位值。本发明频率值、幅值和相位值三者的校正过程相互独立,任一一者的校正误差,不会传递并影响其他两者。
Description
技术领域
本发明涉及数字信号处理领域,具体为一种频谱校正方法。
背景技术
离散频谱分析实现了信号处理从时域到频域的转变,推动了计算机应用技术的发展,在机械、电子、仪器仪表等领域得到了广泛应用。频谱分析准确性对工程应用具有十分重要的意义,然而计算机难以处理实际连续信号,需对信号进行截断和离散,当采样长度不是信号周期整数倍时,将会造成频谱泄漏。
频谱泄漏可分为长程谱泄漏和短程谱泄漏。从连续傅里叶变换到DFT,需要经过时域离散、数据截断以及频域离散过程。信号的时域离散化导致频域周期化,根据Nyquist采样定理,采样频率应大于信号最高频率的两倍,否则会产生假频造成频率混叠现象。计算机处理信号的长度总是有限的,必然要对信号进行截断,若截断信号长度N是信号周期的非整数倍,即非同步采样,在经频域离散造成时域周期化后,截断处会因吉布斯现象产生振荡,这种不连续状态将导致长程谱泄漏,泄漏程度与窗谱的旁瓣特性密切相关。在非同步采样中,频域离散还会使得信号真实频率f0位于两条离散谱线k和k+1之间,造成短程谱泄漏,这种现象称为栅栏效应,相邻谱线间的宽度Δf=fs/N为频率分辨率,直接影响频谱分析精度。
这种频谱泄漏现象会影响频谱分析的准确性,给各类工程应用造成障碍。例如旋转机械振动响应信号包含了转频及其倍频成分,在机械故障诊断中需要研究各倍频轴心轨迹的特征,而短程谱泄漏造成的以DFT频谱分析为基础的频谱相位和幅值的不准确,影响了各频率成分的提取和提纯轴心轨迹的合成。
为减小短程谱泄漏造成的频谱分析误差,研究者提出了多种短程谱泄漏抑制方法,如频谱细化法、插值法、能量重心法、相位差法、三角形法等。上述方法能有效抑制频谱泄漏,进行相对精确的频谱分析和参数估计。但是频谱细化方法存在计算复杂度高或细化频带窄的缺点,而上述其他方法都基于主瓣中心与最近谱线的距离Δx进行频谱校正,Δx的计算误差将影响相位幅值和频率的校正,存在误差传递的问题。
发明内容
针对现有技术中存在的问题,本发明提供一种频谱校正方法,能够实现频谱幅值、频率和相位独立校正,不需要计算主瓣中心与最近谱线的距离Δx,避免了校正误差传递。
本发明是通过以下技术方案来实现:
一种频谱校正方法,包括如下步骤,
步骤1,对采样频率为fs的离散正弦信号S(n)=[s0 s1 s2 … sk],以长度为N的窗函数W0截取信号,获取第一短时信号y0(n)=[s0 s1 … sN-1];其中,k为信号数据点数,1<k<∞,1<N<k;
步骤2,通过FFT算法求解短时信号y0(n),得到其频域离散复数序列通过的最大绝对值序号j0和频率分辨率Δf0计算得到信号的第一近似频率值fa0,并将序号j0处的绝对值记为第一近似幅值Aa0,将序号j0处的虚部与实部反正切值记为第一相位值
步骤4,将m个频率值的均值作为校正的频率值,将m个幅值的最大值最为校正的幅值,将m个相位的均值作为校正的相位值,完成对采样频率为fs的离散正弦信号的频谱校正。
优选的,步骤2中,将最高谱线对应的频率值记为频率真值的近似值,即fa0=(j0-1)Δf0。
与现有技术相比,本发明具有以下有益的技术效果:
本发明利用变窗长方法,求解得到系列幅值、频率和相位。根据主瓣中心与最近谱线的距离随窗长变化的原理,选取不同窗长的窗函数进行截断,获得多组FFT变换结果,则最近谱线相对于主瓣中心的位置将会存在规律性的接近和远离。理论上,只要遍历各种窗函数长度,最近谱线与主瓣中心的距离就可以无限接近。此时,最大幅值对应的谱线位置即与主瓣中心最近,即最接近幅值真值,并且主瓣中心与最近谱线的距离为零均值的分布。因此,可以将系列频率值的均值作为校正的频率值,将系列幅值的最大值作为校正的幅值,将系列相位的均值作为校正的相位值。本发明频率值f、幅值A和相位值三者的校正过程相互独立,任一一者的校正误差,不会传递并影响其他两者。校正过程与窗函数的具体表达式无关,可以适应多种加窗信号。
附图说明
图1为本发明实例中所述的频谱校正方法流程图。
图2为本发明实例中所述的加窗获取短时信号示意图。
图3为本发明实例中所述的局部频带内最高谱线位置示意图。
图4为本发明实例中所述的实例信号时域波形。
图5为本发明实例中所述实例信号的短时信号频谱图。
图6为本发明实例中所述的第一频率成分幅值序列的分布图。
图7为本发明实例中所述的第一频率成分频率序列的分布图。
图8为本发明实例中所述的第一频率成分相位序列的分布图。
图9为本发明实例中所述的第二频率成分幅值序列的分布图。
图10为本发明实例中所述的第二频率成分频率序列的分布图。
图11为本发明实例中所述的第二频率成分相位序列的分布图。
具体实施方式
下面结合具体的实施例对本发明做进一步的详细说明,所述是对本发明的解释而不是限定。
本发明一种频谱校正方法,能够独立校正频谱幅值、频率和相位独立校正,包括如下步骤,
步骤1,对采样频率为fs的离散正弦信号S(n)=[s0 s1 s2 … sk],(k为信号数据点数,1<k<∞)以长度为N(1<N<k)的窗函数W0截取信号,获取第一短时信号y0(n)=[s0 s1 …sN-1];
步骤2,通过FFT算法求解短时信号y0(n),得到其频域离散复数序列进而,通过的最大绝对值序号j0和频率分辨率Δf0计算得到信号的第一近似频率值fa0,并将j0序号处的绝对值记为第一近似幅值Aa0,将j0序号处的虚部与实部反正切值记为第一相位值
步骤4,将m个频率值的均值作为校正的频率值,将m个幅值的最大值最为校正的幅值,将m个相位的均值作为校正的相位值,完成对采样频率为fs的离散正弦信号的频谱校正。
具体的,包括以下步骤:
1)以长度为N的窗函数W0截取采样频率为fs的离散正弦信号S(n)=[s0 s1 s2 …sk],获取短时信号y0(n)=[s0 s1 … sN-1];
本发明所述的方法实施时的基本流程如图1所示,改变窗长为N+i,i=1,2,…m,求解得到系列幅值、频率和相位,利用其分布特征,将系列频率值的均值作为校正的频率值,将系列幅值的最大值最为校正的幅值,将系列相位的均值作为校正的相位值。
首先,如图2所示,以长度为N的窗函数W0截取采样频率为fs的离散正弦信号S(n)=[s0 s1 s2 … sk],获取短时信号y0(n)=[s0 s1 … sN-1];
本优选实例,以信号为例,其时域波形如图7所示。按照图1所示的基本流程,采样频率为fs=1000Hz,取窗长N=1000,得到第一个短时信号的离散频谱如图8所示,按照上述实施方法,计算得到第一组两个频率成分的幅值、频率和相位值为Aa1=4.92,12.85、fa1=53,79和改变窗长为N+i,i=1,2,…m,m=1000分别获得两个频率成分的m个频率值fa1,fa2,…,fam、幅值Aa1,Aa2,…,Aam和相位值其分布图分别如图6~图11所示。最后,将两个频率成分别分对应的m个频率值的均值作为校正的频率值fc,将m个幅值的最大值最为校正的幅值Ac,将m个相位的均值作为校正的相位值校正值、真值及误差如表1所示。为对比说明,给出了本例中直接FFT计算的最大误差值作为参照。从表1中可以看出,直接FFT变换中,频率成分1的幅值、频率和相位误差依次为0.7519、0.42、1.56749,频率成分2的幅值、频率和相位误差依次为2.2601、0.4047、1.56711。其中幅值和相位误差均较大,难以满足工程应用需求,频率误差较小。而校正后,频率成分1的幅值、频率和相位误差依次为-0.0049、-0.0367、-0.00205,频率成分2的幅值、频率和相位误差依次为-0.0141、-0.0538、-0.00423。校正后,幅值和相位误差显著降低,频率误差也进一步降低。
表1校正结果对比
Claims (1)
1.一种频谱校正方法,其特征在于,包括如下步骤,
步骤1,对采样频率为fs的离散正弦信号S(n)=[s0 s1 s2 … sk],以长度为N的窗函数W0截取信号,获取第一短时信号y0(n)=[s0 s1 … sN-1];其中,k为信号数据点数,1<k<∞,1<N<k;
步骤2,通过FFT算法求解短时信号y0(n),得到其频域离散复数序列通过的最大绝对值序号j0和频率分辨率Δf0计算得到信号的第一近似频率值fa0,并将序号j0处的绝对值记为第一近似幅值Aa0,将序号j0处的虚部与实部反正切值记为第一相位值
步骤4,将m个频率值的均值作为校正的频率值,将m个幅值的最大值最为校正的幅值,将m个相位的均值作为校正的相位值,完成对采样频率为fs的离散正弦信号的频谱校正;
校正的幅值为Ac=max(Aai);
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201910888489.3A CN110471018B (zh) | 2019-09-19 | 2019-09-19 | 一种频谱校正方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201910888489.3A CN110471018B (zh) | 2019-09-19 | 2019-09-19 | 一种频谱校正方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN110471018A CN110471018A (zh) | 2019-11-19 |
CN110471018B true CN110471018B (zh) | 2021-12-24 |
Family
ID=68516332
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201910888489.3A Expired - Fee Related CN110471018B (zh) | 2019-09-19 | 2019-09-19 | 一种频谱校正方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN110471018B (zh) |
Families Citing this family (5)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN111257815B (zh) * | 2020-03-06 | 2022-04-05 | 云南电网有限责任公司电力科学研究院 | 一种高精度频谱校正方法 |
CN111884965A (zh) * | 2020-07-22 | 2020-11-03 | 云南电网有限责任公司电力科学研究院 | 基于全泄漏抑制的频谱校正方法及装置 |
CN112986677B (zh) * | 2021-02-04 | 2022-01-28 | 电子科技大学 | 一种基于SoC动态可配的频谱分析的***及实现方法 |
CN113109622A (zh) * | 2021-04-15 | 2021-07-13 | 南方电网科学研究院有限责任公司 | 一种电网信号频谱的分析方法及*** |
CN116125138B (zh) * | 2023-04-17 | 2023-07-14 | 湖南工商大学 | 基于旋转调节的正弦信号频率快速估计方法及装置 |
Family Cites Families (15)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US4782284A (en) * | 1988-01-12 | 1988-11-01 | Bsr North America Ltd. | Frequency analyzer |
JP3147566B2 (ja) * | 1993-02-01 | 2001-03-19 | 日本精工株式会社 | 周波数スペクトル分析装置 |
WO2005101032A1 (en) * | 2004-04-18 | 2005-10-27 | Elspec Ltd. | Power quality monitoring |
DE102004058255A1 (de) * | 2004-12-03 | 2006-05-04 | Red-Ant Measurement Technologies And Services E.K. | Verfahren zur Bestimmung der Periodendauer eines Messsignals |
JP2006276006A (ja) * | 2005-03-01 | 2006-10-12 | Nagoya Institute Of Technology | 電力系統における高調波解析法 |
US7298129B2 (en) * | 2005-11-04 | 2007-11-20 | Tektronix, Inc. | Time arbitrary signal power statistics measurement device and method |
FR2913769B1 (fr) * | 2007-03-12 | 2009-06-05 | Snecma Sa | Procede de detection d'un endommagement d'un roulement de palier d'un moteur |
CN101136896B (zh) * | 2007-09-18 | 2011-06-29 | 东南大学 | 基于快速傅立叶变换的频域迭代均衡方法 |
CN101852826B (zh) * | 2009-03-30 | 2012-09-05 | 西门子公司 | 一种电力***的谐波分析方法及其装置 |
CN101655519B (zh) * | 2009-09-14 | 2012-06-06 | 国电南京自动化股份有限公司 | 数字化变电站测控装置交流采样数据处理方法 |
CN102752062B (zh) * | 2012-06-20 | 2015-04-29 | 京信通信技术(广州)有限公司 | 一种检测频率校正信道的方法及装置 |
CN103454495B (zh) * | 2013-09-13 | 2016-01-20 | 电子科技大学 | 自适应高精度快速频谱分析方法 |
CN105973999B (zh) * | 2016-04-28 | 2018-10-30 | 西安交通大学 | 基于增强相位瀑布图的转子裂纹微弱分数谐波特征识别方法 |
CN106154035A (zh) * | 2016-06-20 | 2016-11-23 | 哈尔滨工业大学 | 一种快速谐波及间谐波检测方法 |
CN107271774B (zh) * | 2017-07-10 | 2019-06-14 | 河南理工大学 | 一种基于频谱泄漏校正算法的apf谐波检测方法 |
-
2019
- 2019-09-19 CN CN201910888489.3A patent/CN110471018B/zh not_active Expired - Fee Related
Also Published As
Publication number | Publication date |
---|---|
CN110471018A (zh) | 2019-11-19 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN110471018B (zh) | 一种频谱校正方法 | |
CN104375111B (zh) | 对密集频谱进行快速高精度细化校正的方法 | |
Kang et al. | Phase difference correction method for phase and frequency in spectral analysis | |
Vold et al. | Theoretical foundations for high performance order tracking with the Vold-Kalman tracking filter | |
CN101113995A (zh) | 基于Nuttall窗双峰插值FFT的基波与谐波检测方法 | |
CN110837001B (zh) | 一种电力***中谐波和间谐波的分析方法与装置 | |
Poletti et al. | Comparison of methods for calculating the sound field due to a rotating monopole | |
Santamaria-Caballero et al. | Improved procedures for estimating amplitudes and phases of harmonics with application to vibration analysis | |
CN105938508B (zh) | 一种精确计算振动或压力脉动信号频率及幅值的方法 | |
CN107632961B (zh) | 基于全相位谱分析的多频内插迭代频率估计方法及估计器 | |
CN109764897B (zh) | 一种正余弦编码器高速信号采集及细分方法和*** | |
CN110598269A (zh) | 一种在低采样点时的离散频谱参数校正方法 | |
CN109374966A (zh) | 一种电网频率估计方法 | |
US6801873B1 (en) | Analysis of rotating machines | |
CN109541304B (zh) | 基于六项最小旁瓣窗插值的电网高次弱幅值谐波检测方法 | |
CN115265691B (zh) | 一种科氏流量计振动频率跟踪方法及*** | |
JP5632576B2 (ja) | マグニチュード値決定方法及び装置 | |
Nguyen et al. | A fast and accurate method for estimating power systems phasors using DFT with interpolation | |
US8032319B1 (en) | Methods for analyzing streaming composite waveforms | |
CN110017957B (zh) | 一种同步分析方法、装置及*** | |
CN113129912B (zh) | 一种单音信号的检测方法 | |
CN110285881A (zh) | 一种基于全相位滤波的密集谱频率估计法 | |
CN109960843B (zh) | 一种基于正交原理的多普勒频移数值仿真方法 | |
JP2000055949A (ja) | 周波数分析方法及び周波数分析装置 | |
CN114114231A (zh) | 一种基于fft-czt的雷达测距方法 |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
PB01 | Publication | ||
PB01 | Publication | ||
SE01 | Entry into force of request for substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
GR01 | Patent grant | ||
GR01 | Patent grant | ||
CF01 | Termination of patent right due to non-payment of annual fee |
Granted publication date: 20211224 |
|
CF01 | Termination of patent right due to non-payment of annual fee |