本篇博文來自博主Imageshop，打賞或想要查閱更多內(nèi)容可以移步至Imageshop。

轉(zhuǎn)載自：https://www.cnblogs.com/Imageshop/p/6166394.html???侵刪

　暗角圖像是一種在現(xiàn)實(shí)中較為常見的圖像，其主要特征就是在圖像四個(gè)角有較為顯著的亮度下降，比如下面兩幅圖。根據(jù)其形成的成因，主要有3種：natural vignetting, pixel vignetting, 以及mechanic?vignetting，當(dāng)然，不管他的成因如何，如果能夠把暗角消除或者局部消除，則就有很好的工程意義。

???　　

? ? ? 這方面的資料和論文也不是很多，我最早是在2014年看到Y(jié). Zheng等人的論文《Single image vignetting correction》以及同樣有他們撰寫的論文《Single image vignetting correction using radial gradient symmetry》有講這方面的算法，不過其實(shí)現(xiàn)的復(fù)雜度較高，即使能編程實(shí)現(xiàn)，速度估計(jì)也很慢，其實(shí)用性就不高了。

? ? ? 前不久，偶爾的機(jī)會(huì)看到一篇名為《Single-Image Vignetting Correction by Constrained Minimization of log-Intensity Entropy》的論文，并且在github上找到了相關(guān)的一些參考代碼，雖然那個(gè)代碼寫的實(shí)在是惡心加無聊，但是對(duì)于我來說這并不重要，只要稍有參考，在結(jié)合論文那自己來實(shí)現(xiàn)就不是難事了。

? ? ?論文里的算法核心其實(shí)說起來也沒啥難的，我就我的理解來簡單的描述下：

? ? ?第一：去暗角可以說是陰影校正的一種特例，而將整副圖像的熵最小化也被證明為進(jìn)行陰影校正的一種有效方法，但是普通的熵在優(yōu)化過程中會(huì)優(yōu)化到局部最優(yōu)的。因此論文中提出了一種對(duì)數(shù)熵的概念（Log-Intensity Entropy），論文中用數(shù)據(jù)做了說明，假設(shè)一副普通正常的圖像其直方圖是單峰分布，那么如果這幅圖像有暗角，其直方圖必然會(huì)存在另外一個(gè)低明度的分布，如下圖所示：

　　我們校正暗角的過程就是使低明度的分布向原來的正常明度靠近，由上圖第一行的數(shù)據(jù)可以看到，普通的熵計(jì)算直到兩個(gè)直方圖有部分重疊的時(shí)候熵才會(huì)下降，之前熵一直都是增加的，而對(duì)數(shù)熵則在沒有重疊前至少是保持不增的，因此能夠更好的獲取全局最優(yōu)解。

? ? ?那么論文提出的對(duì)數(shù)熵的計(jì)算公式為：

? ? ?首先先將亮度進(jìn)行對(duì)數(shù)映射，映射公式為： ? ??

? ? ? ? ? ? ? ? ? ? ??

　　也就是將[0,255]內(nèi)的像素值映射到[0, N-1]內(nèi)，但不是線性映射，而是基于對(duì)數(shù)關(guān)系的映射，通常N就是取256，這樣映射后的像素范圍還是[0,255]，但是注意這里的i(L)已經(jīng)是浮點(diǎn)數(shù)了。我們繪制出N等于256時(shí)上式的曲線：

? ? ? ? ? ? ? ? ? ??

　　可見，這種操作實(shí)際上把圖像整體調(diào)亮了。

　　由于映射后的色階已經(jīng)是浮點(diǎn)數(shù)了，因此，直方圖信息的統(tǒng)計(jì)就必須換種方式了，論文給出的公式為：

? ? ? ? ? ? ?

? ? ? 公式很復(fù)雜， 其實(shí)就是有點(diǎn)類似線性插值那種意思，不認(rèn)識(shí)了那兩個(gè)數(shù)學(xué)符號(hào)了，就是向上取整和向下取整。

? ? ? 這樣的對(duì)數(shù)熵直方圖信息會(huì)由于巨大的色階調(diào)整，導(dǎo)致很多色階是沒有直方圖信息的，一般可以對(duì)這個(gè)直方圖信息進(jìn)行下高斯平滑的，得到新的直方圖。

