随机工程matlab实验.docx

随机工程matlab实验.docx

《随机工程matlab实验.docx》由会员分享，可在线阅读，更多相关《随机工程matlab实验.docx（19页珍藏版）》请在冰点文库上搜索。

随机工程matlab实验

自动化学院

随机过程课程设计

学号：

S3*******1

专业：

导航、制导与控制

学生姓名：

任课教师：

赵希人教授

2008年12月

第一题、

源程序：

x=rand（1,2000）

EX=mean（x）%均质检验

DX=var（x）%方差检验

subplot（2,1,1）,hist（x,10）;

y=linspace（-10,10,21）;

form=-10:

10

mAbs=abs（m）;

s=0;

forn=1:

2000-mAbs

s=s+（x（n+mAbs）-EX）*（x（n）-EX）;

end

y（m+11）=（1/（2000-mAbs））*s;

end

x2=[-15:

30/20:

15];

subplot（2,1,2）,plot（x2,y）;

1.打印前50个数：

Columns1through12

0.33290.51890.54350.90380.24510.13240.09770.69750.37070.19850.57700.4111

Columns13through24

0.15700.71280.34450.64340.47110.31070.96100.22360.06070.26570.38610.3220

Columns25through36

0.39330.30220.71590.43120.59410.90340.13900.11190.81270.17620.51750.0592

Columns37through48

0.45350.71520.16390.16610.24400.22960.97260.54980.79010.35000.28990.8619

Columns49through50

0.99370.3002

2.分布检验：

3.均值检验：

EX=

0.5038

4.方差检验：

DX=

0.0827

5.计算相关函数：

第二题、

源程序：

clearall;

forn=1:

2000

xt（n）=normrnd（0,1）;%产生2000个N（0,1）分布的独立序列

end

plot（xt）,title（'2000个N（0,1）分布的独立序列'）;

figure

fori=1:

5

forj=1:

10

sc（j,i）=xt（（i-1）*5+j）;

end;

end;

disp（[sc]）;

EX=mean（xt）%求平均数并输出

DX=cov（xt）%求方差并输出

subplot（2,1,1）;

p=hist（xt,20）;%将产生的2000个随机数分为20组

p=p/100;t=-2.85:

0.3:

2.85;%求概率密度

bar（t,p,1）,title（'N（0,1）分布的独立序列的直方图'）;

xlabel（'x'）;ylabel（'f（x）'）;

[tx,i]=xcov（xt,10）;%τ取-10到10

Tx=tx/2000;%求自相关函数Γx（τ）

subplot（2,1,2）

plot（i,Tx,'.-'）,title（'自相关函数Γx（τ）分布图'）;

xlabel（'τ'）;ylabel（'Γx（τ）'）;

1.打印前50个数：

0.35990.75621.00201.37470.0046

0.63071.1306-0.45530.0294-0.5836

-0.6737-0.8743-0.0972-1.62000.5674

0.7104-0.80680.7932-0.39710.4200

0.6455-0.0149-1.1756-0.0226-0.3667

0.75621.00201.37470.00460.5132

1.1306-0.45530.0294-0.58361.6777

-0.8743-0.0972-1.62000.5674-0.5122

-0.80680.7932-0.39710.42000.4997

-0.0149-1.1756-0.0226-0.36670.6410

程序运行产生的2000个N（0，1）分布的独立序列：

2.分布检验：

3.均值检验：

EX=

0.0069

4.方差检验：

DX=

0.9903

5.计算相关函数：

第三题、

源程序：

y=normrnd（0,1,1,2000）;

x=zeros（1,2000）;

fork=2:

2000

x（k）=y（k）+4*y（k-1）;

end;

Ex=mean（x）%求EX（k）

Ex2=mean（x.^2）%EX.^2（k）

Dx=var（x）%DX（k）

Bm=zeros（1,11）;

fori=0:

10

xy=zeros（1,2000-i）;

forj=1:

（2000-i）

xy（j）=（x（j+i）-Ex）*（x（j）-Ex）;

end;

Bm（i+1）=sum（xy）/（2000-i）;%Bx（m）

end;

a=Bm（:

11:

-1:

1）;

a=[a（1:

10）,Bm];

Bm=a

q=-10:

10;

plot（q,Bm）

1.求均值EX：

Ex=

0.1767

2.求二阶矩EX2：

EX2=

16.6688

3.求方差：

Dx=

16.6459

4.BX（m）：

Columns1through12

0.5104-0.1167-0.41100.17650.44230.12370.2300-0.3652-0.44334.065216.63764.0652

Columns13through21

-0.4433-0.36520.23000.12370.44230.1765-0.4110-0.11670.5104

第四题、

源程序：

y=normrnd（0,1,1,2000）;

