噪声谱估计范文(精选7篇)
噪声谱估计 第1篇
文中主要采用非参数法功率谱估计,对信号进行功率谱估计,采用周期图法和自相关函数法来对功率谱进行估计。
1 周期图法与自相关法功率谱估计介绍
功率谱估计的方法有周期图法、自相关函数估计法、自协方差法、最大似然法、最大熵法等,本文主要介绍了采用周期图法和自相关法来对信号进行功率谱估计。
1.1 周期图法功率谱估计的原理
周期图法是为了得到一个功率谱估计值,先要取得信号序列的离散傅里叶变换,然后取其幅频特性的平方再除以序列长度n。由于一般的序列x(n)离散傅里叶变换X(k)具有周期性,所以估计的功率谱也具有周期性,常称为周期图。周期图是信号功率谱的一个有偏估计值;而且,当信号序列的长度n增大到无穷大时,估值的方差反而不趋于零。因此,随着所取的信号序列长度的不同,所得到的周期图也相应的不同,这种现象称为随机起伏。由于随机起伏不稳定,所以使用周期图不能得到比较稳定的估计值。
1.2 自相关法功率谱估计的原理
自相关函数估计法是描述随机信号X(t)在任意两个不同时刻t1,t2的取值之间的相关程度;互相关函数给出了在一定的频域内两个时刻信号是否相关的一个判断指标,把两个测时间点之间信号的互谱与各自的自谱建立起了关系。它能用来表示输出信号有多大程度来自于输入信号,对于修正测量中接入噪声源而产生的误差起了很大的作用。维纳-辛钦定理:随机信号的相关函数与其功率谱是一傅里叶变换对,即相关函数的傅里叶变换为功率谱,而功率谱的逆傅里叶变换为相关函数。
由Wiener-Khintchine 定理得知,一个平稳随机过程的功率谱是自相关函数的傅里叶变换:
其中,rxx(m)为平稳随机过程的自相关函数。
r(l)=∫
可见,随机信号的功率谱与自相关函数互为傅里叶变换的关系,他们的关系是不同的,这两个函数分别从频率域和时间域来表达随机信号的基本特征。按上式计算功率谱估值,其运算量非常大,所以,将引入快速傅里叶变换算法(FFt),使运算更加简洁。 参数估计方法是有偏的谱估计方法,谱分辨率随数据序列长度的增加而提高。非参数估计方法大多是无偏的谱估计方法,从而得到更高、更好的谱分辨率。
2 Matlab程序仿真
2.1 Matlab程序如下
Matlab程序如下:
生成输入信号:
clear;
fs=1024;%设采样频率为1024
n=0:1/fs:1;
%产生含有噪声的序列 xn=cos(2*pi*40*n)+3*cos(2*pi*100*n)+randn(size(n));
N=length(n);
W=100*pi;%因方波频率F=1000HZ所以角频率W=100pi
X1n=square(W*n);%方波信号
X2n=randn(1,N);%白噪声信号
xn=X1n+X2n;
nfft=1024; [Pxx,f]=periodogram(xn,window,nfft,Fs); %直接法 plot(f,10*log10(Pxx));
subplot(3,1,1)
plot(n,xn);
xlabel()
ylabel(‘输入信号’)
%绘输入信号图
2.2 周期图法
fs=2000;
n=0:1/fs:1;
xn=cos(2*pi*40*n)+3*cos(2*pi*100*n)+randn(size(n)); nfft=1024;
N=length(n);
W=100*pi;
x1n=square(W*n);
x2n=randn(1,N);
xn=x1n+x2n;
subplot(3,1,1)
plot(n,xn);
plot_Pxx=10*log10(Pxx); plot_Pxx1=10*log10(Pxx1); plot_Pxx2=10*log10(Pxx2); figure(1) plot(f,plot_Pxx); pause; figure(2) plot(f,plot_Pxx1); pause; figure(3) plot(f,plot_Pxx2);
Nfft=1024;N=256;%傅里叶变换的采样点数1024
Pxx=abs(fft(xn,Nfft).^2)/N;
f=(0:length(Pxx)-1)*fs/length(Pxx);
subplot(3,1,2),
plot(f,10*log10(Pxx)),%转成DB单位
xlabel(频率/HZ),ylabel(功率谱/db),title(周期图法);
周期图法的功率谱密度如图1所示。
2.3 相关函数法
fs=1000;
n=0:1/fs:1;
xn=cos(2*pi*40*n)+3*cos(2*pi*100*n)+randn(size(n)); nfft=1024; window=boxcar(length(n));%矩形窗 noverlap=0; %数据无重叠 p=0.9; %置信概率
N=length(n);
W=100*pi;
x1n=square(W*n);
x2n=randn(1,N);
xn=x1n+x2n;
subplot(3,1,1)
plot(n,xn);%输入信号
plot_Pxx=10*log10(Pxx(index+1)); plot_Pxxc=10*log10(Pxxc(index+1)); figure(1) plot(k,plot_Pxx); pause; figure(2) plot(k,[plot_Pxx plot_Pxx-plot_Pxxc plot_Pxx+plot_Pxxc]);
m=-10:10
[r,lag]=xcorr(xn,100,biased)%求XN的自相关函数R,biased为有偏估计lag为R的序列号
subplot(3,1,2)
hndl=stem(m,r);%绘制离散图,分布点从-10—+10
set(hndl,Marker,.)
set(hndl,MarkerSize,2);
ylabel(自相关函数R(m))
%利用间接法计算功率谱
k=0:1000;%取100个点
w=(pi/50)*k;
M=k/50;
X=r*(exp(-j*pi/50).^(m*k));%对R求傅里叶变换
magX=abs(X);
subplot(3,1,3)
plot(M,10*log10(magX));
xlabel(功率谱的改进直接法估计)
自相关函数法的功率谱密度如图2所示。
实验结论:根据几个实验可知,首先要用Wavread函数读入wav文件中的数据,然后再对读入的数据进行处理和分析计算。读入后,共发现为20945个点。由于随机噪声信号是一个典型的非平稳信号,但在10ms~20ms内可以认为是短时平稳的。对数据进行分帧就显得尤为重要,用buffer函数实现。由于wav的采样频率为16000Hz,所以应有50*16000/100=8000个采样点。又根据傅里叶运算法则原理对每帧进行快速傅里叶变换基二fft算法,所以,选择每帧1024个采样点。又由于每帧之间具有关联,所以对50%的帧重叠以保持帧与帧之间的相关性。分帧处理后,对其中的一个帧进行自相关函数的运算,由于自相关函数与信号的功率谱是一对傅里叶变化对,再进行傅里叶变换,再转化成以噪声单位分贝(dB)为单位的数值,最后就得出了功率谱。
可能用到的matlab 指令:
Wavread()音频信号的读取
buffer( )用于信号分帧
fft() 快速傅里叶变换
xcorr() 求自相关函数
窗函数: hamming(),hanning(),blackman()。
原始语音信号如图3所示。
3 结束语
通过本次实验,掌握了如何用自相关法和周期图法对随机噪声信号进行功率谱的估计,以及所涉及到的知识和Matlab程序语言,在随机噪声信号处理中需要对信号进行估计,估计中还要进行分帧处理,在分帧时,由于随机噪声信号各帧之间的关联性,所以要保留一部分相同的数据,这样才能确保处理数据的完整性和可靠性,在对一帧数据进行处理时,如果采用自相关法,需要进行自相关函数求解,快速傅里叶变换,如果采用周期图法,需要进行非参数函数求解,并进行基二fft算法。了解了自相关函数和噪声信号的功率谱是一对傅里叶变换对的结论,对信号的功率谱进行了很好的估计,更深入理解了其原理及其估计方法。
摘要:在通讯与电子信息工程行业及领域中,大部分问题的解决需要进行估计一个随机信号在频率域上的功率谱分布,诸如此类的问题有很多,比如:设计滤波器消除噪声信号,振动随机信号的回波抵消,随机信号的特征抽取与表示等等。功率谱估计的分类:一般分为两大类,一类是参数法功率谱估计,一类是非参数法功率谱估计。参数法功率谱估计通常对数据进行一种建模,比如把数据建模成滑动平均模型(Moving Average),或者自回归(Autoregressive)模型,而非参数法功率谱估计。除了要求信号满足广义平稳之外,不需要其它的统计假设。与非参数法相比较,参数法的优点是在一个给定的数据集合上能够有较少的误差、偏差与方差。
关键词:周期图法,自相关法,Matlab
参考文献
[1]胡广书.数字信号处理——理论、算法与实现[M].北京:清华大学出版社,1997.
[2]S.M.凯依.现代谱估计原理与应用[M].黄建国,等译.北京:科学出版社,1994.
[3]傅广操,樊明捷.Matlab在现代功率谱估计中的应用[J].电脑学习,2003(6):54.
[4]段沛沛,罗丰.一种有效的短序列功率谱估计算法及其应用[J].雷达科学与技术,2004(6):25.
[5]Marsden J E,West M.Discrete mechanics and variational integrator[Z].Acta Numerica 10,2001:357-514.
噪声谱估计 第2篇
有源干扰是雷达电子对抗的主要手段,有源干扰主要包括压制性干扰和欺骗性干扰两种类型。在压制性干扰中,噪声调频信号具有宽的干扰带宽和较大的噪声功率,是目前雷达电子对抗中最常用的干扰形式。为了能让雷达等设备避开或工作在干扰较小的频段上,对噪声干扰信号的实时检测非常重要。随着数字信号处理技术的发展和数字处理芯片速度的提高,采用谱估计的方法来实时检测噪声信号已成为可能。本文针对噪声调频信号,分析了其信号频谱特点,研究通过估计其功率谱的方法来测量噪声调频信号的带宽和载波频率,并进行了仿真验证。
1 噪声调频信号的功率谱分析
由于噪声调频干扰信号具有高度随机性,因此无法像确定性信号那样用数学表达式精确地表述它,而只能用统计平均量来表征它。其中,自相关函数最能完整地表征它的统计平均值。而一个随机信号的功率谱密度正是自相关函数的傅里叶变换。因此可以用功率谱密度(Power spectral density,PSD)来表征它的统计平均谱特性。在实际工程环境中,常见的噪声调频干扰信号形式包括基带噪声调频、基带噪声叠加正弦波调制、基带噪声叠加锯齿波调制等几种形式。
1.1 基带噪声调频干扰信号[1,2]
基带噪声调频干扰信号是瞬时频率随基带调制噪声n(t)的幅度变化而变化的信号,其时域表达式为
式中,u0为噪声调频信号的幅度;kfm为调频指数;ω0为噪声调频信号的载波频率;ξ(t)为零均值、广义平稳随机过程,且具有如下功率带限功率谱:
其相关函数为:
其中:
式中ΔΨn=2πΔfn为调制噪声的谱宽,mfe=Kfmσ/Δfn=fde/Δfn为有效调频指数,fde为有效调频带宽。可求得噪声调频信号功率谱的表示式为:
上式中的积分只有当mfe>>1和mfe<<1时才能近似求解。当mfe>>1时,可以近似得到噪声调频信号的功率谱:
噪声调频干扰信号的频带宽Δfj,可按半功率电平对应的两个频率之间的差值来定,则带宽:
从上式可看出,噪声调频信号的干扰带宽与调制噪声带宽Δfn无关,而决定于调制噪声的功率σ2和调频指数kfm。通过控制fde即可控制噪声调频干扰信号的有效带宽。
1.2 基带噪声加正弦波调制的噪声调频干扰信号
设单音正弦信号表示为:
则基带噪声加上单音正弦波产生的噪声调频信号可表示为:
其傅里叶变换可以近似表示为:
从上式可看出,噪声加上单音正弦波产生的噪声调频信号的功率谱大致相当于单音正弦波调频的频谱与噪声调频信号的功率谱的卷积。
1.3 基带噪声加锯齿波调制的噪声调频干扰信号
锯齿波信号可以表示为:
则基带噪声加上锯齿波产生的噪声调频信号可表示为:
其傅里叶变换可以近似表示为:
从上式可看出,噪声加上锯齿波产生的噪声调频信号的功率谱大致相当于锯齿波调频(即线性调频信号)的频谱与噪声调频信号的功率谱的卷积。
2 噪声调频信号的谱估计方法
在检测噪声调频干扰信号中,主要测量其载波频率f0和带宽Δfj,从而获得干扰信号的频谱信息,进行有效的频谱规避或反干扰。从噪声调频信号的功率谱公式中可看出当f=f0时,Gj(f)取最大值,因此通过估计噪声调频信号频谱峰值来测量其载波频率f0,而其带宽可以用测量其频谱半功率电平对应的频率来计算。
功率谱估计[3,4]是信号分析与处理的一种非常重要的方法,它被广泛应用于雷达、声纳、生物医学、地震、机电等各种信号的处理场所。经典谱估计是基于传统傅里叶变换的一类方法,虽然存在固有缺陷,比如存在泄漏误差和混迭误差,分辨率低。但由于可以采用FFT算法,软硬件实现方便,而且经过改进的周期图谱估计可以减少谱估计的方差。
其中Bartlett法是把长度为N的数据xN(n)分成L段,分别求每一段的功率谱,然后加以平均。Welch法是对Bartlett法的改进,一是在对xN(n)分段时,可允许每一段数据有部分的交叠,二是每一段的数据窗口可以不是矩形窗口,这样可以改善由于矩形窗边瓣较大所产生的谱失真。
通过MATLAB软件的simulink来对噪声调频信号的检测进行建模和仿真。仿真模型如图1所示。系统采样率为200MHz,噪声调频信号产生器采用带宽为1MHz的窄带高斯白噪声信号,调频载波为50MHz,产生的噪声调频信号带宽为20MHz。有效调频指数为
由噪声调频信号的功率谱公式得到的近似频谱图如图2所示。
3 仿真和分析
在仿真中,对噪声调频信号源输出的信号采用不同参数下的Bartlett谱估计和Welch谱估计来测量其载波频率和带宽,并计算仿真测量的STD标准方差和RMS均方误差(用理论值来代替真实值)。
(1)采用Bartlett谱估计,FFT点数取1024点,取100次FFT求平均。采用Welch谱估计,FFT点数取1024点,512点交叠,窗函数采用Hann窗,取200次FFT求平均。结果如表1所示。
(单位:MHz)
(2)分别采用Bartlett谱估计和Welch谱估计,取不同的平均次数。结果如表2和表3所示。
(3)采用Welch谱估计,取200次FFT求平均,采用不同的窗函数。结果如表4所示。
(单位:MHz)
从仿真结果分析,可以得出如下结论:
(1)在对噪声调频信号的谱估计中,采用Bartlett谱估计和Welch谱估计相差不大。
(2)不同的窗函数下,采用Welch谱估计略有差异。
(3)不同的平均次数对测量结果影响很大,平均次数越多,测量的方差越小。说明频谱的积累时间越长,频谱的稳定性越好,对提高测量精度越有利。
因此,在实际应用时,可以根据噪声调频信号长度、测量精度、处理时间等要求,选择合适的谱估计方法和参数设置进行测量系统的设计。
4 结束语
本文分析了电子对抗压制式干扰中常见的几种噪声调频干扰信号的频谱特征讨论了检测噪声调频干扰信号的功率谱特征参数的谱估计方法。通过MAT-LAB/Simulink建模和仿真了不同谱估计方法的性能。在实际应用中,谱估计方法的精度和信号环境、噪声调频信号的时频域特征都有很大关系,需要根据实际情况进一步选择合适的算法和调整参数设置
摘要:对干扰信号的频谱特征进行精确测量和分析是实施有效的电子反干扰(ECCM)的重要步骤。通过对噪声调频干扰信号的频谱特征分析,提出采用合适的功率谱估计方法来检测噪声调频信号的主要频谱特征参数。利用MATLAB/Simulink仿真软件分别建立了噪声调频信号源和谱估计方法的系统模型,研究了系统采用不同的谱估计方法和参数设置的情况下所能够达到的检测精度,并根据仿真结果分析比较了不同谱估计方法的性能。
关键词:电子反干扰,噪声调频信号,信号谱估计
参考文献
[1]肖天南.噪声干扰机中噪声调相调频后干扰信号带宽的计算[J].舰船电子对抗,1998(2):17-20.
[2]夏红娟,陈潜.噪声调频干扰信号仿真及应用[J].上海铁道大学学报,2000(6):22-28.
[3]姚武川,姚天任.经典谱估计方法的MATLAB分析[J].华中理工大学学报,2000(4):45-47.
视频噪声方差估计算法研究 第3篇
关键词:视频噪声方差估计算法,视频去噪,图像处理
0引言
在视频监控系统中,固定安装的摄像机始终保持着场景不变,而带云台和镜头控制的摄像机部分时间会发生场景改变。因此,在监控视频图像序列中通常有大范围的静止背景区域,而运动前景区域的帧差图像不仅包含噪声,还包含视频真实信号。传统的时间差分视频噪声估计算法,往往容易把视频图像的实际信号分量也当成噪声,使得噪声估计结果偏大,尤其是在处理运动前景变化较大情况时,会令图像去噪后 视频序列 失掉一些 非常重要 的信息。
针对上述不足,本文提出了一种改进的噪声方差估计算法。该方法首先用加性高斯白噪声模型来模拟低光照下的摄像机噪声,然后对帧差图像进行高通滤波去除视频真实信号的低频分量,即高通滤波器对帧差图像的像素值依次与相邻像素作横向差分和竖向差分,最后再自适应地选取图像中大部分像素的差分值用来估计噪声,从而使得噪声方差估计更加准确和稳定,达到更好的去噪效果。
1改进的视频噪声方差估计算法
1.1计算帧差图像
假设观察到的视 频图像为Y(i,j,n)= X(i,j,n)+N(i,j,n),其中X(i,j,n)为没有噪声的视频信号,N(i,j,n)是加性高斯白噪声,i,j是像素的横坐标和纵坐标,n是时间坐标。
如果前后两帧完全保持静止,则有X(i,j,n)= X(i,j,n-1),此时的帧差图像D(i,j,n)可表示为:
假定原始图像的噪声方差为σ2N,由概率统计知识,可得出差分图像的噪声方差为:
因此,如果前后两帧图像存在静止块,则可以用静止块的差分图像来计算原始图像的噪声方差,即由式(2)可获得噪声方差的估计值为:
1.2图像分块及场景改变检测
对帧差图像分块[1],设块尺寸为W×W (本文采用4×4块),且这些块不可重叠,每个分块Bkl可表示为:
其中,Bkl表示左下顶点为 (k×W ,l×W),右上顶点为 (k×(W +1)-1,l×(W +1)-1)的矩形窗内像素帧差值的集合。
设分块Bkl的均值为MBkl,则
在静止背景处,帧差图像满足式(1),则MBkl均值如下式:
MBkl 的方差如下式:
因此,在静止背景处,MBkl服从均值为0,方差为2σ2N/W2的高斯分布。
根据高斯分布的性质,在区间 (u-σ,u+σ)的概率是0.6826,则在静止背景区域的块满足如下条件:
本文对多个场景静止序列的实验测试,块均值BMkl满足如下概率:
对多个场景改变序列的实验测试,块均值BMkl满足如下概率:
所以,本文判定块均值MBkl满足式(9)的图像为场景静止。
对于场景不变的图像帧,本文用差分值去计算噪声方差的估计值;对于场景变化的图像帧,本文以参考噪声方差作为当前帧的噪声方差估计值。其中,第一帧图像的参考噪声方差使用Lee和Hoppel[2]提出的基于分块的噪声估计方法计算得到(该方法计算的复杂度是各类噪声方差估计方法中最小的),其它帧的参考噪声方差则采用以前帧的噪声方差估计值。
1.3计算差分值
对帧差图像D(i,j,n)作横向差分和竖向差分,用于除去帧差图像中的原始视频信号的低频分量[3],即利用了视频图像的高频分量估计噪声,如下所示:
其中,D2(i,j,n)表示帧差图像的横向差分值,D3(i,j,n)表示D2(i,j,n)的竖向差分值。
根据高斯加性白噪声独立的性质,可知
则可以根据式(14)获得噪声方差的估计值,即
1.4自适应选取差分值并估计噪声方差
结合式(15),考虑到还有一些其它方向的纹理使得差分后的值包含了原始视频信号分量,从而影响了噪声方差估计的准确性。因此,自适应[6]选取差分值估计噪声方差的计算式如下:
其中,Avail(i,j,n)表示像素点的差分值D3(i,j,n)是否被选用。M和N分别表示图像的长和宽。
若帧差图像的块均值MBkl满足式(16),则使用当前块全部像素的差分值用于估计,即Avail(i,j,n)=1;否则不使用当前块全部像素的 差分值用 于估计,即Avail(i,j,n)=0。
2结语
基于DCT变换的噪声估计算法 第4篇
噪声无处不在, 在视频信号中同样存在噪声。很多视频图像处理都是在没有考虑噪声情况下得出的, 如果存在噪声, 必将影响这些算法的有效性。由于噪声对视频图像处理算法的影响如此大, 因此在开发图像视频处理算法时必须考虑噪声。由此必须有一种可靠的方法能准确地估计噪声, 以使很多图像视频处理算法 (比如边缘检测、图像分割和滤波等) 在噪声存在的情况下依旧保持非常良好的性能。
视频信号中的噪声来源非常广泛, 在视频获取、视频记录、视频传输以及视频处理等这些过程中均可产生噪声。在视频获取过程中, 摄像机在获取信号时会产生噪生, 胶片颗粒的存在也可产生噪声, 尤其在不良的光线条件下, 产生的噪声会更大。电荷藕荷器件摄像机由于高度灵敏性和较长的暴光时间也可能产生噪声。而且, 视频信号在模拟信道传输过程中也会产生噪声, 比如卫星信道和地面广播。此外, 视频记录设备也会产生噪声。
正是由于视频中噪声普遍存在, 噪声方差估计便成为了视频处理系统中的一个必不可少的任务。由此, 准确的噪声估计显得尤其重要。为了进行准确的噪声估计, 必须先对噪声进行建模。视频图像的噪声信号可建模成加性噪声或乘性噪声, 也可以被建模成独立于原图像或不独立于原图像信号。图像噪声也可分为白色噪声和有色噪声。在最普遍的情况下, 视频图像中的噪声信号通常被视为独立的、均匀分布的、平稳的零均值的白噪声。
1 相关的工作
大多数噪声方差估计方法的主要问题是, 在高噪声, 高质量的图像中, 噪声估计精度会大大下降, 尤其是在图像含精细结构或纹理的情况下, 更是如此。因此, 这些噪声估计算法不稳定, 即对一系列视频序列图片进行噪声估计时, 估计出的噪声方差会有较大的波动。对于许多视频处理系统来讲, 噪声估计的稳定性相当重要, 因此有必要开发一种噪声估计算法, 使得它在各种图像中依旧能给出稳定的噪声估计值, 也就是说其估计误差的方差会比较小, 同时计算量不应太大。
噪声估计有两种方法:基于平滑滤波处理和基于块处理。
在基于平滑滤波处理的方法中, 目标图像先经平滑滤波处理, 比如可采用平均值滤波, 于是原来的有噪声的图像与经过滤波后的图像的差值并可认为是近似噪声。然后根据一定的准则由这些近似噪声估计出最后图像的噪声方差。在基于块处理的方法中, 首先计算出图像中各个块的方差, 然后根据一定的准则由这些块的方差估计出最后图像的噪声方差。基于块处理的噪声估计算法比基于平滑滤波处理的噪声估计方法要简单些, 而且也要快些, 但是基于块处理的噪声估计算法在估计噪声时强烈地依赖于输入图像和噪声大小。
一种基于平滑滤波处理的噪声估计算法如下:首先由原来的有噪声的图像与经过滤波后的图像求出近似噪声, 然后由那些强度梯度小于给定门限值的像素的近似噪声估计出最后的噪声方差。这种方法在图像有精细纹理的情况下会失效, 而且往往会夸大地估计了噪声方差。一种基于块处理的噪声估计算法如下:首先计算出图像中各个块的方差, 然后把方差最小的那些块的方差取平均就得到了噪声方差的估计。这种方法存在一些问题, 最小方差的块并不意味着最平滑的块, 因此可能存在误估计的问题。其他的一些不同基于块的噪声估计算法趋向于在质量较好的图片中高估了噪声方差, 而在高噪声的图片中低估了噪声方差, 或者即使效果较好, 但计算量会太大。
提出一种新的基于块处理的噪声估计算法, 这种算法不仅在高质量图像和大噪声情况下给出可靠的噪声估计, 而且也大大减小了计算复杂度。这种方法考虑了图像的频谱特性, 而且并不是采用块的方差作为块均匀度的度量。
2 提出的算法
本文提出的算法是从一组最均匀的块中估计出全局图像噪声方差σ
本文提出的算法由两步构成:①亮度均匀块的检测, ②自适应平均这些块的方差。
2.1 亮度均匀块的检测
本文提出的算法使用一种低复杂度的均匀度度量εBk1来判定一个块是否拥有均匀亮度, 当然这里的均匀度量是按块度量的。块的定义由式 (1) 给出。
假定块中的像素是独立和均匀分布的, 但非零均值的。这里I (i, j) 表示带有噪声的信号, WK1是W×W大小的矩形块, 其左上角的像素位于 (k, 1) 。这些块中均匀样本的方差σ
文献[4]中给出了各种结构检测器的比较, 结果发现精确的结构与边缘检测需要很大的计算开销, 而计算量小些的检测器要么需要手动调整, 要么用于特定模型, 要么就是不够精确。本文提出的算法先将整幅图像分割成固定的W×W大小的块{Bk1} (这些块之间可以相互重叠, 这依赖于像素跳转因子s) 。在每一个块中, 使用局部均匀分析器来计算其均匀度度量εBk1。
均匀分析器是基于块DCT系数中的AC系数的。对于每一个这样的块, 按公式 (3) 计算出其DCT系数:
对u, v不同时为0的那些Fk1 (u, v) 构成块{Bk1}的AC系数。块的DC系数反映的是块中样本的平均值, 块的AC系数反映了某一块是否亮度均匀。众所周知, DCT变换具有压缩能量特性, 能使DCT系数值集中于低频DCT系数上, 也就是说原信号经DCT变换后, 其高频DCT系数很小或接近0。此特性同样适用于二维DCT变换, 只是其系数从低频到高频应按“之”字形排列, 如图1所示。
无噪声均匀亮度块的DCT DC系数即为此块像素的亮度值, 而其余DCT AC系数全为0, 这就是说一个亮度块越均匀, 其AC系数会越小。因此就可以根据W×W亮度块的DCT AC系数作为判断一个块均匀度的度量, 可将所有DCT AC系数的绝对值的和作为均匀度的度量, 但实际上这没有必要。上文提到, DCT变换具有压缩能量的特性, 即DCT系数值通常集中于低频DCT系数上, 因此不需要将所有DCT AC系数的幅值的绝对值相加, 只需取其低频的L-1个DCT AC系数即可, 这样做能减少计算量而极小地影响实际算法的效果。具体说就是将一个W×W亮度块的前L-1个 (按图1“之”字形排序) DCT AC系数的绝对值相加, 得到的和就作为均匀度度量, 即:
其中, Hk1 (n) 是代表块的第n个DCT系数 (按图1“之”字形排列) 。
2.2 自适应方差平均
①方差平均
为了估计全局图像噪声方差, 将m个最均匀块的样本方差取平均, 由此给出估计的噪声方差, σ
其中, V代表m个最小εBK1值的块的索引集。由于噪声是平稳的, m个最均匀块的方差的平均值就代表了整幅图片的噪声。为了达到更快的噪声估计, 并不是对所有的 (k, 1) 计算εBk1, 而是每隔s个像素 (同时水平和垂直两个方向) 去计算这些均匀度度量值。可以使用不同的s去评估算法的性能。选取s=5和W=5能够在效率 (即计算复杂度) 和有效性 (即估计精确度) 上达到很好的折中。值得注意的是, s的选取通常依赖于W的选取, 小的W通常意味着小的s。
②自适应平均
首先按如下方式计算出参考方差。通过上述均匀块检测公式已经计算出了所有块的均匀度度量, 接着找出最均匀的三个亮度块, 最后将这三个亮度块的方差取平均就得出了参考方差。由于此参考方差由最均匀的三个亮度的方差取平均而得, 因此它们最能代表图像的噪声方差。其他更好的高阶滤波方法可以用于得到参考方差。除了平均方法, 可以取中值。仿真结果表明, 取中值的效果更好。因为在某些情况下, 前三个最均匀的亮度块的方差值之间相差较大, 中值滤波能给出参考方差的良好估计值。
公式 (5) 平均过程中的块可能会表现出不同的均匀性, 这将引起对噪声方差的过高估计。为了避免对噪声的过高估计, 仅仅将那些均匀度相似的块纳入平均计算过程。如果一个块的方差σ
其中σ
3 实验与分析
为了评估噪声估计算法的性能, 应先计算出估计误差Ek=|σ
其中Ek (i) 为第i幅噪声方差为σ
本实验用一系列标准灰度图像对本算法进行测试, 图像灰度值在[0, 255]的范围内, 在人为加噪声前, 先将灰度值进行幅度变换, 将其取值范围变换成在[16, 235]范围内。先人为加入某一幅度的白噪声, 再用本文提出的算法算出个每幅测试图片的估计噪声方差, 然后再由式 (7) 求出μEk和σEk。实验结果如表1所示, 其中演示了各种噪声幅度的实验结果。经过实验和分析, 选取s=5和W=5能够在效率 (即计算复杂度) 和有效性 (即估计精确度) 上达到很好的折中。
从表1可以看出, 尽管噪声越高, 其估计准确性会变差, 但是总体上讲, 本文提出的噪声估计算法表现出很好的稳定性。此外, 经过实验表明, 其在不同质量的图像下, 也表现出良好的一致性。
4 结束语
本文提出了一种基于DCT变换的自适应噪声估计算法。理论分析及实验结果表明本算法不仅在低噪声的图片中表现出良好的性能, 而且在高噪声的图片中依旧有效。而且, 还能适应不同质量的图像。此外, 本文提出的算法复杂度也不大, 能适用于各种实时视频图像处理系统。
参考文献
[1] Canny J. A computational approach to edge detection[J]. IEEE Trans. Pattern Anal. Machine Intell., Nov. 1986, 9 (6) :679-698.
[2] Rosin P L. Thresholding for change detection[C]// Proc. IEEE Int. Conf. Computer Vision, Bombay, India, Jan. 1998:274-279.
[3]Amer A, Schr H.der, Anewvideo noise reduction algorithmusingspatial sub-bands[C]//Proc.IEEE Int.Conf.Electron., Cir-cuits, and Syst., Rodos, Greece:Oct.1996.1:45-48.
[4] Lyvers E P. , Mitchell O R. Precision edge contrast and orientation estimation[J], IEEE Trans. Pattern Anal. Machine Intell., Nov. 1988, 10 (6) :927-937.
[5]Rank K, Lendl M, Unbehauen R.Estimation of image noise vari-ance[J].Proc.IEE Vis.Image Signal Process., Apr.1999, 146 (2) :80-84.
[6]Aishy Amer.Eric Dubois Fast and Reliable Structure-Oriented Vide-o Noise Estimation[J].IEEE Transactions on Circuits And SystemsFor Video Technology, Jan.2005, 15 (1) :113-114.
直扩/跳频信号功率谱密度估计 第5篇
1 直扩/跳频通信系统的组成
1.1 直扩/跳频信号的产生
直扩/跳频信号的发射系统如图1所示。信源产生的信号,经编码器后转化成压缩了的数字信号,再经过调制器的相应调制,输出的已调波信号载波频率达到射频通带的要求,通过直扩模块和跳频模块,对其进行了扩频和跳频处理,获得直扩/跳频信号,然后经过功率放大器后,至天线发射出去。
1.2 直扩/跳频信号的接收
直扩/跳频信号的接收系统如图2所示。在接收端,接收到的信号通过功率放大器处理后,送至混频器,再与频率合成器产生的载波信号混频,这些混合信号通过解直扩模块和解调器,可以消除组合波频率成分,恢复出发送的信号。在接收过程中,接收机中的频率合成器必须受同步指令的控制,这样可以有效地抑制干扰信号,不会对直扩/跳频系统产生干扰[2,3]。
2 直扩/跳频信号功率谱密度估计分析
直扩/跳频信号是随机信号,因此无法像确定信号那样用数学表达式来精确地描述它,而只能用它的各种统计平均量来表征。自相关函数最能完整地表征它的特定统计平均量值,而一个随机信号的功率谱密度正是自相关函数的傅里叶变换,可以用功率谱密度来表征它的统计平均谱特性。因此,在统计意义下描述一个随机信号,就需要估计它的功率谱密度[4,5]。
2.1 功率谱密度的估计
对于功率信号,其功率谱密度可定义如下:把f(t)在间隔|t|>T 2以外的部分截去,得到截短函数:
只要T为有限值,则sT(t)的能量ET也是有限值。设ST(ω)为sT(t)的频谱函数,这样sT(t)的能量ET为:
因此平均功率P为:
当T增加时,sT(t)的能量和|ST(ω)|2也增加。当T趋于无穷时,|ST(ω)|2T的极限可能存在,令Tl i→m∞|ST(ω)|2T=Ps(ω),此极限为功率谱密度。
2.2 功率谱估计方法
功率谱估计分为经典谱估计和现代谱估计两大类方法。经典功率谱估计又分为直接法和间接法。直接法又称周期图法,它是直接由傅里叶变换得到的;间接法又称自相关法或BT法,它是通过自相关函数间接得到的。
用周期图法仿真时,直扩/跳频信号的参数设置为:频率集F=[0.5,1.0,1.5,2.0,2.5,3.0,3.5,4.0]MHz,跳频间隔Δf=0.5MHz,扩频因子a=12,采样频率fs=90 MHz。在信噪比SNR=10 dB时,取FFT的长度分别为128和1 024,运用周期图法对直扩/跳频信号进行仿真分析,如图3,图4所示。
由图3和图4可以看出,周期图法谱估计中,数据长度N越大,谱分辨率越大,但N太大会导致方差性能变差,谱线起伏加剧。谱分辨率和谱线起伏成为了一个不可调和的矛盾。这时就需要对直接法进行改进,可采用Welch法进行改进。
改进一,对数据进行分段处理,分段时允许每一段的数据有部分的交叠,这样段数越多估计结果的方差也就越小。但是,由于重叠的段会使各段之间具有统计相依性,反而会导致方差增大,所以在分段数目与重叠之间选择上存在一个折衷。
改进二,每一段的数据窗口可以用海宁窗或海明窗等窗函数代替矩形窗,这样可以改善由于矩形窗旁瓣较大所产生的谱失真[6]。表1给出了几种常用的窗函数的主要性能指标参数值(其中频率变化两倍的区间称为一个倍频程)。
主瓣带宽决定了被截短以后所得序列的频率分辨率,而旁瓣峰值有可能湮没信号频谱分量中较小的成分。选择窗函数时,希望频谱的主瓣尽量窄,旁瓣峰值尽量小且衰减尽量快,使频域的能量尽量集中在主瓣内,减少频谱“泄露”,并且希望选取其频谱恒为正的窗函数。比较表1中的5种窗函数可以看到,矩形窗具有最窄的主瓣,但也有最大的旁瓣峰值和最慢的衰减速度。海宁窗的主瓣较宽,同时有较小的旁瓣和较大的衰减速度,是较为理想的窗函数,因此本文选取海宁窗进行仿真。图5和图6用韦尔奇法对比矩形窗和海宁窗的效果。
从仿真图可以看出,海宁窗对减少“旁瓣效应”,即功率谱泄露,能起到一定的效果,也可以使峰值的宽带增大。多次实践表明,取合适的海宁窗和一半段长度的重叠率,可以最有效地降低估计的偏差。
2.3 功率动态不平衡度测试
直扩/跳频信号带内载波功率动态不平衡度是指跳频过程中发射机带内频点功率变化最大值占平均功率的百分比。仿真时设定每个跳频点的信号幅度相等,因此实验所得的信号功率谱幅度理论上应该是一致的,但是实际情况并非如此。特别说明一下,由于Welch法的功率谱要除以fs,因此所得结果为负值。
仿真条件直扩/跳频信号频率集F=[0.5,1.0,1.5,2.0,2.5,3.0,3.5,4.0]MHz,跳频间隔Δf=0.5 MHz,扩频因子a=12,采样频率fs=90 MHz。在信噪比SNR=10 dB的情况下,根据以上两种方法测试直扩/跳频信号带内载波功率动态不平衡度的结果见表2。
由表2可以看出,无论是周期图法还是Welch法,功率动态不平衡度测试结果均随着窗长的变长而增大,表明带内各个跳频点的功率谱方差变大即一致性变差。在相同的仿真条件下,Welch法比周期图法测试得到的动态功率不平衡度要小,即各个跳频点的功率谱更加接近一致,表明Welch法的测试功率谱的性能要明显优于周期图法。综合以上分析,Welch法能够获得效果良好的跳频信号功率谱,并能够较为准确地测试出跳频信号带内载波功率动态不平衡度。
3 结语
本文主要对直扩/跳频信号的功率谱做了分析,指出了周期图法的不足,再从改进周期图法的思路出发,结果表明韦尔奇法能够获得性能良好的直扩/跳频信号功率谱,并能够较为准确地测试出直扩/跳频信号带内载波功率动态不平衡度。
参考文献
[1]朱攀.直接序列/跳频混合扩频信号检测与参数估计方法研究[D].成都:电子科技大学,2007.
[2]葛哲学,陈仲生.Matlab时频分析技术及其应用[M].北京:人民邮电出版社,2006.
[3]陆根源,陈孝帧.信号检测与参数估计[M].北京:科学出版社,2004.
[4]张明照,何玉彬,刘政波.应用Matlab实现信号分析和处理[M].北京:科学出版社,2006.
[5]赵宏伟.基于时频分布的跳频信号参数检测方法研究[D].西安:西北工业大学,2006.
[6]张峰,石现峰,张学智.Welch功率谱估计算法仿真及分析[J].西安工业大学学报,2009,29(4):353-356.
[7]于训峰,马大玮,魏琳.改进周期图法功率谱估计中窗函数仿真分析[J].计算机仿真,2008,25(3):111-114.
[8]王晓峰,王炳和.周期图及其改进方法中谱分析率的Matlab分析[J].武警工程学院学报,2003(6):75-77.
[9]姚武川,姚天任.经典谱估计方法的Matlab分析[J].华中理工大学学报,2000,28(4):45-47.
一种快速且准确的噪声估计算法 第6篇
I(i,j)=S(i,j)+N(i,j) (1)
其中,S(i,j)为原始图像,I(i,j)为有噪声的图像,N(i,j)为与图像信号相互独立的噪声。
现有的噪声估计算法可分为两大类:帧内方法和帧间方法。帧内方法又可以分为基于块、基于滤波和基于小波的方法[1,2,3]。本文提出的算法只考虑零均值加性高斯白噪声,它是基于块的帧内方法,但同时运用帧间信息来避免排序,进而加快噪声估计的速度,有利于快速实现。
去噪算法普遍要求噪声估计能够做到足够准确[4,5],本文将在测试噪声估计的准确性对去噪算法影响的基础上,提出噪声估计的快速实现算法。在文献[5]中,Kostadin Dabov等人提出了非常出色的BM3D去噪算法,本文以此算法来作为研究对象。
1 噪声估计对去噪算法的影响
假定估计的噪声分别为sigma=15、20、25,对去噪前后的结果进行PSNR分析,表1为分析结果。图1(a)为加了sigma=20高斯白噪声的图像的局部,图1(b)、1(c)、1(d)分别为三种不同噪声估计值对应的BM3D去噪结果。
从表1中可以看出,当估计噪声与真实噪声存在差别时,BM3D的去噪结果会受到影响。当噪声估计值为15时,可以明显看出图像中还存在噪声。当噪声估计值为25时,虽然在图像中看不到明显的噪声,但会有过度平滑的现象出现,图像的细节被掩盖了,所以,准确的噪声估计是良好的去噪算法的保证。接下来,本文介绍噪声估计算法。
2 基于块的噪声估计算法
文献[1]介绍了一种基于块的噪声估计算法。假设噪声为加性零均值高斯白噪声,σundefined是图像的噪声的方差。为了估计σundefined,将图像分割为尺寸相同的块。参考Liwei Guo等人的经验[1]以及本文的实验测试结果,对采用CIF或更低分辨率格式的视频此算法采用8×8尺寸的块,而大一些的视频采用16×16尺寸的块。假设在当前帧共有TBN个块,选择其中方差值最小的3%[1]的块,然后计算它们方差的平均值undefined。再按照如下公式估计当前帧噪声的方差σundefined:
undefined
其中,α是校正因子,通过实验测得的[1]:
undefined
在选择方差值最小的3%的块时难免要用到排序方法,而目前常用的排序方法,其复杂度是o(n2)或者o(nlon2n),无法实现快速的噪声估计。无论是硬件还是软件处理,大部分处理的视频数据是从左往右,从上往下输入的,基于这个特点本文提出一种改进的快速且准确的噪声估计算法。
3 改进的快速估计算法
基于文献[1]中提出的噪声估计算法,本文提出了一种改进算法。对于每帧视频输入,都需要按照以下步骤计算相应的参数。
第一步,假设一帧图像的尺寸为M×N,将当前帧分割成若干个N1×N1大小的像素块,当图像大于CIF时,取N1=16,当图像小于或等于CIF时,取N1=8。
第二步,采用N1×N1尺寸的滑动窗口在当前帧的第一行像素块从左往右滑动,计算每个像素块的方差值,然后计算出此行像素块方差的均值undefined,其中每行的块数由N除以N1后向上取整获得。
第三步,参照步骤第二步,求出其余各行像素块的方差均值undefined,其中Y是每列的块数,由M除以N1后向上取整获得。
第四步,当前帧像素块方差的均值为:
undefined
改进后的算法流程如下:
①首先要对视频序列的第一帧进行中值滤波,然后按照第一步至第四步步骤计算滤波后视频帧像素块的方差均值undefined。
②下一帧数据不进行中值滤波,求出像素块方差的均值undefined,利用①求得的undefined计算出调整系数α:
undefined
其中,50是实验得的数据,然后根据α以及undefined来计算初始参考阈值Vth:
undefined
利用该阈值,在下一帧中选择出的块的数目大约是一帧图像中总像素块数的3%。
③将下一帧作为当前帧,分割成若干个N1×N1大小的像素块,计算每个像素块的方差值,记录小于初始参考阈值Vth像素块的方差值,并同时记录方差小于Vth的像素块的总个数Cv,求出所有记录下来的像素块的方差均值undefined,再根据(2)式求出当前帧噪声的方差σundefined,然后对阈值Vth进行一次调整:
当undefined时,细调阈值Vth,取步长为Lstep=ceil(TBN×0.005),中心为L=ceil(TBN×0.03);
当undefined时,粗调阈值Vth,取步长为Lstep=ceil(TBN×0.01),中心为L=ceil(TBN×0.01)。
当Cv在(L-(i+1)Lstep,L-i*Lstep)范围内时,取参考阈值为Vth+Vstep;当Cv在(L+i*Lstep,L+(i+1)*Lstep)范围内时,取参考阈值为Vth-Vstep。其中,TBN是当前帧N1×N1像素块的数目,Vstep将在后面的参数计算中给出。
重复执行步骤③,经过若干帧的调整后,按不断调整的阈值Vth取出的像素块的数目占总像素块的数目的百分比会收敛于3%,此时求出的σundefined是稳定的噪声估计值。如果undefined大于15%,则回到②重新设定初始参考阈值。
4 参数计算
对于加性零均值高斯白噪声,本文假设其满足X~N(0,σ),那么位于平坦区域的块的方差符合的分布为:
undefined
表2为按不同的百分比选择像素块时所对应的真实阈值,即按特定百分比选出的块中的最大的方差值。
对多幅图像进行处理后拟合得到Vstep的估计方法为:undefined。输入视频的噪声通常变化不大,每帧的Vth变化也不大,Cv一般满足undefined。Vth会随着每帧Cv的值做调整。
视频噪声突变时,Cv的值可能会超出undefined范围,这时对Vth的调整将难以把握,本算法选择:当undefined时,下一帧重新计算方差均值undefined,进而得到新的阈值Vth,当前帧沿用上一帧估计出的噪声方差值,在连续的视频中这种处理是可以接受的;当undefined时,取阈值Vth为Vth+10Vstep,然后重复步骤③。
此算法不需要太大的存储空间,最多需要L+ceil(TBN+0.22)个存储空间。运算过程中不需要处理已经被处理过的像素块,因此复杂度只有o(n)。噪声估计的准确程度很大程度上依赖于块的数目,本文采用总块数的3%[1],块数太多,容易照成估计值偏大;块数太少,估计值也会不准确。
5 实验结果和讨论
本文给出两个测试实验结果。对于尺寸为352×288的图像:TBN=1584, L=48(1%≤Cv/TBN≤5%), L=16(5%≤Cv/TBN≤15%), Lstep=8(1%≤Cv/TBN≤5%), Lstep=16(5%≤Cv/TBN≤15%)。对于尺寸为的图像:TBN=6336, L=192(1%≤Cv/TBN≤5%), L=64(5%≤Cv/TBN≤15%), Lstep=32(1%≤Cv/TBN≤5%), Lstep=64(5%≤Cv/TBN≤15%)。
为了模拟噪声的突变,试验中修改图像中添加的噪声(第8帧sigma由20变为15,第15帧sigma由15变为20)。由图2可知,该算法可以在4帧内达到平稳,在真实的视频监控系统中这种速度是完全可以接受的。
由图3可知,两种算法的准确性都很高,且趋于一致。但由表3可知,本文提出的算法要比原算法快,帧尺寸较小时涉及到的排序次数也较少,效果不明显,但是当尺寸较 大时,速度的增加还是很明显的,由排序法的o(n2)或者o(nlog2)降为o(n)。
6 结束语
提出了一种基于块的快速且可靠的噪声估计算法。此算法根据视频帧的信息自适应地调整部分参数,能够提供可靠的估计结果,同时巧妙地利用阈值比较法避免了排序过程,显著降低有排序过程的噪声估计算法的复杂度,从而实现噪声的快速可靠估计。
参考文献
[1] Olsen S I. Estimation of Noise in images: An Evalution[J]. Graphi-cal models and image processing, 1993, 55(4):319-323.
[2] Dong-Hyuk Shin, Rae-Hong Park, Seungjoon Yang, et al. Block-Based Noise Estimation Using Adaptive Gaussian Filtering[J]. IEEE Trans. Consumer Electronics., 2005, 51(1):218-226.
[3] Donoho D L, Johnstone I M. Ideal Spatial Adaptation by Wavelet Shrinkage[J]. Biometrika, 1994, 81: 425-455.
[4]Liwei Guo,Oscar C Au,Mengyao Ma,et al.Fast Multi-HypothesisMotion Compensated Filter for Video Denoising[J].Journal of Sig-nal Processing Systems,2010,60(3):273-290.
基于DSP的AR谱估计的实现 第7篇
信号处理器(DSP)是一种可编程的高性能处理器,近年来发展很快,不仅适用于数字信号处理,且在图像处理、语音处理、通信等领域都得到了广泛应用[1]。DSP结构设计提供了低成本、低功耗和高性能的处理能力,以DSP作为核心是智能化制造技术中智能控制器的发展方向。AR谱分析是现代谱估计中的一种重要的方法,它相对于经典谱估计法有很多的优势:1)AR谱的分辨率很高,不用加窗函数,不存在数据截取的问题,理论上没有泄漏误差[2]。2)对于长时间序列可以取短序列进行计算,减少计算工作量,能够在很短的时间内完成功率谱估计、趋势预测。目前,AR谱已在机械故障诊断中得到了广泛应用,而DSP的出现使AR谱的计算更加方便快捷。
1 DPS芯片的基本结构和特点
TMS320F2812芯片是美国TI公司最新推出的C2000平台上的定点DSP芯片。F2812芯片具有低成本、低功耗和高性能的处理能力,特别适用于有大量数据处理的监控场合。其具有以下特点[3]:
1)采用高性能静态CMOS技术,I/O供电电压及Flash编程电压为3.3V,内核供电电压降为1.8V(135MHZ)或1.9V(150MHZ),减小了控制器的功耗。150MI/S(百万条指令/S)的执行速度使得指令周期减小到6.67ns,从而提高了控制器的实时控制能力。
2)片内高达128K 16位的Flash存储器,包括4个8K 16位和6个16K 16位的程序存储器、1k 16位的一次性可编程存储器、两个1K 16位的存储器M0和M1、两个4K 16位的存储器L0和L1、一个8K 16位的存储器H0。
3)完善的开发环境CCStudio v3.3,方便开发者进行迅速开发。开发工具包括ANSI C/C++编译器、汇编器、连接器,支持C24X/F240x指令,CCS IDE,DSP/BIOS,JTAG接口。
此外,还有看门狗定时器模块(WDT),控制器局域网络(CAN)2.0模块,串行通信接口(SCI)模块,16位串行外设(SPI)接口模块,40个可单独编程或复用的通用输入/输出引脚(GPIO),12位A/D转换模块等。
2 AR时序模型的建立、模型参数递推估计及定阶
2.1 AR模型的原理
时间序列分析方法作为一种现代数据处理方法,通过所建立的数学模型,不仅研究观测数据本身的统计特性,而且还研究产生此观测数据的系统的动态特性,从而将数据、模型和系统联系起来进行分析处理。正因为这样,对机械工程领域而言,它可以用于机械系统的系统辨识、系统分析、模态参数识别、故障监视与诊断、状态预测与控制等[4]。特别是由于时序方法具有一系列的较优越的方面,具有较强的诊断能力,因此已在机械故障诊断中受到越来越广泛的重视,而且取得了很好的效果。本文根据实验室齿轮箱所观测的数据建立AR模型,以此模型为基础采用burg算法来实现齿轮箱运行状态的谱估计。
对于实验室所采集的齿轮箱运行状态信号,经过滤波、平滑、去除趋势项等处理后,得到平稳、正态、零均值的时间序列{xn},其AR(p)模型的结构为:
式中:xn—时间序列{xn}在n时刻的元素
ak—自回归系数,k=1,2,…,p
en—预测误差
2.2 模型的定阶准则
AR模型的定阶是一个核心问题,它关系到所建立的模型的精确性、分辨率及运算时间,过高的阶次虽然会较大的提高谱的分辨率,但是它所需要的运算时间较长,过低的阶次虽然运算所需时间缩短,但是估计出的AR谱,其精确性和分辨率较差。
本文采用日本学者Akaike提出的最终预测误差准则(FPE)作为AR模型的定阶准则,其表达式为:
式中:N—数据的长度
P—模型阶次
ep—P阶时刻AR(P)模型的预测误差的均方值
当N固定,P增大时,预测误差的均方值ep减小,而增大,取FPE(P)最小值或变化趋势缓慢时,确定最佳合适阶次[5]。
2.3 AR模型的burg算法及定阶实现
1)根据采集的离散时间序列数据,初始化模型0阶时刻向前、向后误差及0阶预测误差:
2)求AR(P)模型P阶时刻反射系数:
3)根据L-D递推公式估计AR(P)模型的参数和预测误差的均方值:
4)递推更高一阶的向前、向后预测误差序列:
5)根据式(3)、(4)、(5)及(6),求出1阶至P阶各阶次的预测误差均方值e(p),带入式(2)求出各阶次{AIC(P)},选择最小{AIC(P)}以确定模型阶次。
6)最后根据所求的AR参数估计值ak(p)和e(p),直接计算burg功率谱的密度为:
3 实验对比
根据以上递推算法和定阶准则,在CCS3.3开发环境中编写AR建模的burg算法和AR谱C程序。用B&K公司Pulse Labshop分析系统对旋转机械振动分析及故障实验平台的齿轮箱进行实测,测得振动时域信号,如图1所示。然后利用B&K公司Pulse Labshop分析系统自带的FFT自功率谱计算分析功能得出所测信号的自功率谱,如图2所示。
在CCS3.3开发环境中利用程序读取Pulse Labshop分析系统所采集的1024个数据,建立AR模型。模型最高搜索阶数定为100,多次的计算分析结果表明,最高阶数为100是可行且具有较好的效果。得到被测信号的AR谱如图3所示,模型参数如表1所示。
由图3与图2的对比可以看出,AR谱中不仅包含了FFT谱中所含有的频率成份,而且它比FFT谱更加的平滑,没有毛刺,谱峰更加突出。从运算量大小来看,对所采集的1024个数据,只需建立63阶AR模型就能较好的估计其功率谱。
4 结束语
AR谱和FFT谱估计所得出的结果是一致的,由于AR谱不存在加窗函数的影响,所以它具有较好的分辨率,没有毛刺并且谱峰突出等优势。AR谱估计的计算量与p2成正比关系,使得他它相对于FFT谱估计的计算量大大的简化。同时,DSP芯片TMS320F2812的高性能的运算处理能力使AR谱估计的实现更加的快捷。这些都能为机械设备的实时在线谱估计提供十分有效的帮助。
参考文献
[1]杜艳生,谢克明,杨斌虎.基于DSP的数据采集及FFT实现[J].太原理工大学学报,2004(5):279-281.
[2]John G.Proakis.统计信号处理算法[M].北京:清华大学出版社,2006.
[3]刘和平,邓力,等.数字信号处理器原理、结构及应用基础—TMS320F28X[M].北京:机械工业出版社,2007.
[4]汪源源.现代信号处理理论和方法[M].上海:复旦大学出版社,2003.