﻿ Low frequency Butterworth and optimal Wiener ECG filters_Circuit Diagram World ﻿
Position: Index > Basic Circuit >

# Low frequency Butterworth and optimal Wiener ECG filters

2015-01-18 23:33
Declaration:We aim to transmit more information by carrying articles . We will delete it soon, if we are involved in the problems of article content ,copyright or other problems.

Regular ad hoc filters don’t guarantee optimal signal filtering as there is no any criteria that evaluates filter characteristics. Usually filter parameters are calculated empirically and filtering is done by best results.
In order to avoid such shortage there are optimal filters used where parameters are optimized by some criteria. The main idea of optimal filtering is to give bigger weight coefficients to signal spectra parts where signal noise has less power and true signal spectral components has bigger power.
Lets project simple Butterworth filter that will be used as comparative filter to optimal Wiener DSP filter.
Butterworth filter transfer characteristics:

Where N a€“ indicates filter Tap number. I will skip description of butterworth filter for now as main idea is constructing optimal wiener filter. Butterworth filter characteristics are pretty plain:

Butterworth characteristics with bandpass of 70Hz

Main disadvantage of Butterworth filter is that signal is distorted on filter output. If you want minimal signal distortions it is better to use optimal Wiener filter. Filter chart looks as follows:

As you can see to make this filter functional wee need additional conditions like signal model reference. Filter coefficients are calculated usingMMSE.rootmeansquareerror. Lets see how it works.
Filter input is ideal signal plus noise:

and filter output depends on filter weights of coefficients wk

Further wee need to optimize filter coefficients wkin order to meet minimal signal error e(n):

As we can see minimalrootmeansquareerror depends on ideal signal model construction.
In frequency domain optimal Wiener filter coefficients are calculated by using Wiener – Hoph equation:

Where S?·(w) a€“ noise power spectral density;

Sd(w) a€“ signal model power spectral density.

Where R(m)- noise self correlation function.

Lets see what results wee get. First of all we have made Butterworth low frequency filter. Signal sampling frequency were used a€“ 1000Hz. We have made couple examples when filter tap numbers are 4 and 8 and cutoff frequencies are 30 and 70Hz:

butterworth_tap4_cut30Hz.jpg

butterworth_tap4_cut70Hz.jpg

butterworth_tap8_cut30Hz.jpg

butterworth_tap8_cut70Hz.jpg

It is really hard to assemble good parameters for Butterwoth filter because there is no good criteria to value filter accuracy. So following step is optimal filter construction.

First of all we need signal model reference:

According to previous formulas we calculate filter coefficients depending on signal and signal model spectral densities:

a) Signal noise spectral density; b) real signal spectral density; c) signal modelspectral density; d) Filter coefficients dependency on frequency.

When we have Wiener-Hoph equation with all required coefficients we can filter ECG signal. First of all filtering is done in frequency domain. Lets calculate signal Fourier transformation, then filtering is nothing more than multiplying signal Fourier transformation by filter coefficients. Results are in following picture:

a) real signalo Fourier transformation; b) filtered signalFourier transformation; c) real signal; d) filtered signal

Of course you can filter signal in time domain. For this you have to make reverse Fourier transformation of filter coefficients. If wee take 30 of filter coefficients we get:

Signal filtered in time domain using 30 of filter coefficients.

Now we can compare filtered signals by using Butterworth and Wiener filters. One of comparative methods can be calculation of group delay. In following results we can see that using Butterworth filter group delay is biggest at 70Hz while Wiener filters resulting group delay oscillates around 0.

a) Butterworth filter group delay; b) Wiener filter group delay

According to results it is really hard to say which filtering method is better. We can see, that using tap 8 and 30Hz Butterworth filter it gives more clear Q segment while using Wiener filter this part of signal is hardly noticeable. Of course signal model is quite rough, this results in less precise filter coefficients.

%%%%%%%%Matlabsource%%%%%%%%%%%%%%%%%%%%%

%Optimal ECG signal filtering using Wiener filter”

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

close all;

clear all;

clc;

%1. Projecting Butterworth filter

Fp=70;

Fd=1000;

eil=8;

%Cut-Off

ecg=ecg/max(ecg(2100:2800));

figure

[B,A]=butter(eil/2,Fp/(Fd/2));

N=length(ecg);

timeall=N/Fd;

tecg=0:1/Fd:timeall;

plot(tecg(1:700),ecg(2101:2800));

grid on;

title(‘measured signal and filtered with Butterworth filter’);

xlabel(‘t, s’);

ylabel(‘Amplitude, mV’);

hold on;

