wte_six_signal_compare.m

六个信号的小波时间熵进行比较

%% load signal
clear;
load AGnzbx;
s1=AGnzbx';
x1=s1(1,1:83000);
load AHtcbx;
s2=AHtcbx';
x2=s2(1,1:83000);
load AIxjbx;
s3=AIxjbx';
x3=s3(1,1:83000);
load ALpbpx;
s4=ALpbpx';
x4=s4(1,1:83000);
load AJqlbw;
s5=AJqlbw';
x5=s5(1,1:83000);
load AKcgtx;
s6=AKcgtx';
x6=s6(1,1:83000);

%%
%%开始小波多尺度分解
[c1,l1]=wavedec(x1,5,'db5');
x1_d5=wrcoef('d',c1,l1,'db5',5);
[c2,l2]=wavedec(x2,5,'db5');
x2_d5=wrcoef('d',c2,l2,'db5',5);
[c3,l3]=wavedec(x3,5,'db5');
x3_d5=wrcoef('d',c3,l3,'db5',5);
[c4,l4]=wavedec(x4,5,'db5');
x4_d5=wrcoef('d',c4,l4,'db5',5);
[c5,l5]=wavedec(x5,5,'db5');
x5_d5=wrcoef('d',c5,l5,'db5',5);
[c6,l6]=wavedec(x6,5,'db5');
x6_d5=wrcoef('d',c6,l6,'db5',5);




win=13000;  %%窗宽
len=length(x1_d5);
dis=1000;
M=(len-win)/dis;
pin_number=10;%%将最大值与最小值之间平均分为10个值那么一共就有9段
duan_number=pin_number-1;
w(M,win)=0;
min_w(1,M)=0;
max_w(1,M)=0;
%%制造滑动窗
for i=1:M
    w1(i,:)=x1_d5(1+i*dis:win+i*dis);
    w2(i,:)=x2_d5(1+i*dis:win+i*dis);
    w3(i,:)=x3_d5(1+i*dis:win+i*dis);
    w4(i,:)=x4_d5(1+i*dis:win+i*dis);
    w5(i,:)=x5_d5(1+i*dis:win+i*dis);
    w6(i,:)=x6_d5(1+i*dis:win+i*dis);
end
Zl_w1(M,pin_number)=0;
Zl_w2(M,pin_number)=0;
Zl_w3(M,pin_number)=0;
Zl_w4(M,pin_number)=0;
Zl_w5(M,pin_number)=0;
Zl_w6(M,pin_number)=0;
%%制造四个不相交的区间
for i=1:M
    min_w1(i)=min(w1(i,:));  %s0
    max_w1(i)=max(w1(i,:));  %sL
    Zl_w1(i,:)=linspace(min_w1(i),max_w1(i),pin_number);
    min_w2(i)=min(w2(i,:));
    max_w2(i)=max(w2(i,:));
    Zl_w2(i,:)=linspace(min_w2(i),max_w2(i),pin_number);
    min_w3(i)=min(w3(i,:));
    max_w3(i)=max(w3(i,:));
    Zl_w3(i,:)=linspace(min_w3(i),max_w3(i),pin_number);
    min_w4(i)=min(w4(i,:));
    max_w4(i)=max(w4(i,:));
    Zl_w4(i,:)=linspace(min_w4(i),max_w4(i),pin_number);
    min_w5(i)=min(w5(i,:));
    max_w5(i)=max(w5(i,:));
    Zl_w5(i,:)=linspace(min_w5(i),max_w5(i),pin_number);
    min_w6(i)=min(w6(i,:));
    max_w6(i)=max(w6(i,:));
    Zl_w6(i,:)=linspace(min_w6(i),max_w6(i),pin_number);
end
%%计算落于区间的概率(一共有9段)
P1(M,duan_number)=0;
P2(M,duan_number)=0;
P3(M,duan_number)=0;
P4(M,duan_number)=0;
P5(M,duan_number)=0;
P6(M,duan_number)=0;
for i=1:M
    for j=1:duan_number
        maxtri_low=ones(1,win).*Zl_w1(i,j);
        maxtri_high=ones(1,win).*Zl_w1(i,j+1);
        compare_result_low=w1(i,:)>=maxtri_low;
        compare_result_high=w1(i,:)<=maxtri_high;
        compare_result=compare_result_low.*compare_result_high;
        number=sum(compare_result);
        P1(i,j)=number/win;
    end
end
for i=1:M
    for j=1:duan_number
        maxtri_low=ones(1,win).*Zl_w2(i,j);
        maxtri_high=ones(1,win).*Zl_w2(i,j+1);
        compare_result_low=w2(i,:)>=maxtri_low;
        compare_result_high=w2(i,:)<maxtri_high;
        compare_result=compare_result_low.*compare_result_high;
        number=sum(compare_result);
        P2(i,j)=number/win;
    end
