局部均值分解（local mean decomposition , LMD）方法同经验模态分解方法（EMD）一样，也是一种自适应信号处理方法。LMD通过改变信号分解过程能有效改进EMD方法存在的包络拟合不准确、边界处发散等问题.
代码如下：
clear all; clc; % 这里是仿真信号，可改成自己的信号 fs=2000; N=2048; n=0:N-1; t=n/fs; x=(1+0.5*t).*sin(2*pi.*20*t)+2*cos(2*pi.*3*t); %绘制仿真信号和其频谱图 figure(1) subplot(211) plot(t,x) subplot(212) y2=x; L=length(y2); NFFT = 2^nextpow2(L); Y = fft(y2,NFFT)/L; f = fs/2*linspace(0,1,NFFT/2); plot(f,2*abs(Y(1:NFFT/2))) % 局域均值分析 x=x'; c = x'; N = length(x); A = ones(1,N); PF = []; aii = 2*A; while(1) si = c; a = 1; while(1) h = si; maxVec = []; minVec = []; % 寻找极大值点和极小值点 for i = 2: N - 1 if h (i - 1) < h (i) & h (i) > h (i + 1) maxVec = [maxVec i]; end if h (i - 1) > h (i) & h (i) < h (i + 1) minVec = [minVec i]; end end % 检查是否有残余 if (length (maxVec) + length (minVec)) < 2 break; end % 原始信号中的两边两个点的判断 lenmax=length(maxVec); lenmin=length(minVec); %先是左边这个点 if h(1)>0 if(maxVec(1)<minVec(1)) yleft_max=h(maxVec(1)); yleft_min=-h(1); else yleft_max=h(1); yleft_min=h(minVec(1)); end else if (maxVec(1)<minVec(1)) yleft_max=h(maxVec(1)); yleft_min=h(1); else yleft_max=-h(1); yleft_min=h(minVec(1)); end end %然后判断右边这个点 if h(N)>0 if(maxVec(lenmax)<minVec(lenmin)) yright_max=h(N); yright_min=h(minVec(lenmin)); else yright_max=h(maxVec(lenmax)); yright_min=-h(N); end else if(maxVec(lenmax)<minVec(lenmin)) yright_max=-h(N); yright_min=h(minVec(lenmin)); else yright_max=h(maxVec(lenmax)); yright_min=h(N); end end %使用三次样条插值方法，对极大值向量和极小值向量进行插值 %spline interpolate maxEnv=spline([1 maxVec N],[yleft_max h(maxVec) yright_max],1:N); minEnv=spline([1 minVec N],[yleft_min h(minVec) yright_min],1:N); mm = (maxEnv + minEnv)/2;%得到局部均值函数 aa = abs(maxEnv - minEnv)/2;%得到包络函数 mmm = mm; aaa = aa; preh = h; h = h-mmm;%从原始信号中分离处局部均值函数 si = h./aaa;%对分离出的信号进行解调 a = a.*aaa; aii = aaa; B = length(aii); C = ones(1,B); bb = norm(aii-C);%返回aii-C的最大奇异值，aii就是那个包络函数 if(bb < 1000)%如果bb<1000，就得到了纯调频函数 break; end end %分解1个Pf分量在这结束 pf = a.*si;%包络函数和纯调频函数相乘，得到PF分量 PF = [PF; pf]; bbb = length (maxVec) + length (minVec); % 简单的一个结束分解的条件 if (length (maxVec) + length (minVec)) < 20 break; end c = c-pf; end m=x'-PF(1,:)-PF(2,:);%如果分解出2个，从原始信号中减去，得到残余分量 figure(2); subplot(4,1,1),plot(n,x),ylabel('X(t)'); subplot(4,1,2),plot(n,PF(1,:)),ylabel('PF_1(t)'); subplot(4,1,3),plot(n,PF(2,:)),ylabel('PF_2(t)'); subplot(4,1,4),plot(n,m),ylabel('u(t)'); xlabel('Time / S');最后分解出来的图：
