FAST DUAL MINIMIZATION OF THE VECTORIAL TOTAL VARIATION NORM AND APPLICATIONS TO COLOR IMAGE PROCESSING

. We propose a regularization algorithm for color/vectorial images which is fast, easy to code and mathematically well-posed. More precisely, the regularization model is based on the dual formulation of the vectorial Total Variation (VTV) norm and it may be regarded as the vectorial extension of the dual approach deﬁned by Chambolle in [13] for gray-scale/scalar images. The proposed model oﬀers several advantages. First, it minimizes the exact VTV norm whereas standard approaches use a regularized norm. Then, the numerical scheme of minimization is straightforward to implement and ﬁnally, the number of iterations to reach the solution is low, which gives a fast regularization algorithm. Finally, and maybe more importantly, the proposed VTV minimization scheme can be easily extended to many standard applications. We apply this L 1 vectorial regularization algorithm to the following problems: color inverse scale space, color denoising with the chromaticity-brightness color representation, color image inpainting, color wavelet shrinkage, color image decomposition, color image deblurring, and color denoising on manifolds. Generally speaking, this VTV minimization scheme can be used in problems that required vector ﬁeld (color, other feature vector) regularization while preserving discontinuities.


1.
Introduction. This paper is devoted to the regularization of multidimensional/ vectorial signals, such as color images or other vector fields (vector of features, normal of level sets, etc.), based on the vectorial TV norm. The regularization of vectorial images in the context of continuous models have already been proposed in the literature. In [33], Sapiro and Ringach proposed an anisotropic diffusion model for vectorial images based on the computation of the eigenvalues of the M -D image defined as a 2-D manifold. The eigenvectors determine the directions of maximal and minimal change of the vectorial image, and their associated eigenvalues give their rate of change. In [8], Blomgren and Chan defined a definition of the VTV norm which satisfies several good properties. Sochen, Kimmel and Malladi in [35] defined a vectorial image denoising models based on a 2-D manifold embedded in a (M + 2)-D space and the Polyakov energy. We propose in this paper a vectorial regularization algorithm that is as efficient as [35] and improve the regularization algorithms [33,8].
The proposed VTV minimization model is based on the dual formulation of the vectorial TV norm. The idea of dualizing the (scalar) TV norm was introduced by Chan, Golub and Mulet (CGM) in [15] to solve the denoising Rudin-Osher-Fatemi (ROF) model [32] for gray-scale images in a faster way than the previous models. CGM introduced in the TV minimization model an extra maximization problem w.r.t. a new (dual) variable. The new minimization-maximization problem can be solved using the Sion's theorem [34]. Sion's theorem is a generalization of John Von Neumann's minimax theorem [40] (see also [2,22]) introduced in game theory for minimizing the maximum possible loss.
Based on the CGM model, Chambolle in [13] developed an efficient dual approach to minimize the scalar ROF model. Chambolle's algorithm is faster than CGM even if the convergence of Chambolle's scheme is linear and the CGM's scheme is quadratic. Chambolle's algorithm is faster because the cost per iteration to use CGM is higher (CGM needs to solve a linear system at each iteration). A recent fast minimization algorithm for the scalar ROF model was proposed by Darbon and Sigelle (DS) in [20] based on graph cuts. Although Chambolle's algorithm is not as fast as the model of DS to solve the variational scalar ROF model, it is still fast and presents some advantages compared with CGM and DS. First, Chambolle's model use the exact scalar TV norm whereas CGM model regularizes it to minimize it. Then, the numerical scheme of [13] is straightforward to implement unlike the CGM and DS algorithms. Besides, the TV norm of DS is anisotropic whereas the TV norm of Chambolle is isotropic. Finally, we will see that the Chambolle's model extends nicely to color/vector images whereas the question of extension is open for the CGM model and the generalization of DS model to color images is not as efficient as in the scalar case (see [19]). This paper proposes these two following contributions: • Extension of Chambolle's model [13] to multidimensional/vectorial images. We will see that the extended algorithm presents the same nice properties than the scalar algorithm since the VTV minimization algorithm is fast, easy to code and well-posed. Unlike [8,4,16], the proposed vectorial scheme does not regularize the VTV to minimize it. Finally, the numerical solution converges to the continuous minimizing solution in the vectorial BV space.
• Extension of the proposed VTV minimization scheme to several standard applications such as deblurring, inpainting, decomposition, denoising on manifolds. In fact, this vectorial regularization scheme can be applied to any problems that require a L 1 regularization process for vectorial components.
The paper is organized in two parts. The first part is the introduction of the VTV minimization algorithm. We introduce some notations and the dual formulation of the VTV. We then introduce the VTV minimization algorithm and prove its convergence to the continuous minimizing solution. Then, we compare our model to related denoising models. The second part of this paper aims at showing that this vectorial regularization process can be easily extended to any problems that require L 1 vector field regularization. We will focus on color image processing problems such as color inverse scale space, color denoising with the chromaticity-brightness color representation, color inpainting, color wavelet shrinkage, color decomposition, color deblurring or color denoising on manifold. But we can also extend this algorithm to other problems such as the vector field regularization of diffusion tensors in medical imaging.
2.1. Notation and definition. In this paper, we are interested in the standard dual definition of the VTV in order to extend the Chambolle's scalar algorithm [13] to vectorial images such as color images or any vector fields. This section introduces some notations and definitions regarding the vector-valued TV norm used throughout the paper. Let us consider a vectorial (or M -dimensional or multichannel) function u, such as a color image or a vector field, defined on a bounded open domain Ω ⊂ R N , as follows: where each scalar function u i : Ω → R, 1 ≤ i ≤ M is one of the M components (or channels) of the vector valued function u.
Definition 2.1. For a given vector valued function u : Ω → R M , the vectorial TV norm is denoted the finite positive measure: , ∇· is the divergence operator such that ∇ · q := (∇ · q 1 , ..., ∇ · q N ) : Ω → R M , ∀q : Ω → R M×N , ∇ · q i := N j=1 ∂ xj q xj i : Ω → R, ∀i ∈ [1, M ], the product < ., . > is the Euclidean scalar product defined as < v, w >: Thus, the VTV norm (1) can be defined of different ways, depending on the set P of functions of the dual variable p. Let us consider two sets: If we assume that u ∈ C 1 (Ω; R M ), i.e. u are smooth functions, and p ∈ C 1 c (Ω; R M×N ) are with compact support, then integration by parts gives: Then if p ∈ P 1 then the supremum of (1) is obtained for: which defines the vectorial TV norm for smooth functions u as: i.e. the sum of the TV of each channel. Then if p ∈ P 2 then the supremum is obtained for: where the vectorial gradient of u is defined as ∇u : Thus, the vectorial TV norm given by P 2 for smooth functions u is: Hence, the VTV defined in (1) have different expressions for smooth functions u for different sets P 1 or P 2 . The difference between sets P 1 or P 2 is important in the denoising context. Indeed, the vectorial TV based on P 1 is defined as the sum of the TV of each channel. This means that channels are considered as independent in the denoising process, which is not true on real-world images. The vectorial TV based on P 2 introduces a coupling between channels. Each channel use information coming from other channels to improve the denoising model as we will see in Section 2.5.
Finally, we want to emphasize that the vector valued TV norm defined by the set P 2 , i.e. the L 2 -norm, is the most standard definition of the VTV norm as introduced in the book of Ambrosio, Fusco and Pallara [1]. From this point, we will use the L 2 -norm for the VTV norm throughout this paper.
We now present standard results regarding the VTV. These results can be found e.g. in [1,27,24].
We define the space BV (Ω; R M ) of vector valued functions as the set of functions u ∈ L 1 (Ω; R M ) such that Ω |Du| < ∞, where the vectorial TV norm Ω |Du| is defined in Definition 2.1 with the set of functions p ∈ P 2 . The space BV (Ω; R M ) endowed with the following norm: where u BV (Ω;R M ) := Ω |Du|, is a Banach space.  Functional Ω |Du| can be decomposed as follows:   [32] is one of the most influential variational and PDE-based image denoising models in image processing. This denoising model removes noise in gray-scale images while preserving main features such as edges. In the case of vector-valued images, the ROF model is naturally extended as follows: where f ∈ L ∞ (Ω; R M ), f := (f 1 , ..., f M ) is the given (noisy) vector valued image and f − u 2 L 2 (Ω;R M ) = Ω |f − u| 2 dx is the L 2 fidelity norm.
Existence of a minimizer for the variational vectorial ROF model (4) can be proved by standard arguments of the vectorial BV space. We include the existence proof in this section for the sake of completeness. First, let us introduce truncated vectorial functions: Let f ∈ L ∞ (Ω; R M ) and u : Ω → R M . We define the truncated functionû : Proof. The proof is, for example, based on [23]. The definition of the truncated function implies that |f −û| ≤ |f − u| from which follows (j).

