Digital Bessel filter problem

I'm getting strange results under certain circumstances when using a digital Bessel filter. A minimal working example is shown below. I get similarly weird results using slight variations, such as using bilinear instead of impinvar, or by using the lower level functions rather than besself. The example raw data is simply a line with slope one (black). I expect the filtered signal (red) to be a line with the same slope, slightly delayed from the original. I can achieve this result with a higher cutoff frequency, or by using a different filter, such as Butterworth. Any ideas what it happening here? Thanks! R2015b

numpoints=10000;
data=1:numpoints;
sampling_freq=10;
order=8;
cutoff=0.03;
[b, a] = besself(order,cutoff*2*pi);
[bz, az] = impinvar(b,a,sampling_freq);
filtered=filter(bz,az,data);
plot(data,'k');
hold on
plot(filtered,'r');

Answer 1

[Scott Webster]

I seem to have solved my problem. It appears that the issue is related to numerical limitations in the transfer function form of the filter. My solution is to use sosfilt instead: a second order sections representation. It appears that perhaps the Bessel transfer function coefficients become problematic at lower orders than for some other filters (for me with an eighth order lowpass filter with "low" cutoff), so it wasn't an immediately obvious concern.

Answer 2

[Star Strider]

Bessel filters are excellent analogue (hardware) designs because of their lack of phase distortion. However the MATLAB documentation says:

• besself does not support the design of digital Bessel filters.

and my reference (Proakis and Manolakis, Digital Signal Processing, 2007 p. 727) says ‘However, we should emphasize that the linear-phase characteristics of the analog [Bessel] filter are destroyed in the process of converting the filter into the digital domain by means of the transformations described previously.’

So if you want linear-phase filtering, use filtfilt with any filter design.

[Scott Webster]

I recognize that there are some of these subtle points regarding the discretization of Bessel filters, but the problem here seems to be something more like a numerical precision issue. The result is highly dependent on small changes in the input parameters near these values, but away from these values the filter seems to work as expected. I guess I don't expect imperfections in the digital conversion to cause the result to be so far wrong. Any specific thoughts about what exactly is going wrong in this case? Thanks!

[Star Strider]

My pleasure.

I don’t use Bessel filters except in hardware, and then principally as anti-aliasing filters, so I have no experience with them as discrete filters. Based on what I’ve read, there is no reliable way of discretising them. I haven’t analysed the reasons they cannot be reliably discretised, since I never needed to do so. I will leave you to explore this if you want to.

Based on my experience and what I’ve read, I would abandon the effort. I would use another filter design that can be reliably discretised, and filtfilt to produce a phase-neutral response.

[Scott Webster]

One issue with filtfilt is that it is non-causal.

[Star Strider]

The second-order-section implementation is superior to others for all filters.

Using it and using filtfilt will produce the sort of filter response you want, but so would any filter design using second-order-section implementation and the filtfilt filter implementation.

[Scott Webster]

Yes, I would recommend using SOS for any filter type, especially for higher orders. I'm not sure if there are any fine points about why you might NOT want to use it if you don't have numerical issues.

However, different types of filters will have different types of responses, even when used with filtfilt in SOS mode, so you still need to choose the right filter for the job, Butterworth, Bessel etc.

For example, in the following image the blue is the original data, green is Bessel filtered, red is Butterworth, same parameters. Ringing in Butterworth.

[Star Strider]

The Butterworth filter likely has a larger band-pass region and sharper roll-off than the Bessel filter (that has probably the most gradual roll-off of any filter design).

That said, this isn’t actually ‘ringing’ (Gibbs phenomenon) but the filter’s frequency selection of the input signal. The Butterworth design is passing more frequencies than the Bessel design, so the approximation of the input waveform is more precise. The Bessel design is approximating the square wave with a sine wave, while the Butterworth design is transmitting an additional harmonic to approximate the square wave. I suspect the freqz function would demonstrate a corresponding difference the the two filter pass-bands.

True Gibbs phenomenon would have ringing only at the sharp transitions that would decay by the end of the ‘constant’ section of the waveform. (I don’t know if that would be true of filtfilt though, since I never did that experiment.)

[Scott Webster]

Thanks, interesting. Yeah, I am not trying to imply that the oscillations are a "problem" with the Butterworth filter, just that the two filters are, of course, different. And, as you point out, it is perhaps not trivial to try to define two different filter types to have the "same parameters" in the first place, so hard to compare.

[Star Strider]

Actually, they don’t have the same parameters. The transfer functions of each are entirely different, as are their pass-bands and stop-bands.

The ‘correct’ filter design depends on what you want to do, and if you want to do it in software (discrete) or hardware (continuous).

Signal processing in any domain is quite definitely not a trivial exercise!