Efficient High-Quality Demosaicing Algorithm for Color Filter Array

(1)

國立台灣科技大學

資訊工程系

National Taiwan University of Science and Technology

Department of Computer Science and Information Engineering

Efficient High-Quality Demosaicing Algorithm for Color Filter Array

鍾國亮 Kuo-Liang Chung 楊偉楨 Wei-Jen Yang

顏文明 Wen-Ming Yan 王中周 Chung-Chou Wang

Technical Report: NTUST-CSIE-07-04

(2)

Abstract Most modern digital cameras use a signal CCD/CMOS sensor with a color filter array (CFA) to capture a color image in which each pixel has only one primary color, and the captured image is called the mosaic image. Without demosaicing processing, this paper first proposes a new approach to capture more accurate color gradient/edge information on mosaic images directly. Next, based on spectral–spatial correlation, a novel adaptive heterogeneity–projection for mosaic images is presented. Combining the captured color gradient/edge information and the adaptive heterogeneity–projection values, a new edge–sensing demosaicing algorithm is presented. Based on twenty–four popular testing mosaic images, experimental results demonstrated that our proposed high–quality demosaicing algorithm has the best image quality performance when compared to cur- rently published four algorithms by Pei and Tam, Lu and Tan, Dubois, and Tsai and Song.

Index Terms: Adaptive heterogeneity projection, Color filter array (CFA), Demosaicing algorithm, Digital cameras, Mosaic images, CPSNR, Sobel operator.

I. INTRODUCTION

Recently, digital cameras have become more and more popular in consumer electronics market. In order to economize the hardware cost, instead of using three sensors, most digital cameras capture a color image with a single sensor array based on the well–known Bayer color filter array (CFA) [1], where each pixel in the captured image has only one measured color and this kind of images is called mosaic images. Fig. 1 depicts the Bayer CFA structure. Because G (green) color channel is the most important factor to determine the luminance of the color image, half of the pixels in Bayer CFA structure are assigned to G channel. R (red) and B (blue) color channels, which share the other half pixels in the Bayer CFA structure, are considered as the chrominance signals. Fig. 2 shows the full color Lighthouse image. After filtering Fig. 2 via the Bayer CFA structure, the resultant mosaic image is shown in Fig. 3.

(3)

Fig. 1: The Bayer CFA structure.

Fig. 2: The full color Lighthouse image.

In order to recover the full color image from the input mosaic image, the demosaicing process is used to estimate the other two color channels for each pixel [9]. Bilinear interpolation (BI) [21] is the simplest demosaicing algorithm in which the unknown two color channels of each pixel are obtained by averaging its proper adjacent pixels. In [8], Gun- turk et al. presented an efficient demosaicing algorithm by using alternating projections.

Based on color correlation concept, Pei and Tam [20] presented an efficient demosaicing algorithm. In [14], Lu and Tan presented an efficient demosaicing algorithm based on exploiting spatial and spectral correlations among adjacent pixels and they also presented a quality measure. Based on primary–consistent soft–decision framework, Wu and Zhang

(4)

Fig. 3: The mosaic Lighthouse image.

[28] presented a multiple estimation to estimate unknown color values according to different edge directions. Based on projection–onto–convex–set approach, Li [13] presented the first iterative demosaicing algorithm. In [4], Dubois presented a novel demosaicing algorithm based on frequency–domain representation. Su [23] presented an improved iterative demosaicing algorithm using weighted–edge and color-difference interpolations. In [11], Lai and Liaw presented an modified mean–removed vector quantization algorithm to improve the quality image performance of the previous algorithm by Pei and Tam [20].

Based on the spatial correlation and the edge–directional information of the neighboring pixels, Lee et al. [12] presented a weighted edge–sensing demosaicing algorithm. Cur- rently, based on the Nth–order directional finite derivative of spectral–special correlation, empirically N = 11, Tsai and Song [27] presented a lined–based demosaicing algorithm and their algorithm has the best image quality performance in average when compared to all other existing demosaicing algorithms. The digital zooming issue on CFA model has been discussed in [3, 16, 17, 18, 29].

After examining most of these previously published edge–sensing demosaicing algorithms, we find that the quality of demosaiced images is heavily dependent on the captured gradient/edge information from input mosaic images, but usually the captured gra-

(5)

dient/edge information on mosaic images is not so accurate. Since each pixel in the mosaic image only has one color channel, the previous color edge detectors [2, 5, 6, 24, 25, 26] can not work well on mosaic images directly. The motivations of this research are two–fold: (1) proposing a new approach to capture more accurate color gradient/edge information on mosaic images directly and (2) developing a new high–quality edge–sensing demosaicing algorithm.

