Spectral preservation fusion for remote sensing images using focus measure operators based on fast discrete curvelet transform and hyperspherical color space

Abstract. How to preserve the spectral information when enhancing the spatial details is a key issue of remote sensing image fusion. The component substitution (CS)-based fusion methods can effectively enhance the spatial details while suffering spectral distortion, and multiresolution analysis (MRA)-based methods have advantages in preserving spectral information but are not satisfactory in terms of spatial details. This paper proposes a hybrid method to integrate the advantages of CS- and MRA-based approaches. The intensity image is first obtained from an original multispectral (MS) image by hyperspherical color space (HCS) transform; then, the intensity image and original panchromatic (PAN) image are decomposed by fast discrete curvelet transform (FDCT). The focus measure operators are introduced to fused low-frequency, middle-high frequency, and fine scale subband coefficients in curvelet domain with specific fusion strategies. The final fusion image is achieved by inverse FDCT and inverse HCS transform. From the aspects of subjective and objective quality assessments, the experimental results on various types of remote sensing images including IKONOS, QuickBird, and WorldView-2 indicate that compared with existing well-known algorithms and commercial softwares, the proposed method shows obvious advantages in preserving spectral information and maintaining the spatial details.


Introduction
With the continuous development of satellite remote sensing technology, remote sensing images have greatly improved in spatial resolution and spectral resolution. However, due to the limitation of technical conditions, the spatial resolution and spectral resolution of the image are conflicting for the same imaging system. Most commercial satellites such as SPOT, IKONOS, QUICKBIRD, and WorldView provide high-resolution panchromatic (PAN) images and low-resolution multispectral (MS) images. For example, IKONOS offers 1-m high-resolution PAN images and 4-m low-resolution MS images. The spatial resolution of QuickBird's PAN and MS images is 0.61 and 2.44 m, respectively. The high-resolution PAN image reflects the spatial structure information and can fully express the detailed features of the ground objects. The low-resolution MS image has rich spectral information and is conducive to the identification and interpretation of the ground targets. More and more remote sensing applications need to combine the advantages of these two types of images to generate MS images with high-spatial resolution to meet the deeper application. [1][2][3] Image fusion is an effective way to solve this problem, which is to spatially register the image data of the same area acquired by different types of sensors, and then use certain methods to combine the advantages or complementarity of each image to generate new images. In general, the spectral resolution of MS images is high, but the spatial resolution is relatively low, that is, the spatial detail expression ability is poor; the PAN image has high spatial resolution but the spectral resolution is low. Therefore, the fused image has high spatial detail performance while preserving the spectral characteristics of the MS image.
Remote sensing image fusion can be divided into three levels, i.e., pixel level, feature level, and decision level. Pixel level is the lowest level. The image is first preprocessed such as spatial registration, and then weighted and summed. The obtained value is the pixel value of the new image on the coordinate. It is mainly to increase the useful information content in the image with higher accuracy, but the amount of information that needs to be processed is too large. Feature level fusion is to fuse at the feature extraction stage. Feature extraction is performed on different images, and fusion is performed on the same type of features on each image. After such fusion, the required thematic image features can be extracted from the fused image with a high degree of confidence. The fused image can provide the feature information required for decision analysis to a large extent. The disadvantage is that its fusion accuracy is relatively poor. Decision fusion is of the highest level. Each type of sensor image is first classified, and the characteristic images in each category are determined. These classification decisions are then combined into a decision tree. Because a complete decision tree requires more than two classification results, for a target, at least two sensors are needed to detect and classify targets at the same time before classification. This level is very fault-tolerant and open, but with low accuracy.
At present, most fusion methods are concentrated at the pixel level. There are four main types: (1) Component substitution (CS). The basic idea of this kind of method is to transform the original MS image into another feature space, then replace one of the components in this feature space with PAN image, and then transform back to the original space to obtain a fused image. The most typical method is the IHS and PCA-based fusion method. The IHS transform algorithm is simple and can obtain high spatial resolution, but it can only be applied to three bands. It is incapable of multiband images. Moreover, since the difference between the PAN image and the intensity image causes the change of saturation, it will cause serious spectral distortion. 4 PCA transform can fuse any multiband image, but the fused image will lose the original physical characteristics, and the fusion effect depends on the correlation between the two replaced main components. 5 In order to reduce these adverse influences, various new CS-based methods are proposed, such as Gram-Schmidt (GS) method, 6 generalized IHS (GIHS), 7 Tu's HIS, 8 hyperspherical color space (HCS), 9 the Amélioration de la Résolution Spatiale par Injection de Structures (ARSIS) concept. 10 In these methods, the PAN image plays an important role. Therefore, the fusion image has more or less spectral distortions.
(2) Multiresolution Analysis (MRA). The introduction of MRA based on wavelet transform (WT) has opened up a new way to reduce spectral distortion, because the WT can describe and extract spatial detail information and spectral information of images at different scales. [11][12][13][14] In recent years, more and more MRA methods have been proposed, such as the discrete wavelet transform (DWT), 13 à trous wavelet transform (ATWT), 15 nonsubsampled Contourlet transform, 16 and nonsubsampled Shearlet transform, 17 fast discrete curvelet transform (FDCT), 18 smoothing filter-based intensity modulation, 19 MTF-GLP with high-pass modulation, 20 and MTF-generalized Laplacian pyramid (MTF-GLP). 21 The original MS and PAN images are decomposed into several scales and the spatial and spectral information are expressed as coefficients in scales. Then, the fused image is obtained by integrating coefficients at every scale using different strategies and the inversing MRA. Compared with CS-based approaches, the MRA-based methods can improve the performance of spectral preservation, but may not show many advantages in maintaining spatial information.
(3) Hybrid fusion. In order to make full use of the advantages of the CS-and MRA-based approaches, many hybrid methods have been proposed to combine the two together. First, the spatial and spectral information of the original MS image are separated by CS, and then the spatial details are extracted by the MRA analysis method from PAN image. Spatial details are then injected into the MS image to achieve a balance between spatial information and spectral information. 22 Such approaches include a large number of research results, such as additive wavelet-based proportioning, 11 curvelet and ICA, 23 discrete contourlet transform and PCA, 24 guided filter PCA (GFPCA), 25 and so on.
(4) Model-based fusion. Model-based fusion methods could use more prior information and other useful data from original MS and PAN images, hence in some cases they may have better performances. A hierarchical Bayesian model is proposed in Ref. 26 by employing an appropriate prior distribution obtained from Hamiltonian Monte Carlo algorithm. Using the learned spectral dictionary and the known high spatial resolution data, a sparse representation-based fusion model is used to desire high quality MS images with high resolution. 27 In order to solve the problem of uncertain result with high reconstruction error, a compressive sensingbased fusion model with fewer nonzero coefficients 28 is proposed to reduce the processing time and achieve high quality fused images. Some other model-based fusion methods such as total variation-based approach, 29,30 dictionary learning-based approach, 31 and Markov random field model-based approach 32 also could guarantee the fusion quality in many applications.
Belonging to the hybrid fusion type, this paper combines the information extraction capabilities of FDCT and HCS together with the focus measure operators, to improve the spatial resolution while preserve better spectral information.
The outline of the paper is as follows: Sec. 2 briefly describes the related work including FDCT, HCS and focus measure operators. The proposed image fusion method with the FDCT and HCS is discussed in Sec. 3. In Sec. 4, the experimental results on various remote sensing images are presented. Finally, Sec. 5 discusses the conclusion of the proposed method.

