Posts

Although the design of FIR filters with linear phase is an easy task. This is certainly not true for IIR filters that usually have a highly non-linear phase response, especially around the filter’s cut-off frequencies. This article discusses the characteristics needed for a digital filter to have linear phase, and how an IIR filter’s passband phase can be modified in order to achieve linear phase using all-pass equalisation filters.

Why do we need linear phase filters?

Digital filters with linear phase have the advantage of delaying all frequency components by the same amount, i.e. they preserve the input signal’s phase relationships. This preservation of phase means that the filtered signal retains the shape of the original input signal. This characteristic is essential for audio applications as the signal shape is paramount for maintaining high fidelity in the filtered audio. Yet another application area that requires this, is ECG biomedical waveform analysis, as any artefacts introduced by the filter may be misinterpreted as heart anomalies.

The following plot shows the filtering performance of a Chebyshev type I lowpass IIR on ECG data – input waveform (shown in blue) shifted by 10 samples (\(\small \Delta=10\)) to approximately compensate for the filter’s group delay. Notice that the filtered signal (shown in red) has attenuated, broadened and added oscillations around the ECG peak, which is undesirable.

Figure 1: IIR lowpass filtering result with phase distortion

In order for a digital filter to have linear phase, its impulse response must have conjugate-even or conjugate-odd symmetry about its midpoint. This is readily seen for an FIR filter,

\(\displaystyle H(z)=\sum\limits_{k=0}^{L-1} b_k z^{-k}\tag{1} \)

With the following constraint on its coefficients,

\(\displaystyle b_k=\pm\, b^{\ast}_{L-1-k}\tag{2} \)

which leads to,

\(\displaystyle z^{L-1}H(z) = \pm\, H^\ast (1/z^\ast)\tag{3} \)

Analysing Eqn. 3, we see that roots (zeros) of \(\small H(z)\) must also be the zeros of  \(\small H^\ast (1/z^\ast)\). This means that the roots of \(\small H(z)\) must occur in conjugate reciprocal pairs, i.e.  if \(\small z_k\) is a zero of \(\small H(z)\), then \(\small H^\ast (1/z^\ast)\) must also be a zero.

Why IIR filters do not have linear phase

A digital filter is said to be bounded input, bounded output stable, or BIBO stable, if every bounded input gives rise to a bounded output. All IIR filters have either poles or both poles and zeros, and must be BIBO stable, i.e.

\(\displaystyle \sum_{k=0}^{\infty}\left|h(k)\right|<\infty \tag{4}\)

Where, \(\small h(k)\) is the filter’s impulse response. Analyzing Eqn. 4, it should be clear that the BIBO stability criterion will only be satisfied if the system’s poles lie inside the unit circle, since the system’s ROC (region of convergence) must include the unit circle. Consequently, it is sufficient to say that a bounded input signal will always produce a bounded output signal if all the poles lie inside the unit circle.

The zeros on the other hand, are not constrained by this requirement, and as a consequence may lie anywhere on z-plane, since they do not directly affect system stability. Therefore, a system stability analysis may be undertaken by firstly calculating the roots of the transfer function (i.e., roots of the numerator and denominator polynomials) and then plotting the corresponding poles and zeros upon the z-plane.

Applying the developed logic to the poles of an IIR filter, we now arrive at a very important conclusion on why IIR filters cannot have linear phase.

A BIBO stable filter must have its poles within the unit circle, and as such in order to get linear phase, an IIR would need conjugate reciprocal poles outside of the unit circle, making it BIBO unstable.

Based upon this statement, it would seem that it’s not possible to design an IIR to have linear phase. However, a discussed below, phase equalisation filters can be used to linearise the passband phase response.

Phase linearisation with all-pass filters

All-pass phase linearisation filters (equalisers) are a well-established method of altering a filter’s phase response while not affecting its magnitude response. A second order (Biquad) all-pass filter is defined as:

\( A(z)=\Large\frac{r^2-2rcos \left( \frac{2\pi f_c}{fs}\right) z^{-1}+z^{-2}}{1-2rcos \left( \frac{2\pi f_c}{fs}\right)z^{-1}+r^2 z^{-2}}\tag{5} \)

