白噪声.docx
《白噪声.docx》由会员分享,可在线阅读,更多相关《白噪声.docx(26页珍藏版)》请在冰豆网上搜索。
白噪声
白噪声
如果一个零均值、平稳随机过程的谱密度为常数,我们称之为白噪声(由白色光联想而得)。
白噪声有以下特点:
1
2
,频谱宽度无限。
∞,τ=0
0,τ≠0
3
,其中,
为Dirac函数,即
=
且
4无记忆性,即t时刻的数值与t时刻以前的过去值无关,也不影响t时刻以后的将来值。
从另一意义上说,即不同时刻的随机信号互不相关。
白噪声的用途:
1作为系统输入时,有
,τ=0,1,2,…,即为系统的单位脉冲响应。
2作为被辨识系统输入时,可以激发系统的所有模态,可对系统充分激励;
3作为被辨识系统输入时,可防止数据病态,保证辨识精度。
4在辨识过程中,以输出估计误差是否具有白色性来判断辨识方法的优劣,也可用来判断模型的结构和参数是否合适。
5产生有色噪声。
白噪声的产生方法:
1(0,1)均匀分布白噪声:
{
,i=1,2,3…},初值可取为:
2正态分布白噪声
:
,
其中ξ为服从(0,1)均匀分布的白噪声。
例2.1用乘同余法产生随机数
选A=6,k=8,M=28=256,递推100次.用MATLAB编程
①编程如下:
A=6;N=100;%初始化;
x0=1;M=255;
fork=1:
N%乘同余法递推100次;
x2=A*x0;%x2和x0分别表示xi和xi-1;
x1=mod(x2,M);%将x2存储器的数除以M,取余数放x1(xi)中;
v1=x1/256;%将x1存储器的数除以256得到小于1的随机数放v1中;
v(:
k)=v1;%将v1中的数(
)存放在矩阵存储器v的第k列中,v(:
k)
%表示行不变、列随递推循环次数变化;
x0=x1;%xi-1=xi;
v0=v1;
end%递推100次结束;
v2=v%该语句末无‘;’,实现矩阵存储器v中随机数放在v2中,%且可直接显示在MATLAB的window中;
k1=k;
%grapher%以下是绘图程序;
k=1:
k1;
plot(k,v,k,v,'r');
xlabel('k'),ylabel('v');title('(0-1)均匀分布的随机序列')
②程序运行结果如图所示
③产生的(0-1)均匀分布的随机序列
在程序运行结束后,产生的(0,1)均匀分布的随机序列,直接从MATLAB的window界面中copy出来如下(v2中每行存6个随机数):
v2=
0.02340.14060.84380.08200.49220.9609
0.78520.72660.37500.25780.55080.3164
0.90230.43360.60940.66800.02340.1406
0.84380.08200.49220.96090.78520.7266
0.37500.25780.55080.31640.90230.4336
0.60940.66800.02340.14060.84380.0820
0.49220.96090.78520.72660.37500.2578
0.55080.31640.90230.43360.60940.6680
0.02340.14060.84380.08200.49220.9609
0.78520.72660.37500.25780.55080.3164
0.90230.43360.60940.66800.02340.1406
0.84380.08200.49220.96090.78520.7266
0.37500.25780.55080.31640.90230.4336
0.60940.66800.02340.14060.84380.0820
0.49220.96090.78520.72660.37500.2578
0.55080.31640.90230.43360.60940.6680
0.02340.14060.84380.0820
。
1编程如下:
A=6;x0=1;M=255;f=2;N=100;%初始化;
x0=1;M=255;
fork=1:
N%乘同余法递推100次;
x2=A*x0;%分别用x2和x0表示xi+1和xi-1;
x1=mod(x2,M);%取x2存储器的数除以M的余数放x1(xi)中;
v1=x1/256;%将x1存储器中的数除以256得到小于1的随机数放v1中;
v(:
k)=(v1-0.5)*f;%将v1中的数(
)减去0.5再乘以存储器f中的系数,存放在矩阵存储器v的第k列中,v(:
k)表示行不变、列随递推循环次数变化;
x0=x1;%xi-1=xi;
v0=v1;
end%递推100次结束;
v2=v%该语句后无‘;’,实现矩阵存储器v中随机数放在v2中,且可直接显示在MATLAB的window中;
k1=k;
%grapher%以下是绘图程序;
k=1:
k1;
plot(k,v,k,v,'r');
xlabel('k'),ylabel('v');title('(-1,+1)均匀分布的白噪声')
②程序运行结果如图2.6所示。
图2.6采用MATLAB产生的(-1,+1)均匀分布的白噪声序列
③产生的(-1,1)均匀分布的白噪声序列
在程序运行结束后,产生的(-1,1)均匀分布的白噪声序列,直接从MATLAB的window界面中copy出来如下(v2中每行存6个随机数):
v2=
-0.9531-0.71880.6875-0.8359-0.01560.9219
0.57030.4531-0.2500-0.48440.1016-0.3672
0.8047-0.13280.21880.3359-0.9531-0.7188
0.6875-0.8359-0.01560.92190.57030.4531
-0.2500-0.48440.1016-0.36720.8047-0.1328
0.21880.3359-0.9531-0.71880.6875-0.8359
-0.01560.92190.57030.4531-0.2500-0.4844
0.1016-0.36720.8047-0.13280.21880.3359
-0.9531-0.71880.6875-0.8359-0.01560.9219
0.57030.4531-0.2500-0.48440.1016-0.3672
0.8047-0.13280.21880.3359-0.9531-0.7188
0.6875-0.8359-0.01560.92190.57030.4531
-0.2500-0.48440.1016-0.36720.8047-0.1328
0.21880.3359-0.9531-0.71880.6875-0.8359
-0.01560.92190.57030.4531-0.2500-0.4844
0.1016-0.36720.8047-0.13280.21880.3359
-0.9531-0.71880.6875-0.8359
例2.3用移位寄存器产生M序列的MATLAB软件实现
①编程如下:
X1=1;X2=0;X3=1;X4=0;%移位寄存器输入Xi初T态(0101),Yi为移位寄存器各级输出
m=60;%置M序列总长度
fori=1:
m%1#
Y4=X4;Y3=X3;Y2=X2;Y1=X1;
X4=Y3;X3=Y2;X2=Y1;
X1=xor(Y3,Y4);%异或运算
ifY4==0
U(i)=-1;
else
U(i)=Y4;
end
end
M=U
%绘图
i1=i
k=1:
1:
i1;
plot(k,U,k,U,'rx')
xlabel('k')
ylabel('M序列')
title('移位寄存器产生的M序列')
②程序运行结果如图2.8所示。
图2.8软件实现的移位寄存器产生的M序列图
.
③'四级移位寄存器产生的M序列
M=
Columns1through10
-11-11111-1-1-1
Columns11through20
1-1-111-11-111
Columns21through30
11-1-1-11-1-111
Columns31through40
-11-11111-1-1-1
Columns41through50
1-1-111-11-111
Columns51through60
11-1-1-11-1-111
i1=
60
最小二乘一次完成算法的MATLAB仿真
考虑仿真对象
式中
是服从正态分布的白噪声N(0,1),输入信号采用4阶M序列,幅度为1。
选择如下形式的辨识模型
构造
和
;数据长度取L=14;加权阵取
;利用
计算参数估计值
。
设输入信号的取值是从k=1到k=16的M序列。
待辨识参数
,观测矩阵
和
的表达式为:
,
,
仿真程序流程图:
%二阶系统的最小二乘一次完成算法辨识程序。
u=[-1,1,-1,1,1,1,1,-1,-1,-1,1,-1,-1,1,1];%系统辨识的输入信号为一个周期的M序列
z=zeros(1,16);%定义输出观测值的长度
fork=3:
16
z(k)=1.5*z(k-1)-0.7*z(k-2)+u(k-1)+0.5*u(k-2);%用理想输出值作为观测值
end
subplot(3,1,1)%画三行一列图形窗口中的第一个图形
stem(u)%画输入信号u的径线图形
subplot(3,1,2)%画三行一列图形窗口中的第二个图形
i=1:
1:
16;%横坐标范围是1到16,步长为1
plot(i,z)%图形的横坐标是采样时刻i,纵坐标是输出观测值z,图形格式为连续曲线
subplot(3,1,3)%画三行一列图形窗口中的第三个图形
stem(z),gridon%画出输出观测值z的径线图形,并显示坐标网格
u,z%显示输入信号和输出观测信号
%L=14%数据长度
HL=[-z
(2)-z
(1)u
(2)u
(1);-z(3)-z
(2)u(3)u
(2);-z(4)-z(3)u(4)u(3);-z(5)-z(4)u(5)u(4);-z(6)-z(5)u(6)u(5);-z(7)-z(6)u(7)u(6);-z(8)-z(7)u(8)u(7);-z(9)-z(8)u(9)u(8);-z(10)-z(9)u(10)u(9);-z(11)-z(10)u(11)u(10);-z(12)-z(11)u(12)u(11);-z(13)-z(12)u(13)u(12);-z(14)-z(13)u(14)u(13);-z(15)-z(14)u(15)u(14)]%给样本矩阵HL赋值
ZL=[z(3);z(4);z(5);z(6);z(7);z(8);z(9);z(10);z(11);z(12);z(13);z(14);z(15);z(16)]%给样本矩阵zL赋值
%CalculatingParameters
c1=HL'*HL;c2=inv(c1);c3=HL'*ZL;c=c2*c3%计算并显示
%DisplayParameters
a1=c
(1),a2=c
(2),b1=c(3),b2=c(4)%从
中分离出并显示a1、a2、b1、b2
%End
程序运行结果:
>>
u=[-1,1,-1,1,1,1,1,-1,-1,-1,1,-1,-1,1,1]
z=[0,0,0.5000,0.2500,0.5250,2.1125,4.3012,6.4731,6.1988,3.2670,-0.9386,-3.1949,-4.6352,6.2165,-5.5800,-2.5185]
HL=
001.0000-1.0000
-0.50000-1.00001.0000
-0.2500-0.50001.0000-1.0000
-0.5250-0.25001.00001.0000
-2.1125-0.52501.00001.0000
-4.3012-2.11251.00001.0000
-6.4731-4.3012-1.00001.0000
-6.1988-6.4731-1.0000-1.0000
-3.2670-6.1988-1.0000-1.0000
0.9386-3.26701.0000-1.0000
3.19490.9386-1.00001.0000
4.63523.1949-1.0000-1.0000
6.21654.63521.0000-1.0000
5.58006.21651.00001.0000
ZL=[0.5000,0.2500,0.5250,2.1125,4.3012,6.4731,6.1988,3.2670,-0.9386,-3.1949,-4.6352,-6.2165,-5.5800,-2.5185]T
c=[-1.5000,0.7000,1.0000,0.5000]T
a1=-1.5000
a2=0.7000
b1=1.0000
b2=0.5000
>>
从仿真结果表可以看出,由于所用的输出观测值没有任何噪声成分,所以辨识结果也无任何误差。
考虑图3.所示的仿真对象,图中,
是服从N
分布的不相关随机噪声。
且
,
,
选择图所示的辨识模型。
仿真对象选择如下的模型结构
其中,
是服从正态分布的白噪声N
。
输入信号采用4位移位寄存器产生的M序
列,幅度为0.03。
按式
构造
;加权阵取单位阵
;利用式递推公式计算l(k)、
和P(k),计算各次参数辨识的相对误差,精度满足要求式适当小的数(当所有参数估计值变化不大时)后停机。
下面给出具体程序。
%最小二乘递推算法辨识程序,
clear%清理工作间变量
L=15;%M序列的周期
y1=1;y2=1;y3=1;y4=0;%四个移位寄存器的输出初始值
fori=1:
L;%开始循环,长度为L
x1=xor(y3,y4);%第一个移位寄存器的输入是第三个与第四个移位寄存器的输出的“或”
x2=y1;%第二个移位寄存器的输入是第一个移位寄存器的输出
x3=y2;%第三个移位寄存器的输入是第二个移位寄存器的输出
x4=y3;%第四个移位寄存器的输入是第三个移位寄存器的输出
y(i)=y4;%取出第四个移位寄存器的幅值为"0"和"1"的输出信号,即M序列
ify(i)>0.5,u(i)=-0.03;%如果M序列的值为"1",辨识的输入信号取“-0.03”
elseu(i)=0.03;%如果M序列的值为"0",辨识的输入信号取“0.03”
end%小循环结束
y1=x1;y2=x2;y3=x3;y4=x4;%为下一次的输入信号做准备
end%大循环结束,产生输入信号u
figure
(1);%第一个图形
stem(u),gridon%显示出输入信号径线图并给图形加上网格
z
(2)=0;z
(1)=0;%设z的前两个初始值为零
fork=3:
15;%循环变量从3到15
z(k)=1.5*z(k-1)-0.7*z(k-2)+u(k-1)+0.5*u(k-2);%输出采样信号
end
%RLS递推最小二乘辨识
c0=[0.0010.0010.0010.001]';%直接给出被辨识参数的初始值,即一个充分小的实向量
p0=10^6*eye(4,4);%直接给出初始状态P0,即一个充分大的实数单位矩阵
E=0.000000005;%取相对误差E=0.000000005
c=[c0,zeros(4,14)];%被辨识参数矩阵的初始值及大小
e=zeros(4,15);%相对误差的初始值及大小
fork=3:
15;%开始求K
h1=[-z(k-1),-z(k-2),u(k-1),u(k-2)]';x=h1'*p0*h1+1;x1=inv(x);%开始求l(k)
k1=p0*h1*x1;%求出l的值
d1=z(k)-h1'*c0;c1=c0+k1*d1;%求被辨识参数c
e1=c1-c0;%求参数当前值与上一次的值的差值
e2=e1./c0;%求参数的相对变化
e(:
k)=e2;%把当前相对变化的列向量加入误差矩阵的最后一列
c0=c1;%新获得的参数作为下一次递推的旧参数
c(:
k)=c1;%把辨识参数c列向量加入辨识参数矩阵的最后一列
p1=p0-k1*k1'*[h1'*p0*h1+1];%求出p(k)的值
p0=p1;%给下次用
ife2<=Ebreak;%如果参数收敛情况满足要求,终止计算
end%小循环结束
end%大循环结束
c,e%显示被辨识参数及其误差(收敛)情况
%分离参数
a1=c(1,:
);a2=c(2,:
);b1=c(3,:
);b2=c(4,:
);ea1=e(1,:
);ea2=e(2,:
);eb1=e(3,:
);eb2=e(4,:
);
figure
(2);%第二个图形
i=1:
15;%横坐标从1到15
plot(i,a1,'r',i,a2,':
',i,b1,'g',i,b2,':
')%画出a1,a2,b1,b2的各次辨识结果
title('ParameterIdentificationwithRecursiveLeastSquaresMethod')%图形标题
figure(3);%第三个图形
i=1:
15;%横坐标从1到15
plot(i,ea1,'r',i,ea2,'g',i,eb1,'b',i,eb2,'r:
')%画出a1,a2,b1,b2的各次辨识结果的收敛情况
title('IdentificationPrecision')%图形标题
程序运行结果:
>>
c=
0.001000.0010-0.4984-1.2328-1.4951-1.4962-1.4991-1.4998-1.4999
0.001000.00100.0010-0.23500.69130.69410.69900.69980.6999
0.001000.25091.24971.06651.00171.00201.00020.99990.9998
0.00100-0.24890.75000.56680.50200.50160.50080.50020.5002
图3.8最小二乘递推算法的参数辨识仿真
-1.5000-1.5000-1.5000-1.4999-1.4999
0.70000.70000.70000.7000-0.7000
0.99990.99990.99990.99990.9999
0.50000.50000.50000.50000.5000
e=
000-499.42001.47340.21280.00070.00200.00040.0000
0000-235.9916-3.94160.00420.00700.00120.0001
00249.86123.9816-0.1466-0.06070.0003-0.0018-0.0003-0.0001
00-249.8612-4.0136-0.2443-0.1143-0.0007-0.0016-0.0012-0.0001
0.00010.0000-0.0000-0.00000.0000
0.0001-0.00000.00000.00000.0000
0.00010.00000.00000.0000-0.0000
-0.00040.0000-0.00000.0000-0.0000
表3.2最小二乘递推算法的辨识结果
参数a1a2b1b2
真值-1.50.71.00.5
估计值-1.49990.70.99990.5000
仿真结果表明,大约递推到第十步时,参数辨识的结果基本达到稳定状态,即a1=-1.4999,a2=0.7000,b1=0.9999,b2=0.5000。
此时,参数的相对变化量E
0.000000005。
从整个辨识过程来看,精度的要求直接影响辨识的速度。
虽然最终的精度可以达到很小,但开始阶段的相对误差会很大,从图3.8(3)图形中可看出,参数的最大相对误差会达到三维数