CN102940489A - Filter designing method of weak electrophysiology signal and filtering method - Google Patents

Abstract

The invention discloses a filter designing method of a weak electrophysiology signal and a filtering method. The filter designing method comprises the following steps of supposing that d is equal to (N-1)/2 according to an interpolation factor D and an prolongation period number N; initializing a filter coefficient h(n) to be a matrix with the length of (2d-1) and all the elements being 1, and carrying out upper sampling of a D point on a filter h, namely interpolating (D-1) zero points behind each coefficient h(n) of the filter h; calculating all nonzero coefficients of the filter coefficient h(n) according to the formula; and modifying and calculating a coefficient of a central point so as to design the filter h. The filtering method comprises the following steps of (1), inputting n sample point sets, namely adjacent two sample point interval D sampling periods of the sample point sets; (2) carrying out accumulation on the n sample points of the sample point sets to acquire the average; and (3) multiplying the central sample point xn by a (1-3/N) value and subtracting an average value calculated by the step (2), namely filter output y(n) at the existing moment, wherein the N is the prolongation period number. Compared with the prior art, the filter designing method of the weak electrophysiology signal and the filtering method, disclosed by the invention, have the advantages of reduced utilization of multiplying units, low calculation complexity, short filtering time and flexible parameter configuration.

Description

A kind of filter design method of faint electricity physiological signal and filtering method

Technical field

The present invention relates to medical signals digital processing technology field, be specifically related to remove the filtering technique that baseline drift interference that faint electricity physiological signal often is subject to and power frequency are disturbed with China.

Background technology

Clinic system carries out long-term electrocardiogram (ECG) monitoring by computer-aided tool to high risk cardiac, and requiring has quick, real-time faint electricity physiological signal analytical method, and high-precision heart disease automatic diagnosis [1] is provided.Like this, faint electricity physiological signal is carried out pretreatment filtering and just seem most important.

The most commonly power frequency interference and base band are disturbed in the multiple external disturbance that faint electricity physiological signal is subject to.Power frequency is (China is 50Hz) of fixing, however baseline drift change indefinitely, baseline drift is common external disturbance, people's breathing, human body slight mobile etc. all can be so that faint electricity physiological signal produces drift.The characteristics of baseline drift are low frequency characteristics, need to eliminate by the means of high-pass filtering.

The removal of baseline drift both can realize by analog filtering, also can realize by digital filtering.Analog filtering is traditional implementation, and its advantage is that time delay is little, cost is low, but filtering accuracy is low, very flexible; Recently, the removal of baseline drift progressively tends to the digital filtering realization, and the digital filtering degree of accuracy is high, motility is high, but computation complexity and cost are high.Thereby the cost (being multiplier and adder resource particularly) and the computation complexity that how to reduce digital filtering are the key problems of eliminating the electrocardio drift.

And the cut-off frequency of baseline drift wave filter really fixes on academia and engineering circles and has dispute and contradiction always, if high pass cut off frequency is too high, means that then baseline drift eliminates too thoroughly, and this easily also filters out T ripple low-frequency component and affects diagnosis; If cut-off frequency is too low, mean that then baseline drift eliminates very slightly, easily the baseline drift that should remove is not filtered out, kept.

Nineteen ninety U.S. AHA group is pointed out, the longest RR interval,, the low-limit frequency (typically referring to the heart rate under the per minute heartbeat 40) with faint electricity physiological signal was corresponding, this just means that high pass cut off frequency is the digital filter of 0.67Hz, can satisfy most baseline drift to eliminating [2].

There are a lot of scholars that the filtering method of the noise of removing faint electricity physiological signal is studied all the time both at home and abroad, proposed a lot of filtering methods.Such as foreign scholar H.Govaerts a kind of method that moving average filter satisfies the requirement of passband and stopband of repeatedly using is proposed.Document [3] then proposes the wave filter [4] that a kind of FIR wave filter with Kaiser window design and rolling average combine; Domestic also have a lot of relevant documents [5-9].Although there has been so many wave filter to meet the demands, these methods have shortcoming separately.Some amounts of calculation are large, are unsuitable for the realization of mini system; Some parameters arrange dumb, such as the selection of the wavelet basis function in the document [6], can only choose with experience by experiment, easily affect the effect of filtering; What have need to carry out pretreatment to faint electricity physiological signal, has caused the delay of signal analysis.