Fast Discrete Curvelet Transform
The digital implementation of the first-generation curvelet is more complex and requires a series of steps such as subband decomposition, smooth partitioning, normalization, and Ridgelet analysis, and the decomposition of the curvelet pyramid also brings a huge amount of data redundancy. So Candes et al. proposed a quick and easy-to-understand curvelet transform algorithm, namely the second-generation curvelet. 33,34 The second-generation curvelet is completely different from the first-generation. The first generation curvelet's construction idea is to approximately treat the curve in each block as a straight line through enough small blocks, and then use the local ridgelet to analyze its characteristics. The second-generation curvelet and ridgelet theory have no relation, and the implementation process also does not need to use ridgelet. The only difference between the two theories lies in the abstract mathematical meaning such as tight support and framework.
The discrete curvelet transform is E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 0 1 ; 1 1 6 ; 2 5 2 cði; j; kÞ ¼ Since the shear block S T θ l ðk 1 × 2 −j ; k 2 × 2 −j∕2 Þ is not a standard rectangle, in order to use the fast Fourier algorithm, the above equation could be rewritten as E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 0 2 ; 1 1 6 ; 1 8 9 cði; j; kÞ ¼ Then, the FDCT could be implemented by the Wrap algorithm in local Fourier Transform. The wrap around the origin is the core idea of the curvelet based on the Wrap algorithm. Any region is mapped to the affine region of the origin through the periodic technique. This mapping is a one-to-one correspondence. Specific steps are as follows: Step 1: Obtaining f½n 1 ; n 2 , −n∕2 ≤ n 1 , n 2 < n∕2 through two-dimensional (2-D) FFT on f½t 1 ; t 2 ∈ L 2 ðRÞ.
Step 2: Multiplying f½n 1 ; n 2 by parabolic window U j;l ½n 1 ; n 2 for each scale and angle parameter set ðj; lÞ to localized f.
Step 4: Implementing 2-D IFFT to localized f 0 j;l and finally obtaining the discrete curvelet coefficients c D ðj; l; kÞ.

