白噪声通过LTI的仿真.docx
- 文档编号:1669554
- 上传时间:2022-10-23
- 格式:DOCX
- 页数:16
- 大小:180.25KB
白噪声通过LTI的仿真.docx
《白噪声通过LTI的仿真.docx》由会员分享,可在线阅读,更多相关《白噪声通过LTI的仿真.docx(16页珍藏版)》请在冰豆网上搜索。
白噪声通过LTI的仿真
实验2白噪声通过LTI的仿真
1、实验目的
了解白噪声通过LTI系统的原理与处理方法,学会运用Matlab函数对随机过程进展均值、相关函数和功率谱的估计,并且通过实验分析理论分析与实验结果之间的差异。
2、实验原理
假定一具有单位方差的抽样序列{X(n)}的白噪声随机过程X(t)通过一脉冲响应为
的线性滤波器,绘出输入输出信号的均值、方差、相关函数及功率谱密度。
设系统冲激响应为h(n),传递函数
,
或者用Z变换,结果为
。
输入为X(n),输出为
,
均值关系:
,
假设平稳有,
自相关函数关系,
,
当是平稳时候,有
题目中假设为白噪声,可以根据白噪声的性质进展理论计算。
白噪声的自相关函数,
这里,假设的是零均值和单位方差,于是,而
对应的功率谱,
在这里,由于,,a=0.95,可以算出输出信号的方差为,
可以用留数法简单计算出来。
下面对输入输出信号的均值、方差、相关函数及功率谱密度分别进展讨论。
均值变化
输入为白噪声,并且均值为0,按照理论公式,可得到
下面对实际值进展分析:
输入的随机序列,服从标准正态分布。
可以用下面的语句产生
x=randn(1,500);%产生题设的随机序列,长度为500点
系统的冲激响应为,可以用下面的语句产生这个冲激信号:
b=[1];a=[1,-0.5];%设置滤波器的参数,b为分子系数,a为分母系数
h=impz(b,a,20);%得到这个系统的冲激响应,就是题设中的h〔n〕
输入信号通过线性系统,可以通过卷积的方法,或者用filter函数,
y1=filter(b,a,x);%用滤波器的方法,点数为500点
y2=conv(x,h);%通过卷积方法得到,点数为519点
实现的MATLAB代码如下:
clearall;
x=randn(1,500);%产生题设的随机序列,长度为500点
b=[1];a=[1,-0.5];%设置滤波器的参数,b为分子系数,a为分母系数
h=impz(b,a,20);%得到这个系统的冲激响应,就是题设中的h〔n〕
y1=filter(b,a,x);%用滤波器的方法,点数为500点
subplot(2,1,1);
plot(y1,'r');
Title('邹先雄——用滤波器的方法,点数为500点');
x=randn(1,500);
y2=conv(x,h);%通过卷积方法得到,点数为519点
subplot(2,1,2);
plot(y2,'b');
title('邹先雄——通过卷积方法得到,点数为519点');
gridon;
下面画出两者得到波形的区别:
〔为了保持一致,对y2的输出取前500点〕
两者的输出波形近似一致,可以采用任意一个进展分析。
就采用y1进展讨论,
输出均值为:
y1_mean=mean(y1);%进展时间平均,求均值
最终值为-0.0973,与理论的零值有一定误差,考虑到输入随机序列的均值不是0,
m_x=mean(x)=-0.0485,按照上面式子,得到m_y=m_xH(0)=2m_x=-0.0970理论值和实际值是非常吻合的。
附运行结果图:
*因为是随机序列,所以每次运行得到y1和m_x的值也是随机的,但是它们始终满足y1=2m_x
方差变化
输入信号方差的理论值就是1,按照公式,输出的功率谱为
下面对实际值进展分析,用y1_var=var(y1);求得输出均值为1.3598,与理论值的1.3333有差距。
如图:
自相关函数的理论与实际值
理论值为:
在题设中,为白噪声,所以
所以,输出的自相关函数理论值为
可以得到,在零点的值就是1.3333,也就是输出信号的平均功率。
由MATLAb计算的结果为1.3608,这和计算结果非常接近,实际的自相关函数曲线为:
clearall;
x=randn(1,500);%产生题设的随机序列,长度为500点
b=[1];a=[1,-0.5];%设置滤波器的参数,b为分子系数,a为分母系数
h=impz(b,a,20);%得到这个系统的冲激响应,就是题设中的h〔n〕
y1=filter(b,a,x);%用滤波器的方法,点数为500点
y2=conv(x,h);%通过卷积方法得到,点数为519点
Y3=var(y1)
title('自相关函数');
Ry=xcorr(y1,20,'coeff');%进展归一化的自相关函数估计,相关长度为20
n=-20:
1:
20;
stem(n,Ry,'MarkerFaceColor','red');
title('邹先雄——实际的自相关函数曲线');
功率谱密度函数的理论与实际值
对于理论的功率谱密度,可以表示为,
而对于观测数据,可以用功率谱估计的方法得到功率谱密度。
首先,采用Welch法估计信号的功率谱。
它的原理是将数据分成等长度的小段,并且允
许数据的重叠,对每段进展估计,再进展平均,得到信号的功率谱。
在Matlab中有专用
函数pwelch,它的用法是:
[Px,f]=pwelch(X,WINDOW,NOVERLAP,NFFT,Fs,'onesided');%window是采用的数据
窗,NOVERLAP是重叠的数目,NFFT是做FFT的点数,Fs是采样频率,onesided是频率取值。
针对本例,可以用下面语句实现:
window=hamming(20);%采用hanmming窗,长度为20
noverlap=10;%重叠的点数
Nfft=512;%做FFT的点数
Fs=1000;%采样频率,为1000Hz
x=randn(1,500);%产生题设的随机序列,长度为500点
b=[1];a=[1,-0.5];%设置滤波器的参数,b为分子系数,a为分母系数
h=impz(b,a,20);%得到这个系统的冲激响应,就是题设中的h〔n〕
y1=filter(b,a,x);
y1_mean=mean(y1);%进展时间平均,求均值
y1_var=var(y1);%进展时间平均,求方差
Ry=xcorr(y1,20,'coeff');%进展归一化的自相关函数估计,相关长度为20
[Py,f]=pwelch(y1,window,noverlap,Nfft,Fs,'onesided');%估计功率谱密度
f=[-fliplr(f)f(1:
end)];%构造一个对称的频率,围是[-Fs/2,Fs/2]
Py=[-fliplr(Py)Py(1:
end)];%对称的功率谱
plot(f,10*log10(abs(Py)),'r');
title('邹先雄——实际功率谱密度曲线');
gridon;
最后,得到的估计值为
根据上述值,可以计算出理论的功率,由于
可以用下面的语句实现:
w=2*pi*f/Fs;;%转化到数字域上面
H=(1+0.25-2*0.5*cos(w)).^(-1);%系统函数模平方
Gy=H/max(H);%归一化处理
Gy=10*log10(Gy);%化成dB形式
plot(f,Gy,'b');
title('邹先雄——理论功率谱密度曲线');
gridon;
画出的图形见下列图:
这是理论的功率谱密度。
为了方便显示,将两幅图画在一起,便于比拟。
window=hamming(20);%采用hanmming窗,长度为20
noverlap=10;%重叠的点数
Nfft=512;%做FFT的点数
Fs=1000;%采样频率,为1000Hz
x=randn(1,500);%产生题设的随机序列,长度为500点
b=[1];a=[1,-0.5];%设置滤波器的参数,b为分子系数,a为分母系数
h=impz(b,a,20);%得到这个系统的冲激响应,就是题设中的h〔n〕
y1=filter(b,a,x);
y1_mean=mean(y1);%进展时间平均,求均值
y1_var=var(y1);%进展时间平均,求方差
Ry=xcorr(y1,20,'coeff');%进展归一化的自相关函数估计,相关长度为20
[Py,f]=pwelch(y1,window,noverlap,Nfft,Fs,'onesided');%估计功率谱密度
f=[-fliplr(f)f(1:
end)];%构造一个对称的频率,围是[-Fs/2,Fs/2]
Py=[-fliplr(Py)Py(1:
end)];%对称的功率谱
Py=Py/max(Py);%归一化处理
w=2*pi*f/Fs;;%转化到数字域上面
H=(1+0.25-2*0.5*cos(w)).^(-1);%系统函数模平方
Gy=H/max(H);%归一化处理
Gy=10*log10(Gy);%化成dB形式
plot(f,10*log10(abs(Py)),'r',f,Gy,'b');
title('邹先雄——实际功率谱和理论功率谱拟合');
legend('','实际值','理论值');
gridon;
结果为:
从结果上可以看出来,两者存在着比拟大的差距,这是由于输入随机序列的功率谱并不是常数的缘故,也就是输入不是严格的白噪声,所以会出现波动。
当随着数据值的增加,拟合的程度会有所改善。
3、实验容
假定一具有单位方差的抽样序列{X(n)的白噪声随机过程X(t)通过一脉冲响应为
的线性滤波器,利用matlab工具绘出输入输出信号的均值、方差、相关函数及功率谱密度。
实现的MATLAB代码和结果如下:
clearall;
x=randn(1,500);%产生题设的随机序列,长度为500点
b=[1];a=[1,-0.6];%设置滤波器的参数,b为分子系数,a为分母系数
h=impz(b,a,20);%得到这个系统的冲激响应,就是题设中的h〔n〕
y1=filter(b,a,x);%用滤波器的方法,点数为500点
figure
(1)
subplot(2,1,1);
plot(y1,'r');
title('邹先雄——用滤波器的方法,点数为500点');
y2=conv(x,h);%通过卷积方法得到,点数为519点
subplot(2,1,2);
plot(y2,'b');
title('邹先雄——通过卷积方法得到,点数为519点');
Y2=mean(y1)%进展时间平均,求均值,为理论值
m_x=mean(x)%输出的实际值可以通过2m_x计算
Y3=var(y1)
figure
(2)
Ry=xcorr(y1,20,'coeff');%进展归一化的自相关函数估计,相关长度为20
n=-20:
1:
20;
stem(n,Ry,'MarkerFaceColor','red');
title('邹先雄——实际的自相关函数曲线');
%实际功率谱密度
window=hamming(20);%采用hanmming窗,长度为20
noverlap=10;%重叠的点数
Nfft=512;%做FFT的点数
Fs=1000;%采样频率,为1000Hz
x=randn(1,500);%产生题设的随机序列,长度为500点
b=[1]
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 噪声 通过 LTI 仿真
![提示](https://static.bdocx.com/images/bang_tan.gif)