Implementing FIR Bandpass Filter for 1-3Hz range

Hi

I am working on a bandpass filter for a converted analog signal from a motion sensor. My desired output is in the 1-3Hz range. The method I am using is an FIR filter that uses pointers and arrays. The coefficient values are generated using Matlab.

When tesing my code I do not observe any noticeable change in the output from the input.

I believe my error or mistake is in the pointers or arrays I am using.

Here is my code. Thanks for the help in advance!

float  B[101] = { -9.744167920817e-18, 0.002903896862104, -0.00579258969526, -0.01305549171185,
                  -0.0010171357064113, 0.0047101380987205, 0.001752792620389, -0.00195890640258,
                  0.007283946726629,  0.01468903320012, 0.0045910141444584, -0.007529394806697,
                  -0.00397974278453, 0.0010101384228858, -0.008734495973165, -0.01631292916464,
                  -0.003423027374945,   0.0113942560321, 0.0067810101635966, -0.0003219402137161,
                  0.0100857719628,  0.01798933486121, 0.001428125389994, -0.01687062505692,
                  -0.010388272310135, -1.745528550173e-18, -0.01128149112134, -0.01990166100632,
                  0.001830986824138,  0.02508643079021,  0.01531465858609, -0.0005330310083937,
                  0.01227059686935,  0.022551016970208, -0.007435810019908,  -0.0390060782669,
                  -0.02307284880661, 0.003081590662441, -0.01301006979061, -0.02762628158982,
                  0.01906859024139,  0.06959286537594,  0.03986831947415, -0.01257647738927,
                  0.013467321013046,  0.04491461919818, -0.06261887829046,   -0.215899673614,
                  -0.1471116759618,   0.1509465615928,    0.326929071269,   0.1509465615928,
                  -0.1471116759618,   -0.215899673614, -0.06261887829046,  0.04491461919818,
                  0.013467321013046, -0.01257647738927,  0.03986831947415,  0.06959286537594,
                  0.01906859024139, -0.02762628158982, -0.01301006979061, 0.003081590662441,
                  -0.02307284880661,  -0.0390060782669, -0.007435810019908,  0.022551016970208,
                  0.01227059686935, -0.0005330310083937,  0.01531465858609,  0.02508643079021,
                  0.001830986824138, -0.01990166100632, -0.01128149112134, -1.745528550173e-18,
                  -0.010388272310135, -0.01687062505692, 0.001428125389994,  0.01798933486121,
                  0.0100857719628, -0.0003219402137161, 0.0067810101635966,   0.0113942560321,
                  -0.003423027374945, -0.01631292916464, -0.008734495973165, 0.0010101384228858,
                  -0.00397974278453, -0.007529394806697, 0.0045910141444584,  0.01468903320012,
                  0.007283946726629, -0.00195890640258, 0.001752792620389, 0.0047101380987205,
                  -0.0010171357064113, -0.01305549171185, -0.00579258969526, 0.002903896862104,
                  -9.744167920817e-18
                };




int N = 101;

float analog[101], *p_vector = analog; //ADC
float fir_output[101];
int i = 0;
int currentPos = 0;

void setup() {
  Serial.begin(115200);
 

}

void loop() {

 

  
  float output, *p;

  *p_vector = analogRead(A0);   // store adc output value
  output = 0;
  p = p_vector;

  if (++p_vector > &analog[N]) {      // update pointer, wrap if necessary
    p_vector = analog;                // and store
  }
  for (i = 0; i <= N; i++) {         // do FIR
    output += *p-- * B[i];      // multiply and accumulate
    
    Serial.println(output);
    delay(5);

    if (p < &analog[0]) {   // check for pointer wrap around
      p = &analog[N];
    }
  }
 }

I suggest to start with a 10 tap filter.

When you multiply an analog input term (in the range of 0-1023) by -9.744167920817e-18 and add to other terms, do you expect to see much effect?

Please get rid of all the commented out junk and repost the code, so that you and other people can begin to make sense of the code.

Yes, sorry for the extra junk. by “10 tap” you are referring to the filter length and not order correct?