　　最后圖像的對(duì)數(shù)熵，計(jì)算方法如下：

? ? ? ? ? ? ? ? ? ? ? ? ? ? ??

　　其中：? ??

? ? ? ? ? ? ? ? ? ? ? ? ? ?

? ? ?第二：論文根據(jù)資料《Radiometric alignment and vignetting calibration》提出了一個(gè)暗角圖像亮度下降的關(guān)系式，而我去看《Radiometric alignment and vignetting calibration》這篇論文時(shí)，他的公式又是從《Vignette and exposure calibration and compensation》中獲取的，所以這個(gè)論文的作者寫得文章還不夠嚴(yán)謹(jǐn)。這個(gè)公式是一個(gè)擁有五個(gè)未知參數(shù)的算式，如下所示：

? ? ? ? ? ? ? ? ???

　　其中：

? ? ? ? ? ? ??

? ? ?其中，x和y是圖像每一點(diǎn)的坐標(biāo)，而則表示暗角的中心位置，他們和a、b、c均為未知量。

　 ?我們可以看到，當(dāng)r=0時(shí)，校正系數(shù)為1，即無需校正。當(dāng)r=1時(shí)，校正系數(shù)為1+a+b+c。

? ? ?那么經(jīng)過暗角校正后的圖像就為：

? ? ? ? ? ? ?

　　按照我們的常識(shí)，暗角圖像從暗角的中心點(diǎn)想四周應(yīng)該是逐漸變暗的，根據(jù)上式函數(shù)g應(yīng)該是隨著r單調(diào)遞增的（因?yàn)槲覀兪切Ｕ到菆D像，所以越靠近邊緣上式的乘法中g(shù)值也就應(yīng)該越大)，因此函數(shù)g的一階導(dǎo)數(shù)應(yīng)該大于0，即：

? ? ? ? ? ? ?

? ? ?同時(shí)，我們注意到參數(shù)r的范圍很明顯應(yīng)該在[0,1]之間，這樣上式則可以轉(zhuǎn)換為：

? ? ? ? ? ? ? ??

　　如果令，則上式變?yōu)?#xff1a;

? ? ? ? ? ?

　　根據(jù)二次不等式相關(guān)知識(shí)，令：

? ? ? ? ?

? ? ?則論文總結(jié)了滿足下述關(guān)系式的a,b,c就能滿足上述要求了：

? ? ? ??

　　這個(gè)我也沒有去驗(yàn)證了。

? ? ??第三:?上面描述了校正暗角圖像的公式（帶參數(shù)）以及評(píng)價(jià)一副圖像是否有暗角的指標(biāo)，那么最后一步就是用這個(gè)指標(biāo)來確定公式的參數(shù)。我們未知的參數(shù)有5個(gè)，即a、b、c以及暗角的中心點(diǎn)。解這種受限的最優(yōu)問題是有專門的算法的，也是非常計(jì)算耗時(shí)的。因此，作者提出了一種快速的算法：Hill climbing with rejection of invalid solutions.

? ? ? 我稍微看了下這個(gè)算法，確實(shí)是個(gè)不錯(cuò)的想法，不過我并沒有去實(shí)踐，我采用了另外一種粗略的優(yōu)化方式。

? ? ? 首先，很明顯，為了計(jì)算這些最優(yōu)參數(shù)，我們沒有必要直接在原圖大小上直接計(jì)算，這點(diǎn)在原論文也有說明，我們即使把他們的寬高分別縮小到原圖的1/5甚至1/10計(jì)算出來的結(jié)果也不會(huì)有太大的差異，而這些參數(shù)的差異對(duì)最終的的結(jié)果影響也不大，但是計(jì)算量就能減少到原來的1/25和1/100。

? ? ?接著，我們觀察到a、b以及c的最優(yōu)結(jié)果范圍一般都在-2和2之間，并且從g的計(jì)算公式中知道，由于r是屬于0和1之間的正數(shù)，r^2, r^4, r^6在數(shù)值遞減的非常快，比如r=0.8,則三者對(duì)應(yīng)的結(jié)果就分別為0.64、0.4096、0.2621，因此，a和b及c在公式中的對(duì)最后結(jié)果的影響也依次越來越小。

