How can ı get transfer function.
显示 更早的评论
Dear sir.
I have input data on time and I have output data on time.
Have can I find transfer function. I attached data.csv and graphic.
Best regard.
采纳的回答
Star Strider
2018-2-27
This is an interesting problem, and eminently solvable with the System Identification Toolbox.
First, the Fourier transform indicates that there are 3 poles, one at the origin, one at about 550 KHz, and one at infinity (actually the Nyquist frequency, that is for all practical purposes infinity), and 2 zeros, between each pair of the poles. That is all we need to estimate the transfer function.
The Code —
[D,S,R] = xlsread('data.csv');
tv = D(:,1);
vi = D(:,2);
vi(1) = 0;
vo = D(:,3);
[ti,ia,ic] = uniquetol(tv, 4.8E-7); % Get Unique Times
vi = vi(ia); % Corresponding Input Voltages
vo = vo(ia); % Corresponding Output Voltages
Ts = mean(diff(ti)); % Sampling Interval (sec)
ti = [(0 : 5E-7 : 1E-6)'; ti+1.5E-6]; % Pad All Vectors With 3 Initial Time Samples
vi = [zeros(3,1); vi]; % Pad All Vectors With 3 Initial Time Samples
vo = [zeros(3,1); vo]; % Pad All Vectors With 3 Initial Time Samples
Fs = 1/Ts; % Sampling Frequency (Hz)
Fn = Fs/2; % Nyquist Frequency
L = numel(ti); % Vector LEngth
FTvi = fft(vi)/L; % Input Fourier Transform
FTvo = fft(vo)/L; % Output Fourier Transform
FTtf = FTvo./FTvi; % Transfer Function Fourier Transform
Fv = linspace(0,1,fix(L/2)+1)*Fn; % Frequency Vector
Iv = 1:numel(Fv); % Index Vector
figure(1)
plot(Fv, 20*log10(abs(FTtf(Iv))*2))
xlabel('Frequency (Hz)')
ylabel('Amplitude (dB)')
grid
dataobj = iddata(vo,vi,Ts); % Prepare Data For Identification
tfobj = tfest(dataobj, 3, 2); % Create Transfer Function Object
NumTF = tfobj.Numerator; % Transfer Function Numerator
DenTF = tfobj.Denominator; % Transfer Function Numerator
TF = tf(tfobj) % Transfer Function
figure(2)
plot(ti, vi, ti, vo) % Plot Input & Output
xlabel('Time (sec)')
ylabel('Amplitude (V)')
grid
[mag,phase,wout] = bode(tfobj);
figure(3)
subplot(2,1,1)
semilogx(wout, 20*log10(squeeze(mag)), '-b', 'LineWidth',1) % Bode Plot: Magnitude
ylabel('H(f) (dB)')
grid
subplot(2,1,2)
semilogx(wout, squeeze(phase), '-b', 'LineWidth',1) % Bode Plot: Phase
xlabel('Frequency (rad/sec)')
ylabel('Phase (°)')
grid
opt = stepDataOptions;
opt.StepAmplitude = vi(end);
[y,t] = step(tfobj, ti(end), opt); % Calculate Step Response
figure(4)
plot(t,squeeze(y)) % Step Response
title('Step Response of Estimated Transfer Function')
xlabel('Time (sec)')
ylabel('Amplitude (V)')
grid
produces:
TF =
From input "u1" to output "y1":
99.51 s^2 + 1.448e06 s + 3.853e09
--------------------------------------
s^3 + 4612 s^2 + 1.002e07 s + 1.076e10
Continuous-time transfer function.
The transfer function polynomial coefficients with full precision are in ‘NumTF’ and ‘DenTF’.
and this plot:
I will let you refine this, since you know your system better than I do.
6 个评论
bahadir safak
2018-2-28
Star Strider
2018-2-28
An .m file would just contain the code I posted.
‘Because the some code not working.’
That gives me no useful information. What parts of my code are ‘not working’?
NOTE — My code requires the System Identification Toolbox (link). You must have access to it for my code to work.
I am using R2017b, however all the functions I use here have been part of the System Identification Toolbox (link) for the last several releases.
bahadir safak
2018-3-4
编辑：bahadir safak
2018-3-4
Star Strider
2018-3-4
It worked when I ran it, as my code demonstrates. (I am using R2017b.)
I have no idea what the problem could be. Only you and your colleagues can troubleshoot the xlsread failure.
The csvread function will not work with your file. (I tried that first.) The textscan function could work. I never tried it.
Anyway, the easiest way to solve this, is for me to save your data to a .mat file and attach it here, so I did. Download the file and save it to your MATLAB path. I saved it as:
save('bahadir safak data.mat', 'D', 'S', 'R');
bahadir safak
2018-3-6
编辑：bahadir safak
2018-3-6
Star Strider
2018-3-6
As always, my pleasure.
更多回答（0 个）
类别
在 帮助中心 和 File Exchange 中查找有关 Linear Model Identification 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!
Translated by 
