Apply a digital filter (direct form II transposed)

Digital filters are used to apply a mathematical operation to a sampled signal. One of the common formulations is the "direct form II transposed" which can represent both infinite impulse response (IIR) and finite impulse response (FIR) filters, as well as being more numerically stable than other forms. ^[1]

Task

Filter a signal using an order 3 lowpass butterworth filter. The coefficients for the filter are a=[1.00000000, -2.77555756e-16, 3.33333333e-01, -1.85037171e-17] and b = [0.16666667, 0.5, 0.5, 0.16666667]

The signal the needs filtering is the following vector: [-0.917843918645,0.141984778794,1.20536903482,0.190286794412,-0.662370894973,-1.00700480494,-0.404707073677,0.800482325044,0.743500089861,1.01090520172,0.741527555207,0.277841675195,0.400833448236,-0.2085993586,-0.172842103641,-0.134316096293,0.0259303398477,0.490105989562,0.549391221511,0.9047198589]

C++

This uses the C++11 method of initializing vectors. In g++, use the -std=c++0x compiler switch.

<lang cpp>#include <vector>

include <iostream>

using namespace std;

void Filter(const vector<float> &b, const vector<float> &a, const vector<float> &in, vector<float> &out) {

out.resize(0); out.resize(in.size());

for(int i=0; i < in.size(); i++) { float tmp = 0.; int j=0; out[i] = 0.f; for(j=0; j < b.size(); j++) { if(i - j < 0) continue; tmp += b[j] * in[i-j]; }

for(j=1; j < a.size(); j++) { if(i - j < 0) continue; tmp -= a[j]*out[i-j]; }

tmp /= a[0]; out[i] = tmp; } }

int main() { vector<float> sig = {-0.917843918645,0.141984778794,1.20536903482,0.190286794412,-0.662370894973,-1.00700480494,\ -0.404707073677,0.800482325044,0.743500089861,1.01090520172,0.741527555207,\ 0.277841675195,0.400833448236,-0.2085993586,-0.172842103641,-0.134316096293,\ 0.0259303398477,0.490105989562,0.549391221511,0.9047198589};

vector<float> a = {1.00000000, -2.77555756e-16, 3.33333333e-01, -1.85037171e-17}; vector<float> b = {0.16666667, 0.5, 0.5, 0.16666667};

vector<float> result; Filter(b, a, sig, result);

for(size_t i=0;i<result.size();i++) cout << result[i] << ","; cout << endl;

return 0; }</lang>

Output:

-0.152974,-0.435258,-0.136043,0.697503,0.656445,-0.435483,-1.08924,-0.537677,0.51705,1.05225,0.961854,0.69569,0.424356,0.196262,-0.0278351,-0.211722,-0.174746,0.0692584,0.385446,0.651771,

Python

<lang python>#!/bin/python from __future__ import print_function from scipy import signal import matplotlib.pyplot as plt

if __name__=="__main__": sig = [-0.917843918645,0.141984778794,1.20536903482,0.190286794412,-0.662370894973,-1.00700480494, -0.404707073677,0.800482325044,0.743500089861,1.01090520172,0.741527555207, 0.277841675195,0.400833448236,-0.2085993586,-0.172842103641,-0.134316096293, 0.0259303398477,0.490105989562,0.549391221511,0.9047198589]

#Create an order 3 lowpass butterworth filter #Generated using b, a = signal.butter(3, 0.5) a = [1.00000000, -2.77555756e-16, 3.33333333e-01, -1.85037171e-17] b = [0.16666667, 0.5, 0.5, 0.16666667]

#Apply the filter to xn. Use lfilter_zi to choose the initial condition of the filter: filt = signal.lfilter(b, a, sig) print (filt)

plt.plot(sig, 'b') plt.plot(filt, 'r--') plt.show()</lang>

Output:

[-0.15297399 -0.43525783 -0.1360434   0.69750333  0.65644469 -0.43548245
 -1.08923946 -0.53767655  0.51704999  1.05224975  0.9618543   0.69569009
  0.4243563   0.19626223 -0.02783512 -0.21172192 -0.17474556  0.06925841
  0.38544587  0.65177084]

References

↑ [1]

[1] [1]

[1]