Retinex图像增强:从SSR到MSRCR的MATLAB实现
写这个系列写到第九篇终于要把图像增强里相对“玄”的一支拿出来聊了——Retinex。前八篇我们把灰度变换、直方图均衡化、空域滤波、频域滤波这些基础操作挨个过了一遍日常的对比度增强、噪声抑制、边缘锐化基本都够用。但如果你手里的图是夜间监控截图、逆光人像、或者室内暗光拍出的产品照片直方图均衡化很容易把噪声一起放大拉普拉斯锐化只会让阴影更黑。光在像素值上做文章已经解决不了“光照不均”这种结构性问题。Retinex消掉光照、保留反射正是这类场景需要的解法。这篇的核心内容围绕单尺度RetinexSSR、多尺度RetinexMSR、带颜色恢复的多尺度RetinexMSRCR展开附带MATLAB实现和调参经验适合已经掌握基础图像处理、想处理低照度或复杂光照图像的读者。1. 为什么做到第九篇才聊Retinex先搭好图像增强的知识坐标系1.1 前八篇走的是“像素值路线”Retinex换了一条赛道前八篇的思路本质上都是把像素值当成最终处理对象——直方图均衡化在重新分配灰度级拉普拉斯算子在做差分中值滤波在取邻域统计值。这些方法在均匀光照下效果稳定但光照一旦不均衡像素值本身就不是单纯的物体属性了。举个例子同一块白布暗处拍出来可能是灰度50亮处能到200直接做直方图均衡化暗处的白布和灰布很可能分不开。Retinex的思路则完全不同。它把图像分解成照度分量和反射分量认为物体颜色由反射分量决定照度是外在干扰。先估计光照再还原反射这就等于直接把光照不均这个“病因”去掉了而不是在“症状”上反复调整。这也是为什么我把这个主题放到第九篇——先掌握基础操作才能理解Retinex是在哪个环节替代了它们。1.2 Retinex理论的一句话拆解与数学表达Retinex模型极其简洁S(x,y) L(x,y) × R(x,y)。S是传感器看到的图像L是入射光照R是物体反射率取值范围在0到1之间。人眼真正感知的主要是R。Retinex处理就是三步从S中估计出L用S除以L得到R。为什么通常在对数域计算一方面乘除运算变成加减简单很多另一方面人眼对亮度的感知本身接近对数响应log域运算更符合视觉模型还能压缩动态范围。公式写出来就是r(x,y) log(S(x,y)) - log(L_est(x,y))其中L_est是环绕函数F与S的卷积结果。这里有一个容易踩的坑log后面通常要加1防止S为0时出现负无穷。如果后续还要还原到0-255这个偏移量会影响映射需要配合拉伸步骤处理别忽略。1.3 SSR、MSR、MSRCR选型对比这三个版本的关系是层层叠加不是互相替代。算法核心特点优点缺点适用场景SSR单一尺度估计光照计算量最小边缘增强明显光晕明显动态范围压缩不够偏色灰度图、快速预览MSR多个尺度加权平均亮暗细节兼顾颜色容易偏灰灰度/彩色增强基础版本MSRCR引入颜色恢复因子颜色自然接近人眼感知参数多计算最重彩色图像视觉质量优先尺度本质是高斯核的宽度。小尺度保留细节但容易产生光晕大尺度把光照信息压得平整。MSR把多个尺度加权原理上类似用多个窗口看同一个场景——小窗口看细节大窗口看全局光照。单独做MSR时RGB三个通道分别算同样运算通道间的相对关系会被改变颜色就漂了。所以MSRCR加上一个颜色恢复因子效果会自然很多。2. 核心原理拆解环绕函数、对数域运算与每一步的“为什么”2.1 高斯环绕函数为什么必须归一化环绕函数F(x,y) λ·exp(-(x²y²)/c²)它模拟的是一个像素感受到的周围光照。参数c越小周围像素参与估计的范围越小估计出的光照越接近细节层c越大光照越平坦更像整体光照。很多人写代码时忽略了归一化导致结果整体偏亮或偏暗。因为卷积后得到的L_est应该是S的一个加权平均如果核的总能量不等于1L_est就会整体偏大或偏小。log差分之后相当于给图像加了一个常数偏移虽然最终的拉伸步骤能抵消一部分但在MSR多尺度叠加时各尺度的常数偏移不统一结果会发灰。所以我在MATLAB里用fspecial生成高斯核之后一定会补一行g g / sum(g(:))确保核和为1。2.2 对数域先减后加到底在干什么很多第一次写Retinex的人会疑惑log(S) - log(L)之后为什么不直接exp回去确实在理想情况下r≈log(R)exp(r)就是反射率R。但实际估计的L不可能是真实光照的完美拟合log域差分后得到的值有正有负动态范围被压缩直接exp回去会得到对比度很低的灰图。这一步的处理思路通常是把r当作一个增强后的对数域图像做线性映射或分位数截断再映射回0到1。换句话说Retinex到这里更像一种特征变换而不是严格的物理还原。我在实现里用prctile取1%和99%的分位数作为上下界做截断比min-max全幅度拉伸更能对抗极端亮暗像素。2.3 从SSR到MSR多尺度的意义SSR只用一个尺度选了小的细节清楚但光照估计不准选了大的光照平整但细节丢失。MSR的做法是同时跑几个尺度然后对输出加权求和R_MSR Σ w_n × r_n常用三个尺度是15、80、250权重各取1/3。为什么是这三个数15个像素的窗口看局部纹理80个像素看中等范围的照明变化250个像素看全局光照分布。这个组合在很多实践里都验证过适合自然图像。如果你的图像分辨率特别高或特别低等比缩放即可不要死守这三个值。2.4 Matlab里构建高斯核的细节与边界处理MATLAB里常见的是fspecial(gaussian, hsize, sigma)但Retinex环绕函数用的是exp(-r²/c²)和fspecial的exp(-r²/(2σ²))差了一个2。所以如果用fspecialsigma要取c/sqrt(2)。这是非常容易出问题的地方参数对不上效果天差地别。边界处理也是坑。卷积时边界附近可用像素少默认的zero padding会让图像四周发黑。我建议imfilter指定replicate用边缘值向外复制能明显减少黑边。还有一个更彻底的办法先用padarray向外扩展一圈处理完再裁掉效果更干净代价是多写两行代码。3. 实操在MATLAB里从零实现SSR与MSRCR附完整代码3.1 读图、归一化、转double预处理三步MATLAB读图进来可能是uint8直接做log和卷积会出问题。第一步用imread第二步转double第三步归一化到[0,1]。顺序不要反先转double再除以255。为什么非要转double一是imfilter对double类型能正确处理小数卷积二是log(0)的问题可以在double数据上做偏移处理uint8的取整误差在暗图像上会被放大。img imread(low_light.jpg); img im2double(img); % 直接把uint8转成[0,1]的double如果图像是灰度图三通道处理的逻辑也同样适用只是循环只跑一次。3.2 SSR函数完整代码与逐段说明这里给出的SSR实现输入是double类型的RGB图像输出是处理后的对数域结果。注意内部用c/sqrt(2)做了高斯核换算避免前面说的那个坑。function out ssr(img, c) % img: double型范围[0,1] % c: 尺度参数例如15、80、250 [h, w, ~] size(img); out zeros(size(img)); for ch 1:3 channel img(:, :, ch); % 构建全图尺寸的高斯环绕核 % Retinex环绕函数是 exp(-r^2/c^2)而fspecial用的是 exp(-r^2/(2*sigma^2)) % 所以 sigma c / sqrt(2) g fspecial(gaussian, [h w], c / sqrt(2)); g g / sum(g(:)); % 归一化保证核能量为1 % 估计光照分量空域卷积 L imfilter(channel, g, replicate, conv); % 对数域差分log(S) - log(L) r log(channel 1) - log(L 1); out(:, :, ch) r; end end关键点有两个。一个是核大小直接取[h, w]因为Retinex理论要求全图卷积不能像普通高斯平滑那样取一个小窗口。小窗口卷积只考虑局部邻域无法估计出整体光照。另一个是imfilter的边界选项replicate它能减少黑边效果比默认的零填充好得多。3.3 MSRCR完整代码与颜色恢复参数MSRCR就是在MSR的基础上每个通道乘一个颜色恢复因子C。这个因子的核心思路大致是如果某通道像素值在当前像素的总亮度中占比高就保留甚至增强这个通道占比低就适当抑制。经典参数α125、β46但这两个参数是面向普通场景的经验值。我自己实测时发现偏色严重的图需要把β调小到20到30颜色恢复会更平缓暗光图适当调小α也有明显改善。function out msrcr(img, scales, alpha, beta) % img: double型RGB范围[0,1] % scales: 三个尺度的数组例如 [15 80 250] % alpha: 颜色恢复调节参数经典取125 % beta: 颜色恢复增益经典取46 n length(scales); [h, w, ~] size(img); msr_out zeros(size(img)); % 第一步多尺度Retinex for ch 1:3 channel img(:, :, ch); msr zeros(h, w); for k 1:n g fspecial(gaussian, [h w], scales(k) / sqrt(2)); g g / sum(g(:)); L imfilter(channel, g, replicate, conv); msr msr (log(channel 1) - log(L 1)) / n; end msr_out(:, :, ch) msr; end % 第二步颜色恢复因子 out zeros(size(img)); sum_img img(:, :, 1) img(:, :, 2) img(:, :, 3) 3; for ch 1:3 C beta * (log(alpha * img(:, :, ch) 1) - log(sum_img)); out(:, :, ch) C .* msr_out(:, :, ch); end % 第三步分位数截断拉伸 for ch 1:3 ch_data out(:, :, ch); lo prctile(ch_data(:), 1); hi prctile(ch_data(:), 99); out(:, :, ch) (ch_data - lo) / max(hi - lo, 1e-6); out(:, :, ch) max(0, min(1, out(:, :, ch))); end end第三步的拉伸极其重要。如果不做分位数截断直接imshow你大概率看到一张灰蒙蒙的图。1%和99%的分位数截断能拉出有效对比度又不让极端亮暗像素主导整个映射范围。3.4 输出图像处理截断拉伸与对比度矫正Retinex输出经常是灰蒙蒙的直接imshow几乎看不清。我的习惯是先做1%到99%分位数截断再线性映射到[0,1]最后用imadjust收尾。注意不要上来就histeqRetinex输出已经改变了直方图形状再用直方图均衡化会把刚才保留的细节再次推平。如果是偏暗的结果我一般会再补一步Gamma校正gamma取0.8到0.9能把暗部细节再提一提而不会把亮部直接打爆。这个组合我用了很多次比单纯调Retinex参数要可控。3.5 实测一张低照度夜景的增强效果观察说一次实测体会。我拿了一张走廊监控角度的夜景截图原图暗部接近纯黑角落里的杂物轮廓完全看不清。直方图均衡化之后噪声被放大得很厉害画面发灰。改用MSRCR之后墙壁纹理和角落轮廓清楚了很多颜色也没有明显偏色。对比SSR和MSRCRSSR出图偏灰且细节略生硬MSRCR色彩自然很多但暗部细节反而比MSR少一点。这符合原理——颜色恢复因子在还原色彩的同时会压缩一部分亮度差异。所以如果只要看细节MSR就够了如果要交付给人看MSRCR更合适。4. 同态滤波另一个攻光照不均的经典思路含MATLAB代码4.1 频域怎么看待光照分量Retinex在空域估计光照同态滤波则从频率角度切分光照分量变化缓慢对应低频反射分量包含细节和边缘对应高频。用一个高通形态的滤波器把低频压下去、高频抬起来就能抑制光照影响同时增强细节。这里要理解一个关键点降低低频不是完全滤除而是一种非均匀放大。如果直接把低频清零图像会丢失整体的明暗关系看起来非常不自然。同态滤波的“高”和“低”是相对的概念靠两个参数控制放大倍数。4.2 同态滤波器的构造参数与MATLAB实现滤波器形式为H(u,v) (γH - γL)·[1 - exp(-c·D²/D0²)] γL参数含义γH 1控制高频放大倍数γL 1控制低频压制D0是截止频率c是过渡陡峭度。典型参数γH2.5、γL0.5、D080、c4但D0要根据图像分辨率调整分辨率翻倍时D0也要相应增加。function out homomorphic_filter(img, gamma_hi, gamma_lo, c, D0) % 灰度图同态滤波彩色图建议先在HSV的V通道上使用 if size(img, 3) 3 img rgb2gray(img); end img im2double(img); % 对数变换进入频域前先压缩动态范围 z log(img 1); [h, w] size(z); % 生成频率坐标注意要中心化 u (0:h-1) - floor(h/2); v (0:w-1) - floor(w/2); [V, U] meshgrid(v, u); D sqrt(U.^2 V.^2); % 构造同态滤波器 Hf (gamma_hi - gamma_lo) * (1 - exp(-c * (D.^2) / (D0^2))) gamma_lo; % FFT、滤波、逆FFT F fftshift(fft2(z)); G F .* Hf; s real(ifft2(ifftshift(G))); % 指数还原再做对比度拉伸 out exp(s) - 1; out mat2gray(out); end代码里有两个细节要注意。fftshift和ifftshift要配对使用不能用错否则频率中心对不上结果会出现“棋盘”状的亮暗条纹。另外逆变换后会有微小的复数残留误差取real是必要的直接用复数结果会报错。4.3 Retinex与同态滤波的适用范围对比对比项RetinexSSR/MSR/MSRCR同态滤波光照估计视角空域卷积局部加权频域滤波全局频率对偏色的影响不做颜色恢复时偏色明显彩色图需转HSV处理参数数量尺度、权重、颜色恢复多组4个参数运算效率多尺度大核卷积较慢两次FFT较快光晕伪影亮暗交界处容易出现过渡平滑光晕少典型场景低照度细节恢复均匀光照校正、老照片提亮光照特别不均匀、阴影边界明显的图Retinex通常效果更直接。同态滤波胜在平滑自然几乎没有光晕但参数不熟时容易把图像搞得发灰发暗。两张图都值得跑一遍挑效果好的用。5. 常见问题排查与调参实测记录5.1 结果发灰暗部还是看不清最常见的原因是没做分位数截断让极端亮暗像素拉低了整体对比度。另一个原因是log偏移量过大r的范围整体偏向负值。解决方式是确保卷积核归一化做1%到99%截断。我测试时还发现尺度调太大会让结果更灰——比如把尺度从15改成30光照估计更平滑细节增强更弱画面对比度就下来了。5.2 边缘发黑、出现halo伪影黑边基本是卷积边界处理不当用replicate能解决。halo伪影则是Retinex的通病尤其在亮暗交界明显的地方会出现一圈光晕。我的经验是尺度不要超过图像短边的四分之一可以明显改善。如果想精雕细节可以用引导滤波替代大尺度高斯卷积做光照估计光晕能大幅减少但代码复杂度会上升不少。5.3 颜色偏紫偏灰MSRCR参数怎么救先检查三通道是否在同一个log偏移下处理。如果每个通道单独拉伸后合并极容易偏色。正确做法是三通道共享同一个拉伸映射或者MSRCR输出后再一起拉伸。颜色恢复因子β调小、α调大可以减少“灰罩子”的感觉。还有一招彩色图先转到HSV只对V通道做RetinexH和S保持不动颜色保真度很好。缺点是这种方法本质是亮度增强对色偏修正帮助不大。5.4 大图跑得慢三种提速方案方案做法优点缺点先缩放再放大把图缩小到256或512处理再插值放大实现最简单细节有损频域实现高斯滤波用fft2构建频域高斯低通替代imfilter速度最快代码改动中等只处理V通道HSV转换后只增强V通道速度接近单通道对颜色校正不生效MATLAB里imfilter虽然会自动选择算法但全尺寸高斯核还是慢。用频域乘法做高斯滤波和imfilter结果几乎一致大图强烈推荐。另外三通道循环体用single类型跑内存和速度都有明显改善double虽然精确但真没必要。5.5 排查速查表整理一个速查表遇到问题直接对号入座。现象可能原因处理建议输出全灰没做分位数截断补1%-99%截断图像整体偏暗log偏移不够或Gamma问题加Gamma0.85边缘黑框卷积边界用零填充改用replicate亮暗交界光晕尺度太大尺度降到短边1/4颜色严重偏紫通道独立拉伸共享拉伸映射颜色灰蒙蒙颜色恢复因子参数不合适β降到20-30运算极慢全尺寸核空域卷积频域实现或先缩图这篇从Retinex的理论模型讲到MATLAB实现再对比了同态滤波核心目的就是让大家在遇到光照不均时不再只会套直方图均衡化。我做了这么多组对比最深的体会是Retinex不是万能的但对低照度、逆光这类图确实比传统方法强一个量级代价是参数调起来确实麻烦。如果你只是为了快速看效果建议先把SSR跑通再上MSRCR如果项目里要批量处理同态滤波的稳定性和速度反而更适合做成标准流程。最后分享一个小技巧Retinex输出后再做一次Gamma校正通常把gamma设到0.8到0.9高光细节能保下来不少这是我在处理过曝区域时屡试不爽的收尾步骤。