Frequently Asked Questions on the Frequency Domain System Identification Toolbox

• How can I identify a stable system?
• How can I remove the effect of a known subsystem (some known poles and/or zeros) from the data?
• How shall I generate time domain data objects from my measured data?
• How shall I deal with frequency response function (FRF) data, that is, magnitude (gain) and phase?
• How shall I average several experiments with changing (maybe random) phases?
• How can I generate a suggestion or bug report which can be quickly processed?
• My fdident copy does not run under Matlab 6.1. What shall I do?
• The GUI behaves oddly under Windows XP. What shall I do?

Back to the Developers' Page

How can I identify a stable system?

Making a good fit with stable poles is sometimes rather difficult. First of all, it is quite likely that the data themselves contain details which drive the identification routine to an unstable fit which is otherwise the best available fit it LS sense. In other words: within the given quality of the data, an unstable model is the best approximation. This means that the data are not good enough to guarantee stability. E.g. in case of slight nonlinear distortions, the best linear approximation might be non-stable indeed.

In summary, possible reasons for unstability:

• nonlinearities
• special pattern of noise
• too simple model for a complex system
• unstability of the system itself (an unstable system can be measured within a stabilizing feedback loop)
• local minimum found (very rare)

Therefore, making repeated and improved mesurements (better SNR, odd multisine etc) is the first advice.

In the identification of high-order systems there is a good chance that some poles will be driven to the unstable region by the noise or by slight nonlinearities. A natural instinct of a researcher is to try yo eliminate these poles. Unfortunately, this is usually not a remedy. Eliminated poles may leave a bad fit which cannot be corrected any more.

The Frequency Domain System Identification Toolbox offers some (artificial) tools to force stable solutions. However, we would like to emphasize here that these are artifical solutions: even reaching the best stable fit is not guaranteed. However, weakly defined poles (overmodeling) may be effectively driven to the stable region.

• Command line
• request in elis that it stabilizes all the poles by reflection or contraction.See "help elis", "elis runmod"

Warning! elis will allow the increase of the cost function when imposing stability. This means that iteration may easily diverge. You may want to observe the evolution of the cost function (see elis fitinfo), and set itmax to the desired value in a second iteration cycle.
• in order to avoid that stabilization destroys the fit, if you are looking for minimum phase zeros, also request the reflection/contraction of the zeros
• along with stabilization, allow the delay be free
• experiment with different starting values of the delay - usually there is a delay value which makes the fit stable
• after an unstable fit with elis (use the new form, see "help elis"), execute  stabmod=stable(model,'with-zeros');
• Graphical User interface
• switch UserLevel to Advanced, and request stabilization in elis, with the above options
• write a vector of starting values into the delay field (works both with fixed or variable delay)

In high-order mechanical systems, the system equations may be badly conditioned. This may cause unstability again. A possible remedy is to switch to orthogonal polynomials.

• GUI: choose "Iimproved numerical stability" in the "Estimate Plant Model" window
• command line: see "elis runmod".

How can I remove the effect of a known subsystem from the data?

Sometimes a part (maybe some poles) of the system are almost perfectly known, and it is necessary to identify only the rest of the dynamics. In this case, we can remove the effect of these poles/zeros from the data.
Warning! We can remove the effect of exactly given poles/zeros only. If they are imprecise, a part of the dynamics they cause can remain in the data.

We enclose a self-explanatory example.

%1. PREPARATIONS
%Close all windows first, and clear previous data:
close all, clear, figure(1)
%Define full system model:
model=fidmodel('s',real(poly([1+j,1-j])),real(poly([-1-j,-1+j,-0.1,-0.5])));
plot(model) %plot original model
%Define helper data object for simulation:
Fdat=fiddata(NaN*ones(10,1),ones(10,1),[0.05:0.05:0.5]);
Fdat.variance=[0.001,0.001];
%Simulate noisy measurement:
data=simfou(model,Fdat);
%Plot simulated data (simulated from original model)
figure(2), plot(data), figure(3)
%
%2. MAIN PART
model_est=elis(data,'s',2,4); %estimate original model
model_known=fidmodel('s',1,real(poly([-0.1]))); %known subsystem
datam=data/model_known; %Modify data, remove effect of one pole
%Plot modified data:
figure(4), plot(datam), figure(5)
model_partial=elis(datam,'s',2,3); %estimate remaining model