end
for i=1:M
    for j=1:duan_number
        maxtri_low=ones(1,win).*Zl_w3(i,j);
        maxtri_high=ones(1,win).*Zl_w3(i,j+1);
        compare_result_low=w3(i,:)>=maxtri_low;
        compare_result_high=w3(i,:)<=maxtri_high;
        compare_result=compare_result_low.*compare_result_high;
        number=sum(compare_result);
        P3(i,j)=number/win;
    end
end
for i=1:M
    for j=1:duan_number
        maxtri_low=ones(1,win).*Zl_w4(i,j);
        maxtri_high=ones(1,win).*Zl_w4(i,j+1);
        compare_result_low=w4(i,:)>=maxtri_low;
        compare_result_high=w4(i,:)<=maxtri_high;
        compare_result=compare_result_low.*compare_result_high;
        number=sum(compare_result);
        P4(i,j)=number/win;
    end
end
for i=1:M
    for j=1:duan_number
        maxtri_low=ones(1,win).*Zl_w5(i,j);
        maxtri_high=ones(1,win).*Zl_w5(i,j+1);
        compare_result_low=w5(i,:)>=maxtri_low;
        compare_result_high=w5(i,:)<=maxtri_high;
        compare_result=compare_result_low.*compare_result_high;
        number=sum(compare_result);
        P5(i,j)=number/win;
    end
end
for i=1:M
    for j=1:duan_number
        maxtri_low=ones(1,win).*Zl_w6(i,j);
        maxtri_high=ones(1,win).*Zl_w6(i,j+1);
        compare_result_low=w6(i,:)>=maxtri_low;
        compare_result_high=w6(i,:)<=maxtri_high;
        compare_result=compare_result_low.*compare_result_high;
        number=sum(compare_result);
        P6(i,j)=number/win;
    end
end

%%计算小波时间熵

result1(M,duan_number)=0;
result2(M,duan_number)=0;
result3(M,duan_number)=0;
result4(M,duan_number)=0;
result5(M,duan_number)=0;
result6(M,duan_number)=0;
for i=1:M
    for j=1:duan_number
        if P1(i,j)==0
            result1(i,j)=0;
        else 
            result1(i,j)=P1(i,j)*log(P1(i,j));
        end
    end
end
for i=1:M
    for j=1:duan_number
        if P2(i,j)==0
            result2(i,j)=0;
        else 
            result2(i,j)=P2(i,j)*log(P2(i,j));
        end
    end
end
for i=1:M
    for j=1:duan_number
        if P3(i,j)==0
            result3(i,j)=0;
        else 
            result3(i,j)=P3(i,j)*log(P3(i,j));
        end
    end
end
for i=1:M
    for j=1:duan_number
        if P4(i,j)==0
            result4(i,j)=0;
        else 
            result4(i,j)=P4(i,j)*log(P4(i,j));
        end
    end
end
for i=1:M
    for j=1:duan_number
        if P5(i,j)==0
            result5(i,j)=0;
        else 
            result5(i,j)=P5(i,j)*log(P5(i,j));
        end
    end
end
for i=1:M
    for j=1:duan_number
        if P6(i,j)==0
            result6(i,j)=0;
        else 
            result6(i,j)=P6(i,j)*log(P6(i,j));
        end
    end
end

WTE_nz(1,M)=0;
WTE_tc(1,M)=0;
WTE_xj(1,M)=0;
WTE_pb(1,M)=0;
WTE_ql(1,M)=0;
WTE_cg(1,M)=0;

for i=1:M
    WTE_nz(1,i)=-(sum(result1(i,1:duan_number)));
end
for i=1:M
    WTE_tc(1,i)=-(sum(result2(i,1:duan_number)));
end
for i=1:M
    WTE_xj(1,i)=-(sum(result3(i,1:duan_number)));
end
for i=1:M
    WTE_pb(1,i)=-(sum(result4(i,1:duan_number)));
end
for i=1:M
    WTE_ql(1,i)=-(sum(result5(i,1:duan_number)));
end
for i=1:M
    WTE_cg(1,i)=-(sum(result6(i,1:duan_number)));
end
m=1:M;
t=win/2+m*dis;
figure(1);
hold on;
plot(t,WTE_nz,'m-','LineWidth',2);

plot(t,WTE_tc,'b--','LineWidth',2);

plot(t,WTE_xj,'r-.','LineWidth',2);

plot(t,WTE_pb,'k:','LineWidth',2);

plot(t,WTE_ql,'k--','LineWidth',2);

plot(t,WTE_cg,'b-.','LineWidth',2);

legend('AGnzbx','AHtcbx','AIxjbx','ALpbpx','AJqlbw','AKcgtx')