## SciPy FIR Filter

#### Mbed OS 2 and Mbed OS 5

This is the handbook for Mbed OS 2. If you’re working with Mbed OS 5, please see the new handbook.

In this paragraph we see how to design a FIR filter using SciPy, NumPy and MatplotLib.

You can install SciPy, NumPy and MatplotLib as normal Python libraries, otherwise you can simply install the Enthought Python distribution.

The Enthought Python distribution is packaging for you a Python interpreter and all the above libraries in a single convenient installer: http://enthought.com/products/epd_free.php

## FIR Filter Design¶

#### Python FIR Filter Design

```from numpy import sin, arange, pi
from scipy.signal import lfilter, firwin
from pylab import figure, plot, grid, show

#------------------------------------------------
# Create a signal for demonstration.
#------------------------------------------------
# 320 samples of (1000Hz + 15000 Hz) at 48 kHz
sample_rate = 48000.
nsamples = 320

F_1KHz = 1000.
A_1KHz = 1.0

F_15KHz = 15000.
A_15KHz = 0.5

t = arange(nsamples) / sample_rate
signal = A_1KHz * sin(2*pi*F_1KHz*t) + A_15KHz*sin(2*pi*F_15KHz*t)

#------------------------------------------------
# Create a FIR filter and apply it to signal.
#------------------------------------------------
# The Nyquist rate of the signal.
nyq_rate = sample_rate / 2.

# The cutoff frequency of the filter: 6KHz
cutoff_hz = 6000.0

# Length of the filter (number of coefficients, i.e. the filter order + 1)
numtaps = 29

# Use firwin to create a lowpass FIR filter
fir_coeff = firwin(numtaps, cutoff_hz/nyq_rate)

# Use lfilter to filter the signal with the FIR filter
filtered_signal = lfilter(fir_coeff, 1.0, signal)

#------------------------------------------------
# Plot the original and filtered signals.
#------------------------------------------------

# The first N-1 samples are "corrupted" by the initial conditions
warmup = numtaps - 1

# The phase delay of the filtered signal
delay = (warmup / 2) / sample_rate

figure(1)
# Plot the original signal
plot(t, signal)

# Plot the filtered signal, shifted to compensate for the phase delay
plot(t-delay, filtered_signal, 'r-')

# Plot just the "good" part of the filtered signal.  The first N-1
# samples are "corrupted" by the initial conditions.
plot(t[warmup:]-delay, filtered_signal[warmup:], 'g', linewidth=4)

grid(True)

show()

#------------------------------------------------
# Print values
#------------------------------------------------
def print_values(label, values):
var = "float32_t %s[%d]" % (label, len(values))
print "%-30s = {%s}" % (var, ', '.join(["%+.10f" % x for x in values]))

print_values('signal', signal)
print_values('fir_coeff', fir_coeff)
print_values('filtered_signal', filtered_signal)
```

#### coefficients

```float32_t fir_coeff[29]        = {-0.0018225230, -0.0015879294, +0.0000000000, +0.0036977508, +0.0080754303, +0.0085302217, -0.0000000000, -0.0173976984, -0.0341458607, -0.0333591565, +0.0000000000, +0.0676308395, +0.1522061835, +0.2229246956, +0.2504960933, +0.2229246956, +0.1522061835, +0.0676308395, +0.0000000000, -0.0333591565, -0.0341458607, -0.0173976984, -0.0000000000, +0.0085302217, +0.0080754303, +0.0036977508, +0.0000000000, -0.0015879294, -0.0018225230}
```

## mbed-dsp library¶

#### Import librarymbed-dsp

CMSIS DSP library

### CMSIS DSP API¶

#### Import programcmsis_dsp_fir - main.cpp