The truncated function is given byû
Using the chain rule, we haveû ∈ BV (Ω; R M ).
Then, it is easy to prove the existence of a minimizer for the VROF. The vectorial ROF model  Therefore, according to Theorem 2.4, there exists a subsequence {û n k } converging to a function u ∈ BV (Ω; R M ). Since T V is lower semicontinuous in BV (Ω; R M ) by definition, we have: and from Fatou's Lemma: which implies that It follows that u ∈ BV (Ω; R M ) is a minimizer of F and the vectorial ROF model (4).
Finally, the VROF model (4) is convex, which guaranties the uniqueness of the minimizer.

2.3.
Proposed algorithm for VROF minimization. In this section, we present the solution minimizing the vectorial ROF model (4). Using the dual definition of the vectorial TV norm given in (1), the variational model (4) can be written as Using the minimax theorem [2,22], inf and sup can be swapped since (5) is convex in u and concave in p and the set {|p| ≤ 1} is bounded and convex. Thus (5) is first minimized w.r.t. u using the Euler-Lagrange's technique, which gives the minimizing solution of the VROF model: and for each channel: Substituting solution (6) in (5), we get the constrained maximization problem: which is equivalent to this minimization problem: The Euler-Lagrange's technique gives the necessary condition for optimality at each x: where α(x) is the Lagrange multiplier associates with the constraint |p| ≤ 1. As noticed by Chambolle in [13] for the scalar case, it can also be eliminated in the vectorial case as follows. If |p| < 1, then the Lagrange multiplier is not active, i.e. α = 0 and ∇ λ∇ · p − f = 0 and if |p| = 1, then the Lagrange multiplier becomes active and α = |∇ λ∇ · p − f | > 0. In any case, the Lagrange multiplier is equal to: Introducing (9) in (8), we get: which can be solved as in [13] with a semi-implicit gradient descent scheme as follows: p n=0 = 0, such that: which means for each channel: Inverse Problems and Imaging Volume 2, No. 4 (2008),  As in [13], the convergence of the iterative scheme (11) also holds since we have the following theorem.

