A Variation-based Nighttime Image Dehazing Flow with a Physically Valid Illumination Estimator and a Luminance-Guided Coloring Model

Image quality is often reduced in hazy weather, especially during the nighttime when image visibility can be further degraded. In this paper, we propose a robust variation-based nighttime image dehazing flow with a physically valid illumination estimator, a luminance-guided coloring model and a transmission refinement procedure to effectively address this problem.We design a new illumination model to better address the non-global air-light issue in nighttime scenes. Then, we introduce a structure-preserving optimization flow based on Retinex theory to obtain ambient illumination. Color consistency is guaranteed because we use the input image as the initial guess of illumination in our coloring model. A variational procedure is developed to smoothen the estimated transmission map, where the block effect and the halos can be eliminated through the procedure. The proposed luminance-based correction mechanism further improves visual image quality in the presence of a large sky region. Our experiments are implemented based on actual hazy images. The user study indicates that the proposed method can effectively provide color consistency, preserve details, and reduce halo artifacts and noise in the resulting images compared to other state-of-the-art algorithms. When tested on real-world nighttime haze images, our dehazing flow quantitatively achieves 3.06 for NIQE and 2.99 for NR-CDIQA.


I. INTRODUCTION
Images captured in hazy weather usually degrade the performance of several vision-based applications, such as surveillance systems [22], remote sensing [20], and driver assistant systems [5]. Therefore, it is crucial to solve the image dehazing problem. Early dehazing studies were based on the physical dichromatic Middleton model [18], where hazy images could be expressed in terms of two light source types: direct transmission and air-light scattering. In the daytime, sunlight was the most dominant light source, making the ambient illumination globally consistent. Various priors have been proposed to recover haze-free images, such as the dark channel prior (DCP) [6], the boundary constraint [17], and the haze line [1]. New techniques, such as intensity projection and wavelet transform, are also introduced to recover degraded scenes. Liu et al. [14] proposed a rankone transmission prior to estimate the transmission. Khmag et al. [9] used the mean of the vector L2-norm to estimate the transmission map and enhance it by a second-generation wavelet transform. Recently, a multitude of learning-based methods have also been developed. For instance, DehazeNet [3] incorporates the Middleton model, the DCP, and image priors into a deep neural network. The AOD-Net [11] directly restores images from hazy input images using a lightweight convolution neural network. When dealing with daytime haze, these methods perform well. However, they are not suitable for nighttime conditions, which are more complicated and challenging. The major difficulty facing nighttime dehazing is an ambient illumination estimation. Images captured at night usually contain the scattering of multiple-colored light sources, such as road lamps, buildings, and vehicle lights. Therefore, illumination should be modeled as a spatially variant map instead of a single vector in the daytime. Furthermore, the most commonly used DCP is not suitable for nighttime scenes due to the large value of minimal intensity in the light source regions. For the reasons outlined above, several studies adopt new priors, such as a glow removal [13], a maximum reflectance prior (MRP) [27], or a bright channel prior (BCP) [26] to address the problems associated with nighttime scenes. However, the aforementioned methods are limited by several problems, including detail losses, noise amplification, and halo effects due to inaccurate ambient illumination and transmission calculations.
In this study, we propose a new dehazing flow with a variational-model approach based on Retinex theory [10]. The dehazed results can achieve better halo suppression and color preservation in different nighttime scenes when compared with the competing methods. To apply the Middletontype model to nighttime dehazing, we have to design new methods to generate correct ambient illumination and transmission maps. First, we estimate our ambient illumination maps by designing a new variational model based on Retinex theory, where the Sobel operator with a pixelwise weighting strategy is introduced as the regularization term. The optimization scheme can preserve the structure of scenes and smoothen textures in the images. This novel approach can handle various light source conditions and reduce thet halo effects in the resulting image.
After we obtain the illumination map, three distinct priors, the dark channel, the bright channel, and the proposed sky region priors, are adopted to calculate the transmission maps. Subsequently, three transmission maps are fused with brightness and the proposed gradient-based sky map to obtain the final estimation. This fusion process can reduce artifacts in the nighttime sky. To preserve the boundary of the original scene structure, the weighted total variation regularization term is introduced to constrain the transmission and the hazefree image, which can effectively smoothen the unwanted textures, except on sharp edges. We can calculate the final refined transmission map using the alternating direction method of multipliers (ADMM) [2]. Moreover, we can obtain the dehazed image using the calculated illumination and transmission maps.
Because there is no benchmark dataset for nighttime image dehazing, most researchers qualitatively measure the visual quality of the dehazed results. In addition to the qualitative analysis, we also perform comprehensive experiments, including a user study and an image quality assessment (IQA). The user study indicates that the proposed method can effectively maintain color consistency, preserve detail, and reduce halo artifacts and noise compared to state-ofthe-art algorithms. In addition, IQA metrics demonstrate the superiority of our method.
In this paper, we propose a robust variation-based nighttime image dehazing flow with a physically valid illumination estimator, a luminance-guided coloring model and a transmission refinement procedure. These three main contributions can be summarized as follows: The Physically Valid Illumination Estimator For nighttime scenes, the estimation of the illumination map is crucial for achieving satisfactory results. Existing methods often use the local maximum pixel [13], [27] or the low-pass filter [26] to estimate ambient illumination maps, resulting in color distortion or artifacts in the light source. We start from the Retinex theory and adopt a weighted L 1 norm regularization to form the objective function for achieving color consistency. Then, it is solved using the efficient ADMM. The resulting ambient illumination is physically valid, and it can smoothen the high-frequency region and preserve the edge details.
The Luminance-Guided Coloring Model Existing nighttime dehazing methods often generate unsatisfactory results in large nighttime sky regions due to the underestimation caused by the physical-invalid model for the transmission map in the sky. We apply the novel gradient map-based method to detect the sky region, and estimate the transmission with the physically based dual-mode luminance model. This model can handle both dark and bright nighttime skies. The experimental results also indicate that the proposed skycorrection mechanism can preserve sky colors while introducing fewer artifacts than the existing methods.
The Transmission Refinement Procedure Rather than simply using the guided filter [7] as in [13], [26], [27], we propose a new variational regularization model that preserves the edges and mitigates the block effect of the refined transmission map. The halo artifacts are effectively suppressed when compared to other competing methods. Comprehensive experiments, including synthetic datasets, real-world hazy images, and user studies indicate that the proposed method qualitatively and quantitatively achieves better performance than other state-of-the-art algorithms. In addition, because of the robustness of the flow, we can apply the same model parameters to all images demonstrated in this work.

