随机过程的计算机模拟

实验2随机过程的计算机模拟

1、实验目的

给定功率谱（相关函数）和概率分布，通过计算机模拟分析产生相应的随机过程，然后

又通过该随即过程的实际功率谱（相关函数）和概率分布验证该实验的有效性

2、实验原理

।、标准正态分布随机序列的产生方法:

利用随机变量函数变换的方法。

设八，G为两个相互独立的（°，1）均匀分布的随机数，如果要产生服从NO%。?）的

随机数X,则可以按如下变换关系产生：

x=m+（y^-2Inrtcos2itr2（1）

2、正态随机矢量的模拟

设有一N维正态随机矢量X=[X-X2,…XNT，其概率密度为

九口,…/Mg乃y\；K「exp{——⑵

其中x=,M=为随机矢量的均值,

“11%121"卜皿

K一女2I女22%2N⑶

KN\^N2"<£VN_

为协方差矩阵，且是正定的。

对于上述正态随机矢量的模拟，首先可以模拟产生一个零均值、单位方差、且各个分量

相互独立的标准正态随机矢量U,然后作如下变换：

X=AU+M（4）

其中

K=AAT

15

aH0000

a21a220…0

aaa

A二313233…°

_aNlaN2aN3…aNN_

元素均可按列依次算出。

对第一列元素有：

%1=拓？%=尤1/%1(7=2,3,…,N)

算出1,2,…,j-1列各元素后，第/列的主对角元素为：

«=1

当，＜〃时.，主对角线以下各元素为

囱=端[局—-《MJ(i=j+l,…,N)

k=\

3、具有有理谱的正态随机序列的模拟

根据随机过程通过线性系统的理论，白噪声通过线性系统后，输出是正态的，且输出功

率谱只与系统的传递函数有关。利用这一性质，我们可以产生正态随机过程。

W(M)----K%(〃)------X(n)

图2.1白噪声通过线性系统

如上图所示，输入W(〃)为白噪声，假定功率谱密度为Gw(z)=l，通过离散线性系统

后，输出X(〃)是正态随机序列,

GX(Z)=GW(Z)"(Z)"(ZT)(5)

由于要求模拟的随机序列具有有理谱，则Gx(Z)可表示为：

Gx(z)=G：(z)GXz)(6)

其中，G；(z)表示有理谱部分，即所有的零极点在单位圆之内，G](z)表示非有理谱

部分，即所有零极点在单位圆之外。由此可以得到

〃(z)=G；(z)(7)

4、满足一定相关函数的平稳正态随机过程的模拟

16

当已知平稳随机过程的相关函数而要确定该随机过程的模拟算法。很显然，只要我们设

计一个合适的滤波器，使得该白噪声通过滤波器后，输出的功率谱满足上述相关函数的傅里

叶变换，就可以模拟得到该随机过程。

3、实验内容

1、产生两组相互独立的（0,1）均匀分布的随机数（随机数个数：500）；

2、按照实验原理中的方法产生一组均值为1,方差为1的正态分布的随机序列（序列

长度：500）；

3、画出功率谱密度为G,（0）=-----1-----的功率谱图（一个周期内），采用均匀采

1.25+COSG

样方法，采样点数为500；

4、模拟产生一个功率谱为GX<y）=-----------的正态随机序列；

1.25+cos。

5、画出4中产生的正态随机序列的功率谱图；

6、估计4中产生的随机序列的功率谱，与理论的功率谱G,（o）=-----!------进行比

1.25+cosg

较。

4、实验要求

（1）用MATLAB编写程序。

（2）写出详细试验报告。

17

实验3随机过程的特征估计

1、实验目的

了解随机过程特征估计的基本概念和方法，学会运用Matlab函数对随机过程进行特征

估计,并且通过实验了解不同估计方法所估计出来结果之间的差异。

2、实验原理

设随机序列X(n)、Y(n)为各态历经过程，样本分别为x(n)、y(n)(n=0』

1、均值的估计

1N7

mx=—(1)

2、方差的估计

方差估计有两种情况，如果均值mx已知,则

1N—1

比=万%("加，)2⑵

IV〃=0

如果均值未知，那么

1N-T

以=1ZE%)⑶

w—1〃=0

3、相关函数估计

]N

*、(〃?)=A,IIZ相+3"⑷

N-阿M

4、功率谱估计

功率谱的估计有几种方法,

(1)自相关法

先求相关函数的估计，

N-网"=0

然后对估计的相关函数做傅立叶变换,

N+I

4(3)=E%(〃。"加"(6)

(2)周期图法

先对序列x(n)做傅立叶变换,

1

N-\

X((o)=Zx(〃)e-"3(7)

?n=0

则功率谱估计为

2

GA((O)=^|X((O)|(8)

周期图法是一种非参数谱估计方法，另外还有一种修正的周期图方法，也叫Welch法，

MATLAB有周期图和Welch法的谱估计函数。

(3)现代谱估计技术

现代谱估计主要有参数谱估计和子空间谱估计。参数谱估计法是假定待估计功率谱的信

号是白噪声驱动线性系统的输出，常用的方法有基于最大埔估计的伯格算法和Yuler-Walk

自回归(AR)方法，这些方法是估计线性系统的参数，通常会得到比经典谱估计方法更好的估

计。子空间法也称为高分辨率谱估计或超分辨率谱估计，常用的方法有MUSIC法和特征矢

量法，这些方法是根据相关矩阵的特征分析或特征分解得到对信号的频率分量的估计，特别

适合于线谱(即正弦信号)的估计，是低信噪比环境下检测正弦信号的有效方法。

MATLAB有许多估计数字特征的统计函数：

(1)均值与方差

mean(A),返回序列的均值，序列用矢量A表示。

VAR(X),返回序列X的方差。

(2)互相关函数估计xcorr

用法：

c=xcorr(x,y)

c=xcorr(x)

c=xcorr(x,y,'option')

c=xcorr(x,'option')

xcorr(x,y)计算X与Y的互相关，矢量X表示序列x(n),矢量Y表示序列y(n)oxcorr(x)

计算X的自相关。option选项是：

]Af—|w|—I

,biased'：有偏估计Rx(/?:)=一£%〃+,一〃

N〃=o

人।N-帆-1

'unbiased':：无偏估计&(愤=-------;~~rVxn+/nxn

N-\m\制

,coeff：m=0的相关函数值归一化为1。

,none1：不作归,化处理。

(3)功率谱估计

MATLAB提供了许多功率铺估计的函数:

函数名用法功能

periodogram[Pxx,w]=periodogram(x)返回序列X用周期图法的功率谱估计

加参数window是采用修正的周期图

法，window指定窗的系数。

[Pxx,w]=periodogram(x,window)

pyulearPxx=pyulear(x,p)使用Yuler-walker方法的功率谱估

计，p表示AR模型的阶数

2

pburgPxx=pburg(x,p)使用Burg算法谱估计，p是估计功率

谱时指定的信号的自回归模型阶数

3、实验内容

1、产生一组均值为1,方差为4的正态分布的随机序列(1000个样本)，估计该序列的均值

与方差。

2、按如下模型产生一组随机序列：

x(n)=0.8x(n-l)+w(n)

其中w(n)为均值为1,方差为4的正态分布白噪声序列。估计过程的自相关函数与功率谱。

3、设信号为x(n)=sin(27djn)+2cos(27tf2n)+w(n),n=l,2,....,N,其中fi=0.05,f2=0.12,w(n)为

正态白噪声，试在N=356和1024点时，分别产生随机序列x(n),画出x(n)的波形并估计x(n)

的相关函数和功率谱。

4、实验要求

(1)用MATLAB编写程序；

(2)写出详细试验报告。

3

实验4数字图像的直方图均衡

1、实验目的

了解随机信号分析与处理的应用领域--图像处理的基本概念和方法，掌握直方图估计

概率密度的方法及其MATLAB实现。

2、实验原理

2.1、数字图像的基本概念

随机过程既可以是随时间变化的过程，也可以是随空间位置变化的过程。•幅图像可以

定义为一个二维函数g(x,y),其中x,y为平面坐标，在平面上一点(x,y)处，函数g的值称

为该点图像的强度或灰度级，当x、y以及g的值都是有限的离散数据时，我们称此图像为

数字图像。数字图像是由有限个元素组成，每个元素都有特定的位置和值，这些元素称为图

像元素或像素。平面上某一点的灰度可以看作为一个随机变量，因此，一幅图像可以看作为

随位置变化的随机序列。

2.2数字图像的直方图

什么是数字图像灰度级的直方图呢？简单的说，灰度级的直方图就是反映一副图像中灰

度级与出现这种灰度概率之间关系的图形。

设变量R代表图像中像素灰度级，R的取值范围为[0,L-l],L为总的灰度级数，具有

L个灰度级的数字图像直方图是一个离散函数：

h(rk)-nkk=O,l,…L(I)

其中〃是第k个灰度级，阪是在图像中具有灰度级人的像素数，/?(〃)=()代表黑色，

%(〃)=1代表白色。在实际中，我们通常对直方图的每个值除以图像中总的像素数n做归

一化处理，因此，归一化的直方图为

n.

/(〃)=,k=O,l,....L(2)

n

/(〃)可以看作为第％个灰度级出现概率的估计。很显然，对于任意的k,

0</(rJ<l⑶

以及

£/(〃)=1⑷

k=0

直方图提供了图像的统计信息，它可以用于图像压缩、图像增强等图像处理技术中。我

1

们对比图1四幅图像，分别表示低亮度和高亮度以及低对比度和高对比度的图片及其直方

图。

我们注意到，在低亮度图像中，直方图分量集中在灰度级的低端(喑端)，而高亮度图

像的直方图分量主要集中在灰度级的高端；低对比度图像的直方图分量集中在比较窄的灰度

级匕而高对比度图像的直方图分布比较平均。

低对比度图片

高对比度图片

图4.1不同亮度和对比度图片及其直方图对比

2.3直方图均衡

设随机变量R在区间［0,1］上取值，他的概率密度为fR(r),对R做变换定义一个新的变

量S,

S=T(R)(5)

其中变换函数T(.)是一个单调函数，且

0<?(/?)<1(6)

随机变量S的概率密度为

人(5)=%(「)乎(7)

as

一类重要的变换是

(8)

S=T(R)=£fR(r)dr=FR(R)

其中FMr)是随机变量R的分布函数。可以证明，随机变量S在［0,1］区间上服从均匀分布,

即

fs(5)=10<5<1(9)

由于数字图像的像素的灰度级是离散值，因此,我们需要用概率和求和取代概率密度和积分。

设R代表图像的灰度级，图像中的灰度级1^=取的概率为

几

P{R=〃}=%(〃)=1k=0,l,2,...,L-l(10)

n

其中n是图像的总的像素数，取是灰度级为w的像素数，做如下变换，

S=T(R)(11)

2

上式变换意味着对R的每一个灰度级都做变换，既

s*=T(〃)k=0,l,2,...,L-l(12)

对于离散形式，(8)式为

kk〃

s*=T(〃)=Z/R(5)=E°k=O,l,....,L-l(13)

)=0j=0n

通过上面的变换，被处理图像的像素灰变换成了Sk，由于S服从均匀分布，所以(13)称

为图像的直方图均衡。图2给出了原始图象和经过均衡后的图，原始图象灰度集中在100-200

之间，经过均衡后，灰度在0-256之间分布比较均匀。

原始图象

均衡化后图象

3、实验方法

直方图均衡是对一幅图像按(13)式做变换，使处理后的图像直方图近似为均匀分布变

其计算步骤如F：

(1)列出原是图像的灰度级/，^=0,1,­••,£-!,其中L是灰度级的个数。

(2)统计各灰度级的像素数目勺，_/=0,1,…，♦一1。

(3)计算原始图像直方图各灰度级的频数/八5)=)，j=0,l,…,L—1,其中〃为原

始图像总的像素数目。

(4)对原始图像每个像素的灰度级按(13)式进行变换，得到新的图像。

(5)按(1)-(3)步骤计算新图像的直方图。

3

实验5典型时间序列模型分析

1、实验目的

熟悉三种典型的时间序列模型：AR模型，MA模型与ARMA模型，学会运用Matlab

工具对对上述三种模型进行统计特性分析，通过对2阶模型的仿真分析，探讨几种模型的适

用范围，并且通过实验分析理论分析与实验结果之间的差异。

2、实验原理

AR模型分析

设有AR(2)模型,

X(n)=-0.3X(n-l)-0.5X(n-2)+W(n)

W(n)是零均值正态白噪声，方差为4。

(1)用MATLAB模拟产生X(n)的500观测点的样本函数，并绘出波形

(2)用产生的500个观测点估计X(n)的均值和方差

(3)画出理论的功率谱

(4)估计X(n)的相关函数和功率谱

【分析】给定二阶的AR过程，可以用递推公式得出最终的输出序列。或者按照一个白噪声

通过线性系统的方式得到，这个系统的传递函数为：

“(Z)=--------:-------

1+0.3—+0*2

这是一个全极点的滤波器，具有无限长的冲激响应。

对于功率谱，可以这样得到,

2

]

巴⑼=5：

-12

1+a.z+a2z~z=exp（j（o）

可以看出，巴(0)完全由两个极点位置决定。

对于AR模型的自相关函数，有下面的公式:

1

G(0)rjl)…々(p)一'1_-2