Proof.
We follow the proof given by Chambolle in [13] in the case of gray-scale images. Let us denote |.| X,Y , < ., . > X,Y the Euclidean norm and inner product We call η := (p n+1 − p n )/δt and we have: which means that p n+1 = p n in both cases.
It remains to show that δt ≤ 1/8. Let us define the discrete divergence operator applied to a vector p i = (p x i , p y i ), 1 ≤ i ≤ M , M being the number of channels, at the grid point (i x , i y ) (following [13]): Then, the norm of the divergence operator is defined by τ 2 := sup |p|≤1 |∇ · p| 2 . The discrete square norm of p is We have: using the zero boundary condition p 0,iy = p Nx,iy = p ix,0 = p ix,Ny = 0, ∀(i x , i y ), which implies that τ 2 = 8 and so δt = 1/8. We notice that the temporal bound does not depend on the number of channels.

Convex analysis and projection point of view.
The VROF minimization algorithm introduced in the previous section can also be deduced from convex analysis. Indeed, the projection algorithm of Chambolle [13] for gray-scale images can be easily extended to the vectorial case.
is given by where Π λK BV (Ω;R M ) is the orthogonal projection operator onto the closed convex set associated to ||.|| BV (Ω;R M ) and defined as follows: Proof. Let us call J BV (u) = ||u|| BV (Ω;R M ) . The Euler-Lagrange of (13) is equal to: where ∂J BV is called the sub-differential of J BV . The sub-differential can be re-written as f −u λ ∈ ∂J BV (u), which, according to convex analysis [22], is equivalent to where J ⋆ BV is the Legendre-Fenchel (LF) transform (or conjugate function) of J BV . In convex analysis [22], the LF of a function J = ∞ convex, one-homogeneous (J(λu) = λJ(u), λ > 0, ∀u ∈ X) and lower semi-continuous, defined as , is the characteristic function of a closed convex set K J : where K J is defined as: The closed convex set K BV can be defined according to (18): which naturally gives: Finally, since (6) gives u = f − λ∇ · p, then the projection operator in (14) is determined by the iterative scheme λ∇ · p n → λ∇ ·p = Π λKBV (f ), introduced in the previous section.