Where, \(\small f_c\) is the centre frequency, \(\small r\) is radius of the poles and \(\small f_s\) is the sampling frequency. Notice how the numerator and denominator coefficients are arranged as a mirror image pair of one another.  The mirror image property is what gives the all-pass filter its desirable property, namely allowing the designer to alter the phase response while keeping the magnitude response constant or flat over the complete frequency spectrum.

Cascading an APF (all-pass filter) equalisation cascade (comprised of multiple APFs) with an IIR filter, the basic idea is that we only need to linearise the phase response the passband region. The other regions, such as the transition band and stopband may be ignored, as any non-linearities in these regions are of little interest to the overall filtering result.

The challenge

The APF cascade sounds like an ideal compromise for this challenge, but in truth a significant amount of time and very careful fine-tuning of the APF positions is required in order to achieve an acceptable result. Each APF has two variables: \(\small f_c\) and \(\small r\) that need to be optimised, which complicates the solution. This is further complicated by the fact that the more APF stages that are added to the cascade, the higher the overall filter’s group delay (latency) becomes. This latter issue may become problematic for fast real-time closed loop control systems that rely on an IIR’s low latency property.

Nevertheless, despite these challenges, the APF equaliser is a good compromise for linearising an IIRs passband phase characteristics.

The APF equaliser

ASN Filter Designer provides designers with a very simple to use graphical all-phase equaliser interface for linearising the passband phase of IIR filters. As seen below, the interface is very intuitive, and allows designers to quickly place and fine-tune APF filters positions with the mouse. The tool automatically calculates \(\small f_c\) and \(\small r\), based on the marker position.

APF equaliser ASN Filter Designer

Right clicking on the frequency response chart or on an existing all-pass design marker displays an options menu, as shown on the left.

You may add up to 10 biquads (professional version only).

An IIR with linear passband phase

Designing an equaliser composed of three APF pairs, and cascading it with the Chebyshev filter of Figure 1, we obtain a filter waveform that has a much a sharper peak with less attenuation and oscillation than the original IIR – see below. However, this improvement comes at the expense of three extra Biquad filters (the APF cascade) and an increased group delay, which has now risen to 24 samples compared with the original 10 samples.

IIR lowpass filtering result with three APF phase equalisation filters
(minimal phase distortion)
IIR lowpass filtering result with three APF phase equalisation filters
(minimal phase distortion)

The frequency response of both the original IIR and the equalised IIR are shown below, where the group delay (shown in purple) is the average delay of the filter and is a simpler way of assessing linearity.

IIR without equalisation cascade
IIR without equalisation cascade

IIR with equalisation cascade
IIR with equalisation cascade

Notice that the group delay of the equalised IIR passband (shown on the right) is almost flat, confirming that the phase is indeed linear.

Automatic code generation to Arm processor cores via CMSIS-DSP

The ASN Filter Designer’s automatic code generation engine facilitates the export of a designed filter to Cortex-M Arm based processors via the CMSIS-DSP software framework. The tool’s built-in analytics and help functions assist the designer in successfully configuring the design for deployment.

Before generating the code, the IIR and equalisation filters (i.e. H1 and Heq filters) need to be firstly re-optimised (merged) to an H1 filter (main filter) structure for deployment. The options menu can be found under the P-Z tab in the main UI.

All floating point IIR filters designs must be based on Single Precision arithmetic and either a Direct Form I or Direct Form II Transposed filter structure. The Direct Form II Transposed structure is advocated for floating point implementation by virtue of its higher numerically accuracy.

Quantisation and filter structure settings can be found under the Q tab (as shown on the left). Setting Arithmetic to Single Precision and Structure to Direct Form II Transposed and clicking on the Apply button configures the IIR considered herein for the CMSIS-DSP software framework.

Select the Arm CMSIS-DSP framework from the selection box in the filter summary window:

The automatically generated C code based on the CMSIS-DSP framework for direct implementation on an Arm based Cortex-M processor is shown below:

The ASN Filter Designer's automatic code generator generates all initialisation code, scaling and data structures needed to implement the linearised filter IIR filter via Arm's CMSIS-DSP library.

The ASN Filter Designer’s automatic code generator generates all initialisation code, scaling and data structures needed to implement the linearised filter IIR filter via Arm’s CMSIS-DSP library.

What we have learnt

The roots of a linear phase digital filter must occur in conjugate reciprocal pairs. Although this no problem for an FIR filter, it becomes infeasible for an IIR filter, as poles would need to be both inside and outside of the unit circle, making the filter BIBO unstable.

