7-EMD端点效应优化（斜率匹配波形延拓法）（MATLAB视频）

Halcom

7-EMD端点效应优化（斜率匹配波形延拓法）（MATLAB视频）
百度网盘链接：https://pan.baidu.com/s/1eiMNqxuXBea7j6yw1-KZ4g
具体链接在halcom.cn论坛，联系人QQ：3283892722
该论坛是一个学习交流平台，我会逐一的和大家分享学习。
欢迎大家录制视频，你可在论坛进行打赏分享。
视频专用播放器：http://halcom.cn/forum.php?mod=viewthread&tid=258&extra=page%3D1\

主函数如下：

1. % 斜率匹配波形延拓法
   
2. clc,clear,close all % 清理命令区、清理工作区、关闭显示图形
   
3. warning off         % 消除警告
   
4. format short        % 数据类型
   
5. tic                 % 运算计时
   
6. %% 初始化信号
   
7. [t, signal] = gen_signal(5e3, 1e3);
   
8. figure(1), plot(t,signal,'linewidth',1)
   
9. xlabel('t');  ylabel('实际信号');  grid on;
   
10. %% 斜率匹配波形延拓法
    
11. %% 极值点和过零点位置提取
    
12. [indmin, indmax, indzer] = extr(signal, t); % 极值点和过零点位置提取
    
13. figure(2), plot(t,signal,'linewidth',1);
    
14. hold on
    
15. plot(t(indmin), signal(indmin), 'm.-')
    
16. plot(t(indmax), signal(indmax), 'r.-')
    
17. plot(t(indzer), signal(indzer), 'g.-')
    
18. xlabel('t');  ylabel('实际信号');  grid on;
    
19. legend('原始信号','极小值','极大值','过零点')
    
20. axis tight;
    
21. hold off;
    
22. %% 斜率匹配波形延拓法
    
23. if(indmax(1)<indmin(1))
    
24.     % 左端点沿拓
    
25.     index = find_first_indmin( indmax(1), indmin, 0 );
    
26.     index0 = find_first_indmin( indmax(1), indzer, 0 );
    
27.     A01 = (signal(indmax(1))-signal(indmin(index)))./( indmax(1)-indmin(index)+eps ); % M2-E0的斜率
    
28.     k=1;
    
29.     for i=2 : length(indmax)
    
30.         index = find_first_indmin( indmax(i), indmin, 0 );
    
31.         indexi = find_first_indmin( indmax(i), indzer, 0 );
    
32.         if(isempty(index))
    
33.             Scope(k,1)=inf;
    
34.             Amp(k,1) = inf;
    
35.             dt(k,1)=inf;
    
36.         else
    
37.             Scope(k,1)=(signal(indmax(i))-signal(indmin(index)))./( indmax(i)-indmin(index)+eps ); % M2-E0的斜率
    
38.             Amp(k,1) = abs( signal(indmax(1))-signal(indmax(i)) );
    
39.             dt(k,1) = abs( abs( indzer(index0)-indmax(1) ) - abs( indzer(indexi)-indmax(i) ) );
    
40.         end
    
41.     end
    
42.     k=k+1;
    
43.     [minV, index] = min(abs(Scope(:,1)-A01) + Amp + dt);
    
44.     newdata = signal( indmax(index) : indmax(index+1) );
    
45.     welddata_left = weld2( signal, newdata, indmax, 0 );
    
46.     % 右端点沿拓
    
47.     index = find_first_indmin( indmax(end), indmin, 1 );
    
48.     index0 = find_first_indmin( indmax(end), indzer, 1 );
    
49.     A01 = (signal(indmax(end))-signal(indmin(index)))./( indmax(end)-indmin(index)+eps ); % M2-E0的斜率
    
50.     k=1;
    
51.     for i=2 : length(indmax)
    
52.         index = find_first_indmin( indmax(i), indmin, 1 );
    
53.         indexi = find_first_indmin( indmax(i), indzer, 1 );
    
54.         if(isempty(index))
    
55.             Scope(k,1)=inf;
    
56.             Amp(k,1) = inf;
    
57.             dt(k,1)=inf;
    
58.         else
    
59.             Scope(k,1)=(signal(indmax(i))-signal(indmin(index)))./( indmax(i)-indmin(index)+eps ); % M2-E0的斜率
    
60.             Amp(k,1) = abs( signal(indmax(end))-signal(indmax(i)) );
    
61.             dt(k,1) = abs( abs( indzer(index0)-indmax(end) ) - abs( indzer(indexi)-indmax(i) ) );
    
62.         end
    
63.     end
    
64.     [minV, index] = min(abs(Scope(:,1)-A01)+Amp+dt);
    
65.     newdata = signal( indmax(index) : indmax(index+1) );
    
66.     welddata_right = weld2( signal, newdata, indmax, 1 );
    
67. 
    
68.     welddata = [welddata_left, signal(indmax(1):indmax(end)), welddata_right];
    
