Finding Upper and Lower Frequency Values in a Bode Plot

13 views (last 30 days)
Hello! I have a problem where I need to find the "stopband" frequencies that are +/- 20dB from the absolute minimum value of the function.
I have found the minimum value of the function to be 6.969 dB, which means that I am trying to find frequencies that are as close to 20.969 dB as possible. The lower frequency of the stopband was found using "q=hmag(abs(hmag-p)==min(abs(hmag-p)))" for a value of wL=70.27, but I seem to be unable to figure out how to isolate the closest magnitude value an wH for the frequency in the higher region of the band.
My upper frequency value should be 18117 (found via manually searching). Does anyone have any suggestions?
%% Problem 1
clear, clc, clf; close all;
format shortG
global w
w=logspace(log10(5),log10(50000),500);
[hmag,hphase]=FTM(w);
figure();
hold on
sp1=subplot(2,1,1);
sp1.Position=sp1.Position + [0.0 0.0 0.0 0.0];
semilogx(w,hmag,'LineStyle','-','color',[0,0,0.8]);
axis([4,55000,0,55]);
title('Bode Plot Magnitude Problem 1'); xlabel('\omega'); ylabel('|X(j\omega)|_{dB}');
%legend('|X(j\omega)|','Location','NorthEast');
grid on
ax=gca;
ax.XAxisLocation = 'origin';
ax.YAxisLocation = 'origin';
sp2=subplot(2,1,2);
sp2.Position=sp2.Position + [0.0 0.0 0.0 0.0];
semilogx(w,hphase,'Linestyle','-','Color',[0,1,0]);
axis([4,55000,(-pi/2)*1.1,(pi/2)*1.1]);
title('Bode Plot Phase Problem 1'); xlabel('\omega'); ylabel('\angleX(j\omega)_{radians}');
%legend('\angleX(j\omega)','Location','NorthEast');
grid on
ax=gca;
ax.XAxisLocation = 'origin';
ax.YAxisLocation = 'origin';
hold off
% Problem 2
hmin=min(hmag);
fprintf('Minimum Value of Bode Plot is: %.3fdB\n',hmin);
% Problem 3
p=hmin+20
q=hmag(abs(hmag-p)==min(abs(hmag-p)))
wmax=[w(hmag==hmin)]
wL=[w(hmag==q)]
% This part is where I manually searched for the values that I should have
hmag
s=hmag(445)
r=w(hmag==s)
% Functions List
function [out1,out2]=FTM(w)
a = (800+1j*w).*(1000+1j*w).*(1500+1j*w).*(2500+1j*w)./(800*(1j*w).*(600+1j*w).*(4000+1j*w).^1);
out1=20*log10(abs(a));
out2=unwrap(angle(a));
end

Accepted Answer

Voss
Voss on 9 Apr 2022
Since the magnitude of the frequency response decreases to its minimum value and then increases again, you can split it into two parts (before and after the minimum), and find the element nearest to 26.969 in each part.
%% Problem 1
clear, clc, clf; close all;
format shortG
global w
w=logspace(log10(5),log10(50000),500);
[hmag,hphase]=FTM(w);
figure();
hold on
sp1=subplot(2,1,1);
sp1.Position=sp1.Position + [0.0 0.0 0.0 0.0];
semilogx(w,hmag,'LineStyle','-','color',[0,0,0.8]);
axis([4,55000,0,55]);
title('Bode Plot Magnitude Problem 1'); xlabel('\omega'); ylabel('|X(j\omega)|_{dB}');
%legend('|X(j\omega)|','Location','NorthEast');
grid on
ax=gca;
ax.XAxisLocation = 'origin';
ax.YAxisLocation = 'origin';
sp2=subplot(2,1,2);
sp2.Position=sp2.Position + [0.0 0.0 0.0 0.0];
semilogx(w,hphase,'Linestyle','-','Color',[0,1,0]);
axis([4,55000,(-pi/2)*1.1,(pi/2)*1.1]);
title('Bode Plot Phase Problem 1'); xlabel('\omega'); ylabel('\angleX(j\omega)_{radians}');
%legend('\angleX(j\omega)','Location','NorthEast');
grid on
ax=gca;
ax.XAxisLocation = 'origin';
ax.YAxisLocation = 'origin';
hold off
% Problem 2
[hmin,idx]=min(hmag);
fprintf('Minimum Value of Bode Plot is: %.3fdB\n',hmin);
Minimum Value of Bode Plot is: 6.969dB
% Problem 3
p=hmin+20;
dmag = abs(hmag-p);
wL = w(find(dmag(1:idx-1) == min(dmag(1:idx-1)),1))
wL =
70.027
wH = w(idx+find(dmag(idx+1:end) == min(dmag(idx+1:end)),1))
wH =
18117
axes(sp1);
hold on
plot([wL wH],[p p],'ro')
% q=hmag(abs(hmag-p)==min(abs(hmag-p)))
% wmax=[w(hmag==hmin)]
% wL=[w(hmag==q)]
% This part is where I manually searched for the values that I should have
% hmag
% s=hmag(445)
% r=w(hmag==s)
% Functions List
function [out1,out2]=FTM(w)
a = (800+1j*w).*(1000+1j*w).*(1500+1j*w).*(2500+1j*w)./(800*(1j*w).*(600+1j*w).*(4000+1j*w).^1);
out1=20*log10(abs(a));
out2=unwrap(angle(a));
end
  2 Comments
Ryan Hampson
Ryan Hampson on 11 Apr 2022
Thank you! Thank works beautifully! Do you have any suggestions on the slope/decade as well?
Voss
Voss on 11 Apr 2022
You're welcome!
For the slope, maybe something like this
w=logspace(log10(5),log10(50000),500);
[hmag,hphase]=FTM(w);
[hmin,idx]=min(hmag);
m1 = median(diff(hmag(1:idx))./diff(log10(w(1:idx))))
m1 = -19.9570
m2 = median(diff(hmag(idx:end))./diff(log10(w(idx:end))))
m2 = 20.1949
semilogx(w,hmag,'LineStyle','-','color',[0,0,0.8]);
axis([4,55000,0,55]);
title('Bode Plot Magnitude'); xlabel('\omega'); ylabel('|X(j\omega)|_{dB}');
grid on
hold on
plot(w([1 idx idx idx end]),[ ...
hmag(1)+[0 m1*log10(w(idx)/w(1))] NaN ...
hmag(end)+[m2*log10(w(idx)/w(end)) 0]], ...
'--','LineWidth',2);
text(sqrt(w(idx)*w(1)),(hmag(1)+hmag(idx))/2,sprintf('%+.2f dB/dec',m1))
text(sqrt(w(idx)*w(end)),(hmag(end)+hmag(idx))/2,sprintf('%+.2f dB/dec',m2),'HorizontalAlignment','right')
function [out1,out2]=FTM(w)
a = (800+1j*w).*(1000+1j*w).*(1500+1j*w).*(2500+1j*w)./(800*(1j*w).*(600+1j*w).*(4000+1j*w).^1);
out1=20*log10(abs(a));
out2=unwrap(angle(a));
end

Sign in to comment.

More Answers (0)

Products


Release

R2021b

Community Treasure Hunt

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

Start Hunting!