matlab 图像去雾

  • 基于非物理模型的去雾方法:不考虑雾导致图像退化的原因,只通过实现对比度增强方法达到去雾目的;例如,子块部分重叠的局部直方图均衡化方法、对比度受限自适应直方图均衡化方法、Retinex增强方法等。
  • 基于物理模型的去雾方法:考虑雾导致图像退化的成因,进行数学建模,补偿退化过程造成的失真,恢复无雾图像。该类方法获得的无雾图像,视觉效果自然,一般无信息损失。

去雾方法

局部直方图均衡化

clc; clear;close all; 
I = imread('sweden_input.jpg');
g1 = GetLocalHisteq(I(:, :, 1));  % 对图像进行局部均衡化增强(分通道进行处理)
g2 = GetLocalHisteq(I(:, :, 2));
g3 = GetLocalHisteq(I(:, :, 3));
In = cat(3, g1, g2, g3);          % 将三个分量整合成三维图像
    figure;
    subplot(2, 2, 1); imshow(I); title('原图像', 'FontWeight', 'Bold');
    subplot(2, 2, 2); imshow(In); title('处理后的图像', 'FontWeight', 'Bold');
    Q = rgb2gray(I);               % 灰度化处理 
    W = rgb2gray(In);
    % 灰度化后进行直方图均衡化处理,并显示原直方图和,均衡化后的直方图  64:为输出图像的灰度级数
    subplot(2, 2, 3); imhist(Q, 64); title('原灰度直方图', 'FontWeight', 'Bold');  
    subplot(2, 2, 4); imhist(W, 64); title('处理后的灰度直方图', 'FontWeight', 'Bold');

GetLocalHisteq()函数

function g = GetLocalHisteq(I)
% 对灰度图像,进行局部直方图均衡化
% 输入参数:
%  I——图像矩阵
% 输出参数:
%  g——结果图像
% 调用库函数adapthisteq,执行局部均衡化增强

% ‘clipLimit’:范围是[0 1]内的标量,用于指定对比度增强的限制。较高的值产生较强的对比度。默认值是0.01
%‘Distribution’:字符串,用于指定图像小片所需的直方图形状: ‘uniform’——平坦的直方图(默认值),‘rayleigh’——钟形直方图,‘exponential’——曲线直方图
g = adapthisteq(I,'clipLimit',0.02,'Distribution','rayleigh');
end

有雾图像仿真 python 图像去雾matlab_直方图

全局直方图均衡化

clc; clear;close all; 
I = imread('sweden_input.jpg');       % 读取图像
R = I(:,:,1);  G = I(:,:,2);  B = I(:,:,3);   % 提取红绿蓝三个通道
%  对三个通道分别进行均衡化处理后整合为三维图像
M = histeq(R);  N = histeq(G);  L = histeq(B);   In = cat(3, M, N, L);
figure;
 subplot(2, 2, 1); imshow(I); title('原图像', 'FontWeight', 'Bold');
 subplot(2, 2, 2); imshow(In); title('处理后的图像', 'FontWeight', 'Bold');
    % 灰度化处理
    Q = rgb2gray(I);      W = rgb2gray(In);
    % 对图形进行均衡化处理,64为输出图像的灰度级数 
    subplot(2, 2, 3);  imhist(Q, 64); title('原灰度直方图', 'FontWeight', 'Bold');
    subplot(2, 2, 4);  imhist(W, 64); title('处理后的灰度直方图', 'FontWeight', 'Bold');

有雾图像仿真 python 图像去雾matlab_直方图_02

基于暗原色先验的图像增强

clc;clear ;close all;
Image=imread( 'scene1.jpg');                % 读取图像
imshow(Image), title('原始图像');           % 展示原图
[height ,width, c]=size(Image);             % size求取矩阵大小,返回其行列值,即hight、width
dark_I = zeros(height,width);               % zeros返回 hight x width 大小的零矩阵
for y=1:height                              
    for x=1:width
        dark_I(y,x) = min(Image(y,x,:));    %计算RGB取最小值后的图像dark_I
    end
