16,003,243 members
Articles / General Programming / Algorithms
Alternative
Tip/Trick

Lowpass, Highpass, and Bandpass Butterworth Filters in C#

Rate me:
9 Jun 2019CPOL1 min read 42K   2.3K   13   7
Implementations of lowpass, highpass, and bandpass nth-order Butterworth filters in C#, with design documentation.

Introduction

This is a C# implementation of digital lowpass, highpass, and bandpass Butterworth filters of arbitrary order (n cascaded 2-pole sections).

Background

A Word document giving the filter design via bilinear z-transformation is included.

The design starts with a continuous-time lowpass Butterworth filter and uses the bilinear z-transform to derive the coefficients of the equivalent digital filter. Then a substitution of variable `s=1/s` transforms the lowpass filter into a highpass filter, and the bilinear z-transform is used to derive the coefficients of the equivalent highpass digital filter.

In the code, the lowpass and highpass filters are implemented according to the coefficient derivation from the design document. A bandpass filter is implemented as a cascade of a lowpass and a highpass filter.

The code includes an implementation of an nth-order FIR filter for the zero (numerator) polynomials and an implementation of an nth-order IIR filter for the pole (denominator) polynomials.

Using the Code

The code is divided into five files:

1. LowpassFilterButterworthImplementation.cs
2. HighpassFilterButterworthImplementation.cs
3. BandpassFilterButterworthImplementation.cs
4. FIRFilterImplementation.cs
5. IIRFilterImplementation.cs

To use the `LowpassFilterButterworthImplementation`, `HighpassFilterButterworthImplementation`, or `BandpassFilterButterworthImplementation`, construct the object with the desired cutoff frequency or frequencies, number of 2-pole sections, and the sample frequency (Fs). Then call the `compute()` method on the object repeatedly, passing in an input sample and receiving a filtered sample as the return value.

C#
```using DSP;
using System;

namespace DSP
{
public class LowpassFilterButterworthImplementation
{
protected LowpassFilterButterworthSection[] section;

public LowpassFilterButterworthImplementation
(double cutoffFrequencyHz, int numSections, double Fs)
{
this.section = new LowpassFilterButterworthSection[numSections];
for (int i = 0; i < numSections; i++)
{
this.section[i] = new LowpassFilterButterworthSection
(cutoffFrequencyHz, i + 1, numSections * 2, Fs);
}
}
public double compute(double input)
{
double output = input;
for (int i = 0; i < this.section.Length; i++)
{
output = this.section[i].compute(output);
}
return output;
}
}
}

public class LowpassFilterButterworthSection
{
protected FIRFilterImplementation firFilter = new FIRFilterImplementation(3);
protected IIRFilterImplementation iirFilter = new IIRFilterImplementation(2);

protected double[] a = new double[3];
protected double[] b = new double[2];
protected double gain;

public LowpassFilterButterworthSection
(double cutoffFrequencyHz, double k, double n, double Fs)
{
// compute the fixed filter coefficients
double omegac = 2.0 * Fs * Math.Tan(Math.PI * cutoffFrequencyHz / Fs);
double zeta = -Math.Cos(Math.PI * (2.0 * k + n - 1.0) / (2.0 * n));

// fir section
this.a[0] = omegac * omegac;
this.a[1] = 2.0 * omegac * omegac;
this.a[2] = omegac * omegac;

//iir section
//normalize coefficients so that b0 = 1,
//and higher-order coefficients are scaled and negated
double b0 = (4.0 * Fs * Fs) + (4.0 * Fs * zeta * omegac) + (omegac * omegac);
this.b[0] = ((2.0 * omegac * omegac) - (8.0 * Fs * Fs)) / (-b0);
this.b[1] = ((4.0 * Fs * Fs) -
(4.0 * Fs * zeta * omegac) + (omegac * omegac)) / (-b0);
this.gain = 1.0 / b0;
}

public double compute(double input)
{
// compute the result as the cascade of the fir and iir filters
return this.iirFilter.compute
(this.firFilter.compute(this.gain * input, this.a), this.b);
}
}```