G⑴G(0)心(。-1)0

%=

G(P)G(P-1)•••G(。)A.0

这称为Yule-Walker方程，当相关长度大于p时，由递推式求出：

rv(r)+tz1(r-l)+---+ap(r-/?)=O

这样，就可以求出理论的AR模型的自相关序列。

1.产生样本函数，并画出波形

题目中的AR过程相当于一个零均值正态白噪声通过线性系统后的输出，可以按照上面

的方法进行描述。

clearalI;

b=[1];a=[10.30.5];%由描述的差分方程，得到系统传递函数

h=impz(b,a,20);%得到系统的单位冲激函数，在20点处已经可以认为值是0

randn('state',0);

w=normrnd(0,2,1,500);%产生题设的白噪声随机序列，标准差为2

x=filter(b,a,w);%通过线形系统，得到输出就是题目中要求的2阶AR过程

plot(x,'r');ylabel('x(n),);ti11e('产生的AR随机序列,);grid

得到的输出序列波形为：

2.估计均值和方差

可以首先计算出理论输出的均值和方差，得到m,=0,对于方差可以先求出理论自相

关输出，然后取零点的值。

rx[m)=*h(-m)*rw(〃?)

并且，(m)=4次加)，带入有

rx(ni)=4[/?(zn)*/z(-m)]