Hyperspherical Color Space-Based Fusion
Converting from the original color space to HCS follows the standard conversion from n-dimensional Cartesian space to n-dimensional hypersphere space. 9 For an image with N bands, an intensity component and N − 1 angular components in hyperspheres are obtained. The equation for the positive transformation from the N-band color space to the HCS is as follows: E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 0 3 ; 1 1 6 ; 5 8 7 8 > > > > > < > > > > > : where χ i is the i'th band of the ogirinal n-bands image. The inverse transform is E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 0 4 ; 1 1 6 ; 4 7 5 8 > > < > > : HCS-based image fusion is similar to the IHS fusion method, which belongs to the CS fusion, and is divided into naïve mode and smart mode. The procedures of the naïve mode are: (1) Calculating the sum of square of MS and PAN images E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 0 5 ; 1 1 6 ; 3 6 5 (2) Calculating the matched PAN E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 0 6 ; 1 1 6 ; 3 0 9 where μ 1 , σ 1 , μ 0 , and σ 0 denote the mean and variance of P 2 and I 2 .
E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 0 7 ; 1 1 6 ; 2 5 8 (3) getting the fused image by substituting I in HCS transform as I adj , and applying the inverse HCS transform.
Step always equals 1. The SML is E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 0 9 ; 1 1 6 ; 7 2 3 where N is the window size used to compute the SML focus measure.
(2) Energy of gradient (EOG) Denoting the horizontal and vertical gradients of an image are f x ¼ fðx þ 1; yÞ − fðx; yÞ and f y ¼ fðx; y þ 1Þ − fðx; yÞ, then the EOG is E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 1 0 ; 1 1 6 ; 6 1 9 (

3) Energy of Laplacian of image (EOL)
For analyzing the high spatial frequencies and the edge sharpness of an image, the EOL is proposed as E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 1 1 ; 1 1 6 ; 5 3 1 where E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 1 2 ; 1 1 6 ; 4 8 0 (4) Tenengrad (TEN) TEN is proposed by Tenenbaum, which employs Sobel operator to describe the gradient energy E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 1 3 ; 1 1 6 ; 3 7 0 where ∇Sðx; yÞ ¼ ½∇S x ðx; yÞ 2 þ ∇S y ðx; yÞ 2 1∕2 , T is the threshold, ∇S x ðx; yÞ and ∇S y ðx; yÞ are the horizontal and vertical Sobel gradients, respectively: E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 1 4 ; 1 1 6 ; 2 9 5 E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 1 5 ; 1 1 6 ; 2 4 4 (5) Spatial frequency (SF) 39 SF is a modification of EOG with the following expression: where E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 1 7 ; 1 1 6 ; 1 0 9 3 Proposed Method

