Introduction
Multicharts and Multitest can calculate a color correction matrix (CCM) from an image of a color test chart. The chart should have a minimum of 9 color patches for 3×3 CCMs (which are generally recommended; 12 for 4×3 CCMs). Excellent results can usually be achieved with the inexpensive, widelyavailable 24patch XRite Colorchecker. The CCM can be applied to images to achieve optimum color reproduction, which is defined as as the minimum mean squared color error between a corrected test chart patches and the corresponding reference values, which can be the standard (published) values, measured with a spectrophotometer, or the saved L*a*b* values from another image (when the intent is to match images).
The default color error parameter is the mean of (DeltaE 94)^{2}, calculated for all patches where 5 ≤ L* ≤ 98 ( where L* is defined in CIELAB color space) and each of the R, G, and B channels are less than 99% of the maximum value (i.e., the patch is unsaturated). Nonlinear optimization is used to calculate the matrix.
The CCM can be used in the Input Device Transform (IDT) of cinema workflows, such as the Academy Color Encoding System. It can be used to achieve consistent color when several different cameras are used to capture a scene.
In setting up color matrix calculations, it is very important that the image be properly linearized because the matrix operates on linear R, G, B values. Several linearization options are available for different image encodings and color spaces (or lack thereof).
The Math
The matrix
Color images are stored in m x n x 3 arrays (m rows (height) x n columns (width) x 3 colors). For the sake of simplicity, we transform the color image to a k x 3 array, where k = m x n. The Original (uncorrected input) pixel data \(O\) can be represented as
\( O = \begin{bmatrix} O_{R1} & O_{G1} & O_{B1} \\ O_{R2} & O_{G2} & O_{B2} \\ & … & \\ O_{Rk} & O_{Gk} & O_{Bk} \\ \end{bmatrix} \)
where the entries of row i, \([O_{Ri}\quad O_{Gi}\quad O_{Bi}]\) , represent the normalized R, G, and B levels of pixel i. The transformed (corrected) array is called \(P\), which is achieved by matrix multiplication with \(A\), the color correction matrix. Imatest allows you to choose two different forms of \(A\), either 3×3 or 4×3.
Case 1: \(P=O \, A\) (\(A\) is a 3x3 matrix)
Each of the R, G, and B values of each output (corrected) pixel is a linear combination of the three input color channels of that pixel.
\(P = \begin{bmatrix} P_{R1} & P_{G1} & P_{B1} \\ P_{R2} & P_{G2} & P_{B2} \\ & … & \\ P_{Rk} & P_{Gk} & P_{Bk} \end{bmatrix} = \begin{bmatrix} O_{R1} & O_{G1} & O_{B1} \\ O_{R2} & O_{G2} & O_{B2} \\ & … & \\ O_{Rk} & O_{Gk} & O_{Bk} \end{bmatrix} \begin{bmatrix} A_{11} & A_{12} & A_{13} \\ A_{21} & A_{22} & A_{23} \\ A_{31} & A_{32} & A_{33} \end{bmatrix}\)
— or —
Case 2: \(P=[O \;\mathbf{1}] \, A\) (\(A\) is a 4x3 matrix)
In this case, a column of 1’s is appended to the A matrix in order to also determine the best offsets for each color channel as well, indicated by the \(A_{41}\), \(A_{42}\), and \(A_{43}\) values. This makes the corrected color values an affine transformation of the original inputs.
\(P = \begin{bmatrix} P_{R1} & P_{G1} & P_{B1} \\ P_{R2} & P_{G2} & P_{B2} \\ & … & \\ P_{Rk} & P_{Gk} & P_{Bk} \end{bmatrix} = \begin{bmatrix} O_{R1} & O_{G1} & O_{B1} & 1 \\ O_{R2} & O_{G2} & O_{B2} & 1 \\ & … & \\ O_{Rk} & O_{Gk} & O_{Bk} & 1 \end{bmatrix} \begin{bmatrix} A_{11} & A_{12} & A_{13} \\ A_{21} & A_{22} & A_{23} \\ A_{31} & A_{32} & A_{33} \\ A_{41} & A_{42} & A_{43} \end{bmatrix}\)
Note: These equations assume that original (input) array \(O\) is a linear representation of the image brightness. It it is not linear (e.g., in sRGB or Adobe RGB color space), \(O\) must be linearized prior to the matrix equation and may need to have a gamma curve applied (i.e., be delinearized) afterwards. See Linearization, below, for details. 
The goal of the calculation is to minimize the difference (the mean square error metric) between \(P\) and the reference array (the ideal chart values) \(R\). The initial values of \(A\) (the starting point for optimization) for the 3x3 and 4x3 cases are, respectively,
\[ A_{(3\times3)} = \begin{bmatrix} k_R & 0 & 0 \\ 0 & k_G & 0 \\ 0 & 0 & k_B \end{bmatrix} \; ; \qquad A_{(4\times3)} = \begin{bmatrix} k_R & 0 & 0 \\ 0 & k_G & 0 \\ 0 & 0 & k_B \\ 0 & 0 & 0 \end{bmatrix} \]
where
\( k_R = \mathrm{mean}(R_{Ri})/\mathrm{mean}(O_{Ri}) \)
\( k_G = \mathrm{mean}(R_{Gi})/\mathrm{mean}(O_{Gi})\)
\( k_B = \mathrm{mean}(R_{Bi})/\mathrm{mean}(O_{Bi})\)
for reference array \(R\) and original array \(O\) where the mean is taken over all \(i\).
These starting values are closer to the final values (have less mean square error) than the identity matrix \((k_R = k_G = k_B = 1)\). They tend to converge slightly better.
Linearization
Although most digital image sensors are linear up to the point where they saturate, many image files are highly nonlinear. Imatest deals with two broad types of image file.
 Files designed to be interchangeable and to be displayed at a specified gamma ( γ ), where display luminance = pixel level^{γ}. These files generally have a color space associated with them (which may or may not be embedded in the file). Gamma = 2.2 for the most commonly used color spaces (sRGB (approximately 2.2) and Adobe RGB (1998)), although some less popular color spaces are designed for display at gamma = 1.8 (ProPhoto, ColorMatch; all RGB).