The passband phase response of an IIR filter may be linearised by using an APF equalisation cascade. The ASN Filter Designer provides designers with everything they need via a very simple to use, graphical all-pass phase equaliser interface, in order to design a suitable APF cascade by just using the mouse!

The linearised IIR filter may be exported via the automatic code generator using Arm’s optimised CMSIS-DSP library functions for deployment on any Cortex-M microcontroller.

 

 

Download demo now

Licencing information

A digital filter is a mathematical algorithm that operates on a digital dataset (e.g. sensor data) in order extract information of interest and remove any unwanted information. Applications of this type of technology, include removing glitches from sensor data or even cleaning up noise on a measured signal for easier data analysis. But how do we choose the best type of digital filter for our application? And what are the differences between a IIR filter and a FIR filter?

Digital filters are divided into the following two categories:

  • Infinite impulse response (IIR)
  • Finite impulse response (FIR)

As the names suggest, each type of filter is categorised by the length of its impulse response. However, before beginning with a detailed mathematical analysis, it is prudent to appreciate the differences in performance and characteristics of each type of filter.

Example

In order to illustrate the differences between an IIR and FIR, the frequency response of a 14th order FIR (solid line), and a 4th order Chebyshev Type I IIR (dashed line) is shown below in Figure 1.  Notice that although the magnitude spectra have a similar degree of attenuation, the phase spectrum of the IIR filter is non-linear in the passband (\(\small 0\rightarrow7.5Hz\)), and becomes very non-linear at the cut-off frequency, \(\small f_c=7.5Hz\). Also notice that the FIR requires a higher number of coefficients (15 vs the IIR’s 10) to match the attenuation characteristics of the IIR.

FIR vs IIR: frequency response of a 14th order FIR (solid line), and a 4th order Chebyshev Type I IIR (dashed line); Fir Filter, IIR Filter
Figure 1: FIR vs IIR: frequency response of a 14th order FIR (solid line), and a 4th order Chebyshev Type I IIR (dashed line)

These are just some of the differences between the two types of filters. A detailed summary of the main advantages and disadvantages of each type of filter will now follow.

IIR filters

IIR (infinite impulse response) filters are generally chosen for applications where linear phase is not too important and memory is limited. They have been widely deployed in audio equalisation, biomedical sensor signal processing, IoT/IIoT smart sensors and high-speed telecommunication/RF applications.

Advantages

  • Low implementation cost: requires less coefficients and memory than FIR filters in order to satisfy a similar set of specifications, i.e., cut-off frequency and stopband attenuation.
  • Low latency: suitable for real-time control and very high-speed RF applications by virtue of the low number of coefficients.
  • Analog equivalent: May be used for mimicking the characteristics of analog filters using s-z plane mapping transforms.

Disadvantages

  • Non-linear phase characteristics: The phase charactersitics of an IIR filter are generally nonlinear, especially near the cut-off frequencies. All-pass equalisation filters can be used in order to improve the passband phase characteristics.
  • More detailed analysis: Requires more scaling and numeric overflow analysis when implemented in fixed point. The Direct form II filter structure is especially sensitive to the effects of quantisation, and requires special care during the design phase.
  • Numerical stability: Less numerically stable than their FIR (finite impulse response) counterparts, due to the feedback paths.

FIR filters

FIR (finite impulse response) filters are generally chosen for applications where linear phase is important and a decent amount of memory and computational performance are available. They have a widely deployed in audio and biomedical signal enhancement applications. Their all-zero structure (discussed below) ensures that they never become unstable for any type of input signal, which gives them a distinct advantage over the IIR.

Advantages

  • Linear phase: FIRs can be easily designed to have linear phase. This means that no phase distortion is introduced into the signal to be filtered, as all frequencies are shifted in time by the same amount – thus maintaining their relative harmonic relationships (i.e. constant group and phase delay). This is certainly not case with IIR filters, that have a non-linear phase characteristic.   
  • Stability: As FIRs do not use previous output values to compute their present output, i.e. they have no feedback, they can never become unstable for any type of input signal, which is gives them a distinct advantage over IIR filters.
  • Arbitrary frequency response: The Parks-McClellan and ASN FilterScript’s firarb() function allow for the design of an FIR with an arbitrary magnitude response. This means that an FIR can be customised more easily than an IIR.
  • Fixed point performance: the effects of quantisation are less severe than that of an IIR.

