1.

Introduction

Weather conditions can severely limit visibility in outdoor scenes. In such cases, atmospheric phenomena such as fog and haze will significantly degrade visibility in the captured scene. Since visibility is dependent on the air, the amount of particles in the air will affect image visibility. This phenomenon is generally composed of water droplets or particles and cannot be ignored. Both the absorption and scattering of light by particles and gases in the atmosphere cause the visibility to decrease, whereas the scattering of particulate in visibility causes more serious damage than the absorption of light. As a result, distant object and part scenes are not visible. That is, the image loses contrast and color fidelity, and the visual quality of the scene is reduced. In a visual sense, the quality of the degraded image is unacceptable. Therefore, a simple and effective image scene recovery method is essential. Image dehazing is a challenging problem, and image recovery technology has attracted the attention of many researchers. The low visibility in hazy images affects the accuracy of computer vision techniques, such as object detection, face tracking, license plate recognition, satellite imaging, and so on, as well as multimedia devices, such as surveillance systems and advanced driver assistance systems. Hence, haze removal techniques are important for improving the visibility of images. Restoring hazy images is a particularly challenging case that requires specific strategies. Therefore, widely varying methods have emerged to solve this problem. In recent years, enhancing images represents a fundamental task in many image processing and vision applications. Proposed strategies for enhancing the visibility of a degraded image include the following.

The first type is the nonmodule method, such as histogram equalization,¹ Retinex theory,² wavelet transform,³ and gamma correction curve.⁴ However, the shortcomings of these methods are that they seriously affect the clear region and also keep color fidelity less effectively.

The second type is the module method, which depends on the physical mode. Compared to the nonmodule method, these methods achieve better dehazing results by modeling scattering and absorption and by using multiple different atmospheric conditions in input images, such as scene depth,^5,6 multiple images,^7–9 polarization angles,^10,11 and geometry models.^12,13 Narasimhan and Nayar^7,8 developed an interactive depth map for removing weather effects, but their method had limited effectiveness. Kopf et al.¹³ presented a novel deep photo system for using prior knowledge of the scene geometry when browsing and enhancing photos. However, the method required multiple images or additional information to get a better estimate of scattering and absorption, which limited its applications. Hautière et al.¹² designed a method using weather conditions and a priori structure of a scene to restore the image contrast for vehicle vision systems.

A novel technique developed in Refs. 10, 14, and 15 exploited the partially polarized properties of airlight. The haze effect was estimated by using different angles of polarized filters to analyze the resulting images of the same scene. In other words, calculating the difference among these images enabled the use of the magnitude of polarization to estimate haze light components. Because the polarization light is not the major degradation factor, these methods have less robustness for scenes with dense haze.

Another recently developed strategy used a module and a single hazy image as input information. This approach has recently become a popular way of eliminating image haze by different strategies.^16–20 Roughly, these methods can be categorized as contrast-based and statistical approaches. An example of a contrast-based approach is the Tan¹⁷ method. In this case, the image restoration maximizes the local contrast while limiting the image intensity to be less than the global atmospheric light value. Tarel and Hautière¹⁹ combined a computationally effective technique with a contrast-based technique. Their method assumed that the depth map must be smooth except along edges with large depth jumps. The second category of statistical approaches includes the technique presented in Fattal,¹⁶ which employs a graphical model to solve ambiguous atmospheric light color and assumes the image shading and scene transmission are partially uncorrelated. According to this assumption, mathematical statistics were utilized to estimate the albedo of a scene and infer the transmission medium. The method provides a physically consistent estimation. However, because the variation of the two functions in Ref. 16 is not obvious, this method requires substantial fluctuation of color information and luminance in the hazy scene. He et al.¹⁸ developed a statistical approach for observing the dark channel and for roughly estimating the transmission map. Then, they refined the final depth map by using a relatively computationally expensive matting strategy.²¹ In this approach, pixels must be found through the entire image, which requires a long computation time. Nishino et al.²⁰ used a Bayesian probabilistic concept by fully leveraging their latent statistical structures to estimate the scene albedo and depth from a single degraded image. A recent study by Gibson and Nguyen²² proposed a new image dehazing method based on the dark channel concept. Unlike the previous dark channel method, their method finds the average of the darkest pixels in each ellipsoid. However, this assumption in Ref. 22 may find several inaccurate pixels for those corresponding to bright objects. Fattal²³ derived a local formation model that explains color lines in the context of hazy scenes and used the model to offset lines for recovering the scene transmission. In addition, Ancuti and Ancuti²⁴ also proposed a fusion-based strategy for enhancing white balance and contrast in two original hazy image inputs. In other words, in order to keep the most significant detected features, the inputs in the fusion process are weighted by the specific calculation maps.

