Color Filter Arrays for Quanta Image Sensors
Abstract
Quanta image sensor (QIS) is to be the next generation image sensor after CCD and CMOS. To enable such technology, significant progress was made over the past five years to advance both the device and image reconstruction algorithms. In this paper, we discuss color imaging using QIS, in particular how to design color filter arrays. Designing color filter arrays for QIS is challenging because at the pixel pitch of 1.1m, maximizing the light efficiency while suppressing aliasing and crosstalk are conflicting tasks. We present an optimizationbased framework to design color filter arrays for very small pixels. The new framework unifies several mainstream color filter array design frameworks by offering generality and flexibility. Compared to the existing frameworks which can only handle one or two design criteria, the new framework can simultaneously handle luminance gain, chrominance gain, crosstalk, antialiasing, manufacturability and orthogonality. Extensive experimental comparisons demonstrate the effectiveness and generality of the framework.
28
I Introduction
Ia Quanta Image Sensor
Quanta Image Sensor (QIS) is a new type of image sensors proposed by E. Fossum in 2005 [1, 2, 3] as a candidate for the third generation digital image sensor after CCD and CMOS. The sensor comprises a massive array of subdiffraction limit singlephoton detectors, called “jots”, with pixel pitch of m. Having readout noise of e and dark current of 0.16e/sec, QIS can count every incoming photon to produce digital output of bit depth in the range of bits without saturation. As reported in [2] and [3], the latest QIS sensor can achieve a resolution of jots and frame rate up to fps by using the commercial 45/65nm 3Dstacked backside illumination CMOS process. Figure 1 illustrates how QIS acquires images by first detecting a stream of binary bitmaps, and then reconstruct the image using algorithms.
Despite the rapid advancement in QIS hardware [1, 2, 3] and algorithms [4, 5, 6, 7] of QIS, all reported findings, todate, are based on monochromatic data. The first color QIS imaging is only recently proposed by Gnanasambandam et al. [8], where they demonstrated a new sensor and new reconstruction algorithms. In this paper, we discuss how to design the color filter array.
IB Color Filter Arrays for QIS
A color filter array (CFA) is a mask placed on top of the sensor to select (filter) wavelengths. As light passes through the color filter array, the resulting image is a mosaic pattern of the three tristimulus RGB colors. Traditionally, CFA is organized as a periodic replica of a 2D kernel called the color atom. The defacto color atom used in the industry is the Bayer pattern [9] because of its simplicity. More sophisticated CFAs have been proposed [10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20] to improve the Bayer CFA.
How to design a good CFA? There are fundamentally three factors that should be taken into consideration:
Aliasing: Since color filtering is a sampling process, aliasing happens when the sampling rate is less than Nyquist. Aliasing causes false color artifacts at color edges, called the Moirè artifacts [11]. Color filters that are susceptible to aliasing, such as the Bayer CFA, require sophisticated demosaicking algorithms to suppress the Moirè artifacts. In contrast, a robust CFA can use simple demosaicking algorithms.
Sensitivity: Since CFA is a filter, it blocks part of the incoming light. This reduces the sensor sensitivity and makes the image more susceptible to noise. A good CFA design should maximize the sensitivity by allowing transparent or “panchromatic” color filters that block as few wavelengths as possible.
Crosstalk: Crosstalk can be either optical or electrical [21]. If not treated, crosstalk will make colors look pale or desaturated. Crosstalk desaturation is corrected by pixelwise multiplication of the RGB color vector using a color correction matrix. However, color correction enhances residual noise in the image [22, 21]. The situation is worsen in QIS because of its small size.
The three factors above are conflicting: Optimizing one factor generally degrades the others. For conventional CMOS image sensors, crosstalk is not severe. Therefore, most CFA designs in the literature consider the first two factors only which are inadequate for QIS. The only available work on QIS color filter array design is by Anzagira and Fossum [21]. However, the aliasing factor was not adequately handled and so their CFA is susceptible to the Moirè artifacts. Other mainstream CFA designs for CMOS lack a unified framework that encompasses all these three design factors in a rigorous optimization framework.
IC Related Work
The design framework we propose in this paper is a unification of several mainstream CMOSbased color filter arrays. To put our paper in the proper context in the literature, we here list a few of the better known results.
Spatial CFA Design: By suppressing the Moirè artifacts and crosstalk while keeping the demosaicing algorithm simple, Lukac and Plataniotis [10] proposed a CFA and compared it with other CFAs using a universal demosaicking method. However, their work did not provide a mathematical framework to analyze the CFA optimality.
SpatioSpectral CFA Design: Hirakawa and Wolfe [11] proposed a method through the spatial and spectral domain analysis. Their CFA reduces aliasing in the frequency domain, and possesses high sensitivity and numerical stability. Condat [23] extended the framework by optimizing luminance and chrominance sensitivity. He defined a new form of orthogonality between chrominance channels in frequency domain. Hao et al. [14] and Wang et al. [15] proposed a framework based on symbolic discrete Fourier transform (DFT). Their CFA maximizes the numerical stability of linear demosaicking process under aliasing and physical constraints.
Learningbased CFA Design: By minimizing the average error on a color dataset, Lu and Vetterli [12] used an iterative algorithm to solve a least squares CFA design problem. Chakrabarti [24] and Henz et al. [25] proposed to learn the optimal CFA pattern by using a deep neural network.
Besides these mainstream CFA design frameworks, there are a number of other CFA designs such as [13, 17, 16, 19, 20, 18]. On the hardware side, [26] and [27] took into account that color filter fabrication technology lags the image sensor technology in terms of miniaturization. They proposed a hardwarefriendly CFA assuming the color filter size is pixel size.
ID Scope and Contributions
In this paper, we propose a rigorous optimization framework that encompasses aliasing, sensitivity and crosstalk in a unified model. This is the first work that incorporates a quantitative crosstalk metric in an optimization framework for CFA design. Our framework is generic as it works with any crosstalk kernel depending on the sensor. Moreover, it works with any color atom size by a mild change in the formulation.
The main contribution of this paper is a general and flexible framework for CFA design. Compared to the existing CFA design framework, the new framework is able to simultaneously (Section III):

