Image Dehazing

Image dehazing, also called haze removal, recovers a clear image of a scene from a photograph degraded by haze, fog or smoke. Light from a scene point is partly scattered away on its way to the camera, and light scattered by the particles in the air adds a veil of its own. Both effects grow with distance, so a hazy photograph loses contrast and colour fidelity unevenly: near objects look almost clear and far ones fade into a grey or bluish wash. Undoing this from a single image is underdetermined, since the amount of haze at each pixel depends on a depth the camera did not measure. Most current methods, according to Lisani and Hessel (2024), build on the dark channel prior of He, Sun and Tang (CVPR 2009; journal version in TPAMI, 2011), an observation about haze-free outdoor images that turns the problem into a few local minimum filters and one division. The article builds on image formation, which describes the photometric chain from scene to pixel, and on image filtering, whose minimum and box filters dehazing uses throughout.

Problem Definition

He et al. call the standard model the haze imaging model. For a pixel $\mathbf{x}$ with colour channels $c \in \{r, g, b\}$,

$$\mathbf{I}(\mathbf{x}) = \mathbf{J}(\mathbf{x})\, t(\mathbf{x}) + \mathbf{A}\,\bigl(1 - t(\mathbf{x})\bigr),$$

where $\mathbf{I}$ is the observed hazy image, $\mathbf{J}$ the scene radiance (the haze-free image), $\mathbf{A}$ the global atmospheric light and $t \in (0, 1]$ the medium transmission, the fraction of light from the scene point that reaches the camera without being scattered (He et al., 2011, Eq. 1). The first term is the direct attenuation, a multiplicative loss; the second is the airlight, an additive veil, a term He et al. trace to Koschmieder (1924). In a homogeneous atmosphere

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

with $\beta$ the scattering coefficient and $d$ the scene depth (He et al., 2011, Eq. 2), so a recovered transmission map is also a depth map up to the unknown scale $\beta$.

Lisani and Hessel (2024, §2) derive the same equation from a three-term model of light in a scattering medium: direct light, forward scatter and backscatter. They neglect forward scatter, which is small compared with backscatter, and write the backscatter as $\mathbf{A}(1 - t)$, which is the airlight in He et al.’s notation. In their underwater form the exponent uses the sum of absorption and scattering coefficients of water. For fog and dense haze the coefficient can be treated as the same for all wavelengths (Lisani and Hessel, citing Narasimhan and Nayar, 2002).

Dehazing estimates $\mathbf{A}$ and $t$ and then inverts the model:

$$\mathbf{J}(\mathbf{x}) = \frac{\mathbf{I}(\mathbf{x}) - \mathbf{A}}{\max\bigl(t(\mathbf{x}), t_0\bigr)} + \mathbf{A}.$$

The lower bound $t_0$, for which He et al. give 0.1 as a typical value, keeps the division from amplifying noise where the haze is dense; it leaves some haze where the haze is densest (He et al., 2011, Eq. 22). Geometrically, the model says that in RGB space $\mathbf{I}(\mathbf{x})$ lies on the segment between $\mathbf{J}(\mathbf{x})$ and $\mathbf{A}$, and $t$ is the ratio of the two lengths (He et al., 2011, Eq. 3).

The problem is ill-posed. For an $N$-pixel colour image the model gives $3N$ equations for $3N$ values of $\mathbf{J}$, $N$ values of $t$ and the three components of $\mathbf{A}$, that is $4N + 3$ unknowns (He et al., 2011, §2). Every single-image method therefore adds a prior.

Input and Output

The input is one RGB image of a hazy scene. The dark channel prior also works on gray-scale images, with the minimum over channels dropped, when the scene has enough shadows (He et al., 2011, §5). The outputs are the dehazed image $\mathbf{J}$, the transmission map $t$, the atmospheric light $\mathbf{A}$ and, by $d = -\ln t / \beta$, a depth map up to scale.

Two practical points come from the sources. He et al. note that the scene radiance is usually darker than the atmospheric light, so the dehazed image looks dim and they raise its exposure for display (2011, §4.4). And the model is linear in light, while stored photographs pass through a camera’s tone curve; the dark channel prior’s statistics were gathered on ordinary photographs, and IPOL’s implementation applies the method to 8-bit values as given (Lisani and Hessel, 2024, Algorithm 4 divides by 255). Gamma correction covers the encoding.

Variants

• Single-image dehazing is the subject here: one photograph, no extra information.
• Multi-image and extra-information methods. Before the dark channel prior, He et al. (2011, §1) group earlier work by what it adds: two or more images with different polarization, several images of the same scene under different weather, and depth from user input or known 3D models, such as Kopf et al.’s Deep Photo, which used 3D models and textures of the scene.
• Underwater images. Light in water is absorbed strongly and unevenly across wavelengths. Lisani and Hessel (2024, §5) review the underwater dark channel prior of Drews et al. (2013), which leaves the unreliable red channel out of the dark channel. In their experiments it usually gives better transmission maps than the original prior, but bright objects such as sand still break it, and the dehazed results show colour artefacts that a better atmospheric-light estimate and a smaller $\omega$ reduce.
• Veiling glare. He et al. (2011, Eq. 13) note that the prior estimates only the additive term, so a model with separate transmissions in the two terms, $\mathbf{I} = \mathbf{J}t_1 + \mathbf{A}(1 - t_2)$, can still be split; the additive term is the “veiling luminance” of vision research.
• Atmospheric correction of satellite data is a related physical correction with its own models, described in multispectral satellite imaging. He et al. (2011, §3) say the prior was partly inspired by the dark-object subtraction technique of multispectral remote sensing, which removes spatially uniform haze by subtracting the value of the darkest object.

Challenges

• Unknown depth. Haze varies with depth, so a single global correction cannot remove it. He et al. (2011, Fig. 17) compare global operators: an auto curve, equivalent to dark-object subtraction per channel, can only remove a thin layer of haze matching the nearest objects, and histogram equalization and unsharp masking cannot tell how much to enhance each region, so distant objects stay hazy.
• Ambiguity between haze and colour. A white wall nearby and a distant grey hill can produce the same pixel colour. Any prior that resolves this can be wrong for some surfaces.
• The atmospheric light. It has to be found in the image, usually in the most haze-opaque region, but the brightest pixels may belong to a white car or building (He et al., 2011, §4.3).
• Noise and compression. The inversion divides by $t$, so where $t$ is small it multiplies noise, banding and JPEG artefacts, which Lisani and Hessel (2024, Fig. 9) show in a dehazed sky.
• No ground truth for real haze. A real hazy photograph has no clear twin taken at the same moment, so quantitative evaluation relies on haze synthesized from depth.

Approaches

Earlier single-image priors

He et al. (2011, §2) describe the two single-image methods that preceded theirs, which we have not read. Tan (2008) observed that a haze-free image has more contrast than the hazy one and maximized local contrast under a Markov random field; He et al. report that the results unveil structure but tend to oversaturate and can show halos at depth discontinuities. Fattal (2008) assumed that within a patch of constant albedo the transmission and the surface shading are statistically uncorrelated and estimated them with independent component analysis; He et al. describe it as physically sound but unreliable in dense haze, where colours are faint and variance is low.

The dark channel prior

He et al. define the dark channel of an image $\mathbf{J}$ as

$$J^{\text{dark}}(\mathbf{x}) = \min_{\mathbf{y} \in \Omega(\mathbf{x})} \Bigl( \min_{c} J^c(\mathbf{y}) \Bigr),$$

a per-pixel minimum over the colour channels followed by a minimum filter over a patch $\Omega(\mathbf{x})$, a rank filter in the sense of image filtering (He et al., 2011, Eq. 5). The prior is that, for an outdoor haze-free image outside the sky, $J^{\text{dark}}$ is low and tends to zero. Its causes, in their account, are shadows (of cars, buildings, leaves, rocks), colourful surfaces that are dark in at least one channel (green plants, red or yellow flowers, blue water), and dark objects such as tree trunks.

To test it, they collected outdoor images from Flickr and other search engines, kept haze-free daytime landscapes and cityscapes, picked 5,000 at random, cut out the sky by hand, resized them to a longest side of 500 pixels and computed dark channels with a 15 × 15 patch. About 75% of the dark-channel pixels were zero and 90% were below 25 (on a 0–255 scale) (He et al., 2011, §3, Fig. 5).

The estimate of $t$ follows in four steps (He et al., 2011, Eqs. 7–12). Divide the model by $A^c$ in each channel:

$$\frac{I^c(\mathbf{x})}{A^c} = t(\mathbf{x}) \frac{J^c(\mathbf{x})}{A^c} + 1 - t(\mathbf{x}).$$

Assume $t$ is constant, $\tilde{t}(\mathbf{x})$, over the patch, and take the minimum over the patch and the channels on both sides; $\tilde{t}$ comes out of the minimum:

$$\min_{\mathbf{y} \in \Omega(\mathbf{x})} \min_c \frac{I^c(\mathbf{y})}{A^c} = \tilde{t}(\mathbf{x}) \min_{\mathbf{y} \in \Omega(\mathbf{x})} \min_c \frac{J^c(\mathbf{y})}{A^c} + 1 - \tilde{t}(\mathbf{x}).$$

By the prior the minimum of $J^c / A^c$ is zero, which removes the multiplicative term:

$$\tilde{t}(\mathbf{x}) = 1 - \omega \min_{\mathbf{y} \in \Omega(\mathbf{x})} \min_c \frac{I^c(\mathbf{y})}{A^c},$$

where He et al. add the constant $0 < \omega \le 1$, fixed at 0.95 in their results, to keep a little haze on distant objects: aerial perspective is a depth cue, and an image with all haze removed can look unnatural (2011, Eq. 12). The estimate is the dark channel of the hazy image normalized by $\mathbf{A}$. In the sky, where the prior fails, the hazy colour is close to $\mathbf{A}$, the normalized minimum tends to 1 and $\tilde{t}$ to 0, which is correct for a region at infinite distance, so the sky needs no separate treatment. He et al. stress that, unlike Tan’s and Fattal’s methods, theirs uses only the additive airlight term.

Keeping $D_J(\mathbf{x}) = \min_{\Omega, c} J^c / A^c$ instead of setting it to zero gives, where $t$ is constant over the patch and $\mathbf{A}$ is exact,

$$\tilde{t} = 1 - \omega + \omega\, t\, \bigl(1 - D_J\bigr),$$

pixel by pixel. This is our step, a direct consequence of the equation above; the example below checks it. With $\omega = 1$ the estimate is $t(1 - D_J)$: exact where the prior holds, and too low by $t D_J$ where it does not. With $\omega = 0.95$ and $D_J = 0$ it is too high by $0.05(1 - t)$, the haze He et al. keep on purpose.

Patch size. Larger patches make the prior more likely to hold, since a patch is more likely to contain a dark pixel; with 3 × 3 patches He et al. find grayish surfaces oversaturated. But larger patches violate the constant-$t$ assumption more often, and halos at depth edges grow. They use 15 × 15 for 600 × 400 images (He et al., 2011, §4.5). Lisani and Hessel (2024, §4) add that larger patches reduce the over-correction of whitish objects but also weaken the dehazing of haze trapped between near objects such as tree branches.

Estimating the atmospheric light

The colour of the most haze-opaque region is a natural estimate of $\mathbf{A}$, but the brightest pixel of the image is the most haze-opaque only under overcast sky; with sunlight, a white car or building can be brighter than the atmospheric light (He et al., 2011, §4.3). He et al. therefore use the dark channel of the hazy image, which approximates haze density: they take the 0.1% of pixels with the brightest dark channel and, among them, select “the pixels with highest intensity in the input image”. Lisani and Hessel (2024, §3.1) read this as the single pixel $\mathbf{x}^*$ of highest intensity, with $\mathbf{A} = \mathbf{I}(\mathbf{x}^*)$. Their implementation follows a variant from a third-party repository, He-Zhang/image_dehaze: $\mathbf{A}$ is the average of the 0.1% brightest dark-channel pixels (Algorithm 1), optionally excluding saturated pixels (§5). The two rules give different results in practice, so an implementation should say which it uses; the example below prints both. He et al. do not define “intensity”; the example uses the mean of the three channels.

A wrong $\mathbf{A}$ shifts colours. Lisani and Hessel (2024, Fig. 5) show an image whose dark channel is brightest on a sunlit building near the camera: $\mathbf{A}$ takes on a yellow tint and the dehazed image a blue cast. With $\mathbf{A}$ taken from a hand-picked region the cast disappears, but the image is over-enhanced, because that $\mathbf{A}$ is not computed from the highest values in the image.

Refining the transmission: soft matting and the guided filter

The patch-wise estimate is blocky and too high near depth edges: a patch that straddles a near object and the background takes its minimum from the near object’s dark pixels. He et al. (2011, §4.2) observe that the haze model has the form of the image matting equation, $\mathbf{I} = \mathbf{F}\alpha + \mathbf{B}(1 - \alpha)$, with $t$ in the role of the alpha matte. They refine $\tilde{t}$ by minimizing

$$E(t) = t^{\top} L\, t + \lambda (t - \tilde{t})^{\top} (t - \tilde{t}),$$

where $L$ is the matting Laplacian of Levin, Lischinski and Weiss and $\lambda = 10^{-4}$ (He et al., 2011, Eqs. 15–17). Solving $(L + \lambda U) t = \lambda \tilde{t}$ with preconditioned conjugate gradients took about 10–20 seconds for a 600 × 400 image on a 3.0 GHz Pentium 4, and a bilateral filter then smooths small-scale texture (§§4.2, 5).

The guided filter of the same authors replaces the linear solve with an explicit filter (He et al., 2013). With a guidance image $I$ (for dehazing, the hazy image or its intensity) and an input $p$ (the raw transmission), it assumes that the output is a linear function of the guide in every window $\omega_k$ of radius $r$:

$$q_i = a_k I_i + b_k, \qquad i \in \omega_k,$$

and fits $a_k, b_k$ by ridge regression of $p$ on $I$ in the window, $\sum_{i \in \omega_k} \bigl( (a_k I_i + b_k - p_i)^2 + \epsilon a_k^2 \bigr)$:

$$a_k = \frac{\frac{1}{|\omega|} \sum_{i \in \omega_k} I_i p_i - \mu_k \bar{p}_k}{\sigma_k^2 + \epsilon}, \qquad b_k = \bar{p}_k - a_k \mu_k,$$

where $\mu_k$ and $\sigma_k^2$ are the mean and variance of $I$ in the window and $\bar{p}_k$ the mean of $p$ (He et al., 2013, Eqs. 6–8). Each pixel lies in many windows, so the output averages their predictions: $q_i = \bar{a}_i I_i + \bar{b}_i$, with $\bar{a}_i, \bar{b}_i$ the means of $a_k, b_k$ over the windows that contain $i$ (Eqs. 9–10). Every step is a box mean, which an integral image computes in constant time per pixel, so the filter is $O(N)$ regardless of radius (He et al., 2013, Algorithm 1 and §4). Lisani and Hessel (2024, Algorithms 2–4) write out the integral-image implementation with a mirror boundary.

Two properties explain its use here. Where $p$ is constant over a window, the covariance is zero, $a_k = 0$ and $b_k = \bar{p}_k$, so a constant transmission passes through unchanged. Where $p$ is a linear function of $I$ and $\epsilon = 0$, the regression recovers it exactly. Between those cases, $\epsilon$ sets which edges of the guide are transferred: a guide variance well above $\epsilon$ makes $a_k$ the plain regression slope, and one well below makes the window average $p$. In the guided-filter paper, He et al. (2013, §3.7) show that the filter is one Jacobi step of the matting-Laplacian system with a particular weighting, which is why it can stand in for soft matting; their Fig. 16 refines a dark-channel transmission with $r = 20$ and $\epsilon = 10^{-3}$. IPOL’s defaults are a 61 × 61 window ($r = 30$) and $\epsilon = 10^{-4}$, with the intensity channel (the RGB average) as guide (Lisani and Hessel, 2024, §3.5). The bilateral filter article compares the guided filter with bilateral filtering as an edge-preserving smoother; here it serves to transfer the hazy image’s edges into $t$.

Lisani and Hessel (2024, §4.2) report how the guided filter’s parameters trade off: a larger radius makes halos less visible but weakens the dehazing, and a larger $\epsilon$ blurs the transmission map, lowering its highest values and raising its lowest, which strengthens the dehazing (the caption of their Fig. 6 says instead that the map becomes more blurry as $\epsilon$ decreases). Their profiles across a tree trunk show that the halo, a gradient on a region of nearly constant depth, survives the default refinement.

Two slips in IPOL’s printed algorithms matter to anyone implementing from them. Algorithm 4, line 11, averages $a$ where it should average $b$ (the text’s Eq. 19 averages $b$); the example below averages $b$. And Algorithm 5, line 6, divides by $\max(t, t_0)$ with the unrefined $t$, although line 5 has just computed the refined $t'$; in our reading the refined map is meant, but the example prints the recovery with each map.

The full algorithm and its costs

Put together, IPOL’s Algorithm 5 computes the dark channel of $\mathbf{I}$, estimates $\mathbf{A}$, computes the raw transmission with $\omega$, refines it with the guided filter and recovers $\mathbf{J}$ with $t_0$. The defaults are a 15 × 15 patch, $p = 0.1\%$ of pixels for $\mathbf{A}$, $\omega = 0.95$, a 61 × 61 guided-filter window, $\epsilon = 10^{-4}$ and $t_0 = 0.1$ (Lisani and Hessel, 2024, §3.5). The minimum filters can be computed in time linear in the image size with van Herk’s algorithm, as He et al. (2011, §5) did; IPOL’s implementation computes them directly, at a cost proportional to the patch area (§3.7).

Failure modes

The sources agree on where the dark channel prior fails:

• Bright objects without shadows. Where objects are inherently similar to the atmospheric light and no shadow falls on them (snowy ground, a white wall, white marble), the prior is invalid, the transmission is underestimated and the haze overestimated (He et al., 2011, §6, Fig. 18). Lisani and Hessel (2024, §4.1, Fig. 4) show a sunlit grey wall whose transmission comes out below the sky’s.
• Halos near objects at contrasting depths, where dehazing is weaker around foreground objects (Lisani and Hessel, 2024, §4.2, Figs. 7–8).
• Wrong atmospheric light, which gives colour casts (above).
• Model failures. He et al. (2011, Fig. 19) show the atmospheric light varying across the sky (bright on one side, dim on the other) and wavelength-dependent transmission in thin haze over kilometres, which leaves the horizon bluish; their model assumes a constant $\mathbf{A}$ and one $t$ for all channels.
• Over-enhancement of noise, banding and JPEG blocks in the sky, reduced by a smaller $\omega$ or a larger $t_0$ (Lisani and Hessel, 2024, Fig. 9).

Other priors: haze-lines

Berman, Treibitz and Avidan (CVPR 2016) replaced patches with a non-local prior. They assume that the colours of a haze-free image are well approximated by a few hundred distinct colours forming tight clusters in RGB space (in their analysis, K-means to 500 colours), and that the pixels of a cluster are spread over the image at different depths. Under haze, each cluster becomes a haze-line in RGB space passing through $\mathbf{A}$, and a pixel’s distance from $\mathbf{A}$ along its line measures its transmission. Assuming that the pixel farthest from $\mathbf{A}$ on each line is haze-free gives $\tilde{t}(\mathbf{x}) = r(\mathbf{x}) / \hat{r}_{\max}(\mathbf{x})$; they then impose the lower bound $t_{LB} = 1 - \min_c I^c / A^c$, which follows from $J \ge 0$ and which they note is, eroded, the basis of the dark channel estimate, and regularize. The algorithm is linear in the image size, needs no training, and takes $\mathbf{A}$ from a separate estimator. They note that it may fail when the airlight is much brighter than the scene, since most pixels then point in the same direction from $\mathbf{A}$.

Learned dehazing

AOD-Net and DCPDN train on hazy images synthesized by the haze model from clear images with depth.

• AOD-Net (Li, Peng, Wang, Xu and Feng, ICCV 2017) rewrites the inversion so that $1/t$ and $\mathbf{A}$ enter one input-dependent variable, $\mathbf{J} = K(\mathbf{x}) \mathbf{I}(\mathbf{x}) - K(\mathbf{x}) + b$ with a constant bias $b$ (default 1), and estimates $K$ with a five-layer convolutional network with three filters per layer, trained under a mean squared error loss. The authors argue that estimating $t$ and $\mathbf{A}$ separately lets their errors compound. They synthesized training data from the indoor NYU2 depth dataset with each channel of $\mathbf{A}$ uniform in [0.6, 1.0] and $\beta \in \{0.4, 0.6, \dots, 1.6\}$, and tested on 3,170 NYU2 images (TestSet A) and 800 synthetic images from the Middlebury stereo database (TestSet B). On TestSet B they report an average PSNR of 21.5412 dB and SSIM of 0.9272 for AOD-Net against 18.9781 dB and 0.8584 for the dark channel prior (their Table 2). They also placed AOD-Net in front of Faster R-CNN to study detection on hazy images.
• DCPDN (Zhang and Patel, CVPR 2018) estimates the transmission with a densely connected encoder–decoder with multi-level pyramid pooling and an edge-preserving loss, estimates $\mathbf{A}$ with a U-Net, embeds the inversion in the network, and adds a joint discriminator, in the manner of a generative adversarial network, that judges pairs of transmission map and dehazed image. They trained on 4,000 images synthesized from 1,000 NYU-depth2 images with four atmospheric-light conditions ($\mathbf{A} \in [0.5, 1]$, $\beta \in [0.4, 1.6]$) and tested on 400 NYU-depth2 images (TestA) and 200 images from Middlebury (40) and Sun3D (160) (TestB).
• All-in-one restoration. Later networks treat haze as one of several degradations. AirNet (Li, Liu, Hu, Wu, Lv and Peng, CVPR 2022) handles noise, rain and haze with one network, training its dehazing on the Outdoor Training Set of RESIDE (72,135 hazy–clean pairs) and testing on its Synthetic Objective Testing Set (500 outdoor pairs), as the AirNet authors describe RESIDE (Li et al., 2019, not read).

In our reading, the shift from priors to networks follows the same logic as elsewhere in restoration: a hand-made prior such as the dark channel is one rule that fails on whole classes of scenes, while a network fitted to synthetic pairs learns whatever statistics the training set contains, including its biases (the synthetic haze, the depth source, the range of $\mathbf{A}$ and $\beta$).

Practical Example: the Dark Channel Prior on Known Haze

The example synthesizes haze with known transmission and atmospheric light on two scenes and runs the dark channel prior against the truth.

• A synthetic scene where the prior holds exactly, 400 × 600 pixels. Its depth is piecewise simple: a sky at 12 m, a wall at 6 m, a ground plane rising from 2 m to just under 6 m, three blocks at constant depths and a white object at 3 m. Its non-sky surfaces are coloured tiles with a black pixel every 6 pixels in each direction, so every 15 × 15 patch over them contains a zero; the white object (0.78, 0.80, 0.82) and the sky have no dark pixel. The scene was designed at planning so that the prior holds exactly where intended, so that part illustrates the derivation above rather than testing the prior.
• The bundled Middlebury motorcycle from scikit-image (data.stereo_motorcycle, a down-sampled scene from the 2014 Middlebury stereo datasets of Scharstein et al.), with haze computed from its true disparity. Depth is $Z = fB/(d + \text{doffs})$, as in the disparity article, with $f = 994.978$ px, $B = 0.193001$ m and doffs $= 31.086$ px from the scikit-image 0.25.2 docstring, which gives this calibration for the down-sampled images. The 7.35% of pixels without a disparity are filled with the nearest valid value. It is an indoor scene, and the prior holds only weakly on it.

Both scenes use $\mathbf{A} = (0.85, 0.88, 0.92)$ and $t = e^{-\beta Z}$ with $\beta = 0.2$ and $0.4$ per metre, design values chosen at planning so that the motorcycle’s transmission spans roughly 0.37–0.66 and 0.13–0.43. Images are floats in [0, 1], and the haze model is applied to the stored values, as IPOL’s implementation applies the method to 8-bit values; there is no noise or quantization.

The script compares:

• the oracle inversion with the true $t$ and $\mathbf{A}$;
• the dark channel prior with IPOL’s defaults (15 × 15 patch, $t_0 = 0.1$, a 61 × 61 guided filter with $\epsilon = 10^{-4}$ guided by the RGB mean of the hazy image), with $\omega = 0.95$ and $1$, the raw and the guided-filter-refined transmission, and three atmospheric lights: the true one, IPOL’s mean of the 0.1% brightest dark-channel pixels (“mean”) and He et al.’s brightest of those pixels (“brightest”, with intensity read as the RGB mean);
• a contrast-stretch baseline, a per-channel linear stretch of the 1st and 99th percentiles to [0, 1].

The guided filter follows IPOL’s Algorithm 4 with box means from an integral image (mirror boundary) and with line 11 averaging $b$. Each row recovers $\mathbf{J}$ with its own transmission, so the raw rows also show IPOL’s Algorithm 5 as printed. The measures are the largest channel error of $\hat{\mathbf{A}}$, the mean absolute error and the mean signed error (bias) of $\hat{t}$, and the PSNR of the recovered image (clipped to [0, 1], peak 1). On the synthetic scene the transmission error is also split by region: prior pixels whose true $t$ is constant over their patch (“pr-cst”), other prior pixels, on depth edges and the ground ramp (“pr-oth”), the white object and the sky. “PSNR-ns” leaves out the sky, and “white dJ” is the mean change of the white object in the recovered image.

The script first checks its components. The oracle must reproduce $\mathbf{J}$ to within $10^{-12}$. The identity $\tilde{t} = 1 - \omega + \omega t(1 - D_J)$ must hold to within $10^{-12}$ on patch-constant pixels. The guided filter must return a constant input unchanged (within $10^{-10}$) and a linear function of the guide exactly when $\epsilon = 0$ (within $10^{-9}$). The integral-image box mean must match SciPy’s uniform_filter with the mirror boundary (within $10^{-10}$), and the minimum filter must equal scikit-image’s erosion with a 15 × 15 footprint. Before the run we expected the following:

1. the raw estimate with the true $\mathbf{A}$ and $\omega = 1$ is exact on patch-constant prior pixels, and $\omega = 0.95$ is not;
2. the white object’s transmission is underestimated and the object comes out darker;
3. with the true $\mathbf{A}$ and $\omega = 1$, the error is larger at depth edges and on the ramp than on patch-constant pixels, and the guided filter lowers it there;
4. on the motorcycle the raw transmission is biased low, and $\omega = 0.95$ gives a higher PSNR than $\omega = 1$ for every atmospheric-light rule and transmission;
5. both rules find $\mathbf{A}$ in the synthetic sky to within $\max t \lvert \mathbf{J} - \mathbf{A} \rvert$ over the sky (0.0363 and 0.0033 at the two $\beta$, from the construction), and less accurately on the motorcycle, which has no sky;
6. the guided filter never raises the mean absolute error of $t$;
7. with IPOL’s defaults the dark channel prior beats the hazy input and the stretch on the synthetic scene, and on the motorcycle beats the hazy input at $\beta = 0.4$ but not at $0.2$.

The motorcycle half of the last expectation was a guess from a planning run with a different depth scaling. The script runs in about 2 seconds with NumPy 2.3.5, SciPy 1.16.3 and scikit-image 0.25.2, and was run as planned, without changes after its first run.

# Image dehazing: Koschmieder's haze model, the dark channel prior and the guided filter
# on a synthetic scene where the prior holds exactly (plus a white object and a sky)
# and on scikit-image's bundled Middlebury motorcycle with haze from its true disparity.
import os
for v in ("OPENBLAS_NUM_THREADS", "OMP_NUM_THREADS", "MKL_NUM_THREADS"):
    os.environ[v] = "1"
import numpy as np
import scipy
import skimage
from scipy import ndimage
from skimage import data, morphology

A_TRUE = np.array([0.85, 0.88, 0.92])   # atmospheric light used to make the haze
BETAS = (0.2, 0.4)                      # scattering coefficient, per metre
SR, P, RR, EPS, T0 = 7, 0.001, 30, 1e-4, 0.1   # IPOL defaults: 15x15 patch, 0.1%, 61x61, eps, t0


def synthetic_scene():
    H, W = 400, 600
    rng = np.random.default_rng(0)
    label = np.zeros((H, W), int)              # 0: prior holds, 1: white object, 2: sky
    Z = np.empty((H, W))                       # depth in metres
    Z[:100] = 12.0; label[:100] = 2            # sky
    Z[100:200] = 6.0                           # far wall
    r = np.arange(200, 400)[:, None]
    Z[200:400] = 6.0 - 4.0 * (r - 199) / 200.0 # ground from 6 m to 2 m
    Z[220:330, 40:160] = 2.5                   # three blocks
    Z[150:260, 220:340] = 4.0
    Z[120:190, 470:570] = 5.0
    Z[240:340, 400:500] = 3.0; label[240:340, 400:500] = 1   # white object
    tiles = rng.uniform(0.15, 0.85, (H // 10, W // 10, 3))   # coloured 10x10 tiles
    J = np.kron(tiles, np.ones((10, 10, 1)))
    lattice = np.zeros((H, W), bool); lattice[::6, ::6] = True
    J[lattice & (label == 0)] = 0.0            # a black "shadow" pixel every 6 px
    top, hor = np.array([0.45, 0.60, 0.90]), np.array([0.70, 0.78, 0.90])
    s = (np.arange(100) / 99.0)[:, None, None]
    J[:100] = top + s * (hor - top)            # sky gradient
    J[label == 1] = [0.78, 0.80, 0.82]         # white object: no dark pixel
    return J, Z, label


def motorcycle_scene():
    L, _, disp = data.stereo_motorcycle()
    bad = ~np.isfinite(disp)
    idx = ndimage.distance_transform_edt(bad, return_distances=False, return_indices=True)
    d = disp[tuple(idx)].astype(float)                 # fill holes with the nearest valid disparity
    Z = 994.978 * 0.193001 / (d + 31.086)              # Z = fB/(d + doffs), metres
    return L / 255.0, Z, np.zeros(Z.shape, int), bad.mean()


def box_mean(x, rr):                                   # IPOL Algorithms 2-3: integral image
    p = np.pad(x, rr + 1, mode="reflect")             # dcb|abcd...|...: NumPy 'reflect'
    ii = p.cumsum(0).cumsum(1)
    k = 2 * rr + 1
    s = ii[k:, k:] - ii[:-k, k:] - ii[k:, :-k] + ii[:-k, :-k]
    return s[: x.shape[0], : x.shape[1]] / k ** 2


def guided_filter(I, p, rr, eps):                      # IPOL Algorithm 4 (line 11 averages b)
    mI, mp = box_mean(I, rr), box_mean(p, rr)
    var = box_mean(I * I, rr) - mI * mI
    cov = box_mean(I * p, rr) - mI * mp
    a = cov / (var + eps)
    b = mp - a * mI
    return box_mean(a, rr) * I + box_mean(b, rr)


def patch_min(x):                                      # min over channels, then a 15x15 minimum filter
    return ndimage.minimum_filter(x.min(2), size=2 * SR + 1, mode="nearest")


def atmospheric_light(I, rule):
    dc = patch_min(I).ravel()
    M = int(np.floor(P * dc.size))
    top = np.argsort(dc, kind="stable")[-M:]           # the 0.1% brightest dark-channel pixels
    pix = I.reshape(-1, 3)[top]
    if rule == "mean":                                 # IPOL Algorithm 1: their average
        return pix.mean(0)
    return pix[np.argmax(pix.mean(1))]                 # He et al.: the brightest of them (RGB mean)


def recover(I, A, t):
    return (I - A) / np.maximum(t, T0)[..., None] + A


def psnr(x, ref, mask=None):
    x = np.clip(x, 0, 1)
    if mask is not None:
        x, ref = x[mask], ref[mask]
    return 10 * np.log10(1.0 / np.mean((x - ref) ** 2))


def stretch(I):                                        # baseline: per-channel 1-99% linear stretch
    lo, hi = np.percentile(I.reshape(-1, 3), [1, 99], axis=0)
    return np.clip((I - lo) / (hi - lo), 0, 1)


print("NumPy %s, SciPy %s, scikit-image %s" % (np.__version__, scipy.__version__, skimage.__version__))
Js, Zs, Ls = synthetic_scene()
Jm, Zm, Lm, holes = motorcycle_scene()
scenes = (("synthetic", Js, Zs, Ls), ("motorcycle", Jm, Zm, Lm))

print("\n[data] synthetic 400x600: depth %.1f..%.1f m" % (Zs.min(), Zs.max()))
DJs = ndimage.minimum_filter((Js / A_TRUE).min(2), 2 * SR + 1, mode="nearest")
print("  max of min_patch,c J/A over prior-region pixels: %.1e" % DJs[Ls == 0].max())
print("  same on white object: median %.4f; on sky: median %.4f" % (np.median(DJs[Ls == 1]), np.median(DJs[Ls == 2])))
print("[data] motorcycle 500x741: holes filled %.4f, depth %.3f..%.3f m" % (holes, Zm.min(), Zm.max()))
Dm = patch_min(Jm)
print("  dark channel (15x15) of the clean image: median %.4f, share = 0: %.4f, share <= 25/255: %.4f"
      % (np.median(Dm), (Dm == 0).mean(), (Dm <= 25 / 255).mean()))

for name, J, Z, lab in scenes:
    DJ = ndimage.minimum_filter((J / A_TRUE).min(2), 2 * SR + 1, mode="nearest")
    for beta in BETAS:
        t = np.exp(-beta * Z)
        I = J * t[..., None] + A_TRUE * (1 - t[..., None])     # Koschmieder: I = J t + A (1 - t)
        gI = I.mean(2)                                         # guidance: intensity in [0, 1]
        const = (ndimage.maximum_filter(t, 2 * SR + 1, mode="nearest") ==
                 ndimage.minimum_filter(t, 2 * SR + 1, mode="nearest"))
        print("\n== %s, beta = %.1f /m: t in [%.4f, %.4f]; pixels with t constant over their patch: %d"
              % (name, beta, t.min(), t.max(), const.sum()))
        print("checks: oracle max|J-J*| %.1e | GF const %.1e | GF linear (eps=0) %.1e | box7 %.1e | box30 %.1e | minfilter-erosion %.1e"
              % (np.abs((I - A_TRUE) / t[..., None] + A_TRUE - J).max(),
                 np.abs(guided_filter(gI, np.full(gI.shape, 0.6), RR, EPS) - 0.6).max(),
                 np.abs(guided_filter(gI, 0.3 * gI + 0.2, RR, 0.0) - (0.3 * gI + 0.2)).max(),
                 np.abs(box_mean(gI, SR) - ndimage.uniform_filter(gI, 2 * SR + 1, mode="mirror")).max(),
                 np.abs(box_mean(gI, RR) - ndimage.uniform_filter(gI, 2 * RR + 1, mode="mirror")).max(),
                 np.abs(patch_min(I) - morphology.erosion(I.min(2), morphology.footprint_rectangle((15, 15)))).max()))
        for om in (0.95, 1.0):
            th = 1 - om * patch_min(I / A_TRUE)
            e = np.abs(th - (1 - om + om * t * (1 - DJ)))
            print("exact t-hat = 1 - w + w t (1 - D_J), true A, w = %.2f: max|err| where t constant %s, over all pixels %.1e"
                  % (om, "%.1e" % e[const].max() if const.any() else "n/a", e.max()))
        sky = lab == 2
        print("PSNR hazy input %.2f dB, contrast stretch %.2f dB" % (psnr(I, J), psnr(stretch(I), J)))
        As = {"true": A_TRUE, "mean": atmospheric_light(I, "mean"), "brightest": atmospheric_light(I, "brightest")}
        print("A estimates: mean %s, brightest %s" % (np.round(As["mean"], 4), np.round(As["brightest"], 4)))
        if name == "synthetic":
            print("%-9s %-4s %-6s %7s %7s %7s %7s %7s %7s %7s %7s %7s %8s"
                  % ("A", "w", "t", "|dA|", "MAE t", "bias t", "pr-cst", "pr-oth", "white", "sky", "PSNR", "PSNR-ns", "white dJ"))
        else:
            print("%-9s %-4s %-6s %7s %7s %7s %7s" % ("A", "w", "t", "|dA|", "MAE t", "bias t", "PSNR"))
        for rule, A in As.items():
            for om in (0.95, 1.0):
                traw = 1 - om * patch_min(I / A)
                for tk, tt in (("raw", traw), ("guided", guided_filter(gI, traw, RR, EPS))):
                    Jr = recover(I, A, tt)
                    err = tt - t
                    dA = np.abs(A - A_TRUE).max()
                    if name == "synthetic":
                        m = lambda k: np.abs(err[k]).mean()
                        print("%-9s %-4.2f %-6s %7.4f %7.4f %7.4f %7.4f %7.4f %7.4f %7.4f %7.2f %7.2f %8.4f"
                              % (rule, om, tk, dA, np.abs(err).mean(), err.mean(), m((lab == 0) & const),
                                 m((lab == 0) & ~const), m(lab == 1), m(sky), psnr(Jr, J), psnr(Jr, J, ~sky),
                                 (np.clip(Jr, 0, 1) - J)[lab == 1].mean()))
                    else:
                        print("%-9s %-4.2f %-6s %7.4f %7.4f %7.4f %7.2f"
                              % (rule, om, tk, dA, np.abs(err).mean(), err.mean(), psnr(Jr, J)))

NumPy 2.3.5, SciPy 1.16.3, scikit-image 0.25.2

[data] synthetic 400x600: depth 2.0..12.0 m
  max of min_patch,c J/A over prior-region pixels: 0.0e+00
  same on white object: median 0.8913; on sky: median 0.6349
[data] motorcycle 500x741: holes filled 0.0735, depth 2.110..5.017 m
  dark channel (15x15) of the clean image: median 0.0824, share = 0: 0.0153, share <= 25/255: 0.5425

== synthetic, beta = 0.2 /m: t in [0.0907, 0.6703]; pixels with t constant over their patch: 124144
checks: oracle max|J-J*| 1.2e-15 | GF const 1.2e-12 | GF linear (eps=0) 4.8e-12 | box7 1.6e-12 | box30 2.6e-13 | minfilter-erosion 0.0e+00
exact t-hat = 1 - w + w t (1 - D_J), true A, w = 0.95: max|err| where t constant 1.7e-16, over all pixels 2.7e-01
exact t-hat = 1 - w + w t (1 - D_J), true A, w = 1.00: max|err| where t constant 1.5e-16, over all pixels 2.9e-01
PSNR hazy input 12.38 dB, contrast stretch 13.57 dB
A estimates: mean [0.8348 0.8698 0.9182], brightest [0.8348 0.8698 0.9182]
A         w    t         |dA|   MAE t  bias t  pr-cst  pr-oth   white     sky    PSNR PSNR-ns white dJ
true      0.95 raw     0.0000  0.0478  0.0128  0.0309  0.0407  0.3635  0.0254   22.63   21.54  -0.2652
true      0.95 guided  0.0000  0.0464  0.0128  0.0328  0.0445  0.3300  0.0163   24.57   23.57  -0.1371
true      1.00 raw     0.0000  0.0406 -0.0206  0.0000  0.0148  0.4055  0.0679   22.91   21.81  -0.2900
true      1.00 guided  0.0000  0.0447 -0.0206  0.0204  0.0300  0.3711  0.0416   24.75   23.63  -0.1984
mean      0.95 raw     0.0152  0.0510  0.0080  0.0297  0.0396  0.3657  0.0407   20.97   21.95  -0.2098
mean      0.95 guided  0.0152  0.0474  0.0080  0.0320  0.0438  0.3320  0.0220   21.60   23.08  -0.1042
mean      1.00 raw     0.0152  0.0447 -0.0256  0.0012  0.0141  0.4079  0.0840   21.56   22.99  -0.2255
mean      1.00 guided  0.0152  0.0482 -0.0256  0.0208  0.0297  0.3732  0.0555   22.06   23.99  -0.1548
brightest 0.95 raw     0.0152  0.0510  0.0080  0.0297  0.0396  0.3657  0.0407   20.97   21.95  -0.2098
brightest 0.95 guided  0.0152  0.0474  0.0080  0.0320  0.0438  0.3320  0.0220   21.60   23.08  -0.1042
brightest 1.00 raw     0.0152  0.0447 -0.0256  0.0012  0.0141  0.4079  0.0840   21.56   22.99  -0.2255
brightest 1.00 guided  0.0152  0.0482 -0.0256  0.0208  0.0297  0.3732  0.0555   22.06   23.99  -0.1548

== synthetic, beta = 0.4 /m: t in [0.0082, 0.4493]; pixels with t constant over their patch: 124144
checks: oracle max|J-J*| 9.7e-15 | GF const 2.7e-12 | GF linear (eps=0) 9.6e-11 | box7 1.6e-12 | box30 2.5e-13 | minfilter-erosion 0.0e+00
exact t-hat = 1 - w + w t (1 - D_J), true A, w = 0.95: max|err| where t constant 1.2e-16, over all pixels 2.5e-01
exact t-hat = 1 - w + w t (1 - D_J), true A, w = 1.00: max|err| where t constant 1.1e-16, over all pixels 2.6e-01
PSNR hazy input 9.77 dB, contrast stretch 10.66 dB
A estimates: mean [0.8486 0.8791 0.9198], brightest [0.8486 0.8791 0.9198]
A         w    t         |dA|   MAE t  bias t  pr-cst  pr-oth   white     sky    PSNR PSNR-ns white dJ
true      0.95 raw     0.0000  0.0546  0.0391  0.0420  0.0520  0.1871  0.0500   17.19   18.47  -0.1358
true      0.95 guided  0.0000  0.0543  0.0391  0.0419  0.0520  0.1678  0.0524   16.79   17.77  -0.1124
true      1.00 raw     0.0000  0.0188 -0.0031  0.0000  0.0145  0.2316  0.0103   19.52   24.38  -0.1453
true      1.00 guided  0.0000  0.0251 -0.0031  0.0122  0.0247  0.2134  0.0077   18.90   22.28  -0.1539
mean      0.95 raw     0.0014  0.0542  0.0386  0.0419  0.0519  0.1873  0.0485   16.93   18.32  -0.1295
mean      0.95 guided  0.0014  0.0539  0.0386  0.0417  0.0518  0.1680  0.0509   16.55   17.64  -0.1071
mean      1.00 raw     0.0014  0.0192 -0.0036  0.0002  0.0144  0.2319  0.0118   19.19   24.22  -0.1388
mean      1.00 guided  0.0014  0.0252 -0.0036  0.0122  0.0246  0.2136  0.0081   18.60   22.15  -0.1470
brightest 0.95 raw     0.0014  0.0542  0.0386  0.0419  0.0519  0.1873  0.0485   16.93   18.32  -0.1295
brightest 0.95 guided  0.0014  0.0539  0.0386  0.0417  0.0518  0.1680  0.0509   16.55   17.64  -0.1071
brightest 1.00 raw     0.0014  0.0192 -0.0036  0.0002  0.0144  0.2319  0.0118   19.19   24.22  -0.1388
brightest 1.00 guided  0.0014  0.0252 -0.0036  0.0122  0.0246  0.2136  0.0081   18.60   22.15  -0.1470

== motorcycle, beta = 0.2 /m: t in [0.3666, 0.6557]; pixels with t constant over their patch: 0
checks: oracle max|J-J*| 5.6e-16 | GF const 2.8e-12 | GF linear (eps=0) 1.6e-12 | box7 8.1e-13 | box30 9.9e-14 | minfilter-erosion 0.0e+00
exact t-hat = 1 - w + w t (1 - D_J), true A, w = 0.95: max|err| where t constant n/a, over all pixels 2.1e-01
exact t-hat = 1 - w + w t (1 - D_J), true A, w = 1.00: max|err| where t constant n/a, over all pixels 2.2e-01
PSNR hazy input 11.96 dB, contrast stretch 17.90 dB
A estimates: mean [0.8921 0.8953 0.9196], brightest [0.9049 0.9055 0.9264]
A         w    t         |dA|   MAE t  bias t    PSNR
true      0.95 raw     0.0000  0.1047 -0.0823   17.16
true      0.95 guided  0.0000  0.1007 -0.0823   18.31
true      1.00 raw     0.0000  0.1229 -0.1110   14.86
true      1.00 guided  0.0000  0.1194 -0.1109   16.09
mean      0.95 raw     0.0421  0.1048 -0.0795   15.81
mean      0.95 guided  0.0421  0.1002 -0.0795   16.93
mean      1.00 raw     0.0421  0.1213 -0.1080   13.58
mean      1.00 guided  0.0421  0.1179 -0.1080   14.70
brightest 0.95 raw     0.0549  0.1034 -0.0745   15.37
brightest 0.95 guided  0.0549  0.0981 -0.0745   16.46
brightest 1.00 raw     0.0549  0.1178 -0.1028   13.15
brightest 1.00 guided  0.0549  0.1142 -0.1027   14.23

== motorcycle, beta = 0.4 /m: t in [0.1344, 0.4299]; pixels with t constant over their patch: 0
checks: oracle max|J-J*| 1.2e-15 | GF const 3.6e-12 | GF linear (eps=0) 1.6e-12 | box7 1.4e-12 | box30 1.6e-13 | minfilter-erosion 0.0e+00
exact t-hat = 1 - w + w t (1 - D_J), true A, w = 0.95: max|err| where t constant n/a, over all pixels 2.2e-01
exact t-hat = 1 - w + w t (1 - D_J), true A, w = 1.00: max|err| where t constant n/a, over all pixels 2.4e-01
PSNR hazy input 8.45 dB, contrast stretch 13.34 dB
A estimates: mean [0.8722 0.8895 0.9211], brightest [0.8763 0.8922 0.9231]
A         w    t         |dA|   MAE t  bias t    PSNR
true      0.95 raw     0.0000  0.0645 -0.0195   17.67
true      0.95 guided  0.0000  0.0606 -0.0194   18.15
true      1.00 raw     0.0000  0.0750 -0.0575   14.81
true      1.00 guided  0.0000  0.0717 -0.0575   15.97
mean      0.95 raw     0.0222  0.0659 -0.0157   16.77
mean      0.95 guided  0.0222  0.0616 -0.0157   17.36
mean      1.00 raw     0.0222  0.0734 -0.0535   13.61
mean      1.00 guided  0.0222  0.0702 -0.0535   14.70
brightest 0.95 raw     0.0263  0.0664 -0.0136   16.55
brightest 0.95 guided  0.0263  0.0619 -0.0135   17.15
brightest 1.00 raw     0.0263  0.0723 -0.0513   13.40
brightest 1.00 guided  0.0263  0.0692 -0.0512   14.47

Checks. Every component check held. The oracle reproduces both scenes to within $10^{-14}$. The derived identity holds to within $2 \times 10^{-16}$ on the 124,144 synthetic pixels whose transmission is constant over their patch, for both $\omega$ and both $\beta$. The motorcycle has none (the count prints 0). Over all pixels the identity misses by 0.2–0.3 on both scenes, which shows how much the pixels where $t$ varies within the patch depart from it. The guided filter passes constants and, with $\epsilon = 0$, linear functions of its guide; the integral-image box means match SciPy; the minimum filter equals the erosion; and the prior holds exactly on the synthetic prior region. The motorcycle’s own dark channel has a median of 0.0824, and 54% of its pixels are at or below 25/255, against He et al.’s 90% for 5,000 outdoor images without sky: the prior holds weakly on this indoor scene.

Where the prior holds. With the true $\mathbf{A}$ and $\omega = 1$, the raw transmission is exact on patch-constant prior pixels at both $\beta$; with $\omega = 0.95$ it is too high there, by the $0.05(1 - t)$ of the derivation. Expectation 1 held. On the white object, the transmission error is the largest of any region, and the recovered object is darker than the truth: the recovery divides its colour difference from $\mathbf{A}$ by an estimated transmission (bounded by $t_0$) below the true one, which pushes it away from $\mathbf{A}$. Expectation 2 held at both $\beta$. This is He et al.’s white-marble failure, reproduced by construction.

Depth edges and the guided filter. The raw error on prior pixels at edges and on the ramp is above the patch-constant error, so the first half of expectation 3 held. Its second half failed at both $\beta$: the guided filter raised the edge-and-ramp error instead of lowering it, and it also moved the patch-constant pixels, where the raw estimate was exact, away from the truth. Expectation 6 held with $\omega = 0.95$ on the synthetic scene and in every case on the motorcycle, and failed with $\omega = 1$ on the synthetic scene for every atmospheric-light rule at both $\beta$. In our reading, the construction explains both failures: the 61 × 61 window reaches up to 30 pixels across depth edges, much farther than the 15 × 15 patch, so it mixes transmissions from both sides, also at pixels where the raw estimate was exact. A lower transmission error did not always mean a higher PSNR: with the true $\mathbf{A}$ and $\omega = 1$, the refined map gives the higher PSNR at $\beta = 0.2$ and the lower at $\beta = 0.4$. In our reading this is because the recovery divides by $t$, so the error in $\mathbf{J}$ depends on where the transmission error falls, not only on its mean.

The atmospheric light. On the synthetic scene both rules return the same $\hat{\mathbf{A}}$: each sky row is constant, so its dark channel is too, and the selected pixels are identical sky pixels; the error is within the bound from the construction at both $\beta$. On the motorcycle, which has no sky, the error is larger for both rules, and larger for He et al.’s single brightest pixel than for IPOL’s mean at both $\beta$. Expectation 5 held. In every printed setting an estimated $\mathbf{A}$ gives a lower full-image PSNR than the true one (not always without the sky).

The motorcycle. The raw transmission with the true $\mathbf{A}$ is biased low at both $\beta$, the direction the identity gives for $D_J > 0$ (although no motorcycle pixel has a patch-constant $t$), and $\omega = 0.95$ beats $\omega = 1$ in PSNR in all twelve comparisons, so expectation 4 held. In our reading, raising $\tilde{t}$ through $\omega < 1$ partly offsets the low bias that the weak prior causes. With IPOL’s defaults the dark channel prior improves on the hazy input at both $\beta$, so expectation 7 held on the synthetic scene and on the motorcycle at $\beta = 0.4$, and failed on the motorcycle at $\beta = 0.2$, where we had predicted that it would not beat the hazy input. Not hypothesized: at $\beta = 0.2$ the global contrast stretch scores higher than the dark channel prior with IPOL’s defaults on the motorcycle (17.90 against 16.93 dB), while at $\beta = 0.4$ it scores far lower.

Datasets and Benchmarks

Real hazy photographs have no ground truth, so the papers above evaluate on synthetic haze:

• NYU Depth v2 (indoor RGB-D) is the source of the training and first test sets of AOD-Net (27,256 training and 3,170 test images) and DCPDN (4,000 training images from 1,000 originals, 400 test images), each sampling $\mathbf{A}$ and $\beta$.
• Middlebury stereo images with their ground-truth disparities give AOD-Net’s TestSet B (800 images) and 40 of DCPDN’s TestB images; the example above uses the same idea on one bundled scene.
• RESIDE (Li et al., 2019, not read) provides, as the AirNet authors describe it, an Outdoor Training Set of 72,135 hazy–clean pairs and a Synthetic Objective Testing Set (SOTS) of 500 outdoor pairs.

Two cautions follow from the construction. The synthetic haze obeys the model exactly, so a method that inverts the model is tested only on its estimates of $t$ and $\mathbf{A}$, not on the model’s failures (non-uniform airlight, wavelength-dependent transmission, noise). And, in our reading, results depend on the ranges of $\beta$, $\mathbf{A}$ and depth a benchmark samples: indoor scenes span a few metres, so their transmission varies less across the image than in a landscape.

Evaluation Metrics

On synthetic data, dehazed images are scored against the clean image by PSNR and SSIM, as in AOD-Net’s and DCPDN’s tables; DCPDN also scores the estimated transmission maps by SSIM. The transmission or depth can be scored against ground truth, and the atmospheric light by its error, as in the example. AOD-Net’s authors split each image’s mean squared error into the error of its mean colour and the error of the residual; on TestSet B their residual error is similar to that of DehazeNet and CAP while the mean-colour error is much lower than theirs, which they relate to a better estimate of $\mathbf{A}$ (their Table 3). On real hazy images without ground truth, the papers above compare results visually; image quality assessment covers no-reference measures in general.

Applications

• Photography. He et al. (2011, §1) name consumer and computational photography: higher visibility and corrected colour shift.
• Input to other vision tasks. They also argue that most vision algorithms assume the input is the scene radiance, so haze biases feature detection, filtering and photometric analysis. AOD-Net’s authors measured object detection on hazy images with and without dehazing in front of Faster R-CNN. Most multi-object trackers on the UAVDT drone benchmark did worse in fog (drone vision), and autonomous driving perception lists weather among the limits of its sensors.
• Depth from haze. The transmission is a depth map up to scale, which He et al. present as a by-product for scene understanding and image editing; depth estimation covers other routes to depth.
• Underwater imaging, with the variants above.

Open Problems

• Priors that fail on whole classes of scenes. The dark channel prior fails on bright surfaces, and haze-lines when the airlight dominates the scene; learned methods inherit the biases of their synthetic training data.
• The model itself. Non-uniform atmospheric light and wavelength-dependent transmission fall outside a model with one $\mathbf{A}$ and one $t$ for all channels (He et al., 2011, §6).
• Evaluation. Synthetic benchmarks test inversion of the model that made the haze, and real hazy images lack ground truth.
• Noise, compression artefacts and halos, which parameter choices trade against dehazing strength (Lisani and Hessel, 2024, §4).