How to find the Transfer function of a Simulink output plot?

1 次查看(过去 30 天)
Hello My friends
I have a lot of figures and I want to find the transfer functions as I=Imax cos(ωt+ϴ) and V=Vmax cos(ωt+ϴ)
here is the values of V and I
Can you Advise?
I attached one of my plots.
I really appreciate your assistance.

采纳的回答

Sam Chak
Sam Chak 2023-9-28
Think you want to apply the curve-fitting method to fit the sinusoidal model to the time-series data.
load('V_I.mat')
Vtime = Vabc3.Time;
Vdat1 = Vabc3.Data(:,1);
ta = Vtime(1:1001);
Va = Vdat1(1:1001);
tb = Vtime(1002:2100);
Vb = Vdat1(1002:2100);
tc = Vtime(2101:3001);
Vc = Vdat1(2101:3001);
figure(1)
plot(Vtime, Vdat1), ylim([-1.5 1.5]), grid on
xlabel('Time (seconds)'), ylabel('V_{abc}')
[Vafit, gof1] = fit(ta, Va, 'sin1')
Vafit =
General model Sin1: Vafit(x) = a1*sin(b1*x+c1) Coefficients (with 95% confidence bounds): a1 = 0.998 (0.998, 0.998) b1 = 314.2 (314.2, 314.2) c1 = 1.031 (1.031, 1.031)
gof1 = struct with fields:
sse: 4.5659e-13 rsquare: 1.0000 dfe: 998 adjrsquare: 1.0000 rmse: 2.1389e-08
figure(2)
plot(Vafit, ta, Va), ylim([-1.5 1.5]), grid on
xlabel('Time (seconds)'), ylabel('V_{a}')
[Vbfit, gof1] = fit(tb, Vb, 'sin1')
Vbfit =
General model Sin1: Vbfit(x) = a1*sin(b1*x+c1) Coefficients (with 95% confidence bounds): a1 = 0.1984 (0.1983, 0.1985) b1 = 314.2 (314.1, 314.2) c1 = 7.268 (7.266, 7.27)
gof1 = struct with fields:
sse: 8.6937e-04 rsquare: 1.0000 dfe: 1096 adjrsquare: 1.0000 rmse: 8.9063e-04
figure(3)
plot(Vbfit, tb, Vb), ylim([-0.3 0.3]), grid on
xlabel('Time (seconds)'), ylabel('V_{b}')
[Vcfit, gof1] = fit(tc, Vc, 'sin1')
Vcfit =
General model Sin1: Vcfit(x) = a1*sin(b1*x+c1) Coefficients (with 95% confidence bounds): a1 = 0.9971 (0.9958, 0.9984) b1 = 314.2 (314.1, 314.2) c1 = 7.304 (7.291, 7.317)
gof1 = struct with fields:
sse: 0.1770 rsquare: 0.9996 dfe: 898 adjrsquare: 0.9996 rmse: 0.0140
figure(4)
plot(Vcfit, tc, Vc), ylim([-1.5 1.5]), grid on
xlabel('Time (seconds)'), ylabel('V_{c}')
  3 个评论
Sam Chak
Sam Chak 2023-11-2
Do you mean a single mathematical function that describes all 3 sine waves?

请先登录,再进行评论。

更多回答(2 个)

Sam Chak
Sam Chak 2023-11-2
I'm not an expert, but I have friends who are, and one of them said this is possible with this math function:
where
, and .
If you find the proposed math function useful, don't forget voting 👍 for the answer.
load('V_I.mat')
Vtime = Vabc3.Time;
Vdat1 = Vabc3.Data(:,1);
ta = Vtime(1:1001);
Va = Vdat1(1:1001);
tb = Vtime(1002:2100);
Vb = Vdat1(1002:2100);
tc = Vtime(2101:3001);
Vc = Vdat1(2101:3001);
subplot(211)
plot(Vtime, Vdat1, 'color', '#528AFA'), ylim([-1.5 1.5]), grid on
title('Plot from Data')
xlabel('Time (seconds)'), ylabel('V_{abc}')
[Vafit, gof1] = fit(ta, Va, 'sin1');
vaAmp = Vafit.a1; vaFreq = Vafit.b1; vaPhi = Vafit.c1;
[Vbfit, gof1] = fit(tb, Vb, 'sin1');
vbAmp = Vbfit.a1; vbFreq = Vbfit.b1; vbPhi = Vbfit.c1;
[Vcfit, gof1] = fit(tc, Vc, 'sin1');
vcAmp = Vcfit.a1; vcFreq = Vcfit.b1; vcPhi = Vcfit.c1;
t = linspace(0, 0.3, 3001);
t1 = 0.1;
t2 = 0.21;
Vaf = vaAmp*sin(vaFreq*t + vaPhi);
Vbf = vbAmp*sin(vbFreq*t + vbPhi);
Vcf = vcAmp*sin(vcFreq*t + vcPhi);
fun1 = sign(t - t1).*(Vaf/4 - Vbf/4); % sigma1
fun2 = Vaf/4 + Vbf/4 + Vcf/2 - fun1;
fun3 = Vaf/4 + Vbf/4 - Vcf/2 - fun1;
f = fun2 - sign(t - t2).*fun3; % <-- this math function
subplot(212)
plot(t, f, 'color', '#FA477A'), ylim([-1.5 1.5]), grid on
title('Plot using Math function')
xlabel('Time (seconds)'), ylabel('V')
  4 个评论