Recently, artificial neural networks (ANNs) have been widely used in many different fields. Research topics related to ANNs have proved suitable for many areas, such as control,^25,26 identification,^27,28 pattern recognition,^29,30 equalization,^31,32 and image processing.^33,34 The cerebellar model articulation controller (CMAC) model proposed by Albus^35,36 is usually applied in ANNs. The CMAC model imitates the structure and function of the cerebellum of a human brain and it is similar to a local network. The CMAC model can be viewed as a basis function network that uses plateau basis functions to compute the output of the model for a given input data point. Therefore, only the basis functions assigned to the hypercube covering the input data are needed. In other words, for a given input vector, only a few of the network nodes (or hypercube cells) are active and will effectively contribute to the corresponding network output. Thus, the CMAC has good learning and generalization capabilities. However, the CMAC requires a large amount of memory for solving the problem of the high dimension,^37,38 is ineffective for online learning systems,³⁹and has relatively poor function approximation ability.^40,41 Another problem is that it is difficult to determine the memory structure, e.g., to adaptively select structural parameters, in the CMAC model.^42,43 Recently, several researchers have proposed various solutions for the above problems, including fuzzy membership functions,⁴⁴ selection of learning parameters,⁴⁵ topology structure,⁴⁶ spline functions,⁴⁷ and fuzzy C-means.⁴⁸ Fuzzy theory embedded in the CMAC model has been widely discussed. Thus, a fuzzy CMAC called FCMAC⁴⁹ was proposed. It takes full advantage of the concept of fuzzy theory and combines it with the local generalization feature of the CMAC model.^49,50 A recurrent network is embedded in the CMAC model by adding feedback connections with a receptive field cell to the CMAC model,⁵¹ which has the advantage of dynamic characteristics (considering past output network information). However, the above-mentioned methods have several drawbacks. For example, the mapping capability of local approximation by hyper-planes is not good enough, and more hypercube cells (rules) are required.

Therefore, this study developed a recurrent fuzzy cerebellar model articulation controller (RFCMAC) model to solve the above problems and to enable applications in widely various fields. A hybrid of the recurrent fuzzy CMAC and weighted strategy is used to process the image dehazing problem. The proposed method provides high-quality images and effectively suppresses halo artifacts. The advantages of the proposed method are as follows:

The rest of this paper is structured as follows. Section 2 discusses the theoretical background of light propagation in such environments. In Sec. 3, we introduce the proposed RFCMAC and weighted strategy for image dehazing. Section 4 presents the experimental results and compares the proposed approach with other state-of-the-art methods. Finally, conclusions are drawn in Sec. 5.

2.

Theory of Light Propagation

Generally, a camera being used to take outdoor photographs obtains an image by the light of the receiving environment, such as the illumination of sunlight, reflecting light from a surface as shown in Fig. 1. Due to absorption and scattering, the light crossing the atmosphere is attenuated and dispersed. In physical terms, the number of suspended particles is low in sunny weather. Thus, the image quality is clear. In contrast, dust and water particles in the air during volatile weather scatter light, which severely degrades image quality. In such degraded circumstances, only 1% of the reflected light reaches the observer, and it causes poor visibility.⁵⁴ McCartney⁵⁵ also noted that haze is an atmospheric phenomenon. That is, the clear sky is obscured by dust, smoke, and other dry particles. In the captured image, the haze generates a distinctive gray hue, reducing visibility for the image. Based on the above, the physical theory of a hazy model can be expressed as

Fig. 1

Haze imaging model.

3.

Proposed Method

This section presents in detail our proposed method, which uses the RFCMAC model and a weighted strategy to recover scenes from the removal of a hazy image. Figure 2 shows the flowchart of the proposed method, and the details are presented in the following sections.

Fig. 2

Flow diagram of the proposed dehazing algorithm.