? ? ?那么，我們可以參考以前的對(duì)比度保留之彩色圖像去色算法---基礎(chǔ)算法也可以上檔次一文中的優(yōu)化方式，把a(bǔ), b ,c 三個(gè)參數(shù)分別在[-2,2]之間離散化，考慮到參數(shù)稍微差異不會(huì)對(duì)結(jié)果有太大的影響，以及a、b、c的重要性，我們可以設(shè)置a、b、c三者的離散間隔分別為0.2、0.3、0.4，然后綜合上述判斷a、b、c是否為合理組合的函數(shù)，離散取樣的計(jì)算量組合大概有300種可能，對(duì)小圖計(jì)算著300種可能性的耗時(shí)是完全可以接受的，甚至考慮極端一點(diǎn)，把c的計(jì)算放到循環(huán)外側(cè)，即C取固定值0，然后計(jì)算出優(yōu)選的a和b值后，在計(jì)算C值。

? ? ? 上述計(jì)算過程并未考慮暗角中心點(diǎn)的范圍，我們是固定把暗角的中心點(diǎn)放置在圖像的正中心位置的，即 (Width/2, Height /2)，實(shí)際上，對(duì)于大部分拍攝的圖來說，暗角就是位于中心位置的，因此這種假設(shè)也無可厚非，因?yàn)榘到侵行挠?jì)算的增加必然會(huì)嚴(yán)重增加計(jì)算量， 為了求出暗角中心的合理位置，我們?cè)谟?jì)算出上述a、b、c后，在小圖中以一定步長按照公式計(jì)算出粗略的中心位置，再放大到原圖中去。

? ? ? 計(jì)算出上述a、b、c以及中心點(diǎn)后，就可以再次按照校正公式來進(jìn)行校正了，注意暗角的影響對(duì)每個(gè)通道都是等同的，因此，每個(gè)通道都應(yīng)該乘以相同的值。

　　下面貼出一些用論文中的算法處理的結(jié)果圖：

? ? ??? ?

? ? ??? ?

? ? ??? ?

? ? ??? ?

? ? ??? ?

? ? ?注意到上面最后一副圖的結(jié)果，那個(gè)女的婚紗以及衣服那些地方已經(jīng)嚴(yán)重的過曝了，我不清楚理論上造成整個(gè)原因的是什么，但是如果把計(jì)算i(L)的公式中的N修改為小一點(diǎn)的值，比如64，則可以避免這個(gè)結(jié)果。

? ? ?github上的那個(gè)代碼則對(duì)這個(gè)對(duì)數(shù)熵的過程做了一點(diǎn)改造，這個(gè)改造相當(dāng)暴力，就是什么呢，他把原來的[0，255]直接量化為8個(gè)等級(jí)，量化的依據(jù)是整形LOG2函數(shù)，即0->0,[1, 2]->1，[3, 4]->3,[5,8]->4,[9,16]->5,[17,32]->6,[33,64]->7,[65,128]->8,[129,255]->9, 原來的一條曲線映射函數(shù)變成了階躍函數(shù)了。這樣直方圖實(shí)際上只有9個(gè)值了，那么也不需要什么直方圖插值和高斯模糊了，直方圖則可以用整形表示，相對(duì)來說速度也能有很大的提升，并且也能克服上述最后一張圖片出現(xiàn)的那個(gè)瑕疵，其結(jié)果如下:

? ????

? ? ? 最后我們貼一些代碼對(duì)上述過程予以解釋：

　　 第一個(gè)是判斷a、b、c是否為合理值的函數(shù)：

// 按論文中公式18得條件判斷是否是合理的參數(shù) bool IsValidABC_Original(float A, float B, float C) {const int MAX_BRIGHTNESS_MULTIPLICATION = 3;if ((1 + A + B + C) > MAX_BRIGHTNESS_MULTIPLICATION) return false; // 當(dāng)r==1時(shí)，出現(xiàn)最大的亮度調(diào)整if (C == 0){if (A <= 0 || (A + 2 * B <= 0)) // 如果C==0，則根據(jù)公式(15)知，當(dāng)r==0時(shí)，A必須大于0，而當(dāng)r==1時(shí)，A+2B必須大于0return false;}else{float D = 4 * B * B - 12 * A * C;float qMins = (-2 * B - sqrtf(D)) / (6 * C); // 公式(17)float qPlus = (-2 * B + sqrtf(D)) / (6 * C);if (C < 0){if (D >= 0){if (qMins > 0 || qPlus < 1)return false;}elsereturn false;}else{if (D >= 0){if (!((qMins <= 0 && qPlus <= 0) || (qMins >= 1 && qPlus >= 1)))return false;}}}return true; }

