An Effective Directional Residual Interpolation Algorithm for Color Image Demosaicking

In this paper, we propose an effective directional Bayer color filter array (CFA) demosaicking algorithm based on residual interpolation (RI). The proposed directional interpolation algorithm aims to reduce computational complexity and get more accurate interpolated pixel values in the complex edge areas. We use the horizontal and vertical weights to combine and smooth color difference estimations. Compared with four directional weights in minimized Laplacian residual interpolation, the proposed algorithm not only guarantees the quality of color images but also reduces the computational complexity. Generally, the directional estimations may be inaccurately calculated because of the false edge information in irregular edges. We alleviate it by using a new method to calculate the directional color difference estimations. Experimental results show that the proposed algorithm provides outstanding performance compared with some previous algorithms, especially in the complex edge areas. In addition, it has lower computational complexity and better visual effect.


Introduction
In recent years, many more people choose digital cameras to take pictures.Each digital color image needs at least three color samples at each pixel.To achieve the requirement, cameras have to set three sensors to get red (R), green (G), and blue (B) components at each pixel.However, in view of equipment costs, there is only one signal sensor set at each pixel.The sensors compose a color filter array (CFA).The CFA needs only one color to be observed and the other two missing colors can be estimated at each pixel.The process of estimation is called demosaicking [1].The method effectively reduces the density of the sensors.
Bayer pattern CFA [2] is the most common pattern for CFA, as shown in Figure 1.The G pixels are estimated on a quincunx grid, while R or B pixels are estimated on a rectangular grid in Bayer pattern CFA.The density of G pixels is as twice as that of R or that of B in Bayer pattern CFA, because the peak sensitivity of human visual system (HVS) lies in the wavelength of G light.
For the original Bayer pattern CFA, the commonly used demosaicking algorithms include the nearest neighbor interpolation [3], bilinear interpolation [4], and bicubic interpolation [5].These demosaicking algorithms generate the other two missing colors at each pixel with linear averaging color pixel values.They can be widely applied in a real-time system because of the lower computational complexity.However, they have significant zipper artifacts around line edges in demosaicking images.To solve this problem, several demosaicking algorithms [6][7][8] have been proposed.Malvar [6] proposed a high-quality linear interpolation algorithm for CFA demosaicking.It presents a simple linear demosaicking filter that has better performance and lower computational complexity.Many previous demosaicking algorithms are based on inter-channel color difference interpolation.For instance, Yu [7] proposed a novel demosaicking algorithm exploiting joint intra-channel color correlation and inter-channel color difference, but it does not perform well enough when images have sharp color transition areas.Normally, G pixel values are first interpolated in most Bayer pattern CFA demosaicking algorithms, since the density of G pixel is as twice as that of R or that of B. On the contrary, Zhang et al. [8] proposed a novel demosaicking algorithm that estimated two independent chrominance components (R and B).Then they estimated the luminance component (G).As a result, the demosaicking algorithm can avoid the associated estimation error.
Although the above mentioned nonadaptive algorithms [6][7][8] have lower computational complexity for the CFA, most of them perform color directional interpolation based on the estimated gradient.However, the gradient estimations of the methods may not be robust in texture areas and edges.Based on the directional linear minimum mean square-error estimation (DLMMSE) framework, Zhang and Wu [9] proposed an adaptive algorithm to improve the problem.They assumed that the primary difference signals between G and R or G and B channels are low-pass.In both horizontal and vertical directions, the missing G pixel values are adaptively estimated by DLMMSE.Some demosaicking algorithms provide novel contributions such as square-on-point neighborhood and lattice variables [10].In [10], a corrective term and localized polling were used to improve the performance of the reconstructed images.However, the corrective term is empirically derived.Shi et al. [11] demonstrated that the input image is divided into two different area types and each type adopts different interpolation respectively.In edge areas, they recovered edge information with weights by using multidirectional components.The demosaicking algorithm has the advantages of region adaptation.Chen and Chang [12] detected edge characteristics by comparing color difference in horizontal and vertical directions.Then they used different weighted interpolation to effectively reduce the color artifacts at the edge of the image and enhance the image quality.Pekkucuksen and Altunbasak [13] exploited an orientation-free edge strength filter to estimate edge information, and applied it to reconstruction process.However, the demosaicking algorithm failed in utilizing spectral correlation.
To solve this problem, Tsai and Song [14] proposed a novel edge-adaptive demosaicking algorithm, which reduces the aliasing error in R and B channels by exploiting high-frequency information of the G channel.Adams and Hamilton [15] exploited chrominance and luminance with 5 × 5 neighborhoods interpolation equation.The advantage of method is computationally efficient in both execution time and memory storage requirements.On the basis of this interpolation method, several demosaicking algorithms [16][17][18] have been proposed.Pekkucuksen and Altunbasak [16] designed a gradient-based threshold free color filter.It combined the estimations of color difference from each direction into a gradient without setting the thresholds.To estimate different directions, in [17], the color differences were combined by using multiscale color gradients adaptively.This algorithm does not need any threshold because it makes no hard decision.Chen et al. [18] proposed an algorithm by combining voting-based edge direction detection with weighting-based interpolation, because some other algorithms may use false edge information in irregular edges.
The residual interpolation (RI) method has been widely used recently [19][20][21][22][23][24][25].The RI algorithm exploited the characteristics of the guided filter [26] and showed outstanding performance [19].Kiku et al. [20] incorporated the RI into the gradient-based threshold free algorithm.The better the interpolation results of images are, the smaller the Laplacian energies of images will have.In this way, the RI algorithm combined with minimizing the Laplacian energies of the residuals, which is termed as minimized Laplacian RI (MLRI) [21].Ye and Ma [22] exploited iterative RI process and reconstructed a highly accurate G channel to reduce reconstruction error for the first time.Wang and Jeon [23] reconstructed only R and B channels by using RI algorithm, while the missing color components of G channel were estimated by eight-direction weighted color difference interpolation.The demosaicking algorithm effectively combined the advantages of color difference interpolation and residual interpolation.Kim and Jeong [24] proposed a four-directional RI algorithm to reduce severe demosaicking artifacts.The adaptive residual interpolation (ARI) algorithm [25] adaptively combined two RI-based algorithms and selected the appropriate number of iterations at each pixel to improve the existing RI-based algorithm.However, the weakness of this algorithm is that the algorithm has high computational complexity.
Furthermore, a demosaicking algorithm was proposed for the multispectral filter array about RI in [27].There are some demosaicking algorithms based on other patterns such as RGB-white color filter arrays [28].The image reconstructed in [28] contains much important information that cannot be perceived in the color image reconstructed by conventional CFA.However, they still have other artifacts, and further improvements are needed.
RI can be used as an alternative to the widely used color different interpolation and it performs better in some image datasets.MLRI has better interpolation results because the images have smaller Laplacian energies.Furthermore, MLRI also produces other artifacts in some complex edge areas.However, we find that the directional color difference estimations are not completely the same as the actual values sometimes.With the process of research, we think the problem can be solved by reducing the weights of directional color difference estimations.This paper proposes an effective directional residual interpolation algorithm for color image demosaicking based on MLRI.We use the two weights, horizontal and vertical weights, while there are four directional weights in MLRI.In this way, our method reduces the computational complexity.Afterwards, we use a novel method to calculate the directional color difference estimations, and to combine and smooth the color difference estimations.The proposed algorithm can effectively reduce the influence of false edge information in irregular edges.In the following experiments, our proposed algorithm outperforms previous algorithms and provides better color fidelity.
The rest of this paper is organized as follows.Section 2 starts with a brief of MLRI and related work.The proposed method is described in Section 3. Experimental results of the proposed method are presented in Section 4. Conclusions and remarks on possible future work are given in Section 5.