end
kenlRatio = .03;                % 定义一个值为0.03的常量
krnlsz = floor(max([3, width*kenlRatio, height*kenlRatio]));    %  先计算max中三个值的最大值; 然后floor:四舍五入最大值
dark_I2 = minfilt2(dark_I, [krnlsz,krnlsz]);    %调用minfilt2计算最小值滤波
dark_I2(height,width)=0;        % 将最大行最大列置零
dark_channel=dark_I2;           % 将dark_I2赋值给dark_channel
hh=floor(width*height*0.001);   
bb=reshape(dark_channel,1,[]);  % reshape就是把指定的矩阵改变形状,但是元素个数不变,将dark_channel转为一维数组
bb=sort(bb,'descend');          % sort函数默认Mode为'ascend'为升序,sort(X,'descend')为降序排列
cc=find(bb>0,hh);       % find:bb中是否有hh个大于0的数
dd=bb(cc(hh));          % 数组bb[cc[hh]]的值
ee=dark_channel(dark_channel>dd);   
AA=(ee);
num=length(find(dark_channel>dd));
sum=0;
for k=1:num
    sum=sum+AA(k);
end
meandc=floor(sum/num);
minAtomsLight = 240;
A= min([minAtomsLight, meandc]);%计算大气光A
w0=0.9;  t0=0.3;%设置调节参数w0值
t=1-w0*(dark_channel/A);%计算透射率t
t=max(t,t0);
img_d = double(Image);
NewImage = zeros(height,width,3);
NewImage(:,:,1) = (img_d(:,:,1) - (1-t)*A)./t;%计算去雾后的R通道
NewImage(:,:,2) = (img_d(:,:,2) - (1-t)*A)./t;%计算去雾后的G通道
NewImage(:,:,3) = (img_d(:,:,3) - (1-t)*A)./t;%计算去雾后的B通道
subplot(121),imshow(uint8(Image)), title('去雾前图像I');%去雾图像
subplot(122),imshow(uint8(NewImage)), title('去雾图像J');%去雾图像
imwrite(uint8(NewImage),'scene2.bmp')           % 存储图像

minfilt2()最小值滤波

function Y = minfilt2(X,varargin)
%  MINFILT2    Two-dimensional min filter
%
%     Y = MINFILT2(X,[M N]) performs two-dimensional minimum
%     filtering on the image X using an M-by-N window. The result
%     Y contains the minimun value in the M-by-N neighborhood around
%     each pixel in the original image. 
%     This function uses the van Herk algorithm for min filters.
%
%     Y = MINFILT2(X,M) is the same as Y = MINFILT2(X,[M M])
%
%     Y = MINFILT2(X) uses a 3-by-3 neighborhood.
%
%     Y = MINFILT2(..., 'shape') returns a subsection of the 2D
%     filtering specified by 'shape' :
%        'full'  - Returns the full filtering result,
%        'same'  - (default) Returns the central filter area that is the
%                   same size as X,
%        'valid' - Returns only the area where no filter elements are outside
%                  the image.
%
%     See also : MAXFILT2, VANHERK
%

% Initialization
[S, shape] = parse_inputs(varargin{:});

% filtering
Y = vanherk(X,S(1),'min',shape);
Y = vanherk(Y,S(2),'min','col',shape);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [S, shape] = parse_inputs(varargin)
shape = 'same';
flag = [0 0]; % size shape

for i = 1 : nargin
   t = varargin{i};
   if strcmp(t,'full') && flag(2) == 0
      shape = 'full';
      flag(2) = 1;
   elseif strcmp(t,'same') && flag(2) == 0
      shape = 'same';
      flag(2) = 1;
   elseif strcmp(t,'valid') && flag(2) == 0
      shape = 'valid';
      flag(2) = 1;
   elseif flag(1) == 0
      S = t;
      flag(1) = 1;
   else
      error(['Too many / Unkown parameter : ' t ])
   end
end

if flag(1) == 0
   S = [3 3];
end
if length(S) == 1;
   S(2) = S(1);
end
if length(S) ~= 2
   error('Wrong window size parameter.')
end

有雾图像仿真 python 图像去雾matlab_matlab_03

基于暗原色先验的低照度图像增强

低照度图像,通常是指在光照较暗或者夜间所拍摄得到的图像。如果将低照度图像反转后,其图像表征及直方图表征与雾天图像具有很高的相似性,因此,可以采用图像去雾技术来实现低照度图像的增强。该算法的基本思想是:对输入的低照度图像进行反转,对反转图像进行去雾处理,然后对去雾结果反转得到低照度图像增强效果。

