Haze removal for a single visible remote sensing image

Signal Processing 137 (2017) 33–43

Contents lists available at ScienceDirect

Signal Processing journal homepage: www.elsevier.com/locate/sigpro

Haze removal for a single visible remote sensing image Qi Liu a, Xinbo Gao b,1,∗, Lihuo He a, Wen Lu a a b

School of Electronic Engineering, Xidian University, Xi’an 710071, China State Key Laboratory of Integrated Services Networks, School of Electronic Engineering, Xidian University, Xi’an 710071, China

a r t i c l e

i n f o

Article history: Received 15 July 2016 Revised 25 January 2017 Accepted 31 January 2017 Available online 1 February 2017 Keywords: Haze removal Visible remote sensing image Haze thickness map Total variation inpainting

a b s t r a c t Satellite remote sensing image often suffers from haze degradation, which deteriorates significantly the effect of data intelligibility and interpretability. Hence, haze removal technique is essential for inferring high quality images with clear visibility to improve the application value of the satellite images. Thin haze removal is a challenging task since the estimation of haze component is easily affected by ground features. To solve the problem, this paper develops an effective haze removal method for a single visible remote sensing image. Firstly, haze is considered as an additive contamination and can be represented by a haze thickness map (HTM). A ground radiance suppressed HTM (GRS-HTM) is then proposed for a more precise estimation of haze distribution. The haze component for each band is calculated via GRS-HTM and can be removed to recover the clear image. Several visible satellite images with different resolutions were tested to validate the effectiveness of the proposed method. The evaluation results with qualitative and quantitative assessments demonstrate that the proposed method is superior to the traditional methods, and can recover a haze-free image with high quality. © 2017 Elsevier B.V. All rights reserved.

1. Introduction The development of remote sensing technology in recent years has led to great increase in the quantity of satellite images. Many remote sensing applications benefit from high resolution images due to the better representability of ground features [1–4]. Unfortunately, the remote sensing imaging is sensitive to atmospheric conditions during the data acquisition process. In many cases, satellite images are affected by the haze which degrades the image quality by obscuring the clarity and reducing the visibility of the scene. However the ground feature is still recognizable, which leaves an opportunity for ground surface radiance recovery [5]. The visible images contaminated by haze lose contrast and color fidelity, which seriously affects their interpretation and usage for humans. Moreover, the performance of many image processing algorithms will inevitably suffer from the degraded scene radiance. Haze removal not only recovers a high quality of the visible images, but also improves the capability for practical usage. It is highly desired in remote sensing image applications. Haze removal for visible remote sensing images is a challenging task since images in haze regions contain both earth surface features and haze information. The ground radiation is mixed with haze scattering light, which makes the problem generally difficult ∗

1

Corresponding author. E-mail addresses: [email protected], [email protected] (X. Gao). Xinbo Gao is a member of Editorial board for EURASIP Signal Processing.

http://dx.doi.org/10.1016/j.sigpro.2017.01.036 0165-1684/© 2017 Elsevier B.V. All rights reserved.

to solve given a single haze image. Some single-image-dehazing techniques can handle out-door hazy images successfully [6–8]. However, these methods have not been proved to be effective for remote sensing images. Recently, some progresses have been made in the research on remote sensing image dehazing, which can be divided into two categories: image enhancement based methods and image restoration based methods. The image enhancement based methods use some assumptions to recover clear scenes from hazy images. Richter [9] proposed a histogram matching method based on the assumption that histograms of clear and hazy regions are the same. In [10], haze effects are assumed to reside in low-frequency parts of the image. Wavelet transform fusion is used to replace the data in this part with a reference haze-free image. Homomorphic filter based methods [11] use the same assumption. However, the low-frequency layer might also contain ground feature. Most of these methods are mainly concentrated on low-resolution MS data such as TM and ETM+ images. The assumption used in these methods may not valid for the high resolution images. Besides, the global operations can lead to the distortion of the clear regions of hazy image inevitably. Image restoration based methods aim at recovering true ground information by correcting hazy pixels based on a imaging model. The dark-object subtraction (DOS) technique is a famous atmospheric correction method that can remove homogeneous haze [12]. Zhang et al. [13] proposed a haze optimized transformation (HOT) to detect spatially varying haze. It combined with DOS

34

Q. Liu et al. / Signal Processing 137 (2017) 33–43