Overall Workflow
The fusion scheme based on FDCT and HCS transform proposed in this paper is shown in the figure below. It is mainly divided into the following steps: first, the original MS image is transformed by HCS to obtain the intensity component and multiple angular components. Then, the intensity component and original PAN image are subjected to FDCT to obtain low-frequency, midhigh frequency, and fine-scale subband coefficients. Different fusion strategies are adopted for each subband coefficient to obtain new low-frequency, midhigh frequency, and fine-scale subband coefficients. According to the characteristic that low-frequency energy concentrates in the intensity image, low-frequency subbands are first calculated by focus measure operators, and then weighted by variance. According to the characteristics of rich high-frequency information in the midhigh frequency and fine-scale subbands of the PAN image, based on the calculation of focus measure operators, the midhigh frequency subbands are fused using the correlation coefficient matching strategy, and the fine-scale subband directly adopts the largest absolute value fusion strategy. Then, the new low frequency, midhigh frequency, and fine-scale subband coefficients are inversely transformed by FDCT to obtain a new intensity component. Finally, the HCS inverse transformation is performed on the new intensity component and the angular components obtained from the HCS transformation of the original MS image in the first step, and the final fused image is obtained (Fig. 1).

HCS and FDCT Transform
First, the original MS image is resampled to the same pixel size as the original PAN image, and an intensity component and multiple angular components could be obtained by applying HCS transform to the resampled MS image. Then, the FDCT is used to the intensity component and original PAN image, respectively, and low-frequency, midhigh frequency, and fine-scale subbands coefficients for intensity component and PAN image are obtained.
Based on the above procedure, the coefficients fusion is performed using different strategies for low-frequency, midhigh frequency, and fine-scale subbands.

Fusion of Low Frequency
The low-frequency coefficient of the image reflects the energy distribution of the image. According to the imaging principle of the MS sensor and the panchromatic sensor, the MS image reflects the spectral characteristics of the target, and the PAN image reflects the spatial details of the target. The gray-scale distribution characteristics of the same target in both images are often very different. Therefore, the traditional average weighting strategy for low-frequency coefficients does not consider the physical characteristics of the sensor, which will increase the spectral distortion of the fused image and lose some useful information. Therefore, this paper first calculates the metric of the low-frequency coefficients by focus evaluation operator, and then uses the local variance weighted strategy to calculate the low-frequency coefficient of the fused image, as E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 1 9 ; 1 1 6 ; 4 6 4 C l fus ðx; yÞ ¼ C l I ðx; yÞ þ σðW pan Þ σðW pan Þ þ σðW I Þ × ½C l pan ðx; yÞ − minðW pan ; W I Þ; where C l I , C l pan , and C l fus represent the low-frequency subband coefficients of the HCS transformed intensity image, PAN image, and fused image, respectively. W represents the local window, σ represents the variance in the local window of the focus measure operator metric value, W represents the mean value in the local window of the original low-frequency coefficient, and minðW pan ; W I Þ are the minimum value of the local window mean value of the original lowfrequency coefficient which represents the common low-frequency characteristic of the MS intensity image and the PAN image. C l pan ðx; yÞ − minðW pan ; W I Þ represents the characteristic low-frequency characteristic of the PAN image. σðW pan Þ σðW pan ÞþσðW I Þ is a weighting coefficient and is used to inject the low-frequency characteristic specific to the PAN image to the MS intensity image. The variance σ is a measure of the difference in energy distribution. After the calculation of the focus measure operator, the variance can better reflect the energy distribution of the lowfrequency coefficients. It can be seen from the equation that the low-frequency coefficients of the fused image consist of the low-frequency coefficients of the MS image and the unique lowfrequency coefficients of the PAN image, so the entire information of the MS image can be completely preserved and the unique information of the PAN image can be fully utilized. This strategy can achieve maximum preservation of the low-frequency information of original images.