3.1.

Estimation of Transmission Map Features Using RFCMAC Model

The transmission map and atmospheric light have important roles in haze removal. Therefore, a good dehazing method with estimation of both the transmission map and the atmospheric light can appropriately process the recovery of a hazy image. Haze, which is generated by light attenuation, depends on the distribution of the number of particles in the air. According to Eq. (1), both the transmission map and the atmospheric light are important factors. Thus, the transmission factor and atmospheric lightness must be improved. This study proposes an RFCMAC model for estimating the transmission map more accurately. The RFCMAC model combines the traditional CMAC model, an interactive feedback mechanism, and a Takagi—Sugeno—Kang (TSK)-type linear function to obtain better solutions. The proposed model also adopts an interactive feedback mechanism, which has the ability to capture critical information from other hypercube cells. The structure of the RFCMAC and associated learning algorithm are presented as follows.

3.1.1.

Structure of the RFCMAC model

The performance of the proposed RFCMAC model is enhanced by using an interactive feedback mechanism in the temporal layer and a TSK-type linear function in the subsequent layer. Figure 3 shows the six-layered structure of the RFCMAC model. The structure realizes a similar fuzzy IF–THEN rule (hypercube cell).

Fig. 3

Structure of the RFCMAC model.

Rule $j$:

$$\begin{gathered} \text{IF}{x}_1\text{is}{A}_{1j}\text{and}{x}_2\text{is}{A}_{2j}\dots \text{and}{x}_i\text{is}{A}_{ij}\dots \text{and}{x}_{N_D}\text{is}{A}_{N_{D}j} \\ \text{THEN}{y}_j=\sum _{j=1}{O}_j^{(4)}({\alpha }_{0j}+\sum _{i=1}^{N_D}{\alpha }_{ij}{x}_i), \end{gathered}$$

where $x_i$ represents the $i$’th input variables, $y_j$ denotes the local output variables, $A_{ij}$ is the linguistic term using the Gaussian membership function in the antecedent part, $O_j^{(4)}$ is the output of the interactive feedback, and ${\alpha }_{0j}+\sum _{i=1}^{N_D}{\alpha }_{ij}{x}_i$ is the basis TSK-type linear function of input variables. The operation functions of the nodes in each layer of the RFCMAC model are described as follows. For the following description, $O^{(l)}$ represents the output of a node in the $l$’th layer.

Layer 1 (input layer): The layer is used as an input feature vector $\stackrel{\rightharpoonup }{\mathbf{x}}=(x_1,x_2,\dots ,x_{N_D})$, and the inputs are crisp values. This layer does not require adjustments of weight parameters. Each node need only directly transmit input values to the next layer. The corresponding outputs are calculated as

Eq. (2)

$$O_i^{(1)}=u_i^{(1)},\text{and}{u}_i^{(1)}=x_i.$$

Layer 2 (fuzzification layer): The layer performs a fuzzification operation and uses a Gaussian membership function to calculate the firing degree of each dimension. The Gaussian membership function is defined as follows:

Eq. (3)

$$O_{ij}^{(2)}=\exp [-\frac{{(u_i^{(1)}-m_{ij})}^2}{{\sigma }_{ij}^2}],\text{and}{u}_i^{(2)}=O_i^{(1)},$$

where $m_{ij}$ and ${\sigma }_{ij}$ denote the mean and variance of the Gaussian membership function, respectively.

Layer 3 (spatial firing layer): Each node of this layer receives the firing strength of the associated hypercube cell by the node of a fuzzy set in layer 2. All layer two outputs are collected in layer three. Specifically, each node performs an algebraic product operation on inputs to generate spatial firing strength ${\alpha }_j$. The layer determines the number of hypercube cells in the current iteration. For each inference node, the output function can be computed as follows:

Eq. (4)

$$O_j^{(3)}=\prod _i^{N_D}{u}_{ij}^{(3)},\text{and}{u}_{ij}^{(3)}=O_{ij}^{(2)},$$

where $\mathrm{\Pi }$ denotes product operation.

