matlab数信号处理周期图程序.docx
《matlab数信号处理周期图程序.docx》由会员分享，可在线阅读，更多相关《matlab数信号处理周期图程序.docx（15页珍藏版）》请在冰豆网上搜索。
matlab数信号处理周期图程序
1、假设一平稳随机信号为
，其中
是均值为0，方差为1的白噪声，数据长度为1024。
（1）、产生符合要求的
和
；
（2）、给出信号x（n）的理想功率谱；
（3）、编写周期图谱估计函数，估计数据长度N=1024及256时信号功率谱，分析估计效果。
（4）、编写Bartlett平均周期图函数，估计当数据长度N=1024及256时，分段数L分别为2和8时信号
的功率谱，分析估计效果。
解:
思路
在matlab中提供的有randn（m.n）函数，其为均值为零，方差为1的函数，所以w（n）可以通过随机序列randn（1,N）来产生，x（n）可以通过对w（n）滤波产生，也可以直接由递推式迭代产生。
由于线性系统的输出功率谱等于输入功率谱乘以传递函数模的平方，X（n）可以看做w（n）
通过一线性系统的输出，
H（z）=1/（1-0.8z）
所以x（n）的理想功率谱
直接产生
Matlab程序：
>>clear;closeall;
F=500;%采样率
N=1024;%观测数据
>>subplot（2,1,1）;
w=sqrt
（1）+randn（1,N）;
>>plot（w）;
>>xlabel（'观察次数'）；ylabel（'功率db'）；title（'白噪声w的分布情况'）;
>>subplot（2,1,2）;
>>x=[w
（1）zeros（1,N-1）];%初始化x（n）,长度1024,x
（1）=w
（1）
fori=2:
N
x（i）=0.8*x（i-1）+w（i）;%迭代产生观测数据x（n）
end
>>plot（x）;
>>xlabel（'观测次数'）;ylabel（'x的分布情况'）;title（'x的分布情况'）;
%%理想功率谱
cn=xcorr（x,x）;%自相关函数
Cn=fft（cn,N）;%快速傅里叶变化
Pxx=abs（Cn）;%取绝对值
index=0:
round（N/2-1）;
k=index*Fs/N;
pxx=10*log10（Pxx（index+1））;%转换成db
figure;
plot（k,pxx）;
>>xlabel（'频率'）;ylabel（'功率db'）;
>>title（'理想功率谱'）；
%%周期图谱
%1024
Pxx1=abs（fft（x））.^2/N;
>>pxx1=10*log10（Pxx1（index+1））;
>>figure;
>>plot（k,pxx1）;
>>title（'周期图1024个点'）;
>>xlabel（'频率'）;ylabel（'功率db'）；
%周期图256个观测点
>>x1=x（1:
4:
N）;
>>Pxx2=abs（fft（x1,1024））.^2/N;
>>pxx2=10*log10（Pxx2（index+1））;
>>figure;
>>plot（k,pxx2）;
>>title（'周期图256'）；xlabel（'频率'）;ylabel（'功率db'）;
%%L=2N=1024
>>L=2;
>>x1=x（1:
L:
N）;
>>x2=x（2:
L:
N）;
>>pxx2_1=abs（fft（x1,1024））.^2/length（x1）;
>>pxx2_2=abs（fft（x2,1024））.^2/length（x2）;
>>pxx_2=（pxx2_1+pxx2_2）/L;
>>figure;
>>subplot（2,1,1）;
>>plot（k,10*log10（pxx_2（index+1）））;
>>title（'N=1024,L=2时的周期图'）;
>>xlabel（'频率'）;
>>ylabel（'功率db'）;
>>L1=8;%当数据长度为8时
x3=zeros（L1,N/L1）;
fori=1:
L1
x3（i,:
）=x（i:
L1:
N）;
end
>>pxx3=zeros（L1,N）;
fori=1:
L1
pxx3（i,:
）=abs（fft（x3（i,:
）,1024））.^2/length（x3（i,:
））;
end
>>fori=1:
1024
Pxx3（i）=sum（pxx3（:
i））/L1;
end
>>subplot（2,1,2）;
plot（k,10*log10（Pxx3（index+1）））;
title（'N=1024,L=8时的周期图'）;
xlabel（'频率'）;
ylabel（'功率db'）;
当长度为256时，指针函数发生变化
clear;closeall;
Fs=500;%采样率
N=256;%观测数据
w=sqrt
（1）+randn（1,N）;%
x=[w
（1）zeros（1,N-1）];
fori=2:
N
x（i）=0.8*x（i-1）+w（i）;
end
index=0:
round（N/2-1）;
k=index*Fs/N;
>>L=2;
>>x1=x（1:
L:
N）;
>>x2=x（2:
L:
N）;
>>pxx2_1=abs（fft（x1,1024））.^2/length（x1）;
>>pxx2_2=abs（fft（x2,1024））.^2/length（x2）;
>>pxx_2=（pxx2_1+pxx2_2）/L;
>>figure;
>>subplot（2,1,1）;
>>plot（k,10*log10（pxx_2（index+1）））;
>>title（'N=256,L=2时的周期图'）;
>>xlabel（'频率'）;
>>ylabel（'功率db'）;
>>L1=8;%当数据长度为8时
x3=zeros（L1,N/L1）;
fori=1:
L1
x3（i,:
）=x（i:
L1:
N）;
end
>>pxx3=zeros（L1,N）;
fori=1:
L1
pxx3（i,:
）=abs（fft（x3（i,:
）,1024））.^2/length（x3（i,:
））;
end
>>fori=1:
1024
Pxx3（i）=sum（pxx3（:
i））/L1;
end
>>subplot（2,1,2）;
plot（k,10*log10（Pxx3（index+1）））;
title（'N=256,L=8时的周期图'）;
xlabel（'频率'）;
ylabel（'功率db'）;
结果:
理想功率谱：
长度为256跟长度为1024时候的周期图谱：
N=1024时，当L=2与L=8时的周期图
当N的长度为256时，L=2与L=8时的周期图谱
从上面的图像可以看出，周期图法得到的功率谱估计，谱线的起伏较大，即估计所得的均方误差较大。
当N加时，摆动的频率加快，而摆动的幅度变化不大。
且N=1024时，谱的分辨率较N=256时大。
采用平均处理后，谱线上下摆动的幅度减小（即均方误差有所降低），
曲线的平滑性也较周期图法好。
N相同时，分段越多，方差越小，曲线越平滑。
这是因为，
N一定时L加大，每一段的数据量就会相应减少，因此估计方差减小，偏移加大，从而分
辨率降低。
N一定时，L与每段数据量相互矛盾，需择中选取。
综上，Bartlett法相对于周期图法来说，较好地减小了估计误差。
2、假设均值为0，方差为1的白噪声
中混有两个正弦信号，该正弦信号的频率分别为100Hz和110Hz，信噪比分别为10dB和30dB，初始相位都为0，采样频率为1000Hz。
（1）、采用自相关法、Burg法、协方差法、修正协方差法估计功率谱，分析数据长度和模型阶次对估计结果的影响（可采用MATLAB自带的功率谱分析函数）。
（2）、调整正弦信号信噪比，分析信噪比的降低对估计效果的影响。
Matlab源程序
clc;
clear;closeall;
fs=1000;%采样率1000
N=1;%改变数据长度
p=50;%AR模型阶数
nfft=512;%fft长度
t=0:
1/fs:
N;
wn=sqrt
（1）+randn（1,N*fs+1）;%白噪声,均值0,方差1
s1=sqrt（20）*sin（2*pi*100*t）;%正弦信号1,信噪比10db
s2=sqrt（2000）*sin（2*pi*110*t）;%正弦信号2,信噪比30db
x=s1+s2+wn;%观测数据
%figure,plot（t,x）;
x1=xcorr（x,'biased'）;
[Pxx,f]=pyulear（x1,p,nfft,fs）;%Yule-Walker方程
figure,plot（f,10*log10（Pxx））;gridon;
title（'自相关法'）;
[Pxx1,f1]=pcov（x,p,nfft,fs）;
figure,plot（f1,10*log10（Pxx1））;gridon;
title（'协方差法'）;
[Pxx2,f2]=pmcov（x,p,nfft,fs）;
figure,plot（f2,10*log10（Pxx2））;gridon;
title（'修正协方差法'）;
[Pxx3,f3]=pburg（x,p,nfft,fs）;figure,plot（f3,10*log10（Pxx3））;gridon;
title（'Burg法'）;
结果：
长度为N=1，p=50时的图像，
长度N=5，p=50是的图片
N=1，p=100时的图谱
N=1，p=80的图谱
分析数据长度和模型阶次对估计结果的影响：
从上面的图像可以看出，当数据的长度变长时，图像的分辨率降到，平滑性争强，模型的阶次加大时，图形的分变率加大，平滑性减弱！
降低信噪比前：
信噪比下降一些：
信噪比再次下降：
通过上面三副图的比较我们可以看到，。
所以相同阶次时，信噪比越高，分辨率越高，信噪比降低则会导致分辨率降低。
对于已知信噪比的观测信号来说，若信噪比低则选用较高的阶次来得到较好地分辨率；信噪比高，则相应地可以减小阶次，以得到较小的误差。
升级会员 