69. 
    
70.     figure(3),hold on
    
71.     plot(1:length(welddata_left), welddata_left,'r.-');
    
72.     plot(length(welddata_left)+1: length(welddata_left)+length(signal(indmax(1):indmax(end))),signal(indmax(1):indmax(end)),'b-');
    
73.     plot(length(welddata)-length(welddata_right)+1:length(welddata),welddata_right,'r.-');
    
74.     xlabel('t');  ylabel('实际信号');  grid on;
    
75.     axis tight;
    
76.     hold off;
    
77.    
    
78. else
    
79.     signal = -1*signal;
    
80.     [indmin, indmax, indzer] = extr(signal, t); % 极值点和过零点位置提取
    
81. %     figure(22), plot(t,signal,'linewidth',1);
    
82. %     hold on
    
83. %     plot(t(indmin), signal(indmin), 'm.-')
    
84. %     plot(t(indmax), signal(indmax), 'r.-')
    
85. %     plot(t(indzer), signal(indzer), 'g.-')
    
86. %     xlabel('t');  ylabel('实际信号');  grid on;
    
87. %     legend('原始信号','极小值','极大值','过零点')
    
88. %     axis tight;
    
89. %     hold off;
    
90.    
    
91.     % 左端点沿拓
    
92.     index = find_first_indmin( indmax(1), indmin, 0 );
    
93.     index0 = find_first_indmin( indmax(1), indzer, 0 );
    
94.     A01 = (signal(indmax(1))-signal(indmin(index)))./( indmax(1)-indmin(index)+eps ); % M2-E0的斜率
    
95.     k=1;
    
96.     for i=2 : length(indmax)
    
97.         index = find_first_indmin( indmax(i), indmin, 0 );
    
98.         indexi = find_first_indmin( indmax(i), indzer, 0 );
    
99.         if(isempty(index))
    
100.             Scope(k,1)=inf;
     
101.             Amp(k,1) = inf;
     
102.             dt(k,1)=inf;
     
103.         else
     
104.             Scope(k,1)=(signal(indmax(i))-signal(indmin(index)))./( indmax(i)-indmin(index)+eps ); % M2-E0的斜率
     
105.             Amp(k,1) = abs( signal(indmax(1))-signal(indmax(i)) );
     
106.             dt(k,1) = abs( abs( indzer(index0)-indmax(1) ) - abs( indzer(indexi)-indmax(i) ) );
     
107.         end
     
108.     end
     
109.     k=k+1;
     
110.     [minV, index] = min(abs(Scope(:,1)-A01) + Amp + dt);
     
111.     newdata = signal( indmax(index) : indmax(index+1) );
     
112.     welddata_left = weld2( signal, newdata, indmax, 0 );
     
113.     % 右端点沿拓
     
114.     index = find_first_indmin( indmax(end), indmin, 1 );
     
115.     index0 = find_first_indmin( indmax(end), indzer, 1 );
     
116.     A01 = (signal(indmax(end))-signal(indmin(index)))./( indmax(end)-indmin(index)+eps ); % M2-E0的斜率
     
117.     k=1;
     
118.     for i=2 : length(indmax)
     
119.         index = find_first_indmin( indmax(i), indmin, 1 );
     
120.         indexi = find_first_indmin( indmax(i), indzer, 1 );
     
121.         if(isempty(index))
     
122.             Scope(k,1)=inf;
     
123.             Amp(k,1) = inf;
     
124.             dt(k,1)=inf;
     
125.         else
     
126.             Scope(k,1)=(signal(indmax(i))-signal(indmin(index)))./( indmax(i)-indmin(index)+eps ); % M2-E0的斜率
     
127.             Amp(k,1) = abs( signal(indmax(end))-signal(indmax(i)) );
     
128.             dt(k,1) = abs( abs( indzer(index0)-indmax(end) ) - abs( indzer(indexi)-indmax(i) ) );
     
129.         end
     
130.     end
     
131.     [minV, index] = min(abs(Scope(:,1)-A01)+Amp+dt);
     
132.     newdata = signal( indmax(index) : indmax(index+1) );
     
133.     welddata_right = weld2( signal, newdata, indmax, 1 );
     
134. 
     
135.     welddata = [welddata_left, signal(indmax(1):indmax(end)), welddata_right];
     
136. 
     
137. %     welddata = -1*welddata;
     
138.    
     
139.     figure(3),hold on
     
140.     plot(1:length(welddata_left), -welddata_left,'r.-');
     
141.     plot(length(welddata_left)+1: length(welddata_left)+length(signal(indmax(1):indmax(end))),-signal(indmax(1):indmax(end)),'b-');
     
142.     plot(length(welddata)-length(welddata_right)+1:length(welddata),-welddata_right,'r.-');
     
143.     xlabel('t');  ylabel('实际信号');  grid on;
     
144.     axis tight;
     
145.     hold off;
     
146.    
     
147. end
     

复制代码