When cameras encode images (a part of the RAW conversion process), they apply a gamma that is an approximate inverse of the display gamma: it may vary considerably from 1/γ, and it often includes a “shoulder” (an area of reduced contrast) in the tonal response curve to minimize highlight burnout. The shoulder improves pictorial quality by making the response more “filmlike”.
 “Raw” files, which are converted from the sensor output without applying a gamma curve (and often little or no signal processing, i.e., sharpening, noise reduction, white balance, etc.) These files are linear with gamma close to 1 (often somewhat below 1 because of lens flare and occasional onsensor signal processing).
If the input image has been gammaencoded it must be linearized prior to applying the correction matrix. Imatest has several linearization options, described below. 1, 5, 6, and 7 linearize the input image. 24 assume the input image is linear. You can find the gamma of the input image for Multicharts or Multitest (for charts that include a grayscale pattern) by viewing the B&W density & White Balance display, shown on the right. This curve has a shoulder.
Optimization steps
 Check the gamma of the input image and linearize it if necessary. See Linearization options, below. If Color space gamma linearization is selected, O_{L} = O^{gamma}. (where gamma is specified by the color space, e.g., 1/2.2 for Adobe RGB).
 Call the optimizer, which
 calculates a (temporary) corrected image T_{L} = O_{L} A.
 (Optionally) removes the linearization: T = T_{L}^{(1/γ)}
 Finds the mean of squares of errors between T and the reference (ideal) array R. The error is one of the standard error measurements: ΔE*_{ab}, ΔC*_{ab}, ΔE_{94}, ΔE_{94}, ΔE_{CMC}, ΔC_{CMC}, ΔE_{00}, or ΔC_{00}, described here. ΔE_{94} is the recommended default. Although the CIEDE2000 color error metrics (ΔE_{00}, …) are more accurate, they contain small discontinuities that can affect optimization, and hence they should be used with caution. See Sharma for details. The ΔC errors are similar to ΔE with luminance (L*) omitted.
 Adjusts A until a minimum value of the sum of squares of the errors is found, using nonlinear optimization.
 Reports the final value of A.