II. RELATED WORK
The dehazing problem faced by daytime scenes can be addressed using various image priors [6], [30], [1] or deep neural networks [3], [11]. However, these approaches are not appropriate for nighttime hazy images due to nonuniform and complicated ambient illumination.
Pei et al. [21] introduced color transfer preprocessing and used a pixelwise DCP to mitigate artifacts. The obtained results were unnatural because of the lack of a physically based model. Zhang et al. [28] introduced an imaging model that included illumination correction and color adjustment. Although these techniques could improve visibility, they also produced halo effects. Li et al. [13] added a glow-associated term to the Middleton model. Although this method reduced the glow problem, it triggered color distortion and artifacts in the light-source region. Zhang et al. [27] proposed a maximum reflectance prior (MRP) estimate of ambient illumination. However, MRP is not appropriate for scenes with distinct colors. Yu et al. [26] used the Retinex theory, and proposed a bright channel prior along with alpha blending to estimate the transmission map. Although it reduced artifacts in the light-source regions, it led to detailed losses   in the dark region. Recently, Zhang et al. [29] proposed a deep neural network using an artificial hazy mechanism to simulate nighttime hazy images. Its performance was substantially related to the training datasets; hence, it could not guarantee satisfactory results under various conditions. These aforementioned methods were limited by either unnatural results or missing details. We form a variational model to estimate varying atmospheric light based on Retinex theory [10] instead of simply using low-pass filters, such as [27], [13], and [26]. Furthermore, we estimate the transmission of the nighttime sky using a luminance-based model. Finally, we refine the transmission using a variational model. By addressing these problems, we infer that our results suppress halos and preserve color consistency. Therefore, the proposed method can obtain more satisfactory haze-free scenes than existing nighttime dehazing methods.