Improve CFA’s luminance and chrominance sensitivity,

Reduce aliasing between luminance and chrominance channels,

Suppress crosstalk between spectral subbands, and

Enforce orthogonality between chrominance channels to permit simple linear demosaicking.
The design framework is presented in the form of optimization. We have two designs: A convex optimization and a nonconvex optimization. In addition to the formulation, we also present an algorithm to solve the nonconvex optimization. (Section IV)
For performance evaluation of different CFAs on QIS images, we propose in Section V a universal demosaicking pipeline. This pipeline comprises iterative ADMM algorithm for removing the CFA masking effect followed by a color correction step for removing the desaturation effect of crosstalk. Experimental evaluation on the Kodak and McMaster color datasets in Section VI shows the robustness of our proposed CFAs compared to other CFAs in literature.
Ii Background and Notations
The design of a robust CFA involves multiple objectives in terms of signal transmittance, color aliasing, crosstalk, and manufacturability. To facilitate readers to understand the design framework, in this section we introduce a few notations and terminologies. We will start in Section IIA by describing the image formation using a CFA, then we discuss CFA in different domains in Sections IIB and IIC. Afterwards, in Section IID, we define the optimization variables to simplify the design framework.
Iia Color Image Formation
Consider a color image of size . We denote the normalized light intensities in the red, green and blue channels for the (th,th) pixel of the color image as
(1) 
where , .
Color Filtering: To obtain color, we place a color filter on top of each jot to collect light for one of the RGB colors. The CFA is a periodic pattern of the same resolution of , defined as
(2) 
where , , are the opacity rates for the red, green and blue pixels, respectively. For example, a red color filter is defined as as it only passes the red color. The light exposure on the QIS after passing through the CFA is denoted as , which is a linear combination of the tristimulus colors:
(3) 
Here, is a positive scalar that represents the sensor gain factor. To write this equation in matrix form, we stack the vectorized color channels of the latent image in long vector, and the mosaic channel in long vector as follows:
(4) 
We define the sampling matrix , where . Then, the color filter sampling can be written as:
(5) 
Photon Arrival. The photon arrival at QIS is modeled as a Poisson process. Let be a vector of nonnegative random integers denoting the number of photons arriving at QIS jots according to the light exposure . Then, the probability of observing a photon count is
(6) 
In this work, we assume singlebit QIS [28] that quantizes the photon counts by QIS jots to a binary values with if and if , where is a threshold. The probability of observing is
(7) 
where is the incomplete Gamma function [7].
Temporal Oversampling. With frame rates that reach fps, QIS is able to catch the scene movement by taking temporal samplings for the same scene. This allows us to utilize multiple independent measurements over time to improve the statistics and decrease noise. Hence, for every jot with light exposure , we have a set of independent binary measurements .
IiB Color Filter Array Analysis in Different Color Spaces
Since the CFA is a periodic pattern, it is sufficient to use a color atom as the optimization variable when designing the CFA. The color atom takes the form
(8) 
where each of , and is an array. For instance, the GRBG Bayer pattern has the following color atom (when ):
because the Bayer pattern has one red pixel and one blue pixel located at two opposite diagonals, and two green pixels located in the remaining two positions. Figure 2 illustrates the idea.
While the primal RGB color is common for making the CFA, it would be more convenient if the colors are decorrelated. To this end, we change the image representation from the canonical RGB basis to an orthornormal basis using a transformation matrix [29, 23]:
(9) 
This transformation maps an RGB image to an image as follows (we drop the spatial indices for simplicity)
where is a luminance (luma) component that contains the high frequency components such as edges and textures, whereas and are chrominance (chroma) components that carry the color information.
Since is orthonormal (i.e., ), we can rewrite the sampling process in (3) in the luma/chroma space:
(10) 
where , and are the luma/chroma representation of the CFA, with
(11) 
The luma/chroma representation of the CFA has a corresponding color atom , and . For instance, the luma/chroma color atom of the GRBG Bayer pattern is
(12) 
where the individual components are
Remark 1.
In principle, there are are infinite choices for the luma/chroma basis . We choose the one in (9) because it makes the components of natural images statistically independent in the first order approximation.
IiC Color Filter Array in Fourier Space
When analyzing the aliasing effects of the CFAs, we need to transform the color atom into the Fourier domain. For simplicity, we represent the Fourier transform of a signal by putting a tilde on top of the symbol, e.g., . The 2D discrete Fourier transform (DFT) of the th color atom is
(13) 
where , .
For example, the discrete Fourier transform of the luma/chroma color atoms in (12) are
Here, we observe that the Fourier transform of the color atom has the same size as the original color atom. The luminance channel has only one baseband components at , whereas the chrominance channel has one baseband component and a component at . The chrominance channel has two components at and . Figure 3 illustrates how these frequency locations are identified from a color atom.
While the Fourier transform of the color atom is useful, for demosaicing we also need to analyze the spectrum of the entire CFA. As shown in by Hao et al. [14], the Fourier transform of the entire CFA can be written in terms of the Fourier transform of the color atoms:
(14) 
where , is the 2D angular frequency, and
(15) 
is the (th,th) 2D angular frequency. It is worth noting that the Fourier transform of the CFA comprises pure sinusoids of amplitudes . These sinusoids are placed at discrete 2D frequencies that divide the 2D frequency plane into equal regions. Therefore, the spectrum of the mosaicked image can be written as
(16) 
where is the standard 2D convolution operator. Having the spectrum of the mosaicked image , we can now discuss the optimization variables in our problem.
IiD Design Variables
To formulate the CFA design problem as an optimization problem, we define the following variables. We denote , and the vectorized representations of the red, green and blue color atoms, respectively. To ensure physical realizability, we require , , , where , and we stack all design variables into one long vector
The design variable is related to the vectorized RGB and luma/chroma color atoms as
where the matrices are defined by (9) as
Given the design variable , we also need to analyze its spectrum. We write the 2D Fourier transform equation (13) as a matrixvector product by using the Fourier transform matrix . Hence, the vectorized spectra of the luma/chroma color atoms can be written in terms of as
(17) 
where , for . The relation between the matrix and the vector forms of the Fourier transform is:
(18) 
where is the Fourier coefficient.
Iii Design Criteria
We now present the design criteria. Our criteria unify the three major approaches in the literature: (i) Sensitivity of luma/chroma channels to noise by Condat [23]; (ii) Aliasing between different color components in the frequency domain by Hirakawa and Wolfe [11]; (iii) Crosstalk between neighboring pixels in the spatial domain by Anzagira and Fossum [21]. Note that the first two criteria were developed for CMOS, whereas the third criterion was developed for QIS. The proposed framework integrates all these criteria into a unified formulation.
Before we go into the details, we refer the readers to Table I for a comprehensive comparison of the literature and design criteria. The proposed optimization framework generalizes all conditions discussed in the previous work and offers a unified formulation. Some criteria are conflicting, e.g., crosstalk and spatial resolution. Thus, no CFA is uniformly superior than the others. The benefit of the proposed optimization is that it allows us to trade off one criterion with another without completely redefining the optimization.
Criterion  Purpose  Regular Pixels  QIS  