fecg1 = filter(B,A,ecg);

for i=1:N

ecg1rev(i)=fecg1(N 1-i);

end

fecg2rev = filter(B,A,ecg1rev);

for i=1:N

fecg(i)=fecg2rev(N 1-i);

end

plot(tecg(1:700),fecg(2101:2800),’red’);

xlabel(‘t, s’);

ylabel(‘Amplitude, mV’);

legend(‘Real signal and filtered signal’);

hold off;

%Butterworth characteristics

figure;

freqz(B,A,700,Fd);

figure

grpdelay(B,A,700,Fd);

%Group delay

clear des_x; clear des_y;

des_x=[]; des_y=[];

%Model Interpolation

des_x=xd-xd(1);

des_y=max(ecg(2100:2800))*yd;

% Aproximation with cubic spine

yi=interp1(des_x,des_y,xi,’cubic’);

figure

plot(des_x,des_y,’o',xi,yi);

grid on;

xlabel(‘t, s’);

ylabel(‘Amplitude, mV’);

des_x = des_x(1:13);

des_y = des_y(1:13);

figure(3);

plot(des_x,des_y,’*',des_x,des_y);

hold on;

grid on;

xlabel(‘Time, ms’);

ylabel(‘Amplitude’);

axis tight;

hold off;

ans = 1;

end

%

x1=trm;

L1=x1(2)-x1(1);

L2=x1(4)-x1(3);

L3=x1(6)-x1(5);

L4=x1(8)-x1(7);

L=[L1; L2; L3; L4];

Lmin=min(L);

ecgtr=[ecg(x1(2):x1(2) Lmin), ecg(x1(3):x1(3) Lmin), ...

ecg((x1(6)-Lmin):x1(6)),ecg(x1(7):x1(7) Lmin)];

Tr=mean(ecgtr,2);

Trvid=mean(Tr,1);

Tr=Tr-Trvid;

%noise power spectral density

[Ptr,f]=periodogram(Tr,[],514,Fd); Q2

%periodogram(Tr,[],512,Fd);

figure

psdplot(Ptr,f,’Hz’);

figure

subplot (2,2,1);psdplot(Ptr,f,’Hz’);

%Real signal spectral density

ecg=ecg-Trvid;

ecg=ecg/max(ecg(2100:2800));

[Psig,f]=periodogram(ecg(2100:2800),[],512,Fd); Q2

subplot(2,2,2);psdplot(Psig,f,’Hz’);

%Model spectral density

[Pmod,f]=periodogram(yi,[],514,Fd); Q2

subplot(2,2,3);psdplot(Pmod,f,’Hz’);

%Wiener coefficients

for ii=1:length(Pmod)

W(ii)=1/(1 (Ptr(ii)/Pmod(ii)));

end

ii=1:length(Pmod);

subplot(2,2,4); plot(ii,W);

grid on;

xlabel(‘f, Hz’);

ylabel(‘W’);

%Low freq filtering

ecgfft = fft(ecg(2100:2800),700);

figure

subplot(2,2,1);plot(abs(ecgfft(1:258)))

grid on;

xlabel(‘f, Hz’);

ylabel(‘Amplitude’);

for ii=1:(length(Pmod))

ecgfd(ii)=W(ii)*ecgfft(ii);

end

subplot(2,2,2);plot(abs(ecgfd));

grid on;

xlabel(‘f, Hz’);

ylabel(‘Amplitude’);

subplot(2,2,3);plot(tecg(1:700),ecg(2101:2800));

grid on;

xlabel(‘t, s’);

ylabel(‘Amplitude, mV’);

ecgiff= (ifft(ecgfd(1:length(Pmod)), 700));

ecgiff=real(ecgiff)/real(max(ecgiff)) Trvid/4;

subplot(2,2,4);plot(tecg(1:700),real(ecgiff));

grid on;

xlabel(‘t, s’);

ylabel(‘Amplitude, mV’);

figure

plot(tecg(1:700),ecg(2101:2800),’red’);

grid on;

xlabel(‘t, s’);

ylabel(‘Amplitude, mV’);

%Filtering in time domain

koef=abs(ifft(W));

ecgtd = filter(koef(1:30),1, ecg);

ecgtd=ecgtd/max(ecgtd);

figure

plot(tecg(1:700),ecgtd(1:700))

grid on;

xlabel(‘t, s’);

ylabel(‘Amplitude, mV’);

figure

subplot(1,2,1);grpdelay(B,A,700,Fd);

subplot(1,2,2);grpdelay(koef(1:30),1,700,Fd);

Files:?? mat1.mat and trm.mat files