cfs27=wpcoef(wpt1,[n,7]);
wpt2=wpdec(s3,n,ji);
rcfs30=wprcoef(wpt2,[n,0]);
rcfs31=wprcoef(wpt2,[n,1]);
rcfs32=wprcoef(wpt2,[n,2]);
rcfs33=wprcoef(wpt2,[n,3]);
for j=1:2^n
text(x(j),p_node(i),num2str(p_node(j),'%0.2f'),...
'HorizontalAlignment','center',...
'VerticalAlignment','bottom')
end
end
Matlab小波包分解后如何求各频带信号的能量值？
for i=1:2^n
pthi(i)= E(i)/E_total;%求每个节点的概率
end
pthi
%%%%以下提取小波包重构信号的故障诊断
clc
clear
fs=12e3;
loadzhengchang.mat;%加载正常信号
load 98.mat;
s1=X098_FE_time(19001:1:20000,1);% 正常情况下驱动端数据
% s3=X283_FE_time(1:1:10000,1);% 滚动体故障情况下驱动端数据
% load fan_out_12.mat;%加载外圈12点位置故障信号
load out_12.mat;%加载外圈12点位置故障信号
load 158.mat;
s4=X158_DE_time(14001:1:15000,1);% 外圈12点位置情况下驱动端数据
cfs15=wpcoef(wpt,[n,5]);
cfs16=wpcoef(wpt,[n,6]);
cfs17=wpcoef(wpt,[n,7]);
wpt1=wpdec(s2,n,ji);
rcfs20=wprcoef(wpt1,[n,0]);
rcfs21=wprcoef(wpt1,[n,1]);
rcfs22=wprcoef(wpt1,[n,2]);
% E
% disp('小波包分解总能量E_total');
E_total=sum(E); %求总能量
for i=1:2^n
pfir(i)= E(i)/E_total;%求每个节点的概率
end
pfir
% load a1_1-2;
load k302_1-4;
% Data=filter(b,a,Data);
end
% E = wenergy(wpt); only get the last layer
figure;
x=1:2^n;
bar(x,p_node);
title(['第',num2str(n),'层']);
axis([0 2^n 0 100]);
xlabel('结点');
ylabel('能量百分比/%');
Data=Data-mean(Data);
Data=detrend(Data);
wpt2=wpdec(Data,n,wpname); %对数据进行小波包分解
for i=1:2^n %wpcoef(wpt1,[n,i-1])是求第n层第i个节点的系数
E(i)=norm(wpcoef(wpt2,[n,i-1]),2);%求第i个节点的范数平方，其实也就是平方和
rex3=wprcoef(wpt,[n 3]);
;
subplot(5,1,1);
plot(x);
title('原始信号');
subplot(5,1,2);
plot(rex0);
title(['重构结点 ',num2str(n) ' 1',' 系数'])
subplot(5,1,3);
plot(rex1);
subplot(5,1,3);
plot(cfs1);
title(['结点 ',num2str(n) ' 2',' 系数'])
subplot(5,1,4);
plot(cfs2);
title(['结点 ',num2str(n) ' 3',' 系数'])
subplot(5,1,5);
plot(cfs3);
wpt3=wpdec(Data,n,wpname); %对数据进行小波包分解
for i=1:2^n %wpcoef(wpt1,[n,i-1])是求第n层第i个节点的系数
E(i)=norm(wpcoef(wpt3,[n,i-1]),2);%求第i个节点的范数平方，其实也就是平方和
end
E_total=sum(E); %求总能量
% s1=X098_FE_time(1:1:10000,1);% 正常情况下驱动端数据
load in.mat;%加载内圈故障信号
% load fan_in.mat;%加载风扇端内圈故障信号
load 106.mat;
s2=X106_DE_time(19001:1:20000,1);% 内圈故障信号情况下驱动端数据
% Plot wavelet packet tree (binary tree)
plot(wpt)
%% wavelet packet coefficients.default：use the front 4.
cfs0=wpcoef(wpt,[n 0]);
cfs1=wpcoef(wpt,[n 1]);
E(1:2^n )=0;
for i=1:2^n
E(i) = norm(wpcoef(wpt,[n,i-1]),2)^2; %20180604更新 原代码：E(i) = norm(wpcoef(wpt,[n,i-1]),2)
end
%求每个节点的概率
E_total=sum(E);
for i=1:2^n
p_node(i)= 100*E(i)/E_total;
小波包分解与重构
function wpt= wavelet_packetdecomposition_reconstruct( x,n,wpname )
%% 对信号进行小波包分解，得到节点的小波包系数。然后对每个节点系数进行重构。
% Decompose x at depth n with wpname wavelet ing Shannon entropy.
wpt1=wpdec(Data,n,wpname); %对数据进行小波包分解
for i=1:2^n %wpcoef(wpt1,[n,i-1])是求第n层第i个节点的系数
E(i)=norm(wpcoef(wpt1,[n,i-1]),2);%求第i个节点的范数平方，其实也就是平方和
end
% disp('每个节点的能量E(i)');
% s4=X305_FE_time(1:1:10000,1);% 外圈12点位置情况下驱动端数据
n=3;
ji='db5';
wpt=wpdec(s1,n,ji);%使用db5小波包分解正常信号到底三层，使用shannon熵
% plot(wpt);
rcfs10=wprcoef(wpt,[n,0]);%利用相应的频段进行重构
end
能量谱函数
function E = wavelet_energy_spectrum( wpt,n )
%% 计算每一层每一个节点的能量
% wpt-wavelet packet tree
% n-第n层能量
%
% Author hubery_zhang
% Date 20170714
%%
% 求第n层第i个节点的系数
cfs32=wpcoef(wpt2,[n,2]);
cfs33=wpcoef(wpt2,[n,3]);
cfs34=wpcoef(wpt2,[n,4]);
cfs35=wpcoef(wpt2,[n,5]);
cfs36=wpcoef(wpt2,[n,6]);
cfs37=wpcoef(wpt2,[n,7]);
clear;
clc;
n=3;
wpname='db3';
% [b,a]=butter(8,[5/100 99/100]);
% load a1_1-1;
load k301_1-4;
% Data=filter(b,a,Data);
Data=Data-mean(Data);
Data=detrend(Data);
rcfs11=wprcoef(wpt,[n,1]);
rcfs12=wprcoef(wpt,[n,2]);
rcfs13=wprcoef(wpt,[n,3]);
rcfs14=wprcoef(wpt,[n,4]);
rcfs15=wprcoef(wpt,[n,5]);
rcfs16=wprcoef(wpt,[n,6]);
title(['重构结点 ',num2str(n) ' 2',' 系数'])
subplot(5,1,4);
plot(rex2);
title(['重构结点 ',num2str(n) ' 3',' 系数'])
subplot(5,1,5);
plot(rex3);
title(['重构结点 ',num2str(n) ' 4',' 系数'])