Disadvantages

  • High computational and memory requirement: FIRs usually require many more coefficients for achieving a sharp cut-off than their IIR counterparts. The consequence of this is that they require much more memory and significantly a higher amount of MAC (multiple and accumulate) operations. However, modern microcontroller architectures based on the Arm’s Cortex-M cores now include DSP hardware support via SIMD (signal instruction, multiple data) that expedite the filtering operation significantly.
  • Higher latency: the higher number of coefficients, means that in general an FIR is less suitable than an IIR for fast high throughput applications. This becomes problematic for real-time closed-loop control applications, where an FIR filter may have too much group delay to achieve loop stability.
  • No analog equivalent: using the Bilinear, matched z-transform (s-z mapping), an analog filter can be easily be transformed into an equivalent IIR filter.  However, this is not possible for an FIR as it has no analog equivalent.

Mathematical definitions

As discussed in the introduction, the name IIR and FIR originate from the mathematical definitions of each type of filter, i.e. an IIR filter is categorised by its theoretically infinite impulse response,

\(\displaystyle
y(n)=\sum_{k=0}^{\infty}h(k)x(n-k)
\)

and an FIR categorised by its finite impulse response,

\(\displaystyle
y(n)=\sum_{k=0}^{N-1}h(k)x(n-k)
\)

We will now analyse the mathematical properties of each type of filter in turn.

IIR definition

As seen above, an IIR filter is categorised by its theoretically infinite impulse response,

\(\displaystyle y(n)=\sum_{k=0}^{\infty}h(k)x(n-k) \)

Practically speaking, it is not possible to compute the output of an IIR using this equation. Therefore, the equation may be re-written in terms of a finite number of poles \(\small p\) and zeros \(\small q\), as defined by the linear constant coefficient difference equation given by:

\(\displaystyle
y(n)=\sum_{k=0}^{q}b_k x(n-k)-\sum_{k=1}^{p}a_ky(n-k)
\)

where, \(\small a_k\) and \(\small b_k\) are the filter’s denominator and numerator polynomial coefficients, who’s roots are equal to the filter’s poles and zeros respectively. Thus, a relationship between the difference equation and the z-transform (transfer function) may therefore be defined by using the z-transform delay property such that,

\(\displaystyle
\sum_{k=0}^{q}b_kx(n-k)-\sum_{k=1}^{p}a_ky(n-k)\quad\stackrel{\displaystyle\mathcal{Z}}{\longleftrightarrow}\quad\frac{\sum\limits_{k=0}^q b_kz^{-k}}{1+\sum\limits_{k=1}^p a_kz^{-k}}
\)

As seen, the transfer function is a frequency domain representation of the filter. Notice also that the poles act on the output data, and the zeros on the input data. Since the poles act on the output data, and affect stability, it is essential that their radii remain inside the unit circle (i.e. <1) for BIBO (bounded input, bounded output) stability. The radii of the zeros are less critical, as they do not affect filter stability. This is the primary reason why all-zero FIR (finite impulse response) filters are always stable.

BIBO stability

A linear time invariant (LTI) system (such as a digital filter) is said to be bounded input, bounded output stable, or BIBO stable, if every bounded input gives rise to a bounded output, as

\(\displaystyle \sum_{k=0}^{\infty}\left|h(k)\right|<\infty \)

Where, \(\small h(k)\) is the LTI system’s impulse response. Analyzing this equation, it should be clear that the BIBO stability criterion will only be satisfied if the system’s poles lie inside the unit circle, since the system’s ROC (region of convergence) must include the unit circle. Consequently, it is sufficient to say that a bounded input signal will always produce a bounded output signal if all the poles lie inside the unit circle.

The zeros on the other hand, are not constrained by this requirement, and as a consequence may lie anywhere on z-plane, since they do not directly affect system stability. Therefore, a system stability analysis may be undertaken by firstly calculating the roots of the transfer function (i.e., roots of the numerator and denominator polynomials) and then plotting the corresponding poles and zeros upon the z-plane.

An interesting situation arises if any poles lie on the unit circle, since the system is said to be marginally stable, as it is neither stable or unstable. Although marginally stable systems are not BIBO stable, they have been exploited by digital oscillator designers, since their impulse response provides a simple method of generating sine waves, which have proved to be invaluable in the field of telecommunications.

