Variational Anisotropic Gradient-Domain Image Processing

Gradient-domain image processing is a technique where, instead of operating directly on the image pixel values, the gradient of the image is computed and processed. The resulting image is obtained by reintegrating the processed gradient. This is normally done by solving the Poisson equation, most often by means of a finite difference implementation of the gradient descent method. However, this technique in some cases lead to severe haloing artefacts in the resulting image. To deal with this, local or anisotropic diffusion has been added as an ad hoc modification of the Poisson equation. In this paper, we show that a version of anisotropic gradient-domain image processing can result from a more general variational formulation through the minimisation of a functional formulated in terms of the eigenvalues of the structure tensor of the differences between the processed gradient and the gradient of the original image. Example applications of linear and nonlinear local contrast enhancement and colour image Daltonisation illustrate the behaviour of the method.


Introduction and Background
In 2002, Fattal et al. [1] introduced the method of gradient-domain high-dynamicrange compression. The technique consisted of first computing the gradient field ∇u 0 of a high-dynamic-range image u 0 : Ω → C, where Ω ⊂ R 2 is the image domain and C ⊂ R 3 is the colour space. The gradient was then rescaled nonlinearly as G = f (∇u 0 ). The resulting tensor field G, which is no longer necessarily a gradient field, was then reintegrated by solving the Poisson equation ∇ 2 u = ∇ · G (1) to obtain the compressed image. This is often solved by a gradient descent, The Poisson equations can be obtained through a variational approach by minimising the functional where Ω ⊂ R 2 is the image domain and the index F indicates the Frobenius norm taken over both image and colour coordinates. This is done by solving the corresponding Euler-Lagrange equations, leading to Equation (1). Perez et al. [2] soon generalised this into the technique called Poisson Image Editing, and showed that it allowed for a broad range of image processing applications such as image inpainting, seamless cloning, texture transfer, feature exchange, insertion of transparent objects through gradient mixing, texture flattening, local illumination correction, and seamless tiling. All this is obtained by various processing of the gradient or gradients of the original images to obtain the tensor field G used in Equation (1).
One significant challenge of Poisson Image Editing is that it tends to produce visual haloing or blurring artefacts in the reconstructed images. The case of G = 0 in Equation (1) is often used for image denoising. For greyscale images, in order to stop the diffusion at edges, and thus reduce the blurring, Perona and Malik [3] introduced an image-dependent, local, nonlinear diffusion method described by the equation where s = |∇u| 2 describing the image 'structure' has been introduced. (Despite the fact that the method is isotropic, but nonlinear and local, the authors termed it 'anisotropic diffusion' in the original publication [3]. This has caused considerable confusion in the terminology in the following literature.) They proposed two different diffusion coefficients with different properties, Instead of designing partial differential equations (PDEs) directly, Rudin et al. [4] used a variational approach like Equation (3) and introduced the concept of total variation. They showed that the minimisation of the functional E = Ω |∇u| dΩ (8) leads to the PDE ∂u ∂t = ∇ · ∇u |∇u| (9) Comparing with the Perona-Malik diffusion, Equation (5), we see that it can be written in the same form choosing For total variation and Perona-Malik diffusion, the extension to colour images is not that straight forward. Blomgren and Chan [5] introduced the concept of colour total variation by using a functional as an 2 norm of the functionals in Equation (8) for each of the colour channels resulting in where E i is the total variation for each colour channel according to Equation (8).
Later approaches have been based on the structure tensor by Di Zenzo [6] and Bigun and Granlund [7], with components S = ∇u · ∇u (12) where the dot product is taken over the colour coordinates. Sapiro and Ringach [8] proposed to use the eigenvalues (13) and the corresponding eigenvectors θ ± of the structure tensor as a basis for constructing the diffusion equations. In terms of these eigenvalues, an alternative to the colour total variation by Blomgren and Chan [5] can be obtained by the functional For greyscale images, where λ + = |∇u| 2 and λ − = 0, this reduces to total variation. For colour images, the corresponding Euler-Lagrange equations become ∂u ∂t which again is on the form of Equation (5) with Notice the coupling between the colour channels introduced by the sum in the denominator of Equation (15).
Tschumperlé and Deriche [9] extended this approach to a anisotropic diffusion by introducing the general Lagrangian density ψ(λ + , λ − ) in the functional The corresponding Euler-Lagrange equations are where D is the diffusion tensor where θ ± are the eigenvectors of the structure sensor S. It should be noted that this actually encompasses all previously presented diffusion methods as follows: gives the solution of Sapiro and Ringach [8] for colour images and total variation of Rudin et al. [4] for greyscale images, give the two equations of Perona and Malik [3], and gives the classical linear diffusion, Equation (2).
in general leads to isotropic equations, since then ∂ψ/∂λ + = ∂ψ/∂λ − = φ (s) and the diffusion tensor in Equation (18) reduces to the scalar diffusion coefficient D = 2φ (s). With the diffusion in the mutually orthogonal directions of maximal and minimal change can be controlled independently, like used by, e.g., Farup [10].
For other application of Poission image editing, i.e., applications of Equation (2) where G = 0, extensions to anisotropic and edge-preserving methods have been obtained by ad hoc modifications of Equation (2). This has been done for e.g., colour gamut mapping [11], colour image demosaicing [12], colour-to-greyscale conversion [13], and colour image Daltonisation [10]. Common to these is that they solve an equation on the form where D is a diffusion tensor constructed from the structure tensor like, e.g., Equation (18). However, a unifying variational formulation of anisotropic gradient-domain image processing has not yet been given. In this paper, we provide this by introducing the difference structure tensor based on the difference of the original image gradient ∇u and the tensor field G and follow the process of Tschumperlé and Deriche [9]. The resulting PDE is similar, but not identical, to the ones obtained by the ad hoc modification of Poisson Image editing, Equation (25). For the applications to linear and nonlinear local contrast enhancement and colour image Daltonisation, we show that the two approaches gives very similar results.