Fusion of Middle-High Frequency
The middle-high frequency coefficients of the image reflect the distribution of the details of the image. The PAN image has rich high-frequency detail information. Relatively speaking, the detailed features of MS intensity images are not obvious. Therefore, this paper first calculates the metrics of the focus measure operator of a middle-high frequency subband, and then uses local correlation coefficients to calculate the matching degree of the metrics, and accordingly chooses the appropriate middle-high frequency coefficient. The local correlation coefficient is defined as where A and B represent two vectors, μ A and μ B denote the mean of A and B, respectively. The value of CC is the matching degree of the middle-high frequency coefficients between the MS intensity image and the PAN image in a certain neighborhood. When selecting the appropriate middle-high frequency coefficients, a threshold T is first set, which is defined as the average of the correlation coefficients of all subbands, i.e., where n is the number of middle-high frequency subbands, C mh I ðiÞ and C mh pan ðiÞ are the i'th subband of the MS intensity image and the PAN image, respectively. CC½C mh I ðiÞ; C mh pan ðiÞ represents the correlation coefficient of the i'th subband. The CC value of a particular subband is compared with the threshold T. CCðiÞ < T indicates that the matching degree between the two subbands is very low, hence it is necessary to retain both information at the same time. At this time, the middle-high frequency coefficient of the fused image is defined as E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 2 2 ; 1 1 6 ; 5 0 0 C mh fus ðx; yÞ ¼ C mh I ðx; yÞ þ C mh pan ðx; yÞ: When CCðiÞ ≥ T, the matching degree is high. In order to reduce the redundancy of information, only high-middle frequency coefficients with relatively large absolute values are preserved. At this time, the middle-high frequency coefficients of the fused image are defined as

Fusion of Fine Scale Frequency
The fine-scale subband coefficients reflect the highest level of detail in the image. Since the spatial resolution of an MS image is much lower than that of a PAN image, the energy of its subband coefficients is much lower than that of a PAN image. Of course, it is not ruled out that in some subbands, the energy is higher than that of a PAN image. Therefore, for the fine-scale subband coefficients, this paper only needs to use the approach with the largest absolute value to keep the details as much as possible. Its fusion strategy is E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 2 4 ; 1 1 6 ; 2 6 3 C fine fus ðx; yÞ ¼ C fine pan ðx; yÞ; if jC fine pan ðx; yÞj > jC fine I ðx; yÞj C fine I ðx; yÞ; otherwise :

Fusion Result
Through the above subband fusion strategies, the fused subbands of the low-frequency subband, the middle-high frequency subband, and the fine-scale subband are obtained. All of the fused subbands were subjected to FDCT inverse transformation to obtain a fused MS intensity image. The HCS inverse transformation is performed on each angular component obtained by HCS transformation of original MS image and the newly generated intensity image, and the final fused image can be obtained.

Experimental Setup
Three types of remote sensing images, i.e., IKONOS, QuickBird, and WorldView-2, are used to test the performance of the proposed method and the results are compared with existing wellknown fusion algorithms and algorithms from commercial softwares. The basic information of these three types of images are listed as following ( Table 1): The size of images used in this paper is all 2048 pixels × 2048 pixels and the MS images are all resampled to the same pixel size with PAN images using bicubic interpolation algorithm. Some subsets images with size of 512 pixels × 512 pixels are shown in Fig. 2.
The GIHS algorithm is selected since it is the generalized representative CS algorithm. To compare the performance between with/without focus measure operators, the same framework named "FDCT + HCS" ("C + HCS" for short) is used without applying focus measure operators before low-, midhigh frequency, and fine scale fusion procedure. The GS algorithm from ENVI software, HCS algorithm from ERDAS IMAGINE software, and UNB algorithm from PCI GEOMATICA software are used to compare the fusion performance since they are three famous and mainstream commercial softwares in remote sensing field and can provide fair fusion results to be compared. Meanwhile, for a compressive comparison, six methods with outstanding Table 1 The basic information of three types of images used in this paper. PAN   (a)

