matlab How to design lpfbpfhpf without builtin functions

related to Filtering voice signal from randomly added frequencies

I'm trying to create a simple low pass filter but i'm not sure what is wrong.

my code:

[y,Fs]=wavread('sound.wav');
m = length(y); %30586
Fc = 500;               %% frequency cut
n = [0:m-1];            %% kernel length / samples
Fs = Fs;    %22000      %% sample frequency / (nyquist - %2*m)

Tc = 2*pi*(Fc/Fs);      %% theta c
Hm = sin(Tc*n) ./ (pi*n);
Hm(1) = Tc/pi;
figure,plot(Hm);

1.

LPF: How to get to the 0-500 range?

We used this code in the class but with Fc=45,m=100,Fs=200 and it looks ok. but when I try this on this wav file it's not working as expected. I thought it will pass from 0-500 but instead it returns some crowled plot.

only when I set Fc to be a low number I can see the sin wave. so i've try to manipulate by setting Fc=real(ceil(m/500)). it better than nothing but it's not it.

2.

hpf = 1-lpf, right? well is seems not to be. my guess is that my lpf is wrong

3.

I will deal the ripple after i'll manage to do this

Answer

You seem to be confused with the time and frequency domain. Your code for an LP filter is the inverse Fourier transform of an Ideal low-pass filter. So by getting a frequency response for the function you wrote you'll get an ideal LP response you want.

However, there are two mistakes:

1) Your frequency for the ideal LP filter needs to be $[-\frac{m-1}{2} : \frac{m-1}{2}]$

2) You need to change the value for the 0-frequency only if M is an odd number.

Now, you can use the function freqz to get your freqency response, so your final code would be like this:

%[y,Fs]=wavread('sound.wav');

%m = length(y); %30586
m = 30586;
Fc = 500;               %% frequency cut
n = -(m-1)/2 : (m-1)/2;            %% kernel length / samples
%Fs = Fs;    %22000      %% sample frequency / (nyquist - %2*m)
Fs = 22000;
Tc = 2*pi*(Fc/Fs);      %% theta c
Hm = sin(Tc*n) ./ (pi*n);
if (rem(m, 2) == 1)
    index = (m+1)/2;

    Hm(index) = Tc/pi;
end

[H w] = freqz(Hm, 1, 1024); % 1024 - or how many points you want for your Fourier Transform

% Now your H is your frequency response, and w is digital frequency in [0,pi].
% You can now plot the amplitude response:
figure,plot(w, abs(H));

And you have your Low-Pass filter. However, I'd like to suggest using the many ways to create filters in Matlab, they were made for a reason, and allow for much easier filter design :)

---

One more remark, I have to suggest that you know the difference between analog and digital frequency, as they have different notation. Analog frequencies you write in uppercase, e.g. $F, \Omega$, and digital frequencies in lower-case, $f, \omega$. Here's a version of the code which is a little prettier:

M = 30586;
Fc = 500;
Fs = 22000;
wc = 2*pi*Fc/Fs;

n = -(M-1)/2:(M-1)/2;    
h = (sin(n*wc))./(n*pi);
if ( M-2*fix(M/2) ) > 0

        index = (M+1)/2;
        h(index) = wc/pi;
        display(index);
end

[H w] = freqz(h, 1, 1024);
plot(w, abs(H));

Cheers :)

EDIT:

When you use freqz, the digital frequency w can be somewhat hard to control, because it's not created in the same way, or more clearly, the points are not selected like a normal frequency we create with, for example, F = 0 : Fs/(M-1) : Fs, even though it should be the "same" frequency. The result is that it's not really possible to get the output filtered signal with Y = H .* X, and y = ifft(Y).

However, I would suggest using the function filter(b, a, x) for filtering your signal, which works good with freqz. You can proceed to do just that:

% ------
% Create a random signal x
n = 0:M-1;
F = -Fs/2:Fs/(M-1):Fs/2;
Fx = [2000 7000 10000];
x = cos(2*pi*Fx(1)/Fs*n) + cos(2*pi*Fx(2)/Fs*n) + cos(2*pi*Fx(3)/Fs*n);
X = abs(fftshift(fft(x, Nfft))); % fftshift because we need domain [-Fs/2, Fs/2]
% ---------


% -------
% Now filter the signal
y = filter(H, 1, x);
Y = abs(fftshift(fft(y)));
% -------

That's it :)