III. THE PROPOSED NIGHTTIME DEHAZING METHOD
The flow chart of the proposed nighttime dehazing method is illustrated in Fig. 1. We first compute the ambient illumination map L using Algorithm 1 before computing the transmission map. Then, we calculate three transmission maps using image priors and Algorithm 2. Subsequently, we fuse the transmission maps together to obtain an initial estimation, and refine it using Algorithm 3. Finally, we restore the hazefree image based on the Middleton model. The output results for each step are presented in Fig. 2. The details of each step are described in the following subsections.

A. THE NIGHTTIME IMAGING MODEL
Based on Koschmieder's law [18], the hazy imaging process can be illustrated by the physically based dichromatic model as follows: where •, I, A, and J represent an elementwise multiplication, the observed hazy image, a global uniform vector representing the air light, and the haze-free radiance of the scene, respectively. The scene transmission t is distance-dependent as follows: where η denotes the extinction coefficient of the haze, and d represents the distance from the camera to the scene. In daytime scenes, sunlight is the dominant light source, which provides globally uniform illumination. However, in the case of nighttime hazy images, the scene usually contains multiple-colored light sources, such as buildings, road lamps, and vehicle lights. The ambient illumination A should be replaced by a spatially variant map L. Hence, the image can be modeled as follows: The objective of nighttime image dehazing is to restore a haze-free image J from an observed hazy image I. Because I is solely known in (3), we need to first estimate L and t using the steps depicted in the following sections.

B. THE PHYSICALLY VALID ILLUMINATION ESTIMATOR
Based on the Retinex theory [10], we can represent a hazefree image as a product of illumination L and intrinsic reflectance ρ. With J = ρ • L, (2) can be reformulated as follows: Ambient illumination is preferred to contain the structure of a scene and it is spatially smooth. To solve L, we design a variational model approach with a Sobel filter-based regularization term. The objective function of the problem in the matrix form can be expressed as follows: where ∥·∥ F and ∥·∥ 1 are the Frobenious and L 1 norms, respectively. • and ⊗ denote the pixelwise multiplication and convolution operators, respectively. Compared with the commonly used total variation (TV) regularizer, we employ high-order operators in (5) to make it more flexible when applying contextual constraints. Furthermore, we select the proper weighting function W(x) based on the following analysis.
Pixels with a similar color should have similar illuminance values in local patches, and illuminance jumps often appear Algorithm 1: illumination estimation initialization: (14); at sharp edges. Accordingly, we can compute the color difference of local pixels to construct the weighting function as follows: The optimization problem (5) can be effectively solved using the alternating direction method of multipliers (ADMM). An auxiliary variable V is introduced to replace S j ⊗ L to simplify the problem for easy solving. Consequently, we have the equivalent optimization problem as follows: The augmented Lagrangian function of (7) can be rewritten as follows: where ⟨·, ·⟩ denotes the matrix inner product, and µ is a penalty scalar. When solving (8), the ADMM iteratively updates one variable at a time by fixing the others. We provide the solutions to the following subproblems: a: The Subproblem of L: Dropping the terms unrelated to L provides the problem as follows: This is a least-squares problem. Hence, the solution can be computed by differentiating with respect to L and setting it to zero.  [13], [26], [27] and our method, respectively.
Then, we can apply the 2-D FFT technique to the above problem with circular boundary conditions, which enables us to compute the solution quickly. Then, we have the equation as follows: where F(·) and F(·) are the 2D FFT operator and its complex conjugate, respectively. F −1 (·) represents the 2D inverse FFT.
b: The Subproblem of V : Collecting the terms involving V from (8) leads to the optimization problem as follows: The closed-form solution of (12) can be obtained by the shrinkage operation as follows: where Ψ ϵ [x] = sign(x)max(|x| − ϵ, 0) represents the shrinkage operator.
c: The Subproblems of ω and µ: The updating of ω and µ can be performed as follows: Fig. 3 shows the comparison of the illumination results for different λ. The smoothness of L increases when parameter λ decreases. Thus, we choose the λ value, which can generate an illumination map with only sufficient contour information, but without high-frequency details.
The entire procedure is summarized in Algorithm 1. Fig.  4 presents a comparison of the ambient illumination map estimated by [13], [26], [27] and the proposed method. Evidently, the illumination map estimated by [13] has color distortions due to the excessive glow removal process. [26], [27] has unnatural dark pixel estimation in the center of the light source, which triggers halo effects in the resulting image. Unlike the above artifacts, our method can achieve color consistency and preserve the edges and structures of the original scene. First, we adopt the dark channel prior (DCP) [6] and the bright channel prior (BCP) [26]. The DCP assumes that in most non-sky patches, at least one color channel has a very low intensity in a local patch. BCP assumes that in most light source regions, at least one color channel has an intensity close to one in a local patch. The equation can be calculated as follows:

C. THE LUMINANCE-GUIDED COLORING MODEL
x denotes the spatial coordinate of each pixel and Ω represents a local patch centered at x. After we obtain L using Algorithm 1, we can compute the transmission maps t DCP and t BCP by applying the DCP and the BCP as follows: Because the DCP is valid for nonlight source regions, it is only valid for light source regions. Therefore, we need to combine them using a weight map P light [26]. This can be calculated as follows: P light is higher when the pixel is bright, and vice versa. Then, the transmission map t is computed as follows: However, the DCP and the BCP assumptions are solely valid for the non-sky region. Hence, to address the nighttime sky region problem, we propose a novel sky transmission correction mechanism with a modified luminance-based model as follows: where β denotes the atmospheric scattering coefficient and L(x) represents the normalized luminance. Notably, we employ the maximum operator to make it suitable for both bright and dark nighttime skies. After we model the transmission map for the nighttime sky, we detect the sky region based on the following analysis. Considering that the sky region usually has a small gradient, we first calculate the gradient map of the luminance channel in the CIE LAB color space. Then, we generate a binary map with a threshold value that indicates whether the gradient in a pixel is sufficiently small. Next, we apply a median filter to eliminate the noise. Finally, the largest connected component region is considered the sky region. Fig. 5 presents examples of the sky detection results. Different light sources and sky conditions are taken into consideration, such as dark skies and bright skies, as well as vari-colored artificial light sources. As shown in Fig. 5, the detection results are close to those of the human visual system. The sky detection procedure is summarized in Algorithm 2. After obtaining P sky , we can combine t L (x) and t BDCP (x) to obtain the final transmission estimation as follows: Fig. 6 presents the dehazed results using t DCP , t BCP , t BDCP , t L , and our proposed t est , respectively. It can be  observed that t DCP tends to over-dehaze the light source regions while preserving details in the nonlight source region. t BCP optimally handles the light source region but overenhances the dark areas. t BDCP strikes a balance between t DCP and t BCP ; however, it still generates halo regions around the building. t L handles the sky problem, but the scenes tend to under-dehaze. The result of the proposed t est benefits from the fusion process with all the advantages of t DCP , t BCP , and t L .

D. THE TRANSMISSION REFINEMENT PROCEDURE
If we directly employ t est for dehazing, the results contain some block effects because the transmission is not always constant in a patch. To further enhance the dehazing performance, a regularization-based variational model for transmis-sion map refinement is presented as follows: arg min where t 0 = t est . In addition, we select the relative total variation [24] as the weighting function as follows: where K σ represents the Gaussian kernel with a variance of σ 2 . Fig. 7 shows the example for the parameter selection. When λ 2 is too small, it causes artifacts around the light source region. However, it might contain a block effect in the resulting images when λ 2 is large. Thus, we choose λ 2 = 1.5 in our model to handle both issues.
The auxiliary variables, G and H, are introduced to replace ∇J and ∇t, and simplify the problem for easy solving. Accordingly, ∇J = G and ∇t = H are introduced as constraints. Similarly to (8), we have the augmented Lagrangian function as follows: The ADMM solver is again applied to solve the problem (22). The setups of the subproblems are provided below. First, we transform the data term F to solve it efficiently using FFT [15]. Then, we drop the terms unrelated to J as follows: This is a least-squares problem. Similarly to (10), we can apply 2-D FFT techniques to the above problem. Then, we have the equation as follows: max(t (k) ,ϵ) + L, ϵ is a small positive value to ensure that the denominator is not zero.
. Then, dropping the terms unrelated to t yields the problem as follows: Similarly to (25), we have a closed form solution as follows: c: The Subproblems of G and H: Similarly to (12), the closed-form solution for the G and H subproblems can be obtained by the shrinkage operation as follows: ). (29) d: The Subproblems of ϕ, δ, and µ : ϕ, δ, and µ can be updated as follows: The entire procedure is summarized in Algorithm 3. Accordingly, we can obtain the final transmission estimation, t ref .