Fusion Quality Assessment Indices
Quality assessment indices are used to objectively compare the performance of all employed fusion algorithms. There are many types of indices which can be divided into three categories, i.e., spatial quality indices, spectral quality indices, and overall quality indices. According to Wald's protocol, two modes with and without reference images are used in this paper. The first mode treats the original MS image as reference image. The original MS and PAN images are first downscaled using a bicubic resampling algorithm, for example, the original MS image with 4-m spatial resolution and PAN image with 1-m spatial resolution are downscaled to the MS image with 16-m spatial resolution and PAN image with 4-m spatial resolution. Then, the downscaled MS and PAN images are fused using some fusion algorithm, and finally the fused result is compared with the original MS image. The basic assumption and principle of this mode is that the fusion algorithm should keep consistency for image fusion procedure. For the second mode, the original MS and PAN images are directly fused in their original spatial resolution and the quality assessment could be conducted without references. The indices for first mode include four spectral indices, i.e., Erreur Relative Golbale Adimensionnelle de Synthese (ERGAS, or relative global dimensional synthesis error), 40,41 spectral angle mapper (SAM), 42 spectral correlation coefficient (CC), 43 Universal Image Quality Index (UIQI), 44 two spatial indices, i.e., spatial frequency index (SFI) and peak signal-to-noise ratio (PSNR), two overall indices, i.e., qquality AB/F (QABF) 45 and structure similarity index (SSIM). 46 The spatial quality indices describe the capacity that the fusion algorithm maintain the spatial information of the PAN image, hence the greater values of SFI and PSNR indicate that the fused image have high spatial quality. The spectral quality measures the performance that the fusion algorithm transmits spectral information of the original MS image to the fused result, therefore, the lower ERGAS, SAM, higher CC, and UIQI show that more spectral information are kept. However, there is no one algorithm that could get higher spatial quality and higher spectral Note: The value in bold denotes the best value of the method. The GS result shows much higher ERGAS and SAM, much lower CC and UIQI indicating that it has the worst spectral information, which is also obvious as seen in Figs. 3 and 4. The spectral quality of GIHS result is worse than the other methods since its spectral quality indices are worse. The C + HCS and HCS show similar quality indices values indicating that they show the similar fusion performance. BDSD has the highest PSNR while PRACS shows worse performance than MTF-based methods due to the over sharpened details. There are a little differences between five methods in the proposed framework, and they have better spectral quality than the other existing fusion methods. Benefiting from the better preservation of spectral information, the overall qualities of five methods in the proposed framework are also better than the other methods, which can be seen from higher values in SSIM, QNR, and GD.
(a) quality. Hence, the overall quality is used to strike a balance between spatial and spectral quality, and the higher the QABF and SSIM are, the better quality the fused image shows. The second mode includes indices named quality with no reference (QNR) 47 and Gaussian scale space-based index (GD). 48 The higher QNR and GD values indicate that the better fused image is achieved.

Experiment on IKONOS Dataset
The original MS/PAN and 10 fused results are shown in Fig. 3. The original size of all images are 2048 pixels × 2048 pixels and here only 512 pixels × 512 pixels subset images are shown. Also the 128 pixels × 64 pixels subset images labeled in the red rectangle in all original and fused images are enlarged in Fig. 4  3(k), 3(m) and 4(g), 4(i), 4(k), 4(m). The objective quality assessment results using 10 indices are shown in Table 2, where "SML" represents "FDCT + HCS + SML," and the other "EOG, EOL, TEN, and SF" are similar as "SML." Figure 5 shows the fusion results on QuickBird dataset, where the size of original images is 2048 pixels × 2048 pixels and only 512 pixels × 512 pixels subset images are shown. Figure 6 is the enlarged areas along the road in Fig. 5 in red rectangles. There are four bands in the original MS and fused images, and the red, green, and blue bands are shown as true color image. Similar to Figs. 3 and 4, the GIHS results in 5(c) and 6(c) show lighter green and the GS results in 5(d) and 6(d) show yellow green in grass and tree areas. Some areas are blurred in the C + HCS results especially along the road which can been seen obviously in Fig. 6(e). The color of UNB and the five methods of the proposed framework are similar to each other, except for some sharper edges in UNB results. The ATWT and MTF_GLP_HPM results are a little over sharpened. The objective quality assessment results are shown in Table 3.