The present invention proposes a kind of digital filtering method based on interpolation, and than the conventional digital filtering method, its computation complexity is very low, and the computational resource that expends considerably less (mainly being multiplier resources) is suitable for portable electrocardiograph and realizes.

List of references:

[1]Kligfeld,P.,L.Gettes,J.Bailey,R.Childers,B.Deal,W.Hancock,G.Herpen,J.Kors,P.Macfarlane,D.Mirvis,O.Pahlm,P.Rautaharju,and?G.Wagner.Recommendations?for?the?standardization?andinterpretation?of?the?electrocardiogram,Part?I.The?electrocardiogram?and?its?technology.J.Am.Coll.Cardiol.49:1109–1127,2007.doi:10.1016/j.jacc.2007.01.024.

[2]Serafim?Tabakov,Ivo?Iliev,Vessela?Krasteva.Online?Digital?Filter?and?QRS?Detector?Applicable?inLow?Resource?ECG?Monitoring?Systems[J].Annals?of?Biomedical?Engineering,2008,36(11):1805-1815.

[3]Faes,Th.,H.Govaerts,B.Tenvoorde,and?O.Rompelman.Frequency?synthesis?of?digital?flters?basedon?repeatedly?applied?unweighed?moving?average?operations.Med.Biol.Eng.Comp.32:698-701,1994.doi:10.1007/BF02524254.

[4]Bai,Y.,W.Chu,Ch.Chen,Y.Lee,Y.Tsai,and?Ch.Tsai.The?combination?of?Kaiser?window?andmoving?average?for?the?low-pass?fltering?of?the?remote?ECG?signals.In:Proc.17th?IEEE?Symposium?onComputer-Based?Medical?Systems,273,2004.

[5] Zou Bo, week is far away. the digital filtering algorithm [J] that a kind of filtering baseline drift and power frequency are disturbed. and Journal of Shenzhen Polytechnic, 2005,4 (4): 3-5

[6] Li An, Li Dong. the application [J] of wavelet analysis method in faint electricity physiological signal digital filtering. Computer Simulation, 2001,18 (6): 70-73

[7] Chen Jing, Qian Lifeng, to thoroughfare, Zhou Linhai. real-time signal-processing method research [J] in the portable cardiac monitoring. biomedical engineering's progress, 2009,30 (4): 200-204

[8] Chen Tianhua, Chen Qian. a kind of digital filtering method [J] of eliminating faint electricity physiological signal noise. World Science technology, 2005,7 (1): 123-126

[9] in tender, bent ripple. the design [J] of the digital filter during faint electricity physiological signal is processed. information technology, 2009, the fifth phase: 61-63

Summary of the invention

Based on above-mentioned technical problem, the present invention proposes a kind of filter design method and filtering method of faint electricity physiological signal, disturb the method for designing proposed wave filter for the baseline drift of removing faint electricity physiological signal.

The filter design method of a kind of faint electricity physiological signal of the present invention may further comprise the steps:

Step 1 according to interpolation factor D and continuation periodicity N, makes d=(N-1)/2, and d represents the group delay of wave filter;

Step 2 is initialized as all 1's matrix that length is 2d-1 with filter coefficient h (n), and wave filter h is carried out the up-sampling that D is ordered, and namely inserts D-1 zero point in each coefficient h (n) back of wave filter h; And

Step 3, according to above-mentioned parameter, the expression formula of filter coefficient h (n) is

N, m represent integer,

Calculate all nonzero coefficients of filter coefficient h (n) according to above-mentioned formula; Step 4, center sampling point x _nThe value that multiply by (1-3/N) deducts the meansigma methods that step 2 is calculated, and is the filtering output y (n) of current time, and N is the continuation periodicity, thereby designs this wave filter h.

