随机信号的数字特征分析

四、实验结果与分析
% 估计两个相似信号间的时间延迟 load ecgdata;load eegdata;load icpdata;load respdata; np = 0:99; p = ones(size(np)); %方波 p2 = exp(-0.06*np); %指数衰减 p3 = sin(pi/5*np); %正弦 p4 = sin(pi/5*np).*exp(-0.06*np); %指数衰减正弦 p5 = (ecgdata (1:100))';%心电 p6 =(eegdata(1:100))'; %脑电 p7 =(icpdata(1:100))'; %颅内压 p8 =(respdata(1:100))'; %呼吸频率
-200
0
200
400
600
800
1000
-0.4 -0.6 -0.8 -1
5 0 -5
Rpx of p(n) and x(n) 200 0 -200 -1000
0
10
20
30
40
50
60
70
80
90
100
0
100
200
300
400
500
600
700
800
900
1000
-800
-600
-400
-200
第 2 页共 13 页
《生物医学信号处理》实验报告
%绘制信号图 figure; subplot(2,2,1); subplot(2,2,2); subplot(2,2,3); subplot(2,2,4); figure; subplot(2,2,1); subplot(2,2,2); subplot(2,2,3); subplot(2,2,4);
0
200
400
600
800
1000
图 3：模板为正弦信号，A=3，噪声均值为 0，方差为 1
指数衰减正弦 1 0.8 0.6 0.4
Rpw of p(n) and w(n) 10
5 Noise
0 -10 -1000 20
0
-800
-600
-400
-200
0
200
400
600
800
1000
-5
Rps of p(n) and s(n)
x x
n nm
m 0, 1, 2,
求和项总数不是 N 而是 N-|m|，因为当 n=N-|m|-1 时，n+|m|=N-1。此时 xn+m 已经到了数据 边沿。这种估计是渐进无偏估计和一致估计。 计算中，只要将其中一个序列反转，就可以用计算线性卷积的程序计算线性相关
rxy m xn y n
400
500
600
700
800
900
1000
-800
-600
-400
-200
0
200
400
600
800
1000
Signal with Noise
Rpx of p(n) and x(n)
0
10
20
30
40
50
60
70
80
90
100
0
100
200
300
400
500
600
700
800
900
1000
-800
因此可以用 FFT 来加速相关运算，即对序列补零后，用循环相关计算线形相关，然后用循 环卷积的快速算法计算循环相关，得到最终结果。 已知发射波形，利用相关技术，在有强背景噪声的情况下检测回波的延时和强度。 首先使用已知信号模版及其若干次衰减延迟生成仿真回波波形，然后与白噪声背景叠加， 构造仿真信号。然后计算模版与仿真信号的相关函数，判断回波位置及相对强度。
plot(np,p);title('方波'); plot(np,p2);title('指数衰减'); plot(np,p3);title('正弦'); plot(np,p4);title('指数衰减正弦'); plot(np,p5);title('心电信号'); plot(np,p6);title('脑电信号'); plot(np,p7);title('颅内压信号'); plot(np,p8);title('呼吸频率信号');
900
1000
Rpx of p(n) and x(n)
n = 0:1000; w = randn(size(n)); %产生白噪声 %%%改变噪声类型 %%A=fspecial('average'); %生成均值滤波器 %%w=filter2(A,w); %用生成的滤波器进行滤波 %%w=medfilt1(w); %中值滤波 %%%改变噪声强度 %%x =sqrt()*w s = zeros(size(n)); A = 3; %衰减系数 %构造仿真回波信号 s(100:199) = s(100:199)+A*(p); s(500:599) = s(500:599)+A/3*(p); s(800:899) = s(800:899)+A/3/3*(p); %构造仿真信号 x = s+w; figure; subplot(3,1,1); plot(n,w); title('Noise'); subplot(3,1,2); plot(n,s); title('Signal'); subplot(3,1,3); plot(n,x); title('Signal with Noise'); p = [p,zeros(1,length(x)-length(p))]; %如果要求归一化相关系数（相干系数） %计算相关系数 %两个序列要同样长 Rps = xcorr(s,p);%加'coeff'求相干系数 Rpw = xcorr(w,p); Rpx = xcorr(x,p); n2 = (n(1)-n(end)):(n(end)-n(1)); %绘制线性相关函数图 figure; subplot(3,1,1); plot(n2,Rpw); title('Rpw of p(n) and w(n)'); subplot(3,1,2); plot(n2,Rps); title('Rps of p(n) and s(n)'); 第 3 页共 13 页
-800
-600
-400
-200
0
200
400
600
800
1000
Rps of p(n) and s(n) 200 0
5
0
0
-0.2
-5
-200 -1000
0 100 200 300 400 500 600 700 800 900 1000 Signal with Noise
-800
-600
-400
-200
0
200
400
600
800
1000
0
100
200
300
400
500
600
700
800
900
1000
图 1：模板为方波信号，A=3，噪声均值为 0，方差为 1
指数衰减 1
Noise
Rpw of p(n) and w(n) 10 0 -10 -1000 40 20 0 -1000 50 0 -50 -1000
四、思考题
尝试修改程序，包括改变仿真信号中模版的形状，噪声的强弱，噪声的类型（对白噪声滤 波可以获得各种有色噪声） ，哪些因素会影响相关函数的结果？
五、实验小结：
一、改变模板的形状
第 4 页共 13 页
《生物医学信号处理》实验报告
方波 2
Rpw of p(n) and w(n) 0.1
Noise 5
1.8 1.6 1.4 1.2
0
0
-0.1 -1000 1
-800
-600
-400
-200
0
200
400600800来自1000-5Rps of p(n) and s(n)
4 2 0
0
100
200
300
400
500 Signal
600
700
800
900
1000
1
0
0.8 0.6 0.4 0.2
-1 -1000 1 0
5 0 -5
0
100
200
300
400
500 Signal
600
700
800
900
1000
0.2
0
0 -0.2 -0.4 -0.6 -0.8
-20 -1000 20 0 -20 -1000
-800
-600
-400
-200
0
200
400
600
800
1000
0
100
200
300
400
500
600
700
800
-800
-600
-400
-200
0
200
400
600
800
1000
0
100
200
300
400
500
600
700
800
900
1000
Rpx of p(n) and x(n)
10 0 -10
Signal with Noise