function [tmin,tmax,zmin,zmax] = boundary_matchwave(indmin,indmax,t,x)
% 一种自适应波形匹配端点延拓法
% 输入参数说明: indmin 极小值地址序列
% indmax 极大值地址序列
% t 时间序列(时序)1:length(x)
% x 时间序列(值)
%
% 输出参数说明: tmin 极小值地址序列
% tmax 极大值地址序列
% zmin 极小值
% zmax 极大值
%
lx = length(x);
indmax=find(diff(sign(diff(x)))==-2)+1;
indmin=find(diff(sign(diff(x)))==2)+1;
% 判断极值点个数
if (length(indmin) + length(indmax) < 3)
error('not enough extrema')
end
% 左端处理
leftmax = indmax(1);
leftmin = indmin(1);
if leftmin > leftmax
left = leftmin; % 左边处理长度
leftindex = indmin;
leftminind = indmax;
leftdisc = leftmax; % 末端长度
else
left = leftmax;
leftindex = indmax;
leftminind = indmin;
leftdisc = leftmin;
end;
leftvalue = x(1:left);
indexlen = length(leftindex);
mintol = 99999;
for i = 2:indexlen
len = leftindex(i);
leftvalue1 = x(len-left+1:len);
[tol,err,extrx1] = SelfAdapMatchWave(leftvalue,leftvalue1);
if tol < mintol
mintol = tol;
%avg=(sum(abs(extrx1))+sum(abs(rightvalue)))/(length(extrx1)*2-2);%有人说此处要改为leftvalue
avg=(sum(abs(extrx1))+sum(abs(leftvalue)))/(length(extrx1)*2-2);
erroo = leftvalue - extrx1;
avgerr=sum(abs(erroo))/(length(erroo)-1);
minextrx1 = extrx1;
newerr = err;
minindex = len;
indmini = i;
end
end;
%minextrx1
if (mintol < avg/(10*2) || mintol < avgerr*2) % 没有找到匹配子波,左端处理,取原始信号最左端的两个相邻极大值点的均值作为左端点的极大值，取信号最左端的两个相邻极小值点的均值作为左端点的极小值
lindmin = indmin(1)-indmin(2);
lvmin = (x(indmin(1))+x(indmin(2)))/2;
lindmax = indmax(1)-indmax(2);
lvmax = (x(indmax(1))+x(indmax(2)))/2;
else % 找到匹配子波,左端处理
leftindex = leftminind(indmini);
leftx = x(leftindex-left-1:leftindex)+newerr;
leftx = [leftx,x(1:2)];
[lhindmin,lhindmax,lhindzer] = extr(leftx);
lindmin = lhindmin(end)-(length(leftx)-2);
lindmax = lhindmax(end)-(length(leftx)-2);
lvmin = leftx(lhindmin(end));
lvmax = leftx(lhindmax(end));
end;
% 右端处理
rightmax = indmax(end);
rightmin = indmin(end);
if rightmin < rightmax
right = lx - rightmin; % 右边处理长度
rightindex = indmin;
rightmaxind = indmax;
rightdisc = lx - rightmax; % 末端长度
else
right = lx - rightmax;
rightindex = indmax;
rightmaxind = indmin;
rightdisc = lx - rightmin;
end;
rightvalue = x(end-right+1:end);
indexlen = length(rightindex)-1;
mintol = 99999;
for i = 1:indexlen
len = rightindex(i);
rightvalue1 = x(len:len+right-1);
[tol,err,extrx1] = SelfAdapMatchWave(rightvalue,rightvalue1);
if tol < mintol
mintol = tol;
avg=(sum(abs(extrx1))+sum(abs(rightvalue)))/(length(extrx1)*2-2);
erroo = rightvalue - extrx1;
avgerr=sum(abs(erroo))/(length(erroo)-1);
minextrx1 = extrx1;
newerr = err;
minindex = len;
indmini = i;
end
end;
%minextrx1
if (mintol < avg/(10*2) || mintol < avgerr*2) % 没有找到匹配子波,右端处理,取原始信号最右端的两个相邻极大值点的均值作为右端点的极大值，取信号最右端的两个相邻极小值点的均值作为右端点的极小值
rindmin = lx + (indmin(end)-indmin(end-1));
rvmin = (x(indmin(end-1))+x(indmin(end)))/2;
rindmax = lx + (indmax(end)-indmax(end-1));
rvmax = (x(indmax(end-1))+x(indmax(end)))/2;
else % 找到匹配子波,右端处理
rightindex = minindex; % rightmaxind(indmini)
rightx = x(rightindex+right:rightindex+2*right+5)+newerr;
rightx =[x(rightmaxind(end)+1:end),rightx];
[rhindmin,rhindmax,rhindzer] = extr(rightx);
rvmin = rightx(rhindmin(1));
rvmax = rightx(rhindmax(1));
rindmin = lx + (rhindmin(1)-length(x(rightmaxind(end):end)));
rindmax = lx + (rhindmax(1)-length(x(rightmaxind(end):end)));
end;
% 取拓延的点,拓延到两边(极大值与极小值)
tlmin = lindmin;
tlmax = lindmax;
trmin = rindmin;
trmax = rindmax;
% 延拓点上的取值(极大值与极小值)
zlmax = lvmax;
zlmin = lvmin;
zrmax = rvmax;
zrmin = rvmin;
% 完成延拓
tmin = [tlmin t(indmin) trmin];
tmax = [tlmax t(indmax) trmax];
zmin = [zlmin x(indmin) zrmin];
zmax = [zlmax x(indmax) zrmax];
end
tstar=min(tlmin tlmax);
tend=max(tlmin tlmax);
%---------------------------------------------------------------------------
function [tol,err,extrx1] = SelfAdapMatchWave(x0,x1)
% 一种自适应波形匹配方法
% 输入参数说明: x0 需要匹配的波形数据,长度为l
% x1 被匹配的波形数据,长度为l
%
% 输出参数说明: tol x0,x1的波形匹配度
% extrx1 x1的极值(extremum)序列,平移后的
%
len0 = length(x0);
len1 = length(x1);
if (len0 ~= len1)
error('x0,x1数据的长度不一致');
end;
% 数据x1平移操作,使x0,x1的极值重合
Max0 = max(x0);
Max1 = max(x1);
err = Max0 - Max1;
extrx1 = x1 + err;
x2 = extrx1;
% 对极值做匹配度的计算
tol = 0;
for i = 1:len0
erroo1(i) = x0(i)-x2(i);
tol = tol+(x0(i)-x2(i))*(x0(i)-x2(i));
end;
avg = (sum(abs(x0))+sum(abs(x2)))/(len0*2-2);
erravg = sum(abs(erroo1))/(len0-1);