Layer 4 (temporal firing layer): Each node is a recurrent hypercube cell node, including the internal feedback (self-loop) and external interactive feedback loop. The output of the recurrent hypercube cell node depends on both the current spatial and previous temporal firing strengths. That is, each node refers to relative information from itself and other nodes. Because the self-feedback of the hypercube cell node is not sufficient to represent the all necessary information, the proposed model refers to relative information not only from the local source (node’s feedback from itself) but also from the global source (feedback from other nodes). The linear combination function of the temporal firing strength is described as follows:

Eq. (5)

$$O_j^{(4)}=\sum _{k=1}[{\lambda }_{kj}^q·O_k^{(4)}(t-1)]+(1-{\gamma }_j^q)·u_j^{(4)},\text{and}{u}_j^{(4)}=O_j^{(3)},$$

where ${\lambda }_{kj}^q$ represents recurrent weights and determines the compromise ratio between the current and previous inputs to the network outputs. ${\gamma }_j^q=\sum _{k=1}^{N_A}{\lambda }_{kj}^q$ and ${\lambda }_{kj}^q=\frac{R_{kj}^q}{N_A}(0\le R_{kj}^q\le 1)$ denote the interactive weights of the hypercube cells from itself and other nodes. $R_{kj}^q$ is a connection weight from the $k$’th node to the $j$’th node and is a random value between 0 and 1. $N_A$ is the number of hypercube cells. Therefore, the compromise ratio between the current and previous inputs is between 0 and 1.

Layer 5 (consequent layer): Each node is a function of a linear combination of input variables in this layer. The equation is expressed as

Eq. (6)

$$O_j^{(5)}=O_j^{(4)}(a_{0j}+\sum _{i=1}^{N_D}{a}_{ij}{x}_i).$$

Layer 6 (output layer): This layer uses the centroid of area (COA) approach to defuzzify a fuzzy output into a scalar output. Then the actual output $y$ is derived as follows:

Eq. (7)

$$y=\frac{{\sum }_{j=1}^{N_A}{O}_j^{(4)}(a_{0j}+{\sum }_{i=1}^{N_D}{a}_{ij}{x}_i)}{{\sum }_{j=1}^{N_A}{O}_j^{(4)}}.$$

3.1.2.

Learning algorithm of the RFCMAC model

The proposed learning algorithm combines structure learning and parameter learning when constructing the RFCMAC model. Figure 4 shows a flowchart of the proposed learning algorithm. First, the self-constructing input space partition in structure learning is based on the degree measure used to appropriately determine the various distributions of the input training data. In other words, the firing strength in structure learning is used to determine whether a new fuzzy hypercube cell (rule) should be added to satisfy the fuzzy partitioning of input variables. Second, the parameter learning procedure performs the back propagation algorithm by minimizing a given cost function to adjust parameters. The RFCMAC model initially has no hypercube cell nodes except the input–output nodes. According to the reception of online incoming training data in the structure and parameter learning processes, the nodes from layer 2 to layer 5 are created automatically. Parameters $R_{kj}^q$ and $a_{ij}$ in the initial model are randomly generated between 0 and 1.

Fig. 4

Flowchart of the proposed structure and parameter learning.

Structure learning algorithm

Generally, the main purpose of structured learning is to determine whether a new hypercube cell should be extracted from the training data. For each incoming pattern $x_i$, the firing strength of the spatial firing layer can be defined as the degree to which the incoming pattern belongs to the corresponding cluster. The entropy measure is used to estimate the distance between each data point and each membership function. Entropy values between data points and current membership functions were calculated to determine whether to add a new hypercube cell. The entropy measure can be calculated using the firing strength from $u_{ij}^{(3)}$ as follows:

Eq. (8)

$${\mathrm{EM}}_j=-\sum _{i=1}^N{D}_{ij}{\log }_2{D}_{ij},$$

where $D_{ij}=\exp (u_{ij}^{{(2)}^{-1}})$ and ${\mathrm{EM}}_j\in [\mathrm{0,1}]$. Based on Eq. (9), the criterion for the degree measure is used to estimate and generate a new hypercube cell of new incoming data $\stackrel{\rightharpoonup }{\mathbf{x}}=(x_1,x_2,\dots ,x_{N_D})$. The maximum entropy measure is calculated as follows:

Eq. (9)

$${\mathrm{EM}}_{\max }=\underset{1\le j\le N_L}{\max }{\mathrm{EM}}_j,$$