which has been shown as an effect way for TM/ETM+ data. However, the method usually fails over snow, water, and bright surfaces. It can induce high values that result in over-correction of these targets [14]. Several strategies have been developed to make HOT method more suitable for high resolution images [15–17]. Since more-detailed confusing land-cover can be distinguished, the good sample selection for HOT and adjustments for confused landcover are tedious tasks. Liu et al. [18] proposed a thin cloud removal method based on a cloud physical model [19]. The model combined with the dark channel prior[8] has also been used in single image dehazing recently [20–22]. However, the atmospheric light is hard to calculate for satellite images and the model is not so effective for uneven haze. Makarau et al. [5] developed a haze thickness map (HTM) to represent uneven haze. The HTM is a development of the DOS technique that searches dark objects locally in whole image. The major problem of the HTM calculation is the overestimation caused by bright ground objects. The ground radiance upraises the HTM value which can lead to over-dehazing of the image. In [5], The bright objects are detected using the threshold segmentation and the HTM values in the regions of bright objects are interpolated using a triangulation. However, the threshold segmentation cannot distinguish between the thick haze regions and the bright ground objects. It will leads to residual haze or overcorrection of ground objects in the results. In this paper, a new haze removal method for visible remote sensing images is presented. The main purpose is to construct a more precise HTM that suffers from a minimal influence of ground brightness. An initial haze detection is given by a HTM. Then a ground radiance suppressed HTM (GRS-HTM) is proposed to calculate a fine haze distribution. Due to the different influences of haze on visible bands, the haze component for each band is calculated via GRS-HTM and finally removed to recover the haze-free image. The proposed GRS-HTM for precise haze detection is the most important step which is also the main contribution of this method. It is based on a physical observation that the distribution of haze thickness is different from that of the ground radiance. The land surface has clear outline whereas the haze is without obvious boundaries. It can be assumed that the haze thickness has a smooth change at edges of ground features. Based on this idea, we propose to refine the HTM by suppressing ground radiance in a local patch centered at the edge point of the ground. The proposed GRS-HTM can eliminate the impact of the ground brightness hence avoid the over-dehazing in traditional HTM method [5]. Several visible images with different resolutions were collected and tested to validate the effectiveness of the proposed method. The evaluation results with qualitative and quantitative assessments demonstrate that the proposed method is superior to the traditional methods, and can recover a haze-free image with high quality. The rest of this paper is organized as follows. Section 2 introduces some background knowledge, including the haze model and the HTM calculation. Section 3 presents our proposed method for haze removal. Section 4 gives a detailed description of test images, results and evaluations. Finally, Section 5 concludes the paper and the further work. 2. Background 2.1. Haze imaging model In the literatures on dehazing remote sensing image, haze is considered as an additive component to the acquired radiance signal at the sensor of the remote sensing systems. Moreover, the additive scatter light depends on the channel wavelength [23]. The atmospheric scattering in the blue band is the most severe, and the scatter light decreases with the increase in the wavelength. Taking a WorldView-2 image for example in Fig. 1, we can find that the

Fig. 1. An example of haze effects on different bands of visible image. (a) RGB color band combination. (b) Blue band. (c) Green band. (d) Red band. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

influence of the haze is different in the visible bands. So this imaging model can be described as

Ic (x, y ) = Jc (x, y ) + Hc (x, y ),

(1)

where (x, y) are the coordinates in an image, I stands for the observed hazy image, J is the haze-free image, H is the haze component which depends on the channel wavelength, and c denotes the index of bands. The clear scene radiance can be restored if the haze components Hc are known. 2.2. Dehazing using HTM Many approaches have been proposed to estimate the haze components. In the DOS method [12], the haze signal in each band is assumed to be a constant value corresponding to the dark pixels searched in the whole scene. The values of these dark pixels are usually very low and close to zero, in the case of hazy condition the values of these dark pixels are mainly contributed by the haze scatter light. So the dark pixels can be used to estimate the haze thickness. Instead of global searching, the HTM method [5] searches dark objects within a local window in the image. It provides a better representation for the spatially uneven haze and can be used to build precise haze components. In this paper, we also use the HTM suggested in [5] as an initial estimation of haze distribution. It is calculated by searching the darkest pixels using a local non-overlapping window in the whole image, and then resized to the original band size

H (x, y ) =

min

(x ,y )∈W (x,y )

(Ib (x , y )),

(2)

where W(x, y) is a window of size w × w centered at (x, y), Ib (x, y) is the blue spectral band which is most affected by haze. The employment of a blue band to calculate HTM can lead to overdehazing of this band. In [5], the blue band is replaced by a synthetic band which is defined by linear extrapolation of two bands:

Q. Liu et al. / Signal Processing 137 (2017) 33–43

35

Fig. 2. An example of the HTM calculation. (a) Original hazy image. (b) HTM calculated on the original blue band. (c) New synthetic band extrapolated using blue and green bands. (d) HTM calculated on the synthetic band in (c). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