In applying A (generally outside of Imatest), a similar linearization should be used. A may be applied without linearization during the RAW conversion process, prior to the application of the gamma + tonal response curve.
There is no guarantee that A is a global minimum. Its final value depends to some extent on its starting value.
The Color correction matrix in Multicharts and Multitest
Options: Color matrix calculation options can be set by clicking Settings, Color matrix in the Multicharts window or by pressing the button in the Multitest settings window. This brings up the dialog box shown below. The options are
 4×3 or 3×3 matrix: Color correction matrix size. 3×3 is generally recommended. 4×3 matrices include a dcoffset (constant) term that slows computations and rarely improves the results.
 Optimize: Select the color error parameter whose mean of squares over patches with L*>10 (nearly black for 95) is to be minimized. Choices include ΔE ab, ΔC ab, ΔE 94, ΔE 94, ΔE CMC, ΔC CMC, ΔE 00, and ΔC 00, described here. ΔE 94 or ΔE 2000 are recommended because they give less weight to chroma differences between highly chromatic (saturated; large a*^{2} + b*^{2 }) colors, i.e., they are closer to the eye’s perception than ΔE ab (the standard ΔE calculation: the geometrical distance between colors in L*a*b* space).
 Weighting: Set the weighting of the patches for optimization. Choices:
 1. Equal weighting [default] Patches are equally weighted. May not be the optimum setting because highlight colors may be more visually prominent than shadow colors.
 2. Emphasize highlights: weight according to CIELAB L* Give more weight to visuallyprominent highlight patches. This is the default.
 3. Strongly emphasize highlights: weight according to L*^2.
 4. Read weighting from a CSV file If selected, the Copy weights, browse, and file boxes are activated. The CSV file should have the same number of weights (typically in the range of 0 to 1) as patches in the chart. The file may one value per line or m values per line on n lines, as long as m*n = the number of patches in the chart.
Here is an example of weights for the 140patch (10 rows x 14 columns) Colorchecker SG.
0,0,0,0,0,0,0,0,0,0,0,0,0,1 
Pressing the Numbers checkbox (on the right of the Multicharts window) when a weight file has been entered and the CCM has been calculated displays both the patch number (the Colorchecker 24 replica section of the SG is shown) and the weights. This can help verify that the correct weights have been entered. 
 Linearization: Choose the method of linearizing the image prior to the matrix calculation and encoding (delinearizing) it after the calculation. You can assume that the input image is linear, gammaencoded (using the gamma ( γ ) of the selected color space (sRGB, Adobe RGB, etc.), a data fit, or a manuallyentered gamma), or logarithmically encoded. The equation for linearization of gammaencoded color spaces is O_{L} = O^{1/γ}. The equation for gamma encoding O = O_{L}^{γ}.
 Linearize input (CS gamma), apply matrix, then apply CS gamma encoding. Usually a good choice for images encoded in a standard color space (CS).
 No linearization: apply matrix to input pixels. For experimentation and for images that start and remain linear. Not a good choice for gammaencoded images.
 Assume linear input, apply matrix, then apply CS gamma encoding. A good choice for images that are not gammaencoded, but need to be converted into a color space. Color space (CS) gamma is not applied in converting the input RGB image to L*a*b*.
 Assume linear input & output (gamma = 1 throughout). No linearization or gamma encoding. Color space gamma is not applied in converting between RGB and L*a*b* spaces. (Same as 2)
 Linearize input (data fit), apply matrix, apply CS gamma encoding (most accurate) Uses a polynomial fit to the grayscale patches to linearize the image. Optimizes linearized image, then gammaencodes it in the Color Space gamma. This often gives good results, but sometimes fails when the inverse of the polynomial fit (which is used for linearization) is becomes unstable (diverges at the extremes).
 Linearize input (data fit), apply matrix, then delinearize (data fit) Same as 5, except that after optimization the image is restored to the original tonal curve.
 Linearize input using gamma (below), apply matrix, then apply CS gamma encoding. Manually enter the gamma, which can be obtained from the B&W density plot if the color chart has a grayscale pattern. This produces excellent results in the frequent cases where the measured gamma differs from the expected color space gamma or 1 (for linear raw conversion). This setting, introduced in Imatest 4.0, is particularly useful when gamma of the input image differs from the standard color space gamma (often 1/2.2) or linear (1).
 Linearize using log slope (for Log ColorSpace, below), apply matrix, then gamma encode (CS) (Imatest 4.0+). For logarithmic color spaces (sometimes used in cinema; never in still photography), which are characterized by the approximate encoding equation, O = max(1+A log_{10}(O_{L}), 0) and inverse (for linearization) O_{L} = 10^((O1)/A). A is Log slope, which is displayed in the B&W density & White Balance display in the Multicharts module in Imatest 4.0+. Logencoded color spaces have approximately straight line response when the Pixels vs. Input Density option (new in 4.0) is selected, while gammaencoded spaces have approximately straight line response when the more traditional Log(Pixels) vs. Input Density option is selected.
