生物医学信号处理第一次上机.docx
《生物医学信号处理第一次上机.docx》由会员分享，可在线阅读，更多相关《生物医学信号处理第一次上机.docx（13页珍藏版）》请在冰豆网上搜索。
生物医学信号处理第一次上机
第一次上机
1,Signalobtainandinput:
EEGineyescloseandopen（eegclose.matandeegopen.mat）.Fs=250Hz
Question1,howlong,andhowmanychannelswereusedwhenwerecordedthesetwodata,respectively?
load（'E:
\1111\lab1\eegclose.mat'）;
load（'E:
\1111\lab1\eegopen.mat'）;
my_close=U（:
11）;
my_open=U（:
22）;
fs=250;
close_length=length（my_close）/fs
open_length=length（my_open）/fs
Myclose=U（:
11）,绘图得到eeclose长度为4500，采样率为：
250hz,open_length=30
Myopen=U（:
22）绘图得到eeopen长度为4500，采样率为：
250hz,close_length=30
2,Pickuponechannelsignalfrombothdatatoanalyze.
max（my_open）=46.1769
max（my_close）=46.1769
min（my_open）=46.1769
min（my_close）=46.1769
Question2,statethemax.andmin.ofbothtimecourses.Inyourlabreport,includeapictureofbothsignals.UnitsaremsanduV.
3,Designtwofilterstoget4-8Hzand8-12Hzsignals.YoucanuseFIRorIIRfilterstofinishthisrequest.（ifyouhavetroubleinthisstep,seecue1）
3.1,fs=250;
Rp=1;Rs=50;
Wp=[2*4/fs2*8/fs];Ws=[2*2/fs2*10/fs];
[N,Wn]=ellipord（Wp,Ws,Rp,Rs）;
[B,A]=ellip（N,Rp,Rs,Wn）;
freqz（B,A）
fs=250;
Rp=1;Rs=50;
Wp=[2*8/fs2*12/fs];Ws=[2*6/fs2*14/fs];
[N,Wn]=ellipord（Wp,Ws,Rp,Rs）;
[D,C]=ellip（N,Rp,Rs,Wn）;
freqz（D,C）
Question3,statethefeaturesofyourfilters.Inyourlabreport,includetwopicturesofbothfilters’frequencyresponse.
•4,Filtertherawsignaltoobtaintwonewsignalsbybothmethods.（filter.mandfiltfilt.m）
loadmy_open
y=filter（B,A,my_open）;
plot（my_open）;
holdon
plot（y,'y'）
y=filter（B,A,my_close）;
plot（my_close）;
holdon
plot（y,'y'）
title（'filtermy_close’）
•
•y=filter（D,C,my_open）;
•plot（my_open）;
•holdon
•plot（y,'y'）
•title（'filtfilt.my_open'）
y=filter（D,C,my_close）;
plot（my_close）;
holdon
plot（y,'y'）
Question4,statethedifferencebetweentwofiltermethods.Inyourlabreport,includeapictureofbothnewsignals.
Nowyouhavetwonewsignalsforbothconditions,respectively.
•5,PowerSpectralDensity（PSD）estimateviaperiodogrammethodandWelch'smethodtoanalyzefournewsignals.（ifyouhavetroubleinthisstep,seecue3）
•fs=250;
•Rp=1;Rs=50;
•Wp=[2*4/fs2*8/fs];Ws=[2*2/fs2*10/fs];
•[N,Wn]=ellipord（Wp,Ws,Rp,Rs）;
•[B,A]=ellip（N,Rp,Rs,Wn）;
•freqz（B,A）
•y=filter（B,A,my_open）;
•plot（my_open）;
•holdon
•plot（y,'r'）
•Fs=1000;
•t=0:
1/Fs:
1;
•xn=y
•figure
（1）
•periodogram（xn,[],[],Fs）;%周期图法
•figure
（2）
•pwelch（xn,[],[],[],Fs）;%welch方法
•
代码2：
4~8hz滤波器实现
fs=250;
Rp=1;Rs=50;
Wp=[2*4/fs2*8/fs];Ws=[2*2/fs2*10/fs];
[N,Wn]=ellipord（Wp,Ws,Rp,Rs）;
[B,A]=ellip（N,Rp,Rs,Wn）;
freqz（B,A）
y=filter（B,A,my_close）;
plot（my_close）;
holdon
plot（y,'r'）
Fs=1000;
t=0:
1/Fs:
1;
xn=y
figure
（1）
periodogram（xn,[],[],Fs）;%周期图法
figure
（2）
pwelch（xn,[],[],[],Fs）;%welch方法
Question5,statethefeaturesofthePSDforeyesopenandclosedata,andpointthedifferencesbetweenthem.
Question6,DescribethedifferencesofbothPSDmethods.
Inyourlabreport,includePSDfigures.
•6,Repeatstep5fortherawsignaltoseethedifferenceinbroadbandsbetweenbothdata.
•
•
•
Question7,statethedifferencebetweentworawdata.
Inyourlabreport,includeapictureofPSDforbothdata.
my_open信号0.3-0.4频域内幅度减少20dB，my_close信号幅度减少35dB
Question8,Couldyouidentifyeyes-openandeyes-closedstatebyPSDobservation?
Andhowdoyoudothat?
figure
psd（my_open）
title（'my_open'）
figure
psd（my_close）
title（'my_close'）
升级会员 