北工大通信系统实验

  1. 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
  2. 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
  3. 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。

信号处理与分析课程设计训练报告

---------- ECG信号高频噪声消除

学号:

姓名:

指导老师:张新峰

完成日期:2016年6月5日

绪论

本实验是信号处理相关方法在心电图ECG(electocardiogram)信号处理中的一个具体应用。通过心电图进行疾病诊断是心血管疾病诊断的一种重要方法。该方法可靠简便,在临床中应用很广泛。近年来给予数字信号处理和模式识别的ECG信号自动分析成为一个研究这点。ECG信号在采集过程中收到工频信号、肌电信号以及基线漂移等干扰,由此带来的噪声对自动分析的结果影响很大。因此在对ECG信号分析之前,必须进行噪声消除研究。本实验内容给定ECG信号,合理设计滤波器,消除高频噪声的影响。

一、实验内容

心电信号是微弱低频人体生理电信号,通常频率在0.05~100Hz,幅值不超过4mv,通过安装在皮肤表面的电极来获取。由于人体是一个复杂的生命系统,存在50H 工频干扰及基线漂移等其他生理电信号的干扰。噪声可能会影响到医生的临床诊断,因此,需对心电信号进行滤波,即必须做好前端数据采集的软硬件设计以保证心电数据的可靠和准确。

(本次实验的心电数据是老师提供,为txt文档)

传统医疗设备分别采用50Hz 带阻滤波器和RC 高通滤波器滤除工频干扰和基线漂移。但带阻滤波器电路复杂,其特性对元器件的精度敏感,而基线漂移本质上是一种缓慢变化的低频信号,采用RC 滤波器很难将高通滤波器的过渡带做得十分陡峭,基线漂移补偿效果不理想。因此,模拟方法往往不太容易获得很好的特性。

在本实验中,我们以C语言为基础,来实现高频噪声的消除。

二、实验原理

1. 低通巴特沃兹模拟滤波器的实现

1)概述

巴特沃兹滤波器的特点是通频带内的频率响应曲线最大限度平坦,没有起伏,而在阻频带则逐渐下降为零。在振幅的对数对角频率的波得图上,从某一边界角频率开始,振幅随着角频率的增加而逐步减少,趋向负无穷大。

一阶巴特沃兹滤波器的衰减率为每倍频6分贝,每十倍频20分贝。二阶巴特沃兹滤波器的衰减率为每倍频12分贝、三阶巴特沃兹滤波器的衰减率为每倍频18分贝、如此类推。巴特沃兹滤波器的振幅对角频率单调下降,并且也是唯一的无论阶数,振幅对角频率曲线都保持同样的形状的滤波器。只不过滤波器阶数越高,在阻频带振幅衰减速度越快。其他滤波器高阶的振幅对角频率图和低级数的振幅对角频率有不同的形状。

2)设计巴特沃兹模拟低通滤波器

1确定参数:

p w ,s w ,p Ω,s Ω

2求滤波器阶次N :

sp λ=s p

ΩΩ sp k

N=

lg lg sp sp k λ

4确定系统函数

N 为奇数时:

∏-=+=2

/)1(1

)(11)(N k k p G p p G N 为偶数时:

∏==

2

/1

)()(N k k p G p G 双线性Z 法把模拟滤波器转换成数字滤波器

3. FFT 的算法

按时间抽选(DIT )的基-2 FFT 算法。当输入序列点数为(L 为整数),先将序列x(n)按n 的奇偶分成以下两组:

则可以将DFT 化为:

k=0,1,2…,

其中

分别为点DFT 。又根据系数的周期性可得到X(k)的

完整表达:

三、实现方法流程图

1. 读取心电信号数据的实现

float agTemp2[1024+2];//定义一个数组用来存放数据

FILE *p;//定义指针

p=fopen("119_heart1.txt","r"); // 打开心电图数据文件

for(i=0;i<1024+2;i++) // 通过循环放入数组中

{

fscanf(p,"%f",&agTemp2[i]);

}

fclose(p);

2. FFT变换的实现

a. 倒位序处理,先将x[n]进行倒位序后,才能进行fft变换

void change(){

complex temp;

unsigned short i=0,j=0,k=0;

double t;

for(i=0;i

k=i;j=0;

t=(log(Num)/log(2));

while( (t--)>0 ){

j=j<<1;

j|=(k & 1);

k=k>>1;

}

if(j>i){

temp=x[i];

x[i]=x[j];

x[j]=temp;

}

}

}

b. 进行蝶形运算,将保存在全局复数数组x[n]中的时域信号转为频域信号:

void fft(){

int i=0,j=0,k=0,l=0;

complex up,down,product;

change();

for(i=0;i< log(Num)/log(2) ;i++){ /*一级蝶形运算*/

l=1<

for(j=0;j

for(k=0;k

mul(x[j+k+l],W[Num*k/2/l],&product);

add(x[j+k],product,&up);

sub(x[j+k],product,&down);

x[j+k]=up;

x[j+k+l]=down;

}

}

}

}

c. 初始化, 复数类全局数组W[k]

void initW(){

int i;

W=(complex *)malloc(sizeof(complex) * Num);

for(i=0;i

{

W[i].real=cos(2*PI/Num*i);

W[i].img=-1*sin(2*PI/Num*i);

}

}

3. BTWZ滤波器和双线性Z转化为数字滤波器的实现

设计巴特沃兹低通模拟滤波器时,我严格按照上学期的方法进行思考,思路是先做出2阶的滤波器,如果阶次高的话,我就把二阶滤波器做级联完成高阶滤波器的滤波。以下是二阶滤波器实现过程。

完成BTWZ滤波器的系数计算

void BTWZ(float fp, float fs,float ap, float as)

{

int Ni,i;

double N,P[1000],Z[3],z[3],Wp;

Wp=2*PI*fp;

N=log(sqrt((pow(10,(as/10))-1)/(pow(10,(ap/10))-1)))/log(2);

if(N==(double)Ni)

Ni=(int)N;

else

Ni=(int)N+1;

if(Ni%2==0)

for(i=0;i<(Ni/2);i++)

P[i]=2*cos((2*i+Ni-1)*PI/(2*N));

else

{

P[0]=0;

for(i=0;i<((Ni-1)/2);i++)

相关文档
最新文档