MATLAB---非线性复扩散滤波

286 篇文章 33 订阅
236 篇文章 15 订阅

在这里插入图片描述

%

非线性复扩散滤波
clc,clear,close all  % 清理命令区、清理工作区、关闭显示图形
warning off       % 消除警告
feature jit off      % 加速代码运行
tic
[filename ,pathname]=...
    uigetfile({'*.bmp';'*.tif';'*.jpg';},'选择图片'); %选择图片路径
str=[pathname filename]; % 合成路径+文件名
im = imread(str);        % 读图
im = imnoise(im,'gaussian',0,1e-3); % 原图像 + 白噪声

%非线性复扩散滤波参数设置
TMAX = 0.80;   % 扩散时间 
[im_e, nIter, dTT] = NCD_filter(im, TMAX);  % 非线性复扩散滤波


figure,
subplot(121),imshow(im);title('原始图像')
colormap(jet)  % 颜色
shading interp % 消隐
subplot(122),imshow(im_e,[]);title('非线性复扩散滤波图像')
colormap(jet)  % 颜色
shading interp % 消隐
toc
function [imgOutput, nIter, dtt] = NCD_filter(imgInput, Tmax)
% 非线性复扩散滤波
% 函数输入:
%        imgInput - 含噪声的图像
%        Tmax  - 扩散时间
% 函数输出:
%        imgOutput - f非线性复扩散滤波图像
%        nIter  - 滤波迭代处理次数
%        dtt    - 每一次迭代的时间步长
% 设定默认的迭代扩散时间
if nargin < 2 ;   % 输入个数小于2
    Tmax  = 3.0 ; % Tmax默认赋值
end 

% 初始化操作
theta     = pi/30;           % 初始化
j         = sqrt(-1);        % 初始化
e_jxtheta = exp(j*theta);    % 初始化
kappamin  =   2.0;           % 初始化
kappamax  =  28.0;           % 初始化

% 高斯滤波器
wsize     = 3;  % 窗口大小 3 x 3
wsigma    = 10; % 方差
filter_gaussian    = fspecial('gaussian', wsize, wsigma);  % 滤波掩模
% fspecial('gaussian', 3, 0.001)
% ans =
%      0     0     0
%      0     1     0
%      0     0     0
% fspecial('gaussian', 5, 0.001)
% ans =
%      0     0     0     0     0
%      0     0     0     0     0
%      0     0     1     0     0
%      0     0     0     0     0
%      0     0     0     0     0
% 扩散滤波器系数
wsigma2    = 0.5; % 方差
filter_gaussian2  = fspecial('gaussian', wsize, wsigma2);  % 滤波掩模

lapKernel      = [0,1,0; 1,-4,1; 0,1,0]; % laplacian kernel
gradKernel     = [-1/2,0,1/2];          % 梯度kernel
[xSize, ySize] = size(imgInput);         % 图像维数

Border = 2; % 图像边界2个像素点之间不进行梯度计算
indexX = 1+Border:xSize+Border; % (1+Border):(xSize+Border)
indexY = 1+Border:ySize+Border; % (1+Border):(xSize+Border)

if ~isa(imgInput,'double')
    yNCDF = double(imgInput); % 图像数据类型转换
end

t_iter = 0;  % 迭代时间累加和
nIter = 0;  % 迭代次数

while (Tmax - t_iter) > eps % 扩散时间
    nIter = nIter + 1;
    ryNCDF = real(yNCDF);  % 实部
    iyNCDF = imag(yNCDF);  % 虚部
    
    % 滤波,见表达式(10.32)
    localAvg = filter2(filter_gaussian, ryNCDF, 'same');
    minlocalAvg = min(localAvg(:));   % 最小值
    
    % 自适应系数 k,见表达式(10.31)
    k_mod = kappamax + (kappamin-kappamax)* ...
        (localAvg - minlocalAvg) / (max(localAvg(:)) - minlocalAvg);
    % 非线性复扩散函数
    coefDif = e_jxtheta ./ (1 + ( (iyNCDF/theta) ./ k_mod ).^2);
    coefDif = filter2(filter_gaussian, coefDif, 'same');

    % 计算 laplacian and gradient
    % lap(yNCDF)
    yAux  = zeros(xSize+2*Border, ySize+2*Border);
    yAux(indexX, indexY) = yNCDF;
    del2Y = conv2(yAux, lapKernel, 'same');
    del2Y = del2Y(indexX, indexY);
    
    % grad(yNCDF)
    dAux  = conv2(yAux, gradKernel, 'same');  % 卷积
    dYx   = dAux(indexX, indexY);
    dAux  = conv2(yAux, gradKernel', 'same'); % 卷积
    dYy   = dAux(indexX, indexY);
    
    % grad(coefDif)
    dDx   = conv2(coefDif, gradKernel, 'same');  % 卷积
    dDy   = conv2(coefDif, gradKernel', 'same'); % 卷积
    
    dIdt  = coefDif.*del2Y + dDx.*dYx + dDy.*dYy;
    
    % 计算自适应时间步长
    dt = (1/4)*( 0.25 + 0.75*exp(- max(max( abs(real(dIdt)) ./ ryNCDF ))) );
    
    dtt(nIter) = dt;
    
    % 约束每一步的最大步长
    if t_iter + dt > Tmax ;
        dt = Tmax - t_iter ; 
    end
    % 更新已经处理已经扩散的时间消耗
    t_iter = t_iter + dt;
    
    % 更新 yNCDF
    yNCDF = yNCDF + dt.*dIdt;
    
end % 结束,对应while (Tmax - t_iter) > eps

imgOutput = real(yNCDF);  % 实部,图像输出

end

  • 0
    点赞
  • 1
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值