本文主要為了通過對代碼進行分析，幫助更進一步了解《Single Image Haze Removal Using Dark Channel Prior》的操作步驟。在看本篇文章前，需要對《Single Image Haze Removal Using Dark Channel Prior》有大致的了解。提供一個文章解析參考的連接：https://www.cnblogs.com/Imageshop/p/3281703.html。本文的代碼分析也是基於上述文章中的代碼資源進行分析的，本文主要在於結合代碼梳理出更加清晰的思路，並且對代碼進行了較為詳細的注釋。
一、主程序
首先在這給出主程序及其相關注釋
clear clc close all kenlRatio = 0.01; minAtomsLight = 240; %設定的大氣光的上限值 % image_name = 'test images\21.bmp'; image_name = 'D:\1.png'; %原圖像 img=imread(image_name); %讀取原始圖片文件中的數據 figure,imshow(uint8(img)), title('src'); %顯示原始圖像 並且顯示標題 sz=size(img); %返回一個行向量sz，該行向量的第一個元素時矩陣的行數，第二個元素是矩陣的列數。 w=sz(2); %讀取所提取矩陣的列數 h=sz(1); %讀取所提取矩陣的行數 dc = zeros(h,w); %返回一個h×w的0矩陣給dc，實際上這個矩陣的size=原輸入圖像size for y=1:h %行遍歷 for x=1:w %列遍歷 dc(y,x) = min(img(y,x,:)); %對原圖像的每一個像素求RGB3個通道中的最小亮度值,然后賦值給dc矩陣的y行x列 end end figure,imshow(uint8(dc)), title('Min(R,G,B)'); %顯示原圖像的灰度圖 krnlsz = floor(max([3, w*kenlRatio, h*kenlRatio])) %krnlsz為最小濾波窗口的邊長尺寸大小 %kenlRatio是一個預設的值，求尺寸時拿kenlRatio和灰度圖矩陣長和寬（即w和h）乘上這個預設值和3進行比較，選擇最大數作為窗口尺寸 dc2 = minfilt2(dc, [krnlsz,krnlsz]); %對得到的灰度圖進行最小值濾波獲得暗通道圖像 dc2(h,w)=0; %dc2矩陣的h行w列賦值0 figure,imshow(uint8(dc2)), title('After filter '); %顯示濾波后的圖像 t = 255 - dc2; %下面4行對應公式11，即求出傳輸特性t（x)的估計值 figure,imshow(uint8(t)),title('t'); %顯示傳輸特性t的圖像(原始的方法求t） t_d=double(t)/255; sum(sum(t_d))/(h*w) %t(x)表示的是x像素出的t的估計值，而此處將t_d的矩陣進行平均處理求出一個平均的t值 A = min([minAtomsLight, max(max(dc2))]) %根據論文提示的方法求A值，不過為了防止A過高，設置了一個上限值minAtomsLight %max(max(dc2))表示從暗通道圖中按照亮度的大小取前0.1%的像素。在這些位置中，在原始有霧圖像I中尋找對應的具有最高亮度的點的值，作為A值 J = zeros(h,w,3); %生成一個0矩陣J 大小為h×w 此處的3表示rgb3個通道 img_d = double(img); %原始圖像類型轉換成double J(:,:,1) = (img_d(:,:,1) - (1-t_d)*A)./t_d; %根據公式22求出圖像矩陣所有元素的R通道的值，下面兩行為G和B的值 J(:,:,2) = (img_d(:,:,2) - (1-t_d)*A)./t_d; J(:,:,3) = (img_d(:,:,3) - (1-t_d)*A)./t_d; figure,imshow(uint8(J)), title('J'); %顯示去霧圖像J %---------------------------------- %由於透射率圖過於粗糙，而通過導向濾波的方式來獲得較好的透射率圖，下面是導向濾波計算t的方式 r = krnlsz*4 eps = 10^-6; % filtered = guidedfilter_color(double(img)/255, t_d, r, eps); filtered = guidedfilter(double(rgb2gray(img))/255, t_d, r, eps);%導向濾波計算t t_d = filtered;%導向濾波計算的值賦值給t_d矩陣 figure,imshow(t_d,[]),title('filtered t'); %顯示導向濾波后求出的t值 J(:,:,1) = (img_d(:,:,1) - (1-t_d)*A)./t_d; %根據公式22求出圖像矩陣所有元素的R通道的值，下面兩行為G和B的值 J(:,:,2) = (img_d(:,:,2) - (1-t_d)*A)./t_d; J(:,:,3) = (img_d(:,:,3) - (1-t_d)*A)./t_d; imwrite(uint8(J),'D:\11.bmp');%將去霧后的圖像數據寫入到文件中 figure,imshow(uint8(J)), title('J_guild_filter');%顯示導向濾波后求出的去霧圖像J
二、主程序拆解分析
1.首先，根據文章我們可以得知：基於暗通道先驗的思想，我們可以先求出t的預估值：
在暗通道先驗的原理下，求t預估值的公式轉化為：
1.首先，根據文章我們可以得知：基於暗通道先驗的思想，我們可以先求出t的預估值：
在暗通道先驗的原理下，求t預估值的公式轉化為：
對於一張確定的圖，我們認為Ac大氣光是一個定值。
我們首先計算出原輸入圖像的灰度圖，即括號內的部分minc Ic(y)。大概思路是對原圖像的所有像素x求RGB通道中亮度最小的值，從而得到灰度圖。
我們首先計算出原輸入圖像的灰度圖，即括號內的部分minc Ic(y)。大概思路是對原圖像的所有像素x求RGB通道中亮度最小的值，從而得到灰度圖。
image_name = 'D:\1.png'; %原圖像 img=imread(image_name); %讀取原始圖片文件中的數據 figure,imshow(uint8(img)), title('src'); %顯示原始圖像 並且顯示標題 sz=size(img); %返回一個行向量sz，該行向量的第一個元素時矩陣的行數，第二個元素是矩陣的列數。 w=sz(2); %讀取所提取矩陣的列數 h=sz(1); %讀取所提取矩陣的行數 dc = zeros(h,w); %返回一個h×w的0矩陣給dc，實際上這個矩陣的size=原輸入圖像size for y=1:h %行遍歷 for x=1:w %列遍歷 dc(y,x) = min(img(y,x,:)); %對原圖像的每一個像素求RGB3個通道中的最小亮度值,然后賦值給dc矩陣的y行x列 end end figure,imshow(uint8(dc)), title('Min(R,G,B)'); %顯示原圖像的灰度圖
對於公式11減號后面內容的實現：對上面得到的灰度圖像進行最小值濾波，minfilt2是一個最小值濾波的函數。
kenlRatio = 0.01; krnlsz = floor(max([3, w*kenlRatio, h*kenlRatio])) %krnlsz為最小濾波窗口的邊長尺寸大小 %kenlRatio是一個預設的值，求尺寸時拿kenlRatio和灰度圖矩陣長和寬（即w和h）乘上這個預設值和3進行比較，選擇最大數作為窗口尺寸 dc2 = minfilt2(dc, [krnlsz,krnlsz]); %對得到的灰度圖進行最小值濾波獲得暗通道圖像
最后進行歸一化、均值處理實現公式11的效果：
dc2(h,w)=0; %dc2矩陣的h行w列賦值0 figure,imshow(uint8(dc2)), title('After filter '); %濾波后的圖像結果 t = 255 - dc2; %下面4行對應公式11，即求出傳輸特性t（x)的估計值 figure,imshow(uint8(t)),title('t'); %顯示傳輸特性t的圖像(原始的方法求t） t_d=double(t)/255; sum(sum(t_d))/(h*w) %t(x)表示的是x像素出的t的估計值，而此處將t_d的矩陣進行平均處理求出一個平均的t值
2.預估了傳輸特性t值后，現在根據文章求出大氣光A的估值。
minAtomsLight = 240; %設定的大氣光的上限值 A = min([minAtomsLight, max(max(dc2))]) %根據論文提示的方法求A值，不過為了防止A過高，設置了一個上限值minAtomsLight
3.根據最終的求去霧圖像J的公式來還原去霧圖像，根據的公式如下：
代碼中的平均估值t_d替換的是公式中的max(t(x),t0)部分。
J = zeros(h,w,3); %生成一個0矩陣J 大小為h×w 此處的3表示rgb3個通道 img_d = double(img); %原始圖像類型轉換成double J(:,:,1) = (img_d(:,:,1) - (1-t_d)*A)./t_d; %根據公式22求出圖像矩陣所有元素的R通道的值，下面兩行為G和B的值 J(:,:,2) = (img_d(:,:,2) - (1-t_d)*A)./t_d; J(:,:,3) = (img_d(:,:,3) - (1-t_d)*A)./t_d; figure,imshow(uint8(J)), title('J'); %顯示去霧圖像J
4.由於透射率圖過於粗糙，而通過導向濾波的方式來獲得較好的透射率圖，下面是導向濾波計算t的方式。
r = krnlsz*4 eps = 10^-6; filtered = guidedfilter(double(rgb2gray(img))/255, t_d, r, eps);%導向濾波計算t t_d = filtered;%導向濾波計算的值賦值給t_d矩陣 figure,imshow(t_d,[]),title('filtered t'); %顯示導向濾波后求出的t值
5.利用導向濾波計算的t_d求出最后的去霧圖像
J(:,:,1) = (img_d(:,:,1) - (1-t_d)*A)./t_d; %根據公式22求出圖像矩陣所有元素的R通道的值，下面兩行為G和B的值 J(:,:,2) = (img_d(:,:,2) - (1-t_d)*A)./t_d; J(:,:,3) = (img_d(:,:,3) - (1-t_d)*A)./t_d; imwrite(uint8(J),'D:\11.bmp');%將去霧后的圖像數據寫入到文件中 figure,imshow(uint8(J)), title('J_guild_filter');%顯示導向濾波后求出的去霧圖像J
三、程序運行結果
1.顯示原始輸入圖像
2.顯示原始圖像的灰度圖
3.顯示對原始圖像灰度圖進行最小濾波后的結果
4.求出利用暗通道先驗計算的t
5.根據求出的t和A計算出去霧圖像J
6.為了細化t，此處求導向濾波后的t
7.利用導向濾波后求出的去霧圖像
四、說明
這里值分析了主程序的思路，對於最小濾波和導向濾波的程序未分析，具體的代碼見文章開頭給出的鏈接。