1 Introduction

Image contrast enhancement is a fundamental pre-processing technique for many applications, such as surveillance system [1], medical image processing [2], analyzing images from satellites [3, 4]. Conventional histogram equalization (HE) can globally enhance the target image. However, HE trends to over-enhance those images with a large proportion of similar regions and results in intensity saturation or noise artifacts effect. Therefore, the approaches of partitioning the input histogram to several sub-histograms and enhancing them separately [5,6,7,8,9,10,11,12] have long been attempts to overcome aforementioned shortcomings and to enhance the input image locally and globally. Brightness-preserving bi-histogram equalization (BBHE) [5] is the earliest work to preserve mean brightness while improving the contrast. Based on the mean brightness, BBHE divides the input histogram into two sub-histograms and applies HE on each sub-histogram independently. Dualistic sub-image histogram equalization (DSIHE) [6] proposed by Wang et al. divides the input histogram into two sub-histograms based on the median value instead of the mean value. Other improvements of BBHE can refer to minimum mean brightness error bi-HE (MMBEBHE) [7], recursive mean-separate HE (RMSHE) [8], and brightness-preserving dynamic HE (BPDHE) [9], etc.

HE performs contrast enhancement based on the accumulative distribution function (CDF). Let \(I_1\) represents an input image of size \(M\times N\) in gray scale\([0,L-1]\). The CDF of the image is defined by

$$\begin{aligned} C(k)=\sum \limits _{k=0}^{L-1}\frac{h(k)}{M\times N},~~k=0,1,\ldots ,L-1, \end{aligned}$$
(1)

where h(k) is the input pixel intensity frequency with gray-level k. Based on the CDF defined by (1), HE maps input gray-level k into output gray-level T(k) by the following transform function

$$\begin{aligned} T(k)=\lfloor (L-1)\times C(k)+0.5 \rfloor . \end{aligned}$$
(2)

The operator \(\lfloor * \rfloor \) means rounding down. From (1) and (2), we can infer the increment in the output gray-level T(k) as

$$\begin{aligned} \bigtriangleup T(k)=(L-1)\times \frac{h(k)}{M\times N}. \end{aligned}$$
(3)

To overcome the shortcomings of HE method, techniques of splitting the input histogram into sub-histograms and applying HE on them become rational choices [5,6,7,8,9]. Some practices utilize nonlinear transform as gamma correction [13, 14], or clip histogram spikes as in exposure based sub-image histogram equalization (ESIHE) [10] and median-mean-based sub-image-clipped histogram equalization (MMSICHE) both by Singh and Kapoor [12], or use multiple sub-histograms equalization as adaptive image enhancement based on bi-histogram equalization (AIEBHE) by Tang and Isa [11]. Some other researchers fall back on using the 2D-histogram as two-dimensional histogram equalization (TDHE) [15], and 2D histogram equalization (2DHE) by Kim [16], etc. They all made their way in avoiding some of the shortcomings.

But above methods do not provide solutions to tackle with the problem causing by histogram spikes and pits, and cannot keep balance on enhancing the image globally and locally. By absorbing and integrating previous research results, we propose a new image contrast enhancement algorithm combing adaptive gamma transform, proportional histogram splitting, and standard deviation-based histogram addition. Gamma transform can adaptively restrain histogram spikes, histogram addition can fill histogram pits, and proportional histogram splitting can preserve mean brightness. The new algorithm can keep balance on contrast enhancement locally and globally, feature preserving, and overall quality. The rest of this paper is organized as follows: Sect. 2 presents the new algorithm, Sect. 3 provides the experimental results and discussions, and Sect. 4 concludes the paper.

Fig. 1
figure 1

From left to right, it is the original image and its histogram, results obtained by HE, adaptive gamma transform, and histogram addition, respectively

2 The proposed algorithm