Described step 2 also comprises the step that the point with the last interpolation of ripple device coefficient h (n) deletes.

The filtering sampling point of the method choose employing: with x _nAnd the x that is adjacent _N-1And x _N+1Centered by sampling point, in the sampling period take D as data decimation, select symmetrically respectively first group of sampling point data x _n, x _N-D, x _N+DSecond group of sampling point data x _N-1, x _N-D-1, x _N+D-1The 3rd group of sampling point data x _N+1, x _N-D+1, x _N+D+1Rule choose.

Described interpolation factor D is set to sample frequency f _sWith power frequency interfering frequency f _pRatio D=f _s/ f _p

Described continuation periodicity N adopts enumerative technique to try to achieve, and specifically may further comprise the steps:

Make ω=2 π f/f _s/ _F=0.67Be updated to $H\left(e^{\mathrm{j\ω}}\right)=\frac{N-2}{N}-\frac{\sin \left[(\frac{N}{2}-1)\mathrm{D\ω}\right]}{N\sin (\mathrm{D\ω}/2)}$ In, w represents frequency;

Change the N value of enumerating, the computer Aided Design is drawn out H (e ^{J ω}) with the change curve of cycle N;

From change curve, find out H (e ^{J ω}) the most approaching

Corresponding odd number N value is the continuation periodicity N that will arrange.

The filtering method of a kind of faint electricity physiological signal of the present invention may further comprise the steps:

Step 1, n sampling point collection of input, namely $\{x_{n-\frac{N-1}{2}D},\·\·\·,x_{n-D},x_n,x_{n+D},\·\·\·,x_{n+\frac{N-1}{2}D}\},$ Two adjacent sampling point interval D sampling periods of this sampling point collection;

Step 2, n the sampling point that sampling point is concentrated add up and are averaging;

Step 3, with center sampling point x _nThe value that multiply by (1-3/N) deducts the meansigma methods that step 2 is calculated, and is the filtering output y (n) of current time, and N is the continuation periodicity.

Compared with prior art, the efficient filter of the faint electricity physiological signal baseline drift of removal that the present invention proposes if be applied to the practical medical engineering field, can produce following beneficial effect:

1, usefulness multiplier seldom can be realized the filtering to electrocardiosignal.

Can find out from the filtering flow process of the faint electricity physiological signal of front Fig. 3, much no matter periodicity N gets, all only need to carry out twice multiplication, this has greatly reduced the number of multiplier.

2, computation complexity is low.

Amount of calculation of the present invention mainly is the calculating of addition.Because the multiplier number is few, adds calculative N-1 sub-addition, amount of calculation reduces greatly.

3, the filtering time is short

Because used multiplier is few, computation complexity is low, therefore than conventional filter, the calculating that filtering spends only has twice multiplying, therefore the filter filtering time that the present invention invents is very short.

4, flexible setting for parameters.

Power frequency interfering frequency f among the present invention _p=50Hz is known constant, generally sample frequency f in practical problem _sAlso be known, therefore interpolation factor can be according to D=f _s/ f _pCalculate.Continuation periodicity N then realizes that according to foregoing computer Aided Design enumerative technique arranges, thereby D and the N of coupling are set so that parameter arrange very flexibly.

Stronger to the noise capacity of resisting disturbance below the 0.67Hz.

The example 2 of front is verified, and this wave filter is also can filtering clean to the noise below the 0.67Hz.

Description of drawings

Fig. 1 is the method for designing flow chart of wave filter of the present invention;

Data chose when Fig. 2 was the different sampling point filtering of the present invention;

Fig. 3 is the filtering flow process of the faint electricity physiological signal of the present invention;

Fig. 4 is the relation curve of transfer function H of the present invention and periodicity N;

Fig. 5 is the output characteristic curve 5 of wave filter of the present invention: wherein (a) is the transfer curve of wave filter for the coefficient of wave filter (b);

Fig. 6 is that amplify the coefficient of wave filter of the present invention and the part of transmission curve;