In this paper, without demosaicing processing, a new approach to capture more accurate color gradient/edge information on mosaic images directly is first proposed. Next, based on spectral–spatial correlation [27], a novel adaptive heterogeneity–projection for mosaic images is presented. Combining the captured color gradient/edge information and the adaptive heterogeneity–projection values, we present a new high–quality edge–sensing demosaicing algorithm. Based on twenty–four popular testing mosaic images, experimental results demonstrated that our proposed new demosaicing algorithm has the best quality performance when compared to current published four demosaicing algorithms by Pei and Tam [20], Lu and Tan [14], Dubois [4], and Tsai and Song [27].

The remainder of this paper is organized as follows. In Section II, a new approach to capture more accurate color gradient/edge information on mosaic images is presented.

In Section III, a novel adaptive heterogeneity–projection for mosaic images is first presented and then combining the captured color gradient/edge information and the adaptive heterogeneity–projection values, our proposed new edge–sensing demosaicing algorithm is presented. In Section IV, some experimental results are carried out to illustrate the quality advantage of our proposed demosaicing algorithm when compared to four previous demosaicing algorithms. Finally, some conclusions are addressed in Section V.

II. NEW APPROACH TO CAPTURE MORE ACCURATE COLOR GRADIENT INFORMATION ON MOSAIC IMAGES

This section presents a new approach to capture more accurate color gradient in-

(6)

formation on mosaic images directly. Although previously, several efficient demosaicing algorithms have been developed [13, 14, 20, 21, 28], we plug the bilinear interpolation technique [21] into our proposed approach on mosaic images due to its simplicity and computational efficiency. In what follows, we first survey the bilinear interpolation technique for demosaicing mosaic images. Next, our proposed new approach to capture more accurate color gradient information on mosaic images directly is presented.

2.1. Bilinear interpolation–based demosaicing process

Based on bilinear interpolation estimation, the other two unknown color values of each pixel on the mosaic image can be estimated by considering the color information of its adjacent pixels. For a demosaiced full color image, suppose the R, G, and B color values of the pixel at position (i, j) are denoted by I_dm^r (i, j), I_dm^g (i, j), and I_dm^b (i, j), respectively;

for a mosaic image, the R, G, and B color pixels located at position (i, j) are denoted by I_mo^r (i, j), I_mo^g (i, j), and I_mo^b (i, j), respectively, where i denotes the location on the vertical axis and j denotes the location on the horizontal axis (see Fig. 1). The full color image can be recovered from the mosaic image by the following rules:

Case 1: Compensate red and blue colors for the pixel at position (i, j), i, j ∈ odd, with green color by computing

I_dm^r (i, j) =1

2[I_mo^r (i, j − 1) + I_mo^r (i, j + 1)]

I_dm^g (i, j) =I_mo^g (i, j) I_dm^b (i, j) =1

2[I_mo^b (i − 1, j) + I_mo^b (i + 1, j)].

(1)

Case 2: Compensate red and blue colors for the pixel at position (i, j), i, j ∈ even, with green color by computing

I_dm^r (i, j) =1

2[I_mo^r (i − 1, j) + I_mo^r (i + 1, j)]

I_dm^g (i, j) =I_mo^g (i, j) I_dm^b (i, j) =1

2[I_mo^b (i, j − 1) + I_mo^b (i, j + 1)].

(2)

(7)

Case 3: Compensate green and blue colors for the pixel at position (i, j), i ∈ odd and j ∈ even, with red color by computing

I_dm^r (i, j) =I_mo^r (i, j) I_dm^g (i, j) =1

4

· I_mo^g (i − 1, j) + I_mo^g (i + 1, j) +I_mo^g (i, j − 1) + I_mo^g (i, j + 1)

¸

I_dm^b (i, j) =1 4

· I_mo^b (i − 1, j − 1) + I_mo^b (i − 1, j + 1) +I_mo^b (i + 1, j − 1) + I_mo^b (i + 1, j + 1)

¸ .

(3)

Case 4: Compensate red and green colors for the pixel at position (i, j), i ∈ even and j ∈ odd, with blue color by computing

I_dm^r (i, j) =1 4

· I_mo^r (i − 1, j − 1) + I_mo^r (i − 1, j + 1) +I_mo^r (i + 1, j − 1) + I_mo^r (i + 1, j + 1)

¸

I_dm^g (i, j) =1 4

· I_mo^g (i − 1, j) + I_mo^g (i + 1, j) +I_mo^g (i, j − 1) + I_mo^g (i, j + 1)

¸