Generating time domain data objects

When the data are available as workspace variables or ASCII files or variables in MAT-files, the "Compose" command in the "Read Time Domain Data" block allows to easily assemble them. The fields of the Compose window can be filled in with valid MATLAB commands, like variables, subscripted variables like data(:,2), or special file references (use "Browse" at the appropriate places). When assembling the objects, the function checks the data if they are consistent indeed, and give warnings if they are not.

A special call allows the independent use of the Compose window:

     fdtool('callback','compinp','init1') %compose time domain objects

In time domain, automatic determination of the period length, the excitation frequencies, and the Fourier coefficients is provided. However, if you can set the first two or at least the second one properly, this can improve speed and calculation accuracy. So, if you can, do not forget to fill in the "Frequencies" and/or the "PeriodLength" property, e.g.

     data.frequencies = 1/0.03*[1:2:15];
     data.periodlength = 0.03;

If the data of more than one experiments are used, and measurements are properly synchronized, set also the property "synchronization", e.g.

     data.synchronization=’on’;

otherwise unnecessarily long processing time may be wasted later, during the execution of the "Variance Analysis" block.

Sample programs for preparing load

Simple example

Let us assume first that we have measured 4096 time domain samples, which are stored in the vectors inp and outp. The sampling frequency is fs = 40960 Hz, the period length in each signal is Tp = 1/80 s (there are 8 periods). The frequency content is fv = [1:2:15]/Tp.

%Make the object:
obj = tiddata(outp,inp,1/40960);
%Set some properties:
obj.periodlength = 1/80;
obj.frequencies = [1:2:15]/80;
get(obj) %check properties

Now the variable obj can be imported to the block "Read Time Domain Data".

Multiple experiments

As a more complex example (multiple experiments), let us assume that we have measured five times 1024 time domain samples, which are stored in the 1024 x 5 arrays inp and outp. The sampling frequency is fs = 40960 Hz, the period length in each signal is 1/160 s (there are 4 periods in each record).

%First prepare the data in 1 x 5 cell arrays:
u = num2cell(inp,1); y = num2cell(outp,1);
%Now make the object:
obj = tiddata(y,u,1/40960);
%Set some properties:
obj.periodlength = 1/160;
obj.synchronization=’on’;
get(obj) %check properties

Now the variable obj can be imported to the block "Read Time Domain Data".

How shall I deal with frequency response function data?

It may happen that we need to start identification from frequency response function (FRF) data: magnitude (gain) and phase. The fdident toolbox is able to deal with these. Here are some pssibilities to illustrate how it works.

When the data are available as workspace variables or ASCII files or variables in MAT-files, the "Compose" command in the "Read Frequency Domain Data" block allows to easily assemble them. The fields of the Compose window can be filled in with valid MATLAB commands, like variables, subscripted variables like data(:,2), or special file references (use "Browse" at the appropriate places). When assembling the objects, the function checks the data if they are consistent indeed, and give warnings if they are not.

Especially for FRF data, please set the "Data type" menu item to the proper value. The amplitudes are either in dB-degrees vectors, or in complex form.

A special call allows the independent use of the Compose window:

     fdtool('callback','compinp','init2') %compose frequency domain objects

