The story so far… in When off-the-shelf won’t do: the step-by-step birth of a custom digital filter, Part 1 I created a continuous-time filter transfer function for my special-purpose 10th order decimation filter. It had a final stopband slope of -18 dB per octave in the far stopband at frequencies much higher than the transmission zeros (something that ‘standard ‘catalogue’ filters can’t do) and a reasonably linear phase response across much of the audio band.
Now, there’s much more to making a practical filter than just getting a transfer function that calculates out nicely. I see a lot of younger engineers do wonders with “filter design packages”, producing beautiful plots of fine filters that work fine as long as you implement them in double-precision floating point. In the real world we have limited resources, whether analog or digital. Most filters I see, analog or digital, are not limited by the approximation stage (finding a transfer function that fits the requirements) but by the synthesis or implementation stages (finding coefficient or component values for a topology that can apply the filtering magic to your signal).
In the case of a filter that’s made from a cascade of 2nd order sections (“biquads”), the two critical factors to consider are the choice of biquad topology, and the way in which you distribute the filter poles and zeros (the roots of the denominator and numerator polynomials respectively). This is just as true in the digital domain as in the analog domain.
In this case, the topology decision was made for me because the filter needed to be implemented straight away using the Direct Form code we had available at the time for the Digital Filter Block in our PSoC device. This topology is definitely not the best that can be done for signal to noise ratio when the ratio of sample frequency to pole frequency is high. But it’s easy to design and understand, and economical in resource terms. It’s the digital equivalent of the “single amplifier biquad”, the second order analog filter made with a single op-amp. These often have poorer performance metrics than more complicated filter sections. But sometimes it’s hard to persuade people that it’s better to use two or even three opamps rather than one. People will thing you’re just trying to sell them op-amps. Anyway, back to the plot…
Step 4: Optimize pole-zero ordering for best passband SNR
The transfer functions of the 7th order filter and the 3rd order equalizer are cascaded to produce a 10th order transfer function. This is the product of five 2nd order rational functions, each of which has a numerator and a denominator with some frequency-dependent behaviour.
Because they are all multiplied together, the numerators can all be swapped round without affecting the overall response. Likewise the denominators. What this means is that there are 5! ways of selecting an order for the numerators, and also 5! ways of independently selecting an order for the denominators. This means that there are up to 14400 ways of making up this transfer function (depending on how many duplicated factors there are in the transfer function). That naturally begs the question: which one of these is ‘best’, in some meaningful way?
This “pole/zero ordering problem” is an important one in filter design. Note that it’s quite distinct from the internal overload problem resolved by calculating scaling factors. I solved this problem for several specific analogue filter design approaches a long time ago, but my old solutions weren’t quite right for the digital biquad case. So I spent some time analyzing the process from scratch and building it into my normal transfer function management tools (also ancient, and with direct lineage back to the BBC Micro days) so that some sort of ‘best’ solution could be automated.
The basic approach works like this (it’s done in an outer loop that executes for a count equal to the number of biquads): look exhaustively at every combination of an available numerator and an available denominator, and work out how much SNR you are “throwing away”. This is calculated by looking at the effective noise bandwidth of all the remaining biquad pieces that are going to come after the section you’re looking at. The higher the noise bandwidth, the more of a contribution of noise from this biquad will ‘punch through’ to the output. It’s straightforward to pick the combination that degrades the overall SNR the least.
The calculation is implicitly scaled by ensuring that the peak passband gain of the tested combination is unity, to prevent simple sinewave overload. This inherently combines the scaling part of the filter design process into the pole/zero ordering process, and no subsequent scaling is needed.
With a well-formed transfer function, this process will always converge to an answer. I ran the process over the 10th order transfer function obtained by combining the raw outputs from the filter design and equalizer design phases, to get a final transfer function complete with implicit scaling factors. This is still an analogue transfer function, but all the scaling will be completely preserved when we ‘go digital’.
I did make one change to the automated transfer function process. One of the filter sections is an allpass section with a numerator that’s the complex conjugate of the denominator. These two pieces have the same magnitude for all frequencies, and the combination will always ‘win’ this selection competition, positioning all the pure allpass sections at the beginning of the sequence. On very ‘difficult’ input signals, this is not optimum from a transient overload perspective.
The pragmatic solution is to relocate the allpass network to a position in the filter after the first ‘proper’ lowpass section, which “takes the sting off” the input signal. This doesn’t invalidate the rationale behind the ordering process. Initially, I placed the equalizer section second, but simulations confirmed that certain high frequency squarewave inputs could cause distinct peaking, and therefore potential overload, to occur. Placing the equalizer section third in the cascade gave good behaviour on these difficult out-of-band signals.
Note that although the scaling has been carried out to ensure no sinewave clipping, it’s possible to introduce ‘pathological’ signals to the filter that cause an output amplitude that’s considerably larger than the input one. This is hard to quantify, and is best managed by leaving an additional ‘safety’ scaling factor for the input signal. My standing recommendation for direct-form filters is “2 bits”, or 12 dB of headroom. I have never seen any (sinewave-scaled) IIR filter breach this limit on any pathological signal, whatever response type is used, and feel it’s sound, if a little empirical.
Now we have a cascade of transfer function blocks whose numerators and denominators have been sequenced to give us (by the particular metric used) the least summed contribution of the residual noise that each section contributes at the output. This process could be made considerably more sophisticated, given some more thinking time. But just addressing the cascade ordering at all, versus ignoring it, gets you most of the way towards an optimum solution.
Figure 1 shows the responses simulated at each of the outputs in the final ordered cascade. You can’t see the blue trace for section 2 because section 3 is the phase equalizer section which has a flat frequency response. I really quite like this plot because it tells a story of the teamwork between poles and zeroes, gradually creating that lovely flat final frequency response.
Want to see how the time-domain response varies as you run down the cascade? Thought you did. Figure 2 shows the peaking that you get if you excite a sharp filter like this with an edge. Another reason why you really should leave some signal level margin at the input. A common problem I’ve seen decade in, decade out is the subtle error caused by internal overload in a filter, especially on transients, that it not visible at the final output.
The magenta trace at the top shows the expected pre- and post-ringing characteristic of a filter with close to linear phase. The pre- and post-ring peaks aren’t equal, reflecting the incomplete passband phase equalization.
And that’s all we need to do to the analogue (continuous time) transfer function. Now it’s time to go digital.
Step 5: Bilinear transform to digital domain and the DFB
This is straightforward, you can easily do it in a spreadsheet. The only additional piece of information you need is the sample rate at which the filter will be running, in this case 352.8 ksps. In this older article I talked briefly about the bilinear transform and how to apply it to a continuous time transfer function.
The spreadsheet can also be used to scale the coefficients to suit the capabilities of the platform. In my case, this means fitting the coefficients to suit the 24-bit fixed point arithmetic to PSoC’s little highly-tuned race-car of a filter engine, the Digital Filter Block. You can use the standard PSoC Creator tool to plot the final result; I’m a belt-and-braces kind of guy and I always like to see something in Excel and something out of a good SPICE simulator (you all know which one I like – it’s the one you all like too…)
There’s a final teaching point right there in figure 4. The continuous time filter (blue) has a frequency response that slides smoothly on down monotonically as the input frequency increases. By contrast, the frequency response of a sampled-domain filter is periodic, and you get copies of it at every integer multiple of the sampling frequency.
Anyway, there you go. We implemented this, and it worked fine. We didn’t pursue the underlying architecture to full production – my awesome idea didn’t work as well as I’d hoped first time round, and so we moved on to something else. But I was really happy that the project ended up with a super filter that I could get so much teaching and mentoring potential from. Hope some of that enthusiasm filtered through to you too! / Kendall