二维经验模态分解的matlab代码，并绘制图像

二维经验模态分解（2D-EMD）是一种用于分解二维信号的方法，可以用于图像处理等领域。下面是一个基于MATLAB的二维经验模态分解的代码示例：

function [imf,residue] = emd2(X)
%二维EMD分解
%输入：X-待分解信号
%输出：imf-分解后的各个本征模态函数（IMF）
%      residue-剩余部分

imf = [];
while 1
    h = X;
    SD = 1;
    while SD > 0.1
        %计算X的2D-Hilbert变换
        H = hilbert2(X);
        %计算X的局部平均
        Havg = (circshift(H,[0 -1])+circshift(H,[0 1])+circshift(H,[-1 0])+circshift(H,[1 0])+circshift(H,[-1 -1])+circshift(H,[1 -1])+circshift(H,[-1 1])+circshift(H,[1 1]))/8;
        %计算X的局部极值
        Hdiff = H-Havg;
        [maxtab,~] = peakdet2(Hdiff,0);
        [mintab,~] = peakdet2(-Hdiff,0);
        %计算X的局部均值
        mmax = mean2(maxtab(:,3));
        mmin = mean2(mintab(:,3));
        %计算SD值
        SD = abs(mmax-mmin);
        %更新X
        X = X - Havg;
    end
    imf = cat(3,imf,X);
    %计算剩余部分
    residue = h - X;
    if sum(sum(abs(residue))) < 1e-6*sum(sum(abs(h)))
        break;
    end
    X = residue;
end
imf = cat(3,imf,residue);
end

function [H] = hilbert2(X)
%计算二维信号的Hilbert变换
%输入：X-待变换信号
%输出：H-变换后的信号

[m,n] = size(X);
H1 = hilbert(X')';
H2 = hilbert(H1);
H = zeros(m,n);
for i=1:m
    for j=1:n
        H(i,j) = H2(i,j) - 1i*H1(i,j);
    end
end
end

function [maxtab,mintab]=peakdet2(v, delta, x)
%找到二维信号中的局部极值
%输入：v-待寻找的信号
%      delta-峰值的阈值
%      x-信号的x轴值（可选）
%输出：maxtab-峰值的位置和幅值
%      mintab-谷值的位置和幅值

maxtab = [];
mintab = [];

if nargin < 3
    x = 1:length(v);
end

v(isnan(v)) = 0;
v = v(:); x = x(:);

if length(delta) > 1
    delta = delta(2);
end

if delta <= 0
    error('Delta must be positive.');
end

mn = Inf; mx = -Inf;
mnpos = NaN; mxpos = NaN;

lookformax = 1;
for i=1:length(v)
    this = v(i);
    if this > mx
        mx = this; mxpos = x(i);
    end
    if this < mn
        mn = this; mnpos = x(i);
    end

    if lookformax
        if this < mx-delta
            maxtab = [maxtab ; mxpos mx]; mn = this; mnpos = x(i); lookformax = 0;
        end
    else
        if this > mn+delta
            mintab = [mintab ; mnpos mn]; mx = this; mxpos = x(i); lookformax = 1;
        end
    end
end

if ~isempty(mintab) && ~isempty(maxtab)
    if mintab(1,2) < maxtab(1,2)
        mintab = mintab(2:end,:);
    end
end
end

下面是一个使用示例：

%读取图像
im = imread('test.jpg');
%将图像转换为灰度图
im = rgb2gray(im);
%对图像进行二维经验模态分解
[imf,residue] = emd2(double(im));
%绘制分解后的各个IMF
figure;
for i=1:size(imf,3)
    subplot(3,3,i);
    imagesc(imf(:,:,i)); colormap(gray); axis off; axis equal;
end
%绘制剩余部分
figure;
imagesc(residue); colormap(gray); axis off; axis equal;

运行上述代码后，将会得到分解后的各个IMF和剩余部分的图像，如下图所示：