1
$\begingroup$

I'm trying to improve my understanding of FIR filters. As an experiment, I've manually created an FIR filter, whose coefficients follow exactly one period of a sine wave. I'm wondering what is the peak gain frequency of such a filter.

For example, I'm filtering a 44.1 kHz audio signal with filter coefficients corresponding to a 440 Hz sine wave1. Intuitively I would expect that the convolution reinforces this frequency. However, I noticed that the frequencies of slightly lower frequencies (roughly a few semitones, ~365 Hz) actually have a stronger gain. This is also confirmed by running the filter coefficients through scipy.signal.freqz (the red line corresponds to 440 Hz, the gray lines indicate frequencies of neighboring semitones):

frequency response

Most likely I've encountered a basic rule of digital filter design that I'm not aware of. I'm wondering:

  • Why is the peak shifted w.r.t. the frequency that is used in the convolution?
  • Is there a formula that describes the relationship between the two, i.e., is it possible to compute the frequency that the sine convolution should have so that the frequency response peaks at 440 Hz?

1 With a sampling frequency of 44.1 kHz, the period of a 440 Hz sine has a fractional length 100.23 samples, which I'm rounding down to 100 filter coefficients. Could this imperfection play a role in the effect I'm seeing?

$\endgroup$
4
  • $\begingroup$ yep, the effect of rectangularly windowing the sine with a rectangle whose length is not a multiple of the sine's period. Try a length that actually fits a whole number of periods and see the analysis become a line spectrum. Or – use a non-rectangular window and suppress the "wideness" of these lobes and get a spectrum that is more closely centered around the target frequency. $\endgroup$ Oct 5, 2019 at 16:36
  • $\begingroup$ @MarcusMüller Hm, I tried it with 441 Hz, so that one sine period matches exactly a filter length of 100, but I'm seeing a similar shift in peak frequency. $\endgroup$
    – bluenote10
    Oct 5, 2019 at 17:04
  • $\begingroup$ did you also set the number of points for freqz to a multiple of the filter length? $\endgroup$ Oct 5, 2019 at 17:06
  • $\begingroup$ @MarcusMüller Yes, that is set to the sampling frequency of 44100. $\endgroup$
    – bluenote10
    Oct 5, 2019 at 18:03

1 Answer 1

3
$\begingroup$

If your signal was a windowed complex exponential, then the maximum of the magnitude of its discrete-time Fourier transform (DTFT) would equal the frequency of the complex exponential:

$$x[n]=\begin{cases}e^{jn\omega_0},\quad&n\in[0,N-1]\\0,\quad& \textrm{otherwise}\end{cases}\tag{1}$$

The DTFT of $(1)$ is easily computed from the convolution of the DTFT of a rectangular window $w[n]=1$ for $n\in[0,N-1]$ and zero otherwise, and the DTFT of $e^{jn\omega_0}$:

$$\textrm{DTFT}\{w[n]\}=W(e^{j\omega})=e^{-j(N-1)\omega/2}\frac{\sin\left(N\omega/2\right)}{\sin(\omega/2)}\tag{2}$$

$$\textrm{DTFT}\{e^{jn\omega_0}\}=2\pi\delta(\omega-\omega_0),\qquad |\omega|<\pi,\quad |\omega_0|<\pi\tag{3}$$

Note that the DTFT in $(3)$ is of course $2\pi$-periodic.

With $(2)$ and $(3)$ we can write the DTFT of $(1)$ as follows:

$$\begin{align}X(e^{j\omega})&=W(e^{j\omega})\star\delta(\omega-\omega_0)\\&=W\left(e^{j(\omega-\omega_0)}\right)\\&=e^{-j(N-1)(\omega-\omega_0)/2}\frac{\sin\left(N(\omega-\omega_0)/2\right)}{\sin((\omega-\omega_0)/2)}\tag{4}\end{align}$$

The magnitude of $(4)$ clearly has a maximum at $\omega=\omega_0$ $(+2\pi k)$

However, for a windowed sinusoidal signal we get two contributions of the form $(2)$ centered at $\pm\omega_0$ because $\sin(n\omega_0)=(e^{jn\omega_0}-e^{-jn\omega_0})/2j$:

$$y[n]=\begin{cases}\sin(n\omega_0),\quad&n\in[0,N-1]\\0,\quad& \textrm{otherwise}\end{cases}\tag{5}$$

$$Y(e^{j\omega})=\frac{1}{2j}\left[e^{-j(N-1)(\omega-\omega_0)/2}\frac{\sin\left(N(\omega-\omega_0)/2\right)}{\sin((\omega-\omega_0)/2)}-\\\quad e^{-j(N-1)(\omega+\omega_0)/2}\frac{\sin\left(N(\omega+\omega_0)/2\right)}{\sin((\omega+\omega_0)/2)}\right]\tag{6}$$

I don't think that there is an analytical formula for the location of the maximum of the magnitude of $(6)$, but it is clear that due to the combination of the two terms in $(6)$, the maximum shifts away from $\omega_0$. That shift can become especially large for frequencies either close to DC or close to Nyquist. It is smallest for frequencies around half the Nyquist frequency (i.e., $f_s/4$). For small frequencies $\omega_0$, the actual frequency where the maximum occurs is smaller than $\omega_0$, whereas for large frequencies it is greater than $\omega_0$.

$\endgroup$
0

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.