《MATLAB之经验模态分解EMD》由会员分享，可在线阅读，更多相关《MATLAB之经验模态分解EMD（18页珍藏版）》请在金锄头文库上搜索。
1、MATLAB之经验模态分解之经验模态分解EMDfunction imf,ort,nbits = emd(varargin)x,t,sd,sd2,tol,MODE_COMPLEX,ndirs,display_sifting,sdt,sd2t,r,imf,k,nbit,NbIt,MAXITERATIONS,FIXE,FIXE_H,MAXMODES,INTERP,mask= init(varargin:);if display_sifting fig_h = figure;endwhile stop_EMD(r,MODE_COMPLEX,ndirs) & (k m = r; mp = m; if F
2、IXE % 如果设定了迭代次数 stop_sift,moyenne = stop_sifting_fixe(t,m,INTERP,MODE_COMPLEX,ndirs); elseif FIXE_H stop_count = 0; stop_sift,moyenne = stop_sifting_fixe_h(t,m,INTERP,stop_count,FIXE_H,MODE_COMPLEX,ndirs); else stop_sift,moyenne = stop_sifting(m,t,sd,sd2,tol,INTERP,MODE_COMPLEX,ndirs); end if (max(a
3、bs(m) if stop_sift warning(emd:warning,forced stop of EMD : too small amplitude) else disp(forced stop of EMD : too small amplitude) end break end while stop_sift & nbit if(MODE_COMPLEX & nbitMAXITERATIONS/5 & mod(nbit,floor(MAXITERATIONS/10)=0 & FIXE & nbit 100) disp(mode ,int2str(k), iteration ,in
4、t2str(nbit) if exist(s,var) disp(stop parameter mean value : ,num2str(s) end im,iM = extr(m); disp(int2str(sum(m(im) 0), minima 0; ,int2str(sum(m(iM) end m = m - moyenne; if FIXE stop_sift,moyenne = stop_sifting_fixe(t,m,INTERP,MODE_COMPLEX,ndirs); elseif FIXE_H stop_sift,moyenne,stop_count = stop_s
5、ifting_fixe_h(t,m,INTERP,stop_count,FIXE_H,MODE_COMPLEX,ndirs); else stop_sift,moyenne,s = stop_sifting(m,t,sd,sd2,tol,INTERP,MODE_COMPLEX,ndirs); end if display_sifting & MODE_COMPLEX NBSYM = 2; indmin,indmax = extr(mp); tmin,tmax,mmin,mmax = boundary_conditions(indmin,indmax,t,mp,mp,NBSYM); envmin
6、p = interp1(tmin,mmin,t,INTERP); envmaxp = interp1(tmax,mmax,t,INTERP); envmoyp = (envminp+envmaxp)/2; if FIXE | FIXE_H display_emd_fixe(t,m,mp,r,envminp,envmaxp,envmoyp,nbit,k,display_sifting) else sxp = 2*(abs(envmoyp)./(abs(envmaxp-envminp); sp = mean(sxp); display_emd(t,m,mp,r,envminp,envmaxp,en
7、vmoyp,s,sp,sxp,sdt,sd2t,nbit,k,display_sifting,stop_sift) end end mp = m; nbit = nbit+1; NbIt = NbIt+1; if (nbit=(MAXITERATIONS-1) & FIXE & nbit 100) if exist(s,var) warning(emd:warning,forced stop of sifting : too many iterations. mode ,int2str(k),. stop parameter mean value : ,num2str(s) else warn
8、ing(emd:warning,forced stop of sifting : too many iterations. mode ,int2str(k),.) end end end imf(k,:) = m; if display_sifting disp(mode ,int2str(k), stored) end nbits(k) = nbit; k = k+1; r = r - m; nbit = 0;end if any(r) & any(mask) imf(k,:) = r;endort = io(x,imf);if display_sifting closeendendfunc
9、tion stop = stop_EMD(r,MODE_COMPLEX,ndirs)if MODE_COMPLEX for k = 1:ndirs phi = (k-1)*pi/ndirs; indmin,indmax = extr(real(exp(i*phi)*r); ner(k) = length(indmin) + length(indmax); end stop = any(ner else indmin,indmax = extr(r); ner = length(indmin) + length(indmax); stop = ner endendfunction envmoy,
10、nem,nzm,amp = mean_and_amplitude(m,t,INTERP,MODE_COMPLEX,ndirs)NBSYM = 2;if MODE_COMPLEX switch MODE_COMPLEX case 1 for k = 1:ndirs phi = (k-1)*pi/ndirs; y = real(exp(-i*phi)*m); indmin,indmax,indzer = extr(y); nem(k) = length(indmin)+length(indmax); nzm(k) = length(indzer); tmin,tmax,zmin,zmax = bo
11、undary_conditions(indmin,indmax,t,y,m,NBSYM); envmin(k,:) = interp1(tmin,zmin,t,INTERP); envmax(k,:) = interp1(tmax,zmax,t,INTERP); end envmoy = mean(envmin+envmax)/2,1); if nargout 3 amp = mean(abs(envmax-envmin),1)/2; end case 2 for k = 1:ndirs phi = (k-1)*pi/ndirs; y = real(exp(-i*phi)*m); indmin
12、,indmax,indzer = extr(y); nem(k) = length(indmin)+length(indmax); nzm(k) = length(indzer); tmin,tmax,zmin,zmax = boundary_conditions(indmin,indmax,t,y,y,NBSYM); envmin(k,:) = exp(i*phi)*interp1(tmin,zmin,t,INTERP); envmax(k,:) = exp(i*phi)*interp1(tmax,zmax,t,INTERP); end envmoy = mean(envmin+envmax
13、),1); if nargout 3 amp = mean(abs(envmax-envmin),1)/2; end endelse indmin,indmax,indzer = extr(m); nem = length(indmin)+length(indmax); nzm = length(indzer); tmin,tmax,mmin,mmax = boundary_conditions(indmin,indmax,t,m,m,NBSYM); envmin = interp1(tmin,mmin,t,INTERP); envmax = interp1(tmax,mmax,t,INTER
14、P); envmoy = (envmin+envmax)/2; if nargout 3 amp = mean(abs(envmax-envmin),1)/2; % 计算包络幅度 endendendfunction stop,envmoy,s = stop_sifting(m,t,sd,sd2,tol,INTERP,MODE_COMPLEX,ndirs)try envmoy,nem,nzm,amp = mean_and_amplitude(m,t,INTERP,MODE_COMPLEX,ndirs); sx = abs(envmoy)./amp; s = mean(sx); stop = (m
15、ean(sx sd) tol | any(sx sd2) & (all(nem 2); if MODE_COMPLEX stop = stop & (abs(nzm-nem)1); endcatch stop = 1; envmoy = zeros(1,length(m); s = NaN;endendfunction stop,moyenne= stop_sifting_fixe(t,m,INTERP,MODE_COMPLEX,ndirs)try moyenne = mean_and_amplitude(m,t,INTERP,MODE_COMPLEX,ndirs); stop = 0;cat
16、ch moyenne = zeros(1,length(m); stop = 1;endendfunction stop,moyenne,stop_count= stop_sifting_fixe_h(t,m,INTERP,stop_count,FIXE_H,MODE_COMPLEX,ndirs)try moyenne,nem,nzm = mean_and_amplitude(m,t,INTERP,MODE_COMPLEX,ndirs); if (all(abs(nzm-nem)1) stop = 0; stop_count = 0; else stop_count = stop_count+1; stop = (stop_count = FIXE_H); endcatch moyenne = zeros(1,length(m); stop = 1;endendfunction display_emd(t,m,mp,r,envmin,envmax,envmoy,s,sb,sx,sdt,sd2t,nbit,k,display_sifting,stop_sift)subplot(4,1,1
