Retinex是一种常用的建立在科学实验和科学分析基础上的图像增强方法。就跟Matlab是由Matrix和Laboratory合成的一样，Retinex也是由两个单词合成的一个词语，他们分别是retina 和cortex，即：视网膜和皮层。Land的retinex模式是建立在以下三个假设之上的：

真实世界是无颜色的，我们所感知的颜色是光与物质的相互作用的结果。我们见到的水是无色的，但是水膜—肥皂膜却是显现五彩缤纷，那是薄膜表面光干涉的结果。
每一颜色区域由给定波长的红、绿、蓝三原色构成的；
三原色决定了每个单位区域的颜色。
Retinex理论的基础理论是物体的颜色是由物体对长波（红色）、中波（绿色）、短波（蓝色）光线的反射能力来决定的，而不是由反射光强度的绝对值来决定的，物体的色彩不受光照非均匀性的影响，具有一致性，即retinex是以色感一致性（颜色恒常性）为基础的。不同于传统的线性、非线性的只能增强图像某一类特征的方法，Retinex可以在动态范围压缩、边缘增强和颜色恒常三个方面达到平衡，因此可以对各种不同类型的图像进行自适应的增强。

原理：

Retinex 理论的基本内容

物体的颜色是由物体对长波（红）、中波（绿）和短波（蓝）光线的反射能力决定的，而不是由反射光强度的绝对值决定的；
物体的色彩不受光照非均性的影响，具有一致性，
即Retinex理论是以色感一致性（颜色恒常性）为基础的。
如下图所示，观察者所看到的物体的图像S是由物体表面对入射光L反射得到的，反射率R由物体本身决定，不受入射光L变化。

基于Retinex的图像增强的目的就是从原始图像S中估计出光照L，从而分解出R，消除光照不均的影响，以改善图像的视觉效果,正如人类视觉系统那样。

Retinex理论，与降噪类似，该理论的关键就是合理地假设了图像的构成。如果将观察者看到的图像看成是一幅带有乘性噪声的图像，那么入射光的分量就是一种乘性的，相对均匀，且变换缓慢的噪声。Retinex算法所做的就是合理地估计图像中各个位置的噪声，并除去它。

在极端情况下，我们大可以认为整幅图像中的分量都是均匀的，那么最简单的估计照度L的方式就是在将图像变换到对数域后对整幅图像求均值。

                                                  Log[R(x,y)] = Log[S(x,y)]-Log[L(x,y)];　

要求得R(x,y)只需得到L(x,y)即可，但是根据数学的理论，L(x,y)是不能够求得的，只能近似求出。我们用S(x,y)和一个高斯核的卷积来近似表示L(x,y)。所以R(x,y)可用下式表示：

                                                  Log(R(x,y)）=(Log(S(x,y))-Log(S(x,y)*G(x,y))) 

由于R是对数域的输出，要转换为数字图像，必须将他们量化为[0,255]的数字图像范畴，所以要进行线性量化，于是就有了：（G和b)为经验参数

 再通过Simplest Color Balance方式，他把数据按照一定的百分比去除最小和最大的部分，然后中间的部分重新线性量化到0和255之间，因为如果直接对MSR处理的结果进行量化，得到的图像往往整体偏灰度，这是由于原始的彩色值经过log处理后的数据范围就比较小了，这样各通道之间的差异也很小，而之后的线性量化比log曲线要平滑很多，因此整体就丧失了彩色。 

        为此，研究者又开发出一种称之为带色彩恢复的多尺度视网膜增强算法（MSRCR,Multi-Scale Retinex with Color Restoration) ，具体讨论的过程详见 <A Multiscale Retinex for Bridging the Gap Between Color Images and the Human Observation of Scenes>这篇论文，但是论文里为了这个又引入了太多的可调参数，增加了算法的复杂性，不利于自动化实现。

    直接从量化的方式上入手，引入了均值和均方差的概念，再加上一个控制图像动态的参数来实现无色偏的调节过程，简要描述如下。

    （1）分别计算出 Log[R(x,y)]中R/G/B各通道数据的均值Mean和均方差Var（注意是均方差）。

    （2）利用类似下述公式计算各通道的Min和Max值。

            Min = Mean – Dynamic * Var;  

            Max = Mean + Dynamic * Var;
    (3)  对Log[R(x,y)]的每一个值Value，进行线性映射： 

           R(x,y) = ( Value – Min ) / (Max – Min) * (255-0) ，同时要注意增加一个溢出判断,即：

           if (R(x,y) > 255)  R(x,y) =255; else if (R(x,y) < 0) R(x,y)=0;

