In two previous articles (1, 2) I've spoken about how to perform image color correction using correction matrix. In this article I explore another idea: using Bezier surfaces instead.
When using correction matrices, we calculate the corrected color of a pixel by multiplying its original color value by a given matrix. It means that the color correction is a linear transformation of the colors in the image. However the color distortion is not necessarily linear and this method may be insufficient. This was a motivation to try to apply the correction after converting to L*a*b* color space for example.
Turning the problem on its head, instead of changing the data to match the linearity of the correction, why not modify the correction to match their non-linearity ? A matrix can't be used anymore, instead I'd like to try Bezier surfaces which I have already implemented in LibCapy as described here. Honestly, I don't know if Bezier surface is the appropriate name: it's a Bezier curve of \(\mathbb{R}^3\) to \(\mathbb{R}^3\). "Multivariate Bezier" ? "Bezier hyper-surface" ? ... If you know, tell me please. Anyway, whatever the name, here is how it works.
We are looking for a function which maps color values in the image \(\vec{I}\in\mathbb{R}^3\) to reference color values \(\vec{R}\in\mathbb{R}^3\). The reference colors chart and its corresponding colors in the image give us a few pairs (\(\vec{I},\vec{R}\)), from which we want to interpolate for any value \(\vec{I}\). Finding the Bezier surface \(\vec{S}\) that fits the best those pairs is easy: we remember that \(\vec{S}(\vec{I})=\vec{W}(\vec{I}).\vec{C}_i=\vec{O}_i,i\in\{0,1,2\}\), where \(\vec{W}()\) and \(\vec{C}\) are respectively the weights and control points of the Bezier surface. Note that even if \(\vec{S}()\) is not a linear function, \(\vec{W}(\vec{I}).\vec{C}_i\) is a linear product of known (\(\vec{W}()\)) and unknown (\(\vec{C}_i\)) variables. This is where using Bezier surface helps us with the non-linearity problem of color correction. Replacing \(\vec{O}\) with the \(\vec{R}\) values from the reference color chart, and calculating \(\vec{W}(\vec{I})\) for the corresponding colors in the image, we obtain a system of linear equations which we can solve to obtain the control points \(\vec{C}\) of \(\vec{S}\). Once we know the control points of the surface, we can calculate the color correction for any other color.
A few things need to be precised. First, the input \(\vec{I}\) should be in the range [0,1]. If we are using sRGB values, either they are already in that range, either they are in [0,255] and a simple division solve the problem. If we are using L*a*b*, values are not strictly bounded but usually within, respectively, [0,100], [-100,100] and [-100,100]. Dividing L* by 100, and applying \(c'=c/200+0.5\) to a* and b* will make us happy.
To solve the linear system of equations, using the pseudo-inverse will do, as long as we choose an appropriate order for the Bezier surface. If the order is equal to 1 the surface is linear and we are almost doing the same thing as we were doing with the correction matrices. So we want a higher order, and possibly the highest possible as it allows, a priori, for better matching the color distortion. We can't however choose an arbitrary high order: the number of control points grows with it, and to keep the system solvable it must stay below the number of equations (which comes from the number of swatches in the reference colors chart). An appropriate order is then: \(o=\lfloor \sqrt[3]{n}\rfloor -1\) where \(n\) is the number of swatches.
To have at least \(o=2\) we need at least \(n=27\). If we are accounting for the reference swatches out of the sRGB gamut in the QP203 color chart, we have only \(n=26\). Not good. And there is another problem: The domain of the Bezier surface is \([0,1]^3\) but most of the known points are well inside that domain (especially if we use brightness matching as explained in my previous article). It means the control points at the boundary of the domain are very loosely constrained, hence the color correction will run wild in the boundary areas. The image below illustrates this (obtained by temporarily using all the reference swatches).
Restraining the color correction to reference colors chart with tons of swatches, including swatches on the domain boundary, is not realistic. Reducing the order makes this solution useless. What can be done to get out of that problem ? I like to see that Bezier surface as a piece of fabric we would like to shape by pushing/pulling it here and there with our finger. Imagine the piece of fabric isn't attached to anything, any attempt to shape it will give nothing but a mess. Now, imagine you attach the four corners of the fabric to put it into tension. Then you have something you can work on.
What does it mean for our Bezier surface ? We need to add artificial constraints to "put it into tension". How should these artificial constraints be created. First I don't want them to interfer too much with the colors to be corrected, i.e. I want them out of the domain. Given that the domain is \([0,1]^3\), I choose arbitrarily to lay them onto the surface defined by \(\{x,y,z\}\) given \(x\in\{-1,2\}\) or \(y\in\{-1,2\}\) or \(z\in\{-1,2\}\).
How many constraints do we add controls how much tension we add: more constraints gives less divergence but less correction, less constraints gives more divergence but more correction. Lets call \(n\ge2\) the number of constraints per axis, with default value equal to 3. The total number of constraints is equal to \(6n^2-12n+8\).
We also want the tension to be equally distributed on the domain. Then the constraints should be located at \(\{3i/(n-1)-1,3j/(n-1)-1,3k/(n-1)-1\}\) with \((i,j,k)\in[0,n-1]^3\) and \((i=0)\lor(i=n-1)\lor(j=0)\lor(j=n-1)\lor(k=0)\lor(k=n-1)\). The image below illustrates how the constraints avoid the divergence at boundaries of the domain.
And we're done with the explanation of how to use Bezier surface to perform color correction. Time to check howit actually performs. I'll reuse the same dataset as the one in my previous article for comparison. Without brightness correction:
| image | initial similarity | final similarity (sRGB) | final similarity (L*a*b*) |
|---|---|---|---|
| #00 | 0.9507 | 0.9588 | 0.9580 |
| #01 | 0.9239 | 0.9330 | 0.9291 |
| #02 | 0.9413 | 0.9490 | 0.9480 |
| #03 | 0.8479 | 0.9225 | 0.9189 |
| #04 | 0.9249 | 0.9324 | 0.9287 |
| #05 | 0.9457 | 0.9588 | 0.9535 |
| #06 | 0.6866 | 0.9238 | 0.8788 |
| #07 | 0.7824 | 0.9467 | 0.9061 |
| #08 | 0.8999 | 0.9553 | 0.9546 |
With brightness correction:
| image | initial similarity | final similarity (sRGB) | final similarity (L*a*b*) |
|---|---|---|---|
| #00 | 0.8979 | 0.9284 | 0.9264 |
| #01 | 0.7666 | 0.9150 | 0.9079 |
| #02 | 0.7287 | 0.9361 | 0.9316 |
| #03 | 0.7429 | 0.9065 | 0.8907 |
| #04 | 0.5236 | 0.9313 | 0.9240 |
| #05 | 0.9007 | 0.9288 | 0.9241 |
| #06 | 0.6593 | 0.9110 | 0.8497 |
| #07 | 0.7720 | 0.9213 | 0.8815 |
| #08 | 0.8742 | 0.9219 | 0.9211 |
The result shows that using Bezier surface does correct the color distortion but less effectively than correction matrix. If applied to sRGB color space the results are generally better than L*a*b* color space. A possible explanation to the worst performance of Bezier surface is that the order is forced to stay low, hence there are few control points to fit the correction to the reference. And these control points are uniformly spread on the domain while the reference points are more concentrated to the center of the domain. This makes it difficult for the Bezier to fit correctly the reference points, hence the lower scores. That's a bit disappointing, but that was worth trying, and I think the idea is worth exploring further with other interpolant.