November – Digital Filtering

Let's construct a simple digital band-pass filter based in the sinc filter.

As everyone know, the sinc function is the impulse response of a low-pass filter. That is, the convolution of a signal $S(t)$ and the sinc function $t \mapsto \frac{\sin(2 \pi B_L t)}{\pi t}$ will produce a low-pass filtered signal $S(t)$ without frequencies higher than $B_L$. Similarly, the function $t \mapsto \delta(t) - \frac{\sin(2 \pi B_H t)}{\pi t},$ where $\delta(t)$ represents the Dirac delta function, is the frequency response of the ideal high-pass filter of frequency $B_H$.

Heuristically, we can right away derive the discrete versions of the impulse responses:

$\displaystyle \mathbb{Z} \ni i \mapsto \frac{\sin(2 \pi \frac{B_L}{f_c} i)}{\pi i},
$

for the low-pass filter, and

$\displaystyle i \mapsto \delta[i] - \frac{\sin(2 \pi \frac{B_H}{f_c} i)}{\pi i},
$

for the high-pass filter. The variable $f_c$ is the sampling frequency and $i \mapsto \delta[i] := \delta_{0i} $ is the Kronecker delta.

Now, having a discrete signal at hand, band-pass filtering is (in the simplest approach) just a matter of discrete convolution of the signal with the functions above. The following octave code does the job.

##Sinc band-pass filter. f0 = B_L/fs, f1 = B_H/fs

function filtered = sincfilter(signal, f0, f1)
  if(mod(length(signal), 2) != 0) #Check if the signal length is even or odd.
    signal = signal(1 : length(signal) - 1);
  end
  
  M = 10000; #Increase this to increase accuracy.
  sincF = zeros(1, M);
  for m =  -M/2 + 1 : 1 : M/2 
    if(!(m  == 0))
      sincF(m + M/2) = sin(2*pi*f1*m)./(pi*m);
    else
      sincF(m + M/2) = 2*f1;
    end
  end

  sincH = zeros(1,M);
  for m =  -M/2 + 1 : 1 : M/2 
    if(!(m  == 0))
      sincH(m + M/2) = -sin(2*pi*f0*m)./(pi*m);
    else
      sincH(m + M/2) = -2*f0 + 1;
    end
  end
  
  ##Plot stuff.
  figure(1)
  plot(sincF)
  figure(2)
  plot(signal)

  tic
  filtered = conv(signal, sincF, "same");
  filtered = conv(filtered, sincH, "same");
  toc

  figure(3)
  plot(filtered)
end

Figure: Signal of form $\cos(20 \pi t) + \sin(2 \pi t)$
Image signal

Figure: The same signal with the first term filtered away
Image filtered

For more accuracy and effectiveness, one should use windowed or recursive filters. You can check in the references for more sophisticated methods!

References: