论文：基于反照系数分解的Fattal去雾算法

Fattal 去雾算法是一类经典的物理模型去雾方法。它从大气散射模型出发，把图像颜色、光照和透射率之间的关系拆开，再利用局部统计假设估计透射率。

📚 参考资料

• 《图像画质算法与底层视觉技术》

✍️ 前言

🌫️ 有雾图像的退化模型

单幅图像去雾通常从大气散射模型开始：

eq 1

$$\begin{array}{}\text{(1)}& I(x)=J(x)t(x)+A[1-t(x)]\end{array}$$

其中：

• $I(x)$ ：相机拍摄到的有雾图像。
• $J(x)$ ：场景辐射（Scene Radiance），也就是希望恢复的无雾清晰图像。
• $t(x)$ ：透射率（Medium Transmission），描述场景光线到达相机时保留下来的比例。
• $A$ ：全局大气光（Global Atmospheric Light）。

如 公式(1) 所示，去雾问题中已知的是 $I(x)$ ，需要求解 $J(x)$ 、 $A$ 和 $t(x)$ 。未知量多于约束量，所以这是一个欠定问题，必须引入合适的先验（prior）才能继续求解。这也是基于物理模型的去雾算法最核心的一步。

透射率 $t(x)$ 与场景深度有关：

$$t(x)=e^{-\beta d(x)}$$

这里 $d(x)$ 表示该点的场景深度， $\beta$ 表示大气散射系数。直观地看，距离越远、雾越浓， $t(x)$ 越小，图像中保留下来的真实场景信息也越少。

---

🧩 基于反照系数分解的 Fattal 去雾算法

该方法来自 Raanan Fattal 的论文 Single Image Dehazing。它的关键想法是把场景辐射 $J(x)$ 拆成反射分量和光照分量：

$$J(x)=R(x)\cdot l(x)$$

这个分解与 Retinex 理论一致，可以把图像形成过程理解成两个部分：

• $R(x)$ ：反射分量（Reflectance）。它代表物体表面的固有属性，比如颜色、纹理、材质和细节，可以理解为物体本身的颜色。
• $l(x)$ ：光照分量（Illumination）。它代表照射到场景中的环境光分布，主要影响整体明暗、亮度和对比度，通常变化更平滑。

如果只看一个局部区域 $\mathrm{\Omega }$ ，可以进一步做两个近似：

• 在同一小区域内， $R(x)$ 可以近似为三维常量 $R$ 。
• $t(x)$ 主要由场景深度和雾的浓度决定， $l(x)$ 主要由光照和表面反射性质决定，因此二者在局部区域内可以看作统计无关。

于是有：

$${\mathrm{Cov}}_{\mathrm{\Omega }}(l,t)=0$$

协方差定义为：

$${\mathrm{Cov}}_{\mathrm{\Omega }}(X,Y)=E[X\cdot Y]-E[X]E[Y]=E[(X-E[X])\cdot (Y-E[Y])]$$

把 $J(x)=R\cdot l(x)$ 代入 公式(1)，可以得到：

eq 2

$$\begin{array}{}\text{(2)}& I(x)=t(x)l(x)R+A[1-t(x)]\end{array}$$

🔍 分解反射分量 $R$

接下来，将 $R$ 分解为与大气光 $A$ 平行和垂直的两个部分：

$$\begin{array}{rl}R& =R^{{}^{\prime }}+R_A\\ R_A& =\frac{⟨R,A⟩}{\parallel A\parallel }\cdot \frac{A}{\parallel A\parallel }\end{array}$$

其中：

• $R^{{}^{\prime }}$ ：与大气光正交的分量。
• $R_A$ ：与大气光平行的分量。
• 二者满足 $R^{{}^{\prime }}\perp R_A$ 。

代入 公式(2) 并整理，可以得到：

eq 3

$$\begin{array}{}\text{(3)}& I(x)=t(x)l^{{}^{\prime }}(x)(\frac{R^{{}^{\prime }}}{\parallel R^{\prime }\parallel }+\eta \frac{A}{\parallel A\parallel })+A[1-t(x)]\end{array}$$