The Outline of MLRI
The proposed algorithm is based on MLRI.First, we introduce the outline of MLRI.The interpolation process of the G pixel can be divided into three steps in MLRI.
Step 1.The G pixel values are calculated through residuals at the location of R and B pixels from horizontal and vertical directions.R and B pixel values also are calculated at the location of G pixels.The calculation of the residuals replaces the calculation of the color difference in Adams and Hamilton's interpolation equation [15].
Step 2. MLRI calculates both horizontal and vertical color difference estimations based on Step 1 at each pixel, then MLRI combines and smooths the color difference estimations.Step 3. The color difference estimations are added to the observed R or B pixel values.It aims to interpolate G pixel values.
After interpolating the G pixel, MLRI interpolates the R and B pixel values, and the method simply replaces the traditional method of using color difference as shown in Figure 2.
In the past decades, color difference interpolation has been used gradually in image demosaicking.However, residual interpolation has been used recently in MLRI.MLRI uses a novel filter called guided filter.We compare the process of MLRI with that of the color difference interpolation in calculating R pixel.The process of calculating B pixel is similar.
We introduce the process of color difference interpolation to help us understand the process of MLRI.Firstly, the missing pixel in G image is the average of four observed G pixels around it.By this way, it can get G pixel values at each pixel to combine the interpolated G image.Secondly, the color differences between R and G are calculated at the locations of R pixel by subtracting the interpolated G pixel values from the observed R pixel values (R − G).Then color differences are interpolated.Finally, an interpolated R image can be calculated by adding G pixel values to the interpolated color image.
Figure 2 shows the interpolation of the R pixel values by using MLRI in detail according to [21].First, it can get the interpolated G image as the same method as the color difference interpolation.To minimize the Laplacian energies of the residuals, the guided image and the input image are first passed through a sparse Laplacian filter.The sparse Laplacian filter with MLRI is shown in Table 1

