系统辨识程序(matlab编写)

%---------------------------理想系统的阶跃响应--------------------------------- A = [1 -1.5 0.7]; B = [0 1 0.5]; Model = idpoly(A,B); %理想系统模型 Step(Model, [0 100]); %理想系统模型的阶越响应 Grid; %---------------------------相关法辨识脉冲函数--------------------------------- %获得辨识所需数据 %1、白噪声 N=1000; %设置实验长度 A = [1 -1.5 0.7]; %建立系统模型 B= [0 1 0.5]; Model= idpoly(A,B); %理想系统模型 Model_white_error= idpoly(A,B,1); %模型中加入白噪声 U= iddata([],idinput(N,'prbs')); %设置输入信号 E= iddata([],idinput(N,'rgs')); %设置噪声信号为白噪声 Y= sim(Model_white_error,[U E]); %获得输出数据 %2、有色噪声 C = [1 0.5]; D = [1 0.5]; Model_1_error= idpoly(A,B,C); %模型中的噪声为有色噪声1 Model_2_error= idpoly(A,B,1,D); %模型中的噪声为有色噪声2 Y1 = sim(Model_1_error,[U E]); %获得输出数据 Y2 = sim(Model_2_error,[U E]); %辨识 %1、白噪声 Z=[Y U]; %得到输入和输出组 Z=dtrend(Z); %滤波处理 [ir,r,cl] = cra(Z,100,0,1); %采用相关分析法估计对象的脉冲响应 %其中ir为对象脉冲响应的估计 plot(ir(2:50)); %绘制辨识结果曲线 grid on; %2、有色噪声 Z1=[Y1 U]; %得到输入和输出组 Z1=dtrend(Z1); %滤波处理 [ir1,r1,cl1] = cra(Z1,100,0,1); %采用相关分析法估计对象的脉冲响应 %其中ir1为对象脉冲响应的估计 plot(ir1(2:50)); %绘制辨识结果曲线 grid on; Z2=[Y2 U]; %得到输入和输出组 Z2=dtrend(Z2); %滤波处理 [ir2,r2,cl2] = cra(Z2,100,0,1); %采用相关分析法估计对象的脉冲响应 %其中ir1为对象脉冲响应的估计 plot(ir2(2:50)); %绘制辨识结果曲线 grid on; %模型检验 plot(ir(2:50)); %绘制辨识结果曲线 grid on; hold on; impulse(Model,'r',50) %绘制实际系统的脉冲响应曲线 plot(ir1(2:50)); %绘制辨识结果曲线 grid on; hold on; impulse(Model,'r',50) %绘制实际系统的脉冲响应曲线 plot(ir2(2:50)); %绘制辨识结果曲线 grid on; hold on; impulse(Model,'r',50) %绘制实际系统的脉冲响应曲线 %---------------------------最小二乘法辨识系统的参数--------------------------------- %辨识系统的阶次 data=iddata(Y,U) NN = struc(1:2,1:2,1:2); Loss_functions = arxstruc(data,data,NN); order = selstruc(Loss_functions,'aic') % 辨识系统的阶,白噪声 data1=iddata(Y1,U) NN = struc(1:2,1:2,1:2); Loss_functions1 = arxstruc(data1,data1,NN); order1 = selstruc(Loss_functions1,'aic'); % 辨识系统的阶,有色噪声1 order1 = [order1(1),order1(2),1,order1(3)]; data2=iddata(Y2,U) NN = struc(1:2,1:2,1:2); Loss_functions2 = arxstruc(data2,data2,NN); order2 = selstruc(Loss_functions2,'aic'); % 辨识系统的阶,有色噪声2 order2 = [order2(1),order2(2),0,1,0,order2(3)]; %辨识系统的参数 Model_parameter = arx(data,order) % AR模型 Model_1_parameter = armax(data1,order1) % 用ARMAX模型 Model_2_parameter = pem(data2,order2) % 用DA模型 %检验模型 compare(data,Model_parameter); %预测输出与实际输比较 resid(Model_parameter,data); %模型预测误差及相关分析 compare(data1,Model_1_parameter); %预测输出与实际输比较 resid(Model_1_parameter,data1); %模型预测误差及相关分析 compare(data2,Model_2_parameter); %预测输出与实际输比较 resid(Model_2_parameter,data2); %模型预测误差及相关分析 %---------------------------其它辨识方法--------------------------------- %辅助变量法 Model_parameter_iv = iv4(data,order) % AR模型 %递推算法 [parameter1,parameter2] = rarx(data,order,'ff',0.95); %遗忘因子法,ARX模型 count=length(parameter1); A=[1,parameter1(count,1),parameter1(count,2)]; B=[0,parameter1(count,3),parameter1(count,4)]; Model_parameter_r=poly2th(A,B); Model_parameter_r=sett(Model_parameter_r,1); [parameter1,parameter2] = rarx(data,order,'ng',0.95); %梯度法,ARX模型 count=length(parameter1); A=[1,parameter1(count,1),parameter1(count,2)]; B=[0,parameter1(count,3),parameter1(count,4)]; Model_parameter_r=poly2th(A,B); Model_parameter_r=sett(Model_parameter_r,1) %---------------------------其它辨识模型--------------------------------- %BJ模型 NN=[2 1 2 2 1]; TH=BJ(Z,NN); present(TH); %OE模型 NN=[2 2 1]; TH=oe(Z,NN); present(TH) NN=[2 2 1]; TH=roe(Z,NN,'ff',0.98); %遗忘因子法 present(TH) %通用模型 NN=[2 1 1 1 1 1]; TH=pem(Z2,NN); present(TH); NN=[2 1 1 1 1 1]; TH=rpem(Z2,NN,'ff',0.98); %遗忘因子法 present(TH);