matlab 分解灰度图像emd代码

以下是用 MATLAB 实现的分解灰度图像 EMD（经验模态分解）的代码：

function [imfs, residual] = emd_gray(im, maxiter)
% EMD_GRAY Empirical Mode Decomposition for gray-scale images
%   [IMFS, RESIDUAL] = EMD_GRAY(IM, MAXITER) performs Empirical Mode 
%   Decomposition (EMD) on the input gray-scale image IM. The decomposition
%   generates a set of Intrinsic Mode Functions (IMFs) and a residual image.
%   The maximum number of iterations for EMD can be specified as MAXITER.
%
%   Input arguments:
%   ------------------
%   IM:      Input gray-scale image (2D matrix)
%   MAXITER: Maximum number of iterations for EMD (default: 100)
%
%   Output arguments:
%   ------------------
%   IMFS:     Set of Intrinsic Mode Functions (IMFs) (3D matrix)
%   RESIDUAL: Residual image (2D matrix)

% Set default value for MAXITER
if nargin < 2
    maxiter = 100;
end

% Pre-allocate memory for IMFs and residual
imfs = zeros(size(im,1), size(im,2), maxiter);
residual = im;

% Perform EMD
for i = 1:maxiter
    % Compute the first IMF
    imf = imemd(residual);
    % Save the first IMF
    imfs(:,:,i) = imf;
    % Subtract the first IMF from the original image
    residual = residual - imf;
    % Stop if the residual becomes too small
    if sum(residual(:).^2) < eps
        break;
    end
end

% Remove unused entries from IMFs
imfs(:,:,i+1:end) = [];

end

function imf = imemd(im)
% IMEMD Intrinsic Mode Function (IMF) computation using EMD
%   IMF = IMEMD(IM) computes the first Intrinsic Mode Function (IMF) of the
%   input gray-scale image IM using Empirical Mode Decomposition (EMD).
%
%   Input arguments:
%   ------------------
%   IM: Input gray-scale image (2D matrix)
%
%   Output arguments:
%   ------------------
%   IMF: First Intrinsic Mode Function (IMF) (2D matrix)

% Set stopping criterion
eps = 1e-5;

% Initialize variables
imf = im;
h = imf;

% Perform EMD
while true
    % Find upper and lower envelopes
    u = h;
    l = h;
    for i = 1:10
        % Find local maxima and minima
        maxlocs = (diff(sign(diff(u))) < 0) + 1;
        minlocs = (diff(sign(diff(l))) > 0) + 1;
        % Interpolate to get upper and lower envelopes
        if length(maxlocs) > 2
            u = interp1(maxlocs, h(maxlocs), 1:length(h), 'spline');
        end
        if length(minlocs) > 2
            l = interp1(minlocs, h(minlocs), 1:length(h), 'spline');
        end
    end
    % Compute mean of upper and lower envelopes
    m = (u + l) / 2;
    % Extract IMF
    imf = imf - m;
    % Compute residue
    r = imf - h;
    % Check if residue is small enough
    if sum(r(:).^2) < eps
        break;
    end
    % Update variables for next iteration
    h = m;
end

end

该代码将输入的灰度图像分解成一组内禀模态函数（IMFs）和一个残差图像。该算法使用经验模态分解（EMD）方法进行分解。由于 EMD 算法的非确定性性质，同一幅图像可能会得到不同的 IMFs。因此，使用不同的参数和初始条件可能会导致不同的结果。