Command-line solution
%***************************
%
%1. First define the data to be able to run this sequence...
%In practice, data often arrive in ASCII files. In such a case,
%load these data into MATLAB, and form vectors. Then you are done with this part.
%See 'help load' or 'help loadasc'.
%If you have the data in a Word table, do the following:
%  Select the table
%  Through the menu, Table / Convert / Table to text
%      Select the separator as Tab mark
%  Delete everything else than the data (the numbers) themselves
%  Save it as a txt file
%  Execute in Matlab
%      load data.txt -ascii
%      Now you have an array which corresponds to the table
%
%The following 3 lines are just preparations for simulation, not needed for true data
F=20; freqv=[1:F]'/F*1e3; %frequency vector (NOT radian frequencies), 1/20 kHz steps 
jomega=j*2*pi*freqv; K=2e-4; %Time constant in seconds
tf=1./(1+jomega*K); %Complex amplitudes, 1st order system
%
%Here are a few examples how to deal with the vectors:
tf_mag=abs(tf); %Magnitude vector
%tf_mag_dB=20*log10(tf_mag); %This would be the magnitude in decibels
tf_phase=angle(tf); %Phase, radians
%tf_phase_degrees=tf_phase/(2*pi)*360; %this would be the phase in degrees
U=ones(size(tf_mag)); %the input amplitudes are all equal to one with zero phase
Y=tf_mag.*exp(j*tf_phase); %Complex output Fourier amplitudes
%or, starting from dB-degrees data:
%Y=10.^(tf_mag_dB/20).*exp(j*tf_phase_degrees/360*2*pi);
%
%*****************************
%
%2. Now form an object from the data. We will use here freqv, Y, and U.
dat=fiddata(Y,U,freqv);
%If you have several repeated experiments:
%dat=fiddata({Y1,Y2,...,YN},{U,U,...,U},freqv);
%If you have variance data of the Fourier amplitudes (scalars or vectors):
%dat.variance=[varY,varU]; %variance column vectors, to prepare the fit
%if the variances are not known, you may leave the variances empty
%
%Remark: advanced users may also use the coherence data instead, to obtain 
%approximate variances, assuming that small coherence is caused by noise alone.
%If these are in the vector 'cohvect',
%  vartf=coh2var(cohvect,Y./U);
%  dat.variance=[vartf,zeros(size(U))];
%Warning! Make sure that no measured coherence is over 1 (this might happen
%in some instruments because of roundoff errors).
%
figure(1)
plot(dat,'=')
drawnow
%
%*****************************
%
%3. Make an estimation
figure(2)
model=elis(dat,'s',0,1);
figure(3), plot(model)
%
%*****************************

Alternatively, a quick way when you already have an fiddata object dat,

• start fdtool,
• open the block "Read Frequency Data", and import your data from the workspace, then push "Close",
  • if you are at Automatic Userlevel (see menu "UserLevel" in the main window), the block "Estimate Plant Model" will automatically open up,
  • if you are at Interactive userlevel (see menu "UserLevel" in the main window), open the block "Variances and/or Averaging", push "Execute", push "Close", and open the block "Estimate Plant Model".
• select proper numerator and denominator orders and make a fit.

Averaging several experiments with changing (maybe random) phases

When nonlinearities may be present, it is a good practice to make several multisine experiments with changing (random) phases in each experiment. In this case, varanal cannnot be used without special manipulations, since the changing phases cannot be brought to the same position.

There are two solutions to tackle this problem.

Correcting phases one by one

When generating the excitation signal, you have the phases of each sine. Let us denote the vector of phases in radians by phases. Here is a sample program for each experiment.

%Do this for each experiment
tdati=tiddata(yi,ui,Ts);
tdati.frequencies=freqvect;
Fdati=tim2fou(tdati);
Fdati.output=exp(-j*phases).*Fdati.output;
Fdati.input=exp(-j*phases).*Fdati.input;
%Now the phases are not influenced by the excitation
%
Fdatall=merge(Fdat1,Fdat2,...); %Multiple experiments
Fdatall.synchronization='on'; %Tell varanal to simply average
Fdat=varanal(Fdatall); %This can be used by elis or by the GUI

Convert data to FRF

Another possibility is to convert the data to FRF = Y/U. By this, the Gaussian nature of the noise is somewhat violated (not much for good input SNR), and the variance of Y/U is approximated, but usually this is not a big problem.

This method works properly in fdident toolbox version 3.1b (21-Jun-2002 or later only) - if  fdtool version  returns something earlier, please update from the web, from index.html.

%Do this for each experiment
tdati=tiddata(yi,ui,Ts);
tdati.frequencies=freqvect;
Fdati=tim2fou(tdati);
Fdatfrfi=stdtfm(Fdati);
%
Fdatall=merge(Fdatfrf1,Fdatfrf2,...); %Multiple experiments
Fdatall.synchronization='on'; %Tell varanal to simply average
Fdat=varanal(Fdatall); %This can be used by elis or by the GUI

How can I prepare and send a report on the toolbox?

If you have suggestions or remarks on the fdident toolbox, or if you want to ask something, or if you want to ask for the correction of a bug you seem to have found, the quickest and most effective procedure is as follows.

