The Program

Program Instructions

The above is a JavaScript version of my BiQuadDesigner program — it's interactive and easy to use:

• Choose a Filter Type from the list at the left.
• Set Parameters for your filter including vertical and horizontal graph scale type (logarithmic or linear).
• Note that, when you change the sample rate or the center frequency, the center frequency for your filter is automatically limited to less than or equal to 1/2 of the sample rate if it has a higher value. This is consistent with the Nyquist-Shannon Sampling theorem, which in essence says that the meaningful bandwidth for a digital channel is limited to 1/2 the sample rate. A more useful value is 1/4 the sample rate, in the middle of the filter's usable bandwidth, but this is under the user's control.
• The filter constants list may appear to be missing one value: a₀. But during filter calculations this value is eliminated (i.e. normalized to equal 1) for the sake of efficiency and operating speed (see the JavaScript source for details, and a more in-depth discussion of this issue below).
• Hover your pointing device over the chart to see specific x,y value pairs.
• Spin your mouse wheel over the numerical inputs to change their values in a fast and intuitive way.
• Many practical programs and code sources are described later in this article, all of which accept the parameter sets this application produces.

Introduction

There are times when language fails us, when direct experience and insight aren't adequately represented by the words we use. For example, it is said (inaccurately) that the Inuit people have 20 words for our word "snow". It's also said (accurately) that people have at least 100 words for "drunk". Which brings us to our topic — the everyday word "filter," the most commonly applied term for this article's topic, just doesn't convey sufficient meaning.

In the present context, "filter" (verb) means to apply a technical method to select something, usually on the basis of frequency, and "filter" (noun) refers to such a method. By this reasoning, radio and TV receivers are all examples of specialized filters, our ears are sophisticated biological Fourier filters, and musical taste, granted that it has no universally agreed measure, might stand as a more abstract example.

We begin our topic by describing the most basic kinds of filters, explain their operation using examples, and gradually progress to a very practical embodiment called a "biquadratic" that is easy to implement and use.

The Basics

Present-day software-based filtering methods have their roots in older analog electronic designs, and filter design represents just one example of a fundamental change in electronic technology, a technology that began with analog designs, moved to integrated-circuit-based digital equivalents, and has finally arrived at pure algorithm-based designs. Another way to say it: filters were once conceived using mathematical equations and then realized in analog circuits, but at present, the mathematical equation is the filter — writing an appropriate equation begins and ends the process in one step, in some cases even if that equation has no closed-form solution.

Here's an example that models the simplest kind of filter — a resistor and capacitor arranged as a low-pass filter:

In the early days of electronics, one would analyze an RC circuit from first principles:

