Transfer function estimation from experimental data
64 views (last 30 days)
Show older comments
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
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.
Answers (1)
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.
0 Comments
See Also
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!