I_dm^b (i, j) =I_mo^b (i, j).

(4)

The data dependence of compensation for the above four cases is depicted in Fig. 4.

(a) (b) (c) (d)

Fig. 4: The data dependence in the bilinear interpolation estimation. (a) Compensating red and blue colors for case 1. (b) Compensating red and blue colors for case 2. (c) Compensating blue and green colors for case 3. (d) Compensating red and green colors for case 4.

2.2. Capturing the full color gradient information on mosaic images directly In this subsection, our proposed new approach to capture more accurate color gradient information on mosaic images directly is presented. The captured gradient information will be used in our proposed new edge–sensing demosaicing algorithm in Section III and it would lead to a high–quality advantage.

(8)

Before presenting the proposed new approach, for completeness, we first introduce how to run Sobel operator [7] on full color images. In [30], Sobel operator for gray images is extended for capturing the gradient information from the gray image domain to the full color image domain. The four masks used in 3 × 3 Sobel operator for full color image are illustrated in Fig. 5. Given a demosaiced full color image, it has been known that the

(a) (b) (c) (d)

Fig. 5: The 3 × 3 Sobel operator. (a) The horizontal mask. (b) The vertical mask. (c) The ^π₄–diagonal mask. (d) The ^−π₄ –diagonal mask.

R, G, and B color values of the pixel at position (i, j) are denoted by I_dm^r (i, j), I_dm^g (i, j), and I_dm^b (i, j), respectively. After running Sobel horizontal, vertical, ^π₄–diagonal, and ^−π₄ – diagonal masks as shown in Fig. 5 on the 3 × 3 full color subimage centered at position (i, j), the horizontal response ∆I_dm^H (i, j), the vertical response ∆I_dm^V (i, j), the ^π₄–diagonal response ∆I_dm^π⁴ (i, j), and the ^−π₄ –diagonal response ∆I_dm^−π⁴ (i, j) are calculated by

∆I_dm^H (i, j) = 0.299∆I_dm^H,r(i, j) + 0.587∆I_dm^H,g(i, j) + 0.114∆I_dm^H,b(i, j)

∆I_dm^V (i, j) = 0.299∆I_dm^V,r(i, j) + 0.587∆I_dm^V,g(i, j) + 0.114∆I_dm^V,b(i, j)

∆I_dm^π⁴ (i, j) = 0.299∆I_dm^π⁴^,r(i, j) + 0.587∆I_dm^π⁴^,g(i, j) + 0.114∆I_dm^π⁴^,b(i, j)

∆I_dm^−π⁴ (i, j) = 0.299∆I_dm^−π⁴ ^,r(i, j) + 0.587∆I_dm^−π⁴ ^,g(i, j) + 0.114∆I_dm^−π⁴ ^,b(i, j)

(5)

(9)

where for C ∈ {r, g, b},

∆I_dm^H,C(i, j) =

½ [I_dm^C (i − 1, j + 1) + 2I_dm^C (i, j + 1) + I_dm^C (i + 1, j + 1)]

−[I_dm^C (i − 1, j − 1) + 2I_dm^C (i, j − 1) + I_dm^C (i + 1, j − 1)]

¾

∆I_dm^V,C(i, j) =

½ [I_dm^C (i + 1, j − 1) + 2I_dm^C (i + 1, j) + I_dm^C (i + 1, j + 1)]

−[I_dm^C (i − 1, j − 1) + 2I_dm^C (i − 1, j) + I_dm^C (i − 1, j + 1)]

¾

∆I_dm^π⁴^,C(i, j) =

½ [I_dm^C (i − 1, j) + 2I_dm^C (i − 1, j + 1) + I_dm^C (i, j + 1)]

−[I_dm^C (i + 1, j) + 2I_dm^C (i + 1, j − 1) + I_dm^C (i, j − 1)]

¾

∆I_dm^−π⁴ ^,C(i, j) =

½ [I_dm^C (i − 1, j) + 2I_dm^C (i − 1, j − 1) + I_dm^C (i, j − 1)]

−[I_dm^C (i + 1, j) + 2I_dm^C (i + 1, j + 1) + I_dm^C (i, j + 1)]

¾

(6)

In order to make Sobel operator workable on mosaic images to extract more accurate gradient information, the bilinear interpolation technique could be embedded into the Sobel operator. Combining Eq. (1)–Eq. (4) and Eq. (6), our proposed new approach to capture more accurate color gradient information on mosaic images directly can be followed. Considering all different cases, we first derive two quad–masks for G channel, and then four quad–masks for R and B channels are derived (Detailed derivations are shown in Appendix I and Appendix II). According to Appendix I, the Sobel– and interpolation–

