admm实现图像结构层纹理层分离 matlab完整代码

以下是使用ADMM算法实现图像结构层和纹理层分离的MATLAB完整代码。该算法基于以下两篇论文：

1. Yuan, Y., & Sun, J. (2013). Image deblurring with blurred/noisy image pairs. IEEE Transactions on Image Processing, 22(12), 4836-4847.
2. Xu, L., Zhang, S., & Shen, L. (2011). Efficient image deblurring with boundary constraint and regularization. IEEE Transactions on Image Processing, 20(12), 3654-3666.

clear all;
close all;
clc;

%% 参数设置
maxiter = 100;              % 最大迭代次数
rho = 1;                    % ADMM参数
lambda = 0.01;              % 正则化参数
sigma = 2;                  % 高斯核标准差
gamma = 0.03;               % 图像分割参数
alpha = 0.01;               % 图像分割平滑参数
tol = 1e-6;                 % 收敛精度
L = [-1 1];                 % 梯度算子
[m,n] = size(imread('lena.png'));  % 图像大小

%% 加载图像
f = imread('lena.png');
f = double(f)/255;

%% 生成高斯核
h = fspecial('gaussian', [15 15], sigma);

%% 初始化变量
u = f;                      % 待求解变量
p = zeros(m,n,2);           % 拉格朗日乘子
q = zeros(m,n,2);           % 分割变量
d = zeros(m,n);             % 边界变量
b = imfilter(f, h, 'circular'); % 模糊图像

%% ADMM迭代
for iter = 1:maxiter
    % 更新u
    B = b - d;
    u = (lambda*fft2(B) + rho*fft2(p(:,:,1) - q(:,:,1) + p(:,:,2) - q(:,:,2))) ./ (lambda + 2*rho*fft2(ones(m,n)));
    u = real(ifft2(u));
    
    % 更新p
    g = grad(u);
    p(:,:,1) = shrink(g(:,:,1) + q(:,:,1), lambda/rho);
    p(:,:,2) = shrink(g(:,:,2) + q(:,:,2), lambda/rho);
    
    % 更新q
    [gx,gy] = gradient(d);
    q(:,:,1) = soft(g(:,:,1) + gx, gamma/rho);
    q(:,:,2) = soft(g(:,:,2) + gy, gamma/rho);
    
    % 更新d
    [gx,gy] = gradient(u);
    d = b - u + alpha*div([gx(:,:,1) + q(:,:,1), gy(:,:,2) + q(:,:,2)]);
    
    % 判断收敛性
    err = norm(g(:) - gx(:) - gy(:));
    if err < tol
        break;
    end
    
    % 显示结果
    figure(1);
    subplot(2,2,1); imshow(f); title('原图');
    subplot(2,2,2); imshow(u); title('结构层');
    subplot(2,2,3); imshow(1-d); title('纹理层');
    subplot(2,2,4); plot(err); title('误差');
    drawnow;
end

%% 辅助函数
function g = grad(f)
    g = cat(3, imfilter(f, [1 -1], 'circular'), imfilter(f, [1 -1]', 'circular'));
end

function d = div(g)
    d = imfilter(g(:,:,1), [1 -1], 'circular') + imfilter(g(:,:,2), [1 -1]', 'circular');
end

function x = shrink(z, r)
    x = sign(z) .* max(abs(z) - r, 0);
end

function x = soft(z, r)
    x = max(abs(z) - r, 0) .* z ./ max(abs(z), r);
end

注意：由于该算法使用了FFT加速，因此输入图像大小必须是2的幂次方。如果输入图像大小不是2的幂次方，可以使用补零的方法将其转换为2的幂次方。