可以采用上面介绍的方法，对式中的卷积进行计算。计算出的卷积输出图形为:

2

1.5

0.5-

0t-a•-•«•♦-.，*.-————————,———»-»-<-*-•-

XI.1

•0.5

0510152025303540

n

在最大值处就是输出的功率，也就是方差，为

<7：=rt(0)=5.6

对实际数据进行估计，均值为mean(x)=-0.0703,而方差为var(x)=5.2795,两者和理论

值吻合的比较好。

3.画出理论的功率谱密度曲线

理论的功率谱为，

以(。)=匕(司〃(/呻=4|〃(/)『

用下面的语句产生：

delta=2*pi/1000;w_min=-pi;w_max=pi;Fs=1000;

w=w_min:delta:w_max;%得到数字域上的频率取样点，范围是[-pi,pi]

Gx=4*(abs(1./(1+0.3*exp(-i*w)+0.5*exp(-2*i*w))).A2);%计算出理论值

Gx=Gx/max(Gx);%归一化处理

f=w*Fs/(2*pi);%转化到模拟域上的频率

plottf.Gx/b'l.gridon;

得到的图形为：

可以看出，这个系统是带通系统。

4.估计自相关函数和功率谱密度

用实际数据估计自相关函数和功率谱的方法前面已经讨论过，在这里仅给出最后的仿真

3

图形。

%计算理论和实际的自相关函数序列

Mlag=20;%定义最大自相关长度

Rx=xcorr(x,Mlag,'coeff1);

m=-Mlag:Mlag;

stem(m,Rx,*r.');

最终的值为

可以看出，它和卜.面的理论输出值吻合程度很好。

实际的功率谱密度可以用类似于上面的方法进行估计，

window=hamming(20);%采用hanmming窗，长度为20

noverIap=10;%重叠的点数

Nfft=512;%做FFT的点数

Fs=1000;%采样频率，为1000Hz

[Px,f]=pwelch(x,window,noverlap,Nfft,Fs,1onesided');%估计功率谱密度

f=[-fHplr(f)f(2:end)];%构造一个对称的频率，范围是[•Fs/2,Fs/2]

Py=[-fIiplr(Py)Py(2:end)];%对称的功率谱

plot(f,10*log10(Py),'b');

估计出来的功率谱密度为，

-12

“WcO-400-300―-200-1000WO―200―荷―400500

f(Hz)

将两幅图画在•起，可以看到拟合的情况比较好：

4

0

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)的相关函数和功率谱

【分析】给定(2,2)的ARMA过程，也可以用递推公式得出最终的输出序列。或者按照一

个白噪声通过线性系统的方式得到，这个系统的传递函数为:

l+0.5：0.2z：

1+0.3/—0.21

对于功率谱，可以这样得到,

对于ARMA过程，当模型的所有极点均落在单位圆内时，才是一个渐进平稳的随机过程。

这个过程的自相关函数不能简单地写成Yule-Walker方程形式，它于模型的参数具有高度的

非线性关系。

1.产生样本函数，并画出波形

题目中的ARMA过程相当于一个零均值正态白噪声通过线性系统后的输出，可以按照

上面的方法进行描述。

5

clearalI;

b=[10,5-0.2];a=[10.3-0.2];%由描述的差分方程，得到系统传递函数

h=impz(b,a,10);%得到系统的单位冲激函数，在10点处已经可以认为值是0

randn「state',0);

w=normrnd(0,2,1,500);%产生题设的白噪声随机序列，标准差为2

x=filter(b,a,w);%通过线形系统，得到输出就是题目中要求的(2,2)阶ARMA过程

plot(x,*r1);

得到的输出序列波形为：

2.估计均值和方差

可以首先计算出理论输出的均值和方差，得到根,=0,对于方差可以先求出理论自相

关输出，然后取零点的值。

rx(〃?)=/i(m)*h(-ni)*rK(〃?)

并且，加)=43(加)，带入有

rv(/?7)=4[/2(m)*/z(-?n)]

可以采用上面介绍的方法，对式中的卷积进行计算。计算出的卷积输出图形为:

6

在最大值处就是输出的功率，也就是方差，为

cr；=rv(0)=4.2

对实际数据进行估计，均值为mean(x)=-0.0547,而方差为var(x)=3.8,两者和理论值吻

合的比较好。

3.画出理论的功率谱密度曲线

理论的功率谱为，

.(0)=.3)忸(一).=4忸(/。旷

用下面的语句产生：

delta=2*pi/1000;w_min=-pi;w_max=pi;Fs=1000;

w=w_min:delta:w_max;%得到数字域上的频率取样点,范围是[•pi,pi]

NS=1+0,5*exp(-i*w)-0.2*exp(-2*i*w);%分子

DS=1+0,3*exp(-i*w)-0.2*exp(-2*i*w);%分母

Gx=4*(abs(NS./DS).A2);%计算出理论值

Gx=Gx/max(Gx);f=w*Fs/(2*pi);%转化到模拟域上的频率

plot(f,Gx,'b'),gridon;

0.9

0.8

0.7

0.6

0.5

0.4

0.3

-200-100

f(HZ)

4.估计相关函数和功率谱密度曲线

用实际数据估计自相关函数和功率谱的方法前面已经讨论过，在这里仅给出仿真图形。

%计算理论和实际的自相关函数序列

Mlag=20;%定义最大自相关长度

Rx=xcorr(x,Mlag,'coeff');

m=-Mlag:Mlag;

stem(m,Rx,'r.');

最终的值为

7

1.2

0.8

0.6

0.4

0.2

0-•-

-20-15-10-505101520

n

实际的功率谱密度可以用类似于上面的方法进行估计，

window=hamming(20);%采用hanmming窗，长度为20

noverIap=10;%重叠的点数

Nfft=512;%做FFT的点数

Fs=1000;%采样频率，为1000Hz

[Px,f]=pwelch(x,window,noverlap,Nfft,Fs,'onesided');%估计功率谱密度

f=(-Hiplr(f)f(2:end)];%构造一个对称的频率，范围是[•Fs/2,Fs/2]

Py=[fIiplr(Py)Py(2:end)];%对称的功率谱

plot(f,10*1oglO(Py),1b');

估计出来的功率谱密度为，

-400-300-200-1000100200300400500

f(Hz)

把两幅图画在一起，可以得到下面的图形，可以看出两者的吻合度比较高。

8

3、实验内容

1、熟悉实验原理，将实验原理上的程序应用matlab工具实现；

2、设有MA(2)模型，x(n)=W(n)-0.3W(n-l)+0.2W(n-2)

W(n)是零均值正态白噪声，方差为4。

(1)用MATLAB模拟产生X(n)的500观测点的样本函数，并绘出波形

(2)用产生的500个观测点估计X(n)的均值和方差

(3)画出理论的功率谱

(4)估计X(n)的相关函数和功率谱

4、实验要求

(1)用MATLAB编写程序。

(2)写出详细试验报告(要有自己对实验结果的结论)。

9

实验6随机过程通过线性系统分析

1、实验目的

了解随机信号通过系统的时域分析方法和频域分析方法，了解随机过程通过系统后均

值、相关函数和频谱与输入随机过程的关系，从而掌握随机过程通过系统的分析方法。

2、实验原理

随机信号通过线性系统分析的中心问题是：给定系统的输入函数（或统计特性：均值和

自相关函数）和线性系统的特性，求输出函数。如下图所示，T为线性变换，信号X（f）为

系统输入，y«）为系统的输出，也是随机信号。

X。）-------->H---------►y（z）

图6.1随机信号通过系统示意图

则有:

〃[x（f）]=y（，）（1）

均值和相关函数估计

对.（i）式两边取数学期望，

E{H[X（t）]}^E[Y（t）]（2）

即

H[mx（0]=rnY（r）（3）

相关函数可通过下面方法求出：

RXY=H[Rx）1（4）

Ry（小，2）="八[呢丫G/2）]（5）

我们比较常见的时域分析方法是冲击响应法.

•冲击响应法

均值为：

mY=£mxh（t）dt=mxh（t）dt=（6）

其中H（0）为系统传递函数在G=0时的值。

1

自相关函数为:

勺(4,/2)=/!(/,)0h(t2)0Rx(r,,t2)(7)

如果X(t)是平稳随机信号，贝ij：

RY(r)=h(-r)®/z(r)0Rx(r)

其中//«)为系统的传递函数。

功率谱估计

Gxy3)="*(MGx(M

Gy(co)=H(a))GXY(ey)

GYG)=H*(o)H(MGx(iy)=\HG)『Gx⑼(11)

simulink模拟方法：

图2.2.1是用MATLAB的sumulink模拟白噪声通过图222的RC电路，用示波器观

察输入和输出的波形，改变RC的值，使电路时间常数改变，观察输出波形的变化。

图2.2.2RC电路

【分析】本例演示白噪声通过线性系统后的变化。可以从时域和频域两方面进行分析。

对于上述低通RC滤波器，可以用下面的微分方程描述，

用传递函数描述，令。=—L,有"(5)=—竺

RCs+a

其单位冲激响应为，h(t)=—e-,,RCu(t).

RC

在Similink里，有时域连续系统的传递函数模块，如图所示：

在左图中，三个选项分别为“分母系数”、“分子系数”和“误差限制”，根据上面的分析,

可知

s+a

2

式中a=—,即为滤波器的时间常数，在仿真中，可以根据实际情况自由选择。比如可以

RC

选取a=l,这样便有

"(s)=」7

5+1

所以，上面显示的系数便是正确的，不用修改。

输入的白噪声可以有两种方式产生，•种是直接调用Simulink里面的“带限白噪声”

或者“正态随机数发生器”，另外一种是用M文件编写，通过工作区进行数据传递。简单起

见，采用第一种方法。先对二者进行简单介绍，

上述的带限白噪声模块中Seed参数代表了随机数发生器初始化的状态，可以自行设定。

具体的含义参议参考Matlab中随机数发生器的原理说明文件。

3

上述的随机数发生器中的"InitialSeed”是同样的含义。通过设定均值和方差，可以产生

服从正态分布（高斯分布）的随机数，并且根据采样时间进行零阶保持。

下面是实验的原理框图：

Scope2

原始白噪声信号

4

a=l系统输出信号

a=0.1系统输出信号

同时，我们可以看出，经过不同的系统之后，输入白噪声的变化有快又慢。直观地，可

以看到，如果系统的时间常数大，对应输出信号的变化更加强烈，反之较为平稳。可以进行

简单的分析，分析如下：

该系统的传递函数为：”（。）=——

a+jco

输出的功率谱密度：Gy（M=g|"⑹2吟//

输出的自相关函数为，&«）=等咖

1（12a

可以求出等效噪声能带颂=——[——彳da）=—

e2开小4+疗4

所以，随即过程的相关时间为：r0=-即有，。

比a4甑

相关时间与系统的噪声等效能通能带成反比。白噪声通过该系统，如果af8,会导

致相关时间趋于零，此时RC电路成为了全通系统，输出仍然为白噪声。当a->0时，相关

时间很大，白噪声只有低频部分能通过，输出噪声的变化就慢了。

下面对随机过程的统计特性进行分析，讨论随机过程通过线性系统后均值、方差、自相

关函数、功率谱密度的变化，并对平稳性进行简单讨论。

均值的变化

开始假定随机过程的均值为0,方差为1,这可以通过对实验数据进行验证。在仿真时,

可以选择示波器以选定的矩阵方式来向工作空间保存数据，以便于进行进一步的分析。

我们可以用下面两句指令获取存储的数据和时间，

x=ScopeData.signals.values;

5

t=ScopeData.time;

通过计算，得到输入信号的均值为0.0219,输出信号均值为0.0425。根据公式，可知

mr=jmx//(r)JT=mxH(0)

而〃(o)=—!—=i,实验值和理论值之间存在着一定的误差，但是比较小。

S+J。

自相关函数的变化

按照公式，可知

RXY(T)=h(-T)*Rx(T)

Ry(r)=/?(r)*RXY(r)=A(r)*h(-T)*Rx(r)

功率谱的关系

由相关理论，得知Gy3)=Gx(M|〃3)「

3、实验内容

1.模拟产生如下图所示的包含信号S(f)的随机过程X"),即XQ)=S(f)+NQ)(其

fS⑺

中，N(f)为零均值方差为1的高斯白噪声，假设采样率为1000Hz,采样长度为10

1

秒(0-10))o

0510

t(s)

2.计算1中产生的随机过程的自相关函数，通过比较该随机序列的相关函数和5(/)的

相关函数之间的差别，确认模拟产生的随机过程的有效性。

3.让模拟产生的随机过程通过-RC系统，关系图如下：

R

X(f)。=Y(t)

其中RC=1,对系统进行建模，求出传递函数及其冲击响应，并对冲击响应离散化,

6

采样率为1000Hz,采样长度为10秒(取前10秒)。

4.用冲击响应法求输出的自相关函数。

5.用频谱法计算出输出信号的功率谱，通过傅立叶反变换求出相关函数，和4进行比

较，有何异同？

6.改变1种噪声的方差，重复1-5,判断异同点。

7.用simulink框图方法实现上面的流程。

4、实验要求

(1)用MATLAB编写程序。

(2)画出相关波形，并便于对比。

(3)写出详细试验报告。

7

实验7窄带随机过程的产生

1、实验目的

熟悉窄带随机过程的定义，了解窄带随机过程产生的原理与方法，最后估计实验产生的

窄带随机过程的功率谱。

2、实验原理

本实验模拟产生一个窄带随机过程。首先产生两个相互独立的随机过程A,⑴和A»(t),

并将用两个正交载波cos2W°f和sin2W°f进行调制，如图5.25所示，然后进行抽样得到窄

带过程的抽样。

4，(i)cos

-4(j)sin

带通随机过程的产生

低通过程AKt)和A《)的抽样是通过将两个独立的白噪声序列通过两个相同的低通滤波

器产生的，分别得到&⑴和Ae)的抽样Ac(n)和As(n),低通滤波器的传递函数选择为

I

"(Z)=

1-0.9Z-1

然后将AJn)乘以cos2矶仃，As(n)乘以sin2班"，再通过相减得到窄带随机过程的抽样

值，其中T为抽样间隔，假定T=l,/o=1000/7：,估计窄带过程的相关函数和功率谱并绘

出图形。

3、实验内容

(1)按照实验原理