How can i multiply the transfer function with a vector)

71 ビュー (過去 30 日間)

採用された回答

Teja Muppirala
Teja Muppirala 2011 年 11 月 21 日
You can use the LSIM command.
For example:
s = tf('s');
G = (1.659e010*s^2 + 9.123e012*s + 4.147e014)/ (1.011e-005*s^7 + 0.02896*s^6 + 99.2*s^5 + 1.461e005*s^4+ 2.187e008*s^3 + 1.042e011*s^2 + 1.308e013*s+ 4.147e014)
t = 0:0.001:10;
u = sin(t.^3);
y = lsim(G,u,t);
plot(t,u,t,y);
legend({'Input (u)', 'Output (y)'})
  1 件のコメント
Mohamed Nassef
Mohamed Nassef 2011 年 11 月 21 日
Thanks Teja so much for answering my question. My question is now how can i get the input vector (u) if I have the output vector (y) and the transfer function (G). i mean how can i divide the output vector by the transfer function to get the input again?

サインインしてコメントする。

その他の回答 (2 件)

Teja Muppirala
Teja Muppirala 2011 年 11 月 21 日
You need to invert the transfer function. But simply taking 1/G leaves you with more zeros than poles, and so you cannot exactly solve the inverse problem. But you can get a good estimate if you multiply 1/G with a low pass filter to make it proper (number of poles = number of zeros). Then you can call LSIM once again. You will need a high sampling rate to do this accurately.
%%The forward problem
s = tf('s');
G = (1.659e010*s^2 + 9.123e012*s + 4.147e014)/ (1.011e-005*s^7 + 0.02896*s^6 + 99.2*s^5 + 1.461e005*s^4+ 2.187e008*s^3 + 1.042e011*s^2 + 1.308e013*s+ 4.147e014);
t = 0:0.00001:10;
u = sin(t.^3);
y = lsim(G,u,t);
plot(t,u,t,y);
legend({'Input (u)', 'Output (y)'})
%%The inverse problem:
% Step 1. Invert G
G_inverse = 1/G;
% Step 2. Make the low pass filter
lpf_order = numel(pole(G)) - numel(zero(G));
lpf_freq = 10*max(abs(real([pole(G); zero(G)]))); %Cutoff frequency 10x higher than highest zero/pole in G
[num,den] = butter(lpf_order,lpf_freq,'s');
% Step 3. Multiply the inverse by the low pass filter
G_inverse = G_inverse * tf(num,den)
% Step 4. Use LSIM with 'y' as the input
u_estimate = lsim(G_inverse,y,t);
figure;
plot(t,u,t,y,t,u_estimate);
legend({'Actual Input (u)', 'Output (y)', 'Estimated input (u estimate)'})
  1 件のコメント
Mohamed Nassef
Mohamed Nassef 2011 年 11 月 21 日
thanks so much..really appreciate it ..it works for me

サインインしてコメントする。


Alex
Alex 2011 年 11 月 20 日
What form is your transfer function in?
  1 件のコメント
Mohamed Nassef
Mohamed Nassef 2011 年 11 月 21 日
Dear Alex thanks for your help but it didnt work...the transfer function is in s domain and as follows : TF= (1.659e010 s^2 + 9.123e012 s + 4.147e014)/ (1.011e-005 s^7 + 0.02896 s^6 + 99.2 s^5 + 1.461e005 s^4+ 2.187e008 s^3 + 1.042e011 s^2 + 1.308e013 s+ 4.147e014)
when i apply ilaplace it gives me the following result and inside it there is term called r3 and i dont know what is r3 and still dont know how to multiply the result with a vector in time domain
sum((382493238368367552757760000000000000000*exp(r3*t) + 8414482309222611969638400000000000000*r3*exp(r3*t) + 15301574209142073065472000000000000*r3^2*exp(r3*t))/(65273803904821246875*r3^6 + 160265312512388584439808*r3^5 + 457479253027996880076800000*r3^4 + 539013861833793098219520000000*r3^3 + 605145439338041840762880000000000*r3^2 + 192215073248053527838720000000000000*r3 + 12064170624206046756864000000000000000), r3 in RootOf(s3^7 + (26710885418731430739968*s3^6)/9324829129260178125 + (146393360968959001624576*s3^5)/14919726606816285 + (14373702982234482619187200*s3^4)/994648440454419 + (7172094095858273668300800000*s3^3)/331549480151473 + (30754411719688564454195200000000*s3^2)/2983945321363257 + (1286844866581978320732160000000000*s3)/994648440454419 + 122397836277877616882483200000000000/2983945321363257, s3))

サインインしてコメントする。

Community Treasure Hunt

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

Start Hunting!

Translated by