　　可見除了文中公式的一些限制，我們還增加了幾個(gè)額外的小限制，比如最大亮度調(diào)節(jié)比例為3等等。

??　 第二個(gè)是計(jì)算指定參數(shù)下計(jì)算對(duì)數(shù)熵的過程：

// 計(jì)算不同參數(shù)修復(fù)后的圖像的整體對(duì)數(shù)熵 float CalculateEntropyFromImage(unsigned char *Src, int Width, int Height, float A, float B, float C, int CenterX, int CenterY) {float Histgram[256] = { 0 };float Invert = 1.0f / (CenterX * CenterX + CenterY * CenterY + 1);float Mul_Factor = 256.f / log(256.0f);for (int Y = 0; Y < Height; Y++){unsigned char *LinePS = Src + Y * Width * 4;int SquareY = (Y - CenterY) * (Y - CenterY);for (int X = 0; X < Width; X++){int Intensity = (LinePS[0] + (LinePS[1] << 1) + LinePS[2]) >> 2; // 公式(2)int RadiusSqua2 = (X - CenterX) * (X - CenterX) + SquareY;float R = RadiusSqua2 * Invert; // 公式(12)float Gain = 1 + (R * (A + R * (B + R * C))); // gain = 1 + a * r^2 + b * r^4 + c * r^6 ，公式(11)float Correction = Gain * Intensity; // 直接校正后的結(jié)果值if (Correction >= 255){Correction = 255; // It is possible that, due to local intensity increases applied by devignetting, the corrected image intensity range exceeds 255. Histgram[255]++; // In this case the algorithm simply adds histogram bins at the upper end without rescaling the distribution,}else{float Pos = Mul_Factor * log(Correction + 1); // 公式(6)int Int = int(Pos);Histgram[Int] += 1 - (Pos - Int); // 公式(7)Histgram[Int + 1] += Pos - Int;}LinePS += 4;}}float TempHist[256 + 2 * 4]; // SmoothRadius = 4TempHist[0] = Histgram[4]; TempHist[1] = Histgram[3]; TempHist[2] = Histgram[2]; TempHist[3] = Histgram[1];TempHist[260] = Histgram[254]; TempHist[261] = Histgram[253];TempHist[262] = Histgram[252]; TempHist[263] = Histgram[251];memcpy(TempHist + 4, Histgram, 256 * sizeof(float));for (int X = 0; X < 256; X++) // 公式(8),進(jìn)行一個(gè)平滑操作Histgram[X] = (TempHist[X] + 2 * TempHist[X + 1] + 3 * TempHist[X + 2] + 4 * TempHist[X + 3] + 5 * TempHist[X + 4] + 4 * TempHist[X + 5] + 3 * TempHist[X + 6] + 2 * TempHist[X + 7]) + TempHist[X + 8] / 25.0f;return CalculateEntropyFromHistgram_Original(Histgram, 256); }

　　其中計(jì)算熵的函數(shù)為：

// 從直方圖中計(jì)算熵值,原論文中直方圖肯定是浮點(diǎn)數(shù) float CalculateEntropyFromHistgram_Original(float Histgram[], int Length) {float Sum = 0;for (int X = 0; X < Length; X++){Sum += Histgram[X];}float Entropy = 0;for (int X = 0; X < Length; X++){if (Histgram[X] == 0) continue;float p = (float)Histgram[X] / Sum;Entropy += p * logf(p);}return -Entropy; }

　　其中

int Int = int(Pos);Histgram[Int] += 1 - (Pos - Int); // 公式(7)Histgram[Int + 1] += Pos - Int;

　　就是公式7所描述的過程的實(shí)現(xiàn)。