Isyn (x, y ) = 2 · Ib (x, y ) − 0.95 · Ig (x, y ). Isyn (x, y) is the synthetic band, and Ib (x, y) and Ig (x, y) are the blue and green bands respectively. The use of the synthetic band is to alleviate the influence of ground objects’ brightness. We follow this method to use the synthetic band. As we can see in Fig. 2(b), the HTM built on the original blue band contains many ground textures. While the situation is better on the HTM in Fig. 2(d) calculated using the synthetic band shown in Fig. 2(c). However, there are still some bright objects shown in Fig. 2(d) and HTM values in these areas are higher than real haze thickness. It is difficult to capture the local dark pixels using a fixed window since the ground objects are always varied in size. No matter how large the window is, the local search will suffer from ground radiance inevitably. Large window will breach the structure of haze and leads to an underestimated HTM. A smaller window allows better estimation of haze structure and the influence of the ground radiance can be eliminated in the subsequent processing [5]. So a smaller window is more suitable for creating HTM. The major problem of the initial HTM is overestimation cause by ground reflected light. To correct the ground bright object, the traditional HTM based method [5] labeled them using a threshold segmentation of the red band. However, the simple threshold segmentation cannot handle thick haze well. Since the pixels in thick haze regions also have high intensity value and can easily be mislabeled as bright objects. This will lead to the residual haze in the final results. Increasing the threshold will result in smaller detected region and hence can remove thick haze eventually. But it will cause distortion in the undetected bright regions. To obtain better dehazing results, we should calculate a HTM with a minimal influence of ground radiance. 3. Haze removal based on a precise HTM In this section a reliable haze component detection and estimation method is proposed and the dehazing procedure is described thoroughly. 3.1. Ground radiance suppressed HTM According to Mie theory, haze is mainly generated by large particles which are evenly distributed above the ground in the atmosphere [24]. The spatial variation of haze distribution is smooth since its propagation is slower. So the haze component is mainly located in the low-frequency layer of the image [10]. However, it can be seen from Fig. 2(d) that the original HTM calculated by searching local minima still contains excessive small textures. Some objects with high reflectance cause very high HTM values than the real haze thickness. Besides, the ground objects bring false edges and shapes on the HTM. These features can be treated as noise and should be removed to get a smooth HTM. To avoid the blending of the ground brightness in haze component, we propose a ground radiance suppressed HTM which contains two steps. The first step is to remove the textures and high reflectance objects

in the original HTM. The high reflectance objects are usually small isolate pixels or patches, such as blue roofs, snow, white buildings and etc. Different form the method in [5], we do not intend to segment and remove all bright objects. This is because the moderately bright objects like bare soils shown in Fig. 2(d) have a large size. It is hard to get an accurate correction for them. In this step we only correct small highly bright patches and leave the large object in the following procedure. The next step is to suppress reflected radiance for large objects with moderately high intensity on the corrected HTM. We do not try to divide the bright objects and haze regions into two parts like the method in [5]. Instead, we use edge information to eliminate influence of the ground radiance. 3.1.1. Textures and highly bright patches removal To remove textures and highly bright patches in the original HTM, an image inpainting model [25] based on the total variation regularization is employed. Inpainting is the process of reconstructing lost or deteriorated parts of images, which can be used to remove small defects and replace the corrupted parts of the image data [25]. The highly bright patches can be treated as abnormal HTM values. So they can be reconstructed by using the inpainting technique. Additionally, the total variation regularization has been proved effective in removing unwanted details and textures [26]. According to the total variation inpainting model, a smoothed HTM can be obtained by solving the following minimization problem

hˆ = arg min h

 i

|∇ hi | + μ

  h j − x j 2 ,

(3)

j∈N

where h is the updated HTM, x is the HTM obtained by local minimum search, i is the pixel index, ∇ is the gradient operator, μ is the positive parameter and N is the index set of the pixels except those in highly bright patches. The first term of the objective function in Eq. (3) is the regularization term, which controls the smoothness of the reconstructed HTM. The second item of the objective function is data error term measures the fidelity of the output with the initial HTM. It can be found from Eq. (3) that the minimization problem is equal to the total variation denoising model [26] outside the bright patch areas. The parameter μ can balance the smoothness and data fidelity of the result so as to removal the excessive details on the initial HTM. For the pixels of labeled patches, the data error term in energy function is not is unused and h is only influenced by the smoothness term. In this way the labeled regions can be reconstructed. A first-order method suggested in [27] is used to solve the minimization problem. Fig. 3(c) shows the corrected HTM in which the small textures have been removed and highly bright objects have been corrected. The highly bright patches are detected on the synthetic band using a threshold. The threshold should not be large due to two reasons. First, the pixels in thick haze also have high values. The thick haze regions can easily be mislabeled as bright objects. The large threshold will mislabel the thick haze areas and cause wrong correction for them. Fig. 4 shows the HTM correction results using different

36

Q. Liu et al. / Signal Processing 137 (2017) 33–43

Fig. 5. An illustration of the difference between haze component and ground features. (a) Hazy blue band image. (b) Subset image in the upper red rectangle. (c) Subset image in the blue rectangle underneath. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

Fig. 3. An illustration of textures and bright patches removal for original HTM. (a) Original HTM established on synthetic band. (b) Mask for highly bright patches. (c) Corrected HTM.

patch at edges of ground features. Based on this idea, we propose to suppress ground reflected radiance using edge information. The Canny edge detector [28] is used to detect edges on the synthetic band image. In a local patch around an edge point, the HTM values should be constants or have a small variation. The low values in this patch are closer to the real haze thickness than higher values which should be diminished. Due to the presence of noise, it is not suitable to use the minimum of the patch to represent the haze thickness. Besides, we just hope to depress the higher values. So we calculate the mean value as the upper bound of the patch. The suppression is performed on the smoothed HTM obtained in the previous step. The highly bright patches should be removed in the smoothed HTM because it is hard to reduce the high values using a mean value. The suppression operation can be expressed as follows:

thresholds. It can be seen that the mask obtained by threshold 3% contains haze regions, and the corrected HTM in Fig. 4(e) has a lower value than that in Fig. 4(d). The other reason is that more labeled regions will bring more computational cost for the correction operation. In this paper, the top 1% brightest pixels in the synthetic band are labeled. It can also be set by users according to the land cover type of the image. 3.1.2. Suppress ground reflected radiance for large objects For large flat land surfaces in haze areas, the HTM values are contributed by both ground radiance and haze scattering light. Without prior information, it is hard to separate the two components. Actually the two signals are quite different in pattern and property. The haze component is a low-frequency signal due to its physical property whereas the component from ground radiance contains a variety of textures and structures. The blue band image in Fig. 1(b) is used to show the difference between them as shown in Fig. 5. The enlarged regions from red and blue rectangles contain borders of haze area and ground features respectively. It can be found that the land surface has clear outline whereas the haze is without evident boundaries. Around the edges of the image, the ground surface brightness changes dramatically whereas the haze thickness varies slowly. From both theoretical analysis and practical observation, we can assume that the haze thickness has a smooth change in a local

(a)

hm =

1  Hs (x, y ), N

(4)

(x,y )∈

Hg (x, y ) = min (Hs (x, y ), hm ),

∀(x, y ) ∈ 

(5)

where  is a local patch centered at an edge point, N is the pixel numbers of the patch, Hs is the smoothed HTM, Hg is a GRS-HTM. The whole procedure of the ground radiance suppression is described in Algorithm 1. As we can see from the result of this refine process in Fig. 6, the ground features brightness such as roads, roofs of buildings

(b)

(c)

(d)

(e)

Fig. 4. HTM correction results using different thresholds. (a) Synthetic band. (b) Top 1% brightest pixels in the Synthetic band. (c) Top 3% brightest pixels in the synthetic band. (d) Corrected HTM using the mask in (b). (e) Corrected HTM using the mask in (c).

Q. Liu et al. / Signal Processing 137 (2017) 33–43

37

3.3. Dehazing using GRS-HTM

Algorithm 1 The GRS-HTM calculation. Input: Smoothed HTM Hs , patch size ω × ω, edge point set P = { ( xi , yi )|i = 1, 2, . . . , m} Initialize: Hg = Hs For i = 1 to m a) Select patch centered at point (xi , yi ); b) Compute mean of the patch by Eq. (4); c) Update pixel values for Hg in the patch by Eq. (5);

In practice, even on clear days the atmosphere is not absolutely free of any particle. The HTM includes both the thickness and aerosol thickness. If we directly subtract the HTM from the original image using Eq. (1), the image will be over-dehazed and the spectral profiles in the non-haze area will change [5]. The aerosol thickness is simply excluded by subtracting the minimum value of the HTM. So the dehazing is performed using the following equation



End for Output: GRS-HTM Hg



Jc (x, y ) = Ic (x, y ) − Hc (x, y ) − min Hc .

(10)

Putting Eq. (6) into Eq. (10), we can get and others are suppressed. The main problem of this patch-wise operation is the block artifacts. A guided filtering [29] is applied to remove the block artifacts and get a final smooth GRS-HTM.





Jc (x, y ) = Ic (x, y ) − kc · Hg (x, y ) − min Hg .

(11)

After the GRS-HTM and coefficients for each band are obtained, we can recover the haze-free image using Eq. (11).

3.2. HTM calculation for each band 4. Experiments and evaluation According to the imaging model in Eq. (1), the additive haze components are different for each spectral band. So the dehazing for a visible image is to subtract a band-specific HTM. The haze component for each band is estimated relative to the GRS-HTM based on the linear relationship between them [5]:

Hc (x, y ) = kc · Hg (x, y ),

(6)

where kc is the linear coefficient. The coefficient kc can be calculated through a linear regression. The regression needs input and output pixel values. We use pixels values of GRS-HTM as the independent variable. A temporary HTM calculated using the original band was assumed as the dependent variable. Not all pixel values of the two HTMs are involved in regression. As we discussed above, the temporary HTM calculated using the original band is not accurate due to the influence of the ground radiance. The higher value pixels in the temporary HTM have more chance to be influenced, whereas the lower values are most likely represented for the real haze thickness. So the pixels are selected in a more reliable way. The mean value of the HTM is used as a threshold to decide whether a pixel can be used. Two point sets are established by thresholding the temporary HTM and GRS-HTM. The final valid point set is the intersection of the two sets



  1   P1 = (x, y )Ht (x, y ) < Ht (i, j ) , L 

(7)

(i, j )∈Ht

  1   P2 = (x, y )Hg (x, y ) < Hg (i, j ) , L