where $N_L$ is the number of the hypercube cell and $\overline{\mathrm{EM}}\in [\mathrm{0,1}]$ is a prespecified threshold. In order to limit the number of hypercube cells in the proposed RFCMAC model, the threshold value will decay during the learning process. A low threshold leads to the learning of coarse clusters (i.e., a low number of hypercube cells are generated), whereas a high threshold leads to the learning of fine clusters (i.e., a high number of hypercube cells are generated). Therefore, the selection of the threshold value $\overline{\mathrm{EM}}$ critically affects the simulation results. That is, $\overline{\mathrm{EM}}$ determines whether the proper new hypercube cell is generated. Therefore, if ${\mathrm{EM}}_{\max }\le \overline{\mathrm{EM}}$, then a new hypercube cell is generated. Otherwise, the hypercube cell is not added.

Parameter learning algorithm

Five parameters of the model are entered in the learning algorithm and optimized based on the training data. The parameter learning occurs concurrently with the structure learning. For each piece of incoming data, five parameters (i.e., $m_{ij}$, ${\sigma }_{ij}$, $a_0$, $a_{ij}$, and ${\lambda }_{kj}^q$) are tuned in the RFCMAC model when the hypercube cells are newly generated or originally existed. Here, the gradient descent method is used to adjust the parameters of the receptive field functions and the TSK-type function. To clarify, consider the single-output case. The goal of the minimizing cost function $E$ is described as

Eq. (10)

$$E(t)=\frac{1}{2}{[y^d(t)-y(t)]}^2,$$

where $y^d(t)$ denotes the desired output and $y(t)$ is the model output for each discrete time $t$. In each training cycle, from the starting input variables to the activity of the model output, $y(t)$ are calculated by a feed-forward pass operation. According to Eq. (10), the error is used to regulate the weighted vector of the proposed RFCMAC model in a given number of training cycles. The well-known learning method of the backpropagation algorithm can be simplified as follows:

Eq. (11)

$$W(t+1)=W(t)+W(t)=\mathrm{\Delta }W(t)+[-\eta \frac{\partial E(t)}{\partial W(t)}],$$

where $\eta$ and $W$ represent the learning rate and the free parameters, respectively. $\eta$ denotes the pace factor for the learning rate in the search space. A low value may lead to a local optimal solution, whereas a high value leads to premature convergence that cannot obtain a better optimal solution. Therefore, the initial settings for $\overline{\alpha }$ and $\eta$ are based on experience estimation. According to Eq. (10), with respect to an arbitrary weight vector, $W$ is calculated by

Eq. (12)

$$\frac{\partial E(t)}{\partial W}=e(t)\frac{\partial y(t)}{\partial W}.$$

The corresponding antecedent and consequent parameters of the RFCMAC model are then adjusted using the chain rule to perform the error term recursive operation. With the RFCMAC model and the cost function as defined in Eq. (10), the update rule for $a_{ij}$ can be derived as

Eq. (13)

$$a_{ij}(t+1)=a_{ij}(t)+\mathrm{\Delta }{a}_{ij}(t),$$

where

Eq. (14)

$$\begin{gathered} a_{ij}(t)=-\eta ·\frac{\partial E}{\partial a_{ij}} \\ =-\eta ·\frac{\partial E}{\partial y}·\frac{\partial y}{\partial O_j^{(5)}}·\frac{\partial O_j^{(5)}}{\partial a_{ij}}. \end{gathered}$$

The equations used to update the recurrent weight parameter ${\lambda }_{kj}^q$ cell are

Eq. (15)

$${\lambda }_{kj}^q(t+1)={\lambda }_{kj}^q(t)+\mathrm{\Delta }{\lambda }_{kj}^q(t),$$

where

Eq. (16)

$$\begin{gathered} \mathrm{\Delta }{\lambda }_{kj}^q(t)=-\eta ·\frac{\partial E}{\partial {\lambda }_{kj}^q} \\ =-\eta ·\frac{\partial E}{\partial y}·\frac{\partial y}{\partial O_j^{(5)}}·\frac{\partial O_j^{(5)}}{\partial O_j^{(4)}}·\frac{\partial O_j^{(4)}}{\partial {\lambda }_{kj}^q} \\ =-\eta ·e·\frac{(a_{0j}+\sum _{i=1}^{N_D}{a}_{ij}{x}_i)\sum _{j=1}^{N_L}{O}_j^{(4)}-\sum _{j=1}^{N_L}{O}_j^{(4)}(a_{0j}+\sum _{i=1}^{N_D}{a}_{ij}{x}_i)}{{\left(\sum _j^{N_L}{O}_j^{(4)}\right)}^2}[O_j^{(4)}(t-1)-{\alpha }_j], \end{gathered}$$

