www.gusucode.com > ident_featured 案例代码 matlab源码程序 > ident_featured/iddemo9.m

    %% Dealing with Multi-Variable Systems: Identification and Analysis
% This example shows how to deal with data with several input and output
% channels (MIMO data). Common operations, such as viewing the MIMO data,
% estimating and comparing models, and viewing the corresponding model
% responses are highlighted.

% Copyright 1986-2014 The MathWorks, Inc.

%% The Data Set
% We start by looking at the data set SteamEng.
load SteamEng

%%
% This data set is collected from a laboratory scale steam engine. It has
% the inputs 
% *Pressure* of the steam (actually compressed air) after the control valve, and
% *Magnetization voltage* over the generator connected to the output axis.
%
% The outputs are
% *Generated voltage* in the generator and
% *Rotational speed* of the generator (Frequency of the generated AC
% voltage).The sample time was 50 ms.

%%
% First collect the measured channels into an |iddata| object:
steam = iddata([GenVolt,Speed],[Pressure,MagVolt],0.05);
steam.InputName  = {'Pressure';'MagVolt'};
steam.OutputName = {'GenVolt';'Speed'};

%% 
% Let us have a look at the data
plot(steam(:,1,1))

%%
plot(steam(:,1,2))

%%
plot(steam(:,2,1))

%%  
plot(steam(:,2,2))
  
%% Step and Impulse Responses
% A first step to get a feel for the dynamics is to look at the
% step responses between the different channels estimated directly from
% data:

mi = impulseest(steam,50);
clf, step(mi)

%% Responses with Confidence Regions
% To look at the significance of the responses, the impulse plot
% can be used instead, with confidence regions corresponding to
% 3 standard deviations:
showConfidence(impulseplot(mi),3)

%%
% Clearly the off-diagonal influences dominate (Compare the y-scales!)
% That is, |GenVolt| is primarily affected by |MagVolt| (not much dynamics)
% and |Speed| primarily depends on |Pressure|. Apparently the response from
% |MagVolt| to |Speed| is not very significant. 

%% A Two-Input-Two-Output Model
% A quick first test is also to look a a default continuous time
% state-space prediction error model. Use only the first half of the data 
% for estimation:
mp = ssest(steam(1:250))

%%
% Compare with the step responses estimated directly from data:
h = stepplot(mi,'b',mp,'r',2); % Blue for direct estimate, red for mp
showConfidence(h)

%%
% The agreement is good with the variation permissible within the shown
% confidence bounds.

%%
% To test the quality of the state-space model, simulate it on the part
% of data that was not used for estimation and compare the outputs:
compare(steam(251:450),mp)

%%
% The model is very good at reproducing the Generated Voltage for the
% validation data, and does a reasonable job also for the speed. 
% (Use the pull-down menu to see the fits for the different outputs.)

%% Spectral Analysis
% Similarly, comparisons of the frequency response of |mp| with a spectral
% analysis estimate gives:
msp = spa(steam);

%%    
% |bode(msp,mp)|
clf, bode(msp,'b',mp,'r')

%%
% You can right-click on the plot and select the different I/O pairs for
% close looks. You can also choose 'Characteristics: Confidence Region' for
% a picture of the reliability of the bode plot.

%% 
% As before the response from |MagVolt| to |Speed| is insignificant and
% difficult to estimate.

%% Single-Input-Single-Output (SISO) Models
% This data set quickly gave good models. Otherwise you often have to try
% out sub-models for certain channels, to see significant influences
% The toolbox objects give full support to the necessary bookkeeping in
% such work. The input and output names are central for this.
%
% The step responses indicate that |MagVolt| primarily influences |GenVolt|
% while |Pressure| primarily affects |Speed|. Build two simple SISO model
% for this: Both names and numbers can be used when selecting channels.
m1 = tfest(steam(1:250,'Speed','Pressure'),2,1); % TF model with 2 poles 1 zero
m2 = tfest(steam(1:250,1,2),1,0) % Simple TF model with 1 pole.

%%
% Compare these models with the MIMO model mp:
compare(steam(251:450),m1,m2,mp)

%%
% The SISO models compare well with the full model. Let us now compare the
% Nyquist plots. |m1| is blue, |m2| is green and |mp| is red. Note that the
% sorting is automatic. |mp| describes all input output pairs, while |m1|
% only contains |Pressure| to |Speed| and |m2| only contains |MagVolt| to
% |GenVolt|.
clf
showConfidence(nyquistplot(m1,'b',m2,'g',mp,'r'),3)
 
%%
% The SISO models do a good job to reproduce their respective outputs.
%
% The rule-of-thumb is that the model fitting becomes harder when you add
% more outputs (more to explain!) and simpler when you add more inputs.

%% Two-Input-Single-Output Model
% To do a good job on the output |GenVolt|, both inputs could be used.

m3 = armax(steam(1:250,'GenVolt',:),'na',4,'nb',[4 4],'nc',2,'nk',[1 1]);
m4 = tfest(steam(1:250,'GenVolt',:),2,1);
compare(steam(251:450),mp,m3,m4,m2)

%%
% About 10% improvement was possible by including the input |Pressure| in  
% the models |m3| (discrete time) and |m4| (continuous time),
% compared to |m2| that uses just |MagVolt| as input.

%% Merging SISO Models
% If desired, the two SISO models |m1| and |m2| can be put together as one 
% "Off-Diagonal" model by first creating a zero dummy model:
mdum = idss(zeros(2,2),zeros(2,2),zeros(2,2),zeros(2,2));
mdum.InputName = steam.InputName;
mdum.OutputName = steam.OutputName;
mdum.ts = 0; % Continuous time model
m12 = [idss(m1),mdum('Speed','MagVolt')];    % Adding Inputs. 
                                             % From both inputs to Speed
m22 = [mdum('GenVolt','Pressure'),idss(m2)]; % Adding Inputs. 
                                             % From both inputs to GenVolt

mm = [m12;m22]; % Adding the outputs to a 2-by-2 model.

compare(steam(251:450),mp,mm)

%%
% Clearly the "Off-Diagonal" model |mm| performs like |m1| and |m2| in
% explaining the outputs.

%% Additional Information
% For more information on identification of dynamic systems using System
% Identification Toolbox(TM) visit the
% <http://www.mathworks.com/products/sysid/ System Identification Toolbox>
% product information page.