• The instantaneous current flowing in this simple circuit is proportional to the voltage across the resistor ($I = \frac{E}{R}$, Ohm's Law).
• The rate at which the capacitor's voltage changes is also proportional to current.
• At time zero, the voltage across the capacitor is zero, therefore the resistor's voltage, and the circuit's current, are at a maximum.
• As time passes, the voltage across the capacitor increases, which causes the resistor's voltage to decline. This means the circuit's current declines, which causes the capacitor's rate of change in voltage to decline.
• Another way to say this: the rate at which the capacitor's voltage approaches the source voltage declines in step with the remaining voltage difference between the capacitor and the source.

• This relationship can be expressed as a differential equation. For an unknown function f(t) representing the voltage across the capacitor as a function of time:

  • First term: the capacitor voltage, represented by trace B in the above graphic, is equal to zero at time zero (i.e. the resistor voltage is equal to v): \begin{equation} \displaystyle f(0) = 0 \end{equation}
  • Second term: the sum of the function and its first derivative, multiplied by the resistance and capacitance values, is equal to v: \begin{equation} \displaystyle f(t) + r c f'(t) = v \end{equation}
  • Term (2) above means the capacitor's rate of change in voltage is proportional to the circuit's instantaneous current.

• The above differential equation terms yield this classic solution:

  \begin{equation} \displaystyle f(t) = v e^{\frac{-t}{r c}} \end{equation} (confirm this at Wolfram Alpha)

• For expressions (1,2,3) above, the following definitions apply:

  • t: time, seconds
  • r: circuit resistance, ohms
  • c: circuit capacitance, farads
  • v: circuit voltage, volts
  • e: base of natural logarithms, approximately equal to 2.71828, and exactly equal to: \begin{equation} \displaystyle e = \lim_{n \to \infty} \left(\frac{1}{n} \right)^n \end{equation}
  • f(t): resistor voltage as a function of time

• Here is the same solution, restated to give the voltage across the capacitor:

  \begin{equation} \displaystyle f(t) = v - v e^{\frac{-t}{r c}} \end{equation}
• A real-time measurement of this analytical solution is identical to the modeled result graphed above.

The above analysis produces a result that is accurate for a simple step-function, the closing of a hypothetical switch at time zero. But what about other kinds of signals? How about a sine wave — is there a solution for that case? Yes — here is a diagram for a sinewave-driven RC circuit, and the resulting waveforms:

As before, there is a differential equation expression and solution for this RC circuit with a driving sine wave:

• The unknown function at any given time is equal to: \begin{equation} \displaystyle f(t) + r c f'(t) = \sin(\omega t) \end{equation}

  (where $\omega = 2 \pi f$ and $f$ = frequency)

• The solution for this differential equation is: \begin{equation} \displaystyle f(t) = -\frac{c \omega r \cos\left(\omega t\right) - \sin\left(\omega t\right)}{c^{2} \omega^{2} r^{2} + 1} \end{equation} (confirm this at Wolfram Alpha)

My point in presenting these examples is to show that there are analytical solutions for many similar differential equations, but they quickly become complex, and such solutions can only be applied to particular waveforms. But for the general case, for a filter that can be applied to any waveform or mixture of waveforms, there typically are no analytical solutions and a different approach is used. Here's how modern signal processing works:

• The first step in computer signal processing is to convert real-world signals from analog to digital form.
• The digital versions of the original signals arrive on a regular schedule, one set by the system's "sample rate".
• A sample rate is a sort of clock, a reliable timing standard. In practice one can establish times and frequencies by counting samples.
• Given a digital data stream and a known sample rate, we can write software equivalents for any of the many kinds of filters that have been designed over the years.

For an easily understood example, the following Python code models the above RC filter circuit — pay particular attention to the green code lines:

The result produced by the adjacent code

import math, pylab

# samples per second
sample_rate = 40000.0

# data lists for plotting
ta = []
da = []
db = []

# resistor value ohms
r = 1e3 

# capacitor value farads
c = 0.5e-6 

# time constant
tc =  1.0 / (r * c * sample_rate) 

# integration variable
vi = 0 

# signal frequency
f = 440 

for n in range(int(sample_rate) / 200):
    # t = time, seconds
    t = n / sample_rate
    ta.append(t)
    # generate the test signal
    v = math.sin(2 * math.pi * f * t)
    da.append(v)
    # perform the filtering function
    vi += (v-vi) * tc 
    db.append(vi)

pylab.plot(ta,da)
pylab.plot(ta,db)

pylab.show()
          

To understand this and similar examples, remember that a current of one ampere causes a one-farad capacitor to change voltage at a rate of one volt per second. Now notice this highlighted line from the above listing:

vi += (v-vi) * tc

Technically, this code line represents a difference equation, an equation that operates on the difference between adjacent values. Here is how this code line looks when expressed using formal mathematical notation:

\begin{equation} vi(n) = vi(n-1) + (v - vi(n-1))\ tc \end{equation}

Where $n$ represents a time-based sample index.

The above code line, executed on each data sample, produces a transfer function arbitrarily close to that produced by the analog RC circuit (and differential equation) listed above. The algorithm is written in such a way that the rate of change over time in vi (which represents the capacitor's voltage) is proportional to the remaining distance between vi and v, multiplied by the time constant — just as with the RC circuit being modeled.

The above algorithm represents a stepwise numerical solution to a differential equation that is not amenable to closed-form solution in the general case. The algorithm robustly models the original RC circuit for any applied waveform, however complex. Indeed, all software-based filters represent numerical solutions to differential equations, carried out in small, finite slices of time.

For completeness, here is code that models the earlier, step-function example, using the same methods:

The result produced by the adjacent code

import pylab

sample_rate = 4000.0

ta = []
da = []
db = []

# resistor value, ohms
r = 10e3 

# capacitor value, farads
c = 100e-6 

# time constant
tc =  1.0 / (r * c * sample_rate) 

vi = 0 # integration variable

for n in range(-int(sample_rate),int(sample_rate) * 5):
    t = n / sample_rate
    ta.append(t)
    # generate the test step function
    v = (0,1)[t > 0]
    da.append(v)
    # perform the filtering function
    vi += (v-vi) * tc 
    db.append(vi)

pylab.plot(ta,da)
pylab.plot(ta,db)

pylab.show()
          

As before, the algorithm produces the result one would expect to see from an analog circuit having the same component values and applied waveforms.

Remember about this class of numerical solver that, according to the Nyquist-Shannon sampling theorem, its legitimate bandwidth is limited to 1/2 the sampling rate — to increase bandwidth, one must increase the sampling rate.

This section has introduced the basics — the key ideas behind modern signal processing and software filter design. Next we show how to create and apply more sophisticated filters.

The Biquadratic Filter

In the section above we modeled a simple RC filter by applying a difference equation to a time series of digital data samples. Rather than try to create a possibly nonexistent closed-form solution to a difference equation, such modeling methods produce numerical results by applying the equation to time-series data samples.

The digital biquadratic filter is a more complex example of a difference equation, one that retains, and recursively applies, prior values from the data stream. Here is the biquadratic filter transfer function where the time-ordering of the data samples is represented by the $n$ in $z^n$, and the context of the equation is the Z-transform domain:

\begin{equation} \displaystyle H(z) = \frac{b_0 + b_1 z^{-1} + b_2 z^{-2}}{1 + a_1 z^{-1} + a_2 z^{-2}} \end{equation}

The name "biquadratic" refers to the fact that equation (9) contains two quadratic expressions, one in the numerator and another in the denominator.

Notice about equation (9) that variable a0 has been removed from the classical form of this equation. This results from an intentional optimization step in which the entire constant set {a0, a1, a2, b0, b1, b2} is divided by a0, such that (1) variable a0 is eliminated (it becomes equal to 1) and (2) the resulting code is faster and more efficient in practice. The JavaScript program at the top of this page follows this convention, as do all the other digital biquadratic filter embodiments available from this page.

As it happens, equation (9) has a relatively simple difference-equation form, one that can be easily applied to a time-based data stream:

\begin{equation} \displaystyle y(n) = b_0 x(n) + b_1 x(n-1) + b_2 x(n-2) - a_1 y(n-1) - a_2 y(n-2) \end{equation}

Biquadratic difference equation flow graph
(horizontal = time, vertical = data flow)

Where:

• $y(n)$ = present y result in the time series enumerated by n
• $y(n-1)$ = y result from time slice-1
• $y(n-2)$ = y result from time slice-2
• $x(n)$ = present x argument in the time series enumerated by n
• $x(n-1)$ = x argument from time slice-1
• $x(n-2)$ = x argument from time slice-2
• $b_0, b_1, b_2, a_1, a_2$ = precomputed constants

Here are some general comments about the biquadratic:

• The biquadratic filter scheme has sufficient flexibility to produce many kinds of filters — bandpass, low- and high-pass, peak and notch, among others — and with a wide variety of performance characteristics.
• All these traits are determined by the values chosen for the constant set {b₀,b₁,b₂,a₁,a₂}.
• As it happens, the constant set can be easily generated based on a few user-provided parameters — desired filter type, center (or corner) frequency, Q, and (for just a few filter types) gain.
• In most cases a biquadratic filter can be reconfigured as it runs — it can be modified as it processes data.
• There are obviously more complex filter types, but they tend toward instability in direct proportion to their complexity.
• If a biquadratic filter can't meet a particular filtering requirement, the reader is encouraged to consider cascading biquadratic filters (there is an example below) over choosing a more complex filter type that might be more difficult to understand or control.

Here is an example snippet of code (from this Java source file) that solves the biquadratic difference equation (Java):

// perform one filtering step
public double filter(double x) {
    y = b0 * x + b1 * x1 + b2 * x2 - a1 * y1 - a2 * y2;
    x2 = x1;
    x1 = x;
    y2 = y1;
    y1 = y;
    return (y);
}

            

In practice, the above function is applied to each arriving sample in the input data stream, and the return value represents the current filter state.

The biquadratic filter constants are computed based on user-provided parameters (filter type, center/corner frequency, Q, gain). The process is deterministic and relatively simple (here is a reference), and the function "reconfigure()" in this source file shows a typical embodiment. Again, in most cases this configuration step can be performed while the filter is operating.

For filter design and evaluation purposes, there is an equation that provides frequency-domain response levels (not time-domain results) — it provides the data for the BiQuadDesigner program's real-time graph displays. Here it is, written as a function:

\begin{equation} f(\phi)^{2} = \frac{16 \, b_{0} b_{2} \phi^{2} + {\left(b_{0} + b_{1} + b_{2}\right)}^{2} - 4 \, {\left(b_{0} b_{1} + 4 \, b_{0} b_{2} + b_{1} b_{2}\right)} \phi}{16 \, a_{2} \phi^{2} + {\left(a_{1} + a_{2} + 1\right)}^{2} - 4 \, {\left(a_{1} a_{2} + a_{1} + 4 \, a_{2}\right)} \phi} \end{equation}

Where $\phi = \sin(\frac{2 \pi f}{2\ \text{sample rate}})^2$

Obviously one could model filters by sweeping a test signal across their intended frequency domain (and I have used this method), but this approach has a number of difficulties — if the sweep rate is slow, ths process requires an inordinate amount of time, but if the sweep rate is fast, one ends up testing the filter's response time rather than its spectrum.

This brings up an important point — each biquadratic filter design has a characteristic spectrum (the graph BiQuadDesigner provides) but also a response time. For reasons beyond the scope of this article, filters with high Q or gain settings tend to have slow response times in proportion to their selectivity. This is an important factor in filter design, and it often requires a more complete kind of testing than can be provided by a modeling tool like BiQuadDesigner.

To summarize this section:

• Biquadratic filters are relatively simple to create and use.
• Creating biquadratic constant sets is a well-established practice, requiring only a handful of user-provided parameters.
• When confronted by a filtering task that cannot be accomplished using a single biquadratic filter, readers are encouraged to try cascading filters before considering a more complex filter type.
• GPL source code is readily available in several languages to implements biquadratic filters (examples: Java, Python).
• I have also prepared an interactive Sage worksheet that creates and demonstrates biquadratic filters.

Filter Design Case History

I've recently been writing software radios for various purposes, mostly having to do with boating and collecting weather data from a shortwave radio while underway. Examples are:

• JWX — a java-based weather chart receiver/converter. This program receives and displays daily weather chart broadcasts acquired from the shortwave bands.
• JNX — almost the same name, and a similar purpose. This program receives and decodes a text-based weather service called Navtex.

Both these programs acquire their data from a personal computer's sound card, in essence an A/D converter, connected to a shortwave radio. Both of them need to decode an FM signal to acquire the information of interest. So both these projects are perfect applications for biquadratic filters.

How does one apply a biquadratic filter to the task of decoding an FM signal? Well, there are many approaches to decoding FM, some exotic, some relatively ineffective. After evaluating many approaches and conscious of the fact that both my programs need to be able to deal with relatively high noise levels and weak signals, I began evaluating different arrangements of biquadratic filters.

A perfect FM demodulator would reject amplitude noise and out-of-band signals, then efficiently convert the FM signal into an amplitude with reasonable linearity. After much experimentation, I chose this filter arrangement:

• First, a relatively broad (i.e. low Q) bandpass filter to help attenuate out-of-band signals:
  
• Next, I added two narrower (higher Q) bandpass filters meant to detect the space and mark frequencies — here's how the three filters look when plotted together:
  
• I then applied this transform to the filter combination:
  \begin{equation} v = \left(\text{mark}^2 - \text{space}^2\right) \times \text{bandpass} \end{equation}
• After application of a final low-pass filter that works in a different frequency domain (not shown), the system achieves this result:
  

That's a pretty ideal transfer function — essentially a straight line between the mark and space frequencies. I should add that, in the pre-computer era, such a linear transfer function was rare, even in a radio meant to decode high-fidelity music broadcasts.

Remember about a design like this that one can easily focus on one aspect of the problem at the expense of others. For example, one is tempted to increase the Q setting (the setting that makes the mark and space filters narrower than the bandpass filter) to decrease the system's response to noise. But in doing so one can easily make the system's impulse response too slow to efficiently decode the desired signals. This is why theoretical testing with a tool like BiQuadDesigner must be accompanied by data gathering using realistic signals mixed with noise. It happens I have another free Java application named JSigGen that is easily configured to generate various combinations of modulated signals and noise.

Program Background

The program at the top of this page can be used to explore the properties of biquadratic filters.This table outlines the properties of the available filter types:

Filter Type	Q adj	Gain adj	Comments	Image
Bandpass	Y	N	The most generally useful filter type.
Low-pass	Y	N	For low-pass and high-pass biquadratic filters, one normally sets Q = 0.707 ($\frac{1}{\sqrt{2}}$) to achieve a Butterworth filter transfer function with a 3 DB drop at the specified operating frequency. Higher Q settings produce an often-undesirable peak near the center frequency and dynamic instability in operation.
High-pass	Y	N
Peak	Y	Y	This filter is a bit tricky to set up, because both Q and gain are effective. The idea is that one can use the gain control to set a nonzero base gain level that applies to all frequencies, then use the frequency and Q controls to set a narrow peak to exceed that level. Note also that, with a negative gain setting, the relation between the plateau and peak is reversed.
Notch	Y	N	This filter is more or less the opposite of the "Peak" filter — it creates a narrow rejection band, the width of which is set by the Q control. (But no plateau as with "Peak".)
Lowshelf	N	Y	Lowshelf and highshelf filters provide a sort of "plateau" effect, under control of the gain setting, and not unlike the "Peak" filter described above. Note that negative gain settings reverse the identity of the filter — lowshelf becomes highshelf and the reverse.
Highshelf	N	Y

Here are some points about the JavaScript program:

• The numerical entry controls will accept numbers that are typed in, or (on most browsers) they can be controlled by your mouse — spin your mouse wheel over the numerical input controls to increment/decrement their values.
• The program's "Filter Constants" list shows the properties of the present filter. In typical practice, the data available there is not very important because normal biquadratic filter applications and classes include a way to create any desired filter constants on the fly during instantiation of the class. This means BiQuadDesigner is really meant to allow interactive exploration of biquadratic filters of various type and settings, rather than pin down a design's details.
• Notice that the vertical chart dimension can have a linear or decibel scale. The decibel scale is useful for comparing biquadratic filters to other filter types and to each other.
• Notice about the graph that hovering your pointing device over the graph produces a readout of x,y data pairs, with both linear and decibel values. This allows the user to establish details of the filter's properties.
• While using the program, remember the Nyquist-Shannon sampling theorem — the meaningful domain for filter design and modeling is between 0 Hz and 1/2 the sampling rate. In some cases this fact can be turned to the user's advantage — setting a low sampling rate may uncover a filter design that may be optimal for low-clock-rate environments.

Download

BiQuadDesigner is Copyright 2017, P. Lutus, and is released under the GPL.

Here are the download options:

• First, to run the Java version of BiQuadDesigner, you must have Java installed. It's free, and it is here. Those running Linux may already have an open-source Java runtime engine installed, which is more than satisfactory.
• Java JAR file. Suitable for Linux, Macintosh, essentially everyone not running Windows. To run a Java JAR file, issue this command:

  java -jar BiQuadDesigner.jar
                  

  The above may be of course placed in a shell script.
• Windows executable installer. Suitable for all current Windows environments.
• Python script for configuring and running Biquadratic filters.
• Source archive/NetBeans project. A full source download organized to work with the current Netbeans Java development environment.

Version History

• Version 1.7 10.25.2019. By user request increased the downloadable Java version's maximum frequency entry range to 100 MHz.
• Version 1.6 10.23.2017. Updated the Java application's code to bring it into correspondence with the JavaScript version.
• Version 1.5 10.22.2017. Changed the decibel readout convention to 20 * log10 in all forms of the biquadratic filter, to be consistent with standard voltage/current measurement conventions.
• Version 1.4 10.21.2017. Fixed an error in the function that provides a static filter output for testing and display purposes.
• Version 1.3 10.21.2017. Extensively changed the Java application to better correspond with this JavaScript version and its features.
• Version 1.2 12.09.2013. Added explicit clipboard copy function and control (only works when operated as desktop application).
• Version 1.1 03.24.2011. Improved the mouse cursor readout to display simultaneous unit and decibel values.
• Version 1.0 03.23.2011. Initial public release.

References

• Example biquadratic filter class source file (Java)
• Example biquadratic filter class source file (Python)
• Sage biquadratic filter interactive worksheet
• RC Circuit (Wikipedia)
• Digital Biquad Filter (Wikipedia)
• Cookbook formulae for audio EQ biquad filter coefficients