Note: For the Rec. 709 color space, which has allowable values between 16 and 235 (of 255), linearization maps 16235 to the full 01 (internal) range and delinearization maps 01 back to 16235. The image is displayed full range.
 Gamma (only enabled for Linearization setting 7: Linearize using input gamma…). Manually enter the value of gamma measured from a grayscale stepchart (which is included in most color charts and is displayed in the B&W density plot).
 Chroma multiplier (for reference) Normally 1. Increase it for more saturated colors; decrease it for less saturated colors. Sometimes results are improved by increasing it to 1.05 or even 1.1.
 Optimization constraint
 No constraints. This is the default and should only be changed with good reason.
 Rows sum to 1.
 Columns sum to 1.
 Calculation details Two types of calculation are available: (1) Direct RGB_{in} to RGB_{out}, and (2) RGB_{in} to XYZ (indirectly to RGB_{out}). Both use an optimizer: they start with a first guess, then the optimizer adjusts the values until a minimum value of the selected error metric is found.
 Calculate RGB matrix directly (original algorithm) Start with a guess for the RGB_{in} to RGB_{out} matrix, then find L*a*b*_{out} from RGB_{out}. Under control of the optimizer adjust the matrix so the color error (between L*a*b*_{out} and the reference values (adjusted by the chroma multiplier if ≠ 1) is minimum.
 Calculate XYZ matrix and xy primaries (newer algorithm) Start with a rough guess for the RGB_{in} to XYZ matrix (sRGB coefficients are used). Then convert to XYZ to L*a*b* and adjust the matrix under control of the optimizer until the selected error metric is minimized. When optimization is complete, convert XYZ to RGB_{out}. Both RGB to XYZ and RGB to RGB matrices are calculated and included in the CSV and JSON output files.
 Reference file (in the Multicharts window; not in the Color Matrix window) The reference file contains the ideal chart L*a*b* values, used for the optimization. Normally it contains colorimetricallymeasured L*a*b* values for the chart patches. This is always the case for the Default values and typically the case for reference files, but reference files may contain different values if a special “look” is desired, as described in Matching images, below.
Matching images If you have an image with high quality colors (it might be a “gold standard”) that you want to match, you can save the L*a*b* values in a named file by clicking File, Save L*a*b* results as CSV reference in the Multicharts window. Then you can use this file as a reference file by clicking on the Ref (Reference file) box, selecting LAB data file, and selecting the saved file. 
Acquire an image of the test chart, taking care, if possible, to avoid saturating any of the patches. If patches are saturated, i.e., if any of the channels have an average pixel level above 99% of the maximum for the bit depth (e.g., 255 for bit depth = 8; 65535 for bit depth = 16) they will be omitted from the optimization and they will be labeled as Sat in the split view. It may not be visually obvious when a patch is saturated, especially if only one of the three channels is saturated (as shown below). Saturated patches tend to look very strange (see below) when the matrix is applied. (This may appear to be a bug, but it isn’t, as explained below.)
Saturation
An image with a saturated patch (row 4, left) is shown on the right. It’s not obvious because only the green channel is saturated. This is apparent from the density response curve (from Multicharts or Multitest) shown below: the green channel is saturated at maximum exposure and the slope of the green curve drops just below maximum exposure. As a result of this, the green channel has a lower pixel level (relative to the other patches) than it would have had if it were not saturated. 

The corrected value of the saturated patch on the lowerleft (shown on the bottom– it’s magenta!) is way off because the saturated green channel is lower than it would have been if the patch were not saturated.  
Note: Saturated patches are not included in the matrix calculation. Though we recommend avoiding them, they will not have a major effect on the calculated matrix. 
To calculate the color correction matrix,

read the image into Multicharts or Multitest.

Find the gamma from the density response curve in the B&W density display
 Set the Color matrix settings, paying special attention to linearization.
 Select the reference file. Most of the time the standard default values are used, but you can enter a reference file created from spectrophotometer measurements or by saving the L*a*b* values of a good (“gold standard”) image.
 Press the button, shown on the right. The display will change, as shown below: Ideal (reference) colors remain on the upperleft of each patch, input colors remain on the upper right, and corrected colors are now displayed on the bottom. The input image (from our archive) has relatively poor white balance. The improvement is quite dramatic.
Split view, showing reference, input, and corrected patch colors. The dropdown menu that shows
Display input (on the right) can be changed to Display corrected to get [Corrected – ideal] statistics.
The , highlighted with a pink background. The correction matrix cannot be recalculated until an image property changes (new image, color space, reference file, or color matrix setting). The Display input (or Corrected) dropdown menu, immediately to its left, is enabled. You can choose one of two selections.
button changes toThe EXIF data and Color matrix display contains a summary of results.
The color correction matrix, results summary, and both [input – ideal] and [corrected ideal] color difference summaries are shown. The initial and final error numbers shows how much the selected metric (in this case the sum of squares of DeltaE 94 for all patches with L*>10 and L*
Raw image results
Here is an example of Color matrix optimization for a raw image. The original image was acquired as a raw (RW2) image (a JPEG was also acquired for reference) from a Panasonic GF1 Micro FourThirds (mirrorless interchangeable lens) camera, and converted using dcraw with Manual settings, Normal RAW conversion, Output gamma = 1.0 (Linear), Auto white level off (unchecked), White Balance = None, Output color space = RAW. These settings applied no image processing other than demosaicing.
The file was almost perfectly linear: the tonal response curve was a straight line with gamma = 1.01. Linearization option 3 (Assume linear input, apply matrix, then apply CS gamma encoding.) was selected.
The output is nearly perfect; about as good as it can get for a 3color sensor. This is typical of what you might expect working from wellcharacterized raw images. We recently (May 2014) tested the calculation with four cameras that have RAW output (Panasonic Lumix LX7 and G3, Canon EOS6D, and Nikon D800E) under six different lighting conditions (2700K and 5000K LED, 3500K and 5000K CFL, standard (2700K) incandescent, and 4700K (measured 4000K) SoLux (dichroicfiltered) halogen lamps), and we achieved remarkably consistent colorcorrected images.
Results from a linear raw file (gamma = 1; no color correction)
x,y and X,Y Z primaries Starting with Imatest 4.0, Multicharts and Multitest can calculate the x,y and X,Y Z primaries of the input (camera) color space, display them (the bottom three lines on the right), and save them in CSV and JSON output files. This calculation works best when the input image is linear (converted from raw without a gamma curve) or in ACES color space. It can be useful for calculating signal processing coefficients. Note that the x,y primaries are unrelated to the color gamut of the camera. These primaries are the the x and y values corresponding to R,G,B = {max,0,0}, {0,max,0}, [0,0,max}, which are unattainable in practice because of color crosstalk in the sensor and are often outside the spectral locus of real colors.
Saving the matrix from Multicharts


Apply the matrix – Colorcorrect images

Uncorrected and corrected images 
This method works only with single images, not batches.
2. Apply the matrix in the Imatest Image Processing module
In Imatest 5.0+ the Image Processing module can be used to apply a Color Correction Matrix (saved from Multicharts as an individual CSV file, as described above) to individual images or to batches of images. To do so, check the CCM (color matrix) checkbox on the left of the Image Processing window. The matrix is applied before the highly nonlinear tone mapping and bilateral filter operations.
CCM settings can be entered or viewed by clicking the Set CCM button to the right of the checkbox. The key settings include linearization (described above), Color space, input gamma (for one linearization setting), and the matrix multiplier (which is can lighten or darken the image; it is normally left at 1).
CCM settings for Image Processing
Image Processing can operate on batches of images simply be selecting multiple images instead of a single image.
3. Apply the matrix externally with Matlab code
The CCM can be applied outside of Imatest. We show an example, using Matlab, of how to apply the matrix. You can try this out if you have Matlab (or convert it into another language like C or Python if you don’t).
The 24patch Colorchecker image shown on the right is used as the input to the Color Correction Matrix calculation. It was derived from a raw image, converted to a JPEG file with no color correction, then reduced. You can click (or rightclick) on it to download it.
This image was read into Multicharts, then colorcorrected with the following settings.
Optimize: DeltaE 00
Linearization: Assume linear input, apply matrix, then apply Color Space gamma encoding.
Weighting: 2. Emphasize highlights
Chroma multiplier: 1.06
Calculate XYZ matrix…
Here are the results of the correction. The corrected colors are extremely accurate— better than we usually get when we apply the matrix to an image that has already been processed (ΔE_{ab} mean = 5.32; ΔE_{00} mean = 3.11).
The Color Correction Matrix was obtained by clicking File, Copy color matrix to clipboard, then pasting directly into this page.
1.48663 0.192007 0.033443
0.612018 2.03673 0.796356
0.367863 0.580001 3.04927
Here is a program for applying the CCM. You can select the text, copy it to the clipboard (controlC), paste it into the Matlab (or other) editor, then save as an mfile (ccm_test.m).
function ccm_test
% Test program for applying Color Correction matrix.
% Output for display gamma = 2.2  sRGB and Adobe RGB color spaces.
% Read image file
disp('Enter the file for the Color Correction Matrix (CCM)');
imagefile = input('','s');
gamma = input('Enter encoding gamma (typically 0.5 to 1) ');
imageRGB = imread(imagefile);
imtype = class(imageRGB); % Image type: typically uint8 or uint16.
imageRGB = double(imageRGB)/double(intmax(imtype)); % Normalize to maximum for image type.
figure; image(imageRGB);
title(['Original image: gamma = ' num2str(gamma)]); % Works without Image Processing Toolbox.
linearRGB = imageRGB.^(1/gamma); % Linearize the image (apply inverse of encoding gamma).
% The CCM has 3 lines, which makes it tricky to enter.
disp('(Calculate a 3x3 CCM in Multicharts, then press File, ''Copy Color Matrix to Clipboard'').');
disp('Press ''enter'' when the the CCM is in the clipboard');
inp = input('','s'); % pause fails after copy to clipboard. inp is not used.
ccm = clipboard('paste');
ccm = strrep(ccm,10,';'); % Bring back semicolons.
ccmnum = str2num(ccm), % str2double doesn't work.
% Change to 2D to apply matrix; then change back.
[my, mx, mc] = size(linearRGB); % rows, columns, colors (3)
linearRGB = reshape(linearRGB,my*mx,mc);
correctedRGB = linearRGB*ccmnum;
correctedRGB = min(correctedRGB,1); correctedRGB = max(correctedRGB,0); % Place limits on output.
correctedRGB = correctedRGB.^(1/2.2); % Apply gamma for sRGB, Adobe RGB color space.
correctedRGB = reshape(correctedRGB, my, mx, mc);
figure; image(correctedRGB); title('Corrected image');
Read the following instructions carefully before running the program.
 Download and save the uncorrected image (above).
 Select the program text (above), paste it into an empty (new) Matlab (or text editor) window, save it as ccm_test.m, and run it.
 Enter the full path name of the uncorrected image file.
 Enter the value of the encoding gamma (1 for this file: you can derive this from the density plot in Multicharts).
 Now the tricky part. Copy the CCM into the clipboard. You can copy it from this page (above), or run Multicharts yourself, following the above instructions.
 Click return. The calculations will proceed.
The output file display is shown on the right. You can add code to save it if you wish.
Links
Some of the background for the calculation can be found in Color Correction Matrix for Digital Still and Video Imaging Systems by Stephen Wolf, though the Imatest calculation differs in many respects: there is no issue with outliers and optimization is performed using one of the standard color difference metrics.