[11]  [23]  [14]  [26]  [21]  Ours  
Proposition 1  To minimize noise power after linear demosaicking  
Proposition 2  To simplify denoising of luminance channel  
Proposition 3  To simplify denoising of chrominance channel  
Proposition 4  To maximize spatial resolution  
Proposition 5  To mitigate crosstalk  
Definition 4  Enforce total orthogonality  
Definition 4  Enforce quadrature orthogonality 
Iiia Luminance and Chrominance Gain
Definition 1.
The luminance gain and the chrominance gain of a CFA with color atom of size are defined as
(19) 
where is a normalization factor.
Intuitively, the luminance and chrominance gains are measures of the signal power that can be transmitted through the color filter. A more transparent color filter allows more light to pass through, and hence the signal power is higher. This is reflected by the magnitudes for , which according to Parseval’s Theorem they are equivalent to .
The following proposition shows how can we compute and in terms of the optimization vector .
Proposition 1.
For a CFA with color atoms represented by the vector , the luminance and chrominance gain can be calculated as
(20) 
where , and .
Proof:
See Appendix A. ∎
The luminance gain and the chrominance gain cannot be arbitrarily chosen. One practical consideration is to ensure uniform noise power across the luma channel so that the denoising procedure can be simplified (because the noise will be i.i.d.). Thus, the luminance color atom should be constant, i.e., , where is a positive constant. Taking Fourier transform, this means that comprises only one impulse at baseband , and no impulses at all other frequencies. In vector form, we need
(21) 
where is the standard basis. Putting in terms of the optimization variable , we have a constraint.
Proposition 2 (Uniform Luminance Constraint).
If a CFA has a uniform luminance gain, then needs to satisfy
(22) 
Similarly, we can impose a constraint to the chrominance channel. For chrominance, we require that the color filter passes the same amount of red, green, and blue so that the primary colors have uniform gain [23, 17]. This leads to
Putting into vector form, we have the following constraint.
Proposition 3 (Uniform Chrominance Constraint).
If a CFA has a uniform chrominance gain, then needs to satisfy
(23) 
where , , and .
IiiB AntiAliasing
In the frequency domain, the luminance controls the baseband whereas the chrominance controls the sideband of the spectrum. To minimize spectral interference, aka aliasing, it is necessary to modulate the chrominance as far as possible from the baseband. However, the luminance has approximately a diamond shape since it has large frequency components at and . Therefore, to mitigate aliasing, we should avoid modulating the chrominance on vertical and horizontal axes. Figure 6 shows a CFA obtained by our framework. In this example, no chrominance components are modulated on the vertical and horizontal frequency axes.
The antialiasing requirement can be formulated as forcing the Fourier coefficients of the chrominance color atoms at and to zero for all and . This translates to chrominance color atom whose first column and first row are zeroed out (See Figure 3). In terms of our design variable , we require the following constraint.
Proposition 4 (Antialiasing Constraint).
The chrominance in the vertical and horizontal directions must be set to 0. Hence, must satisfy
(24) 
where and are the matrices formed by choosing the rows in and , respectively, that correspond to the vertical and horizontal frequency components, i.e., rows in the set .
To quantify the amount of aliasing for every CFA, we define the following aliasing criterion.
Definition 2.
For a CFA, aliasing between luminance and chrominance channels is measured by
(25) 
where , and denote the power spectral density of the luminance channel and the two chrominance channels and , respectively.
IiiC Crosstalk
Crosstalk is caused by the leakage of electrical and optical charge from a pixel to its adjacent pixels [30, 21]. Crosstalk leads to saturation of color, which is formally known as color desaturation. To model crosstalk, we follow [21] by defining three scalars , , and representing the proportion of leaked charges to neighboring pixels. These three scalars then form a crosstalk kernel,
(26) 
which can be considered as a spatial lowpass filter of the mosaiced image. Applying the crosstalk kernel to the CFA is equivalent to a spatially invariant convolution
where denotes the effective CFA in the presence of crosstalk.
The effect of crosstalk is more severe when the adjacent colors are different. For example, in Figure 5, the red and blue pixels are surround by 8 neighbors of different colors and the green pixels are surrounded by 4 neighbors of different colors. This is equivalent to saying that there is a red pixel having a value and is surrounded by pixels having the value . Using similar argument, we can see that if the color atoms have more rapid changes of colors, then the resulting CFA is more susceptible to crosstalk.
One way to quantify the variation of the color atoms (and hence crosstalk) is by means of measuring the total variation of the color atom. The total variation is a proxy of the complexity of the color filter array. A color filter array with high total variation means a more complicated pattern and so it is more susceptible to crosstalk. Our total variation is defined as follows.
Definition 3 (Total Variation).
For a CFA defined by the color atoms , and , the weighted total variation is defined as
(27) 
where is an operator that computes the vertical and horizontal derivatives with circular boundary conditions, and is the leakage factor defined in the crosstalk kernel in (26).
(a) Red Channel  (b) Green Channel  (c) Blue Channel 
To control the amount of variations in the CFA (so that we can limit the amount of crosstalk), we upper bound the total variation by a scalar . This leads to the following constraint.
Proposition 5 (Crosstalk Constraint).
The crosstalk is limited by upperbounding the total variation metric :
(28) 
Figure 6 shows two CFAs proposed in literature. The first one, proposed in [14] is more robust to aliasing than the second one proposed in [21]. This is because the chrominance channels are modulated at high frequencies which are far from baseband luminance. However, [21] is more robust to crosstalk than [14] because the color atom have less variation in colors. We can also see this in the total variation values ( compared to ). This tradeoff constitutes a gap in literature: Color filter designs can improve robustness of either aliasing or crosstalk, but not for both. Our proposed design framework allows us to optimize them simultaneously.
[14] 