Experiment on QuickBird Dataset
It can be seen from Table 3 and Fig. 7 that the spectral qualities of GIHS, GS, HCS, ATWT, BDSD, and MTF_GLP_HPM are worse than the other methods since they have higher ERGAS and SAM, and lower CC and UIQI. For spatial quality, the GS, HCS, and UNB results are better than the other methods because they have sharper edges. From the overall quality aspect, the five methods in the proposed framework as well as MTF_GLP show better performance than the other existing methods, since they have better spectral quality and relatively better spatial quality which can be seen from spectral and spatial quality indices.

Experiments on WorldView-2 Datasets
The WorldView-2 images with size of 2048 pixels × 2048 pixels are tested using 10 fusion methods, and the results are shown in Figs. 7 and 8. Similar to the previous experiments on (a) IKONOS and QuickBird dataset, only 512 pixels × 512 pixels subset and 128 pixels × 64 pixels enlarged subset with true color are shown. The GS results in Figs. 7(d) and 8(d) are obviously overenhanced and the color distortion is much more than the other methods. Except the GS result, the other methods show similar spectral information and spatial details. In addition, the UNB, ATWT, and MTF_GLP_HPM results are a little sharper in edge areas. The objective quality assessment results are shown in Fig. 9.
The quality assessment results are consistent with the subjective observations that the GS result shows the worst overall quality which can also be seen from highest ERGAS and SAM values, and the other methods have similar spectral quality. The five methods in the proposed framework still have the best spectral performance resulting in the best overall qualities.
From the above experiments on three typical remote sensing image datasets with various spatial resolutions and different bands, we can see that the proposed framework could effectively preserve spectral information while maintaining spatial details resulting in the outstanding overall qualities. The introduction of focus measure operators could effectively improve the fusion performances compared with the FDCT + HCS framework, because the focus measure operators are helpful in aggregating energies in every frequency subbands of curvelet domain to enhance the spatial details, and finally achieving better overall performances together with the spectral preservation capacity of FDCT + HCS framework.    This paper proposed a hybrid spectral preservation image fusion method combing with CS-and MRA-based approaches to try to achieve a balance between spatial details and spectral information. The HCS transform is used to separate the spatial and spectral information into an intensity component and multiple angular components. In order to inject the spatial details of PAN image into the intensity component, the FDCT is employed to decompose the intensity and PAN images into low-frequency, middle-high frequency, and fine scale subbands. By introducing the focus measure operators to concentrate energies of subband coefficients, each subband is fused to generate a new intensity component according to different strategies, i.e., the low-frequency subbands use the local variance weighted strategy, the middle-high subbands use the local correlation coefficients matching strategy, and the largest absolute value strategy is used for fine-scale subbands. The final fusion results are obtained by performing inverse HCS transform using the newly generated intensity component and original angular components. Various remote sensing images including IKONOS, QuickBird, and WorldView-2 are used to test the performance of the proposed method and existing well known and commercial softwares. According to Wald's protocol, some widely used objective quality assessment indices including ERGAS, UIQI, QNR, and so on are used to compare the performance through the aspects of spatial quality, spectral quality, and overall quality. The experimental results indicate that the proposed method can better preserve the spectral information while enhancing the spatial details.
The experimental results also show that the focus measure operators are helpful to improve the fusion performance. The future work will focus on developing more accurate strategies to further enhance the fusion capacities.