[通信电子]频谱分析完整版.doc
《[通信电子]频谱分析完整版.doc》由会员分享,可在线阅读,更多相关《[通信电子]频谱分析完整版.doc(29页珍藏版)》请在三一办公上搜索。
1、Matlab 信号处理工具箱 帮助文档 谱估计专题 翻译:无名网友 & Lyra频谱分析Spectral estimation(谱估计)的目标是基于一个有限的数据集合描述一个信号的功率(在频率上的)分布。功率谱估计在很多场合下都是有用的,包括对宽带噪声湮没下的信号的检测。从数学上看,一个平稳随机过程的power spectrum(功率谱)和correlation sequence(相关序列)通过discrete-time Fourier transform(离散时间傅立叶变换)构成联系。从normalized frequency(归一化角频率)角度看,有下式注:,其中。其matlab近似为X=
2、fft(x,N)/sqrt(N),在下文中就是指matlab fft函数的计算结果了使用关系可以写成物理频率的函数,其中是采样频率相关序列可以从功率谱用IDFT变换求得:序列在整个Nyquist间隔上的平均功率可以表示为上式中的以及被定义为平稳随机信号的power spectral density (PSD)(功率谱密度)一个信号在频带上的平均功率可以通过对PSD在频带上积分求出从上式中可以看出是一个信号在一个无穷小频带上的功率浓度,这也是为什么它叫做功率谱密度。PSD的单位是功率(e.g 瓦特)每单位频率。在的情况下,这是瓦特/弧度/抽或只是瓦特/弧度。在的情况下单位是瓦特/赫兹。PSD对频
3、率的积分得到的单位是瓦特,正如平均功率所期望的那样。对实信号,PSD是关于直流信号对称的,所以的就足够完整的描述PSD了。然而要获得整个Nyquist间隔上的平均功率,有必要引入单边PSD的概念:信号在频带上的平均功率可以用单边PSD求出频谱估计方法Matlab 信号处理工具箱提供了三种方法Nonparametric methods(非参量类方法)PSD直接从信号本身估计出来。最简单的就是periodogram(周期图法),一种改进的周期图法是Welchs method。更现代的一种方法是multitaper method(多椎体法)。Parametric methods (参量类方法)这类方
4、法是假设信号是一个由白噪声驱动的线性系统的输出。这类方法的例子是Yule-Walker autoregressive (AR) method和Burg method。这些方法先估计假设的产生信号的线性系统的参数。这些方法想要对可用数据相对较少的情况产生优于传统非参数方法的结果。Subspace methods (子空间类)又称为high-resolution methods(高分辨率法)或者super-resolution methods(超分辨率方法)基于对自相关矩阵的特征分析或者特征值分解产生信号的频率分量。代表方法有multiple signal classification (MUSI
5、C) method或eigenvector (EV) method。这类方法对线谱(正弦信号的谱)最合适,对检测噪声下的正弦信号很有效,特别是低信噪比的情况。方法描述函数周期图PSD 估计spectrum.periodogram, periodogramWelch重叠,加窗的信号段的平均周期图spectrum.welch, pwelch, cpsd, tfestimate, mscohere多椎体多个正交窗(称为锥)的组合做谱估计spectrum.mtm, pmtmYule-Walker AR时间序列的估计的自相关函数计算自回归(AR)谱估计spectrum.yulear, pyulearBu
6、rg通过最小化线性预测误差计算自回归(AR)谱估计spectrum.burg, pburgCovariance(协方差)通过最小化前向预测误差做时间序列的自回归(AR)谱估计spectrum.cov, pcov修正协方差通过最小化前向及后向预测误差做时间序列的自回归(AR)谱估计spectrum.mcov, pmcovMUSIC多重信号分类spectrum.music, pmusic特征向量法虚谱估计spectrum.eigenvector, peigNonparametric Methods非参数法下面讨论periodogram, modified periodogram, Welch, 和
7、 multitaper法。同时也讨论CPSD函数,传输函数估计和相关函数。Periodogram周期图法一个估计功率谱的简单方法是直接求随机过程抽样的DFT,然后取结果的幅度的平方。这样的方法叫做周期图法。一个长L的信号的PSD的周期图估计是注:这里运用的是matlab里面的fft的定义不带归一化系数Matlab FFT函数未做归一化,所以要除以L其中实际对的计算可以只在有限的频率点上执行并且使用FFT。实践上大多数周期图法的应用都计算N点PSD估计,其中选择N是大于L的下一个2的幂次是明智的,要计算我们直接对补零到长度为N。假如LN,在计算前,我们必须绕回模N。作为一个例子,考虑下面1001
8、元素信号,它包含了2个正弦信号和噪声randn(state,0);fs = 1000; % Sampling frequencyt = (0:fs)/fs; % One second worth of samplesA = 1 2; % Sinusoid amplitudes (row vector)f = 150;140; % Sinusoid frequencies (column vector)xn = A*sin(2*pi*f*t) + 0.1*randn(size(t);注意:最后三行表明了一个方便的表示正弦之和的方法,它等价于:xn = sin(2*pi*150*t) + 2*si
9、n(2*pi*140*t) + 0.1*randn(size(t); 对这个PSD的周期图估计可以通过产生一个周期图对象(periodogram object)来计算Hs = spectrum.periodogram(Hamming);估计的图形可以用psd函数显示。psd(Hs,xn,Fs,fs,NFFT,1024,SpectrumType,twosided)平均功率通过用下述求和去近似积分 求得Pxx,F = psd(Hs,xn,fs,twosided);Pow = (fs/length(Pxx) * sum(Pxx)Pow = 2.5059你还可以用单边PSD去计算平均功率Pxxo,F
10、= psd(Hs,xn,fs,onesided);Pow = (fs/(2*length(Pxxo) * sum(Pxxo)Pow = 2.5011周期图性能下面从四个角度讨论周期图法估计的性能:泄漏,分辨率,偏差和方差。频谱泄漏考虑有限长信号,把它表示成无限长序列乘以一个有限长矩形窗的乘积的形式经常很有用:因为时域的乘积等效于频域的卷积,所以上式的傅立叶变换是前文中导出的表达式说明卷积对周期图有影响。正弦数据的卷积影响最容易理解。假设是M个复正弦的和其频谱是对一个有限长序列,就变成了所以在有限长信号的频谱中,Dirac函数被替换成了形式为的项,该项对应于矩形窗的中心在的频率响应。一个矩形窗的
11、频率响应形状是一个sinc信号,如下所示该图显示了一个主瓣和若干旁瓣,最大旁瓣大约在主瓣下方13.5dB处。这些旁瓣说明了频谱泄漏效应。无限长信号的功率严格的集中在离散频率点处,而有限长信号在离散频率点附近有连续的功率。因为矩形窗越短,它的频率响应对Dirac冲击的近似性越差,所以数据越短它的频谱泄漏越明显。考虑下面的100个采样的序列randn(state,0)fs = 1000; % Sampling frequencyt = (0:fs/10)/fs; % One-tenth of a second worth of samplesA = 1 2; % Sinusoid amplitud
12、esf = 150;140; % Sinusoid frequenciesxn = A*sin(2*pi*f*t) + 0.1*randn(size(t);Hs = spectrum.periodogram;psd(Hs,xn,Fs,fs,NFFT,1024)注意到频谱泄露只视数据长度而定。周期图确实只对有限数据样本进行计算,但是这和频谱泄露无关。分辨率分辨率指的是区分频谱特征的能力,是分析谱估计性能的关键概念。要区分两个在频率上离得很近的正弦,要求两个频率差大于任何一个信号泄漏频谱的主瓣宽度。主瓣宽度定义为主瓣上峰值功率一半的点间的距离(3dB带宽)。该宽度近似等于两个频率为的正弦信号,可分
13、辨条件是上例中频率间隔10Hz,数据长度要大于100抽才能使得周期图中两个频率可分辨。下图是只有67个数据长度的情况randn(state,0)fs = 1000; % Sampling frequencyt = (0:fs/15)./fs; % 67 samplesA = 1 2; % Sinusoid amplitudesf = 150;140; % Sinusoid frequenciesxn = A*sin(2*pi*f*t) + 0.1*randn(size(t);Hs=spectrum.periodogram;psd(Hs,xn,Fs,fs,NFFT,1024)上述对分辨率的讨论都
14、是在高信噪比的情况进行的,因此没有考虑噪声。当信噪比低的时候,谱特征的分辨更难,而且周期图上会出现一些噪声的伪像,如下所示randn(state,0)fs = 1000; % Sampling frequencyt = (0:fs/10)./fs; % One-tenth of a second worth of samplesA = 1 2; % Sinusoid amplitudesf = 150;140; % Sinusoid frequenciesxn = A*sin(2*pi*f*t) + 2*randn(size(t);Hs=spectrum.periodogram;psd(Hs,
15、xn,Fs,fs,NFFT,1024)估计偏差周期图是对PSD的有偏估计。期望值可以是该式和频谱泄漏中的式相似,除了这里的表达式用的是平均功率而不是幅度。这暗示了周期图产生的估计对应于一个有泄漏的PSD而非真正的PSD。注意本质上是一个三角Bartlett窗(事实是两个矩形脉冲的卷积是三角脉冲。)这导致了最大旁瓣峰值比主瓣峰值低27dB,大致是非平方矩形窗的2倍。周期图估计是渐进无偏的。这从早期的一个观察结果可以明显看出,随着记录数据趋于无穷大,矩形窗对频谱对Dirac函数的近似也就越来越好。然而在某些情况下,周期图法估计很差劲即使数据够长,这是因为周期图法的方差,如下所述。周期图法的方差L趋
16、于无穷大,方差也不趋于0。用统计学术语讲,该估计不是无偏估计。然而周期图在信噪比大的时候仍然是有用的谱估计器,特别是数据够长。Modified Periodogram修正周期图法在fft前先加窗,平滑数据的边缘。可以降低旁瓣的高度。旁瓣是使用矩形窗产生的陡峭的剪切引入的寄生频率,对于非矩形窗,结束点衰减的平滑,所以引入较小的寄生频率。但是,非矩形窗增宽了主瓣,因此降低了频谱分辨率。函数periodogram允许指定对数据加的窗,例如默认的矩形窗和Hamming窗randn(state,0)fs = 1000; % Sampling frequencyt = (0:fs/10)./fs; % O
17、ne-tenth of a second worth of samplesA = 1 2; % Sinusoid amplitudesf = 150;140; % Sinusoid frequenciesxn = A*sin(2*pi*f*t) + 0.1*randn(size(t);Hrect = spectrum.periodogram;psd(Hrect,xn,Fs,fs,NFFT,1024);Hhamm = spectrum.periodogram(Hamming);psd(Hhamm,xn,Fs,fs,NFFT,1024);事实上加Hamming窗后信号的主瓣大约是矩形窗主瓣的2倍。
18、对固定长度信号,Hamming窗能达到的谱估计分辨率大约是矩形窗分辨率的一半。这种冲突可以在某种程度上被变化窗所解决,例如Kaiser窗。非矩形窗会影响信号的功率,因为一些采样被削弱了。为了解决这个问题函数periodogram将窗归一化,有平均单位功率。这样的窗不影响信号的平均功率。修正周期图法估计的PSD是其中U是窗归一化常数假如U保证估计是渐进无偏的。Welch法包括:将数据序列划分为不同的段(可以有重叠),对每段进行改进周期图法估计,再平均。用spectrum.welch对象,或pwelch函数。默认情况下数据划分为4段,50%重叠,应用Hamming窗。取平均的目的是减小方差,重叠会
19、引入冗余但是加Hamming窗可以部分消除这些冗余,因为窗给边缘数据的权重比较小。数据段的缩短和非矩形窗的使用使得频谱分辨率下降。下面的例子展示Welch法的折衷。首先用周期图法估计一个小信噪比下信号的PSD:randn(state,1)fs = 1000; % Sampling frequencyt = (0:0.3*fs)./fs; % 301 samplesA = 2 8; % Sinusoid amplitudes (row vector)f = 150;140; % Sinusoid frequencies (column vector)xn = A*sin(2*pi*f*t) +
20、5*randn(size(t);Hs = spectrum.periodogram(rectangular)psd(Hs,xn,Fs,fs,NFFT,1024);可以看出由于噪声太大,150Hz正弦信号已经无法识别。Hs = spectrum.welch(rectangular,150,50);psd(Hs,xn,Fs,fs,NFFT,512)可以看出两个信号峰,但是如果进一步削减方差,主瓣增宽也使得信号不可识别。Hs = spectrum.welch(rectangular,100,75);psd(Hs,xn,Fs,fs,NFFT,512);Welch法的偏差其中是分段数据的长度,是窗归一化
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 通信电子 通信 电子 频谱 分析 完整版
链接地址:https://www.31ppt.com/p-4152424.html