Report on the GUI

Start the gui (fdtool), and open the Action Recorder. You can record your GUI actions with the use of it. Then, after STOP, you can save the history record into a MAT-file. This can be replayed anywhere under the fdident toolbox, if you carefully follow the guidelines below.

You want to have your actions you recorded be reproducible by the developers. This will only be successful if the GUI is started from exactly the same state. Therefore, at the beginning of each history record, we suggest to do precise GUI status setting. Please push RECORD, and do the following.

1. "File / Clear / Whole GUI" (set standard starting status), and
2. "UserLevel / Interactive", or "UserLevel / Automatic", or "UserLevel / Advanced", as you prefer.

Now you can record your actions. When you are ready, push STOP, then select "File / Save history as...", and give it a file name with the extension .mat .

It is also possible to automatically load your data (e.g. from a MAT-file) to the workspace before your actions. To do this, push STOP after setting the UserLevel (and before recording your actions), select in the recorder "Edit / Insert MATLAB Command", and in the field "Param" type 'load .mat'. Now push PLAY to really execute the load command. After this, and pushing RECORD, you can record your actions.

It is also possible to insert the load command later. Write e.g. number 3 into the Index field of the recorder, push , select in the recorder "Edit / Insert MATLAB Command", and then insert the load command into the field "Param". It will be inserted before action 3.

You can also extend the history file to send more information. If you click on "Pause" at an action, later the history file will stop running before this action. You can even type your comments into the big edit box, and the developers will see them when reproducing. After editing, do not forget to save the history file again.

You might have also introduced some superfluous actions during recording. Please delete these (Edit / Delete menu item - this acts on the current action). After getting more familiar with the recorder, you can even insert actions into a history: go to the action you want to make insertion before, push RECORD, and do your recording. Edit / Cut and Edit / Paste also work for history actions.

To be totally sure that your history can be reproduced elsewhere, please

• save the history into a mat-file, and
• execute the following commands:
    fdtool exit force %close all GUI-related windows
    fdtool quick %this starts the GUI without any local settings
    runhist .mat

Please write an email to fdident@vub.ac.be with a short decription of the problem, and also give the results of the following commands (exact reproduction may depend on using the same toolbox and Matlab versions, moreover, common error sources can often be localized by these informations):

version
computer
fdident version
fdident license
which elis -all
elis date

Attach the history file prepared as above, attach also the MAT-file with the data (if you have some), and send the email.

Report on the command-line toolbox

Making a report on the command-line toolbox is maybe even easier. Please write a short M-file which reproduces the problem or question. Please write an email to fdident@vub.ac.be with a short decription of the problem, and also give the results of the following commands (exact reproduction may depend on using the same toolbox and Matlab versions, moreover, common error sources can often be localized by these informations):

version
computer
fdident version
fdident license
which elis -all
elis date

Attach the M-file, and also the MAT-file with the data (if you have some), and send the email.

Running FDIDENT under MATLAB 6.1 or later

Switching from MATLAB 6.0 to 6.1 or later needs the following steps:

• install MATLAB 6.x
• the passcode for installation of fdident is different under MATLAB 6.0, 6.1 or 6.5. If you have a valid fdident license, please contact fdident@gamax.hu for a new passcode for installation.
• install fdident from the CD, or from the WEB
• fdident versions 3.0f and earlier need an additional new file to be able to run properly under Matlab 6.x. Please go to directory toolboxfdidentfdident@iddat, and create a file numel.m there with the following contents:

  function ne2 = numel(dat,varargin)
  ne2=1;
   
  

  The version on the WEB server already contains this file, so an alternative is to upgrade fdident from the WEB.

Run the GUI by typing fdtool under MATLAB.

The GUI behaves oddly under Windows XP. What shall I do?

If MATLAB 6.1 runs under Windows XP, in some cases a strange phenomenon ccurs. If someone wants to double click on e.g. the arrows of the GUI, the panel that pops-up has only the upper bar, the rest of it being "hidden" - so no actions can be initiated on the arrow data.

The cause is most probably that MATLAB 6.1 was earlier issued than XP, therefore information exchange is not perfect. The easiest remedy: upgrade MATLAB to version 6.5 (or use MATLAB 6.1 under Windows 2000 :-))