clear;clc;close all;
img=imread( 'Lyewan.bmp');
subplot(121),imshow(uint8(img)), title('原始低照度图像');
img(:,:,1)=255-img(:,:,1);          % 各通道取反
img(:,:,2)=255-img(:,:,2);
img(:,:,3)=255-img(:,:,3);
sz=size(img);                        % size求取矩阵大小,返回其行列值,即hight、width
w=sz(2);
h=sz(1);
%计算RGB取最小值后的图像dark_I
dark_I = zeros(h,w);
for y=1:h
    for x=1:w
        dark_I(y,x) = min(img(y,x,:));%计算每个像素点RGB中最小的值
    end
end

kenlRatio = .03;
krnlsz = floor(max([3, w*kenlRatio, h*kenlRatio]));
dark_I2 = minfilt2(dark_I, [krnlsz,krnlsz]);%最小值滤波
dark_I2(h,w)=0;
dark_I2=uint8(dark_I2);

dark_channel=double(dark_I2);       % 增加精度
aa=dark_channel;
hh=floor(w*h*0.001);
bb=reshape(aa,1,[]);%将a转换为行向量
bb=sort(bb,'descend');%将b按从大到小顺序排列。
cc=find(bb>0,hh);%找到前4个大于0(非NaN)的数的位置
dd=bb(cc(hh));%找到第4大的数
ee=aa(aa>dd);%较大的前3个数
[mm,nn]=find(aa>dd);%较大的前3个数对应下标
AA=[ee mm nn];
bw=zeros(h,w);
num=length(find(aa>dd));
sum=0;
for y=1:h
    for x=1:w
        for k=1:num
            if y==AA(k,2) && x==AA(k,3) && dark_channel(y,x)==AA(k,1)
                bw(y,x)=255;
                sum=sum+AA(k,1);
            end
        end
    end
end
meandc=floor(sum/num);
minAtomsLight = 240;
A= min([minAtomsLight, meandc]);%计算大气光A

w0=0.9;  
t=1-w0*(dark_channel/A);%计算透射率t
t0=0.1;
t=max(t,t0);
img_d = double(img);
J = zeros(h,w,3);
J(:,:,1) = (img_d(:,:,1) - (1-t)*A)./t;%计算去雾后的R通道
J(:,:,2) = (img_d(:,:,2) - (1-t)*A)./t;%计算去雾后的G通道
J(:,:,3) = (img_d(:,:,3) - (1-t)*A)./t;%计算去雾后的B通道
J=uint8(J);         % 把大于255的数全部强制置为255,而小于255的部分则保持原样不变。
J(:,:,1)=255-J(:,:,1);
J(:,:,2)=255-J(:,:,2);
J(:,:,3)=255-J(:,:,3);
subplot(122),imshow(J), title('基于暗原色先验的低照度图像增强J');%去雾图像
imwrite(J,'Lyewan.bmp');

vanherk()函数

function Y = vanherk(X,N,TYPE,varargin)
%  VANHERK    Fast max/min 1D filter
%
%    Y = VANHERK(X,N,TYPE) performs the 1D max/min filtering of the row
%    vector X using a N-length filter.
%    The filtering type is defined by TYPE = 'max' or 'min'. This function
%    uses the van Herk algorithm for min/max filters that demands only 3
%    min/max calculations per element, independently of the filter size.
%
%    If X is a 2D matrix, each row will be filtered separately.
%    
%    Y = VANHERK(...,'col') performs the filtering on the columns of X.
%    
%    Y = VANHERK(...,'shape') returns the subset of the filtering specified
%    by 'shape' :
%        'full'  - Returns the full filtering result,
%        'same'  - (default) Returns the central filter area that is the
%                   same size as X,
%        'valid' - Returns only the area where no filter elements are outside
%                  the image.
%
%    X can be uint8 or double. If X is uint8 the processing is quite faster, so
%    dont't use X as double, unless it is really necessary.
%

% Initialization
[direc, shape] = parse_inputs(varargin{:});
if strcmp(direc,'col')
   X = X';
end
if strcmp(TYPE,'max')
   maxfilt = 1;
elseif strcmp(TYPE,'min')
   maxfilt = 0;
else
   error([ 'TYPE must be ' char(39) 'max' char(39) ' or ' char(39) 'min' char(39) '.'])
end

% Correcting X size
fixsize = 0;
addel = 0;
if mod(size(X,2),N) ~= 0
   fixsize = 1;
   addel = N-mod(size(X,2),N);
   if maxfilt
      f = [ X zeros(size(X,1), addel) ];
   else
      f = [X repmat(X(:,end),1,addel)];
   end
else
   f = X;
end
lf = size(f,2);
lx = size(X,2);
clear X