x=zeros（1,2000）;

fork=2:

2000

x（k）=y（k）-0.707*x（k-1）;

end;

M=zeros（1,901）;

fork=101:

1001

x（k）=x（k）+x（k-1）;

M（k-100）=x（k-1）;

end;

Ex=x（k）/900%EX（k）

Ex2=sum（M.^2）/900%EX.^2（k）

Dx=Ex2-（Ex）.^2%DX（k）

Bm=zeros（1,11）;

fori=0:

10

xy=zeros（1,1900-i）;

forj=100:

（2000-i）

xy（j-99）=（x（j+i）-Ex）*（x（j）-Ex）;

end;

Bm（i+1）=sum（xy）/1900%Bx（m）

end;

a=Bm（:

11:

-1:

1）;

a=[a（1:

10）,Bm];

q=-10:

10;

plot（q,a）%Bx（m）二维分布曲线

1.求EX：

Ex=

-0.0025

2.求Ex2：

Ex2=

39.8177

3.求DX：

Dx=

39.8177

4.求Bx（m）：

Bm=

19.859317.743619.024318.024518.457418.015918.086317.838117.772517.591917.4769

第五题、

源程序：

f='sin（x）';%画出符号函数

subplot（3,2,1）;

ezplot（f）

subplot（3,2,2）;

ezplot（f）

n=-20:

20;

y=sin（n*pi/2）;

subplot（3,2,3）;

k=-10:

10;

plot（k,sin（k*pi/2）,'--rs','MarkerEdgeColor','k','MarkerFaceColor','g','MarkerSize',5）

title（'采样信号（±10点）'）;

subplot（3,2,4）;

k=-20:

20;

plot（k,sin（k*pi/2）,'--rs','MarkerEdgeColor','k','MarkerFaceColor','g','MarkerSize',5）

title（'采样信号（±20点）'）;

D=0.05;

z1=1;

fort1=-5*pi:

D:

5*pi;

s1=0;

form1=-10:

10

s1=s1+y（m1+11）*sinc（（1/pi）*（t1-pi*m1/2））;

%sinc函数内插恢复，重建信号过程

end

fa（z1）=s1;z1=z1+1;

end

subplot（3,2,5）

xlab1=linspace（-5*pi,5*pi,length（fa））;

plot（xlab1,fa）;title（'内插恢复信号（±10点）'）;

z2=1;

fort1=-10*pi:

D:

10*pi;

s2=0;

form2=-20:

20

s2=s2+y（m2+21）*sinc（（1/pi）*（t1-pi*m2/2））;

%sinc函数内插恢复，重建信号过程

end

fb（z2）=s2;z2=z2+1;

end

subplot（3,2,6）

xlab2=linspace（-10*pi,10*pi,length（fb））;

plot（xlab2,fb）;title（'内插恢复信号（±20点）'）;

运行结果：

说明：

图中第一行两图为sin函数的符号函数示意图，区间-2π～2π。

图中第二行两图为sin函数的采样函数示意图，左图区间-5π～5π（±10点），右图区间-10π～10π（±20点）。

图中第三行两图为sin函数的恢复函数示意图，右图比左图更加接近原函数，但是由于是有限采样点恢复，并且，内插sinc函数本身是非因果的连续信号，物理上不可实现，故只能用大样本离散值模拟，程序中取步长D=0.05。

所以，可以看到两图在函数边缘有较大过冲。

这是由于sinc函数被截断而引起的。

第六题、

源程序：

A=zeros（40,21）;

k=zeros（1,20）;

fori=1:

21

A（1,i）=2*i-1;

end

forj=1:

20

ifmod（j,2）==1;

A（2,j）=A（1,j+1）;

else

A（2,j）=0;

end;

end;

k

（1）=A（1,1）/A（2,1）;

fori=3:

40

ifmod（i,2）==1;

forj=（i+1）/2:

2:

20;

A（i,j）=A（i-2,j）;

ifj~=20;

A（i,j+1）=A（i-2,j+1）-A（i-1,j+1）*k（（i-1）/2）;

else

A（i,j+1）=41;

end;

end

else

forj=i/2:

2:

20

A（i,j）=A（i-1,j+1）;

end;

k（i/2）=A（i-1,i/2）/A（i,i/2）;

end;

end;

disp（A）

disp（k）

1.奥斯特姆表：

A=

Columns1through12

1.00003.00005.00007.00009.000011.000013.000015.000017.000019.000021.000023.0000

3.000007.0000011.0000015.0000019.0000023.00000

03.00002.66677.00005.333311.00008.000015.000010.666719.000013.333323.0000

02.666705.333308.0000010.6667013.3333016.0000

002.66671.00005.33332.00008.00003.000010.66674.000013.33335.0000

001.000002.000003.000004.000005.00000

