1 /*
2 * Copyright 2016 Google Inc.
3 *
4 * Use of this source code is governed by a BSD-style license that can be
5 * found in the LICENSE file.
6 */
7
8 #include "SkColorLookUpTable.h"
9 #include "SkFloatingPoint.h"
10
interp(float * dst,const float * src) const11 void SkColorLookUpTable::interp(float* dst, const float* src) const {
12 if (fInputChannels == 3) {
13 interp3D(dst, src);
14 } else {
15 SkASSERT(dst != src);
16 // index gets initialized as the algorithm proceeds by interpDimension.
17 // It's just there to store the choice of low/high so far.
18 int index[kMaxColorChannels];
19 for (uint8_t outputDimension = 0; outputDimension < kOutputChannels; ++outputDimension) {
20 dst[outputDimension] = interpDimension(src, fInputChannels - 1, outputDimension,
21 index);
22 }
23 }
24 }
25
interp3D(float * dst,const float * src) const26 void SkColorLookUpTable::interp3D(float* dst, const float* src) const {
27 SkASSERT(3 == kOutputChannels);
28 // Call the src components x, y, and z.
29 const uint8_t maxX = fGridPoints[0] - 1;
30 const uint8_t maxY = fGridPoints[1] - 1;
31 const uint8_t maxZ = fGridPoints[2] - 1;
32
33 // An approximate index into each of the three dimensions of the table.
34 const float x = src[0] * maxX;
35 const float y = src[1] * maxY;
36 const float z = src[2] * maxZ;
37
38 // This gives us the low index for our interpolation.
39 int ix = sk_float_floor2int(x);
40 int iy = sk_float_floor2int(y);
41 int iz = sk_float_floor2int(z);
42
43 // Make sure the low index is not also the max index.
44 ix = (maxX == ix) ? ix - 1 : ix;
45 iy = (maxY == iy) ? iy - 1 : iy;
46 iz = (maxZ == iz) ? iz - 1 : iz;
47
48 // Weighting factors for the interpolation.
49 const float diffX = x - ix;
50 const float diffY = y - iy;
51 const float diffZ = z - iz;
52
53 // Constants to help us navigate the 3D table.
54 // Ex: Assume x = a, y = b, z = c.
55 // table[a * n001 + b * n010 + c * n100] logically equals table[a][b][c].
56 const int n000 = 0;
57 const int n001 = 3 * fGridPoints[1] * fGridPoints[2];
58 const int n010 = 3 * fGridPoints[2];
59 const int n011 = n001 + n010;
60 const int n100 = 3;
61 const int n101 = n100 + n001;
62 const int n110 = n100 + n010;
63 const int n111 = n110 + n001;
64
65 // Base ptr into the table.
66 const float* ptr = &(table()[ix*n001 + iy*n010 + iz*n100]);
67
68 // The code below performs a tetrahedral interpolation for each of the three
69 // dst components. Once the tetrahedron containing the interpolation point is
70 // identified, the interpolation is a weighted sum of grid values at the
71 // vertices of the tetrahedron. The claim is that tetrahedral interpolation
72 // provides a more accurate color conversion.
73 // blogs.mathworks.com/steve/2006/11/24/tetrahedral-interpolation-for-colorspace-conversion/
74 //
75 // I have one test image, and visually I can't tell the difference between
76 // tetrahedral and trilinear interpolation. In terms of computation, the
77 // tetrahedral code requires more branches but less computation. The
78 // SampleICC library provides an option for the client to choose either
79 // tetrahedral or trilinear.
80 for (int i = 0; i < 3; i++) {
81 if (diffZ < diffY) {
82 if (diffZ > diffX) {
83 dst[i] = (ptr[n000] + diffZ * (ptr[n110] - ptr[n010]) +
84 diffY * (ptr[n010] - ptr[n000]) +
85 diffX * (ptr[n111] - ptr[n110]));
86 } else if (diffY < diffX) {
87 dst[i] = (ptr[n000] + diffZ * (ptr[n111] - ptr[n011]) +
88 diffY * (ptr[n011] - ptr[n001]) +
89 diffX * (ptr[n001] - ptr[n000]));
90 } else {
91 dst[i] = (ptr[n000] + diffZ * (ptr[n111] - ptr[n011]) +
92 diffY * (ptr[n010] - ptr[n000]) +
93 diffX * (ptr[n011] - ptr[n010]));
94 }
95 } else {
96 if (diffZ < diffX) {
97 dst[i] = (ptr[n000] + diffZ * (ptr[n101] - ptr[n001]) +
98 diffY * (ptr[n111] - ptr[n101]) +
99 diffX * (ptr[n001] - ptr[n000]));
100 } else if (diffY < diffX) {
101 dst[i] = (ptr[n000] + diffZ * (ptr[n100] - ptr[n000]) +
102 diffY * (ptr[n111] - ptr[n101]) +
103 diffX * (ptr[n101] - ptr[n100]));
104 } else {
105 dst[i] = (ptr[n000] + diffZ * (ptr[n100] - ptr[n000]) +
106 diffY * (ptr[n110] - ptr[n100]) +
107 diffX * (ptr[n111] - ptr[n110]));
108 }
109 }
110
111 // |src| is guaranteed to be in the 0-1 range as are all entries
112 // in the table. For "increasing" tables, outputs will also be
113 // in the 0-1 range. While this property is logical for color
114 // look up tables, we don't check for it.
115 // And for arbitrary, non-increasing tables, it is easy to see how
116 // the output might not be 0-1. So we clamp here.
117 if (dst[i] > 1.f) {
118 dst[i] = 1.f;
119 } else if (dst[i] < 0.f) {
120 dst[i] = 0.f;
121 }
122
123 // Increment the table ptr in order to handle the next component.
124 // Note that this is the how table is designed: all of nXXX
125 // variables are multiples of 3 because there are 3 output
126 // components.
127 ptr++;
128 }
129 }
130
interpDimension(const float * src,int inputDimension,int outputDimension,int index[kMaxColorChannels]) const131 float SkColorLookUpTable::interpDimension(const float* src, int inputDimension,
132 int outputDimension,
133 int index[kMaxColorChannels]) const {
134 // Base case. We've already decided whether to use the low or high point for each dimension
135 // which is stored inside of index[] where index[i] gives the point in the CLUT to use for
136 // input dimension i.
137 if (inputDimension < 0) {
138 // compute index into CLUT and look up the colour
139 int outputIndex = outputDimension;
140 int indexMultiplier = kOutputChannels;
141 for (int i = fInputChannels - 1; i >= 0; --i) {
142 outputIndex += index[i] * indexMultiplier;
143 indexMultiplier *= fGridPoints[i];
144 }
145 return table()[outputIndex];
146 }
147 // for each dimension (input channel), try both the low and high point for it
148 // and then do the same recursively for the later dimensions.
149 // Finally, we need to LERP the results. ie LERP X then LERP Y then LERP Z.
150 const float x = src[inputDimension] * (fGridPoints[inputDimension] - 1);
151 // try the low point for this dimension
152 index[inputDimension] = sk_float_floor2int(x);
153 const float diff = x - index[inputDimension];
154 // and recursively LERP all sub-dimensions with the current dimension fixed to the low point
155 const float lo = interpDimension(src, inputDimension - 1, outputDimension, index);
156 // now try the high point for this dimension
157 index[inputDimension] = sk_float_ceil2int(x);
158 // and recursively LERP all sub-dimensions with the current dimension fixed to the high point
159 const float hi = interpDimension(src, inputDimension - 1, outputDimension, index);
160 // then LERP the results based on the current dimension
161 return (1 - diff) * lo + diff * hi;
162 }
163