C#
```using DSP;
using System;

namespace DSP
{
public class HighpassFilterButterworthImplementation
{
protected HighpassFilterButterworthSection[] section;

public HighpassFilterButterworthImplementation
(double cutoffFrequencyHz, int numSections, double Fs)
{
this.section = new HighpassFilterButterworthSection[numSections];
for (int i = 0; i < numSections; i++)
{
this.section[i] = new HighpassFilterButterworthSection
(cutoffFrequencyHz, i + 1, numSections * 2, Fs);
}
}
public double compute(double input)
{
double output = input;
for (int i = 0; i < this.section.Length; i++)
{
output = this.section[i].compute(output);
}
return output;
}
}
}

public class HighpassFilterButterworthSection
{
protected FIRFilterImplementation firFilter = new FIRFilterImplementation(3);
protected IIRFilterImplementation iirFilter = new IIRFilterImplementation(2);

protected double[] a = new double[3];
protected double[] b = new double[2];
protected double gain;

public HighpassFilterButterworthSection(double cutoffFrequencyHz, double k, double n, double Fs)
{
// pre-warp omegac and invert it
double omegac = 1.0 /(2.0 * Fs * Math.Tan(Math.PI * cutoffFrequencyHz / Fs));

// compute zeta
double zeta = -Math.Cos(Math.PI * (2.0 * k + n - 1.0) / (2.0 * n));

// fir section
this.a[0] = 4.0 * Fs * Fs;
this.a[1] = -8.0 * Fs * Fs;
this.a[2] = 4.0 * Fs * Fs;

//iir section
//normalize coefficients so that b0 = 1
//and higher-order coefficients are scaled and negated
double b0 = (4.0 * Fs * Fs) + (4.0 * Fs * zeta / omegac) + (1.0 / (omegac * omegac));
this.b[0] = ((2.0 / (omegac * omegac)) - (8.0 * Fs * Fs)) / (-b0);
this.b[1] = ((4.0 * Fs * Fs)
- (4.0 * Fs * zeta / omegac) + (1.0 / (omegac * omegac))) / (-b0);
this.gain = 1.0 / b0;
}

public double compute(double input)
{
// compute the result as the cascade of the fir and iir filters
return this.iirFilter.compute
(this.firFilter.compute(this.gain * input, this.a), this.b);
}
}```

C#
```using System;

namespace DSP
{
public class BandpassFilterButterworthImplementation
{
protected LowpassFilterButterworthImplementation lowpassFilter;
protected HighpassFilterButterworthImplementation highpassFilter;

public BandpassFilterButterworthImplementation
(double bottomFrequencyHz, double topFrequencyHz, int numSections, double Fs)
{
this.lowpassFilter = new LowpassFilterButterworthImplementation
(topFrequencyHz, numSections, Fs);
this.highpassFilter = new HighpassFilterButterworthImplementation
(bottomFrequencyHz, numSections, Fs);
}

public double compute(double input)
{
// compute the result as the cascade of the highpass and lowpass filters
return this.highpassFilter.compute(this.lowpassFilter.compute(input));
}
}
}```
C#
```using System;

namespace DSP
{
public class FIRFilterImplementation
{
protected double[] z;
public FIRFilterImplementation(int order)
{
this.z = new double[order];
}

public double compute(double input, double[] a)
{
// computes y(t) = a0*x(t) + a1*x(t-1) + a2*x(t-2) + ... an*x(t-n)
double result = 0;

for (int t = a.Length - 1; t >= 0; t--)
{
if (t > 0)
{
this.z[t] = this.z[t - 1];
}
else
{
this.z[t] = input;
}
result += a[t] * this.z[t];
}
return result;
}
}
}```
C#
```namespace DSP
{
public class IIRFilterImplementation
{
protected double[] z;
public IIRFilterImplementation(int order)
{
this.z = new double[order];
}

public double compute(double input, double[] a)
{
// computes y(t) = x(t) + a1*y(t-1) + a2*y(t-2) + ... an*y(t-n)
// z-transform: H(z) = 1 / (1 - sum(1 to n) [an * y(t-n)])
// a0 is assumed to be 1
// y(t) is not stored, so y(t-1) is stored at z[0],
// and a1 is stored as coefficient[0]

double result = input;

for (int t = 0; t < a.Length; t++)
{
result += a[t] * this.z[t];
}
for (int t = a.Length - 1; t >= 0; t--)
{
if (t > 0)
{
this.z[t] = this.z[t - 1];
}
else
{
this.z[t] = result;
}
}
return result;
}
}
}```

History

• 9th June, 2019: Initial submission

 First Prev Next
 Fs is a rather confusing parameter Member 1014767325-May-23 5:37 Member 10147673 25-May-23 5:37
 FYI - The download is missing the actual source code Lee.W.Spencer25-Feb-22 2:44 Lee.W.Spencer 25-Feb-22 2:44
 Quick and Dirty Frontend Pookalton19-Jan-21 20:47 Pookalton 19-Jan-21 20:47
 A Suggested Method To Add BKrell23-Nov-20 23:30 BKrell 23-Nov-20 23:30
 How can I use this on a real signal? Please add explanation. Hamed Nik12-Sep-19 0:02 Hamed Nik 12-Sep-19 0:02
 Re: How can I use this on a real signal? Please add explanation. HardworkingPearl5-Jan-20 21:10 HardworkingPearl 5-Jan-20 21:10
 Please add a Winform and a couple of nice signal sources for a real demo asiwel11-Jun-19 12:04 asiwel 11-Jun-19 12:04
 Last Visit: 31-Dec-99 18:00     Last Update: 19-Sep-24 14:34 Refresh 1