based (SI–bsaed) two quad–masks for i + j ∈ even and i + j ∈ odd are shown in Fig.

6 and Fig. 7, respectively. For G channel, the four masks can be used to capture the horizontal response, the vertical response, the ^π₄–diagonal response, and the ^−π₄ –diagonal response without demosaicing processing. Because the derivation of four quad–masks for R channel is the same as the B channel, we thus only derive the masks for R channel.

According to the result of Appendix II, the derived four SI–based quad–masks for R/B channel are illustrated in Fig. 8–Fig. 11.

After running the above proper SI–based masks on the 5×5 mosaic subimage centered at position (i, j), the obtained horizontal responses of R, G, and B channels are denoted by ∆I_dm^H,r(i, j), ∆I_dm^H,g(i, j), and ∆I_dm^H,b(i, j), respectively; the vertical responses of R, G, and B channels are denoted by ∆I_dm^V,r(i, j), ∆I_dm^V,g(i, j), and ∆I_dm^V,b(i, j), respectively; the

π

4–diagonal responses of R, G, and B channels are denoted by ∆I_dm^π⁴^,r(i, j), ∆I_dm^π⁴^,g(i, j), and ∆I_dm^π⁴^,b(i, j), respectively, and the ^−π₄ –diagonal responses of R, G, and B channels

(10)

(a) (b) (c) (d)

Fig. 6: For i + j ∈ even, the four SI–based masks for G channel. (a) The horizontal SI–based mask. (b) The vertical SI–based mask. (c) The ^π₄–diagonal SI–based mask. (d) The ^−π₄ –diagonal SI–based mask.

(a) (b) (c) (d)

Fig. 7: For i + j ∈ odd, the four SI–based masks for G channel. (a) The horizontal SI–based mask. (b) The vertical SI–based mask. (c) The ^π₄–diagonal SI–based mask. (d) The ^−π₄ –diagonal SI–based mask.

are denoted by ∆I_dm^−π⁴ ^,r(i, j), ∆I_dm^−π⁴ ^,g(i, j), and ∆I_dm^−π⁴ ^,b(i, j), respectively. By Eq. (5), combining the three horizontal responses, we have the total horizontal response ∆I_dm^H (i, j).

Similarly, we have the total vertical response ∆I_dm^V (i, j), the total ^π₄–diagonal response

∆I_dm^π⁴ (i, j), and the total ^−π₄ –diagonal response ∆I_dm^−π⁴ (i, j).

III. THE PROPOSED NEW EDGE–SENSING DEMOSAICING ALGORITHM Based on the accurate gradient information obtained in last section, this section presents our proposed new high–quality edge–sensing demosaicing algorithm. Our proposed novel adaptive heterogeneity–projection for mosaic images is first presented in subsection 3.1 to capture more accurate horizontal and vertical heterogeneity–projection value. Based on the two obtained heterogeneity–projection values and those obtained gra-

(11)

(a) (b) (c) (d)

Fig. 8: For i ∈ odd, j ∈ even (i ∈ even, j ∈ odd), the four SI–based masks for R (B) channel. (a) The horizontal SI–based mask. (b) The vertical SI–based mask. (c) The

π

4–diagonal SI–based mask. (d) The ^−π₄ –diagonal SI–based mask.

(a) (b) (c) (d)

Fig. 9: For i ∈ even, j ∈ odd (i ∈ odd, j ∈ even), the four SI–based masks for R (B) channel. (a) The horizontal SI–based mask. (b) The vertical SI–based mask. (c) The

π

4–diagonal SI–based mask. (d) The ^−π₄ –diagonal SI–based mask.

dient information mentioned in Section II, the proposed new edge–sensing interpolation estimation is presented. The proposed interpolation estimation consists of two steps, the first step for G channel and the second step for R and B channels.

3.1. Novel adaptive heterogeneity–projection for mosaic images

In this subsection, a novel adaptive heterogeneity–projection for mosaic images is presented. Given an original mosaic image I_mo, its horizontal heterogeneity–projection map HP_H−map and vertical heterogeneity–projection map HP_{V −map} can be obtained by running the following two 1–D Laplacian operations [27]:

HP_H−map = |I_mo∗ (F^{N ×(N −3)}T^{(N −3)×1})^t| HP_{V −map} = |I_mo∗ (F^{N ×(N −3)}T^{(N −3)×1})|

(7)

(12)

(a) (b) (c) (d)

Fig. 10: For i, j ∈ odd (i, j ∈ even), the four SI–based masks for R (B) channel. (a) The horizontal SI–based mask. (b) The vertical SI–based mask. (c) The ^π₄–diagonal SI–based mask. (d) The ^−π₄ –diagonal SI–based mask.