Fig. 7 is the time domain waveform before and after the faint electricity physiological signal filtering of 900Hz;

Fig. 8 is the time domain waveform before and after the Filtering of ECG Signal of making an uproar adding of 900Hz;

Fig. 9 is the human body location drawing that respectively leads.

The specific embodiment

Below in conjunction with accompanying drawing the specific embodiment of the present invention is described in detail.

The present invention proposes following wave filter h method for designing: the supposition sampling rate is f _s, interpolation factor is D, the continuation periodicity is N, make d=(N-1)/2, then designing a length is FIR (FiniteImpulse Response, FIR) the wave filter h of 2 (d-1) D+1=(N-3) D+1, its coefficient h (n) as the formula (1)

Wave filter h coefficient h (n) forming process corresponding with formula (1) as shown in Figure 1.

Design of filter flow process as shown in Figure 1: at first according to sample frequency f _sWith power frequency interfering frequency f _pCalculate interpolation factor D, adopt enumerative technique continuation periodicity N to be set, this seasonal d=(N-1)/2; Filter coefficient h (n) is initialized as all 1's matrix that length is 2d-1, then h is carried out the up-sampling (being interpolation) that D is ordered, insert D-1 zero point in each coefficient h (n) back of h like this, simultaneously in order to reduce the length of h, the point of the last interpolation of filter coefficient is deleted, then the coefficient expression formula according to wave filter calculates all nonzero coefficient h (n), revises at last the coefficient of computer center's point, thereby designs this wave filter.

For example, make D=3, N=5, form flow process according to filter parameter shown in Figure 1, the coefficient that can calculate first this wave filter be h (n)=1/500-1/500-1/5}, and then revise to determine the coefficient of central point, thus obtain h (n)={ 0.2000000.400000-0.2000}.From the coefficient h (n) of wave filter, can find out, insert D-1 null value between per two nonzero coefficients, and except the central point coefficient, other zero coefficient values all equate.

Filter factor h (n) forming process from above-mentioned design of filter flow process can be found out, among the filter coefficient h (n) a lot of null values are arranged, and these zero valued coefficients are inoperative to filtering, only use nonzero coefficient when namely carrying out the filtering of faint electricity physiological signal.Can save the memory data output of (D-1)/D when this means filtering.

Baseline drift and the power frequency removed in the faint electricity physiological signal of human body are disturbed, can be good by the data decimation with sampling point filtering shown in Figure 2 first, then carry out filtering.As shown in Figure 2, if with x _nCentered by sampling point carry out filtering, take D as the data break distance, select symmetrically respectively first group of sampling point data x _n, x _N-D, x _N+D;In like manner, second group with x _N-1Centered by sampling point carry out filtering, sampling point data x _N-1, x _N-D-1, x _N+D-1The 3rd group with x _N+1Centered by sampling point carry out filtering sampling point data x _N+1, x _N-D+1, x _N+D+1

Filtering flow process of the present invention may further comprise the steps:

1, N sampling point collection of input, namely $\{x_{n-\frac{N-1}{2}D},\·\·\·,x_{n-D},x_n,x_{n+D},\·\·\·,x_{n+\frac{N-1}{2}D}\},$ Two adjacent sampling point interval D sampling periods of this sampling point collection;

2, N the sampling point of sampling point being concentrated adds up and is averaging;

3, with center sampling point x _nThe value that multiply by (1-3/N) deducts the meansigma methods that step (ii) is calculated, and is current n filtering output y (n) constantly.

As shown in Figure 3, choose 5 sampling points and illustrate the filtering flow process of faint electricity physiological signal:

As can be seen from Figure 3, very simple of the filtering of signal only need carry out N-1 sub-addition and 2 multiplication, the multiplier of consumption seldom, therefore computation complexity is very low, filtration efficiency improves greatly.

The Transmission characteristics of the wave filter that, is obtained by above-mentioned method for designing

Filter coefficient h (n) is carried out Fourier transformation get final product to get frequency response function H (e ^{J ω}), this is the important indicator of wave filter transmission characteristic, namely formula (1) is asked for its Fourier transformation:

$H\left(e^{\mathrm{j\ω}}\right)=\underset{n=-(d-1)D}{\overset{(d-1)D}{\mathrm{\Σ}}}h\left(n\right)e^{-\mathrm{jn\ω}}$

$=1-3/N-\frac{(e^{\mathrm{jD\&omega;}}+e^{-\mathrm{jD\&omega;}})+({\mathrm{<figure-callout\; id="2"\; label="e\; j"\; filenames="HDA00002494477700041.png,HDA00002494477700051.png"\; state="\{\{state\}\}">e</figure-callout>}}^j+e^{-\mathrm{<figure-callout\; id="2"\; label="j"\; filenames="HDA00002494477700041.png,HDA00002494477700051.png"\; state="\{\{state\}\}">j</figure-callout>}2\mathrm{D\&omega;}})+\&CenterDot;\&CenterDot;\&CenterDot;+(e^{j(d-1)\mathrm{D\&omega;}}+e^{-j(d-1)\mathrm{D\&omega;}})}{N}$

$=1-3/N-\frac{2[\cos \mathrm{D\&omega;}+\mathrm{<figure-callout\; id="2"\; label="cos"\; filenames="HDA00002494477700041.png,HDA00002494477700051.png"\; state="\{\{state\}\}">cos</figure-callout>}2\mathrm{D\&omega;}+\&CenterDot;\&CenterDot;\&CenterDot;\cos (d-1)\mathrm{D\&omega;}]}{N}$

$=1-3/N-\frac{2\sin D\mathrm{\&omega;}\cos \mathrm{D\&omega;}+2\sin D\mathrm{\&omega;}\mathrm{<figure-callout\; id="2"\; label="cos"\; filenames="HDA00002494477700041.png,HDA00002494477700051.png"\; state="\{\{state\}\}">cos</figure-callout>}2\mathrm{D\&omega;}+2\sin D\mathrm{\&omega;}\cos (d-1)\mathrm{D\&omega;}}{N\sin \mathrm{D\&omega;}}$

$=1-3/N-\frac{\mathrm{<figure-callout\; id="2"\; label="sin"\; filenames="HDA00002494477700041.png,HDA00002494477700051.png"\; state="\{\{state\}\}">sin</figure-callout>}2\mathrm{D\&omega;}+(\mathrm{<figure-callout\; id="3"\; label="sin"\; filenames="HDA00002494477700041.png,HDA00002494477700051.png"\; state="\{\{state\}\}">sin</figure-callout>}3\mathrm{D\&omega;}-\sin \mathrm{D\&omega;})+(\mathrm{<figure-callout\; id="4"\; label="sin"\; filenames="HDA00002494477700041.png,HDA00002494477700051.png"\; state="\{\{state\}\}">sin</figure-callout>}4\mathrm{D\&omega;}-\mathrm{<figure-callout\; id="2"\; label="sin"\; filenames="HDA00002494477700041.png,HDA00002494477700051.png"\; state="\{\{state\}\}">sin</figure-callout>}2\mathrm{D\&omega;})+\&CenterDot;\&CenterDot;\&CenterDot;+(\sin \mathrm{dD\&omega;}-\sin (d-2)\mathrm{D\&omega;})}{N\sin \mathrm{D\&omega;}}$

$=\frac{N\sin \mathrm{D\&omega;}-3\sin \mathrm{D\&omega;}-(\sin \mathrm{dD\&omega;}+\sin (d-1)\mathrm{D\&omega;}-\sin \mathrm{D\&omega;})}{N\sin \mathrm{D\&omega;}}$

$=\frac{(N-2)\sin \mathrm{D\&omega;}-\sin \mathrm{dD\&omega;}-\sin (d-1)\mathrm{D\&omega;}}{N\sin \mathrm{D\&omega;}}$