　　論文中的高斯模糊，我這里只是借助了一個(gè)簡單的線性模糊來代替，這個(gè)不會(huì)對(duì)結(jié)果造成本質(zhì)的區(qū)別。

? ? ?最后圖像的校正代碼大概如下：

int Devignetting_Original(unsigned char *Src, unsigned char *Dest, int Width, int Height) {if ((Src == NULL) || (Dest == NULL)) return STATUS_NULLREFRENCE;if ((Width <= 0) || (Height <= 0)) return STATUS_INVALIDPARAMETER;const float Step = 0.2f; //` 粗選A\B\C三個(gè)變量的步長float SmallestEntropy = 1000000000000.0f;float A = 0, B = 0, C = 0; int CenterX = Width / 2, CenterY = Height / 2; // 中心就默認(rèn)為圖片中心for (int X = -10; X <= 10; X++) // 多次測(cè)試，表面最優(yōu)的A\B\C的范圍均在[-2,2]之間{for (int Y = -10; Y <= 10; Y++){for (int Z = -10; Z <= 10; Z++){if (IsValidABC_Original(X * Step, Y * Step, Z * Step) == true) // 判斷這個(gè)組合時(shí)候有效{float Entropy = CalculateEntropyFromImage(Src, Width, Height, X * Step, Y * Step, Z * Step, CenterX, CenterY);if (Entropy < SmallestEntropy) // 取熵值最小的{A = X * Step;B = Y * Step;C = Z * Step;SmallestEntropy = Entropy;}}}}}float Invert = 1.0 / (CenterX * CenterX + CenterY * CenterY + 1);for (int Y = 0; Y < Height; Y++){byte *LinePS = Src + Y * Width * 4;byte *LinePD = Dest + Y * Width * 4;int SquareY = (Y - CenterY) * (Y - CenterY);for (int X = 0; X < Width; X++){int RadiusSqua2 = (X - CenterX) * (X - CenterX) + SquareY;float R2 = RadiusSqua2 * Invert; // 公式12float Gain = 1 + (R2 * (A + R2 * (B + R2 * C))); // 公式11 LinePD[0] = ClampToByte(LinePS[0] * Gain); // 增益LinePD[1] = ClampToByte(LinePS[1] * Gain);LinePD[2] = ClampToByte(LinePS[2] * Gain);LinePD += 4;LinePS += 4;}}return STATUS_OK; }

　 ?上面的代碼是未經(jīng)特別的優(yōu)化的，只是表達(dá)了大概的意思，并且把暗角的中心點(diǎn)默認(rèn)為圖像的中心點(diǎn)，如果暗角的中心點(diǎn)不是圖像中心，要注意計(jì)算r時(shí)可能會(huì)出現(xiàn)r>1的情況，這個(gè)時(shí)候一定要注意重置r=1，否則結(jié)果就會(huì)完全不對(duì)了。

? ??經(jīng)過測(cè)試，對(duì)于沒有暗角的圖像，一般來說該算法不會(huì)對(duì)圖片產(chǎn)生很大的影響，很多圖片基本是無變換的，有些可能會(huì)有點(diǎn)區(qū)別，也就是整體變亮一點(diǎn)而已，因此，還是有比較好的普適性的。

? ? ?由于論文中的暗角強(qiáng)度減弱公式是根據(jù)一些光學(xué)原理獲得的，其可能并不符合一些軟件自己添加的暗角的規(guī)律，所以如果你用這中測(cè)試圖去測(cè)試算法，可能不會(huì)獲得非常滿意的結(jié)果，

? ? ?算法速度上，經(jīng)過我上面的描述的優(yōu)化后的算法，對(duì)于800*600的彩色圖，一般能在15ms左右處理完成，而未優(yōu)化的上述代碼可能要10000ms。?基本上算法的核心代碼為已經(jīng)貼出，這里不共享我優(yōu)化后的快速實(shí)現(xiàn)，有能力的朋友自然能快速搞定他們，共享一個(gè)測(cè)試工程吧。

? ? ??http://files.cnblogs.com/files/Imageshop/Devignetting_Test.rar

總結(jié)

以上是生活随笔為你收集整理的图像增强系列之图像自动去暗角算法的全部內(nèi)容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網(wǎng)站內(nèi)容還不錯(cuò)，歡迎將生活随笔推薦給好友。