RGBCWY [21] 

(a) Red Mask  (b) Green Mask  (c) Blue Mask  (d) CFA Atom  (e) Spectrum 
IiiD Orthogonality of Chrominance Channels
When designing a CFA, one should take into consideration of the complexity of the demosaicking process. Recall (IIC) where we show that
This is a modulation of the signal by a modulating frequency . Therefore, to reconstruct the signal, one approach is to demodulate by shifting the channels to the baseband by multiplying pure sinusoids and then applying a lowpass filter [23]. Demodulation can be done efficiently if there is orthogonality between the channels. Following the literature, our optimization takes into account of two forms of orthogonality.

Quadrature Orthogonality [23]: The idea is to make one chroma channel real and the other imaginary at any , i.e., the two channels are modulated by a frequency but in quadrature phase. Translating the spatial domain, this means that
(29a) (29b) where , , and is the phase angle. In this way, the two channels can be easily separated during the demosaicking process using the orthogonality of cosine and sine functions.
We formulate the orthogonality criteria as a penalty function that is a surrogate of both approaches.
Definition 4 (Orthogonality Penalty).
For a CFA having chrominance channels with spectra and , the orthogonality penalty is defined as
(30) 
which can be written as a function in as follows
(31) 
Looking at the first summation in (4), we notice that for every 2D frequency , the term is the norm of a 2dimensional vector . Therefore, minimizing this norm promotes either one of the components to zero (or both). Similar argument applies for the imaginary components in the second summation. As a result, the total variation can be regarded as a proxy to the orthogonality condition.
Iv Formulation of Optimal CFA Design Problem
Using the variables and constraints defined in the previous section, we present two different optimization formulations of the CFA design problem in this section: (i) A nonconvex formulation that integrates all the above information into a single optimization, and (ii) convex relaxation that makes the problem more tractable.
Iva NonConvex CFA Design
The nonconvex CFA optimization puts all the objectives and constraints defined in the previous section into a single optimization problem. This gives us
(32)  
where and are the regularization parameters controlling the relative weights of the luminance gain and the orthogonality penalty. The penalty function is added to the objective with a negative sign so that it is minimized. By lower bounding with a constant , we can rewrite (32) as
(33)  
In this optimization problem, the objective and constraints are convex except for (33)(f) and (33)(g). This is because these inequalities include convex quadratic form in the “” side, where convexity comes from the fact that and are positive semidefinite matrices. Hence, the optimization problem is nonconvex.
IvB Solving the Optimization
While problem (33) is nonconvex, we can find a local minimum by successive convex approximations [31]. The idea of successive convex approximation is to replace the quadratic terms in the nonconvex constrains (33)(f) and (33)(g) by first order approximations around the initial guess . Since the quadratic form is convex, its first order approximation constitutes a lower bound. Hence, we are replacing the nonconvex constraints (33)(f) and (33)(g) with convex but tighter constraints that limit the feasible set of . The algorithm repeats until converges to a fixedpoint, which is the final chrominance gain.
The overall algorithm is summarized in Algorithm 1. Figure 7 shows the convergence of Algorithm 1 for designing a color atom. We notice the monotonic increase of until it converges to a fixed point. Since the original problem is nonconvex, solution to the problem could be a local minimum depending on how the initialization is done. In practice, we solve the problem for multiple instances with different randomly generated initial guesses which approximately cover the design space (e.g., using the Latin hypercube sampling [32]), and pick the best solution among them.
IvC Convex CFA Design
The relaxation from nonconvex to convex can be done by explicitly forcing part of the chrominance components to zero. Specifically, we modulate the chrominance channels on the same frequency 29). In terms of , these two equations can be written as: using the quadrature orthogonality mentioned in (
(34) 
where and are constant vectors that represent the vectorized version of the cosine and sine signals on the right hand side of (29a) and (29b), respectively, i.e.,
(35a)  
(35b) 
Since we explicitly choose the modulation frequencies of chrominance channels manually, we can drop the aliasing constraint in Proposition 4. However, we still need the uniform luminance and chrominance constrains in Propositions 2 and 3. Moreover, since the luminance and chrominance gains are adversarial, the objective of this formulation is to maximize their weighted sum. To this end, the problem is written as:
(36)  
In our terminology, the optimization problem of [23] is obtained from (36) by removing the crosstalk constraint (36)(d). Hence, our optimization limits the search space of the optimization in [23] to get CFAs that have acceptable crosstalk performance.
Figure 8 shows two color atoms obtained using the convex and nonconvex formalizations. In the convex formulation, we select the modulation frequency to be and the phase that maximizes at this frequency is found to be . Then, we solve the problem to get . As for the nonconvex formulation, we let the optimization to choose modulation frequencies subject to crosstalk and aliasing constraints. Solving the nonconvex formulation yields . We notice that the nonconvex formulation achieves higher chrominance gain because of its flexibility in choosing the modulation frequencies.
Color Atom 


Spectrum 

(a) Convex Formulation  (b) NonConvex Formulation  
V Universal Demosaicking
In this section, we present a universal demosaicking algorithm. The algorithm is called universal because it can be applied to any CFA designed by our framework. This is in big contrast to many existing demosaicing algorithms which are customized for certain types of CFAs, e.g., frequency selection [11]. Our algorithm comprises two main parts: (i) a demosaicking step to remove the color filtering effect (Section VB) and (ii) a color correction step to mitigate the crosstalk effect (Section VC).
Va Special Consideration for QIS.
Before we talk about the demosaicing algorithm, we should briefly discuss the photon statistics of QIS. In CMOS, the measured voltage can be modeled as a nominate value corrupted by i.i.d. Gaussian noise. In QIS, previous work showed that the measured photon counts follow a truncated Poisson process [6]. When averaging over a number of temporal frames, the truncated Poisson becomes a Binomial [7]. If the photon count is sufficiently high, this binomial will approximately approaching to a Gaussian. Applying the law of large numbers on the distribution of in (7), the average is
and the maximumlikelihood estimate of the signal is
As discussed in [7], we can regard this equation as a tonemapping of the photon counts. We regard as the th pixel of the mosaicked image generated by the CFA. The goal of demosaicing is to reconstruct a color image from .
VB Demosaicking Algorithm
By recalling the forward model (5), we can write the inverse problem for obtaining the latent color image im from the light exposure on QIS as follows
(37) 
where is the color filter sampling operator. The first term in the cost function is a datafidelity term that forces to agree with the measurements . The second term is a regularization term to improve the conditioning of our illposed problem. is a positive scalar that controls the amount of regularization.
To solve the inverse problem (37), we may use any optimization toolbox since it is convex. Here, we report our results using the PlugandPlay (PnP) ADMM algorithm [33], which has demonstrated effectiveness in image restoration tasks. Starting from an initial guess , the PnP ADMM algorithm iteratively updates its estimate via two steps:
Demosaicking Module:
(38) 
Denoising Module:
(39) 
and updates the Lagrange multiplier by . For additional details on PnP ADMM, we refer the readers to, e.g., [33]. Here, is an internal parameter that controls the convergence. The operator is an offtheshelf image denoiser, e.g., BM3D in our experiments. The subscript denotes the denoising strength, i.e., the hypothesized “noise variance”. The inversion in the demosaicking module is performed in closed form because exhibits a block diagonal structure.
VC Color Correction
The optimization problem in (37) does not take into account of the crosstalk effect. ^{1}^{1}1In principle we can incorporate the crosstalk kernel into the matrix but then will have a complicated structure which does not allow simple inversion. Like most of the mainstream image and signal processing (ISP) pipelines, we reduce the crosstalk via a color correction step.
Before Color Correction  After Color Correction 
Given the demosaicked color pixel , the color correction multiplies by a matrix such that is the colorcorrected pixel. The matrix is learned by comparing a measured color pixel to a known color chart value. Mathematically, suppose we have measured color pixels forming a matrix , and a corresponding true color values forming another matrix , we can estimate by solving
(40) 
To minimize the noise amplification, it is advised to add regularization when estimating [34]. Since a standard color chart comprises color patches, we can estimate the noise by computing the norm of covariance matrix of every color patch, and get the average value over the color patches. Hence, the optimization problem is rewritten as
(41) 
where represents the pixels of the th color patch, and is a positive scaler that controls the noise amplification effect. By varying , we can draw a tradeoff curve between color reproduction accuracy and noise amplification.
Figure 10 shows reconstructed images before and after color correction. In this figure, the crosstalk parameters are . We can see the effect of color correction in the more saturated red and yellow feathers and in the green leaves in the background.
Vi Experimental Evaluation
In this section, we present CFAs obtained using our optimization framework in Section VIA. Afterwards, we evaluate the performance of different CFAs using the universal demosaicking algorithm proposed in Section V. First, using the Macbeth color chart, we show the noisecolor tradeoff of our robust CFAs comparsed to others in Section VIB. Second, we show in Section VIC a quantitative and qualitative comparison of the reconstruction performance of all CFAs on Kodak [35] and McMaster [36] color datasets.
Via Experiment 1: Proposed Solutions of CFA Design Problem
We focus on the nonconvex formulation (33) since it is more flexible than the convex formulation (36). We set the parameters of (33) as and . We run multiple instances of Algorithm 1 ( instances) using different random initializations for the color atoms . Then, we pick the solution with the highest chrominance gain. To ensure that the initial guess spans the feasible set of , we use uniform Latin hypercube sampling of the domain .
Table II shows our proposed CFAs compared to other CFAs in the literature. For every CFA, we compute 1) the luminance and chrominance gains ( and ) in (19) to measure robustness to noise, 2) the total variation metric (Proposition 5) to measure robustness to crosstalk, and 3) the aliasing metric in (25) to measure robustness to aliasing.
[23]  [11]  [26]  RGBCY [21]  RGBCWY [21]  [15] 
Ours  Ours  Ours  [14]  Ours  Ours 
Size  CFA Pattern  CFA Parameters  CPSNR on McM Dataset  CPSNR on Kodak Dataset  