where $\eta$ represents the learning rate of the recurrent $\lambda$ for the fuzzy weight functions and is set between 0 and 1, and $e$ denotes the error between the desired output and actual output, i.e., $y^d–y$.

$m_{ij}$ and ${\sigma }_{ij}$ represent the mean and variance of the receptive field functions, respectively. The adjustable parameters of the receptive field functions are calculated by

Eq. (17)

$$m_{ij}(t+1)=m_{ij}(t)+m_{ij}(t),$$

and

Eq. (18)

$${\sigma }_{ij}(t+1)={\sigma }_{ij}(t)+{\sigma }_{ij}(t),$$

where

Eq. (19)

$$\begin{gathered} \mathrm{\Delta }{m}_{ij}=-\eta ·\frac{\partial E}{\partial y}·\frac{\partial y}{\partial O_j^{(5)}}·\frac{\partial O_j^{(5)}}{\partial O_j^{(4)}}·\frac{\partial O_j^{(4)}}{\partial O_j^{(3)}}·\frac{\partial O_j^{(3)}}{\partial O_j^{(2)}}·\frac{\partial O_j^{(2)}}{\partial m_{ij}} \\ =-\eta ·e·\frac{(a_{0j}+\sum _{i=1}^{N_D}{a}_{ij}{x}_i)\sum _{j=1}^{N_L}{O}_j^{(4)}-\sum _{j=1}^{N_L}{O}_j^{(4)}(a_{0j}+\sum _{i=1}^{N_D}{a}_{ij}{x}_i)}{{\left(\sum _j^{N_L}{O}_j^{(4)}\right)}^2}·(1-{\gamma }_j^q)·{\alpha }_j·\frac{2(u_i^{(1)}-m_{ij})}{{\sigma }_{ij}^2}, \end{gathered}$$

and

Eq. (20)

$$\begin{gathered} \mathrm{\Delta }{\sigma }_{ij}=-\eta ·\frac{\partial E}{\partial y}·\frac{\partial y}{\partial O_j^{(5)}}·\frac{\partial O_j^{(5)}}{\partial O_j^{(4)}}·\frac{\partial O_j^{(4)}}{\partial O_j^{(3)}}·\frac{\partial O_j^{(3)}}{\partial O_j^{(2)}}·\frac{\partial O_j^{(2)}}{\partial {\sigma }_{ij}} \\ =-\eta ·e·\frac{(a_{0j}+\sum _{i=1}^{N_D}{a}_{ij}{x}_i)\sum _{j=1}^{N_L}{O}_j^{(4)}-\sum _{j=1}^{N_L}{O}_j^{(4)}(a_{0j}+\sum _{i=1}^{N_D}{a}_{ij}{x}_i)}{{\left(\sum _j^{N_L}{O}_j^{(4)}\right)}^2}·(1-{\gamma }_j^q)·{\alpha }_j·\frac{2{(u_i^{(1)}-m_{ij})}^2}{{\sigma }_{ij}^3}, \end{gathered}$$

where $i$ denotes the $i$’th input dimension for $i=\mathrm{1,2},\dots ,n$, and $j$ denotes the $j$’th hypercube cell.

3.2.

Weighted Strategy for Adaptively Refining the Transmission Map

In the real world, the transmission is not always constant within a window, especially around the contour of an object. In the inconstant regions, the image of the recovery scene generates some halos and block artifacts. The proposed solution is to use a pixel-window ratio (PWR) method to detect the possible regions of the halo artifact in the recovered scene and to use an adaptive weighting technique to mitigate the artifact. The PWR is defined as the ratio of the pixel itself and the $7\times 7$ mask of the window. The PWR is derived by

Eq. (21)

$$\mathrm{PWR}=\frac{\mathrm{PTM}}{\mathrm{WTM}},$$