2.5.
Code and results. The iterative scheme (12) is straightforward to implement 1 . We already defined in the previous section the discrete divergence operator. It remains to define the discrete gradient operator applied to a function v at the point (i x , i y ) (following [13]): The iterative process (12) is stopped if max(|u n+1 − u n |) ≤ r, where r is a given residue and we recall that u n+1 = f − λ∇ · p n+1 . In all experiments, initial values are chosen to be p = 0, the residual r = 10 −4 .
Before testing the vectorial denoising ROF model, we come back to (12) where it is interesting to notice that a coupling term between the M channels appears, precisely N j=1 |∇ ∇ · p n j − f j /λ | 2 . Blomgren and Chan have already noticed in [8] the existence of a coupling term in the minimization equation of the vectorial ROF model. They observed that the coupling improves the denoising process compared with a direct application of the scalar ROF model channel-by-channel. Figure 1 gives an illustrative example of comparison between the vectorial ROF model and the scalar ROF model channel-by-channel. Globally speaking the coupling term helps to better restore parts in images where intensities are weak.
As a summary of the first part of this paper, the standard dual definition of the vectorial TV norm, introduced in (1), will be our basic tool to develop different color image processing models such as image denoising, vector field denoising, image decomposition, image deblurring and image inpainting. This definition of the vectorial TV norm uses many advantages. First of all, there is no need to introduce a regularization parameter of the TV norm as it was done in [8,4,16] and many other papers to carry out the minimization process. Indeed, the vectorial TV norm is usually not directly minimized but the regularized version Ω |∇u| 2 + ǫ dx, ǫ being a very small parameter to be faithful to the original norm and useful to avoid numerical instabilities in the minimization flow. Secondly, the minimization process based on the dual approach is much faster than in a standard regularized approach [8] (explicit gradient descent flow). Indeed, the direct consequence of the regularization parameter ǫ is the obligation to use a very small temporal step to ensure a correct and stable minimization process. Thus the minimization process is slow because a large number of iterations is necessary to reach the steady state minimizing solution. Figure 3 presents a comparison between our dual-based ROF model and the standard ǫ-based ROF model used in [8]. Our model converges fast to the denoising result whereas the standard model converges slowly because of the regularization parameter. Besides, the quality of the denoising process is slightly better with our model, still because of the ǫ parameter. Finally, another advantage of our approach is to provide a numerical minimization scheme easy to code.
2.6.1. Sapiro-Ringach and Blomgren-Chan's Models. In [33], Sapiro and Ringach propose an anisotropic diffusion model for vectorial images based on Riemannian geometry. From a Riemannian point of view, the vectorial image u : (x, y) → (u 1 (x, y), ..., u M (x, y)) can be seen as a parametric 2-D manifold embedded in a M -D Euclidean space. The (square) distance between two close points on the 2-D Riemannian manifold is given by ds 2 := µ,ν g µν dµdν, µ, ν = {x, y}, and g µν :=< ∂ µ u, ∂ ν u >, ∂ µ u := ∂u/∂µ, is the first fundamental form or metric tensor [28]. The metric tensor g µν is considered as a "vectorial edge" indicator function, i.e. a detector of discontinuities for vectorial images. More precisely, the eigenvectors of the image metric g µν determine the directions of maximal and minimal change at a give point on the manifold-image, and their associated eigenvalues give their rate of change. Eigenvalues of g µν are provided by: λ ± = 1 2 g xx + g yy ± (g xx − g yy ) 2 + 4g 2 xy and the eigenvectors are (cos θ ± , sin θ ± ), where θ + = 1 2 arctan 2gxy gxx−gyy and θ − = θ + + π 2 . Based on the previous analysis of the image metric g µν , Sapiro and Ringach propose to define the anisotropic diffusion for vectorial images by smoothing in the direction of minimal change, θ − , i.e. parallel to the "vectorial" edge, using the flow: . This approach is consistent with the scalar case presented in [12]. Indeed, in the case of gray-scale images, i.e. M = 1, it is easy to show that (cos θ − , sin θ − ) = ∇u |∇u| which implies that ∂u ∂t = ∂ 2 u ∂θ 2 − = ∇·( ∇u |∇u| ) and the TV norm is equal to T V (u) = Ω λ + dx since λ + = |∇u| 2 and λ − = 0. From the latter scalar approach, Sapiro and Ringach propose in [33] the following vectorial TV norm: However, they do not recommend a specific f . One natural choice, as proposed by Blomgren and Chan in [8], is f := λ + + λ − because it is equivalent to the TV norm for the scalar case and it is a non-decreasing function of λ ± .
We notice that (19) while considering f = λ + + λ − is strictly equal to the definition of the vectorial TV norm used in this paper and presented in Equation Finally, Blomgren and Chan noticed that the norm (20) has the tendency to produce a color smearing contrary to their definition of vectorial TV given by Since our model corresponds to the norm (20) for smooth functions u, our algorithm still gets the color smearing. However, despite of this drawback, as we have already said, the norm (20) offers many important advantages. We remind that no regularization of the TV norm is done. Besides, the number of iterations to converge to the denoised vectorial image is much less important in our approach than in [33,8] as shown on Figure 3 and our numerical scheme is straightforward to implement. Finally, a dual formulation of the vectorial TV norm provides the basic element to prove the existence of a denoising solution in Section 2.2, which is difficult to prove with T V BC in [8].
2.6.2. Sochen-Kimmel-Malladi's Model. In [35], Sochen, Kimmel and Malladi introduce another Riemannian framework to smooth vectorial images such as color or texture images. In their approach, they represent the given vectorial image as a parametric 2-D manifold embedded in a (M + 2)-D space, M being the dimension of the vectorial image. The core of their framework is the Polyakov functional, which measures the weight of a mapping X : Σ → M between an embedded manifold (the image manifold) Σ and an embedding manifold M as follows: where g µν is the first fundamental form of the manifold Σ, dΣ is the integration element w.r.t. the local coordinates on Σ, h ij is the metric tensor of the embedding space M, g µν is the inverse metric of g µν , g is the determinant of g µν , µ, ν = x, y and i, j = 1, ..., M in our paper and ∂ µ X i := ∂X i /∂µ. Finally, the Einstein summation convention is used in (21), which means that when identical indices appear one up and one down, they are summed over. The Polyakov functional is related with harmonic maps which are geodesics or minimal surfaces for curves and surfaces. Indeed, if the metric tensor of the embedded manifold Σ is chosen to be the induced metric tensor, g µν := ∂ µ X i ∂ ν X j h ij , by pullback procedure, then X which minimizes the Polyakov action are called harmonic maps and the Polyakov functional is reduced to the Euler functional: which describes the length/(hyper-)area of a curve/(hyper-)surface Σ and √ g being the invariant-area element on Σ.
In the case of vectorial images, we have X : (x, y) → (x, y, u 1 (x, y), ..., u M (x, y)) and the components of the metric tensor are thus: This implies the denoising functional for vectorial images [35]: where the term M i,j=1 < ∇u i , ∇u j > 2 is a cross-correlation term of orientation between different channels that directs different channels to align together. Minimizing (22) w.r.t. u i gives the Beltrami flow: ∂ui where ∆ g is the Laplace-Beltrami operator.
There exist some connections between the previous model and the dual TV norm. If we define a new vectorial image as u : (x, y) → ( 1 √ 2 x, 1 √ 2 y, u 1 (x, y), ..., u M (x, y)), then the dual vector TV norm is equal to: Functional (23) is equal to (22) when the cross-correlation term is null. The solution of the vectorial ROF model inf .., f M ) as the given image, is given by: and p i is the steady state solution of the following flow: An interesting application of (24) and (25) is for gray-scale images. Indeed, Functionals (22) and (23) are equal for gray-scale images, S(u) = u BV (Ω) = Ω 1 + |∇u| 2 dx, since the cross-correlation term is null. This means that the Beltrami flow can be "dualized" in the case of gray-scale images as proposed by Chambolle in [13]. In other words, the minimizing solution of where u : Ω → R, is given either in the standard Euler-Lagrange approach by: or in the dual approach by: where p is the steady state of: We have noticed that the dual approach (28) is faster than the Euler-Lagrange (EL) approach (27). We mean that the computation of the minimum of (26) is obtained faster with solution (28) than solution (27), even if we consider an explicit steepest gradient descent scheme for both approaches. We think that the explicit computation of the Lagrange multiplier in the dual approach helps to speed up the minimization speed. We have already observed that the dual approach is faster than the EL approach in the ROF denoising model. Indeed, Figure 3 presents a comparison between the dual and the EL approaches, where the dual is clearly faster. Finally, we have also noticed that the semi-implicit scheme does not really speed up the minimization in the dual approach. Similar minimization speeds are obtained with the implementation of the explicit gradient descent scheme.
3. Part II: Extension to standard color image processing. In the first part of this paper, we have introduced a regularization algorithm for color/vectorial images based on the VTV norm. In the second part of this paper, we extend this L 1 vectorial regularization algorithm to standard color image processing problems. We will see that the extension works well for many color problems such as the color inverse scale space, the color denoising with the chromaticity-brightness color representation, the color wavelet shrinkage, the color inpainting, the color decomposition, the color deblurring and finally the color denoising on manifolds. Generally speaking, the proposed VTV minimization algorithm can be used in problems that required vector field (color, other vector features) L 1 regularization while preserving discontinuities.

