Obtain transfer function from (noisy) measurement data and use this transfer function for control in simulink
27 次查看(过去 30 天)
显示 更早的评论
I have a setup with an accelerometer and a speaker. The speaker is used excited with a random noise, and I measure the acceleration in g (after normalization of the output voltage). My question now is what would be the best solution to obtain the transfer function of the given real system in a certain frequency range? Ideally, I would like to avoid to have to guess the number of poles and zeros with a function like tfest.
Currently I have the following code for the estimation of the transfer function. With this solution I was able to obtain the transfer function, but I am at a loss how to implement it into Simulink.
load('datset_y_u.mat')
disp("calculating psd of input/output")
nfft = 32 * samples_per_frame_sensor;
window = nfft;
[pxx,fxx] = pwelch(excitation_signal,window, [], nfft, f_s_audio_sensor);
[pyy,fyy] = pwelch(raw_sensor_signal,window, [], nfft, f_s_audio_sensor);
disp("transfer function estimation")
[Txy,f] = tfestimate(pxx,pyy,window/2,[],nfft/2,f_s_audio_sensor);
figure(1);
subplot(3,1,1);
title("input");
plot(fxx, 10*log10(pxx))
xlabel('Frequency (Hz)')
ylabel('PSD (dB/Hz)')
xlim([0 2000]);
subplot(3,1,2);
title("output");
plot(fyy, mag2db(abs(pyy)))
xlim([0 2000]);
xlabel('Frequency (Hz)')
ylabel('PSD (dB/Hz)')
subplot(3,1,3);
title("transfer function");
plot(f,mag2db(abs(Txy)))
xlim([0 2000]);
xlabel('Frequency (Hz)')
ylabel('PSD (dB/Hz)')
0 个评论
回答(2 个)
Sam Chak
2022-8-10
Hi @tiessen
Do you want to estimate the frequency response the Transfer Function or a real mathematical function-type of the Transfer Function?
If it is the latter, then you should consider using the tfest() function to obtain the numerator and denominator coefficients of estimated model. If you know behavior of the system, then you can try guessing the number of poles and zeros.
Having the numerator and denominator coefficients, in Simulink, you can input the parameters in the Transfer Function block.
% estimate a transfer function Gp with the number of poles and the number of zeros
Gp = tfest(data, np, nz);
% If np < 4, then try using pidtune to autotune the PID gains for the estimated Gp
Gc = pidtune(Gp, 'pid');
Rajiv Singh
2022-8-10
Automatic order determination is not easy/trivial in general. The closest we have to this is automatic selection of state-space model order in ssest/n4sid commands, as in model = ssest(data, 'best')
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Classical Control Design 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!