How to manually replicate the bode() gain plot from a transfer function
28 views (last 30 days)
Show older comments
Hello,
I have been trying to replicate a bode gain plot from a given transfer function (see code below).
I am not seeing where I went wrong, but my manual plots looks like this:
But the actual bode plots should look like this:
I thought that I was setting up my amplitude equation exactly the same as how the bode() function operates, but apparently not. Can you see where I made the mistake?
Here is my code:
syms C R_1 R_2 R_3 L s V_o i omega
%Equations
R_1 = 10;
R_2 = 10;
R_3 = 10;
L = 0.001;
C = 2*10^-6;
Z_1 = R_1+L*s;
Z_2 = 1/(C*s)+R_2;
Z_3 = R_3;
Z_23 = 1/(1/(Z_2)+1/(Z_3));
Z_tot = Z_1+Z_23;
V_i = Z_1*i+V_o;
V_o = (V_i/(Z_1))/((1/((1/(C*s))+R_2))+1/R_3+1/(Z_1));
Transfer_func(s) = vpa(simplify(V_o/V_i),5);
amp = (sqrt((5.903*10^24)^2+(1.1806*10^20*omega)^2))/(sqrt((1.1806*10^25+(-2.3613*10^16*omega^2))^2+(9.4447*10^20*omega)^2))
amp_dB = 20*log(amp)
fplot(amp_dB,[0.001 1*10^6])
figure
subplot(2,1,1)
fplot(20*log(abs(Transfer_func)), [0.001 1*10^6])
legend('|H( j\omega )|')
title('Gain (dB)')
grid
subplot(2,1,2)
fplot(180/pi()*angle(Transfer_func), [0.001 1*10^6])
grid
title('Phase (degrees)')
xlabel('Frequency (rad/s)')
legend('\phi( j\omega )')
H = tf([118059162071741125000 5902958103587056517120000],[23611832414348225 944473296573929026712 11805916207174113034240000]);
figure
bode(H,{0.001,1*10^6})
0 Comments
Accepted Answer
Star Strider
on 4 Jun 2020
‘But the actual bode plots should look like this:’
It does.
Harkening back to yesterday:
syms C R_1 R_2 R_3 L s V_o i omega
%Equations
R_1 = 10;
R_2 = 10;
R_3 = 10;
L = 0.001;
C = 2*10^-6;
Z_1 = R_1+L*s;
Z_2 = 1/(C*s)+R_2;
Z_3 = R_3;
Z_23 = 1/(1/(Z_2)+1/(Z_3));
Z_tot = Z_1+Z_23;
V_i = Z_1*i+V_o;
V_o = (V_i/(Z_1))/((1/((1/(C*s))+R_2))+1/R_3+1/(Z_1));
Transfer_func(s) = vpa(simplify(V_o/V_i), 5);
Transfer_func(omega) = subs(Transfer_func, {s},{1j*omega});
figure
subplot(2,1,1)
fplot(20*log10(abs(Transfer_func)), [0.001 1E6])
set(gca, 'XScale','log')
legend('|H( j\omega )|')
title('Gain (dB)')
grid
subplot(2,1,2)
fplot(180/pi*angle(Transfer_func), [0.001 1E6])
set(gca, 'XScale','log')
grid
title('Phase (degrees)')
xlabel('Frequency (rad/s)')
legend('\phi( j\omega )')
produces:
Note the slight changes between your code and mine.
7 Comments
Star Strider
on 9 Jun 2020
I suspect Shane Palmer is encountering this for the first time and is experimenting in order to understand how it all works. (I did something similar when I first encountered them, so I understand.)
More Answers (0)
See Also
Categories
Find more on Get Started with Control System Toolbox in Help Center and File Exchange
Products
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!