Atefeh
Atefeh 2023-12-13
Dear Sam
Would you mind telling me the reference of the attached formulas?
Sam Chak
Sam Chak 2023-12-14
The provided signum-based formula does describe the behavior of the piecewise function. However, I find it unnecessary to cite the formula, as it involves a clever but unintuitive manipulation of the sign functions to separate the sub-functions at different intervals in the domain.
Moreover, I have provided a more intuitively understandable formula in my Answer below. If you still wish to cite, consider referencing the special function used in the formulation, acknowledging the mathematical properties of that function.

请先登录,再进行评论。


Sam Chak
Sam Chak 2023-12-14
Given the piecewise function defined for
I propose a more intuitively understandable Heaviside-based formula, which I affectionately refer to as the Piecewise Function Put Together (PFPT) formula:
with .
In this formulation, the 1st term represents throughout the entire range of x. The 2nd term is activated by the Heaviside step function when x equals , plotting and canceling out from ​ onwards. Upon reaching , the Heaviside function triggers the 3rd term, plotting and canceling out from onwards.
If you appreciate the explanation of the PFPT formula, don't forget to vote 👍 for the answer as a token of appreciation.
load('V_I.mat')
Vtime = Vabc3.Time;
Vdat1 = Vabc3.Data(:,1);
ta = Vtime(1:1001);
Va = Vdat1(1:1001);
tb = Vtime(1002:2100);
Vb = Vdat1(1002:2100);
tc = Vtime(2101:3001);
Vc = Vdat1(2101:3001);
subplot(211)
plot(Vtime, Vdat1, 'color', '#528AFA'), ylim([-1.5 1.5]), grid on
title('Plot from Data')
xlabel('Time (seconds)'), ylabel('V_{abc}')
[Vafit, gof1] = fit(ta, Va, 'sin1');
vaAmp = Vafit.a1; vaFreq = Vafit.b1; vaPhi = Vafit.c1;
[Vbfit, gof1] = fit(tb, Vb, 'sin1');
vbAmp = Vbfit.a1; vbFreq = Vbfit.b1; vbPhi = Vbfit.c1;
[Vcfit, gof1] = fit(tc, Vc, 'sin1');
vcAmp = Vcfit.a1; vcFreq = Vcfit.b1; vcPhi = Vcfit.c1;
t = linspace(0, 0.3, 3001);
t1 = 0.1;
t2 = 0.21;
Vaf = vaAmp*sin(vaFreq*t + vaPhi); % f1, function at the 1st interval, t < t1
Vbf = vbAmp*sin(vbFreq*t + vbPhi); % f2, function at the 2nd interval, t1 < t < t2
Vcf = vcAmp*sin(vcFreq*t + vcPhi); % f3, function at the 3rd interval, t2 < t
%% Old Piecewise Function formula
% fun1 = sign(t - t1).*(Vaf/4 - Vbf/4); % sigma1
% fun2 = Vaf/4 + Vbf/4 + Vcf/2 - fun1;
% fun3 = Vaf/4 + Vbf/4 - Vcf/2 - fun1;
% fold = fun2 - sign(t - t2).*fun3; % <-- this math function
%% New Piecewise Function Put Together (PFPT) formula
fnew = Vaf + (Vbf - Vaf).*heaviside(t - t1) + (Vcf - Vbf).*heaviside(t - t2);
subplot(212)
plot(t, fnew, 'color', '#FA477A'), ylim([-1.5 1.5]), grid on
title('Plot using PFPT formula')
xlabel('Time (seconds)'), ylabel('V')

类别

Help CenterFile Exchange 中查找有关 Financial Toolbox 的更多信息

Community Treasure Hunt

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

Start Hunting!

Translated by