Total variation (TV) regularization is one of the well-known approach in solving the ill-posed inverse problem of image restoration. The main advantage of the TV is that it preserves edges without compromising the quality of denoising capabilities. The TV regularization is written as minimization of the following, (1) min u TV ( u ) = ∫ Ω | ∇ u ( x ) | d x where u : Ω ⊂ R 2 → R, and the minimization is taken over a suitable space such as the space of functions of bounded variation B V ( Ω ). Although the TV regularization removes noise while preserving salient edges, it can create artificial edges in homogeneous areas which is terms as ’staircasing’ artifacts. This is an undesirable quality while denoising medical imagery as creating new structures can cause confusions in the final medical diagnostic systems. To avoid the staircasing artifacts while maintaining the important edge preserving property of TV has motivated many studies to come up with adaptive solutions (Prasath et al., 2015).

To solve the energy minimization problem with weighted TV regularization (3) we utilize the alternating direction method of multipliers (ADMM) algorithm which is very efficient in terms of computational time. We provide a brief treatment as for the adaptive TV regularization the general ADMM formulation carries over with only little modifications, and we refer to Figueiredo and Biouscas-Dias (2010) for more details. The following iterative formulations constitute the ADMM for our problem, (10) u k + 1 = arg min u ∫ Ω ω ( x ) | ∇ u ( x ) | d x + α 2 ‖ u ( x ) − v k − b k ‖ 2 2 , (11) v k + 1 = arg min v ∫ Ω λ ( v − u 0 log ( v ) ) d v + α 2 ‖ u k + 1 − v − b k ‖ 2 2 , (12) b k + 1 = b k − ( u k + 1 − v k + 1 ) . Here λ is the Lagrangian multiplier, b auxiliary variable, and α > 00$]]> is a parameter (with a relation λ k = − α b k ). To solve the first minimization problem (weighted TV regularization with L 2 fidelity) for u various optimization methods can be utilized. In this work, we use Chambolle’s dual minimization method (Chambolle, 2004), (13) u k + 1 = v k + b k − 1 α ∇ · ( ω ( x ) p k ) , where p is the dual variable and can be obtained by the following formula, (14) p k + 1 = p k + Δ t ( ∇ div ( ω ( x ) p k ) − α ( v k + b k ) ) 1 + Δ t | ∇ div ( ω ( x ) p k ) − α ( v k + b k ) | . To solve for v, a simple solution exists, (15) v k + 1 = u k + 1 − λ α − b k 2 + ( u k + 1 − λ α − b k 2 ) + λ α u 0 .

We rescaled the images used here to the continuous domain [ 0 , 1 ]. The initial parameters of ADMM and dual minimization α = 1, b t = 0 = 0, p t = 0 = 0, v = u 0 , inner iteration of 5 (for dual variable p in Eq. (14)) were set in all our experimental results reported here. For the weights we utilize the pre-smoothing Gaussian kernel with σ = 1, and the neighbourhood size r = 5 (for the local histogram in the proposed generalized inverse gradient weight (9) in our comparison results. The Lagrangian multiplier λ and the contrast parameter k were modified according to optimal results presented below. Typically 20–30 iterations of the ADMM is sufficient to obtain good restoration results on a MATLAB implementation (on a Mac Pro Laptop, 2.3 GHz Intel Core i 7 processor, 8 GB 1600 MHz DDR3 memory), and we show the optimized results in the experiments reported below.

For comparing the quality of image restorations in this paper we utilize the root mean squared error (RMSE), relative error (RE), Pratt’s figure of merit (FOM), and structural similarity (SSIM) (Wang et al., 2004). These are given as, (16) RMSE = ( m n ) − 1 ∑ ∑ ( u − u o ) 2 , where u o is the original (noise-free) image, and the image is of size m × n. (17) RE = ‖ u − u o ‖ 2 ‖ u o ‖ 2 . Pratt’s FOM is calculated as, (18) FOM = 1 max ( E a , E d ) ∑ i = 1 E d 1 1 + γ d i 2 , where E a , E d are the number of actual and detected edge pixels (using the classical Sobel edge detector), γ = 0.1 scaling parameter fixed. The higher value of FOM indicate better quality edge map and value closer to 1 shows the denoising method has good edge preserving property. The SSIM is calculated between two windows ω 1 and ω 2 of common size N × N, and is given by, (19) SSIM ( ω 1 , ω 2 ) = ( 2 μ ω 1 μ ω 2 + c 1 ) ( 2 σ ω 1 ω 2 + c 2 ) ( μ ω 1 2 + μ ω 2 2 + c 1 ) ( σ ω 1 2 + σ ω 2 2 + c 2 ) where μ ω i the average of ω i , σ ω i 2 the variance of ω i , σ ω 1 ω 2 the covariance, and c 1 , c 2 stabilization parameters. The mean value of SSIM (MSSIM) is taken as the final error metric and if the value is closer to 1 indicating better structural similarity with u o the original image.

In our first experiment we consider a synthetic S h a p e s image of size 320 × 320 pixels which consists of various objects with different scales and intensity values. To test different TV regularization filters we add Poisson noise of strength 20 % to the original pristine image as shown in Fig. 2(a). Figure 2(b–d) show optimal restoration results (with respect to the highest MSSIM value) with TV (Le et al., 2007), adaptive TV with (5), and our proposed adaptive TV with (9). As can be seen in the bottom row of concentric circles (from the full S h a p e s image top right corner part), our proposed approach obtains better noise removal with edges preserved in all the objects. Moreover, compared to TV regularization we obtain no staircasing artifacts in flat regions. The slight smoothing of sharp corners is part of the TV regularization formulation and further corner adaptive weight in Eq. (4) is required to preserve them. To see the denoising effects clearly, in Fig. 3 we show a cross section taken in Fig. 2 and their corresponding denoising results with TV regularization filters. As can be seen, our regularization obtains better edge preservation without staircasing artifacts in flat regions.

Table 1 shows error metrics values for restoring corrupted Lena, Cameraman and Boat test images from the USC-SIPI Miscellaneous dataset. We fixed all the common parameters equal for the TV regularization variants and vary only the contrast parameter k, and Lagrangian multiplier λ, see Zhou and Lia (2012).

Figure 4 shows comparison restoration results with TV, ATV, and our proposed regularization filters on different X-ray imagery of hands. Our proposed generalized inverse gradient weight based TV regularization obtains better edge preservation with overall pixel intensity improvement indicating that the signal dependent Poisson noise is removed effectively. To show the denoising clearly in Fig. 5 we show a line taken across from the X-ray image in Fig. 4(a) top and the corresponding results from Fig. 4(b–d). We only show the first finger on the left (300 pixels long). It is clear that our proposed regularization obtains better restoration of improved intensity profile (contrast enhanced) and without staircasing artifacts. In summary, we see that our proposed ATV regularization works well in preserving important details while removing noise effectively when compared with related denoising methods from the literature.