图像复原--Color Balance and Fusion for Underwater Image Enhancement

图像复原–Color Balance and Fusion for Underwater Image Enhancement

参考链接: （建议先阅读原论文和参考链接内容，本文纯属记录性质）

https://blog.csdn.net/qq_43287277/article/details/104901501
https://github.com/fergaletto/Color-Balance-and-fusion-for-underwater-image-enhancement.-.
https://ieeexplore.ieee.org/stamp/stamp.jsp?tp=&arnumber=8058463&isnumber=8071125

1. 摘要:

​ 对水下因为介质散射和吸收导致的图像质量下降问题, 提出了一种融合颜色补偿和白平衡图像的算法，同时定义两幅融合图像及其相关的权重图，强化边缘和颜色对比度向输出图像的传输。为了避免权值变换在重建图像的低频分量中产生伪影，还采用了多尺度融合策略。实验表明增强图像全局对比度和边缘清晰度能够明显改善。

2. 背景

​ 水下成像通常包含三个分量，分别是直接分量(direct composed)、前向散射(forward scattering)、后向散射(back scattering)，前向散射在实际模型中经常被忽略， 所以最终的成像是直接分量和后向散射两个分量的叠加：
$$I_c=J_c{e}^{-\eta d(x)}+B_c^{\infty }(1-e^{-\eta d(x)})\text{\hspace{1em}(1)}$$
$J_c$是未衰减的直接分量，$d$是观察者与物体之间的距离函数，$\eta$是衰减系数，$B_c^{\infty }$是背景光，指数项是通过水下介质的传输公式。sea-Thru算法对$\eta$取同一个值提出了质疑，采用两个不同的系数控制衰减，取得了比较好的结果。

3. 算法实现

​ 文章提出的图像增强方法采用两步策略，分别为白平衡和图像融合。白平衡的目的是补偿颜色随深度的选择性吸收而产生的颜色投射，而图像融合则被考虑用于增强场景的边缘和细节，以减轻由于后向散射而造成的对比度损失。

白平衡

​ 不同光源色温也不同。当拍摄画面出现偏色时，多采用白平衡来解决这一问题。传统的白平衡有很多种方法，比如灰色世界，灰色边缘，灰色阴影等，原文作者对各种白平衡方法进行比较，发现灰色世界法对于水下图像有良好的效果，但经过观察发现，这种方法存在严重的红色伪影。

​ 灰色世界白平衡法假设对于一幅图像来说, RGB 三个分量的平均值趋于同一个灰度K。灰度K可以自己定义，也可以使用三个通道的平均值，之后使用K比上每一个通道的平均值，算出每一个通道的权重，若一个通道在图像中的成分过多，则权值会小于1，反之，权重值会大于1。最后再将各个通道依照此权重进行相加即可。灰色世界白平衡法对水下图像使用会产生的红色伪影。按照其算法原理来分析，**红色分量在水下图像（以蓝色或绿色为主题颜色）中，比重实在太小，导致权值过大，从而被极大地放大。那么调整红色的权重就很关键。**原文作者提出的白色补偿白平衡法的表达式如下:
$$I_{rc}(x)=I_r(x)+a(\overline{I_g}-\overline{I_r})(1-I_r(x))I_g(x)\text{\hspace{1em}(2)}$$
其中，$I_r,I_g$表示图像$I$的红色和绿色通道，每个通道的值要根据本通道的最大最小值范围归一化到[0,1]内，$\overline{I_r}$，$\overline{I_g}$表示$I_r,I_g$的平均值，$a$表示常量参数，在实际测试中，$a=1$的值适合于各种照明条件和采集设置。式(2)具体含义为，$I_{rc}$代表增加补偿项的红色分量，$a$通常为1，$\overline{I_r}$，$\overline{I_g}$相减得到一个系数，这个系数与红色分量成反比，红色分量越小时，系数越大，$1-I_r$代表后面绿色通道的补偿系数，如果红色分量越多，那么绿色通道的补偿则越少，相反，则越大。可以看出这个补偿项是将绿色成分的一部分拿来补偿红色，另外也加入了图像整体红色分量的比重和各个像素点的红色分量比重，考虑的非常完备。

不过有的时候，因为水下过于浑浊或者微生物太多的原因，水下会呈现绿色而非蓝色，这个时候只需要再将蓝色通道进行相同方式的补偿即可。