$=\frac{N-2}{N}-\frac{2\sin \left[(d-0.5)\mathrm{D\&omega;}\right]\cos (\mathrm{D\&omega;}/2)}{N\sin \mathrm{D\&omega;}}=\frac{N-2}{N}-\frac{2\sin \left[(d-0.5)\mathrm{D\&omega;}\right]\cos (\mathrm{D\&omega;}/2)}{2N\sin (\mathrm{D\&omega;}/2)\cos (\mathrm{D\&omega;}/2)}$

$=\frac{N-2}{N}-\frac{\sin \left[(d-0.5)\mathrm{D\&omega;}\right]}{N\sin (\mathrm{D\&omega;}/2)}=\frac{N-2}{N}-\frac{\sin \left[(\frac{N}{2}-1)\mathrm{D\&omega;}\right]}{N\sin (\mathrm{D\&omega;}/2)}---\left(2\right)$

For probing into the filtering characteristic at power frequency place, might as well study specific digital angular frequency, i.e. ω _kFollowing formula value during=k2 π/D, its corresponding analog frequency is f _k=ω _k/ (2 π) f _s=kf _s/ D, then the 2nd limit item that becomes 0/0 type of following formula according to the Luo Bita rule, can obtain ultimate value to this differentiate, namely simultaneously $\underset{\mathrm{\&omega;}\&RightArrow;\mathrm{<figure-callout\; id="2"\; label="k"\; filenames="HDA00002494477700041.png,HDA00002494477700051.png"\; state="\{\{state\}\}">k</figure-callout>}2\mathrm{\&pi;}/D}{\lim }H\left(e^{\mathrm{j\&omega;}}\right)=\underset{\mathrm{\&omega;}\&RightArrow;\mathrm{<figure-callout\; id="2"\; label="k"\; filenames="HDA00002494477700041.png,HDA00002494477700051.png"\; state="\{\{state\}\}">k</figure-callout>}2\mathrm{\&pi;}/D}{\lim }\frac{(N-2)}{N}-\frac{\cos \left[(N/2-1)\mathrm{D\&omega;}\right]\&CenterDot;(N/2-1)D}{N\cos (\mathrm{D\&omega;}/2)\&CenterDot;D/2}$

（3）

$=\frac{(N-2)}{N}-\frac{(N/2-1)D}{N\&CenterDot;D/2}=\frac{(N-2)}{N}-\frac{(N-2)}{N}=0$

Two, filter parameter method to set up

For so that wave filter is removed the interference of the baseline drift composition that faint electricity physiological signal is subject to below the 0.67Hz, and the 50Hz power frequency disturbed remove totally, need reasonably arrange filter parameter.After sampling rate was determined, the parameter setting here comprised determining interpolation factor D, wave filter continuation periodicity N.

1, the setting of interpolation factor D and setting principle thereof

The method of interpolation factor D parameter setting is: interpolation factor D is set to the ratio of sampling rate and power frequency.

Its principle is: interpolation factor D arranges according to being the interference of requirement removal power frequency.Suppose that the power frequency interfering frequency is f _p(China f for example _p=50Hz), then in order to remove the interference at power frequency place fully, be not difficult to know by inference by formula (3): can make k=1, then f _k=f ₁=f _s/ D=f _p, D=f namely satisfies condition _s/ f _pThe time, wave filter can be removed power frequency fully and disturb.

For example suppose f _p=50Hz, f _s=300Hz, then D=f _s/ f _p=6.

2, the computer Aided Design realizes enumerative technique parameters N

Before this wave filter of design, be set as the interpolation factor D that removes the power frequency interference and introduce except calculating, also need calculate continuation periodicity N(in order to prevent the phase shift between input signal and the output signal, N must be odd number).At known sample frequency f _sIn the situation of interpolation factor D, (corresponding 3dB cut-off frequency, namely amplitude-frequency response is at 0.67Hz to make cut-off frequency

Corresponding frequency location), like this periodicity N there is strict requirement.

In order to calculate more accurately the value of N, the present invention proposes a kind of enumerative technique.According to the transmission characteristic of wave filter:

$H\left(e^{\mathrm{j\&omega;}}\right)=\frac{N-2}{N}-\frac{\sin \left[(\frac{N}{2}-1)\mathrm{D\&omega;}\right]}{N\sin (\mathrm{D\&omega;}/2)}---\left(4\right)$

Thereby the method to set up of periodicity N is

(1) makes ω=2 π f/f _s/ _F=0.67Be updated to $H\left(e^{\mathrm{j\&omega;}}\right)=\frac{N-2}{N}-\frac{\sin \left[(\frac{N}{2}-1)\mathrm{D\&omega;}\right]}{N\sin (\mathrm{D\&omega;}/2)}$ In;

(2) change the N value of enumerating, the computer Aided Design is drawn out H (e ^{J ω}) with the change curve of cycle N;

(3) from change curve, find out H (e ^{J ω}) the most approaching

Corresponding odd number N value is the cycle parameter that will arrange.

For example suppose f _p=50Hz, f _s=300Hz, then D=6 according to above-mentioned three steps, can draw out H (e ^{J ω}) with the relation curve of enumerating Integer N as shown in Figure 4, be 0.707 for guaranteeing frequency response values, because the N value is set to 59.

Three, wave filter transmission curve

Suppose f _p=50Hz, sample frequency f _s=300Hz, in the practical application, these two frequencies are normally fixed.Making N=59 again, can obtain the output characteristic curve of this wave filter as shown in Figure 5, can visually see from Fig. 5, (a) is the coefficient h (n) of wave filter, and the length of wave filter is (N-3) D+1, i.e. (59-3) * 6+1=337.And except the coefficients comparison of central point was large, other nonzero coefficient size all equated, this and previous calculations coefficient out are consistent.(b) figure is the transfer curve of wave filter, accurately is 0 at the frequency response H at 50Hz place (j2 π f), and namely the trap effect is very good.Although this wave filter has corresponding trap at the integer multiple frequency place of backward 50Hz, but because the trap scope is very little, be limited within the frequency of 2Hz, and the above faint electricity physiological signal composition of 100Hz is less, therefore little on the impact of faint electricity physiological signal itself.

In order better to analyze and explanation, now Fig. 5 to be carried out the part and amplify as shown in Figure 6.From (a), not only can find out the distribution situation of nonzero coefficient size, and 5 (being D-1) null values are arranged between per two nonzero coefficients.(b) to be illustrated in the frequency response at 0.67Hz place be 0.695 to two red dotted lines among the figure, and error is 0.012, and wave filter of this explanation the present invention design satisfies the requirement that cut-off frequency is 0.67Hz.

Five, experiment is for example:

Example 1: with sampling rate f _s=900Hz, the human ecg signal that primary signal gathers for experiment.Because the sample rate of this design of filter is 300Hz, thus first primary signal is carried out 3 times down-sampling, and then carry out filtering with the wave filter of design.At this moment, each parameter of wave filter is respectively: f _s=300Hz, D=6, N=59, the waveform before and after the filtering is as shown in Figure 7.Can find out from the original electrocardiographicdigital signal waveform of 900Hz, the baseline drift of signal is very large, and the power frequency of 50Hz is disturbed more serious, through behind the wave filter, basically eliminate baseline drift and power frequency and disturbed, obtained comparatively clean electrocardiosignal, so that the further processing of signal.This illustrates that this wave filter can remove the baseline drift of electrocardiosignal really, and the power frequency of 50Hz is disturbed to have well to fall into remove.

Example 2: with sampling rate f _s=900Hz, the human ecg signal that primary signal gathers for experiment.Adding a frequency first in primary signal is 0.2Hz, and amplitude is 10 ⁴Noise signal, the electrocardiosignal of making an uproar that adds that will obtain is again carried out filtering, the time domain waveform before and after the filtering as shown in Figure 8.Can clearly find out, baseline drift can not only well be removed by this wave filter and power frequency is disturbed, and to the noise signal below the 0.67Hz also well filtering, thereby obtains pure electrocardiosignal.