(a) (b) (c) (d)

Fig. 11: For i, j ∈ even (i, j ∈ odd), the four SI–based masks for R (B) channel. (a) The horizontal SI–based mask. (b) The vertical SI–based mask. (c) The ^π₄–diagonal SI–based mask. (d) The ^−π₄ –diagonal SI–based mask.

where N (=11 empirically [27]) denotes the vector length (or the 1–D mask size); F^{N ×(N −3)}= [1 − 1 − 1 1]^t ∗ I(N −3)×(N −3) denotes a N × (N − 3) coefficient matrix; T^{(N −3)×1} = Q_{N −4}

x=1[1 − 1]^t∗ I(N −3−x)×(N −3−x) denotes a (N − 3) × 1 coefficient vector; I^{M ×M} denotes an M × M identity matrix; the symbol “∗” denotes the 2–D convolution operator; | · | denotes the absolute value operator and the operator “t” denotes the transpose operator.

In what follows, our approach can determine the suitable value of N adaptively.

Let us examine some cases of the current pixel on the mosaic image. If the surrounding region of the current pixel is homogenous, the two responses by Eq. (7) are almost the same whether a large mask size or small mask size is adopted. If there is one tiny hori- zontal (vertical) edge passing through the current pixel, a small mask size for HP_{V −map}

(13)

(HP_H−map) is enough rather than a large mask size. In [27], the adopted mask size is fixed and is set to N = 11. According to the above discussion, we now present our pro- posed adaptive heterogeneity–projection for mosaic images such that the used mask size is as small as possible and the computed responses are more accurate than those computed by the mask with size 11. Experimental results reveal that our proposed adaptive heterogeneity–projection approach under proper horizontal and vertical mask sizes has the computation–saving and more accurate advantages.

We utilize the horizontal (vertical) spectral–spatial correlation (SSC) [27] between the current pixel at location (i, j) and its neighboring pixel at location (i, j + 1) ((i + 1, j)) to determine the proper horizontal (vertical) mask size N_H(i, j) (N_V(i, j)). For a horizontal SSC map, the horizontal SSC value at location (i, j) can be obtained by using the following rules:

S_H(i, j) =









|I_mo^g (i, j) − I_mo^r (i, j + 1)| if the current pixel is green and i ∈ odd.

|I_mo^g (i, j) − I_mo^b (i, j + 1)| if the current pixel is green and i ∈ even.

|I_mo^r (i, j) − I_mo^g (i, j + 1)| if the current pixel is red.

|I_mo^b (i, j) − I_mo^g (i, j + 1)| if the current pixel is blue.

Similarly, the vertical SSC value at location (i, j) can be obtained by using the following rules:

SV(i, j) =









|I_mo^g (i, j) − I_mo^r (i + 1, j)| if the current pixel is green and j ∈ even.

|I_mo^g (i, j) − I_mo^b (i + 1, j)| if the current pixel is green and j ∈ odd.

|I_mo^r (i, j) − I_mo^g (i + 1, j)| if the current pixel is red.

|I_mo^b (i, j) − I_mo^g (i + 1, j)| if the current pixel is blue.

The horizontal SSC map and the vertical SSC map of the mosaic Lighthouse image in Fig. 3 are illustrated in Fig. 12(a) and Fig. 12(b), respectively, and it is observed that the SSC values are locally constant in homogeneous regions.

Because the determination of the proper horizontal mask size N_H(i, j) is the same as that for the vertical mask size N_V(i, j), we thus only focus on the determination of horizontal mask size. Fig. 13 depicts our proposed approach to determine the proper horizontal mask size NH(i, j). Centered at location (i, j), Fig. 13 illustrates a 1 × 11 row data extracted from the horizontal SSC map. The procedure to determine the proper horizontal mask size N_H(i, j) consists of the following three steps.

(14)

(a) (b)

Fig. 12: The two SSC maps of the mosaic Lighthouse image. (a) The horizontal SSC map. (b) The vertical SSC map. (In order to show the images more clear, the gray value 80 is used to represent the value 0)

Fig. 13: The depiction of our proposed approach to determine the proper horizontal mask size N_H(i, j).

Step 1: Initially, we compute the two SSC values S_H(i, j − 2) and S_H(i, j + 2), and temporarily set N_H(i, j) = 5, x_{lef t} = j − 2, and x_right = j + 2.

Step 2: If the condition N_H(i, j) = 11 holds, output N_H(i, j) and stop. Here, the maxi- mal horizontal mask size is bounded by eleven. Otherwise, go to Step 3.