'''
白平衡
'''
def white_balance(img):B, G, R = np.double(img[:, :, 0]), np.double(img[:, :, 1]), np.double(img[:, :, 2])# 获取每个通道归一化后[0, 1]的值R_norm = signal_channel_normal(R)G_norm = signal_channel_normal(G)B_norm = signal_channel_normal(B)# 获取归一化后的每个通道的均值mean_R_norm = np.mean(R_norm)mean_G_norm = np.mean(G_norm)mean_B_norm = np.mean(B_norm)# 对红蓝通道补偿rc = R_norm + (mean_G_norm - mean_R_norm)*(1-R_norm)*G_normbc = B_norm + (mean_G_norm - mean_B_norm)*(1-B_norm)*G_norm# 判断补偿参数,大于1的变成1rc[rc > 1] = 1bc[bc > 1] = 1# 现在红蓝通道是补偿后的值newr = rc # 红色通道补偿newb = bc # 蓝色通道补偿newg = G_norm# 进行白平衡B_ave, G_ave, R_ave = np.mean(newb), np.mean(newg), np.mean(newr)# factor = G_ave/B_ave;K = (B_ave + G_ave + R_ave) / 3Kb, Kg, Kr = K / B_ave, K / G_ave, K / R_aveBa = (newb * Kb)Ga = (newg * Kg)Ra = (newr * Kr)for i in range(len(Ba)):for j in range(len(Ba[0])):Ba[i][j] = 255 if Ba[i][j] > 255 else Ba[i][j]Ga[i][j] = 255 if Ga[i][j] > 255 else Ga[i][j]Ra[i][j] = 255 if Ra[i][j] > 255 else Ra[i][j]# print(np.mean(Ba), np.mean(Ga), np.mean(Ra))dst_img = np.double(np.zeros_like(img))dst_img[:, :, 0] = Badst_img[:, :, 1] = Gadst_img[:, :, 2] = Rareturn dst_img

锐化

​ 白平衡的目的是补偿颜色随深度的选择性吸收而产生的颜色投射，而图像融合则用于增强场景的边缘和细节，以减轻由于后向散射而造成的对比度损失。接下来对以上算法处理之后得到的图像进行锐化处理，传统的锐化公式:
$$S=I+\beta (1-G\ast I)$$
​

$I$是要锐化的图像（白平衡图像），$G\ast I$表示对I进行高斯滤波，β是参数。在实际应用中，β的选择并不简单。较小的β不能锐化I，但过大的β会导致过饱和，具有较亮的高光和较暗的阴影。为了避免这个问题，新的锐化图像算法表达式如下：
$$S=\frac{(I+N(1-G\ast I))}{2}$$
N表示线性归一化运算，也是直方图拉伸运算，上述方式称为normalized unsharp masking不需要考虑$\beta$值的选取，此处要注意的是，直方图拉伸并不等于直方图均衡化。

def imageSharp_(img):# 图像线性拉伸锐化# S=(I+N(1-G*I))/2# 高斯滤波sigma = 20kernel_size = (0,0) # 0,0函数会自己计算卷积核Igauss_ = imgN=30for i in range(N):Igauss = cv2.GaussianBlur(img, kernel_size, sigmaX=sigma, borderType=cv2.BORDER_REPLICATE)Igauss_ = np.minimum(img, Igauss)gain = 1 # Norm = img-gain*Igauss_#执行直方图拉伸操作,对三个通道分别执行B, G, R = np.double(Norm[:, :, 0]), np.double(Norm[:, :, 1]), np.double(Norm[:, :, 2])# 计算原图灰度直方图equ_R = cv2.normalize(R, None, 1.0, 0.0, cv2.NORM_MINMAX, dtype=cv2.CV_64F)equ_G = cv2.normalize(G, None, 1.0, 0.0, cv2.NORM_MINMAX, dtype=cv2.CV_64F)equ_B = cv2.normalize(B, None, 1.0, 0.0, cv2.NORM_MINMAX, dtype=cv2.CV_64F)hist_norm_img = np.double(np.zeros_like(img))hist_norm_img[:, :, 0] = equ_Bhist_norm_img[:, :, 1] = equ_Ghist_norm_img[:, :, 2] = equ_Rout = (img + hist_norm_img)/2return out

权重计算

​ 图像融合肯定是要依靠一定的权重的，原文作者提供了如下几个权重：Laplacian contrast weight（拉普拉斯对比度权重）、Saliency weight（显著性权重）、Saturation weight（饱和权重）每个权重都要对白平衡后的图像和锐化后的图像做。

Laplacian contrast weight（拉普拉斯对比度权重）: 主要是针对亮度通道进行拉普拉斯滤波而来的。首先把rgb图像转化为hsv图像，提取亮度层，然后进行拉普拉斯滤波便可。

def getLaplacianContrast_(img):#提取亮通道lab_color = color.rgb2lab(img)R1 = np.double(lab_color[:, :, 0])/255localContrast = img[:, :, 0]- R1localContrast = np.sqrt((1.0/3)*(np.power((lab_color[:, :, 0]-R1),2)+np.power((lab_color[:, :, 1]-R1),2)+np.power((lab_color[:, :, 2]-R1),2)))return localContrast