E. HAZE-FREE IMAGE RECONSTRUCTION
After calculating the ambient illumination L and the refined transmission t ref maps, we can obtain the haze-free image J by applying Middleton's model as follows: where t lb denotes the lower bound of the transmission map. As suggested in [6], selecting t lb as 0.1 is suitable for most practical applications.

IV. EXPERIMENTAL RESULT
In this section, the experimental results are presented. The proposed method is quantitatively and qualitatively compared with existing state-of-the-art methods.

A. EXPERIMENTAL SETTINGS
All experiments were run on MATLAB R2020b with 16 GB RAM and an Intel Core i5-9400F @2.90 GHz. The parameters λ and µ in (9) are empirically fixed to 2 and 1, respectively. In (23), the parameters are λ 1 = 0.01, λ 2 = 1.5, λ 3 = 0.5, λ 4 = 1.0 and µ 1 = µ 2 = 1. In general, these settings performed well. To evaluate the effectiveness of the proposed dehazing flow, we compare the proposed method with stateof-the-art prior-based approaches, i.e., glow removal of light sources (GS) [13], maximum reflectance prior (MRP) [27], and pixelwise alpha blending (PAB) [26]. In addition, a deep learning-based method, optimal-scale fusion-based dehazing (OSFD) [29], is included for comparison. The results of these methods are generated by the code downloaded from the websites maintained by the original authors of the related work, with the recommended experimental settings. Our code is available at https://github.com/chyang913/NightDehazing

C. QUANTITATIVE COMPARISON OF SYNTHETIC DATA
First, we conduct quantitative evaluations of the dehazing methods using synthetic data. We follow the procedures described in [26] to synthesize nighttime hazy images. Twelve images from the Middlebury 2005 and 2006 datasets [23] are selected to generate synthetic data. The transmission map t is computed from the normalized depth map d using t(x) = 0.8d(x). We assume that the light source is in the center of the image, and the illumination map can be regarded as L(x) = 1 − 0.8 × dist(x), where dist(x) represents the distance from the light source to the pixel. Finally, we adopt (4) to obtain synthetic hazy images. Table 1 presents the PSNR and the SSIM values of various methods. We can conclude that our approach remains competitive with stateof-the-art methods. Fig. 9 illustrates visual comparisons for the dehazed results using synthetic hazy images. We can observe that the MRP and the OSFD obtain unnatural shadows around the light source because of their excessive color adjustment procedure. GS generates over-dehazed images, while our method produces the most natural results. (c) GS [13] (d) MRP [27] (e) OSFD [29] (f) Ours

D. A QUALITATIVE COMPARISON ON THE REAL-WORLD DATASET
Because synthetic images cannot reproduce all practical natural nighttime scenes, we further conduct qualitative evaluations of the proposed dehazing flow using actual data. We test 250 real-world nighttime hazy images provided by the authors of [27], [29], [12]. Each contains 20, 130, and 100 images, respectively. Fig. 8 presents comparisons between actual nighttime dehazed results generated by different methods. The proposed method produces a more natural result with fewer artifacts than state-of-the-art methods. Fig. 10 and Fig. 11 present more dehazed results than their close-up versions. In Figs. 10, we can observe that the GS and the PAB over-enhance the lighthouse. The MRP and the OSFD produce unnatural shadows around the light sources. These methods also amplify the noise in the sky region. Unlike the aforementioned methods, the proposed dehazing flow not only provides satisfactory results without artifacts but also preserves details from the input image, such as the bottom of the metal column in the train station. In Figs. 11, we can observe that the MRP and the OSFD trigger color distortion results due to the excessive color adjustment procedures, while our results maintain the color consistency of the original scene. The GS and the PAB overexpose the light source regions, thereby reducing the visibility of hands on the clock. In contrast, we could eliminate the undesirable halo effects, resulting in a more natural-looking appearance, as well as enhancing the visual appearance. From the above visual comparisons, we can conclude that the proposed method can generally better preserve the detailed texture, maintain color consistency, preserve the shape of light sources without a halo, and produce fewer artifacts in sky regions.