Step 3: Examine whether the neighboring horizontal SSC values are locally chosen by using the following testing condition:

Max(DS_H^{lef t}, DS_H^right) < T_h

where DS_H^{lef t} = |S_H(i, x_{lef t}) − S_H(i, x_{lef t}− 1)| + |S_H(i, x_{lef t}) − S_H(i, x_{lef t}+ 1)|, DS_H^right = |S_H(i, x_right) − S_H(i, x_right− 1)| + |S_H(i, x_right) − S_H(i, x_right + 1)|, and the threshold is set to Th = 8 empirically. If the above testing condition

(15)

holds, output the value of N_H(i, j) as the proper horizontal mask size and stop.

Otherwise, perform the operation N_H(i, j) = N_H(i, j) + 2, x_{lef t}= x_{lef t}− 1, and x_right = x_right+ 1, and then go to Step 2.

The above three–step procedure can be easily applied to determine the proper horizontal mask size. In order to normalize the masks for different sizes, the normal factor _Q¹

(N ) is used to normalize the coefficients of the mask. In other words, instead of using F^{N ×(N −3)}T^{(N −3)×1}, we use _Q¹

(N )(F^{N ×(N −3)}T^{(N −3)×1}) to obtain the normalized heterogeneity–projection map. The value of Q_{(N )}is defined as the sum of the positive coef- ficients in the mask. For example, if N = 5, the mask F^5×(5−3)T^{(N −3)×1}= [1 −2 0 2 −1]^t, can be normalized to ^F^5×2_Q^T^2×1

(5) = [1 −2 0 2 −1]^t

3 .

Finally, in order to reduce the estimation error, we use the local mean to tune the heterogeneity–projection maps. For the HP_H−map (HP_{V −map}), the horizontal (vertical) heterogeneity–projection value at location (i, j) is denoted by HP_H(i, j) (HP_V(i, j)).

The tuned horizontal heterogeneity–projection value HP_H⁰ (i, j) and the tuned vertical heterogeneity–projection value HP_V⁰(i, j) can be computed by the following operations:

HP_H⁰ (i, j) = 1 6

X2 k=−2

δkHPH(i, j + k)

HP_V⁰ (i, j) = 1 6

X2 k=−2

δ_kHP_V(i + k, j)

where δ_k= 2 if k = 0; δ_k = 1, otherwise.

After performing the above adaptive heterogeneity–projection for mosaic images, the values HP_H⁰ (i, j) and HP_H⁰ (i, j) can be obtained for each mosaic pixel. In next two subsections, the six responses HP_H⁰ (i, j), HP_H⁰ (i, j), ∆I_dm^H (i, j), ∆I_dm^V (i, j), ∆I_dm^π⁴ (i, j), and ∆I_dm^−π⁴ (i, j) will be used in our proposed new edge–sensing interpolation estimation for demosaicing mosaic images.

3.2. The edge–sensing interpolation estimation for G channel

(16)

Fig. 14: The 7 × 7 mosaic subimage.

In this subsection, we present the proposed high–quality edge–sensing interpolation estimation for G channel. For exposition, let us take a 7 × 7 mosaic subimage as shown in Fig. 14 to explain how to estimate the G channel value located at the center position of Fig. 14. Before performing the interpolation estimation for G channel, assume the gradient information of the current mosaic pixel at position (i, j) and the eight neighboring pixels at positions Ω_n = {(x, y)|(x, y) = (i ± 1, j), (i ± 2, j), (i, j ± 1), (i, j ± 2)} have been captured by using our proposed method described in last section, and the captured nine horizontal gradient magnitudes and nine vertical gradient magnitudes are denoted by

∆I_dm^H (∗, ∗) and ∆I_dm^V (∗, ∗), respectively.

According to the tuned horizontal heterogeneity–projection value HP_H⁰ (i, j) and the tuned vertical heterogeneity–projection value HP_V⁰ (i, j) of the current mosaic pixel at position (i, j), the interpolation estimation scheme in our proposed edge–sensing demo- saicing algorithm considers three cases, namely (1) horizontal variation as shown in Fig.

15(a), (2) vertical variation as shown in Fig. 15(b), and (3) the other variations as shown in Fig. 15(c). The arrows in Fig. 15 denote the data dependence.

In addition, in order to estimate I_dm^g (i, j) more accurately from its four neighbors,

(17)

four proper weights in terms of gradient information are assigned to corresponding four spectral–correlation terms in the interpolation estimation. Considering the neighboring pixel located at location (i−1, j), if there is (is no) a horizontal edge passing through it, i.e.

