Estimating a second order transfer function for a over damped second order system.

I have been tasked with estimating a second order transfer function for a set of data given that has been given to me. The form of this transfer function is;

G(S) = (k * w^2) / (s^2 + 2Cws + w^2)

where k is the steady state gain, w is the undamped natural frequency and C is the dampening ratio.

I have previously used the tfest function to estimate a first order transfer function which seems to be correct, but I am unsure how to create the second order transfer function in the overdamped shape that I need.

Im not looking for the complete answer as I still want to solve the problem myself, but if anybody could suggest some points where I may be going wrong that would be appreciated.

[Mathieu NOE]

hello

so If you were successfull on a first order system, what blocks you for a second order system ?

in my 2020b release the function for transfer estimate is : tfestimate

have you time domain data and already computed the complex transfer function Txy = tfestimate(X,Y) ?

all the best

Answer 1

[Nagasai Bharat]

Hi,

In order to do an estimate of a second order transfer function the tfest function can be used with suitable arguments of np,nz (number of poles and zeros in the transfer funcition) along with the data you have.

Do give a clear look into the tfest documentation for more precise information.

PS: tfestimate is a function in the signal processing toolbox which does the similar tasks as tfest in the system identification toolbox

[Gregory Woodhall]

Thankyou for your reply,

I feel like I should of been clearer in my orignal question, I understand how to use the tfest function and have read through its documentation. The issue im having is that when I review the poles of the system using the pole() function, I receive the following values;

-0.4087 + 0.1584i

-0.4087 - 0.1584i

This as far as im aware, is incorrect for what I am trying to acheive, which is an overdamped system. This should produce the same value for the imaginary axis as seen in the image bellow.

What I am trying to find out is if there is any way to adjust my simulation to meet these requirements, as after reviewing different forums and the documentation I am still lost.

[Mathieu NOE]

well

could you share your data and your code and explain why you expect to identify an overdamped system ?

[Gregory Woodhall]

I expect to identify an overdamped system because the question speceifically states before hand that the system is overdamped and that I am to create a simulation and compare the two, that is why I am asking if there is anyway to produce one.

The data given is the output of the actuator, the step input is 5 Amps and the data is taken across 200 seconds starting at 1, recorded in 1 second increments.

Answer 2

[Paul]

I haven't yet tried this myself, but in the System Identification Toolbox idproc looks like it will allow you to specify a second order model with two real poles and procest is the function to estimate the model parameters

doc idproc
doc procest

Answer 3

[Mathieu NOE]

so this is the poor's man (no ID toolbox ) , frequency domain ID

it may help you

all the best

clc,close,clear;

File1 = fopen('Set9_1.txt');

Set9_1 = fscanf(File1,'%f',200);

dt = 1;
T = 1:dt:200;
step_resp = Set9_1;

% Actuator_Exp = plot(T,step_resp,'r');

% convert to impulse response = time derivative of step response
imp_resp = [diff(step_resp)./dt; 0];

% frequency response function (impulse response = FIR filter)
Fs = 1/dt;
freq = logspace(-3,log10(Fs/3),100);
b = imp_resp';
a = [1 zeros(1,length(imp_resp)-1)];

[g,p] = dbode(b,a,dt,2*pi*freq);
p = p-p(1); % phase must start at zero

% from there we can compute all transfer function parameters 
static_gain = g(1);

% at cut off frequency (fc) : transfer function phase is - 90°
fc = interp1(p,freq,-90);      % fc =  0.0555 Hz
wc = 2*pi*fc;

% search for the Q factor : Q is the amplitude of the normalized (modulus = 1 at f = 0 Hz) 
% transfer function at  f = fc 
mod_at_fc = interp1(freq,g,fc);      
Q = mod_at_fc/(static_gain); % Q = 0.61

% let's compare the model frequency response function with the identified model
num = [0 0 static_gain];
den = [1/wc^2 1/(wc*Q) 1];

[g_id,p_id] = bode(num,den,2*pi*freq);

figure(2), 
subplot(2,1,1),loglog(freq,g,'b',freq,g_id,'r');grid
title('Sys TF');
ylabel('magnitude')
legend('measure','model')

subplot(2,1,2),semilogx(freq,p,'b',freq,p_id,'r');grid
ylabel('phase (°)')
xlabel('Frequency (Hz)')
legend('measure','model')
axis([min(freq) max(freq) -360 0]);

% step response of identified model

sr_id = step(num,den,T);

figure(3),
plot(T,step_resp,'b',T,sr_id,'r');grid
title('Sys step response');
legend('measure','model')

[Mathieu NOE]

hello

so , how are we doing here ?