Equation (3) implies that the increment of gray-level \(\bigtriangleup T(k)\) is proportional to its pixel intensity frequency h(k). Therefore, if there exist big histogram spikes in the input histogram, the corresponding output gray levels will occupy broad grayscale bands and squeeze other gray-bins with histogram pits, which causes the intrinsic shortcomings of over-enhancing some regions (regions with high-frequency bins) and contrast loss in the other regions. On the other hand, if the input histogram h is close to a uniformly distributed histogram [that is, h(k) is almost equal to each other for all k], \(\bigtriangleup T(k)\) will be almost equal to each other for all k too, that is, the output histogram will be almost uniform. Therefore, before applying HE, we can modify the input histogram as close to a uniformly distributed histogram as possible to fully exploit the dynamic range [14, 17].

Based on the above considerations, the proposed algorithm consists of three steps, that is, (i) adaptive gamma transform of the input histogram, (ii) histogram splitting and histogram bins redistribution between sub-histograms, and (iii) histogram addition and equalization. Adaptive gamma transform is employed to smooth histogram spikes, histogram addition is applied to fill the histogram pits, and histogram splitting is introduced to preserve the histogram mean. The gamma transform is defined as

$$\begin{aligned} \tilde{h}=ch^\gamma , \end{aligned}$$
(4)

where c and \(\gamma \) are positive constants, h represents the original input histogram, and \(\tilde{h}\) is the corresponding output. Applying gamma transform on the input histogram can smooth histogram spikes and restrain the noise artifacts effect. As the value of \(\gamma \) varying, different levels of smooth and restraining can be achieved. Since different images require different levels of smooth and restraining, in the proposed algorithm, the parameter \(\gamma \) is calculated adaptively according to the image intensity exposure. Thus, before gamma transform, the intensity exposure threshold is obtained by [12]

$$\begin{aligned} \varepsilon =\frac{1}{L}\frac{\sum \nolimits _{k=0}^{L-1}h(k)k}{\sum \nolimits _{k=0}^{L-1}h(k)}. \end{aligned}$$
(5)

For the proposed algorithm, gamma is defined by