(8)

(i, j )∈Hg

Pv = P1 ∩ P2 ,

(9)

where L is the pixel number of the image, Pv is the final valid point set. The points in this set can then be used to calculate a more precise regression coefficient.

To validate the effectiveness of the proposed method, it was tested on several hazy visible images from two kinds of satellite data, including four Landsat-8 images and two Worldview-2 images. The dataset contains various land-cover types in both urban and rural area, including vegetated areas, water bodies, bare soil and urban targets. The performance of the proposed method is also compared with the traditional HTM method [5] and HOT method [13] both qualitatively and quantitatively. The traditional HTM method [5] detected the large bright objects using threshold segmentation on the red band, and the threshold value is equal to the mean of the red band. We consider two aspects to evaluate the dehazing results: the recovery of hazy area and the distortion of non-haze area. However, the evaluation on quality of the dehazed image is difficult since the ground truth data is not available. So the haze removal effects of different methods are compared visually first. For each Landsat-8 image, we also find a corresponding haze-free image as a reference. The haze-free image is taken at the same place of the hazy image with a minimal difference of acquisition time and conditions, which makes it reasonable for the comparison on the quality of the dehazed image and the true color of the dehazed areas. The distortion of the clear region is a negative side effect of the dehazing method. To measure the distortion of the non-haze area quantitatively, we use mean absolute error (MAE) to indicate the dissimilarity between the original pixel values and the dehazed version for different visible bands. A smaller MAE value usually indicates that the result has a higher fidelity in the clear region. The major parameters of the proposed method are set as follows. For all the images, the window size in the initial HTM calculation is set to 3 × 3. When calculating the temporary HTM for each visible band, we use a larger window with the size 11 × 11 which allows easier detection of dark pixels. The parameters used in TV inpainting are set to the default value suggested in [27]. The

Fig. 6. An example of the GRS-HTM construction. (a) Edge points detected in synthetic band. (b) Corrected HTM by removing textures and highly bright objects. (c) Ground radiance suppression for (b). (d) Final GRS-HTM smoothed by guided filter.

38

Q. Liu et al. / Signal Processing 137 (2017) 33–43

Fig. 7. Haze removal results of the Landsat-8 image. (a) Original hazy image. (b) Reference image. (c)–(e) Results of the HOT, the traditional HTM, and the proposed GRSHTM, respectively. (f)–(i) Enlarged subset images of (b)–(e) in the rectangles. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

parameters for GRS-HTM calculation are not same for the images of different resolution. For Landsat-8 images the patch size is set to 19 × 19. The threshold used in Canny edge detector is set to 0.1 in the experiments. The calculation for GRS-HTM needs edge information. If the edges cannot be detected under haze area, the brightness of ground features will not be suppressed and lead to over-estimate of haze thickness. The threshold should be given a proper value to find adequate edges on synthetic band. But too small value brings more operations for the GRS-HTM calculation and cannot further improve the results. So the threshold is tuned manually to a proper value for all the test images. For Worldview2 images the suppression operation is carried out using a larger patch with the size 25 × 25, and the threshold is set to 0.05 for edge detection. All the experiments are conducted on a desktop PC with 3.2 GHz Intel Core i5 CPU and 8GB memory, using MATLAB 2012a for implementation of the algorithms. 4.1. Qualitative evaluation Fig. 7 shows the dehazing results for a Landsat-8 image. The image has a size of 1200 × 1200 pixels with a spatial resolution of 30m. This scene contains vegetated targets and bare soil in both hazy and clear areas. It can be observed that the proposed method provides the best dehazing result compared to the other methods. The HOT method has little dehazing effect on this image, and the result remains hazy as shown in Fig. 7(c). Both of the traditional HTM method and the proposed method improve the quality of the hazy image visually. However, the result of the traditional HTM method still remains a little hazy. For a better visibility, the subset

images in the rectangles of Fig. 7(b)–(e) are shown in Fig. 7(f)–(i), respectively. A further inspection on the enlarged views of the subsets demonstrates that the result of the proposed method is more natural and more consistent with the reference image. Fig. 8 shows the other three Landsat-8 images, the reference images, and the corresponding haze removal results of the different methods, respectively. The three scenes have the same size and spatial resolution with the image in Fig. 7(a). They are all covered by spatially varying haze and contain diverse urban and suburban land cover types under the haze. The quality of the dehazed images can be assessed by comparing them with reference images. The results of HOT method have a serious color cast in the hazy area where the color tends to be white as shown in Fig. 8(g) and (h). Besides, there are also distortions in the thick haze area which can be seen in the red rectangle regions. The results of the traditional HTM method suffer little from color distortion. However, it cannot remove dense haze as shown in the yellow rectangles in Fig. 8(k) and (l). It is due to the underestimate of the HTM values for thick haze that the pixels in thick haze regions also have high intensity value and are mislabeled as bright objects. So the HTM values in thick haze regions will be depressed in the triangulation operation. The proposed method can remove the thick haze effectively while keep the spectral feature of the original non-haze area. The results exhibit the best color and contrast which are more consistent with the reference images compared with the results of the other methods. Figs. 9 and 10 show the dehazing results on two Worldview2 satellite images with a spatial resolution of 0.5m. The size of the original hazy images shown in Figs. 9 and 10 are 800 × 800

Q. Liu et al. / Signal Processing 137 (2017) 33–43

39

Fig. 8. Haze removal results of the other three Landsat-8 images. The first row shows the original hazy images. The second row is the corresponding reference images. (g)–(i) Results of the HOT. (j)–(l) Results of the traditional HTM. (m)–(o) Results of the proposed GRS-HTM.

and 1600 × 1600, respectively. Both the two images contain many bright objects on the ground. Though all the three methods improve the quality of the hazy images, a visual inspection of image details shows that our method gives better results. The results of HOT method and traditional HTM method suffers from color distortions as shown in the red rectangles of Figs. 9(b)(c) and 10(b)(c). Most of the blue roofs become black in the results of the HOT method. This artifact is due to the over-estimation of the HOT levels and subtraction of a higher value in these areas on the blue band images. The proposed method is robust against the ground

brightness thus can avoid the color distortion. So the results are more natural and of higher quality. 4.2. Quantitative evaluation In order to calculate the MAE for the dehazing results, we labeled the clear regions manually in original images as shown surrounded by the red polygons in Fig. 11. The values of the metric for the different bands of the dehazing results are listed in Table 1. Due to the different quantization bits of the Landsat-8 images and

40

Q. Liu et al. / Signal Processing 137 (2017) 33–43

(a)

(b)

(c)

(d)

Fig. 9. Results of the Worldview-2 image. (a) Original hazy image. (b)–(d) Results of the HOT, the traditional HTM, and the proposed GRS-HTM, respectively.

Table 1 MAE values of the different dehazing results. Image bands

HOT [13]

HTM [5]

GRS-HTM

Fig. 7

24.67 15.06 4.83 22.58 15.92 19.18 11.90 7.05 3.68 36.71 33.54 31.45 5.24 6.29 5.23 0.31 0.34 0.12

59.22 54.48 49.75 101.62 87.25 79.12 70.69 70.65 64.99 295.67 284.43 260.78 11.95 10.99 9.52 4.16 3.99 2.89

13.07 11.90 10.72 15.25 12.10 10.39 4.22 4.01 3.55 25.23 22.88 20.54 2.43 1.77 1.45 0.21 0.20 0.19

Fig. 8(a)

Fig. 8(b)

Fig. 8(c)

Fig. 9

Fig. 10

Blue Green Red Blue Green Red Blue Green Red Blue Green Red Blue Green Red Blue Green Red

the Worldview-2 images, the results of two kinds of satellite images have different orders of magnitude. It can be observed that the proposed method has the smallest MAE value except on two

red bands. The traditional HTM results in the largest values of the metric, while the proposed GRS-HTM method has a significant reduction. This indicates that the proposed method can effectively detect the haze thickness and does not breach the spectral properties of clear region data. In a conclusion, the qualitative comparison agrees with the visual analysis, the proposed method provides the superior performance both subjectively and objectively. 4.3. Complexity experiment Suppose the number of bands for dehazing is p, the number of the pixels in each band is N. The computational complexity for obtaining the initial HTM is O(N2 /w2 ) with bicubic interpolation, where w is number of pixels of searching window. The complexity is O(k · N) for correcting HTM using the TV inpainting, where k is the number of iterations which depending on the missing pixels. Since we used a percent to control the number of missing pixels, the number of iterations will not be too large. The efficiency of Algorithm 1 is largely depending on the number of edge points Ne and the size of patch Np . It is a linear computational complexity of O(Ne · Np ). The complexity for calculating coefficient for band-specific HTM is O(p · Ns ), where Ns is the number of valid

Q. Liu et al. / Signal Processing 137 (2017) 33–43

41

Fig. 10. Results of the other Worldview-2 image. (a) Original hazy image. (b)–(d) Results of the HOT, the traditional HTM, and the proposed GRS-HTM, respectively. Table 2 Executive time(s) on different images. Images

HOT [13]

HTM [5]

GRS-HTM

Fig. Fig. Fig. Fig. Fig. Fig.

19.08 45.49 32.27 52.78 2.19 4.28

10.63 11.89 10.59 12.06 4.59 22.41

9.40 10.43 8.22 10.98 3.48 20.69

7 8(a) 8(b) 8(c) 9 10

pixels we selected for linear regression. Taking all the procedures into account, the total computational complexity for the proposed method is O(N2 /w2 ). Furthermore, we compared the executive time of the different methods. The results are shown in Table 2. It can be seen that the proposed method is comparative in most cases. 5. Conclusion In this paper, we presented an effective method for a single visible remote sensing image dehazing. The uneven haze is an additive component to the ground radiance and can be represented by a HTM. The original HTM calculation is prone to overestimation of haze thickness which can lead to over-dehazing of the image. For