3.1.
Color inverse scale space. In this section, we improve the VROF model, introduced in Section 2.2, with the iterative regularization process proposed by Osher, Burger, Goldfarb, Xu and Yin in [31]. This process, also called inverse scale space, enhances the denoising task defined by the ROF model by better preserving the contrast. The proposed iterative regularization process is based on geometry. The basic idea is to match not only the intensity values between the denoising image and the noisy one but also the normals of level sets in each image while applying the TV norm which preserves the location of edges while smoothing out the noise. In this section, we improve the VROF model by extending the iterative regularization scheme [31] to color images.
We consider u ∈ C 1 (Ω; R M ) to explain the model. Given the noisy vectorial image f ∈ L ∞ (Ω, R M ), an "approximation of normals of level sets" can be computed by carrying out the minimization model: where n 1 := ∇u1 |∇u1| is called the approximation of "normals of level sets". The Euler-Lagrange equation of the previous equation is equal to: where v i := f i − u i , ∀i. We now introduce n 1 in the following minimization model: In the previous equation, we have: Using (29), we obtain: which gives: where η is independent of u and therefore useless in the minimization process. Finally, we end up with the modified vectorial ROF model: where v 1 := f − u 1 . Hence, the implementing iterative procedure for the vectorial image case is as follows: 1. Initialize: u 0 = 0 and v 0 = 0. 2. For k = 0, 1, 2, ...: compute u k+1 as a minimizer of the modified vectorial ROF model: and update v k+1 := v k + f − u k+1 . Figure 4 presents an application of the proposed iterative regularization process, which better preserve the color contrast.