Variational Anisotropic Gradient-Domain Formulation
The variational formulation of anisotropic diffusion by Tschumperlé and Deriche [9] was based on the structure tensor S, Equation (18). In analogy, we introduce the difference structure tensor, with eigenvalues and corresponding eigenvectors θ ± . In analogy with Tschumperlé and Deriche [9], we define the functional in terms of the eigenvalues of this difference structure tensor, The corresponding Euler-Lagrange equations are It is clear that ∂ψ/∂u = 0. To derive the explicit form of the Euler-Lagrange equations, we will need to look at the individual components. Using Greek indices for the colour coordinates, Latin indices for the spatial coordinates, letting a comma denote partial differentiation with respect to the following coordinates, and using Einstein's convention of summing over repeated indices, the expression within the parenthesis in the first term can be written ∂ψ ∂u where the index p is used for summing over the two eigenvalues of the difference structure tensor. The first factor, ∂ψ/∂λ p , can be computed directly since the Lagrangian ψ is defined explicitly in terms of the eigenvalues of the difference structure tensor, and will thus depend on the design of the Lagrangian density.
The second factor can be computed implicitly following the method of Tschumperlé and Deriche [9] as follows: Multiplying with θ mk and θ ml , summing over the k and l indexes, using that θ mk θ pk = δ pm , and θ pk (∂θ pk /∂S ij ) = 0, the latter due to the orthonormality of the θ s gives (no sum). The last term, ∂S kl /∂u ρ ,i deviates from the derivation of Tschumperlé and Deriche due to the definition of the difference structure tensor, Equation (26), Inserted into Equation (30) and exploiting the symmetry of Equation (32) gives where the diffusion tensor has been defined. Inserting this into Equation (29) and solving by gradient descent gives the variational anisotropic gradient-domain image processing PDE (switching back to the regular vector notation): The form of this equation is similar to the one based on the ad hoc approach, Equation (25). The only difference is that, in this case, the diffusion tensor is computed from the difference structure tensor, and not the common structure tensor.

Results and Discussion
In order to demonstrate the usefulness of the proposed method, we apply it to three example problems: linear local contrast enhancement, nonlinear local contrast enhancement, and colour image Daltonisation.