where the numerator is the minimum channel by $1\times 1$ mask for RGB color space and the denominator is the window transmission map (WTM) by $7\times 7$ mask. A PWR value very close to 1 means that the transmission within the WTM is nearly constant. Although the halo situation cannot occur, the relative color saturation in the image is very high. In contrast, if the value of PWR is far $>1$, this means that the transmission within the window is inconstant and the halo artifact will occur. However, excessive color saturation is not a problem. Although the halo artifact region can be found by the value of PWR, the main problem is how to mitigate artifacts from these regions. The proposed solution to this problem is to use a weighted strategy approach to improve refinement of the transmission map and to mitigate the halo artifact. The weighted strategy approach is defined as follows:

Eq. (22)

$$t=\{\begin{array}{ll}\omega \times \left[\right(1-\frac{\alpha }{\mathrm{PWR}})\times \mathrm{PTM}+\frac{\alpha }{\mathrm{PWR}}\times \mathrm{WTM}],& \text{if}\mathrm{PWR}>T_{\text{upper}}\\ \omega \times \left[\right(1-\frac{\beta }{\mathrm{PWR}})\times \mathrm{PTM}+\frac{\beta }{\mathrm{PWR}}\times \mathrm{WTM}],& \text{if}{T}_{\text{lower}}<\mathrm{PWR}\le T_{\text{upper}}\\ \omega \times \mathrm{WTM},& \text{otherwise}\end{array}$$

where $\alpha$ and $\beta$ are the weighting factors for mitigating the artifacts. The range of $\alpha$ and $\beta$ is set as $0<\alpha <\beta <1$. In Eq. (22), if the PWR value is greater than $T_{\text{upper}}$, it means that the transmission is greatly different from the WTM, the weighting estimation of the WTM is decreased, and the weighting estimation of the PTM is increased. Therefore, this situation requires a very small weighting factor $\alpha$ to adjust the transmission rapidly so that the halo artifact can be eliminated. If the PWR value is between $T_{\text{upper}}$ and $T_{\text{lower}}$, this means that the transmission is a little different from the WTM. For this situation, the weighting factor $\beta$ is greater than the weighting factor $\alpha$ and it is applied to adjust the transmission smoothly. Otherwise, the WTM value is directly used as an estimation value. Parameter values $\alpha$ and $\beta$ are based on computational analysis of the intensity values associated with the halos. Figure 5(a) shows the original hazy image and Figs. 5(b)–5(j) show the results using the different values of $\alpha$ and $\beta$. Based on the above computational analysis, weighting factors $\alpha$ and $\beta$ are set appropriately to improve the quality of the image dehazing.

Fig. 5

(a) Original haze image; (b)–(j) the results using different $\alpha$ and $\beta$ values, where (b) $\alpha =0.1$, $\beta =0.1$; (c) $\alpha =0.1$, $\beta =0.5$; (d) $\alpha =0.1$, $\beta =0.9$; (e) $\alpha =0.5$, $\beta =0.1$; (f) $\alpha =0.5$, $\beta =0.5$; (g) $\alpha =0.5$, $\beta =0.9$; (h) $\alpha =0.9$, $\beta =0.1$; (i) $\alpha =0.9$, $\beta =0.5$; and (j) $\alpha =0.9$, $\beta =0.9$.

3.3.

Atmospheric Light Estimation

The atmospheric light factor must be carefully selected for effective image dehazing. An incorrectly selected atmospheric light factor will obtain very poor dehazing results. In some situations, many objects are considered atmospheric light, which results in erroneous image restoration. To solve this problem, the proposed solution is to use an average value of the brightest 1% in the transmission $t$ to refine the atmospheric light level. The average value is calculated as follows:

Eq. (23)

$$A^c=\frac{{\sum }_{\text{pixel}\in x}{\text{pixel}}^c}{|x|},$$

where $A$ is the atmospheric light and $c$ is the color channel. Figure 6 shows the results of scene radiance recovery.

Fig. 6

Estimation using an average value: (a) original image; (b) estimate of transmission map; (c) image of atmospheric light; and (d) scene radiance recovery.

4.

Results and Discussion

The experiments were performed in the C language on a Pentium(R) i7-3770 CPU @3.20 GHz. The effectiveness and robustness of the proposed method were verified by testing several hazy images, namely, “New York,” “ny12,” “ny17,” “y01,” and “y16”. The proposed approach was also compared with other well-known haze removal methods.^{13,16,17–20,24} Performance testing was divided into three parts: (1) results of removing the halo, (2) assessment of the visual images, and (3) analysis of the quantitative measurement.

4.1.

Results of Removing the Halo

Figure 7 shows the results of removing the halo for different images. In Fig. 7(a), the estimated transmission map is from an input hazy image using the patch size $7\times 7$. Although the dehazing results are good, some block effects (halo artifacts) exist in the blue blocks of Fig. 7(a). The phenomenon is because the transmission is not always a constant value in a patch. In Fig. 7(b), the halo artifacts are suppressed by the proposed method in the red blocks. Therefore, the halo artifacts do not exist using the proposed method.

Fig. 7

Removal of halo artifacts for different images. (a) Halo artifacts and (b) removal of the halo artifacts.

4.2.

Estimation of the Visual Image

Figure 8 shows the comparison results. This figure shows that the dehazing results obtained by the proposed method are better than those of Fattal,¹⁶ Tarel and Hautière,¹⁹ and Ancuti and Ancuti.²⁴ Additionally, Schechner and Averbuch¹⁴ adopted a multi-image polarization-based dehazing method that employs the worst and the best polarization states among the existing image versions. For a comparison with the method developed in Schechner and Averbuch,¹⁴ we processed only one input used in that study.r14 The dehazing results obtained by the proposed method are superior to those of Schechner and Averbuch.¹⁴

Fig. 8

Comparison of dehazing results using various methods.

Figures 9 and 10 also show the comparison results for the proposed approach and other state-of-the-art methods. Figure 9 shows that, compared with the techniques developed by Tan¹⁷ and by Tarel and Hautière,¹⁹ the proposed method preserves the fine transitions in the hazy regions and does not generate unpleasing artifacts. Moreover, the techniques of Tan¹⁷ and Tarel and Hautière¹⁹ produce oversaturated colors. Although the technique developed by Fattal¹⁶ obtains good dehazing results, its applications are limited in dense haze situations. The poor performance mainly results from the use of a statistical analysis method that needs to estimate the variance of the depth map. The technique of Kopf et al.¹³ obtains a good result in the color contrast, but only a little detailed texture is presented in the image. The technique of He et al.¹⁸ gets an obvious color difference in some regions. Recently, the technique developed by Nishino et al.²⁰ yields aesthetically pleasing results, but some artifacts are introduced in those regions, which are considered at infinite depth. The method developed by Ancuti and Ancuti²⁴ obtains a natural image, but color differences are visible in some regions, such as objects. The proposed method can effectively perform hazing, halation, and color cast.

Fig. 9

Comparison of dehazing techniques for city scene images: (a) ny12 and (b) ny17.

Fig. 10

Comparison of dehazing techniques for mountain scene images: (a) y01 and (b) y16.

An image of a mountain was also used for comparison with other state-of-the-art methods. Figure 10 shows the dehazing results using various methods. Comparisons showed that the Tan method¹⁷ produces oversaturation phenomena and causes color differences and halo artifacts. A good color contrast is obtained by the Fattal¹⁶ method, but some differences in detailed textures and color differences are visible. The results of Kopf¹³ are similar to those of Fattal.¹⁶ Though Tarel and Hautière’s¹⁹ method has a good detailed texture, the color difference problem is generated. Because of color differences caused by oversaturation, the results obtained by the He et al.¹⁸ method are unnatural. The technique developed by Nishino et al.²⁰ obtains a good overall image, but an unnatural phenomenon is visible in clouds in the sky. The technique of Ancuti and Ancuti²⁴ performs well in terms of true color contrast; however, a slight unnatural phenomenon still occurs around the sky. Overall, the results obtained by the proposed method are superior to those of other methods.

5.

Conclusions

The hybrid RFCMAC model and weighted strategy developed in this study effectively solve hazy and foggy images. The proposed RFCMAC model performs estimation of the transmission map and accurately selects the average of the brightest 1% of atmospheric light. An adaptively weighted strategy is applied to generate a refined transmission map for removing the halo effect. Experimental results demonstrate the superiority of the proposed method in enhancing color contrast, balancing color saturation, removing halo artifacts, and reducing computation time.