```00001 #include "arm_math.h"
00002 #include "math_helper.h"
00003 #include <stdio.h>
00004
00005 #define BLOCK_SIZE          32
00006 #define NUM_BLOCKS          10
00007
00008 #define TEST_LENGTH_SAMPLES (BLOCK_SIZE * NUM_BLOCKS)
00009
00010 #define SNR_THRESHOLD_F32   140.0f
00011 #define NUM_TAPS            29
00012
00013 /* -------------------------------------------------------------------
00014  * The input signal and reference output (computed with MATLAB)
00015  * are defined externally in arm_fir_lpf_data.c.
00016  * ------------------------------------------------------------------- */
00017 extern float32_t testInput_f32_1kHz_15kHz[TEST_LENGTH_SAMPLES];
00018 extern float32_t refOutput[TEST_LENGTH_SAMPLES];
00019
00020 /* -------------------------------------------------------------------
00021  * Declare State buffer of size (numTaps + blockSize - 1)
00022  * ------------------------------------------------------------------- */
00023 static float32_t firStateF32[BLOCK_SIZE + NUM_TAPS - 1];
00024
00025 /* ----------------------------------------------------------------------
00026  * FIR Coefficients buffer generated using fir1() MATLAB function.
00027  * fir1(28, 6/24)
00028  * ------------------------------------------------------------------- */
00029 const float32_t firCoeffs32[NUM_TAPS] = {
00030     -0.0018225230f, -0.0015879294f, +0.0000000000f, +0.0036977508f, +0.0080754303f,
00031     +0.0085302217f, -0.0000000000f, -0.0173976984f, -0.0341458607f, -0.0333591565f,
00032     +0.0000000000f, +0.0676308395f, +0.1522061835f, +0.2229246956f, +0.2504960933f,
00033     +0.2229246956f, +0.1522061835f, +0.0676308395f, +0.0000000000f, -0.0333591565f,
00034     -0.0341458607f, -0.0173976984f, -0.0000000000f, +0.0085302217f, +0.0080754303f,
00035     +0.0036977508f, +0.0000000000f, -0.0015879294f, -0.0018225230f
00036 };
00037
00038 /* ----------------------------------------------------------------------
00039  * FIR LPF Example
00040  * ------------------------------------------------------------------- */
00041 int32_t main(void) {
00042     /* Call FIR init function to initialize the instance structure. */
00043     arm_fir_instance_f32 S;
00044     arm_fir_init_f32(&S, NUM_TAPS, (float32_t *)&firCoeffs32[0], &firStateF32[0], BLOCK_SIZE);
00045
00046     /* ----------------------------------------------------------------------
00047      * Call the FIR process function for every blockSize samples
00048      * ------------------------------------------------------------------- */
00049     for (uint32_t i=0; i < NUM_BLOCKS; i++) {
00050         float32_t* signal =  testInput_f32_1kHz_15kHz + (i * BLOCK_SIZE);
00051         arm_fir_f32(&S, signal, signal, BLOCK_SIZE);
00052     }
00053
00054     /* ----------------------------------------------------------------------
00055      * Compare the generated output against the reference output computed
00056      * in MATLAB.
00057      * ------------------------------------------------------------------- */
00058     float32_t snr = arm_snr_f32(refOutput, testInput_f32_1kHz_15kHz, TEST_LENGTH_SAMPLES);
00059     printf("snr: %f\n\r", snr);
00060     if (snr < SNR_THRESHOLD_F32) {
00061         printf("Failed\n\r");
00062     } else {
00063         printf("Success\n\r");
00064     }
00065 }
```

### mbed DSP API¶

#### Import programdsp_fir - main.cpp

```00001 #include "mbed.h"
00002 #include "dsp.h"
00003
00004 #define BLOCK_SIZE              (32)
00005 #define NUM_BLOCKS              (10)
00006 #define TEST_LENGTH_SAMPLES     (BLOCK_SIZE * NUM_BLOCKS)
00007
00008 #define SAMPLE_RATE             (48000)
00009
00010 #define SNR_THRESHOLD_F32       (50.0f)
00011
00012 float32_t expected_output[TEST_LENGTH_SAMPLES];
00013 float32_t          output[TEST_LENGTH_SAMPLES];
00014
00015 /* FIR Coefficients buffer generated using fir1() MATLAB function: fir1(28, 6/24) */
00016 #define NUM_TAPS            29
00017 const float32_t firCoeffs32[NUM_TAPS] = {
00018     -0.0018225230f, -0.0015879294f, +0.0000000000f, +0.0036977508f, +0.0080754303f,
00019     +0.0085302217f, -0.0000000000f, -0.0173976984f, -0.0341458607f, -0.0333591565f,
00020     +0.0000000000f, +0.0676308395f, +0.1522061835f, +0.2229246956f, +0.2504960933f,
00021     +0.2229246956f, +0.1522061835f, +0.0676308395f, +0.0000000000f, -0.0333591565f,
00022     -0.0341458607f, -0.0173976984f, -0.0000000000f, +0.0085302217f, +0.0080754303f,
00023     +0.0036977508f, +0.0000000000f, -0.0015879294f, -0.0018225230f
00024 };
00025 #define WARMUP    (NUM_TAPS-1)
00026 #define DELAY     (WARMUP/2)
00027
00028 int main() {
00029     Sine_f32 sine_1KHz(  1000, SAMPLE_RATE, 1.0);
00030     Sine_f32 sine_15KHz(15000, SAMPLE_RATE, 0.5);
00031     FIR_f32<NUM_TAPS> fir(firCoeffs32);
00032
00033     float32_t buffer_a[BLOCK_SIZE];
00034     float32_t buffer_b[BLOCK_SIZE];
00035     for (float32_t *sgn=output; sgn<(output+TEST_LENGTH_SAMPLES); sgn += BLOCK_SIZE) {
00036         sine_1KHz.generate(buffer_a);           // Generate a 1KHz sine wave
00037         sine_15KHz.process(buffer_a, buffer_b); // Add a 15KHz sine wave
00038         fir.process(buffer_b, sgn);             // FIR low pass filter: 6KHz cutoff
00039     }
00040
00041     sine_1KHz.reset();
00042     for (float32_t *sgn=expected_output; sgn<(expected_output+TEST_LENGTH_SAMPLES); sgn += BLOCK_SIZE) {
00043         sine_1KHz.generate(sgn);        // Generate a 1KHz sine wave
00044     }
00045
00046     float snr = arm_snr_f32(&expected_output[DELAY-1], &output[WARMUP-1], TEST_LENGTH_SAMPLES-WARMUP);
00047     printf("snr: %f\n\r", snr);
00048     if (snr < SNR_THRESHOLD_F32) {
00049         printf("Failed\n\r");
00050     } else {
00051         printf("Success\n\r");
00052     }
00053 }
```