Implementation
The algorithm represented by Equation (36) is implemented using the finite difference method (FDM). For the time derivative, the explicit (forward) Euler method is used. For the spatial derivatives, forward differences are used for the ∇u and G terms, and backward difference is used for the divergence in order to balance the overall resulting numerical scheme. The resulting code is available online. https://github.com/ifarup/variationalanisotropic-gradient-domain (accessed on 28 September 2021) A selection of colourful and detailed test images was downloaded from Pixnio. https://pixnio.com (accessed on 28 September 2021). The images are shown in Figure 1. All images are available under the CC0 licence, https://creativecommons.org/share-yourwork/public-domain/cc0/ (accessed on 28 September 2021).

The Diffusion Tensor
The purpose of the diffusion tensor is to direct the diffusion of the image information. The idea is that the diffusion should be avoided across the edges of the image, but allowed to some degree along the edges. The expression for the actual diffusion tensor is given in Equation (36). A choice has to be made for the function ψ. Here, we will use the Perona-Malik motivated functional The resulting diffusion coefficients along and across the edges, represented by 2∂ψ /∂λ + and 2∂ψ /∂λ − , respectively, are shown for one of the test images in Figure 2, where it is easily seen how the edges influence the directional diffusion.

Linear Local Contrast Enhancement
Local contrast enhancement in the gradient domain is a technique that is known to be particularly prone to haloing problems. First, we follow a simple approach and define the contrast enhancement as a scalar multiplication of the original image gradient where a > 1.
We compare this with the standard Poisson method, Equation (2), and the ad hoc addition with the diffusion tensor derived from the structure tensor instead of the difference structure tensor, Equation (18). We use a = 2 in order to obtain a quite extreme contrast enhancement, and use with K = 10 −3 for the variational approach based on the difference structure tensor, and K = 3 × 10 −4 for the ad hoc solution based on the standard structure tensor. Different constants are needed since the scaling of the structure tensor will be different in the two cases. Example results are shown in Figure 3. We can see that, with this choice of K, the results of the variational and the ad hoc approaches are indistinguishable, whereas the Poisson solutions exhibit severe haloing artefacts, as expected.

Nonlinear Local Contrast Enhancement
A somewhat more sophisticated contrast enhancement can be achieved by using the nonlinear gamma compression in the gradient domain, computed element-wise. Solving as for the linear case with the same choice of parameters as above gives the results shown in Figure 4. We can observe the indistinguishable behaviour between the variational approach and the ad hoc solution, and that the haloing problem of the Poisson solution is solved.

Colour Image Daltonisation
Colour image Daltonisation denotes the recolouring of colour images to increase detail visibility for colour-deficient observers. It is a problem that is known to be particularly prone to the loss of finer image details and textures. It has previously been shown that this process can be performed in the gradient-domain [10] by constructing a gradient field based on the original image as where u 0 is the original colour image, e d is the first principal component of the difference between the original image and a colour-vision-deficiency simulation the same image, and e c is a unit vector of maximum chromatic visibility for the colour-vision-deficient observer (in practice, a vector orthogonal to both the lightness axis and e d ). Figure 5 shows the resulting daltonised images (left) and the corresponding colour-vision-deficient simulations (right). We observe that the variational approach solves the same problem of lost detail visibility as the ad hoc anisotropic solution, and again gives results that are indistinguishable.

Conclusions
In this paper, a variational formulation of anisotropic gradient-domain image processing has been introduced. It generalises previous formulations of anisotropic processing and also introduces the difference structure tensor. The resulting PDE deviates somewhat from the one usually found when anisotropic diffusion is added in and ad hoc manner to do gradient-domain image processing. The difference is that the diffusion tensor is based on the difference structure tensor rather than the conventional structure tensor. The example applications of linear and nonlinear local contrast enhancement and colour image Daltonisation illustrate that the behaviour is very similar to what is obtained using the conventional approach. This indicates that the proposed variational formulation is well suited for rigorous derivations of anisotropic gradient-domain image processing for a broad range of applications.