Guided Filter
As an example, we describe how to get the R image by the guided filter in detail.The guided filter can accurately transfer the structures of the guided image to the filtered output.It has a good performance in edge-preserving smoothing and does not suffer from reversal artifacts.Thus, it can be used in demosaicking area.Furthermore, the guided filter also can be applied in denoising, edge-preserving, dehazing, and guided feathering.
The interpolated G image can be regarded as a guided image, and the observed R image can be regarded as an input image.To minimize the Laplacian energies of the residuals, the guided image and the input image are first passed through a sparse Laplacian filter.After that, the guided image uses a mask to obtain the value at the location of input image.
The guided filter assumption is a local linear model between the guided image I and the filtering output q in a window as shown in Equation (1): where i is the pixel in w k , and w k is a window centered at the pixel k; (a k , b k ), which are assumed to be constant, are linear coefficients in w k .To determine the linear coefficients (a k , b k ), the cost function is given as: where p i is a pixel at w k in the input image.To minimize Equation ( 2), the solution is given by: where µ k and σ 2 k are the mean and variance of I in w k , |w| is the number of pixels in w k , and p k = 1 |w| ∑ i∈w k p i is the mean of p in w k .However, a pixel i is involved in all the overlapping windows w k that covers i.The value of q i in Equation ( 1) is not identical when it is computed in different windows, so we rewrite Equation (1) by: where |w| are the average coefficients of all windows overlapping i.

The Proposed Demosaicking Algorithm
Our proposed demosaicking algorithm, which is a high-performance algorithm, is based on MLRI.We designed a new method that could effectively calculate the horizontal and vertical weights.It can also combine and smooth the horizontal and vertical color difference estimations.It improves the visual effect of the complex texture region.Compared with four directional weights in MLRI, the amount of calculation in the proposed demosaicking algorithm is efficiently reduced.
It is vital to interpolate the missing G pixel values because the peak sensitivity in HVS lies in the wavelength of the green light.Since the missing G pixel values are the key to the performance of the proposed algorithm, we first estimate the G pixel values at each pixel.Then we estimate the R and B pixel values, which are based on the interpolated G pixel values.
Similarly, compared with MLRI, the interpolation processes of the missing G pixel values can be divided into three steps in the same way as that in MLRI [21].
Step 1.The linear color difference interpolation in Adams and Hamilton's interpolation equation [15] can be replaced by MLRI [21].Then the new equation estimates the G, R, and B pixels in horizontal and vertical directions, so we can get G H , G V , R H , R V , B H , and B V .Step 2. The horizontal and vertical color difference estimations are calculated, and we can generate the horizontal and vertical weights at each pixel.Step 3. To get color difference, the horizontal and vertical color difference estimations are combined and smoothed by two directional weights.As a result, the G pixel values at the location of R and B pixels are generated by adding final color difference to the observed R or B pixel values.
To sum up, the outline of the interpolation processes of the missing G pixel values is shown in Figure 3.