Biquad IIR filters

The IIR filter implementation discussed herein is said to be biquad, since it has two poles and two zeros as illustrated below in Figure 2. The biquad implementation is particularly useful for fixed point implementations, as the effects of quantization and numerical stability are minimised. However, the overall success of any biquad implementation is dependent upon the available number precision, which must be sufficient enough in order to ensure that the quantised poles are always inside the unit circle.

Direct Form I (biquad) IIR filter realization and transfer function.; Direct Form; Biquad filter

Figure 2: Direct Form I (biquad) IIR filter realization and transfer function.

Analysing Figure 2, it can be seen that the biquad structure is actually comprised of two feedback paths (scaled by \(\small a_1\) and \(\small a_2\)), three feed forward paths (scaled by \(\small b_0, b_1\) and \(\small b_2\)) and a section gain, \(\small K\). Thus, the filtering operation of Figure 1 can be summarised by the following simple recursive equation:

\(\displaystyle y(n)=K\times\Big[b_0 x(n) + b_1 x(n-1) + b_2 x(n-2)\Big] – a_1 y(n-1)-a_2 y(n-2)\)

Analysing the equation, notice that the biquad implementation only requires four additions (requiring only one accumulator) and five multiplications, which can be easily accommodated on any Cortex-M microcontroller. The section gain, \(\small K\) may also be pre-multiplied with the forward path coefficients before implementation.

A collection of Biquad filters is referred to as a Biquad Cascade, as illustrated below.

Biquad Cascade; Biquad filter

The ASN Filter Designer can design and implement a cascade of up to 50 biquads (Professional edition only).

Floating point implementation

When implementing a filter in floating point (i.e. using double or single precision arithmetic) Direct Form II structures are considered to be a better choice than the Direct Form I structure. The Direct Form II Transposed structure is considered the most numerically accurate for floating point implementation, as the undesirable effects of numerical swamping are minimised as seen by analysing the difference equations.

Direct Form II Transposed strucutre, transfer function and difference equations; IIR Filters; Biquad Filters

Figure 3 – Direct Form II Transposed strucutre, transfer function and difference equations

The filter summary (shown in Figure 4) provides the designer with a detailed overview of the designed filter, including a detailed summary of the technical specifications and the filter coefficients, which presents a quick and simple route to documenting your design.

The ASN Filter Designer supports the design and implementation of both single section and Biquad (default setting) IIR filters.

Biquad filter ASN Filter Designer DSP
Figure 4: detailed specification.

FIR definition

Returning the IIR’s linear constant coefficient difference equation, i.e.

\(\displaystyle
y(n)=\sum_{k=0}^{q}b_kx(n-k)-\sum_{k=1}^{p}a_ky(n-k)
\)

Notice that when we set the \(\small a_k\) coefficients (i.e. the feedback) to zero, the definition reduces to our original the FIR filter definition, meaning that the FIR computation is just based on past and present inputs values, namely:

\(\displaystyle
y(n)=\sum_{k=0}^{q}b_kx(n-k)
\)

Implementation

Although several practical implementations for FIRs exist, the direct form structure and its transposed cousin are perhaps the most commonly used, and as such, all designed filter coefficients are intended for implementation in a Direct form structure.

The Direct form structure and associated difference equation are shown below. The Direct Form is advocated for fixed point implementation by virtue of the single accumulator concept.

\(\displaystyle y(n) = b_0x(n) + b_1x(n-1) + b_2x(n-2) + …. +b_qx(n-q) \)

Direct form; Direct form structure

The recommended (default) structure within the ASN Filter Designer is the Direct Form Transposed structure, as this offers superior numerical accuracy when using floating point arithmetic. This can be readily seen by analysing the difference equations below (used for implementation), as the undesirable effects of numerical swamping are minimised, since floating point addition is performed on numbers of similar magnitude.

\(\displaystyle \begin{eqnarray}y(n) & = &b_0x(n) &+& w_1(n-1) \\ w_1(n)&=&b_1x(n) &+& w_2(n-1) \\ w_2(n)&=&b_2x(n) &+& w_3(n-1) \\ \vdots\quad &=& \quad\vdots &+&\quad\vdots \\ w_q(n)&=&b_qx(n) \end{eqnarray}\)

