Trying to apply Parks-McClellan optimal FIR filter design

Question

My goal is to apply the Remez exchange algorithm in R to design an FIR band pass filter for a specific frequency window. The question is not supposed to be referring to R of course but about how to choose the parameters correctly and why the result looks like as it does.

The toy signal is an additive combinations of four sine waves each adding a single frequency component to the final signal. The frequencies are 1/20, 1/300, 1/5000 and 1/7000. Its FFT plot is shown in figure (s).

Now what I apparently didn't get yet is how I would have to choose the parameters for the filter design - being the pass band and the stop band cutoff points.

The 1st filter uses the cutoff points:

• 0
• 1/10 - 2/1000
• 1/10 - 1/1000
• 1/10 + 1/1000
• 1/10 + 2/1000
• 1

With 1 being the Nyquist-frequency and the valid cutoff point interval being (0,1).

As you can see in figure (r1) this results in an attenuation of all frequency components except for 1/5000. (I would assume that also the frequency component 1/5000 is actually attenuated but relatively less then the other components)

So my first set of questions would be - why is that so? How could I calculate f.x. the frequency located in the center of the pass band? What is the quantitative relation between cutoff points on an interval (0,1) and the pass band?

The 2nd filter uses the cutoff points:

• 0
• 1/10 - 2/500
• 1/10 - 1/500
• 1/10 + 1/500
• 1/10 + 2/500
• 1

And the result is an (relative) attenuation of all components except for 1/300. What I am wondering about in this case is why the pass band is apparently completely different from the above filter even though I just made the filter's pass band wider. I would assume at least that 1/5000 is still located within the pass band and not attenuated. A possible explanation might be that 1/300 is now included within the pass band and relative to 1/5000 less attenuated making 1/5000 look just smaller.

I would like to add that I am - as you can see - quite new to this subject and actually from a different area and trying my best to understand what's going on.

---

The canonical question would be: How do I have to choose the cutoff points in case of a possible interval (0,1=Nyquist-Freq) to attenuate all frequency components except for 1/N and its close neighbourhood?

---

library(signal)
par(mfrow=c(5,1))

# the toy signal
ch <- sin(2*pi*1:100000/floor(100000/20)) + 
  sin(2*pi*1:100000/floor(100000/300)) +
  sin(2*pi*1:100000/floor(100000/5000)) +
  sin(2*pi*1:100000/floor(100000/7000))

# (s)
barplot(abs(fft(ch))[1:10000])

filter <- function(c0,d1,d2) {

  # the filter design with cutoffs specified
  fir <- remez(n=1000,
    f=c(0,c0-d2,c0-d1,c0+d1,c0+d2,1),
    a=c(0,0,1,1,0,0))

  # (fN)
  plot(freq$f,abs(freq$h),type="l")
  freq <- freqz(fir,n=100000)
  y <- signal::filter(as.vector(fir), 1, x=ch)

  # (sN)
  barplot(abs(fft(y))[1:10000])
}

filter(1/10,1/1000,2/1000)
filter(1/10,1/500,2/500)

Answer 1

There are a few problems with your designs. First, the desired bandwidths and transition bands (between passbands and stopbands) are very narrow. Second, some of the center frequencies are very low, especially the ones at 20/100000 and 300/100000. And finally, your chosen filter length is very large. I do not know the implementation of the Remez exchange algorithm that you use, but I've seen several implementations that run into serious numerical problems for such filter lengths. You can see in your designs that around Nyquist there are some small peaks in the filter response. These are due to numerical problems.

So you should try to make the passband and the transition bands a bit wider and reduce the filter length to see if you can get a reasonable design. From there you could try to make the specifications a bit more stringent and increase the filter order and see how far you can go.

In Octave, the following specification gives a reasonable result (center frequency: 0.1, filter length: 300, '1' corresponds to Nyquist):

h=remez(300,[0,.07,.085,.115,.13,1],[0,0,1,1,0,0]);