The Calculation Process of Directionaly Estimated Pixel Value
In Step 1, the Adams and Hamilton's interpolation equation [15] can be considered as a linear model, as shown in Figure 4a.The linear color model can be replaced by MLRI [21].The locations of R and B are symmetrical and the amount of R pixels and B pixels are same.To simplify the explanation, we explain the estimation of the R pixel values in the horizontal direction in detail, as shown in Figure 4b.The linear model in Adams and Hamilton's interpolation equation can be expressed as: where R H is the horizontally estimated R pixel value, the suffix (i, j) represents the target pixel coordinates.This interpolation equation can be considered the linear interpolation, which is shown in Equation ( 7): Here G H is the G pixel value horizontally estimated at the R pixel and G H is calculated as follows: Similarly, we can calculate G H and G V at the R and B pixels.MLRI uses the tentative estimations to replace of these estimated G pixel values as shown in Figure 4b.As a result, Equation ( 6) can be rewritten as: where R H is the tentative estimated R image.To sum up, the horizontally estimated G pixel value at the location of R pixel can be calculated by Equations ( 8) and ( 9) at first.Then the interpolated G H image is used as the guided image.
As a result, the tentative estimation R H is generated by guided filter.The residual can be calculated after subtracting the tentative estimation R H value from the observed R pixel values.After that, residuals are interpolated by Equation ( 10) at each pixel.We use the residuals to add the tentative estimation R H and generate the horizontally estimated R pixel value R H . R V , B H , and B V at G pixel can be similarly generated by Step 1.

The Calculation Process of Directional Weights
We can get an observed pixel value and two directionally estimated values at each pixel after Step 1.In this section, we use the directional color difference estimations and the absolute color difference gradients to calculate the directional weights.The process is shown in Figure 5.The directional color difference estimations can be calculated by Equation (11) in MLRI.
where ∆ H g,r is the horizontal color difference estimation between G and R image, and ∆ V g,r is the vertical color difference estimation between G and R image.The effective combination of the directional color difference estimations is the most important factor of a successful algorithm.However, the horizontally estimated R pixel value R H in Equation ( 11) is not completely the same as the actual value.There are differences between the missing actual value during imaging and the reconstructed value.To reduce the influence of differences and to get estimated values, we reduce the weight of R H in Equation (11).Compared with MLRI, we use Equation (12) to calculate the directional color difference estimations.
Similarly, we can calculate the directional color difference estimations between G and B image.As a result, the observed pixel value and the two directional color difference estimations are generated at each pixel.It provides convenience for the reconstruction of missing G pixel values.Furthermore, the directional estimation may be inaccurately calculated because of the false edge information in irregular edges.However, experimental results show that Equation ( 12) can effectively reduce this influence of it.From the above steps, we use the directional color difference estimations to calculate the color difference gradients.The absolute color difference gradients at pixel coordinates (i, j) are given by: where ∆ H (i, j) is the horizontal color difference estimation in (i, j), and ∆ V (i, j) is the vertical color difference estimation in (i, j).Then we use the absolute color difference gradients to calculate horizontal and vertical weights.As we all know, the average without directional weights usually causes aliasing in the edge and the texture areas, so it is necessary to distinguish direction to avoid aliasing.In a word, the values around the target pixel are averaged according to different weights.It can effectively reduce aliasing which is caused by the abrupt color change.Thus, the weight will become smaller when the abrupt color change becomes larger.As a result, directional color difference estimations will be given a small weight in color difference combination if there is a sharp color transition area.The weights in horizontal and vertical, w h and w v , are given as:

The Calculation Process of Estimated Pixel Values
We can get the weights and the directional color difference estimations at each pixel after Step 2.Then, the directional color difference estimations are combined and smoothed by Equation ( 15): where f is a linear filter f = ( 1/4 2/4 1/4 ).To make the sum of the weights equal to 1, the weights should be allocated in the combination and smoothness.
Once the color difference is estimated, we add it to the available target pixels to obtain the estimated G pixel value: After G image is interpolated, the R and B pixel values are interpolated by guided filter.The interpolation is described in Figure 2.

