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

用MATLAB进行AR模型功率谱分析随机信号序列x(n)是均值为0方差为1的高斯型白噪声经过AR模型43219606.01697.29403.22137.211zzzzzH后的输出，采样长度为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短时傅里叶变换（ShortTimeFourierTransform,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);00.20.40.60.81050100150200PSDestimatedforAR(3)NormalizedFrequency/RelativePowerMagnititude00.20.40.60.810100200300400PSDestimatedforAR(4)NormalizedFrequency/RelativePowerMagnititude00.10.20.30.40.50.60.70.80.9102004006008001000PSDestimatedforAR(5)NormalizedFrequency/RelativePowerMagnititudey=filter(1,[1-2.21372.9403-2.16970.9606],x);[px,w]=pyulear(y,3,512);subplot(2,2,1);plot(w/pi,px);title('PSDestimatedforAR(3)');xlabel('NormalizedFrequency/{\pi}');ylabel('RelativePowerMagnititude');gridon;[px,w]=pyulear(y,4,512);subplot(2,2,2);plot(w/pi,px);title('PSDestimatedforAR(4)');xlabel('NormalizedFrequency/{\pi}');ylabel('RelativePowerMagnititude');gridon;[px,w]=pyulear(y,5,512);subplot(2,1,2);plot(w/pi,px);title('PSDestimatedforAR(5)');xlabel('NormalizedFrequency/{\pi}');ylabel('RelativePowerMagnititude');gridon;