Color Denoising with the Chromaticity-Brightness representation.
In Section 2.2, the color denoising model has been formulated as a vectorial model in the RGB representation of color images. In [16], authors showed that the Chromaticity-Brightness model or the Hue-Saturation-Value (HSV) enable to better restore color images because these models are closer to human perception. In this section, we propose to perform the color image denoising in the CB representation of color images. In the RGB color model, the given image is given by In the CB color model, the image f is separated into the brigthness component B f := |f | : Ω → R and the chromaticity component Thus, the brightness component is a gray-scale image and the chromaticity component stores the color information which takes values on the unit sphere. Following [27], the previous constraint on the sphere is introduced as a penalization term in the variational model: which denoises the chromaticity component. (30) is regularized using the Ginzburg-Landau functional [11] as follows: The Ginzburg-Landau functional ( |D| 2 − 1 2 2 ) is introduced in order to use the fast minimization algorithm proposed in Section 2.3. Indeed, we have two minimization problems. The first one is the VROF: and the second one is defined by: whose the Euler-Lagrange equation is given by D( 2β α (|D| 2 − 1) + 1) + C = 0. We use the following semi-implicit iterative scheme to minimize (31): .
The brightness component is a gray-scale image that can be denoised with the model [13]. Finally, the denoised image is given by u = C/B, where C and B are respectively the denoised chromaticity component and the denoised brightness component. Figure 5 presents the denoised color image with the CB color representation. The SNR is 14.57 with the RGB color representation and 14.72 with the CB representation.

3.3.
Color image inpainting. Image inpainting consists in recovering/interpolating an image in regions where the original information is missing, see [6,17]. Given the image domain Ω, a color image f ∈ L ∞ (Ω, R M ), a region Ω I ⊂ Ω where the image will be inpainted, the variational color image inpainting model is as follows: We propose to solve the previous minimization problem with boundary condition with the following convex variational model: The minimization problem w.r.t. u is the vectorial ROF model and the minimization w.r.t. v gives v |Ω I = u and v |Ω\ΩI = f . 3.4. Color wavelet shrinkage. In this paper, we have considered for vectorial images the space of regularization as the BV space. However, other spaces of smoothness can obviously be used. A popular regularization space is L 2 (Ω; R M ). In this case, the minimizer (not really useful for denoising) of ||u|| 2 L 2 (Ω;R M ) + 1 2λ ||f − u|| 2 L 2 (Ω;R M ) is given by the Tikhonov method [38]. In this section, we consider the Besov space for color images, which will relate variational model to wavelet shrinkage through a denoising process as shown in [14].
is given by where B α p (L p (Ω; R M )) is a Besov space. In the case of p = 1, α = 1, the closed convex set K B α p is given by: Proof. The proof is the straightforward extension of the paper [14]. In this paper, Chambolle, DeVore, Lee and Lucier showed relations between wavelet-based algorithms and variational models in the case of the Besov space and gray-scale images.
Let us express the vectorial image u : Ω → R M in an orthogonal wavelet basis for L 2 (Ω) that we denote ψ j,k , then for f ∈ L 2 (Ω; R M ), we have: where c f j,k,ψ : Ω → R M . The norm of f in L 2 (Ω; R M ) is: where i ∈ [1, M ] and the norm of f in B α p (L p (Ω; R M )) is ( [14]): It can be shown that minimizing (32) is equivalent to the wavelet shrinkage algorithm of Donoho and Johnstone defined in [21]. If we consider f = j,k,ψ c f j,k,ψ ψ j,k and u = j,k,ψ c u j,k,ψ ψ j,k , then (32) is equal to: which implies the minimization of for each i, j, k, ψ. (34) is equivalent to minimize E(s) = |s| p + 1 2λ |s − t| 2 . In the special case of p = 1, the exact minimizer is given by s = t |t| max(|t|−λ, 0). Thus the wavelet coefficients c ui j,k,ψ are shrinked toward zero by an amount of λ to obtain the exact minimizer. This corresponds to the wavelet shrinkage algorithm of Donoho and Johnstone. Thus, wavelet shrinkage can be obtained with the variational model (32). And the minimizer of the minimization problem (32) can be obtained by a projection algorithm, whose the closed convex set K B 1 1 (L1(Ω;R M )) (with p = 1, α = 1) is defined as follows: |c si j,k,ψ | , and so: K B 1 1 (L1(Ω;R M )) = s ∈ L 2 (Ω; R M )| ∀i, j, k, ψ, |c si j,k,ψ | ≤ 1 .

