Transfer function estimation from experimental data

30 views (last 30 days)
I have a particular system's frequency response data from experimental measurements and I would like to estimate a transfer function. Based on what I've read here, I have the magnitude and phase response response stored as an idfrd and I have tried playing around with the system identification toolbox (transfer function and state space models) without much luck. Nothing I've come up with really seems to fit it well. Am I missing something here?
Here's what the Bode plot looks like. Actual data is also attached.
  2 Comments
Star Strider
Star Strider on 1 Jul 2018
It may help to have your original, time-domain data. There doesn’t appear to be any relevant information in the idfrd data you posted. There are no poles, and I cannot even identify any zeros.
david
david on 2 Jul 2018
Don't have the original time domain data that corresponds to this set of idfrd data, but I took another data set using a chirped sinusoid. Sampled at 100 kHz. File has the system input and output response in time domain.

Sign in to comment.

Answers (1)

Star Strider
Star Strider on 2 Jul 2018
This is the best I can do:
d = load('td_data.mat');
Fs = 100E+3; % Hz
tdu = d.input;
tdy = d.output;
L = numel(tdu);
t = linspace(0, L, L)/Fs;
figure
plot(t, tdu, t, tdy)
grid
Fn = Fs/2;
FTy = fft(tdy)/L;
FTu = fft(tdu)/L;
FTH = FTy./FTu; % Empirical Transfer Function
Fv = linspace(0, 1, fix(L/2)+1)*Fn;
Iv = 1:numel(Fv);
figure
semilogy(Fv, abs(FTy(Iv))*2, Fv, abs(FTu(Iv))*2, Fv, abs(FTH(Iv))*2)
grid
xlim([0 5E+3])
legend('Output', 'Input', '‘Transfer Function’', 'Location','E')
figure
plot(Fv, imag(FTH(Iv))) % Poles, Zeros On Imaginary Axis
grid
xlim([0 5E+3])
figure
plot(Fv, abs(FTH(Iv))*2)
grid
xlim([0 5E+3])
D = iddata(tdy, tdu, 1/Fs);
sys = tfest(D, 5, 5);
figure
bode(sys)
An easy way to estimate the number of poles and zeros your system has is to plot the imaginary part of the transfer function as a function of frequency. The poles and zeros will be readily visible. I counted 5 poles before your system degenerated into noise. This is in figure(3).
I would incrementally add more poles and the equivalent number of zeros until it shows the desired high-pass characteristic. It has a zero at the origin, a pole at infinity, and a maximum between about 250-2000 Hz (in figure(4)), so I would experiment to get those characteristics in ‘sys’. I began to see that with 10 poles and 10 zeros.

Categories

Find more on Linear Model Identification in Help Center and File Exchange

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!