Determine -3db gain on Freq- gain plot
44 views (last 30 days)
Show older comments
Moussa Ihab
on 18 Nov 2022
Commented: Star Strider
on 18 Nov 2022
Hello,
I got some freq-gain data from a filter and I am trying to determine -3db gain on the plot and get its frequency.
I am trying "interp1" but it comes with only one point.
How can I get the other one ?
Thank you
freq=t1.f;
noise=t1.noise;
gain=t1.gain;
figure
plot(freq,noise)
xlabel('Frequency')
ylabel('noise')
figure
[Max_gain,index1]=max(gain);
freq_max=freq(index1);
Gain_3db=Max_gain-3;
freq_3db=interp1(gain,freq,Gain_3db)
plot(freq,gain,'b')
ylabel('gain')
xlabel('Frequency')
hold on
plot(freq_3db,Gain_3db,'r*')
0 Comments
Accepted Answer
Star Strider
on 18 Nov 2022
Edited: Star Strider
on 18 Nov 2022
Use:
find(diff(sign(gain-Gain_3db)))
or something similar, to find the indices of all the points where the curve equals the specified value, then include a short range (I use [-1 0 1]) so that interp1 can interpolate over a limited set of values. This prevents problems with functions that can vary significantly, especially those that change frequently.
This approach works regardless of the number of times the transfer function equals the desired decibel value (so it would work for a filter with multiple passbands and stopbands), so long that at least one does.
I do not have your data, however this should work as an illustration —
freq = linspace(0, 10, 150).'; % Assuming Column Vectors
gain = mag2db(exp(-(freq-5).^2)); % Create Function, Copnvert To 'dB'
[Max_gain,index1]=max(gain);
freq_max=freq(index1);
Gain_3db=Max_gain-3;
idxr = find(diff(sign(gain-Gain_3db)))+(-1:1) % Matrix Of Index Ranges
for k = 1:size(idxr,1)
freq_3db(k) = interp1(gain(idxr(k,:)),freq(idxr(k,:)),Gain_3db);
end
freq_3db % Desired Frequencies
figure
hp1 = plot(freq, gain, 'DisplayName','Transfer Function');
hold on
hp2 = plot(freq_3db, Gain_3db, '+r', 'DisplayName','Transfer Function -3 dB Points');
hold off
grid
legend([hp1; hp2(1)], 'Location','best')
ylim([-20, 5])
Make appropriate changes to work with your data.
EDIT — (18 Nov 2022 at 13:15)
Clarified description. Code unchanged.
EDIT — (18 Nov 2022 at 19:48)
With the provided data —
freq = [2400 2470 2540 2610 2680 2750 2820 2890 2960 3030 3100 3170 3240 3310 3380 3450 3520 3590 3660 3730 3800 3870 3940 4010 4080 4150 4220 4290 4360 4430 4500 4570 4640 4710 4780 4850 4920 4990 5060 5130 5200];
gain = [-18.6680 -14.4060 -10.2174 -6.4374 -3.7861 -2.4066 -1.7975 -1.6368 -1.7193 -1.7336 -1.7873 -1.4274 -1.7189 -1.3573 -1.6372 -1.3042 -1.2758 -1.1481 -1.0076 -1.1632 -1.0454 -1.1071 -1.1686 -1.1704 -1.3548 -1.2611 -1.4895 -1.3286 -1.6556 -1.4886 -2.1840 -2.1738 -3.9351 -3.9649 -6.5240 -7.0555 -10.1172 -11.1532 -13.9887 -16.3625 -18.6674];
freq = freq(:);
gain = gain(:);
[Max_gain,index1]=max(gain);
freq_max=freq(index1);
Gain_3db=Max_gain-3;
idxr = find(diff(sign(gain-Gain_3db)))+(-1:1) % Matrix Of Index Ranges
for k = 1:size(idxr,1)
freq_3db(k) = interp1(gain(idxr(k,:)),freq(idxr(k,:)),Gain_3db);
end
fprintf('-3 dB Frequency = %.3f Hz\n',freq_3db)
% freq_3db % Desired Frequencies
figure
hp1 = plot(freq, gain, 'DisplayName','Transfer Function');
hold on
hp2 = plot(freq_3db, Gain_3db, '+r', 'DisplayName','Transfer Function -3 dB Points');
hold off
grid
legend([hp1; hp2(1)], 'Location','best')
ylim([-20, 5])
txtc = compose(' \\leftarrow %9.3f Hz',freq_3db);
text(freq_3db+15, [1 1]*Gain_3db, txtc, 'Horiz','left', 'Vert','middle', 'Rotation',-90)
Make appropriate changes to get different results.
.
2 Comments
More Answers (4)
cr
on 18 Nov 2022
interp1() outputs linearly interpolated values that you are asking for. if Gain_3db is a scalar value the output is a scalar. The output you get this way gives you the freq at which gain fell by 3dB from its max value.
Perhaps the other value you are looking for is the start freq.
This is also not a great way to do what you seem to be looking for. It only works if gain values are monotonic.
0 Comments
Askic V
on 18 Nov 2022
I think this code can help you, it is something I often use. It is not most efficient, but it works.
Please have a look at the example:
clear
clc
close all
num = 1;
den = [1 1];
g = tf(num,den);
[mag, phase, w] = bode(g);
magdb = squeeze(20*log10(mag));
semilogx(w,magdb);
ylim([-15, 0]);
grid on;
ind = find(magdb>=-3,1,'last');
% find two nearest points
P1w = [w(ind), w(ind+1)];
P2m = [magdb(ind), magdb(ind+1)];
% Find slope and intercept
c = [[1; 1] P1w(:)]\P2m(:);
% y = ax+b;
a = c(2);
b = c(1);
% Find angular frequency where magdb is -3 dB
w_3db = (-3-b)/a;
% Check if this is the case
mag_3db = a*w_3db+b;
hold on
plot(w_3db, mag_3db,'ro');
text(w_3db, mag_3db, [' -3 dB point,', newline, 'w = ', num2str(w_3db), 'rad/s'])
0 Comments
Moussa Ihab
on 18 Nov 2022
1 Comment
Star Strider
on 18 Nov 2022
It would help to have your data. Please post the frequency and gain vectors.
See Also
Categories
Find more on Fourier Analysis and Filtering in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!