3.5.
Color image decomposition. In this section, we consider the color image decompsition task introduced by Aujol and Kang introduced in [4] (another relevant work is [39]). Image decomposition aims at splitting an image into a smooth/geometrical component and a textural/oscillating component. This decomposition can be used to enhance image processing tasks such as the denoising [3], segmentation [10] or inpainting [7]. The model of Aujol and Kang considers a regularized version of the vectorial TV norm, Ω |∇u| 2 + ǫ dx, ǫ > 0, which induces a slow decomposition process as we explained in Section 2.1. We will use the vectorial dual approach introduced in this paper to speed up the color image decomposition process.
is given by where J u , J v are the functionals whose describe the smooth and textural parts of f . Since we are interested in the image decomposition algorithm proposed by Aujol and Kang, then J u = ||u|| BV (Ω;R M ) and J v = ||v|| BV ⋆ (Ω;R M ) such that: We observe that the norm ||.|| BV ⋆ (Ω;R M ) was introduced by Aujol and Kang in [4] and it corresponds to the multidimensional version of the G-norm introduced by Meyer in [30] to capture oscillating patterns in gray-scale images.
Proof. It is easy to prove (37) using the following result. The minimizer of J BV (u)+ 1 2 ||u−g|| 2 is also the minimizer of J BV ⋆ (z 1 )+ 1 2 ||z 1 −g|| 2 with z 1 = Π KBV (g). Hence, the minimizer of J BV ⋆ (u) + 1 2λ ||u − g|| 2 is the minimizer of J BV ⋆⋆ (z 2 ) + 1 2 ||z 2 − g|| 2 with z 2 = Π λK BV ⋆ (g). Since we have J BV ⋆⋆ = J BV , this implies that we also have z 2 = (1 − Π λKBV )(g). This means that Π λK BV ⋆ = 1 − Π λKBV , which gives (37). 3.6. Color image deblurring. We consider the fidelity term like ||f −Au|| 2 , where A : X → X is a linear self-adjoint and positive operator satisfying ||A|| < 1. This new fidelity term allows to deal with the image deblurring problem. We remind that image deblurring aims at compensating the motion of the camera during the exposure, see [18,25] for gray-scale TV deblurring. Given a color blurred and noisy image f , we propose the following variational color image deblurring and denoising model: is computed via the following regularized energy: Since (38) is convex, the minimizer is given by: Proof. The proof is immediate using the VROF model and computing the Euler-Lagrange of (39) w.r.t. v. Figure 9 presents a color image deblurring and denoising result. The deblurring process is a non-blind deblurring process since the kernel A is assumed to be known.

3.7.
Generalized projection algorithm of chambolle. Based on the work of Bect, Blanc-Feraud, Aubert and Chambolle [5], the regularization operator can be generalized to any linear operator for color images. Chambolle) In the case of a regularization term defined as J(u) = ||Qu|| 1 , where Q is a linear operator Q : X → Y , the projection Π λK onto the convex set closed set λK can be numerically computed by a fixed point method. We consider vectors p ∈ {p ∈ Y : |p| ≤ 1} such that:
Proof. The proof is similar to the proof of Theorem 2.3.1. We have: inf and sup can be swapped using the minimax theorem. The Euler-Lagrange The constrained Euler-Lagrange equation is −Q(Q ⋆ p−f /λ)+αp = 0. Following Section 2.3, we have (40). Besides, Section 2.3 shows that the temporal step that guarantees the stability is δt ≥ 1/τ 2 with τ = ||Q ⋆ ||. Application 2: VTV on Manifold. Using Theorem 3.4, the VROF model can also be used to denoise color/vectorial images defined on a manifold Σ: where λ > 0, f is the given (noisy) vectorial image on Σ, ||u|| BV (Σ;R M ) = Σ |∇ Σ u|, ∇ Σ is the gradient operator on Σ and ||f − Au|| 2 L 2 (Σ;R M ) = Σ |f − u| 2 dΣ where dΣ is the invariant-area element. It follows from the Stokes' theorem (integration by parts on the manifold Σ) that < ∇ Σ u, p >= − < u, ∇ Σ · p >, where ∇ Σ · is the divergence operator on Σ, which is the adjoint of the gradient operator ∇ Σ . This implies that the solution of (42) is equal to: where p is the steady state of the following flow: given Q = ∇ Σ and Q ⋆ = −∇ Σ ·.
Gradient operator and divergence operator in (43) and (44) can be defined in a two different ways. If the manifold Σ is defined as a parametric manifold, then the gradient of a function h : Σ → R and the divergence of a vector field V = (V x , V y ) : Σ → R 2 in local coordinates are given by [28]: where µ, ν = {x, y}, g µν is the inverse of the first fundamental form g µν of the parametric manifold Σ and √ g is the square root of the determinant of g µν .
If the manifold Σ is represented by an implicit level set function as described in [29], then the gradient and the divergence operator are given by: where φ Σ is a level set function representing the surface Σ := {x : φ Σ (x) = 0} and P ∇φΣ := I − ∇φΣ⊗∇φΣ |∇φΣ| 2 is a projector operator onto the plane orthogonal to ∇φ Σ , which is the normal to the surface Σ.
One interesting application of the previous denoising model is for omnidirectional images, which contain information of a scene in all directions simultaneously without rotating a camera. Figure 10 presents the denoising of a color omnidirectional image. The manifold considered is a paraboloid. Its metric, developed in [9], is: which implies the following components of the gradient operator: and the divergence operator:

Conclusion.
In this paper, we have presented a regularization algorithm for color/vectorial images based on the dual definition of the vectorial TV norm. Unlike the CGM model [15] or the DS model [20], the dual approach of Chambolle [13] extends nicely to a variety of vector models, producing useful results which are better than standard methods in many cases. Indeed, this regularization algorithm offers important advantages. Firstly, it minimizes the exact TV norm unlike usual approaches such as [8,4,16]. Secondly, the numerical scheme to carry out the minimization process is easy to code. Then, the number of iterations to reach the steady state solution is low, which offers a fast color denoising model compared with other related models such as [8,33]. Finally, the dual formulation of the vectorial TV norm offers well-founded mathematical theorems to prove the existence of a minimizing solution.
We have compared our approach to related vector denoising models such as the model of Sapiro and Ringach in [33]. In our case, a fast and easy to implement denoising scheme is proposed. We have also proved the existence of a minimizer in the vector valued BV space. In [8], Blomgren and Chan propose a different version of the vectorial TV norm, which produces less color smearing than our model. However [8] does not minimize the exact TV norm and the minimization scheme is slow as shown in Figure 3. In [35], Sochen, Kimmel and Malladi propose a variational model which shares some connections with our model except the existence of a cross-correlation term in their energy which does not exist in ours. However, the cross-correlation term vanishes for gray-scale images, which induces the equivalence of our model with the Beltrami model. This also means that the Beltrami flow can be "dualized" for gray-scale images. Besides, it is difficult to study the existence of a denoising solution in the Beltrami framework unlike the TV framework.
We have also showed that the proposed fast vectorial regularization scheme can be easily applied to several standard color image processing such as the color inverse scale space, color denoising with the chromaticity-brightness color representation, color image inpainting, color wavelet shrinkage, color image decomposition, color image deblurring, and color denoising on manifolds.
Finally, and maybe more importantly, we would like to emphasize that this proposed L 1 vectorial regularization algorithm can be used in problems that required vector field regularization while preserving discontinuities. A potential interesting application of vector field denoising is in medical imaging where the proposed regularization process can be used to denoise Diffusion Tensor Images (DTI). Figure  11 presents the regularization of a synthetic vector field that has been damaged by noise. Note that the regularization is able to restore the four regions with same orientation although it is very difficult to see them in Figure 11(e). Other applications may be optical flow regularization [26], image colorization [27], image inpainting [37], multiscale hierarchical decomposition of color images [36], etc.    . Comparison between our vectorial denoising ROF model (red plot) and the denoising ROF model (blue plot) proposed by Blomgren and Chan in [8]. Both models are applied on Figure 2. x-axis represents the number of iterations and y-axis the Signal-to-Noise Ratio (SNR). Our numerical scheme, based on the dual approach, converges much faster to the denoising solution than the original model [8]. Our model converges in 60 iterations (1.8sec) whereas the original model needs more than 2000 iterations (∼500sec) to reach the denoising solution. We also notice than the SNR of the original model is slightly smaller than the SNR obtained by the dual approach. This is due to the regularization of the vectorial TV norm in the standard approach [8], which does not exist in the dual approach.         Figures (b,c) present the image which lost 56% of original information and the inpainting result. The inpainting process takes 2.3sec and 60 iterations and the picture size is 250x303. Figures (d,e) show the image which lost 79% of original information and the inpainting result. The inpainting process takes 5.9sec and 150 iterations. (Note: All pictures are in color.)  (c) Noisy Omnidirectional Image on the Parametric Plane.
(d) Zoom In the Noisy Image.
(e) Zoom In the Denoised Image.  Figure 11. Denoising of a vector field with the vectorial ROF model. The denoising process takes less than 0.1sec, 30 iterations for #1 and 100 iterations for #2. The picture size is 30x30. We can notice the good performance of our denoising model with the noisy vector field (e), which has lost the main directions (d). Our model is able to recover the main structure (f).