Practical application example of the present invention explanation: with collection and the filtering of the faint electricity physiological signal among the present invention, and the processing procedure on Matlab7.9 and Visio studio platform is example, main divided data collection, processing, three parts of filtering.

1, data acquisition

File is preserved in being captured in of the faint electricity physiological signal of human body on the processor, conducting wire is connected the human body position of respectively leading.Position shown in Figure 9 connects each conducting wire.

2, date processing

Read the realization with deal with data on the Matlab7.9 platform.

3, filtering

Each data of leading of calculating are above carried out filtering according to filtering flow process shown in Figure 3 can obtain the faint electricity physiological signal of pure human body, faint electricity physiological signal used in the present invention all obtains according to this embodiment.

Claims (6)

1. the filter design method of a faint electricity physiological signal is characterized in that, the method may further comprise the steps:

Step 1 according to interpolation factor D and continuation periodicity N, makes d=(N-1)/2, and d represents the group delay of wave filter;

Step 2 is initialized as all 1's matrix that length is 2d-1 with filter coefficient h (n), and wave filter h is carried out the up-sampling that D is ordered, and namely inserts D-1 zero point in each coefficient h (n) back of wave filter h; And

Step 3, according to above-mentioned parameter, the expression formula of filter coefficient h (n) is

N, m represent integer,

Calculate all nonzero coefficients of filter coefficient h (n) according to above-mentioned formula; Step 4, center sampling point x _nThe value that multiply by (1-3/N) deducts the meansigma methods that step 2 is calculated, and is the filtering output y (n) of current time, and N is the continuation periodicity, thereby designs this wave filter h.

2. the filter design method of faint electricity physiological signal as claimed in claim 1 is characterized in that, described step 2 also comprises the step that the point with the last interpolation of ripple device coefficient h (n) deletes.

3. the filter design method of faint electricity physiological signal as claimed in claim 1 is characterized in that, the filtering sampling point of the method choose employing: with x _nAnd the x that is adjacent _N-1And x _N+1Centered by sampling point, in the sampling period take D as data decimation, select symmetrically respectively first group of sampling point data x _n, x _N-D, x _N+DSecond group of sampling point data x _N-1, x _N-D-1, x _N+D-1The 3rd group of sampling point data x _N+1, x _N-D+1, x _N+D+1Rule choose.

4. the filter design method of faint electricity physiological signal as claimed in claim 1 is characterized in that, described interpolation factor D is set to sample frequency f _sWith power frequency interfering frequency f _pRatio D=f _s/ f _p

5. the filter design method of faint electricity physiological signal as claimed in claim 1 is characterized in that, described continuation periodicity N adopts enumerative technique to try to achieve, and specifically may further comprise the steps:

Make ω=2 π f/f _s/ _F=0.67Be updated to $H\left(e^{\mathrm{j\&omega;}}\right)=\frac{N-2}{N}-\frac{\sin \left[(\frac{N}{2}-1)\mathrm{D\&omega;}\right]}{N\sin (\mathrm{D\&omega;}/2)}$ In, w represents frequency;

Change the N value of enumerating, the computer Aided Design is drawn out H (e ^{J ω}) with the change curve of cycle N;

From change curve, find out H (e ^{J ω}) the most approaching

Corresponding odd number N value is the continuation periodicity N that will arrange.

6. the filtering method of a faint electricity physiological signal is characterized in that, the method may further comprise the steps:

Step 1, n sampling point collection of input, namely $\{x_{n-\frac{N-1}{2}D},\&CenterDot;\&CenterDot;\&CenterDot;,x_{n-D},x_n,x_{n+D},\&CenterDot;\&CenterDot;\&CenterDot;,x_{n+\frac{N-1}{2}D}\},$ Two adjacent sampling point interval D sampling periods of this sampling point collection;

Step 2, n the sampling point that sampling point is concentrated add up and are averaging;

Step 3, with center sampling point x _nThe value that multiply by (1-3/N) deducts the meansigma methods that step 2 is calculated, and is the filtering output y (n) of current time, and N is the continuation periodicity.

Images