Re: Plotting the frequency response of a filter in "real" Hz

[Guilherme Ritter]

2016-04-21 20:09 GMT-03:00 Maynard Wright <address@hidden>:
> On Wednesday, April 20, 2016 09:44:42 PM Guilherme Ritter wrote:
>> Hi everyone.
>>
>> I've just started in filter design at college and I'm learning to use
>> Octave for it. I want to see the frequency response of filters I design.
>> I've managed to find code on the internet, but the output's x axis is in
>> radian frequency. I'd like it to show "actual" Hz. For example, if the
>> cutoff frequency is 5,5 kHz, I'd like for it to be represented in the
>> plot's x axis at 5500 or 5,5.
>>
>> I've searched a lot but couldn't find anything, only some solutions that
>> work in MatLab but not in Octave. At college, I'm using Octave 4.0.0,
>> packages control 3.0.0 and signal 1.3.2, Windows 7 Enterprise x64. At home,
>> all the versions are up to date, Xubuntu 14.04 x64.
>>
>> I've found the code here:
>> https://ccrma.stanford.edu/~jos/fp/Example_LPF_Frequency_Response.html
>>
>> Can I use Octave's functions to get that plot the way I want it?
>>
>> Thanks in advance.
>
>
> If you are plotting using
>
> plot(w,abs(H));
>
> as in the first example in the code you reference,  w is in radians per 
> second.
> You can create another vector that contains corresponding frequencies in Hz
> and then use it in the plot:
>
> fHz = w / (2 * pi);
> plot(fHz,abs(H));
>
> You will be using w in all the calculations and fHz just for the plotting.
>
> Maynard Wright
>

@ Mike Miller

freqz (b, a, [], Fs);
Mostly worked, but the scale of the x axis doesn't seem quite right.
For example, I learned that the Butterworth filter, in the frequency
response, have an y axis value of exactly 0,5 at the cutoff frequency.
Using the code above, for cutoff = 10 kHz and sampling rate = 44,1
kHz, the value at y = 0,5 is 10.382,8 Hz. Just to make sure we're on
the same page, here's the code:

% code BEGIN
format long; pkg load signal;
sampling_frequency = 44100;
Nyquist_frequency = sampling_frequency / 2;
cutoff_frequency = 10000;
[filter_num, filter_den] = butter(10, cutoff_frequency / Nyquist_frequency);
[frequency_H, frequency_w] = freqz(filter_num, filter_den, [],
sampling_frequency);
plot(frequency_w, abs(frequency_H));
% code END

--//--

@ Sergei Steshenko

error: freqresp: second argument 'w' must be a real-valued vector of frequencies

That's the last line in this code:

% code BEGIN
format long; pkg load signal;
sampling_frequency = 44100;
Nyquist_frequency = sampling_frequency / 2;
cutoff_frequency = 10000;
[filter_num, filter_den] = butter(10, cutoff_frequency / Nyquist_frequency);
filter_s = tf(filter_num, filter_den);
filter_z = c2d(filter_s, 1);
number_of_points = 1000; % the bigger number of points, the better is
resolution of the plot
zfrequency_range = exp(pi * i * (0:number_of_points)
/number_of_points); % from 0Hz to Nyquist frequency
frequency_range = 0.5 * sampling_frequency * (0:number_of_points)
/number_of_points); % from 0Hz to Nyquist frequency
semilogx(frequency_range(2:end), 20 *
log10(abs(filter_z(zfrequency_range(2:end))))); % 'filter_z' is your
filter transfer function in 'z' domain
% code END

I used c2d as per the following link to transform the transfer
function from S domain to Z domain. I don't know if it's the right
way.
http://octave.1599824.n4.nabble.com/Converting-from-S-to-Z-domain-td3317272.html

And I added a missing open paren on the line that starts with
"frequency_range = ..." where I imagined it was supposed to be. I
tried with it in different positions but it didn't change the result.

--//--

@ Maynard Wright

That just makes the x axis be between 0 and 0,5. Not what I'm looking
for. I'm looking for Hz in the audio range.