matlab计算周期信号功率谱密度，功率谱密度相关⽅法的 MATLAB实现
功率谱密度相关⽅法的MATLAB实现
1. 基本⽅法 周期图法是直接将信号的采样数据x(n)进⾏Fourier变换求取功率谱密度估计的⽅法。假定有限长随机信号序列为x(n)。它的Fourier变换和功率谱密度估计存在下⾯的关系： 式中，N为随机信号序列x(n)的长度。在离散的频率点f=kΔf，有： 其中，FFT[x(n)]为对序列x(n)的Fourier变换，由于FFT[x(n)]的周期为N，求得的功率谱估计以N为周期，因此这种⽅法称为周期图法。下⾯⽤例⼦说明如何采⽤这种⽅法进⾏功率谱 ⽤有限长样本序列的Fourier变换来表⽰随机序列的功率谱，只是⼀种估计或近似，不可避免存在误差。为了减少误差，使功率谱估计更加平滑，可采⽤分段平均周期图法(Bartlett法)、加窗平均周期图法(Welch法)等⽅法加以改进。 2. 分段平均周期图法(Bartlett法) 将信号序列x(n)，n=0,1,…,N-1，分成互不重叠的P个⼩段，每⼩段由m个采样值，则P*m=N。对每个⼩段信号序列进⾏功率谱估计，然后再取平均作为整个序列x(n)的功率谱估计。 平均周期图法还可以对信号x(n)进⾏重叠分段，如按2:1重叠分段，即前⼀段信号和后⼀段信号有⼀半是重叠的。对每⼀⼩段信号序列进⾏功率谱估计，然后再取平均值作为整个序列x(n)的功率谱估计。这两种⽅法都称为平均周期图法，⼀般后者⽐前者好。程序运⾏结果为图9-5，上图采⽤不重叠分段
法的功率谱估计，下图为2:1重叠分段的功率谱估计，可见后者估计曲线较为平滑。与上例⽐较，平均周期图法功率谱估计具有明显效果(涨落曲线靠近0dB)。
大连理工大学bbs
3.加窗平均周期图法 加窗平均周期图法是对分段平均周期图法的改进。在信号序列x(n)分段后，⽤⾮矩形窗⼝对每⼀⼩段信号序列进⾏预处理，再采⽤前述分段平均周期图法进⾏整个信号序列x(n)的功率谱估计。由窗函数的基本知识(第7章)可知，采⽤合适的⾮矩形窗⼝对信号进⾏处理可减⼩“频谱泄露”，同时可增加频峰的宽度，从⽽提⾼频谱分辨率。 其中上图采⽤⽆重叠数据分段的加窗平均周期图法进⾏功率谱估计，⽽下图采⽤重叠数据分段的加窗平均周期图法进⾏功率谱估计，显然后者是更佳的，信号谱峰加宽，⽽噪声谱均在0dB附近，更为平坦(注意采⽤⽆重叠数据分段噪声的最⼤的下降分贝数⼤于5dB，⽽重叠数据分段周期图法噪声的最⼤下降分贝数⼩于5dB)。 导电浆料
巡回法庭4. Welch法估计及其MATLAB函数Welch功率谱密度就是⽤改进的平均周期图法来求取随机信号的功率谱密度估计的。Welch法采⽤信号重叠分段、加窗函数和FFT算法等计算⼀个信号序列的⾃功率谱估计(PSD如上例中的下半部分的求法)和两个信号序列的互功率谱估计(CSD)。 MATLAB信号处理⼯具箱函数提供了专门的函数PSD和CSD⾃动实现Welch法估计，⽽不需要⾃⼰编程。 (1) 函数psd利⽤Welch法估计⼀个信号⾃功率谱密度，函数调⽤格式为： [Pxx[,f]]=psd(x[,Nfft,Fs,window,Noverlap,’dflag’]) 式中，x为信号序列；Nfft为采⽤的FFT长度。这⼀值决定了功率谱估计速度，当Nfft采⽤2的幂时，程序采⽤快速算法；Fs为采样频率；Window定义窗函
数和x分段序列的长度。窗函数长度必须⼩于或等于Nfft，否则会给出错误信息；Noverlap为分段序列重叠的采样点数(长度)，它应⼩于Nfft;dflag为去除信号趋势分量的选择
项：’linear’，去除线性趋势分量，’mean’去除均值分量，’none’不做去除趋势处理。Pxx为信号x的⾃功率谱密度估计。f为返回的频率向量，它和Pxx对应，并且有相同长度。 在psd函数调⽤格式中，缺省值为：Nfft=min(256,length(x)),Fs=2Hz,
全球紧急卫生事件window=hanning(Nfft),noverlap=0. 若x是实序列，函数psd仅计算频率为正的功率 注意程序前半部分中频率向量f的创建⽅法。它与函数psd的输出Pxx长度的关系如下：若x为实序列，当Nfft为奇数时，f=(0:(Nfft+1)/2-1)/Nfft；当Nfft为偶数时，f=(0:Nfft/2)/Nfft。 函数还有⼀种缺省返回值的调⽤格式，⽤于直接绘制信号序列x的功率谱估计曲线。 函数还可以计算带有置信区间的功率谱估计，调⽤格式为： [Pxx,Pxxc,f]=psd(x,Nfft,Fs,window,Noverlap,p) 式中，p为置信区间，0<=p<=1。 由此可知，滤波器输⼊⽩噪声序列的输出信号的功率谱或⾃相关可以确定滤波器的频率特性。 (2)函数csd利⽤welch法估计两个信号的互功率谱密度，函数调⽤格式为：
w980i
[Pxy[,f]]=csd(x,y,Nfft,Fs,window,Noverlap,’dflag’) [Pxy,Pxyc[,f]]=csd(x,y,Nfft,Fs,window,Noverlap,p) 这⾥，x,y为两个信号序列；Pxy为x,y的互功率谱估计；其他参数的意义同⾃功率谱函数psd。 可以看到，两个⽩噪声信号的互功率谱(上图)杂乱⽆章，看不出周期成分，⼤部分功率谱在-5dB以下。然⽽
小文件存储⽩噪声与带有噪声的周期信号的功率谱在其周期(频率为1000Hz)处有⼀峰值，清楚地表明了周期信号的周期或频率。因此，利⽤未知信号与⽩噪声信号的互功率谱也可以检测未知信号中所含有的频率成分。 5 多 窗 ⼝ 法 多窗⼝法(Multitaper ，简称MTM法)利⽤多个正交窗⼝(Tapers)获得各⾃独⽴的近似功率谱估计，然后综合这些估计得到⼀个序列的功率谱估计。相对于普通的周期图法，这种功率谱估计具有更⼤的⾃由度，并在估计精度和估计波动⽅⾯均有较好的效果。普通的功率谱估计只利⽤单⼀窗⼝，因此在序列始端和末端均会丢失相关信息，⽽且⽆法回。⽽MTM法估计增加窗⼝使得丢失的信息尽量减少。 MTM法简单地采⽤⼀个参数：时间带宽积(Time-bandwidth product)NW，这个参数⽤以定义计算功率谱所⽤窗的数⽬，为2*NW-1。NW越⼤，功率谱计算次数越多，时间域分辨率越⾼，⽽频率域分辨率降低，使得功率谱估计的波动减⼩。随着NW增⼤，每次估计中谱泄漏增多，总功率谱估计的偏差增⼤。对于每⼀个数据组，通常有⼀个最优的NW使得在估计偏差和估计波动两⽅⾯求得折中，这需要在程序中反复调试来获得。MATLAB信号处理⼯具箱中函数PMTM就是采⽤MTM法估计功率谱密度。函数调⽤格式为： [Pxx[,f]]=pmtm(x[,nw,Nfft,Fs]) 式中，x为信号序列；nw为时间带宽积，缺省值为4。通常可取2，5/2，3，7/2；Nfft为FFT长度；Fs为采样频率。 上⾯的函数还可以通过⽆返回值⽽绘出置信区间，如pmtm(x,nw,Nfft,Fs,’option’,p)绘制带置信区间的功率谱密度估计曲线，0<=p<=1。 6 最 ⼤ 熵 法(Maxmum entropy ， MEM
豫公网安备41160202000603
站长QQ:729038198
关于我们
投诉建议