Direct form Transposed

What have we learned?

Digital filters are divided into the following two categories:

  • Infinite impulse response (IIR)
  • Finite impulse response (FIR)

IIR (infinite impulse response) filters are generally chosen for applications where linear phase is not too important and memory is limited. They have been widely deployed in audio equalisation, biomedical sensor signal processing, IoT/IIoT smart sensors and high-speed telecommunication/RF applications.

FIR (finite impulse response) filters are generally chosen for applications where linear phase is important and a decent amount of memory and computational performance are available. They have a widely deployed in audio and biomedical signal enhancement applications.

ASN Filter Designer provides engineers with everything they need to design, experiment and deploy complex IIR and FIR digital filters for a variety of sensor measurement applications. These advantages coupled with automatic documentation and code generation functionality allow engineers to design and validate an IIR/FIR digital filter within minutes rather than hours.

 

 

Download demo now

Licencing information

In ECG signal processing, the Removal of 50/60Hz powerline interference from delicate information rich ECG biomedical waveforms is a challenging task! The challenge is further complicated by adjusting for the effects of EMG, such as a patient limb/torso movement or even breathing. A traditional approach adopted by many is to use a 2nd order IIR notch filter:

\(\displaystyle H(z)=\frac{1-2cosw_oz^{-1}+z^{-2}}{1-2rcosw_oz^{-1}+r^2z^{-2}}\)

where, \(w_o=\frac{2\pi f_o}{fs}\) controls the centre frequency, \(f_o\) of the notch, and \(r=1-\frac{\pi BW}{fs}\) controls the bandwidth (-3dB point) of the notch.

What’s the challenge?

As seen above, \(H(z) \) is simple to implement, but the difficulty lies in finding an optimal value of \(r\), as a desirable sharp notch means that the poles are close to unit circle (see right).

In the presence of stationary interference, e.g. the patient is absolutely still and effects of breathing on the sensor data are minimal this may not be a problem.

However, when considering the effects of EMG on the captured waveform (a much more realistic situation), the IIR filter’s feedback (poles) causes ringing on the filtered waveform, as illustrated below:

Contaminated ECG with non-stationary 50Hz powerline interference (FIR filtering), ECG sigal processing, ECG DSP, ECG measurement

Contaminated ECG with non-stationary 50Hz powerline interference (IIR filtering)

As seen above, although a majority of the 50Hz powerline interference has been removed, there is still significant ringing around the main peaks (filtered output shown in red). This ringing is undesirable for many biomedical applications, as vital cardiac information such as the ST segment cannot be clearly analysed.

The frequency reponse of the IIR used to filter the above ECG data is shown below.

IIR notch filter frequency response, ECG signal processing, ECG DSP, ECG  measurement

IIR notch filter frequency response

Analysing the plot it can be seen that the filter’s group delay (or average delay) is non-linear but almost zero in the passbands, which means no distortion. The group delay at 50Hz rises to 15 samples, which is the source of the ringing – where the closer to poles are to unit circle the greater the group delay.

ASN FilterScript offers designers the notch() function, which is a direct implemention of H(z), as shown below:

ClearH1;  // clear primary filter from cascade
ShowH2DM;   // show DM on chart

interface BW={0.1,10,.1,1};

Main()

F=50;
Hd=notch(F,BW,"symbolic");
Num = getnum(Hd); // define numerator coefficients
Den = getden(Hd); // define denominator coefficients
Gain = getgain(Hd); // define gain

Savitzky-Golay FIR filters

A solution to the aforementioned mentioned ringing as well as noise reduction can be achieved by virtue of a Savitzky-Golay lowpass smoothing filter. These filters are FIR filters, and thus have no feedback coefficients and no ringing!

Savitzky-Golay (polynomial) smoothing filters or least-squares smoothing filters are generalizations of the FIR average filter that can better preserve the high-frequency content of the desired signal, at the expense of not removing as much noise as an FIR average. The particular formulation of Savitzky-Golay filters preserves various moment orders better than other smoothing methods, which tend to preserve peak widths and heights better than Savitzky-Golay. As such, Savitzky-Golay filters are very suitable for biomedical data, such as ECG datasets.

Eliminating the 50Hz powerline component