$$\begin{aligned} \gamma =\left\{ \begin{array}{ll} 1-\varepsilon ,&{}\quad {0.5\le \varepsilon<1}\\ \varepsilon ,&{}\quad {\varepsilon <0.5} \end{array}\right. \end{aligned}$$
(6)

The first column of Fig. 1 is the original image and its histogram. The second column is the results of conventional HE. And the third column presents the results of applying adaptive gamma transform on HE. Comparing the original histogram with that of adaptive gamma transform-based HE and that of convention HE, we observe that the original histogram features are well preserved in the output histogram and the over-enhancement is nonexistence now. We can also observe that applying adaptive gamma transform on HE can lighten the noise artifacts effect as shown on the up right corner of the processed images.

Based on the definition of exposure threshold, the splitting threshold is defined as

$$\begin{aligned} T_s=L\times \varepsilon , \end{aligned}$$
(7)

where \(T_s\) is the threshold for splitting \(\tilde{h}\) into under exposed sub-histogram \(h_u\) and over exposed \(h_o\).

Unlike aforementioned contrast enhancement methods [5,6,7,8,9,10,11,12] that enhance each sub-histograms within the splitting thresholds, we redistribute the input under exposed sub-histogram from range \([0,T_s]\) to output range [0, U], and the over exposed sub-histogram from \([T_s+1,L-1]\) to output range \([U+1,L-1]\) to fully exploit the dynamic range, where U is the output threshold for under exposed and over exposed sub-histograms and is calculated as

$$\begin{aligned} U=L \times \frac{\left( \sum \nolimits _{k=0}^{T_s} h_u\right) ^\gamma }{\left( \sum \nolimits _{k=0}^{T_s} h_u\right) ^\gamma +\left( \sum \nolimits _{k=T_s+1}^{L-1} h_o\right) ^\gamma }. \end{aligned}$$
(8)

The fourth column of Fig. 1 presents the results of splitting the input histogram and enhancing them independently. We can observe that applying splitting histogram alone can preserve the histogram shape, but the noise artifacts effect is still obvious.

To deal with the problem of detail loss causing by histogram pits, we add two terms to \(h_u\) and \(h_o\), respectively, to further smooth the histogram \(\tilde{h}\). The two addition terms are their standard deviations of \(h_u\) and \(h_o\), respectively.

$$\begin{aligned} S_u= & {} \sqrt{\frac{1}{m}\sum \limits _{i}^m(h_{u}(i)-\psi )}, ~~i=0,1,\ldots ,U, \end{aligned}$$
(9)
$$\begin{aligned} S_o= & {} \sqrt{\frac{1}{n}\sum \limits _{j}^n(h_{o}(j)-\nu )} ~~j=U+1,U+2,\ldots ,L-1.\nonumber \\ \end{aligned}$$
(10)

where m and n are the number of gray levels in sub-histogram \(h_u\) and \(h_o\), respectively, and \(\psi \) and \(\nu \) are their corresponding sub-histogram mean. Finally, the sub-histograms for applying HE are defined by

$$\begin{aligned} \tilde{h_u}=S_u+h_u,~~~~ \tilde{h_o}=S_o+h_o. \end{aligned}$$
(11)

What calls for special attention is that the output histogram for \(\tilde{h_u}\) and \(\tilde{h_o}\) are in range \([0,T_s]\), and \([T_s+1,L-1]\), respectively.

The fifth column of Fig. 1 presents the results of adding a term (here is the standard deviation of the input histogram) to the original histogram and applying HE on the modified histogram. We can observe that histogram addition can preserve the histogram shape too. By incorporating adaptive gamma transform, histogram splitting, and histogram addition together, the output histogram will be close to a uniformly distributed histogram as well as close to the input histogram to the most extend. Thus, the output image will be a visually pleasing enhanced image.

3 Performance evaluation and discussion

To extensively evaluate the performance of the proposed algorithm, comparative experiments are conducted on 300 test images, which are (i) 25 reference images from the TID2013 [18]; (ii) 75 images from the miscellaneous, and aerials series of the USC-SIPI Image Database [19]; and (iii) 200 training images form the Berkeley Image Data Set [20]. The proposed algorithm is compared with HE, the weighted histogram approximation method (WAHE) [7], MMSICHE [12], AGCWD [13], 2DHE [16], and ESIHE [10] method. Authors of ESIHE compared their method against the BBHE, DSIHE, MMBEBHE, and RMSHE method [5,6,7,8], and showed their superiority.

We assess the performance of these six methods subjectively and objectively. Subjective assessment focuses on evaluating the visual quality of the enhanced images. Objective assessment involves image details, contrast, level of noise, natural appearance, etc., and is evaluated in terms of discrete entropy (DE) [1, 4, 10,11,12, 14, 15, 17], peak signal-to-noise-ratio (PSNR) [3, 11, 21], edge-based contrast measurement (EBCM) [4], gradient magnitude similarity deviation (GMSD) [16, 21, 22], and multiscale contrast similarity deviation (MCSD) [22].

  1. (i)

    Discrete entropy (DE) [1, 4, 10,11,12, 14, 15, 17]

DE characterizes the information contained in an image. Therefore, no results of any enhancement method can outperform the original image on DE, which is defined (in bits) by

$$\begin{aligned} \mathrm{De}=\sum \limits _{k=0}^{L-1}-p(k)\log _2 p(k). \end{aligned}$$
(12)
  1. (ii)

    Peak signal-to-noise-ratio (PSNR) [3, 11, 14, 21]

PSNR measures the noise level of the result, and a good enhancement method should not amplify the noise level of the origin. For an input image \(I_1\) and its enhanced image \(I_2\) with dimension of \(M\times N\), PSNR is defined by

$$\begin{aligned} \mathrm{PSNR} = 20 \times \log _{10}\left( \frac{\mathrm{MAX}}{\sqrt{\mathrm{MSE}}}\right) , \end{aligned}$$
(13)

where MAX is the maximum intensity value, e.g., 255 for 8-bit grayscale images, and MSE refers to the mean square error defined by

$$\begin{aligned} \mathrm{MSE}=\frac{1}{\mathrm{MN}}\sum \limits _{i=0}^{M-1}\sum \limits _{j=0}^{N-1}(I_1(i,j)-I_2(i,j))^2. \end{aligned}$$
(14)
  1. (iii)

    Gradient magnitude similarity deviation (GMSD) [16, 21, 22]

GMSD is a perceptual image quality index with high prediction accuracy [21]. It is an efficient distortion assessment metric which can measure the perceptual image quality of a distorted image against the reference [21, 22]. It predicts the quality of image by combining pixel-wise gradient magnitude similarity (GMS) with the standard deviation of the GMS map. The horizontal and vertical gradient magnitude images of \(I_{1}\) and \(I_{2}\) are defined as [21, 23]

$$\begin{aligned} \begin{aligned} h_{r}(i)&=\sqrt{(I_{1}\otimes h_{x})^2(i)+(I_{1}\otimes h_{y})^2(i)},\\ V_{d}(i)&=\sqrt{(I_{2}\otimes h_{x})^2(i)+(I_{2}\otimes h_{y})^2(i)}, \end{aligned} \end{aligned}$$
(15)

where symbol \(\otimes \) means convolving, \(h_{x}\) and \(h_{y}\) are the Prewitt filters in the direction of horizontal x and vertical y, respectively. And the gradient magnitude similarity (GMS) map is defined as

$$\begin{aligned} G(i)=\frac{2h_{r}(i)V_{d}(i)+\delta }{h^2_{r}(i)+V^2_{d}(i)+\delta }, \end{aligned}$$
(16)

where \(\delta \) is a positive constant for keeping numerical stability. Then, we compute the gradient magnitude similarity mean (GMSM) as

$$\begin{aligned} \mathrm{GM}=\frac{1}{\chi }\sum _{i=1}^{\chi }G(i). \end{aligned}$$
(17)

where \(\chi \) is the total number of pixels in image. And finally, the gradient magnitude similarity deviation of the GMS map is defined as

$$\begin{aligned} \mathrm{GMSD}=\sqrt{\frac{1}{\chi }\sum _{i=1}^{\chi }(G(i)-\mathrm{GM})^2}. \end{aligned}$$
(18)

Lower GMSD score denotes better image perceptual quality.

Fig. 2
figure 2

From left to right and top to bottom, the images are from the origin, HE, AGCWD, ESIHE, 2DHE, MMSICHE, WAHE, and the proposed, respectively. The corresponding contrast measurement values are listed in Table 1

  1. (iv)

    Edge-based contrast measurement (EBCM) [4]

The EBCM is based on the observation that human perception mechanisms are very sensitive to contours (or edges). An enhancement method should yield bigger EBCM result than the input image. The EBCM for image \(I_1\) of size \(M\times N\) is defined as

$$\begin{aligned} \mathrm{EBCM}=\frac{\sum \nolimits _{i=1}^{M}\sum \nolimits _{j=1}^{N}C(i,j)}{\mathrm{MN}}, \end{aligned}$$
(19)

where C(ij) represents contrast of pixel at location (ij) with intensity \(I_1(i,j)\) and is defined as follows

$$\begin{aligned} C(i,j)=\frac{|I_1(i,j)-f(i,j)|}{|I_1(i,j)+f(i,j)|}, \end{aligned}$$
(20)

where f(ij) is the mean edge gray level defined by

$$\begin{aligned} f(i,j)=\frac{\sum _{(y,z)_{\epsilon \alpha }(i,j)}g(y,z)I_1(y,z)}{\sum _{(y,z)_{\epsilon \alpha }(i,j)}g(y,z)}, \end{aligned}$$
(21)

where \(\alpha (i,j)\) is the neighboring pixel set of \(I_1(i,j)\), and g(yz) is the edge value at location (yz), which is the magnitude of the image gradient by Sobel operators [4]. Higher EBCM value means more edges information. Therefore, we consider the contrast of target image \(I_2\) is enhanced when \(\mathrm{EBCM}(I_2)>\mathrm{EBCM}(I_1)\). But it must be pointed out that higher EBCM value does not necessarily mean better visual quality.

  1. (v)

    Multiscale contrast similarity deviation (MCSD) [22]

MCSD is another perceptual image quality assessment with high correlation with human judgments [22]. The MCSD measures the contrast features resorting to images multiscale representation. First, contrast similarity deviations (CSD) for the reference image \(I_1\) and its distorted version \(I_2\) at three reduced resolutions are computed. Then, the final MCSD score is defined by

$$\begin{aligned} \mathrm{MCSD}=\prod \limits _{j=1}^{n\mathrm{Scales}}\mathrm{CSD}_j^{\beta _i}, \end{aligned}$$
(22)

where \(n\mathrm{Scales}\) represents the total scales, \(\beta _{i}\) is the weight of the \(j\mathrm{th}\) scale and \(\sum \nolimits _{j=1}^{n\mathrm{Scales}}\beta _j=1\). The contrast similarity deviation is defined by

$$\begin{aligned} \mathrm{CSD}=\sqrt{\frac{1}{\mathrm{MN}}\sum \limits _{i=1}^{M}\sum \limits _{j=1}^{N}(\mathrm{CSM}(i,j)-\mathrm{MCS})^2}, \end{aligned}$$
(23)

where \(\mathrm{MCS}\) (mean contrast similarity) is the mean pooling of contrast similarity map

$$\begin{aligned} \mathrm{MCS}=\frac{1}{\mathrm{MN}}\sum \limits _{i=1}^{M}\sum \limits _{j=1}^{N}S(i,j). \end{aligned}$$
(24)

And the contrast similarity map (CSM) between the original image and its enhanced one is

$$\begin{aligned} \mathrm{CS}=\frac{2\mathrm{CM}r.\times \mathrm{CM}d+2}{\mathrm{CM}r.^2+\mathrm{CM}d.^2+a}, \end{aligned}$$
(25)

where ‘\(.\times \)’ refers to element-wise multiplication of two matrices, ‘.2’ indicates element-wise square, and a is a positive constant avoiding divide by zero. CMr and CMd are the contrast maps for the original image and its enhanced one, respectively, and are defined as

$$\begin{aligned} \mathrm{CM} = \sqrt{\sum \limits _{i=1}^{M}\sum \limits _{j}^{N}w_{i,j}(I_{i,j}-\mu _I)^2},~~ \mu _I=\sum \limits _{i=1}^{M}\sum \limits _{j}^{N}w_{i,j}I_{i,j}, \end{aligned}$$
(26)

where \(w_{i,j}\) is the local window centered at (ij), and \(\mu _I\) is the local mean of I.

Fig. 3
figure 3

From left to right and top to bottom, the images are from the origin, HE, AGCWD, ESIHE, 2DHE, MMSICHE, WAHE, and the proposed, respectively. The corresponding contrast measurement values are listed in Table 2

Fig. 4
figure 4

From left to right and top to bottom, the images are from the origin, HE, AGCWD, ESIHE, 2DHE, MMSICHE, WAHE, and the proposed, respectively. The corresponding contrast measurement values are listed in Table 3

3.1 Subjective assessment

On most of the test images, MMSICHE, WAHE, 2DHE, ESIHE, and the proposed method present similar visual quality, the AGCWD method presents darker results, and HE provides over-enhanced bright images. On some of the test images, ESIHE and the proposed method obviously outperform the other four reference methods visually. On quite a few images, results of the proposed method are superior to all the other six reference ones. Due to limited space, only three images as Figs. 2, 3, and 4 are selected from those images that cause different visual quality.

Figure 2 presents a “TANK” on grassland. The processing results of HE, 2DHE, and WAHE are over-enhanced, making the grassland like snowfield. The result of AGCWD method, on the contrary, is too dark to be discernible. All methods except the proposed lose details on turret of the “TANK” as shown on the bottom right corner of the processed images.

Figure 3 shows that the HE, 2DHE, and MMSICHE methods over-enhance the moon, making the margin of the moon scattered on the results. The WAHE method enhances the sky too much making the background sky and the moon become an entire bright sky. The AGCWD method extends the gray levels to both ends, making the bright sky brighter and the dark background darker. The ESIHE method does enhance the background, but it introduces halo effect to the moon making the moon blurred. The HE, ESIHE, 2DHE, MMSICHE, and WAHE methods amplify the noise at different levels. Only the proposed method provides visually enhanced image. The output image of the proposed is more bright than the input, and the moon is not blurred.

In Fig. 4, the HE, ESIHE, 2DHE, MMSICHE, and WAHE methods over-enhance some of the bricks with different levels. The result of AGCWD is too dark and makes the building less clear. Only the proposed produces visually pleasing result.

Table 1 Enhancement comparison for applying the reference methods on Fig. 2
Table 2 Enhancement comparison for applying the reference methods on Fig. 3
Table 3 Enhancement comparison for applying the reference methods on Fig. 4
Table 4 Enhancement comparison on 300 test images by average

3.2 Objective assessment

Tables 1, 2, 3 list the performance values of DE, PSNR, EBCM, GMSD, and MCSD as shown in Figs. 2, 3, 4. Table 4 presents the average values of aforementioned five performance metrics on 300 test images.

Table 1 shows that, of all five measurement metrics, the proposed method outperforms the other six reference methods. On PSNR, GMSD, and MCSD metrics, the proposed method is much better than all the other reference methods. Especially for the GMSD and MCSD measurements, the proposed method is better almost by one order of magnitude even compared with the second best ones.

Table 2 demonstrates that the proposed method outperforms the other six reference methods measured by DE, GMSD, and MCSD. The proposed method ranks the third on PSNR, and the first two are the MMSICHE and AGCWD methods. However, they perform bad on EBCM as inferior to all other methods by one order of magnitude. Since EBCM measures the enhancement level, both the MMSICHE and AGCWD methods do not enhance the target image. The proposed method performs the third on EBCM, and the WAHE and 2DHE are the best two. But they attain high EBCM by augmenting the noise level, which is indicated by their worst PSNR scores and visual quality shown in Fig. 3.

As for Table 3, the proposed method produces the best image quality measured by DE, PSNR, GMSD, and MCSD. One thing that needs to be noted is the proposed method attains the same second best EBCM score as the 2DHE method, but provides much better image quality measured by the other five metrics.

Table 4 presents the average experimental results on 300 test images. The ESIHE and MMSICHE attain better PSNR values than the proposed method. But on average, both method barely enhance the edge information as indicated by their EBCM value of 101.850 and 102.004, respectively, which are almost the same as the original average EBCM value 101.386. The best two average values on EBCM are 117.041 and 116.094 by the WAHE, and HE method, respectively. However, both methods suffer from noise artifacts effect as indicated by their low PSNR values. The proposed method performs the best on DE, GMSD, and MCSD, and the second best on PSNR.

4 Conclusion

we proposed a contrast enhancement method that performs well on many image quality metrics. On average, the proposed algorithm have some improvements of 2.89, 9.83, 28.32, and 26.38% over the second best ESIHE algorithm on DE, EBCM, GMSD, and MCSD. Though the proposed method ranks the third on PSNR, the two methods with the highest PSNR values (ESIHE and MMSICHE) barely improve the contrast as denoted by their EBCM values which are almost the same as the input. The proposed method ranks the third on EBCM too. The best two are WAHE and HE, but they performs much worse on the other four metrics. In a word, the overall performance of the proposed method is superior to state-of-the-art reference methods.