Desktop.rar_EEMD matlab_differenceq5x_pulse shaping

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % 由淘宝同学店铺---面对面数学---郑一原创改进。2017年1月8日 % % 旺旺：正一算法程序；旺旺：面对面数学 % % 只提供拍得《经验模态分解大类算法及其应用与MATlAB程序分析视频》课程的君子学者学习参考% % 请保留这里的文字记录，手下留情。不要上传到网络。 % % %%%%%%%%%%%%%%%% % % PPT课件：EMD10_基于集合经验模态分解(EEMD)建立最优降噪整形算法.ppt % % 参考文献引用信息：郑一.经验模态分解大类算法与应用及MATLAB程序视频. % % https://shop108509999.taobao.com. 2017-01-08. % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% clc clear close all; %%% 导入数据，合成信号，显示信号 danweimaichong %输出有xx,yy(1,4856) [xxr xxc]=size(xx); figure(1) subplot(4,1,1); plot(xx,yy); axis([0 120 -0.5 1.2]); xlabel('时间'); ylabel('单位脉冲'); title('单位脉冲信号'); ysin = (1+0.4*sin(2*pi*5*xx)).*sin(2*pi*25*xx.^2+2*pi*9*xx);%调幅调频信号 %ysin=0.3*sin(2*pi*0.3*xx+21); subplot(4,1,2); plot(xx,ysin); ylim([-0.6 0.6]) ylabel('调幅调频'); title('调幅调频信号'); % % ywgn = wgn(xxr,xxc,0.6); % stdywgn=std(ywgn); % mywgn=mean(ywgn); % a=0; % b=1.2; %b=0.6; % ygs=a+b.*((ywgn-mywgn)/stdywgn); % stdygs=std(ygs); load ('YGS.mat'); subplot(4,1,3); plot(xx,ygs); ylabel('高斯噪声'); title('高斯噪声信号'); xy = yy+ysin+ygs; %合成信号 subplot(4,1,4); plot(xx,xy); xlabel('时间'); ylabel('含噪方波合成信号'); title('含噪方波合成信号'); stdxy=std(xy); %%% EEMD分解，构建10个带通滤波信号 Nstd=0.4; NE=50; TNM= 10; % c1=eemd_ZHYzhushi(xy,Nstd,NE); c1=ceemd(xy,Nstd,NE,TNM); c=c1'; %c第1行是原始数据,最后1行是余项。c是13行4856列矩阵 [m,n]=size(c); CEE_ifm11=c(11,:); CEE_ifm10=c(10,:)+c(11,:); CEE_ifm9=c(9,:)+c(10,:)+c(11,:); CEE_ifm8=c(8,:)+c(9,:)+c(10,:)+c(11,:); CEE_ifm7=c(7,:)+c(8,:)+c(9,:)+c(10,:)+c(11,:); CEE_ifm6=c(6,:)+c(7,:)+c(8,:)+c(9,:)+c(10,:)+c(11,:); CEE_ifm5=c(5,:)+c(6,:)+c(7,:)+c(8,:)+c(9,:)+c(10,:)+c(11,:); CEE_ifm4=c(4,:)+c(5,:)+c(6,:)+c(7,:)+c(8,:)+c(9,:)+c(10,:)+c(11,:); CEE_ifm3=c(3,:)+c(4,:)+c(5,:)+c(6,:)+c(7,:)+c(8,:)+c(9,:)+c(10,:)+c(11,:); CEE_ifm2=c(2,:)+c(3,:)+c(4,:)+c(5,:)+c(6,:)+c(7,:)+c(8,:)+c(9,:)+c(10,:)+c(11,:); % CEE_ifm11=xy-c(2,:)-c(3,:)-c(4,:)-c(5,:)-c(6,:)-c(7,:)-c(8,:)-c(9,:)-c(10,:)-c(11,:); %低通 % CEE_ifm10=xy-c(2,:)-c(3,:)-c(4,:)-c(5,:)-c(6,:)-c(7,:)-c(8,:)-c(9,:)-c(10,:); % CEE_ifm9=xy-c(2,:)-c(3,:)-c(4,:)-c(5,:)-c(6,:)-c(7,:)-c(8,:)-c(9,:); % CEE_ifm8=xy-c(2,:)-c(3,:)-c(4,:)-c(5,:)-c(6,:)-c(7,:)-c(8,:); % CEE_ifm7=xy-c(2,:)-c(3,:)-c(4,:)-c(5,:)-c(6,:)-c(7,:); % CEE_ifm6=xy-c(2,:)-c(3,:)-c(4,:)-c(5,:)-c(6,:); % CEE_ifm5=xy-c(2,:)-c(3,:)-c(4,:)-c(5,:); % CEE_ifm4=xy-c(2,:)-c(3,:)-c(4,:); % CEE_ifm3=xy-c(2,:)-c(3,:); % CEE_ifm2=xy-c(2,:); %% 画出原始信号图形及部分IMF分量及最后一个剩余分量residual的图形 figure(2) subplot(6,1,1) plot(xx,xy) set(gca,'fontname','simher') set(gca,'fontsize',14.0) ylabel(['Signal']) title('固有模态分量分解') %axis([-2,2,0,1.2]) subplot(6,1,2); set(gcf,'color','w') plot(xx,c(2,:),'k') %画IMF1 set(gca,'fontname','simher') set(gca,'fontsize',14.0) ylabel(['IMF ',int2str(1)]) subplot(6,1,3); set(gcf,'color','w') plot(xx,c(4,:),'k') %画IMF3 set(gca,'fontname','simher') set(gca,'fontsize',14.0) ylabel(['IMF ',int2str(3)]) subplot(6,1,4); set(gcf,'color','w') plot(xx,c(6,:),'k') %画IMF5 set(gca,'fontname','simher') set(gca,'fontsize',14.0) ylabel(['IMF ',int2str(5)]) subplot(6,1,5); set(gcf,'color','w') plot(xx,c(5,:),'k') %画IMF4 set(gca,'fontname','simher') set(gca,'fontsize',14.0) ylabel(['IMF ',int2str(4)]) subplot(6,1,6); set(gcf,'color','w') plot(xx,c(9,:),'k') %画IMF8 set(gca,'fontname','simher') set(gca,'fontsize',14.0) ylabel(['IMF ',int2str(8)]) figure(3) subplot(3,1,1); set(gcf,'color','w') plot(xx,c(10,:),'k') %画IMF9 set(gca,'fontname','simher') set(gca,'fontsize',14.0) ylabel(['IMF ',int2str(9)]) subplot(3,1,2); set(gcf,'color','w') plot(xx,c(11,:),'k') %画IMF10 set(gca,'fontname','simher') set(gca,'fontsize',14.0) ylabel(['IMF ',int2str(10)]) subplot(3,1,3); set(gcf,'color','w') plot(xx,c(12,:),'k') %画余项 set(gca,'fontname','simher') set(gca,'fontsize',14.0) ylabel('RS ') %% 画出EE_ifm3，5,6,8构成的带通滤波算法的效果 figure(4); subplot(2,2,1); plot(xx,CEE_ifm3,'b-') % title('算法CEE-ifm3滤波效果'); xlabel('时间');ylabel('信号幅值'); subplot(2,2,2); plot(xx,CEE_ifm4,'b-') % title('算法CEE-ifm4滤波效果'); xlabel('时间');ylabel('信号幅值'); SNR_singlech(CEE_ifm4,xy); subplot(2,2,3); plot(xx,CEE_ifm5,'b-') % title('算法CEE-ifm5滤波效果'); xlabel('时间');ylabel('信号幅值'); SNR_singlech(CEE_ifm5,xy); subplot(2,2,4); plot(xx,CEE_ifm6,'b-') % title('算法CEE-ifm6滤波效果'); xlabel('时间');ylabel('信号幅值'); SNR_singlech(CEE_ifm6,xy); %%% 平顶脉冲整形处理 yuzhi1=0.25; yuzhi2=0.75; Zhengxing2 = Zhengxing(CEE_ifm2,xxc,yuzhi1,yuzhi2); Zhengxing3 = Zhengxing(CEE_ifm3,xxc,yuzhi1,yuzhi2); Zhengxing4 = Zhengxing(CEE_ifm4,xxc,yuzhi1,yuzhi2); Zhengxing5 = Zhengxing(CEE_ifm5,xxc,yuzhi1,yuzhi2); Zhengxing6 = Zhengxing(CEE_ifm6,xxc,yuzhi1,yuzhi2); Zhengxing7 = Zhengxing(CEE_ifm7,xxc,yuzhi1,yuzhi2); Zhengxing8 = Zhengxing(CEE_ifm8,xxc,yuzhi1,yuzhi2); Zhengxing9 = Zhengxing(CEE_ifm9,xxc,yuzhi1,yuzhi2); Zhengxing10 = Zhengxing(CEE_ifm10,xxc,yuzhi1,yuzhi2); Zhengxing11 = Zhengxing(CEE_ifm11,xxc,yuzhi1,yuzhi2); %%% 综合显示部分平顶整形效果 figure(5); subplot(2,3,1); plot(xx,Zhengxing2) axis([0 120 -0.5 1.5]); title('算法EE-ifm2整形结果'); xlabel('时间'); ylabel('信号幅值'); subplot(2,3,2); plot(xx,Zhengxing4) axis([0 120 -0.5 1.5]); title('算法EE-ifm4整形结果'); xlabel('时间'); ylabel('信号幅值') subplot(2,3,3); plot(xx,Zhengxing5) axis([0 120 -0.5 1.5]); title('算法EE-ifm5整形结果'); xlabel('时间'); ylabel('信号幅值') subplot(2,3,4); plot(xx,Zhengxing6) axis([0 120 -0.5 1.5]); title('算法EE-ifm6整形结果'); xlabel('时间'); ylabel('信号幅值'); subplot(2,3,5); plot(xx,Zhengxing7) axis([0 120 -0.5 1.5]); title('算法EE-ifm7整形结果'); xlabel('时间'); ylabel('信号幅值') subplot(2,3,6); plot(xx,Zhengxing9) axis([0 120 -0.5 1.5]); title('算法EE-ifm9整形结果'); xlabel('时间'); ylabel('信号幅值'); %%% 求带通滤波算法后数据与原始信号的均方误差，计算算法逼近度 dt2std=sum((xy-CEE_ifm2).^2)/(xxr*xxc); dt3std=sum((xy-CEE_ifm3).^2)/(xxr*xxc); dt4std=sum((xy-CEE_ifm4).^2)/(xxr*xxc); dt5std=sum((xy-CEE_ifm5).^2)/(xxr*xxc); dt6std=sum((xy-CEE_ifm6).^2)/(xxr*xxc); dt7std=sum((xy-CEE_ifm7).^2)/(xxr*xxc); dt8std=sum((xy-CEE_ifm8).^2)/(xxr*xxc); dt9std=sum((xy-CEE_ifm9).^2)/(xxr*xxc); dt10std=sum((xy-CEE_ifm10).^2)/(xxr*xxc); dt11std=sum((xy-CEE_ifm11).^2)/(xxr*xxc); MSE_2=1/(dt2std+1); MSE_3=1/(dt3std+1); MSE_4=1/(dt4std+1); MSE_5=1/(dt5std+1); MSE_6=1/(dt6std+1); MSE_7=1/(dt7std+1); MSE_8=1/(dt8std+1); MSE_9=1/(dt9std+1); MSE_10=1/(dt10std+1); MSE_11=1/(dt11std+1); %%% 求带通滤波整形后数据与原始单位脉冲信号的相关系数，计算算法线性相关度 %计算不含噪脉冲信号yy和滤波平顶信号maichongZhengxing2的相关性。P小于0.05时，R有意义。 [ROH2,P2]=corr(yy',Zhengxing2'); [ROH3,P3]=corr(yy',Zhengxing3'); [ROH4,P4]=corr(yy',Zhengxing4'); [ROH5,P5]=corr(yy',Zhengxing5'); [ROH6,P6]=corr(yy',Zhengxing6'); [ROH7,P7]=corr(yy',Zhengxing7'); [ROH8,P8]=corr(yy',Zhengxing8'); [ROH9,P9]=corr(yy',Zhengxing9'); [ROH10,P10]=corr(yy',Zhengxing10'); [ROH11,P11]=corr(yy',Zhengxing11'); ROH9=0; ROH10=0; ROH11=0; %%% 计算目标函数值 alpha=0.4; fdt2s=alpha*MSE_2+(1-alpha)*ROH2; fdt3s=alpha*MSE_3+(1-alpha)*ROH3; fdt4s=alpha*MSE_4+(1-alpha)*ROH4; fdt5s=alpha*MSE_5+(1-alpha)*ROH5; fdt6s=alpha*MSE_6+(1-alpha)*ROH6; fdt7s=alpha*MSE_7+(1-alpha)*ROH7; fdt8s=alpha*MSE_8+(1-alpha)*ROH8; fdt9s=alpha*MSE_9+(1-alpha)*ROH9; fdt10s=alpha*MSE_10+(1-alpha)*ROH