1、1随机信号分析实验报告学院:电信学院班级:电子 3 班姓名:王敏学号:082603352实验名称: 随机过程的特征估计实验目的:了解随机过程特征估计的基本概念和方法,学会运用 Matlab 函数对随机过程进行特征估计,并且通过实验了解不同估计方法所估计出来结果之间的差异. 实验原理:设随机序列 X(n),Y(n)为各态历经过程,样本分别为 x(n),y(n) (n=0,1,N-1).1. 均值的估计2. 方差的估计 方差估计有两种情况, 如果均值 已知,则如果均值未知,那么3. 相关函数估计4. 功率谱估计功率谱的估计有几种方法,(1) 自相关法先求相关函数的估计,然后对估计的相关函数做傅立叶
2、变换,3(2) 周期图法先对序列 x(n)做傅立叶变换,则功率谱估计为周期图法是一种非参数谱估计方法,另外还有一种修正的周期图方法,也叫 Yelch 法,MATLAB有周期图和 Welch 法的谱估计函数 .实验内容 : 1、 产生一组均值为1,方差为4 的正态分布的随机序列(1000 个样本),估计该序列的均值与方差。X=-10:0.1:10;Y=normpdf(X,1,2);stem(X,Y);xlabel(X);ylabel(Y);title(;2、按如下模型产生一组随机序列:x(n)=0.8x(n-1)+w(n)其中w(n)为均值为1,方差为4 的正态分布白噪声序列。估计过程的自相关函
3、数与功率谱。自相关函数4clear all; b=1; a=1 -0.8;h=impz(b,a,20); randn(state,0); w=normrnd(0,2,1,500); x=filter(b,a,w);subplot(1,2,1);plot(x,r);ylabel(x(n);title(AR);grid Mlag=20;Rx=xcorr(x,Mlag,coeff); m=-Mlag:Mlag; stem(m,Rx,r.);clear all; delta=2*pi/1000;w_min=-pi;w_max=pi;Fs=1000;w=w_min:delta:w_max;NS=1DS=
4、1-0.8*exp(-i*w); Gx=4*(abs(NS./DS).2); Gx=Gx/max(Gx);f=w*Fs/(2*pi);window=hamming(20);noverlap=10; Nfft=512; Fs=1000;Px,f=pwelch(x,window,noverlap,Nfft,Fs, onesided);plot(f,10*log10(Px),b);53、设信号为x(n)=sin(2f1n)+2cos(2f2n)+w(n),n=1,2,N,其中f1=0.05, f2=0.12,w(n)为正态白噪声,试在N=356 和1024 点时,分别产生随机序列x(n)、画出 x(
5、n)的波形并估计 x(n)的相关函数和功率谱。clear all;n=1:1:356;w = randn(1,356);x=sin(2*pi*0.05*n)+2*cos(2*pi*0.12*n)+w;plot(x,r);Mlag=20; Rx=xcorr(x,Mlag,coeff);m=-Mlag:Mlag;stem(m,Rx,r.)window=hamming(20); noverlap=10; Nfft=512;Fs=1000; Px,f=pwelch(x,window,noverlap,Nfft,Fs, onesided); plot(f,10*log10(Px),b)6clear al
6、l;n=1:1:1024;w = randn(1,1024);x=sin(2*pi*0.05*n)+2*cos(2*pi*0.12*n)+w;plot(x,r);Mlag=20; Rx=xcorr(x,Mlag,coeff);m=-Mlag:Mlag;stem(m,Rx,r.)window=hamming(20); noverlap=10; Nfft=512;7Fs=1000; Px,f=pwelch(x,window,noverlap,Nfft,Fs, onesided); plot(f,10*log10(Px),b)8实验总结通过做随机过程的特征估计实验,我了解了如何运用 MATLAB 函
7、数对随机过程进行特征估计,了解了不同估计方法所估计出来结果之间的差异。9实验名称 典型时间序列模型分析1、实验目的熟悉三种典型的时间序列模型:AR 模型,MA 模型与ARMA 模型,学会运用Matlab工具对对上述三种模型进行统计特性分析,通过对2 阶模型的仿真分析,探讨几种模型的适用范围,并且通过实验分析理论分析与实验结果之间的差异。2、实验原理AR 模型分析设有AR(2)模型,X(n)=-0.3X(n-1)-0.5X(n-2)+W(n)W(n)是零均值正态白噪声,方差为4。(1)用MATLAB 模拟产生X(n)的500 观测点的样本函数,并绘出波形(2)用产生的500 个观测点估计X(n)
8、的均值和方差(3)画出理论的功率谱(4)估计X(n)的相关函数和功率谱1. 产生样本函数,并画出波形题目中的AR 过程相当于一个零均值正态白噪声通过线性系统后的输出,可以按照上面的方法进行描述。clear all;b=1; a=1 0.3 0.5; % 由描述的差分方程,得到系统传递函数h=impz(b,a,20); % 得到系统的单位冲激函数,在20 点处已经可以认为值是0randn(state,0);w=normrnd(0,2,1,500); % 产生题设的白噪声随机序列,标准差为2x=filter(b,a,w); % 通过线形系统,得到输出就是题目中要求的2 阶AR 过程plot(x,r
9、);ylabel(x(n);title(产生的AR 随机序列);grid得到的输出序列波形为:2. 估计均值和方差可以首先计算出理论输出的均值和方差,得到0 x m = ,对于方差可以先求出理论自相关输出,然后取零点的值。10可以采用上面介绍的方法,对式中的卷积进行计算。计算出的卷积输出图形为:在最大值处就是输出的功率,也就是方差,为对实际数据进行估计,均值为mean(x)= -0.0703,而方差为var(x)= 5.2795,两者和理论值吻合的比较好。3. 画出理论的功率谱密度曲线理论的功率谱为,用下面的语句产生:delta=2*pi/1000;w_min=-pi;w_max=pi;Fs=
10、1000;w=w_min:delta:w_max; % 得到数字域上的频率取样点,范围是-pi,piGx=4*(abs(1./(1+0.3*exp(-i*w)+0.5*exp(-2*i*w).2); % 计算出理论值Gx=Gx/max(Gx); % 归一化处理f=w*Fs/(2*pi); % 转化到模拟域上的频率plot(f,Gx,b),grid on;得到的图形为:11可以看出,这个系统是带通系统。4. 估计自相关函数和功率谱密度用实际数据估计自相关函数和功率谱的方法前面已经讨论过,在这里仅给出最后的仿真图形。% 计算理论和实际的自相关函数序列Mlag=20; % 定义最大自相关长度Rx=x
11、corr(x,Mlag,coeff);m=-Mlag:Mlag;stem(m,Rx,r.);最终的值为可以看出,它和上面的理论输出值吻合程度很好。实际的功率谱密度可以用类似于上面的方法进行估计,window=hamming(20); % 采用hanmming 窗,长度为20noverlap=10; % 重叠的点数Nfft=512; % 做FFT 的点数Fs=1000; % 采样频率,为1000HzPx,f=pwelch(x,window,noverlap,Nfft,Fs, onesided);f=-fliplr(f) f(2:end); Py=-fliplr(Py) Py(2:end); pl
12、ot(f,10*log10(Py),b);估计出来的功率谱密度为,12将两幅图画在一起,可以看到拟合的情况比较好:ARMA 模型分析设有ARMA(2,2)模型,X(n)+0.3X(n-1)-0.2X(n-2)=W(n)+0.5W(n-1)-0.2W(n-2)W(n)是零均值正态白噪声,方差为4。(1)用MATLAB 模拟产生X(n)的500 观测点的样本函数,并绘出波形(2)用产生的500 个观测点估计X(n)的均值和方差(3)画出理论的功率谱(4)估计X(n)的相关函数和功率谱1. 产生样本函数,并画出波形题目中的ARMA 过程相当于一个零均值正态白噪声通过线性系统后的输出,可以按照上面的方
13、法进行描述。clear all;b=1 0.5 -0.2; a=1 0.3 -0.2; % 由描述的差分方程,得到系统传递函数h=impz(b,a,10); % 得到系统的单位冲激函数,在10 点处已经可以认为值是0randn(state,0);w=normrnd(0,2,1,500); % 产生题设的白噪声随机序列,标准差为2x=filter(b,a,w); % 通过线形系统,得到输出就是题目中要求的(2,2) 阶ARMA 过程plot(x,r);13得到的输出序列波形为2. 估计均值和方差可以首先计算出理论输出的均值和方差,得到0 x m = ,对于方差可以先求出理论自相关输出,然后取零点
14、的值。可以采用上面介绍的方法,对式中的卷积进行计算。计算出的卷积输出图形为:在最大值处就是输出的功率,也就是方差,为对实际数据进行估计,均值为mean(x)= -0.0547,而方差为var(x)=3.8,两者和理论值吻合的比较好。3. 画出理论的功率谱密度曲线理论的功率谱为,14用下面的语句产生:delta=2*pi/1000;w_min=-pi;w_max=pi;Fs=1000;w=w_min:delta:w_max; % 得到数字域上的频率取样点,范围是-pi,piNS=1+0.5*exp(-i*w)-0.2*exp(-2*i*w); % 分子DS=1+0.3*exp(-i*w)-0.2
15、*exp(-2*i*w); % 分母Gx=4*(abs(NS./DS).2); % 计算出理论值Gx=Gx/max(Gx);f=w*Fs/(2*pi); % 转化到模拟域上的频率plot(f,Gx,b),grid on;4. 估计相关函数和功率谱密度曲线用实际数据估计自相关函数和功率谱的方法前面已经讨论过,在这里仅给出仿真图形。% 计算理论和实际的自相关函数序列Mlag=20; % 定义最大自相关长度Rx=xcorr(x,Mlag,coeff);m=-Mlag:Mlag;stem(m,Rx,r.);最终的值为实际的功率谱密度可以用类似于上面的方法进行估计,window=hamming(20);
16、 % 采用hanmming 窗,长度为20noverlap=10; % 重叠的点数Nfft=512; % 做FFT 的点数15Fs=1000; % 采样频率,为1000HzPx,f=pwelch(x,window,noverlap,Nfft,Fs, onesided); % 估计功率谱密度f=-fliplr(f) f(2:end); % 构造一个对称的频率,范围是-Fs/2, Fs/2Py=fliplr(Py) Py(2:end); % 对称的功率谱plot(f,10*log10(Py),b);估计出来的功率谱密度为,把两幅图画在一起,可以得到下面的图形,可以看出两者的吻合度比较高。3、实验内
17、容1、 熟悉实验原理,将实验原理上的程序应用matlab 工具实现;其实现见试验原理。2、设有MA(2)模型, x(n) =W(n)-0.3W(n-1)+0.2W(n-2)W(n)是零均值正态白噪声,方差为4。(1) 用MATLAB 模拟产生 X(n)的500 观测点的样本函数,并绘出波形(2) 用产生的500 个观测点估计X(n)的均值和方差(3) 画出理论的功率谱(4) 估计X(n)的相关函数和功率谱clear all; b=1 -0.3 0.2; a=1; % 由描述的差分方程,得到系统传递函数h=impz(b,a,10); % 得到系统的单位冲激函数,在 10 点处已经可以认为值是 0
18、 randn(state,0); w=normrnd(0,2,1,500); % 产生题设的白噪声随机序列,标准差为 2 x=filter(b,a,w); % 通过线形系统,得到输出就是题目中要求的(2,2) 阶 ARMA 过程16plot(x,r);mean(w)sigma2=var(w)m = -0.1259sigma2 =3.5852delta=2*pi/1000;w_min=-pi;w_max=pi;Fs=1000; w=w_min:delta:w_max; NS=1+0*exp(-i*w)+0*exp(-2*i*w); DS=1-0.3*exp(-i*w)+0.2*exp(-2*i*
19、w); Gx=4*(abs(NS./DS).2); Gx=Gx/max(Gx);f=w*Fs/(2*pi); plot(f,Gx,b),grid on;clear all; b=1 -0.3 0.2 ; a=1 ; h=impz(b,a,10); randn(state,0); w=normrnd(0,2,1,500); x=filter(b,a,w); Mlag=20; Rx=xcorr(x,Mlag,coeff); m=-Mlag:Mlag; stem(m,Rx,r.); delta=2*pi/1000;w_min=-pi;w_max=pi;Fs=1000;w=w_min:delta:w_
20、max; NS=1-0.3*exp(-i*w)+0.2*exp(-2*i*w); DS=1Gx=4*(abs(NS./DS).2); Gx=Gx/max(Gx);f=w*Fs/(2*pi);window=hamming(20); noverlap=10; Nfft=512;Fs=1000; Px,f=pwelch(x,window,noverlap,Nfft,Fs, onesided); plot(f,10*log10(Px),b);17184、实验总结通过做典型时间序列模型分析,了解了三种典型的时间序列模型:AR 模型,MA 模型与 ARMA 模型,通过对 2 阶模型的仿真分析,知道了 MA
21、 模型的适用范围,并且通过实验观察得出了其差异。19实验名称 白噪声通过LTI 的仿真1、实验目的了解白噪声通过LTI 系统的原理与处理方法,学会运用Matlab 函数对随机过程进行均值、相关函数和功率谱的估计,并且通过实验分析理论分析与实验结果之间的差异。2、实验原理假定一具有单位方差的抽样序列X(n)的白噪声随机过程X(t) 通过一脉冲响应为的线性滤波器,绘出输入输出信号的均值、方差、相关函数及功率谱密度。设系统冲激响应为 h(n),传递函数为题目中假设为白噪声,可以根据白噪声的性质进行理论计算。20均值变化输入为白噪声,并且均值为0,按照理论公式,可得到下面对实际值进行分析:输入的随机序
22、列,服从标准正态分布。可以用下面的语句产生x = randn(1,500); % 产生题设的随机序列,长度为500 点系统的冲激响应为 h(n) =(0.5)nu(n),可以用下面的语句产生这个冲激信号: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 点下面画出两者得到波形的区别
23、:(为了保持一致,对y2 的输出取前500 点)两者的输出波形近似一致,可以采用任意一个进行分析。就采用y1 进行讨论,输出均值为:y1_mean=mean(y1); % 进行时间平均,求均值最终值为-0.1248,与理论的零值有一定误差,考虑到输入随机序列的均值不是0,理论值和实际值是非常吻合的。21方差变化输入信号方差的理论值就是1,按照公式,输出的功率谱为下面对实际值进行分析,y1_var=var(y1); % 进行时间平均,求均值得到,输出均值为1.1745,与理论值的1.3333 有差距。自相关函数的理论与实际值理论值为:在题设中,为白噪声,所以所以,输出的自相关函数理论值为用图形表
24、示,上图分别是 h(n)和 h(n),最终得到的卷积值为:可以看到,在零点的值就是1.3333,也就是输出信号的平均功率。3、实验内容22假定一具有单位方差的抽样序列X(n)的白噪声随机过程X(t) 通过一脉冲响应为的线性滤波器,利用matlab 工具绘出输入输出信号的均值、方差、相关函数及功率谱密度。x = randn(1,500); plot(x);b=1;a=1,-0.6;h=impz(b,a,20);plot(h); y=filter(b,a,x); plot(y);x_mean=mean(x)x_var=var(x)y_mean=mean(y)y_var=var(y) Mlag=20
25、; Rx=xcorr(x,Mlag,coeff); m=-Mlag:Mlag; stem(m,Rx,r.); Mlag=20; Rx=xcorr(y,Mlag,coeff); m=-Mlag:Mlag; stem(m,Rx,r.);window=hamming(20); noverlap=10; Nfft=512; Fs=1000; Px,f=pwelch(x,window,noverlap,Nfft,Fs,onesided); plot(f,10*log10(abs(Px),r),grid on; window=hamming(20);noverlap=10; Nfft=512; Fs=1000; Py,f=pwelch(y,window,noverlap,Nfft,Fs,onesided); plot(f,10*log10(abs(Py),r),grid on; 2324x_mean =-0.0294x_var =1.0716y_mean =-0.0707y_var =1.6139254、实验总结通过做白噪声通过 LTI 的仿真实验,了解了许多,可以用功率谱估计的方法得到功率谱密度,可以采用 Welch 法估计信号的功率谱。它的原理是将数据分成等长度的小段,并且允许数据的重叠,对每段进行估计,再进行平均,得到信号的功率谱。对白噪声的仿真有了更好的了解。