Experimental Results
In this section, 18 images in the MCMaster dataset [29] are used to evaluate the effectiveness of the proposed algorithm.These images have a resolution of 500 × 500 as shown in Figure 6.The MCMaster dataset images are more in line with the characteristics of digital color images.The proposed algorithm is compared with DLMMSE [9], LDI-NAT [29], LDI-NLM [29], VDI [18], RI [19], and MLRI [21].Especially, the guided filter is applied in the proposed algorithm, RI, and MLRI.To simulate the proposed algorithm, we use MATLAB R2015b with Inter Core i7-5500U and 2.4 GHz CPU processors.Each calculation excludes the boundary of 10 pixels around the interpolated image in order to avoid the boundary effect.
Two indexes, peak signal-to-noise ratio (PSNR) and composite PSNR (CPSNR), are used to calculate the quality of the reconstructed demosaicking images, which are defined as: and where M and N are the row and column sizes of image; I in and I out are the input image and the output image, respectively; i and j are the locations of pixels in the color plane, and k represents the color plane.
Table 2 shows the PSNR of different algorithms to G components, since G components are the peak sensitivity of HVS.The highest PSNR of each image is marked in bold.Table 3 shows the CPSNR of different algorithms.The highest CPSNR of each image is marked in bold.
As shown in Table 2, the proposed algorithm performs better than the other six algorithms on 10 out of 18 images.The average PSNR of G components in the proposed algorithm higher than MLRI and RI by 0.31 dB and 0.20 dB, respectively.Moreover, as for CPSNR in Table 3, the proposed algorithm performs better than the other six algorithms on 11 out of 18 images.The average CPSNR of the proposed algorithm is higher than that of MLRI and RI by 0.35 dB and 0.25 dB, respectively.In brief, among the three methods, RI [19], MLRI [21], and the proposed algorithm, those use the guided filter, the proposed algorithm is superior to MLRI and RI.This demonstrates that the proposed algorithm is effective for image demosaicking.Objective measures are not reliable enough to judge the performance of the proposed method sometimes.Thus, we choose Figure 6a,h to show visual comparison of the interpolated images.Figure 7 shows the visual comparison of the yellow rectangular frame in Figure 6a. Figure 8 shows the visual comparison of the yellow rectangular frame in Figure 6h.We can find that the proposed algorithm produces precise image color information in details.In Figure 7, the Figure 7i has fewer artifacts than other subfigures.In Figure 8, the Figure 8i reconstructs the dots more exactly than other subfigures.We use CPU processing time to compare computational complexity.The original source codes of [9,19,21,29] are used to record the CPU processing time under the same test conditions.The average computation time for each method is shown in Table 4.For each image, the proposed algorithm takes 1.50 s while MLRI takes 1.57 s.The proposed algorithm reduces two directional-weight calculation processes.Therefore, it shows lower complexity.LDI-NAT and LDI-NLM take more than 60 s in the same condition because of the nonlocal mean method that is used to refine the interpolation.To sum up, the proposed algorithm provides a relatively low computational complexity.

Conclusions
In this paper, an effective directional residual interpolation algorithm for color image demosaicking is presented.With the use of directional weights by residual interpolation, more precise color components are interpolated.The proposed algorithm can be an alternative to color difference interpolation and MLRI, because it provides lower computational complexity by using horizontal and vertical weights.Besides, based on calculating the color difference estimations, we use a novel method to effectively reduce the influence of false edge information in irregular edges.It provides better color fidelity.Experimental results show that the proposed algorithm has outstanding performance on more complicated edge areas with lower computational complexity.In addition, it shows better subjective quality.For further work, we will continue to focus on the improved algorithm with better visual effect on more types of images.However, this improvement does not necessarily make the estimated values closest to the true values.Therefore, future research efforts will focus on finding an optimal algorithm to remove the difference caused by the directional color difference estimations absolutely.
. Then MLRI uses the guided filter to generate the tentative estimation of the R ( R) image.The interpolated G image is used as the guided image, and the observed R image is used as the input image.After that, tentative estimation of the R image uses a mask to obtain the value at the location of R image.Then minimized Laplacian residuals are calculated at the location of R pixel by subtracting the tentative estimation R H values from the observed R pixel values (R − R).Then MLRI interpolates the residuals.Finally, an interpolated R image can be calculated by adding the R image to the interpolated residual image.

Figure 2 .
Figure 2. The interpolation of the R pixel values by using MLRI.

Figure 3 .
Figure 3.The outline of the proposed algorithm.

Figure 5 .
Figure 5.The calculation process of the directional weights.

Table 1 .
The sparse Laplacian filter used in MLRI.

Table 2 .
The PSNR (dB) of G components for the MCMaster images.

Table 4 .
Average computation time (s) on MCMaster images.