this purpose, we proposed a GRS-HTM to provide a precise detection and estimation of the haze component. Spatially varying haze can be removed based on the GRS-HTM. Several visible satellite images with different resolutions were tested to validate the effectiveness of the proposed method. The performance of the method is also compared with the other representative methods both qualitatively and quantitatively. It turns out that our method can effectively eliminate uneven haze and improve image quality visually. The objective assessment metrics also demonstrate the outstanding performance of our method. The proposed method also has a limitation. The GRS-HTM calculation is performed at the edge point. When the edges cannot be detected around the high bright ground objects, especially in some large flat ground surface, the method will over-estimate the haze thickness values. The parameter in edge detection is set manually for different images to get adequate edge points. The detected points should not be too more since they increase the computational burden in the suppression operation. The window size used in GRS-HTM calculation is determined experimentally, which is also not a constant value for the images with different resolution. Future improvements of the method will deal with these defects and make the algorithm suitable for different types of land cover. In addition, more reasonable evaluation criteria should be considered to judge the dehazing performance quantitatively with the

42

Q. Liu et al. / Signal Processing 137 (2017) 33–43

Fig. 11. Clear regions labeled on the original hazy images.

reference images. The influence of the dehazing operation to the spectral properties of the ground objects is important and should be examined more effectively. Acknowledgements This work was supported in part by the National Natural Science Foundation of China under Grant 61432014, U1605252 and 61501349, in part by Key Industrial Innovation Chain in Industrial Domain under Grant 2016KTZDGY-02, in part by the Fundamental Research Funds for the Central Universities under Grant BDZ021403 and Grant JB149901, in part by National High-Level Talents Special Support Program of China under Grant CS3111720 0 0 01, in part by the Program for Changjiang Scholars and Innovative Research Team in University of China under Grant IRT13088. References [1] M.D. Mura, J.A. Benediktsson, B. Waske, L. Bruzzone, Morphological attribute profiles for the analysis of very high resolution images, IEEE Trans. Geosci. Remote Sensing 48 (10) (2010) 3747–3762. [2] D. Brunner, G. Lemoine, L. Bruzzone, Arthquake damage assessment of buildings using vhr optical and sar imagery, IEEE Trans. Geosci. Remote Sensing 48 (5) (2010) 2403–2420. [3] S. Balla-Arabe, X. Gao, B. Wang, F. Yang, V. Brost, Multi-kernel implicit curve evolution for selected texture region segmentation in vhr satellite images, IEEE Trans. Geosci. Remote Sensing 52 (8) (2014) 5183–5192. [4] F. Pacifici, F. Del Frate, Automatic change detection in very high resolution images with pulse-coupled neural networks, IEEE Geosci. Remote Sens. Lett. 7 (1) (2010) 58–62. [5] A. Makarau, R. Richter, R. Muller, P. Reinartz, Haze detection and removal in remotely sensed multispectral imagery, IEEE Trans. Geosci. Remote Sensing 52 (9) (2014) 5895–5905. [6] R. Fattal, Single image dehazing, ACM Trans. Graph. 27 (3) (2008) 72. [7] C.O. Ancuti, C. Ancuti, Single image dehazing by multi-scale fusion, IEEE Trans. Image Process. 22 (8) (2013) 3271–3282. [8] K. He, J. Sun, X. Tang, Single image haze removal using dark channel prior, IEEE Trans. Pattern Anal. Mach. Intell. 33 (12) (2011) 2341–2353. [9] R. Richter, Atmospheric correction of satellite data with haze removal including a haze/clear transition region, Comput. Geosci. 22 (6) (1996) 675–681. [10] Y. Du, B. Guindon, J. Cihlar, Haze detection and removal in high resolution satellite image with wavelet analysis, IEEE Trans. Geosci. Remote Sensing 40 (1) (2002) 210–217.

