一種可實時處理 O(1)復雜度圖像去霧算法的實現。

本文轉載自查看原文 2013-11-06 12:35 18178 [05] 圖像去霧技術/ 去霧/ 均值模糊

　　在我博文的一系列的文章，有不少算法都於去霧有關，比如限制對比度自適應直方圖均衡化算法原理、實現及效果、局部自適應自動色階/對比度算法在圖像增強上的應用這兩個增強算法都有一定的去霧能力，而最直接的就是《Single Image Haze Removal Using Dark Channel Prior》一文中圖像去霧算法的原理、實現、效果及其他一文，描述了暗通道去霧這一state-of-the-art algorithms的過程和實現，雖幾經優化，對於常用的視頻1024*768大小的圖片，算法處理部分還是需要70MS的時間（I3 筆記本CPU），因此，這一算法用於實時要求時還有一定的難度，並且優化后的算法基本無法並行，而可並行的算法重復計算大，由於不熟悉GPU方面的理念，不曉得使用不優化的算法靠GPU是否能有多大速度的提升。

為此，我一直在找尋相關的論文，這種找尋的蹤跡一般就是看到一篇好論文--》看其參考文獻--》再看參考文獻的參考文獻，這樣循環下去。然后有某種機會或巧合，又看到一篇好論文，重復前面的過程，你就會發現很多交集，慢慢的就會有一些好運向你招手。話說我原本只看英文的文獻，所以一直忽略了國內的文章，前幾日，一個QQ朋友推薦了一篇清華大學的論文，下載后稍微看了下，覺得其描述的結果還是比較吸引人的，於是就實現了下，實時的效果應該說很不錯，這里就簡單的介紹並推薦給大家。

原始論文下載：基於單幅圖像的快速去霧.pdf ，作者劉倩，陳茂銀，周東華，感謝他們（她們）。

算法原理沒有什么復雜的地方，其實說原理，還不如說經驗或實驗，因為論文中可以用理論來推導的公式確實不多。不過這也沒關系，有用的東西就應該拿來用。

算法的執行流程直接貼用原圖的文字來說明吧：

先挑點小瑕疵，比如步驟5中的L有個下標，而其他涉及到L的地方去沒有，這叫前后不一致。論文中也還有一些其他的地方有些小錯誤，所以看啊，即使是清華的文章，在編輯審核這一塊還是相當的不嚴謹。

我們先來看看算法，1、2、3、4步都沒有什么說的，第五步是求大氣透射率的過程，這里ρ是一個用來調節的參數，當ρ值越大時，結果圖像整體越暗，去霧的效果更明顯，ρ較小時，圖像偏白，有明顯的霧氣。第六步求全局大氣光A值，用了很簡單的方式，即求原始圖像的RGB所有像素分量的最大值（這個估計99%都為255了）何暗通道的最大值的平均值，並且注意到RGB三個通道用的A值都為同一個數字。這個和我在何凱明的文章的分析也有相似的地方。第七步用了標准的去霧模型來求結果值。

在來看看算法的效率問題，從算法的初步分析來看，算法的效率取決於第3步和第7步，第三步中，使用了均值模糊，目前已經存在大量的O(1)均值模糊算法，可是O(1)只能表示算法的執行速度和參數的大小無關系，並不表示算法就很快。比如基於積分圖的模糊算法是廣為認知的O(1)算法，但是他也存在很多問題，最嚴重的就是數據的溢出，當圖像較大和偏白時，對圖像積分圖的累加和存在超出int.Maxvalue所能表達的范圍的問題，解決辦法就是積分圖內的數據全部使用long類型表示，這將導致程序多占用Width*Height*4字節的大小的內存，且在32位系統還流行的情況進一步降低程序的速度（32位系統64位整數的計算速度要比32位整數慢)。積分圖的另外一個問題就是計算積分圖的過程難以並行化，因為一個像素的積分值是依賴於其前面一系列像素的相關結果值的。另外一種優化方式就是先計算行方向的平均值，然后再計算列方向的值。這種方式在同一行（列）內，算法依舊必行順序執行，這也是因為前后影響的原因。但是不同行（列）之間的計算是沒有任何關系，因此非常適合GPU這種可大規模並行計算的場合，但不適於CPU這種重量級的線程並發（反而會慢）。這種算法如果為了精度會需要一個和原圖一樣大小，占用字節Width*Height*4字節大小的的中轉區用來保存中間計算的結果。在彩色圖像高速模糊之懶惰算法一文中，我采用了另外一種處理方法，利用列直方圖相關的技術，只需對每個循環的起始位置處的像素做特殊處理，其他位置的利用簡單的一加一簡即可獲得累加和，從而快速的實現模糊，我實際的編碼表明，這種方式比其他的方式都要快。但是有一個缺點，不適合於並行計算，不過在CPU上這個很有優勢。