0001.000002.00000.00003.00000.00004.000005.0000

000000.000000.00000000.0000

00000-Inf0.0000-Inf0.0000NaN0-Inf

0000-Inf0-Inf0NaN0-Inf0

00000-InfNaN-InfNaNNaNNaN-Inf

Columns13through21

25.000027.000029.000031.000033.000035.000037.000039.000041.0000

27.0000031.0000035.0000039.000000

16.000027.000018.666731.000021.333335.000024.000039.000041.0000

018.6667021.3333024.0000041.00000

16.00006.000018.66677.000021.33338.000024.0000-7.12500

6.000007.000008.00000-7.125000

0.00006.00000.00007.000008.000043.0000-7.125041.0000

00.000000043.0000041.00000

0.0000-Inf0.0000NaN0-Inf43.0000-Inf0

-Inf0NaN0-Inf0-Inf00

NaN-InfNaNNaNNaN-InfNaN-Inf41.0000

Ak=0.3333,1.1250,2.6667,NaN,NaN,NaN,NaN,NaN,NaN,NaN,NaN,NaN

B=

Columns1through12

000000000000

3.000007.0000011.0000015.0000019.0000023.00000

000000000000

02.666705.333308.0000010.6667013.3333016.0000

000000000000

001.000002.000003.000004.000005.00000

000000000000

000000.000000.000000.000000.0000

0000NaNNaNNaNNaNNaNNaNNaNNaN

0000-Inf0-Inf0NaN0-Inf0

00000NaNNaNNaNNaNNaNNaNNaN

Columns13through21

00000001.00000

27.0000031.0000035.0000039.000000

00000001.00000

018.6667021.3333024.0000041.00000

00000001.00000

6.000007.000008.00000-7.125000

00000001.00000

00.000000043.0000041.00000

NaNNaNNaNNaNNaNNaNNaNNaNNaN

-Inf0NaN0-Inf0-Inf00

NaNNaNNaNNaNNaNNaNNaNNaNNaN

Bk=0,0,0,NaN,NaNNaN,NaN,NaN,NaN,NaN,NaN,NaN,NaN,NaN,NaN,NaN,NaN

2.由奥斯特姆表可以看出系统不稳定。

3.

第七题、

源程序：

n=7;

N=13;

x=[10.50.58-0.01-0.01210.000050.00006];

y=[10.5500000];

%x=[10.70.5-0.3];

%y=[10.30.20.1];

A=zeros（N,n）;%分配内存空间

B=zeros（N,n）;

a=zeros（n-1,1）;

b=a;

fori=1:

n,%赋表的初始值，前两行及a

（1）,b

（1）;

A（1,i）=x（i）;

B（1,i）=y（i）;

end;

fori=1:

n,

A（2,i）=A（1,n-i+1）;

end;

fori=1:

n,

B（2,i）=A（2,i）;

end;

a

（1）=A（1,n）/A（2,n）;

b

（1）=B（1,n）/B（2,n）;

fori=2:

n,

forj=1:

n-i+1,

A（2*i-1,j）=A（2*i-3,j）-A（2*i-2,j）*a（i-1）;

B（2*i-1,j）=B（2*i-3,j）-B（2*i-2,j）*b（i-1）;

end;

ifi==n,

break;

end;

forj=1:

n-i+1,

A（2*i,j）=A（2*i-1,n-i+2-j）;

end;

forj=1:

n,

B（2*i,j）=A（2*i,j）;

end;

a（i）=A（2*i-1,n-i+1）/A（2*i,n-i+1）;

b（i）=B（2*i-1,n-i+1）/B（2*i,n-i+1）;

end;

A

B

1.奥斯特姆表：

A=

1.00000.50000.5800-0.0100-0.01210.00010.0001

0.00010.0001-0.0121-0.01000.58000.50001.0000

1.00000.50000.5800-0.0100-0.01210.00000

0.0000-0.0121-0.01000.58000.50001.00000

1.00000.50000.5800-0.0100-0.012100

-0.0121-0.01000.58000.50001.000000

0.99990.49990.5870-0.0039000

-0.00390.58700.49990.9999000

0.99980.50220.58900000

0.58900.50220.99980000

0.65280.206300000

0.20630.652800000

0.5876000000

B=

1.00000.550000000

0.00010.0001-0.0121-0.01000.58000.50001.0000

1.00000.550000000

0.0000-0.0121-0.01000.58000.50001.00000

1.00000.550000000

-0.0121-0.01000.58000.50001.000000

1.00000.550000000

-0.00390.58700.49990.9999000

1.00000.550000000

0.58900.50220.99980000

1.00000.550000000

0.20630.652800000

0.8262000000

2.由奥斯特姆表可看出系统稳定。

3.