[11] H. Shen, H. Li, Y. Qian, L. Zhang, Q. Yuan, An effective thin cloud removal procedure for visible remote sensing images, ISPRS J. Photogramm. Remote Sens. 96 (11) (2014) 224–235. [12] P.S. Chavez, An improved dark-object subtraction technique for atmospheric scattering correction of multispectral data, Remote Sens. Environ. 24 (88) (1988) 459–479. [13] Y. Zhang, B. Guindon, J. Cihlar, An image transform to characterize and compensate for spatial variations in thin cloud contamination of landsat images, Remote Sens. Environ. 82 (2) (2002) 173–187. [14] Y. Zhang, B. Guindon, Quantitative assessment of a haze suppression methodology for satellite imagery: effect on land cover classification performance, IEEE Trans. Geosci. Remote Sensing 41 (5) (2003) 1082–1089. [15] G.D. Moro, L. Halounova, Haze removal for high-resolution satellite data: a case study, Int. J. Remote Sens. 28 (10) (2007) 2187–2205. [16] X.Y. He, J.B. Hu, W. Chen, X.Y. Li, Haze removal based on advanced haze optimized transformation (ahot) for multispectral imagery, Int. J. Remote Sens. 31 (20) (2010) 5331–5348. [17] C. Liu, J. Hu, Y. Lin, S. Wu, W. Huang, Haze detection, perfection and removal for high spatial resolution satellite imagery, Int. J. Remote Sens. 32 (23) (2011) 8685–8697. [18] J. Liu, X. Wang, M. Chen, S. Liu, X. Zhou, Z. Shao, P. Liu, Thin cloud removal from single satellite images, OPT. EXPRESS 22 (1) (2014). 618–32 [19] Z. Liu, B.R. Hunt, A new approach to removing cloud cover from satellite imagery, Comput. Vis. Graph. Image Process. 25 (2) (1984) 252–256. [20] X. Lan, L. Zhang, H. Shen, Q. Yuan, H. Li, Single image haze removal considering sensor blur and noise, EURASIP J. Adv. Signal Process. 23 (1) (2013) 1–13. [21] J. Long, Z. Shi, W. Tang, C. Zhang, Single remote sensing image dehazing, IEEE Geosci. Remote Sens. Lett. 11 (1) (2014) 59–63. [22] X. Pan, F. Xie, Z. Jiang, J. Yin, Haze removal for a single remote sensing image based on deformed haze imaging model, IEEE Geosci. Remote Sens. Lett. 22 (10) (2015) 1806–1810. [23] S.G. Narasimhan, S.K. Nayar, Vision and the atmosphere, Int. J. Comput. Vision 48 (3) (2002) 233–254. [24] W.C. Hinds, Aerosol Technology: Properties, Behavior, and Measurement of Airborne Particles, John Wiley, 2012. [25] T.F. Chan, J. Shen, Mathematical models for local nontexture inpaintings, SIAM J. Appl. Math. 62 (3) (2002) 1019–1043. [26] L.I. Rudin, S. Osher, E. Fatemi, Nonlinear total variation based noise removal algorithms, Physica D 60 (1) (1992) 259–268. [27] J. Dahl, P.C. Hansen, S.H. Jensen, T.L. Jensen, Algorithms and software for total variation image reconstruction via first-order methods, Numer. Algorithms 53 (1) (2010) 67–92. [28] J. Canny, A computational approach to edge detection, IEEE Trans. Pattern Anal. Mach. Intell. 8 (6) (1986) 679–698. [29] K. He, J. Sun, X. Tang, Guided image filtering, IEEE Trans. Pattern Anal. Mach. Intell. 35 (6) (2013) 1397–1409.

Q. Liu et al. / Signal Processing 137 (2017) 33–43

43

Qi Liu received the B.Sc. degree in Electronic Information Science and Technology from Huainan Normal University, China, in 2007 and the M.Sc. degrees in Bio-Medical Engineering from Northwestern Polytechnical University, Xi’an, China, in 2010. He is currently pursuing the Ph.D. degree in pattern recognition and intelligent system at Xidian University, Xi’an, China. His research interests include remote sensing image processing, image enhancement and computer vision.

Xinbo Gao (M’02-SM’07) received the B.Eng., M.Sc. and Ph.D. degrees in signal and information processing from Xidian University, China, in 1994, 1997 and 1999 respectively. From 1997 to 1998, he was a research fellow in the Department of Computer Science at Shizuoka University, Japan. From 20 0 0 to 20 01, he was a postdoctoral research fellow in the Department of Information Engineering at the Chinese University of Hong Kong. Since 2001, he joined the School of Electronic Engineering at Xidian University. Currently, he is a Professor of Pattern Recognition and Intelligent System, and Director of the State Key Laboratory of Integrated Service Networks, Xidian University. His research interests are computational intelligence, machine learning, computer vision, and pattern recognition. In these areas, he has published 6 books and around 150 technical articles in refereed journals and proceedings including IEEE TIP, TCSVT, TNN, TSMC, IJCV, Pattern Recognition etc. He is on the editorial boards of several journals including Signal Processing (Elsevier), and Neurocomputing (Elsevier). He served as general chair/co-chair or program committee chair/co-chair or PC member for around 30 major international conferences. Now, he is a Fellow of IET and Senior Member of IEEE.

Lihuo He is currently an Associate Professor at Xidian University. He received the B.Sc. degree in Electronic and Information Engineering and Ph.D. degree in Pattern Recognition and Intelligent Systems from Xidian University, Xi’an, China, in 2008 and 2013 respectively. His research interests focus on image/video quality assessment, cognitive computing, and computational vision.

Wen Lu received the B.Sc., M.Sc. and Ph.D. degrees in Signal and Information Processing from Xidian University, China, in 20 02, 20 06 and 2009 respectively. From 2010 to 2012, he was a postdoctoral fellow in the department of electrical engineering at Stanford University, USA. He is currently an Associate Professor at Xidian University. His research interests include image & video quality metric, human vision system, computational vision. He has published 2 books and around 30 technical articles in refereed journals and proceedings including IEEE TIP, TSMC, Neurocomputing, Signal processing, etc.