E. USER STUDY ON A REAL-WORLD DATASET
Subjective feelings for the visual quality of the actual scenes might differ for each person. We conduct a user study using the two-alternative forced-choice (2-AFC) paradigm to make a qualitative comparison. We adopt 20 input hazy images provided by the authors of [27], [29], [12], such that there are ten paired comparisons for each image. All comparisons are randomly arranged. Subsequently, we recruit 20 participants from Amazon Mechanical Turk. Instructions inform users that a good result should contain the following properties: no under-or overexposed regions, no visible noise, and no color distortion. Table 2 presents the preference matrix of the user study. We can observe that every element in our row is larger than 200, which is half of the total comparisons made by the two methods. In other words, most viewers prefer our dehazed results when compared with the other four methods. Thus, we can conclude that our dehazing flow generally achieves a better performance. To compare the visual quality of the dehazed results while handling various actual scenes, we further analyze the user preference for each image. Table 3 summarizes the ranking distributions of the five methods. The proposed approach outperforms other state-of-the-art methods in 13 out of 20 images, and does not obtain an image with the fifth rank.
The results indicate that the proposed method can obtain (a) Input hazy image (b) GS [13] (c) MRP [27] (d) PAB [26] (e) OSFD [29] (f) Ours (g) Input hazy image (h) GS [13] (i) MRP [27] (j) PAB [26] (k) OSFD [29] (l) Ours sequently, we assess the objective qualities of the datasets provided by [27], [12], and [29]. For the NIQE, CNN, and MEON, smaller values represent better image qualities, and for the NR-CDIQA and the CEIQ, larger values indicate better qualities of image contrast. The results of the different approaches are summarized in Table 4. The NIQE measures the distance between natural scene statistics (NSS) features. The CNN predicts image quality by adopting image patches as the input and applying the CNN in the spatial domain without using handcrafted features. The MEON is a deep neural network comprising distortion identification and quality prediction. The NR-CDIQA employs support vector regression (SVR) to predict the human mean opinion score (MOS) of contrast-distorted images from multiple NSS features. The CEIQ takes features of the SSIM, entropy, and cross-entropy between the original and the enhanced images to predict the quality of the contrast-changed images.
From Table 4, we can deduce that the results of our method exhibit the best quality on these datasets for all objective measurements. The NIQE, CNN, and MEON prefer natural scenes without artifacts. The NR-CDIQA and the CEIQ evaluate the visibility of the resulting images. In summary, the proposed method can obtain natural results and generate fewer artifacts while enhancing visual acuity.

G. LIMITATION
It should be noted that this work still has limitations for images taken in rainy and snowy weather conditions. As shown in Fig. 12, there are some streak artifacts caused by rain. The main reason for this is that the Middleton-typed model does not perform well under rainy conditions. In such cases, more sophisticated models will be necessary in order to overcome the issues.

V. CONCLUSION
In this paper, we propose an effective method for dehazing nighttime images. The key to nighttime image dehazing is estimating the ambient illumination and the transmission maps. Therefore, a structure-aware smoothing model is developed to consistently improve the maps. Furthermore, ADMM is adopted to effectively optimize the objective function. We present a comprehensive experimental analysis of the proposed model using both subjective and objective assessments. The obtained experimental results indicate that the proposed method is better than other state-of-the-art methods. In addition, our nighttime image dehazing technique can be adopted as the preprocessing stage of several computer vision applications, such as driver assistant systems, remote sensing, and object detection, thereby improving their performance. YI-CHANG LU (S'00-M'05-SM'12) received the B.S. degree in electrical engineering from National Taiwan University, Taipei, Taiwan, in 1993, and the M.S. degree in electrical engineering, the M.S. degree in engineering-economic systems, and the Ph.D. degree in electrical engineering from Stanford University, Stanford, CA, USA, in 1997, 1999, and 2005, respectively.
He was an Engineering Officer with the Naval Surveillance and Communication Command Department, Suao, Taiwan, from 1993 to 1995. He was a Post-Doctoral Research Fellow with Stanford University in 2005. Since 2006, he has been with the Graduate Institute of Electronics Engineering and the Department of Electrical Engineering, National Taiwan University, Taipei, Taiwan, where he is currently a Professor. His research interests include VLSI circuit and architecture design, digital signal processing, and high performance scientific computing.