% Declaring aux. mat.
g = f;
h = g;

% Filling g & h (aux. mat.)
ig = 1:N:size(f,2);
ih = ig + N - 1;

g(:,ig) = f(:,ig);
h(:,ih) = f(:,ih);

if maxfilt
   for i = 2 : N
      igold = ig;
      ihold = ih;
      
      ig = ig + 1;
      ih = ih - 1;
      
      g(:,ig) = max(f(:,ig),g(:,igold));
      h(:,ih) = max(f(:,ih),h(:,ihold));
   end
else
   for i = 2 : N
      igold = ig;
      ihold = ih;
      
      ig = ig + 1;
      ih = ih - 1;
      
      g(:,ig) = min(f(:,ig),g(:,igold));
      h(:,ih) = min(f(:,ih),h(:,ihold));
   end
end
clear f

% Comparing g & h
if strcmp(shape,'full')
   ig = (N : 1 : lf );
   ih = ( 1 : 1 : lf-N+1 );
   if fixsize
      if maxfilt
         Y = [ g(:,1:N-1)  max(g(:,ig), h(:,ih))  h(:,end-N+2:end-addel) ];
      else
         Y = [ g(:,1:N-1)  min(g(:,ig), h(:,ih))  h(:,end-N+2:end-addel) ];
      end
   else
      if maxfilt
         Y = [ g(:,1:N-1)  max(g(:,ig), h(:,ih))  h(:,end-N+2:end) ];
      else
         Y = [ g(:,1:N-1)  min(g(:,ig), h(:,ih))  h(:,end-N+2:end) ];
      end
   end
   
elseif strcmp(shape,'same')
   if fixsize
      if addel > (N-1)/2
         disp('hoi')
         ig = ( N : 1 : lf - addel + floor((N-1)/2) );
         ih = ( 1 : 1 : lf-N+1 - addel + floor((N-1)/2));
         if maxfilt
            Y = [ g(:,1+ceil((N-1)/2):N-1)  max(g(:,ig), h(:,ih)) ];
         else
            Y = [ g(:,1+ceil((N-1)/2):N-1)  min(g(:,ig), h(:,ih)) ];
         end
      else   
         ig = ( N : 1 : lf );
         ih = ( 1 : 1 : lf-N+1 );
         if maxfilt
            Y = [ g(:,1+ceil((N-1)/2):N-1)  max(g(:,ig), h(:,ih))  h(:,lf-N+2:lf-N+1+floor((N-1)/2)-addel) ];
         else
            Y = [ g(:,1+ceil((N-1)/2):N-1)  min(g(:,ig), h(:,ih))  h(:,lf-N+2:lf-N+1+floor((N-1)/2)-addel) ];
         end
      end            
   else % not fixsize (addel=0, lf=lx) 
      ig = ( N : 1 : lx );
      ih = ( 1 : 1 : lx-N+1 );
      if maxfilt
         Y = [  g(:,N-ceil((N-1)/2):N-1) max( g(:,ig), h(:,ih) )  h(:,lx-N+2:lx-N+1+floor((N-1)/2)) ];
      else
         Y = [  g(:,N-ceil((N-1)/2):N-1) min( g(:,ig), h(:,ih) )  h(:,lx-N+2:lx-N+1+floor((N-1)/2)) ];
      end
   end      
   
elseif strcmp(shape,'valid')
   ig = ( N : 1 : lx);
   ih = ( 1 : 1: lx-N+1);
   if maxfilt
      Y = ( max( g(:,ig), h(:,ih) ) );
   else
      Y = ( min( g(:,ig), h(:,ih) ) );
   end
end

if strcmp(direc,'col')
   Y = Y';
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [direc, shape] = parse_inputs(varargin)
direc = 'lin';
shape = 'same';
flag = [0 0]; % [dir shape]

for i = 1 : nargin
   t = varargin{i};
   if strcmp(t,'col') && flag(1) == 0
      direc = 'col';
      flag(1) = 1;
   elseif strcmp(t,'full') && flag(2) == 0
      shape = 'full';
      flag(2) = 1;
   elseif strcmp(t,'same') && flag(2) == 0
      shape = 'same';
      flag(2) = 1;
   elseif strcmp(t,'valid') && flag(2) == 0
      shape = 'valid';
      flag(2) = 1;
   else
      error(['Too many / Unkown parameter : ' t ])
   end
end

有雾图像仿真 python 图像去雾matlab_直方图_04