the vertical gradient magnitudes of the pixel at location (i−1, j) is large (small), the color difference assumption [14, 20] reveals that the green component of this pixel makes less (more) contribution to the estimation of green component for the current pixel at location (i, j). On the other hand, if the gradient magnitudes of the pixels at location (i − 2, j) and (i, j) are large (small), the pixel at location (i − 1, j) is located in a nonhomogeneous (homogeneous) region and we claim that the pixel at location (i − 1, j) makes less (more) green contribution to the estimation of green component for the current mosaic pixel at lo- cation (i, j). Combining the above analysis of gradient and direction effects, the weight of pixel at location (i − 1, j) can be given by w_g(i − 1, j) = _1+β[^∆^Vdm(i−2,j)+2∆¹^V_dm(i−1,j)+∆^V_dm(i,j)]. Following the similar discussion, the weights of the four neighbors of the current pixel are expressed by wg(i − 1, j) = _1+β[^P²k=0δ_k¹⁰∆^V_dm(i−k,j)], wg(i + 1, j) = _1+β[^P²k=0δ⁰_k¹∆^V_dm(i+k,j)], w_g(i, j − 1) = ¹

1+β[^P²k=0δ⁰_k∆^H_dm(i,j−k)], and w_g(i, j + 1) = ¹

1+β[^P²k=0δ_k⁰∆^H_dm(i,j+k)] where δ_k⁰ = 2 if k = 1; δ⁰_k= 1, otherwise; the parameter β is set to β = 0.7 empirically.

According to the above description, the value of I_dm^g (i, j) of the current pixel at location (i, j) can be estimated by the following rules:

I_dm^g (i, j) = I_mo^b (i, j) + P

(x,y)∈Ωgw_g(x, y)D_g(x, y) P

(x,y)∈Ωgw_g(x, y) Ωg =





{(i ± 1, j)} if HP_V⁰(i, j) < αHP_H⁰ (i, j) (horizontal variation).

{(i, j ± 1)} if HP_H⁰ (i, j) < αHP_V⁰(i, j) (vertical variation).

{(i ± 1, j), (i, j ± 1)} otherwise (other variations).