import cv2
import numpy as np
import math




def replaceZeroes(data):
    min_nonzero = min(data[np.nonzero(data)])
    data[data == 0] = min_nonzero
    return data


def simple_color_balance(input_img, s1, s2):
    h, w = input_img.shape[:2]
    out_img = np.zeros([h, w])
    sort_img = input_img.copy()
    one_dim_array = sort_img.flatten()  # 转化为一维数组
    sort_array = sorted(one_dim_array)  # 对一维数组按升序排序

    per1 = int((h * w) * s1 / 100)
    minvalue = sort_array[per1]

    per2 = int((h * w) * s2 / 100)
    maxvalue = sort_array[(h * w) - 1 - per2]

    # 实施简单白平衡算法
    if (maxvalue <= minvalue):
        for i in range(h):
            for j in range(w):
                out_img[i, j] = maxvalue
    else:
        scale = 255.0 / (maxvalue - minvalue)
        for m in range(h):
            for n in range(w):
                if (input_img[m, n] < minvalue):
                    out_img[m, n] = 0
                elif (input_img[m, n] > maxvalue):
                    out_img[m, n] = 255
                else:
                    out_img[m, n] = scale * (input_img[m, n] - minvalue)  # 映射中间段的图像像素

    out_img = cv2.convertScaleAbs(out_img)
    return out_img
def MSRCR(img, scales, s1, s2):
    h, w = img.shape[:2]
    # print(h, w)
    scles_size = len(scales)
    img = np.array(img, dtype=np.float64)
    # print(img)
    log_R = np.zeros((h, w), dtype=np.float64)
    img_sum = np.add(img[:,:,0],img[:,:,1],img[:,:,2])
    # print(img_sum)
    # print("11111111111111111111111111111111")
    img_sum = replaceZeroes(img_sum)
    # print(img_sum)
    gray_img = []

    for j in range(3):
        img[:, :, j] = replaceZeroes(img[:, :, j])
        for i in range(0, scles_size):
            L_blur = cv2.GaussianBlur(img[:, :, j], (scales[i], scales[i]), 0)
            L_blur = replaceZeroes(L_blur)

            dst_img = cv2.log(img[:, :, j])
            dst_Lblur = cv2.log(L_blur)
            # dst_ixl = cv2.multiply(dst_img, dst_Lblur)
            log_R += cv2.subtract(dst_img, dst_Lblur)
            # print(i)
        # print(scles_size)
        MSR = log_R / 3.0
        '''
            img_sum_log = np.zeros((h, w))
            for i in range(0, h):
                for k in range(0, w):
                    img_sum_log[i,k] = 125.0*math.log(img[i,k,j]) - math.log(img_sum[i,k])
            MSRCR = MSR * (img_sum_log[:, :])
            print(img_sum)
            # x = cv2.log(img_sum)
            '''
        MSRCR = MSR * (cv2.log(125.0 * img[:, :, j]) - cv2.log(img_sum))
        gray = simple_color_balance(MSRCR, s1, s2)
        gray_img.append(gray)

    return gray_img

if __name__ == '__main__':
    scales = [15, 101, 301]
    s1, s2 = 2,3
    src_img = cv2.imread('4.jpg')
    src_img = cv2.cvtColor(src_img, cv2.COLOR_BGR2RGB)
    cv2.imshow('img', src_img)
    MSRCR_Out = MSRCR(src_img, scales, s1, s2)

    result = cv2.merge([MSRCR_Out[0], MSRCR_Out[1], MSRCR_Out[2]])


    cv2.imshow('MSR_result', result)
    cv2.waitKey(0)
    cv2.destroyAllWindows()

————————————————
版权声明：本文为CSDN博主「SongpingWang」的原创文章，遵循CC 4.0 BY-SA版权协议，转载请附上原文出处链接及本声明。

致谢：
https://blog.csdn.net/wsp_1138886114/article/details/83096109

带色彩恢复的多尺度视网膜增强算法（MSRCR）的原理、实现及应用。 – Imageshop – 博客园 (cnblogs.com)