TV  w/o Ctk  w/ Ctk  w/o Ctk  w/ Ctk  

Hao et al. [14]  0.577  0.109  0.413  0.094  21.69  26.81  27.92  29.68 
RGBCWY [21]  0.577  0.125  0.263  1.354  30.07  29.86  31.14  30.80  
Ours  0.090  0.263  0.147  29.94  29.90  31.25  30.45  

Cheng et al. [26]  0.577  0.167  0.350  1.638  29.39  29.11  29.50  28.52 
Ours  0.088  0.325  0.097  30.78  30.13  31.32  31.00  

Condat [23]  0.250  0.633  0.182  31.13  30.57  33.29  32.59  
Ours  0.187  0.554  0.157  28.37  32.03  33.22  32.68  

Hirakawa et al. [11]  0.866  0.125  0.550  0.173  26.49  30.23  31.59  31.28 
Ours  0.187  0.513  0.181  26.72  30.70  32.04  32.01  

Wang et al. [15]  0.577  0.200  0.443  0.178  30.47  30.00  30.21  29.77 
Ours  0.191  0.440  0.177  31.54  30.81  31.01  30.37 

: Among all CFAs in Table II, [14] is the most robust CFA to aliasing, but the least robust to crosstalk; whereas RGBCWY [21] is the most robust to crosstalk, but it is not as robust to aliasing. Our CFA achieves the best of both worlds by having the same totalvariation like RGBCWY, and good aliasing metric.