这里引入了两个记号：

• $l^{{}^{\prime }}(x)=l(x)\parallel R^{{}^{\prime }}\parallel$
• $\eta =\frac{⟨R,A⟩}{\parallel R^{{}^{\prime }}\parallel \cdot \parallel A\parallel }$

---

🧮 分解观测图像 $I(x)$

同样地，把 $I(x)$ 也投影到与 $A$ 平行和垂直的方向上：

eq 4

$$\begin{array}{}\text{(4)}& \begin{array}{rl}{I}_A(x)& =\frac{⟨I(x),A⟩}{\parallel A\parallel }\\ & =t(x)l^{{}^{\prime }}(x)\eta +[1-t(x)]\parallel A\parallel \\ I_{R^{{}^{\prime }}}(x)& =\frac{⟨I(x),R^{{}^{\prime }}⟩}{\parallel R^{{}^{\prime }}\parallel }\\ & =t(x)l^{{}^{\prime }}(x)\\ & =\sqrt{\parallel I(x){\parallel }^2-I_A(x{)}^2}\end{array}\end{array}$$

这里用到了两个关系：

• $R^{\prime }\cdot A=0$
• $A\cdot A={\parallel A\parallel }^2$

将 $I_{R^{{}^{\prime }}}$ 代入 $I_A(x)$ ，有：

$$I_A(x)=I_{R^{{}^{\prime }}}\eta +[1-t(x)]\parallel A\parallel$$

移项整理后得到透射率表达式：

eq 5

$$\begin{array}{}\text{(5)}& t(x)=1-\frac{I_A(x)-\eta I_{R^{{}^{\prime }}}(x)}{\parallel A\parallel }\end{array}$$

到这里可以看到：如果已经有大气光 $A$ 的估计结果，那么根据 公式(4) 可以分解得到 $I_A(x)$ 和 $I_{R^{{}^{\prime }}}(x)$ 。要求 $t(x)$ ，现在只剩下 $\eta$ 这个未知量。

---

🧠 引入 $h$ 消去 $t(x)$

前面用到的关键假设是 ${\mathrm{Cov}}_{\mathrm{\Omega }}(l,t)=0$ 。为了利用这个假设，可以构造一个与 $t(x)$ 无关的变量 $h$ ，再通过协方差把 $t(x)$ 项消掉，只留下 $\eta$ 。

把 公式(5) 中的 $\parallel A\parallel$ 记为 $a$ ，移项可得：

eq 6

$$\begin{array}{}\text{(6)}& I_A=a-at+\eta I_{R^{{}^{\prime }}}\end{array}$$

等式两边同时对 $h$ 计算协方差：

$$\begin{array}{rl}{\mathrm{Cov}}_{\mathrm{\Omega }}(I_A,h)& =\eta {\mathrm{Cov}}_{\mathrm{\Omega }}(I_{R^{{}^{\prime }}},h)\\ \eta & =\frac{{\mathrm{Cov}}_{\mathrm{\Omega }}(I_A,h)}{{\mathrm{Cov}}_{\mathrm{\Omega }}(I_{R^{{}^{\prime }}},h)}\end{array}$$

其中：

• $\mathrm{Cov}(a,h)=0$ ，因为常数与任意变量的协方差为 0。
• $\mathrm{Cov}(t,h)=0$ ，因为这里构造的 $h$ 与 $t$ 无关。

---

🧭 $h$ 如何取？

由 公式(6) 将 $t$ 移项：

$$at=a-I_A+\eta I_{R^{{}^{\prime }}}$$

等式两边同除 $I_{R^{{}^{\prime }}}$ ：

$$\frac{a}{l^{\prime }}-\eta =\frac{at}{I_{R^{{}^{\prime }}}}-\eta =\frac{a-I_A}{I_{R^{{}^{\prime }}}}$$

因此可以令：

$$h=\frac{a-I_A}{I_{R^{{}^{\prime }}}}=\frac{\parallel A\parallel -I_A(x)}{I_{R^{{}^{\prime }}}(x)}$$