You have some errors in handling the pointers.

  if (++p_vector > &analog[N]) {

The > should be >=

      p = &analog[N];

This should be

      p = &analog[N-1];

The for loop limit is incorrect. It should be

  for (i = 0; i < N; i++) {

Your code will go round the loop 102 times instead of 101.

But you have much bigger problems than that. An FIR filter is designed for a specific sampling rate but you have made no attempt to maintain a specific sampling rate. The 101 multiply and accumulate operations will take many milliseconds (probably hundreds of ms) to execute on most Arduinos which will limit what your sampling rate can be. Also, you will need to ensure that your input signal is properly band-limited so that its frequency components are all less than half the sampling rate otherwise the results will be meaningless.
What sampling rate did you specify to Matlab?

My desired output is in the 1-3Hz range

Do you mean you have a filter that is 1-3Hz wide at some frequency or that you are filtering out everything except the signal from 1-3Hz? Either way, such a narrow filter is going to be very difficult to implement on an Arduino.

Pete

@el_supremo

Thank you for catching the errors.

So, reducing the order filter would shorten the time taken to execute?

The sampling rate I specified in Matlab was 8Hz.

I want to filter out frequencies except 1-3Hz.

Thanks again!

At a sampling rate of 8Hz you would have to filter out everything above 4Hz in the input before it gets to the FIR filter, otherwise you will get aliased frequencies.
You also need to fix your code so that it actually samples at a rate of 8Hz.

I've just spotted another problem in your code. The FIR filter assumes that the signal is centred at 0V whereas the ADC returns values from 0 to 1023 which represents a voltage range from zero to +5V. Make the change shown here to offset the value properly:

    output += (*p-- - 512) * B[i];

Have you offset your input signal so that it is within the range 0 to 5V? A negative voltage on pin A0 could destroy the pin or, in the worst case, the Arduino as a whole.

Pete

Okay, thank you for catching that error.

I did have an offset but, also tested a couple times without. Luckily, nothing was damaged.

For the low pass filter with a cut off of 4Hz, would the low pass filter need to be an analog filter from the sensor to the arduino?

And could I simply add a delay(125) to set the code to sample at 8Hz?

For the low pass filter with a cut off of 4Hz, would the low pass filter need to be an analog filter from the sensor to the arduino?

Yes.
It would also need to be a much better filter than a first order RC filter.

And could I simply add a delay(125) to set the code to sample at 8Hz?

In theory yes but only if all the other instructions took zero time to execute. You are better off with a timer calling an ISR and doing the sampling in the ISR.

Yes. I don't really have the time or space to build a high-order filter. Is there an opamp ic like the LM741 that would do the trick?

Thank you. I think I will got with the ISR approach.

EDIT. Thinking about it further, I don't think an ic would help the sharpness in drop off. Hopefully, I am wrong?

I don't really have the time or space to build a high-order filter.

Then you have to rethink this problem completely.

It is absolutely essential that no frequency higher than half the sampling frequency gets to the input, otherwise you have aliasing, a severe problem pointed out earlier.

With no input filter at all, you must sample at twice the highest frequency in the input, then decimate and filter.

theswellylife:
Yes. I don't really have the time or space to build a high-order filter. Is there an opamp ic like the LM741 that would do the trick?

Thank you. I think I will got with the ISR approach.

EDIT. Thinking about it further, I don't think an ic would help the sharpness in drop off. Hopefully, I am wrong?

You can implement a 2nd order Sallen-Key filter with a single op amp that would roll off twice as quickly as an RC filter. The op amp would also allow you to set the input level to the A/D converter.

Without a very high order filter, you're going to have to sample at a higher rate where the filter has sufficiently attenuated the high frequency content of your source signal.

IIR filters are more computationally efficient than FIR filters, so that would be preferred unless you need linear phase or some other characteristic of the FIR filter. In either case, it will likely be more computationally efficient to implement the digital filter in two or more stages with decimation.

It seems to me that the FIR part of the project is well in hand, with the Matlab design.

For stripping the high-frequencies, I would not bother with anything more than a single-pole filter. Just crank up the sample rate until you can be absolutely sure that the nyquist frequency is completely covered by the analog filter. Even though you are looking for very low frequencies, there's not much penalty for using a high sample frequency. Either just use a large FIR filter to directly find the frequencies of interest or use a digital filter (of any sort) to strip the frequencies above the final sample rate's nyquist frequency before decimating.

My favourite reference on DSP: The Scientist and Engineer's Guide to
Digital Signal Processing
By Steven W. Smith, Ph.D.

My favourite FIR filter design tool: TFilter

Could you tell us the specifications of the filter you wish to implement? Both analogue and digital. Why?

Bandwidth - 1-3 Hz - OK. I presume -3dB points.

attenuation rate outside the passband in dB/octave.

Stopband attenuation.

Acceptable ripple.

Phase response.

Given that I may be able to suggest a solution. An analogue elliptic may be a good approach. What type does your digital filter implement? - given it's high order perhaps a Bessel....

is this an X-Y problem? ie are you asking a detail question on just a part of your perceived solution - which may not be optimal.

What does your overall project hope to achieve? There may be a better/easier way.

Allan

filter. Even though you are looking for very low frequencies, there's not much penalty for using a high sample frequency. Either just use a large FIR filter to directly find the frequencies of interest or use a digit

Well I would say their is a great penalty in using a higher sample rate, you have to compute the filter coefficients between each sample, with a lot like here it will define the maximum sample rate you can use

there's not much penalty for using a high sample frequency.

There is an enormous penalty for high sample rate on a 16 MHz Arduino. As pointed out above, that leaves little time for computing the filter, especially something on the order of 100 taps.

100 FP multiplications alone takes around 1 ms, which means a sample rate of 1 kHz max even if the ADC sampling completely overlaps the filter calculations.

The OP is looking for 1-3Hz. In this context, 100Hz might be a high sample rate. It certainly doesn't need kHz.

It's still an XY problem though.

It certainly doesn't need kHz.

You cannot possibly know that, because you do not know the frequency components of the input signal.

With no input filter, the sample frequency must be at least twice the highest frequency component.

allanhurst:
Could you tell us the specifications of the filter you wish to implement? Both analogue and digital. Why?

Bandwidth - 1-3 Hz - OK. I presume -3dB points.

attenuation rate outside the passband in dB/octave.

Stopband attenuation.

Acceptable ripple.

Phase response.

So I am building a contactless heart monitor, using a parallax x-band motion detector. The detector at its highest sensitivity can detect movement anywhere from the beating of a heart to a person walking. So I am trying to to just detect the heart beating which is typically between 1-2Hz but, I've expanded that window that width to 1-3Hz in case of a very high heart rate.

Attenuation rate would need to be high. I don't want any signals above 4Hz

Stopband attenuation -10dB or less

I am not necessarily worried about ripple.

I am not really worried about phase response either though, I am reading values in real time so maybe I should be.

EDIT. I am using Matlab's Filter Design and Analysis Tool. To design a Bandpass Equiripple filter

I don't want any signals above 4Hz

Stopband attenuation -10dB or less

These are completely inconsistent specifications.

@jremington

Sorry wasn’t paying attention when I wrote the stopband attenuation I meant -100dB