Designing an 18th order Savitzky-Golay filter with a 4th order polynomial fit (see the example code below), we obtain an FIR filter with a zero distribution as shown on the right. However, as we wish to eliminate the 50Hz component completely, the tool’s P-Z editor can be used to nudge a zero pair (shown in green) to exactly 50Hz.

The resulting frequency response is shown below, where it can be seen that there is notch at exactly 50Hz, and the group delay of 9 samples (shown in purple) is constant across the frequency band.

FIR  Savitzky-Golay filter frequency response, ECG signal processing, ECG DSP, ECG measurement

FIR  Savitzky-Golay filter frequency response

Passing the tainted ECG dataset through our tweaked Savitzky-Golay filter, and adjusting for the group delay we obtain:

Contaminated ECG with non-stationary 50Hz powerline interference (FIR filtering), ECG signal processing, ECG digital filter, ECG filter designa

Contaminated ECG with non-stationary 50Hz powerline interference (FIR filtering)

As seen, there are no signs of ringing and the ST segments are now clearly visible for analysis. Notice also how the filter (shown in red) has reduced the measurement noise, emphasising the practicality of Savitzky-Golay filter’s for biomedical signal processing.

A Savitzky-Golay may be designed and optimised in ASN FilterScript via the savgolay() function, as follows:

ClearH1;  // clear primary filter from cascade

interface L = {2, 50,2,24};
interface P = {2, 10,1,4};

Main()

Hd=savgolay(L,P,"numeric");  // Design Savitzky-Golay lowpass
Num=getnum(Hd);
Den={1};
Gain=getgain(Hd);

Deployment

This filter may now be deployed to variety of domains via the tool’s automatic code generator, enabling rapid deployment in Matlab, Python and embedded Arm Cortex-M devices.

Generations have been entertained by the gadgets and future technology portrayed in Sci-fi series, such as Star Trek, but is it all science fiction?


The Tricorder

One device that intrigued me for years was the so called ‘Tricorder’ and the ability for the doctor to read a person’s vital life signs or VLS (e.g. heart beat and respiration) with a handheld device from about a metre away.


Back to the 21st century!

With advances in radar technology over the last few years, a few chip manufacturers are now producing affordable radar devices suitable for biomedical VLS measurement. Radar technology that used to cost thousands of Euros, and was primarily aimed a military technology, is now available for a few hundred Euros, making it viable for home medical products.


Sounds great, but what can UWB (ultrawide band) pulse doppler radar do?

  • millimetre accuracy: allowing for detection of the smallest changes, such as respiration and heart rate from several metres away.
  • Penetrates duvets, blankets and clothes: The virtue of the small wavelength ensures accurate detection of humans in a bed or sitting in chair watching TV or reading a book.
  • Penetrates Walls and doors: Tracking of VLS and movement when radar is mounted within a ceiling or behind a wall – no ugly module on the wall!
  • High sensitivity: able to see the VLS from tiny premature babies.
  • Ultra-safe technology: RF emission is 0.01% the energy typically found in a household WiFi router – meaning that prolonged exposure will have no detrimental effects on human health.

After receiving a request from a client about monitoring the health of an eldery person living alone or in a nursing home, we decided to conduct a few tests of our own to see what was possible with this new technology.


VLS of a subject lying in bed

After building a demonstrator, and placing the radar sensor about 1 metre from a subject sleeping (similar to the doctor in Star Trek), we obtained the following waveform:

Captured VLS data captured from a UWB radar – containing both heart beat and respiration information

Wow! Was our initial reaction to the test data – this is millimetre movement through a duvet! Notice how slow the biomedical signal is, as an average adult’s respiration rate at rest is about 12bpm, requiring relatively long data acquisition times (tens of seconds) for meaningful data analysis.


What’s the respiration and heart beat?

Passing the signal through our algorithm, we could easily estimate the respiration (RR) and heartbeat (HB) from the plot (see the two red squares on the two peaks). However, in order to be objective, we attached a clip-on pulse oximeter to the subject finger, and as seen they matched very well.

What does this all mean for me?
Contactless VLS measurement for home use is closer than you think, and is certainly not science-fiction anymore. This technology opens up many possibilities for monitoring when normal sensors are infeasible, such as premature babies, patients with dementia and even sleep trend analysis. We’ll improve our algorithm in order to make it more robust and faster, but as seen our results are very promising indeed, and open up the possibility of contactless vital life signs (VLS) measurements for many practical applications!