where for (x₁, y₁) ∈ {(i ± 1, j)}, D_g(x₁, y₁) = I_mo^g (x₁, y₁) − ^I^mo^b ^(x¹^+1,y¹^)+I₂ ^mo^b ^(x¹^−1,y¹⁾; for (x₂, y₂) ∈ {(i, j ± 1)}, D_g(x₂, y₂) = I_mo^g (x₂, y₂) −^I^mo^b ^(x²^,y²^+1)+I₂ ^mo^b ^(x²^,y²⁻¹⁾; the parameter α is set to α = 0.5.

After performing the above edge–sensing interpolation estimation for G channel, the G channel of the demosaiced image is fully populated. In next subsection, the fully populated G channel of the image will be used to assist the interpolation of R and B

(18)

channels.

(a) (b) (c)

Fig. 15: The data dependence of our proposed interpolation estimation for G channel.

(a) Horizontal variation (vertical edge). (b) Vertical variation (horizontal edge). (c) The other variations.

3.3. The edge–sensing interpolation estimation for R and B channels

Because the number of R pixels or B pixels is less than G in the mosaic image, the interpolation estimation for R and B channels should be partitioned into two steps: (1) estimating the red (blue) values at blue (red) pixels and then (2) recovering the remains of the red and blue values at green pixels. Because the interpolation estimation for R channel is the same as it for B channel, we thus only present it for R channel. In our proposed interpolation estimation for R and B channels, the fully populated G channel is used to assist the interpolation of R and B channels. For convenience, we still use Fig.

14 to explain how to estimate the R channel value for the current pixel at position (i, j).

Similar to the interpolation estimation for G channel, assume the gradient infor- mation of the current mosaic pixel at position (i, j) and the eight pixels at positions Ω⁰_n = {(x, y)|(x, y) = (i ± 1, j ± 1), (i ± 2, j ± 2)}, respectively, have been captured by using our proposed method. By the same arguments, in order to estimate I_dm^r (i, j) more accurately from its four red neighbors, four proper weights in terms of gradient/direction information are assigned to the corresponding four spectral–correlation terms in the in- terpolation estimation. For estimating I_dm^r (i, j), we consider four diagonal variations of

(19)

the mosaic pixel at position (i, j) to determine the four weights which will be used in the proposed edge–sensing estimation of I_dm^r (i, j). Considering the neighboring pixel at location (i − 1, j − 1), if there is (is no) a ^π₄–diagonal edge passing through it, i.e. the

−π

4 –diagonal gradient magnitude of the pixel at location (i − 1, j − 1) is large (small), and then it claims that the red component of this pixel makes less (more) contribution to the estimation of red component for the current pixel at location (i, j). On the contrary, if the gradient magnitudes of the pixel at location (i − 1, j − 1) and (i − 2, j − 2) are large (small), i.e. the pixel at location (i − 1, j − 1) is inside nonhomogeneous (homogeneous) region, the color difference assumption indicates that the green component of this pixel makes less (more) contribution to the estimation of red component for the current pixel at location (i, j). Following the above analysis, the weights of the four diagonal red neigh- bors of the current pixel can be expressed by w_r(i − 1, j − 1) = ¹

1+β

P2 k=0δ_k⁰∆I

−π4

dm(i−k,j−k)

, w_r(i − 1, j + 1) = ¹

1+β

P2

k=0δ_k⁰∆I_dm^π⁴ (i−k,j+k)

, w_r(i + 1, j − 1) = ¹

1+β

P2

k=0δ_k⁰∆I_dm^π⁴ (i+k,j−k)

, and wr(i + 1, j + 1) = ¹

1+β

P2 k=0δ_k⁰∆I

−π4

dm(i+k,j+k)

where δ⁰_k = 2 if k = 1; δ_k⁰ = 1, other- wise; the parameter set to β = 0.5 empirically. Based on the four weights obtained above and the color difference concept, the demosaiced full red color for the blue pixel in the mosaic image, r_f^d(i, j), can be estimated by

I_dm^r (i, j) = I_dm^g (i, j) + P

(x,y)∈Ωrw_r(x, y)D_r(x, y) P

(x,y)∈Ωrw_r(x, y) , Ωr = {(x, y)|(i ± 1, j ± 1)}

where for (x, y) ∈ {(i ± 1, j ± 1)}, Dr(x, y) = I_mo^r (x, y) − I_dm^g (x, y).

After describing how to estimate demosaiced full red colors for those blue pixels in the mosaic image, we now introduce how to estimate full red colors for those green pixels.

Fig. 16 illustrates the pattern of R channel. Referring to Fig. 16, the full red color for the green pixel, I_dm^r (i, j), can be estimated by the following rules:

I_dm^r (i, j) = I_dm^g (i, j) + P

(x,y)∈Ω⁰rwr(x, y)Dr(x, y) P

(x,y)∈Ω⁰_rw_r(x, y) Ω⁰_r =





{(i ± 1, j)} if HP_V⁰(i, j) < αHP_H⁰ (i, j) (horizontal variation).

{(i, j ± 1)} if HP_H⁰ (i, j) < αHP_V⁰(i, j) (vertical variation).

{(i ± 1, j), (i, j ± 1)} otherwise (other variations).

(20)

Fig. 16: The patterns of the R channel.

where for (x, y) ∈ {(i ± 1, j), (i, j ± 1)}, we perform D_r(x,y) = I_mo^r (x, y) − I_dm^g (x, y). If k = 1, we set δ_k⁰ = 2; otherwise, we set δ_k⁰ = 1. We further perform wr(i − 1, j) =

1

1+β[^P²k=0δ⁰_k∆I_dm^V (i−k,j)], w_r(i + 1, j) = ¹

1+β[^P²k=0δ⁰_k∆I_dm^V (i+k,j)], w_r(i, j − 1) = ¹

1+β[^P²k=0δ_k⁰∆I_dm^H (i,j−k)], and w_r(i, j + 1) = ¹

1+β[^P²k=0δ⁰_k∆I_dm^H (i,j+k)]. Empirically, the parameter β is set to β = 0.7.

After performing the interpolation estimation, the postprocessing proposed in [15]

could be used to improve the demosaiced image quality. After presenting our proposed new edge–sensing demosaicing algorithm, experimental results in next section will illustrate the quality advantage of our proposed new edge–sensing demosaicing algorithm.

IV. EXPERIMENTAL RESULTS

In this section, based on twenty–four popular testing mosaic images, some experimental results are demonstrated to show that our proposed new demosaicing algorithm has better image quality performance when compared to the previous four algorithms by Pei and Tam [20], Lu and Tan [14], Dubois [4], and Tsai and Song [27], respectively. The concerned algorithms are implemented on the IBM compatible computer with Intel Core 2 Duo CPU 1.6GHz and 1GB RAM. The operating system used is MS–Windows XP and the program developing environment is Borland C++ Builder 6.0. Our program has been uploaded in [32].

Fig. 17 illustrates the twenty–four testing images from Kodak PhotoCD. In our ex-