Saliency weight (显著性权重): 首先要将图像从RGB转化到LAB，L是明度，A, B是两种颜色。然后计算三个通道各个的平均值，记为$mL,mA,mB$再计算各个通道的均方差，最后进行加和，表达式如下：
$$W_{sal}=(L-mL{)}^2+(A-mA{)}^2+(B-mB{)}^2$$

def getSaliencyWeight(img):lab_color = color.rgb2lab(img)l = np.double(lab_color[:, :, 0])a = np.double(lab_color[:, :, 1])b = np.double(lab_color[:, :, 2])lm = np.mean(l)am = np.mean(a)bm = np.mean(b)sm_weight = np.power((l-lm),2) +np.power((a-am),2)+ np.power((b-bm),2)return sm_weight

Saturation weight (饱和权重):主要是利用RGB三个通道与亮度的偏差来计算的，表达式如下:
$$W_{sat}=\sqrt{\frac{1}{3}[(R-I{)}^2+(G-I{)}^2+(B-I{)}^2]}$$
其中$RGB$分别表示的是红绿蓝三个通道，$I$代表亮度通道。

def getSaturationWeight_(img):lab_color = color.rgb2lab(img)R1 = np.double(lab_color[:, :, 1]) / 255localContrast = np.sqrt((1.0/3)*(np.power((lab_color[:, :, 0]-R1),2)+np.power((lab_color[:, :, 1]-R1),2)+np.power((lab_color[:, :, 2]-R1),2)))return localContrast

融合

​ 图像融合采用的是高斯金字塔和拉普拉斯金字塔来实现，图像金字塔其实就是由一组不同尺寸的图像组成的，把尺寸小的图像放在最上层，把尺寸大的图像放在最下层，这样看起来就像一座金字塔，这是最直观的解释。但是这一组不同尺寸的图片可不能是随意组成的，他们都是由最底层的原始图像，经过滤波然后隔行隔列采样而来，所以上一层的图像的尺寸是下一层图像尺寸的四分之一。之所以叫高斯金字塔，是因为他的滤波过程选择的是高斯低通滤波器，这样构成的图像金字塔就是高斯金字塔了，而拉普拉斯金字塔是用来重建高斯金字塔的，就是高斯金字塔的逆过程。

​ 作者将白平衡图像的三个权重和锐化图像的三个权重看成两个整体，再将二者进行归一化处理，用两个归一化权重来构建高斯金字塔，再使用锐化后的白平衡图像和进行过gamma变换的白平衡图像建立拉普拉斯金字塔。这里使用gamma变换是将图像做一个非线性的变换（幂函数映射），是一种简单的图像增强方法，为的是提升图像的视觉效果。接下来就是将权重金字塔的对应层与图像金字塔的对应层相乘，形成一个新的加入权重的图像金字塔，具体表达式如下:
$$N_k=Gw{b}_k\times Lw{b}_k+G{s}_k\times L{s}_k$$
其中$k$代表金字塔的层数，$Gwb$是对白平衡图像归一化权重的高斯金字塔，$Lwb$是白平衡图像的拉普拉斯金字塔，$G_s$为对锐化图像归一化权重的高斯金字塔，$L_s$为锐化图像的拉普拉斯金字塔。

def gaussian_pyramid2(img, level):h = 1/16 * np.array([1, 4, 6, 4, 1])filt = np.outer(h, h)out = [convolve2d(img, filt, mode='same', boundary='wrap')]temp_img = img.copy()for i in range(2, level+1):temp_img = temp_img[::2, ::2]out.append(convolve2d(temp_img, filt, mode='same', boundary='wrap'))return outdef laplacian_pyramid2(img, level):h = 1/16 * np.array([1, 4, 6, 4, 1])out = [img]temp_img = img.copy()for i in range(2, level + 1):temp_img = temp_img[::2, ::2]out.append(temp_img.copy())# Calculate the DoG (Difference of Gaussians)for i in range(level - 1):m, n = out[i].shapeout[i] = out[i] - cv2.resize(out[i + 1], (n, m))return outdef pyramid_reconstruct2(pyramid):level = len(pyramid)for i in range(level - 1, 0, -1):m, n = pyramid[i - 1].shapepyramid[i - 1] = pyramid[i - 1] + cv2.resize(pyramid[i], (n, m))out = pyramid[0]return out

4. 结果图

上面图是对应的naive_funsion结果，中间图是对应的金字塔多尺度融合结果， 最后三联图是对应RGB三通道多尺度融合后每个通道结果，native_funsion明显有局部地方出现偏白的问题：

5. python复现中出现的问题: 使用cv2.imshow/plt.imshow显示图片正常，但保存下来的图片不正常或者全黑

最后输出结果保存的时候，用opencv或matplotlib都可以正常显示，但是保存下来的图像是黑的，感谢以下博主提供的解决方案。

参考链接: https://blog.csdn.net/qq_40682833/article/details/124049695