再觀察下第七步。第七步存在兩個需要優化的地方，第一，存在除法；第二，有浮點運算。如果直接編碼必然會帶來性能損失，但是，觀察下在第七步的公式中，只有兩個自變量，H(X)和L(X),並且自變量的取值都為[0,255]之間的整數，因此，如果事先建議一個查找表，由於這個查找表的計算量只有 256*256次，要遠遠的小於直接計算的次數，必然能提高程序的速度。256這個數字還有個好處，就是可以用移位來輔助計算查表表的下標，部分參考代碼如下所示：

unsigned char * Table = (unsigned char  *) malloc (256*256*sizeof(unsigned char));
for (Y = 0; Y < 256; Y++)
{
    Index=Y<<8;
    for (X = 0; X < 256; X++)
    {
        Value = (Y-X) /(1-X*InvA);
        if (Value > 255)
            Value = 255;
        else if (Value < 0)
            Value = 0;
        Table[Index++]= Value;
    }
}

　　其中InvA = 1 / A;A為全局大氣光值。 Value需為double類型的變量。

當我們需要計算F(x)時，查表的方式為 F(X)=Table[(H[X]<<8 )+ L[X]];

實際的效果表明，這樣的方式對於1024*768的圖，可以提速10ms。

那么對於其他步驟也有很多優化的注意事項，比如計算M(X)中所有元素的平均值M_av這一塊，完全沒有必要在開一個循環，而是可以在進行步驟3的時候同步進行，大家知道，循環楚了要計算循環體內部的東西外，還要有個循環計數器的更新的，何必浪費這個時間呢。

for (Y = 0, DarkPt = DarkChannel; Y < Height; Y++)
{
    ImgPt = Src + Y * Stride;
    for (X = 0; X < Width; X++)
    {
        Min = *ImgPt;
        if (Min > *(ImgPt + 1)) Min = *(ImgPt + 1);
        if (Min > *(ImgPt + 2)) Min = *(ImgPt + 2);
        *DarkPt = Min;                                        //    三通道的最小值
        Sum    +=    Min;                                        //  累積以方便后面求平均值
        ImgPt += 3;
        DarkPt++;
    }
}
Mean =(double)Sum /(Width*Height*255);

還有比如第6步中，分別求原圖RGB三像素最大值以及安通道中的最大值的過程，傳統的過程如下代碼：

int Max1 =0 ;
for (Y = 0; Y < Height; Y++)
{
    ImgPt = Src + Y * Stride;
    for (X = 0; X < Width; X++)
    {
        if (Max1 < *(ImgPt)) Max1 = *(ImgPt);
        if (Max1 < *(ImgPt + 1)) Max1 = *(ImgPt + 1);
        if (Max1 < *(ImgPt + 2)) Max1 = *(ImgPt + 2);
        ImgPt += 3;
    }
}

　　特別是對於求原圖的最大值，實際上很多情況下這個值都為255，因此如果Max1變量已經是255，則循環完全沒有必要進行下去了，因此，如果改為下述代碼，必然可以減少計算量：

for (Y = 0; Y < Height; Y++)
{
    ImgPt = Src + Y * Stride;
    for (X = 0; X < Width; X++)
    {
        if (Max1 < *(ImgPt)) Max1 = *(ImgPt);
        if (Max1 < *(ImgPt + 1)) Max1 = *(ImgPt + 1);
        if (Max1 < *(ImgPt + 2)) Max1 = *(ImgPt + 2);
        ImgPt += 3;
    }
    if (Max1==255) break;
}

注意，這個break語句必須放在Y循環中，如果放在X循環中，雖然提前退出循環的可能性會增加，但是判斷的工作量帶來的損失更多。
綜合上述優化，我用C++寫了個DLL，對於1024*768的圖像在我的I3的機器上平均能達到18ms每副圖像的計算速度，相當於56fps,只占用了單核的資源，考慮解碼、顯示等等其他過程所占用的時間，應該是能夠靠CPU實現20fps的實時速度的。

在內存占用上，約需要> 3*Width*Height+256*256字節的空間（不包括圖像本身的），如果用在連續的視頻處理上，這部分內存就不需要頻繁的分配和釋放，可能也對速度的保證有好處。

程序下載地址： http://files.cnblogs.com/Imageshop/FastHazeRemovalTest.rar

　參數的選取：

　為了獲得好的效果，該算法需要選擇恰當的參數，我們為此做了一些測試。對於半徑參數，我的個人建議是取值不要小於50或圖像寬度和高度最大值的的1/20，比如對下面的圖像（原圖大小不是下圖中的大小，這里是為了方便瀏覽縮小顯示的），ρ取1.28時（對於上圖中的去霧程序選擇為128），半徑取不同參數時的效果：