Image Fusion Based on the ∆ − 1 − TV 0 Energy Function

This article proposes a ∆−1− TV0 energy function to fuse a multi-spectral image with a panchromatic image. The proposed energy function consists of two components, a TV0 component and a ∆−1 component. The TV0 term uses the sparse priority to increase the detailed spatial information; while the ∆−1 term removes the block effect of the multi-spectral image. Furthermore, as the proposed energy function is non-convex, we also adopt an alternative minimization algorithm and the L0 gradient minimization to solve it. Experimental results demonstrate the improved performance of the proposed method over existing methods.


Introduction
With the application of multi-source imaging in many areas, such as remote sensing [1,2], medical imaging [3] and quality and defect detection [4], image fusion has become an attractive and important research area in image processing.More specifically, in remote sensing, there is a significant interest in the fusion of high spatial and high spectral resolution images to provide a better description and visual representation of a scene.However, typically, both high spatial and spectral resolution images are not simultaneously presented in a single image, due to commercial constraints.More specifically, panchromatic images contain high spatial resolution with limited spectral resolution, while multi-spectral images contain high spectral resolution with low spatial resolution.
Consequently, many researchers have expressed interest in the fusion of multi-spectral images and panchromatic images.The most common methods are based on the intensity-hue-saturation (IHS) transform [5,6].However, using the IHS method results in spectral degradation.To solve this problem, some algorithms based on the wavelet transform are proposed.Wavelet transform is an important tool in signal processing [7] and image processing [8,9], especially for image fusion.Núñez et al. [10] fused a high spatial resolution panchromatic image (Satellite pour l' observation de la Terre (SPOT)) with a low spatial resolution multi-spectral image (Landsat Thematic Mapper (TM)) using the additive wavelet algorithm (AW).The Atrous wavelet approximation of the SPOT panchromatic image is substituted by the bands of TM image.Li et al. [11] proposed a choose-max wavelet fusion algorithm (CMW) based on the wavelet transform.An algorithm based on the multi-scale first fundamental form (MFF) was presented by Scheunders [12], which used a multi-valued image wavelet representation method to fuse images.However, this method is prone to large wavelet coefficients.Chen et al. [13] addressed this issue by proposing a weighted multi-scale first fundamental form method (WMFF).The image fusion methods based on wavelet transform, while reporting good results, have three limitations: empirical selection of the predefined basis; selection of the decomposition level; the wavelet method focuses on preserving spectral information.In order to solve these issues, techniques, such as the T V (total variation)−L 1 model [14][15][16], were proposed as a balanced spatial and spectral information model, which use the primal-dual algorithm and Bregman-splitting algorithm to maximize the matching measurement.In the T V − L 1 model, the T V regularization preserves the discontinuity of the magnitudes of gradient difference; and the L 1 data term tends to measure the low frequency information by calculating the difference between the result and the multi-spectral image.Another class of literature utilizes population-based optimization to fuse multiple images.Saeedi [17] utilized the population-based optimization to select the weights of the low-frequency wavelet coefficients to improve the infrared and visible image fusion.Saeedi [18] found the optimal pan-sharpened image by applying the multi-objective particle swarm optimization algorithm and using the two initial pan-sharpened results generated by two different fusion rules.Finally, Lacewell [19] used genetic algorithms to get a more optimized combined image for visual and thermal satellite image fusion.While the reported literature produces good results, the multi-spectral images are prone to the block effect, and the existing literature does not address this problem.In this article, we propose to address the block effect within the multi-spectral images by using the T V 0 − ∆ −1 energy function.The ∆ −1 component not only deletes the block effect, but also retains the spectral information.Additionally, the L 0 gradient term of the proposed energy function increases the spatial resolution within the final fused result.
The contributions of our paper in image fusion can be summarized as follows: (1) The ∆ −1 term within the fusion energy function is proposed to remove the block effect in multi-spectral images without affecting the spatial information.
(2) To improve the image fusion accuracy, an alternative minimization algorithm using a non-convex L 0 regularization term is proposed.
The rest of the paper is organized as follows: Section 2 presents a detailed description of the proposed fusion algorithm.Section 3 provides experimental results and discussions on IKONOS and QUICKBIRD images.Finally, in Section 4, our conclusions are discussed.

Flowchart of the Proposed Algorithm
The flowchart of the proposed algorithm is presented in Figure 1.In the proposed algorithm, the color multi-spectral image is transformed into the IHS color model [20].The intensity component is selected as the primary fusion variable.The method uses an energy minimization process solved by an iterative process, which is based on the inverse Laplace transform and sparse norm of the gradient term.The ∆ −1 operator is proposed by [21,22] and also applied to image fusion quality assessment [16,23].Additionally, [22] computes ∆ −1 by discrete Fourier transform.The discrete ∆ operator is given by: where m and n are the pixel position of image f with size [M, N ].Applying the discrete Fourier transform on (1), we obtain: where p and q are discrete Fourier sample points and F presents the discrete Fourier transform.In order to avoid singularity, we input a small constant and obtain the ∆ −1 as follows: As shown in Figure 2, the weighting function of Equation ( 3) focuses on the low frequency component of the image f .An example of the application ∆ −1 on the IKONOS's image is shown in Figure 3.It can be seen that ∆ −1 removes the block effect, while preserving the general information.

( )
remove the block effect w hil e pr eser ving t he general information

Functional Form
Xie et al. [14] proposed a T V − L 1 fusion method as follows: where J is an energy function about the fusion result R, G is a panchromatic image, T is the intensity part of IHS transformation of a multi-spectral image, λ is a weighted parameter and ∇• 1 stands for a discrete total variation norm.If we use the anisotropic total variation form, we can transform Equation (4) as follows: min where R x is the partial derivative of R with respect to x, R y is the partial derivative of R with respect to y, G x is the partial derivative of G with respect to x and G y is the partial derivative of G with respect to y.In recent years, many studies have empirically demonstrated through experiments that non-convex potential functions are suitable for optimization-based computer vision problems.Krishnan [24] used a normalized sparsity measure similar to L 0 for the image blind deconvolution model.A L p norm [25] is proposed to solve the optical flow and image denoising model.Fu and Zhang [26] proposed an adaptive non-convex model to solve the image restoration model.Consequently, in our article, we investigate the use of a sparse L 0 norm to increase the image details of the fused result.A non-convex L 0 norm is substituted for the L 1 norm of anisotropic total variation of the T V − L 1 model as follows: where R y − G y 0 = ⊕ p, |R y,p − G y,p | = 0 , p is the pixels and ⊕ is the counting operator.We use the L 0 norm of the gradient to transform the spatial information, which can retain the spectral information, whose tendency is inconsistent with the spatial information.While the T V − L 1 model demonstrates good results for multi-spectral fusion, the final fused images contain the block effect, as discussed earlier.To address this issue, we apply ∆ −1 on the data fitting term and substitute where R is the fused result, which retains the spatial information and the spectral information.Moreover, we employ the ∆ −1 operator, within T V − L 1 as ∆ −1 R − ∆ −1 T 1 , resulting in the function focusing on the low frequency information of R − T and neglecting the high frequency information.Finally, as Model ( 7) is a non-convex model, we solve the problem using an alternative minimization algorithm, as detailed in the next subsection.

Alternative Minimization Algorithm
We use the splitting scheme similar to [20,27] to solve Equation (7).It is well-known that the L 0 problem is strongly NP-hard to solve, but we solve it based on the conclusion reported in [26].Firstly, we introduce two auxiliary variables p 1 and p 2 to transform Model (7) as follows: The alternative minimization algorithm [28] is employed to solve Model (8).The most important part of the alternative minimization is to solve: It is well known that Model ( 9) is strongly NP-hard.We use the technique of [27] to decompose it into the following model: where i is the pixel position of the image and H is defined as: Next, Model ( 9) is rewritten by Model (10).Additionally, solve Model (10) by computing each term as follows: To solve Model ( 12), we need Proposition 1.
Proposition 1. Model ( 12) obtains the minimum under the conditions: The proof is similar to [27].In order for this paper to be self-contained, we present the proof of Proposition 1 as follows. Proof.
Combining ( 14) and ( 15), we conclude that if ) Combining ( 16) and ( 17), we conclude that if 12) obtains the minimum λ β .Thus, we complete the proof of Proposition 1.Because every point is considered to be a single point individual, we simplify Equation ( 13) as follows: Although Proposition 1 can solve the L 0 approximately, another sub-problem of the alternate minimization is to solve as follows: Based on the definition of ∆ −1 and the discrete Fourier transform, we transform (19) as follows: min where • represents the element-wise multiplication.Based on (3), Equation ( 20) can be transformed as follows: Then, we obtain the solution: where † represents the conjugate transpose operator, ./ is the element-wise division, F −1 is the inverse Fourier transform and (∂ x , ∂ y ) are differentiated operators.
After we solve the above important sub-optimization, we present our framework by Algorithms 1 and 2 as follows:

Algorithm 1 Image fusion algorithm
Input: a multispectral image M , a panchromatic image G, parameters λ, β 0 , β max and v Initialization: where i is the pixel position of the image.

Experimental Results
Some experiments on IKONOS and QUICKBIRD images are used to validate the proposed algorithm.We use two metrics to quantify the performance.The correlation measure (CM) [10] is defined as follows: where A and B are the images in the lexicographic order.CM computes correlation coefficients of the red, the green and the blue channels between the multispectral image and the fused result; this computation can be used to assess the preservation of the spectral information of the result.
A feature-based metric Q AB/F [29] is adopted.The Q AB/F metric is defined as follows: where s x A (m, n) and s y A (m, n) are the outputs of the horizontal and vertical Sobel templates centered on pixel (m, n) and convolved with the corresponding pixels of image A. Similarly, g B (m, n) and g F (m, n) fundamental forms fusion method (MFF) [12] and the weighted multi-scale fundamental forms fusion method (WMFF) [13].In order to illustrate the superiority of the proposed L 0 gradient minimization, we perform a comparative experimentation with L 1 gradient minimization-based methods [14,15], and the results are shown in Tables 1 and 2.
From the enlarged portion of the fused result, Figure 5 shows that our proposed algorithm has removed the block effect.We selected the parameter to make CM (r) metric, as close as possible, which means the retained spectral information is as close as possible.We applied the above metric to quantify the results in Table 1.In Table 1, CM (r), CM (g) and CM (b) correspond to the correlation measures of the red channel, the green channel and the blue channel.This shows that our fused result preserves more spectral information than the other methods.The fourth column of Table 1 illustrates that the image obtained by the proposed method gets more feature information from the original image than other methods.The sixth column of Table 1 obtained by the proposed method indicates that our method obtains clearer results than the other methods.5. Fused IKONOS images using AW, CMW, MFF, WMFF, method of [14], method of [15] and the proposed method.
AW CMW MFF WMFF Method [14] Method [15] Proposed Method Figure 6 shows the fused QUICKBIRD's images that are obtained by CMW, MMF, WMMF, the methods of [14,15] and the proposed method.The above metric is applied to quantify the results in Table 2.In Table 2, CM (r), CM (g) and CM (b) correspond to the correlation measure of the red channel, the green channel and the blue channel; this shows that our fused result preserves more spectral information than other methods.The fourth column of Table 2 illustrates that the image obtained by the proposed method gets more feature information from the original images than other methods.The sixth column of Table 1 obtained by the proposed method indicates that our method obtains clearer results than the other methods.Figure 6.Fused QUICKBIRD images using CMW, MFF, WMFF, method of [14], method of [15] and the proposed method.

Conclusion
The proposed algorithm removes the block effect efficiently by using the ∆ −1 operator.In addition, we use a non-convex functional form to improve our fusion result.Simultaneously, our proposed method presents a new energy minimization algorithm to retain the spectral and high frequency information from the original images.We conduct some experiments on the IKONOS and QUICKBIRD images, and the proposed algorithm obtains a better fusion result than other fusion algorithms, in terms of not only visual quality, but also some fusion indicators.

Figure 4 .
Figure 4. Generation of the parameter λ = A.

Table 2 .
Quantitative comparison results of Figure6.