: Compared to [23], our CFA is more robust to crosstalk and aliasing.

: Compared to [11], our CFA has higher chrominance gain and it is more robust to crosstalk.

: Compared to [15], our CFA is more robust to crosstalk and aliasing.
To calculate the aliasing metric for [14], [15] and [11], we use the transformation matrices proposed by the authors to ensure that the chrominance channel are orthogonal in these bases.
ViB Experiment 2: ColorNoise Tradeoff
In this experiment, we compare the tradeoff between noise amplification and color accuracy of our proposed CFAs and other CFAs in literature. To do so, we use the Macbeth color chart that comprises 24 color patches. The forward model consists of mosaicking using a CFA and crosstalk using the crosstalk kernels (26) with . QIS parameters are , and . We use the PlugandPlay ADMM algorithm for image demosiacking with TV denoiser and . Since the ground truth color values of Macbeth color chart are known, we compute the color correction matrix by solving (41).
To draw the noisecolor tradeoff curve, we vary the parameter in (41) from to on the logscale. Color error is obtained by calculating the mean square error of the color channels after converting to CIELAB color space. Noise is obtained by calculating the norm of each color patch covariance matrix as in (41). Since both noise and color error should be minimized, the tradeoff curve is better when the total area under it is smaller.
Figure 12 shows the tradeoff curves for the proposed CFAs and other CFAs. Our , and CFAs are better than other CFAs for nearly all values of . This happens because they are more robust to crosstalk by design. As for CFA, if we restrict to small color error, then our CFA is better. However, if we allow larger color error, then [11] is better.
Figure 13 shows the ground truth and reconstructed Macbeth charts using CFA proposed in [11] and our robust CFA. The first and second rows show the reconstructed images before and after doing the color correction. Our CFA achieves higher quality (measured in PSNR) than the others.
(a)  (b) 
(c)  (d) 
ViC Experiment 3: Color Image Reconstruction
In this experiment, we perform color image reconstruction using the 24 and 18 color images in Kodak and McMaster datasets, respectively. QIS parameters are