Finding Upper and Lower Frequency Values in a Bode Plot

12 Ansichten (letzte 30 Tage)
Ryan Hampson
Ryan Hampson am 9 Apr. 2022
Kommentiert: Voss am 11 Apr. 2022
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

Akzeptierte Antwort

Voss
Voss am 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 Kommentare
Ryan Hampson
Ryan Hampson am 11 Apr. 2022
Thank you! Thank works beautifully! Do you have any suggestions on the slope/decade as well?
Voss
Voss am 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

Melden Sie sich an, um zu kommentieren.

Weitere Antworten (0)

Produkte


Version

R2021b

Community Treasure Hunt

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

Start Hunting!

Translated by