这样 $h$ 与 $t$ 无关，同时也能表达 $l^{\prime }$ 的相关信息。接下来就可以用局部协方差估计 $\eta$ 。

---

🧾 最后求解 $J(x)$

整体流程可以概括为：

1. 给定或估计大气光 $A$ 。
2. 根据 公式(4) 分解得到 $I_A$ 和 $I_R$ 。
3. 构造 $h$ ，再根据 公式(6) 估计 $\eta$ 。
4. 将 $\eta$ 代入 公式(5)，得到透射率 $t(x)$ 。
5. 最后把 $t(x)$ 带回 公式(1)，恢复无雾图像 $J(x)$ 。

---

📌 总结

✅ 优点

优点	说明
物理意义明确	推导基于大气散射模型，不是单纯的经验增强。颜色线假设在自然场景中有一定统计合理性。
色彩还原自然	相比早期直方图拉伸、对比度增强等方法，颜色线模型更容易保留物体固有颜色，色偏相对较小。
推导较紧凑	通过局部统计关系估计透射率，不需要像部分变分方法那样手动调大量平滑参数。

⚠️ 局限

局限	说明
缺乏鲁棒的大气光自动估计机制	算法本身不擅长自动解耦 $A$ 和 $t$ ，实际效果容易依赖人工指定或启发式估计。
天空和高亮区域容易出问题	天空、高亮反射区域不一定满足颜色线假设，去雾后可能出现变暗、断层或颜色异常。
对噪声和压缩敏感	JPEG 压缩块、传感器噪声会破坏局部颜色分布，导致协方差估计不稳定。

个人理解上，Fattal 方法的价值不在于工程上足够鲁棒，而在于它提供了一个很清晰的视角：去雾不只是增强对比度，而是在估计图像形成过程里的透射率和大气光。

---

💻 Python 程序

import cv2
import matplotlib.pyplot as plt
import numpy as np


def estimate_trans(I, A, tmin=0.4):
    A = np.expand_dims(A, axis=[1])  # [3, 1]
    norm_A = np.linalg.norm(A)

    IA = np.dot(I, A) / norm_A  # [N, 1]
    norm_I2 = np.linalg.norm(I, axis=1, keepdims=True) ** 2
    IR = norm_I2 - IA**2
    IR = np.maximum(IR, 0)
    IR = np.sqrt(IR)

    h = (norm_A - IA) / (IR + 1e-8)
    cov_IA = np.cov(IA.flatten(), h.flatten())
    cov_IR = np.cov(IR.flatten(), h.flatten())
    eta = (cov_IA / cov_IR)[0, 1]

    t_map = 1 - (IA - eta * IR) / norm_A
    t_map = (t_map - t_map.min()) / (t_map.max() - t_map.min())
    t_map = np.clip(t_map, a_min=tmin, a_max=1.0)
    return t_map


def dehaze_fattal(img, atmo, tmin=0.3):
    H, W, C = img.shape
    I = img.reshape((H * W, C))

    t_map = estimate_trans(I, atmo, tmin)
    J = (I - (1 - t_map) * np.expand_dims(atmo, axis=[0])) / (t_map + 1e-16)

    dehazed = J.reshape((H, W, C))
    t_map = t_map.reshape((H, W))
    return dehazed, t_map


if __name__ == "__main__":
    img_path = "Imgs/Fog1 (102).png"
    img = cv2.imread(img_path)[:, :, ::-1] / 255.0
    atmo = np.array([210.0, 217.0, 233.0]) / 255.0

    out, tmap = dehaze_fattal(img, atmo)
    out = np.clip(out, 0, 1)

    fig, axes = plt.subplots(1, 3, figsize=(15, 5))

    axes[0].imshow(img)
    axes[0].set_title("hazy input")
    axes[0].axis("off")

    axes[1].imshow(out)
    axes[1].set_title("fattal dehazy output")
    axes[1].axis("off")

    axes[2].imshow(tmap, cmap="gray")
    axes[2].set_title("transmission map")
    axes[2].axis("off")

    plt.tight_layout()
    plt.savefig("comparison_result.png", dpi=300)