《[2017年整理]脉搏信号处理y105》由会员分享，可在线阅读，更多相关《[2017年整理]脉搏信号处理y105（18页珍藏版）》请在金锄头文库上搜索。
1、 数字信号处理综合训练说 明 书题目：脉搏信号间期序列的时域、频域分析学院：电气工程与信息工程学院班级：电子信息科学与技术（1）班姓名：马艳霞学号：082601372011年7月15日8目录1 脉搏信号处理的基本思路.12 信号预处理.23 脉搏信号的时域分析.54 脉搏信号频域分析.85 程序清单.96 心得体会.237 参考文献.2495 程序清单%加载采样数据；y=y10 y11 y12 y13 y14 y15;t=1:6000; %时域范围；figure(1);subplot(2,1,1);y=-y; plot(t,y); grid on;title(脉搏信号); xlabel(时间/
2、ms); ylabel(相对幅度);%截取整数倍周期data1=min(y(1:500);data2=min(y(5000:5500);for i=1:500if y(i)=data1m1=i;endendfor j=5000:5500if y(j)=data2m2=j;endend10mb=m2-m1; y=y(m1:m2); t=1:mb+1;subplot(2,1,2); plot(t,y); title(脉搏信号); xlabel(时间/ms); ylabel(相对幅度); grid on;原始的脉搏信号和经过整数倍提取后的信号如下图所示图1%带通滤波器,虑除工频.基线漂移和肌电干扰；
3、figure(2);fs=1500;%采样频率；wp=0.9,30/(fs/2);%设置通带截止频率;ws=0.1,200/(fs/2);%设置阻带截止频率;rp=0.5;%通带波纹系数;11rs=40;%阻带波纹系数;N,wc=buttord(wp,ws,rp,rs);num,den=butter(N,wc);%滤波器分子分母系数向量；H,W=freqz(num,den);subplot(2,1,1);plot(fs*W/(2*pi),abs(H);title(butterworth带通滤波器幅频响应);xlabel(频率/Hz); ylabel(幅度); grid on;y1=filtfi
4、lt(num,den,y);%脉搏信号通过零相位滤波器；subplot(2,1,2); plot(y1);title(脉搏信号通过带通滤波器);xlabel(时间/ms); ylabel(幅度); grid on;脉搏信号通过零相位带通滤波器后，高频信号和基线漂移信号有所衰减，信号变得清晰。图212%小波包消噪；%进行一维小波消噪figure(3);z=wden(y1,heursure,s,mln,3,sym8);%XD,CXD,LXD=wden(X,tptr,sorh,scal,n,wavename)使用小波系数阈值，返回输入信号X除噪后的信号XD，输出参数CXD,LXD表示XD的小波分解结
5、构。 %输入参数中，tptr同thselect（）函数；sorh为 s或 h表示软硬阈值；n表示在n层上的小波分解；wavename指定小波名称；scal定义阈值调整比例：%one不设定比例； %sln使用的基于第单层系数噪声估计，设置比例； %mln用噪声层的层相关估计，调整比例。 subplot(2,1,1);plot(z); title(小波消噪后的信号);xlabel(时间/ms); ylabel(幅度); grid on;%小波包消噪，固定阈值nn=length(y1);13thr=sqrt(2*log(nn*log(nn)/log(2);y2=wpdencmp(y1,s,4,db4
6、,sure,thr,1);subplot(2,1,2); plot(y2);title(小波包消噪后的信号y2);xlabel(时间/ms); ylabel(幅度); grid on;图3小波包消除信号中夹杂的突变信号，使得脉搏信号能准确反映人体的生理状况。经过消噪后的脉搏信号如图3所示，比未经处理的信号更清晰，准确。figure(4);s1=fft(y,fs);%原始脉搏信号傅里叶变换；subplot(2,1,1);plot(abs(s1);%幅频特性； axis(0,400,0,40);%定坐标范围；title(原脉搏信号幅频响应);14xlabel(频率/Hz); ylabel(幅度);
7、 grid on;s2=fft(y2,fs);subplot(2,1,2); plot(abs(s2);axis(0,400,0,15);%定坐标范围；title(消噪后脉搏信号幅频特性);xlabel(频率/Hz); ylabel(幅度); grid on;图4原始的脉搏信号经过带通滤波器以及小波消噪后，主要频率集中在040HZ之间，极低频段的基线漂移信号被滤除，高频段的工频干扰和肌电干扰也被滤除。%时域分析；%求1s中出现的最大值，取其0.6倍作为阈值，计算y的平均幅值pRv13=y1; v13=v13-mean(v13); %去基频直流分量;15v13=v13 zeros(1,lengt
8、h(y)-length(v13);m=ceil(length(v13)/1000)-1;for j=1:ms(j+1),t(j)=max(v13(1+1000*(j-1):1000*j);s5(j)=max(s(j+1);%s4,t=max(x3(2000:2890);endpR=sum(s5)/m;%平均幅值； thr1=0.6*pR;%阈值；%求0.5s中出现的最大值，幅值超过阈值的认为是R波，计算脉率率Prate；rr=0; z=round(length(v13)/500)-1;for q=1:zs0(q+1),t(q)=max(v13(1+500*(q-1):500*q);s8(q)=
9、max(s0(q+1);if s8(q)=thr1rr=rr+1;%峰值点的个数；endendPrate=20*rr %脉率；%标出y主波的位置；for l=1:length(v13)%范围；if v13(l)thr1 %大于阈值时，v13=v14;16v14(l)=v13(l);elsev14(l)=0;%小于阈值时v14=0;endendz1=round(length(v13)/1000)-1;for j=1:z1s(j+1),t(j)=max(v14(1+1000*(j-1):1000*j);s9(j)=max(s(j+1);%s4,t=max(x3(2000:2890);for l=1
10、+1000*(j-1):1000*jif v14(l)=s9(j)v14(l)=v14(l);elsev14(l)=0;endendendR1=0;for l=1:length(v14)if v14(l)0R1=R1+1; %峰值个数；17endendR1=R1-1;figure(5);subplot(2,1,1);plot(v13); %去基频直流分量后的脉搏波;xlabel(时间/ms); ylabel(幅值); grid on;subplot(2,1,2); plot(v14); title(峰值检测点); xlabel(时间/ms); ylabel(幅值); grid on; 图5虽然
11、每个脉搏波的波峰值大小并不相同，但他们总在一个范围内波动，波动的范围基本上不超过最大波形高度的0.3倍，为了更为可靠的检测波峰，以0.6倍为参考。因此可以认为波峰点是在每个脉18搏周期中波形的最大值附近，大于其邻域内所有点的点。%提取一个周期波形figure(6);data1=min(y2(750:1500);data2=min(y2(1500:2500);for i=750:1500if y2(i)=data1N11=i;%单周期起点；endendfor j=1500:2500if y2(j)=data2N21=j;%单周期终点；endendNc=N21-N11;%单周期长度；y4=y2(N
12、11:N21);subplot(2,1,1);plot(y4);%绘制单周期波形； grid on;title(小波包消噪后的信号y2的单周期波形);xlabel(时间/ms); ylabel(幅值);s4=fft(y4,fs);%傅里叶变换；19subplot(2,1,2);plot(abs(s2); axis(0,400,0,10);%定坐标范围；title(单周期脉搏信号幅频特性);xlabel(频率/Hz); ylabel(幅度);grid on;图6上图为脉搏信号的单周期时域波形和频域波形图%求R点位置,计算PULSE间期MEAN;R=zeros(1,rr);u3=0;for r=1
13、:(length(v14)-2)v15(r)=v14(r+1)-v14(r);v15(r+1)=v14(r+2)-v14(r+1);if v15(r)0&v15(r+1)=2RR(u3-1)=R(u3)-R(u3-1); %间期序列；endendendMEAN=sum(RR)/(R1) %脉搏信号间期平均值；%计算PRR间期总体标准差差SDHR、心率变异指数CV;for i1=1:length(RR)HR(i1)=(RR(i1)-MEAN).2; endSDRR=sqrt(sum(HR)/(length(HR)-1) %脉搏信号间期标准差；for i1=1:length(RR)-1HR1(i1
14、)=HR(i1+1)-HR(i1);HR2(i1)=HR1(i1).2;HR3(i1)=mean(HR1);endRMSSD=sqrt(sum(HR2)/(length(HR)-1)%相邻PR间期差的均方根；SDSD=sqrt(sum(HR1-HR3).2)/(length(HR)-1)%全部PR间期21的标准差；CV=SDRR/MEAN %心率变异系数；根据公式计算的脉搏信号间期序列的平均值，标准差，均方根，差值标准差以及心率的变异系数。%频域分析（功率谱估计）;figure(7);y1=y1-mean(y1); %去直流,去趋势；y1=detrend(y1); m1=mb;nfft=256
15、;%设定FFT算法的长度;n1=0:m1-1; t=n1/fs; %数据长度、时间序列window=hanning(256); %选用的窗口noverlap=128; %分段序列重叠的采样点数（长度）dflag=none; %不做趋势处理subplot(2,1,1);pxx=psd(y,nfft,fs,window,noverlap,dflag); %功率谱估计f=(0:nfft/2)*fs/nfft; %求得对应的频率向量plot(f,10*log10(pxx); %绘制功率谱xlabel(频率/Hz); ylabel(功率谱/dB);title(滤波前的功率谱图：PSDWelch方法); grid on;subplot(2,1,2);pxx=psd(y1,nfft,fs,window,nove
