用MATLAB进行AR模型功率谱的分析.docxVIP

用MATLAB进行AR模型功率谱分析随机信号序列x(n)是均值为0方差为1的高斯型白噪声经过AR模型后的输出，采样长度为512，AR模型阶次取3，4，5，用L-D算法估计功率谱密度。分析：MATLAB函数pyulear()的用法pyulear()是基于自相关法、利用Levesion-Durbin算法估计功率谱密度。[px,w]=pyulear(x,p,[nfft],’range’)x为随即信号序列，是由白噪声经AR模型产生的，在MATLAB中可以由白噪声序列u经过表示AR模型的数字滤波器后得到，使用的是filter函数；p为AR模型阶次；nfft为由模型参数计算频谱时的频域采样点数，默认为256；range用于选择输出是为单边[0,π]，还是双边[0,2π]；w的范围[0,π]，还是 [0,2π]由range确定或由nfft的奇偶性确定；该函数返回实际频率w下的功率谱密度向量，w的单位即为rad/sample，默认sample为1Hz，若要转化为归一化频率，只需用w/π即可。实验结果如图三.1（对应程序为shiyan3.m）：图.1短时傅里叶变换（Short Time Fourier Transform, STFT）法，在MATLAB中做短时傅里叶变换的函数为spectrogram：spectrogram(x,window,overlap,f,fs)[s,f,t,p]=spectrogram(x,window,overlap,f,fs)x为被分析序列，window为窗函数及长度，默认为hamming窗，overlap为相邻两个短时序列之间重叠的数据点数，f为一向量，确定在某一个频率范围内做短时傅里叶变换，fs为采样频率。附程序：%%------------------------------------------------------------------------%%功能：随机信号AR模型功率谱分析方法%%------------------------------------------------------------------------x=randn(1,512);y=filter(1,[1 -2.2137 2.9403 -2.1697 0.9606],x);[px,w]=pyulear(y,3,512);subplot(2,2,1);plot(w/pi,px);title(PSD estimated for AR(3));xlabel(Normalized Frequency/{\pi});ylabel(Relative Power Magnititude);grid on;[px,w]=pyulear(y,4,512);subplot(2,2,2);plot(w/pi,px);title(PSD estimated for AR(4));xlabel(Normalized Frequency/{\pi});ylabel(Relative Power Magnititude);grid on;[px,w]=pyulear(y,5,512);subplot(2,1,2);plot(w/pi,px);title(PSD estimated for AR(5));xlabel(Normalized Frequency/{\pi});ylabel(Relative Power Magnititude);grid on;