OLD | NEW |
| (Empty) |
1 #include <cmath> | |
2 | |
3 #include "SkBitmap.h" | |
4 #include "skpdiff_util.h" | |
5 #include "SkPMetric.h" | |
6 #include "SkPMetricUtil_generated.h" | |
7 | |
8 struct RGB { | |
9 float r, g, b; | |
10 }; | |
11 | |
12 struct LAB { | |
13 float l, a, b; | |
14 }; | |
15 | |
16 template<class T> | |
17 struct Image2D { | |
18 int width; | |
19 int height; | |
20 T* image; | |
21 | |
22 Image2D(int w, int h) | |
23 : width(w), | |
24 height(h) { | |
25 SkASSERT(w > 0); | |
26 SkASSERT(h > 0); | |
27 image = SkNEW_ARRAY(T, w * h); | |
28 } | |
29 | |
30 ~Image2D() { | |
31 SkDELETE_ARRAY(image); | |
32 } | |
33 | |
34 void readPixel(int x, int y, T* pixel) const { | |
35 SkASSERT(x >= 0); | |
36 SkASSERT(y >= 0); | |
37 SkASSERT(x < width); | |
38 SkASSERT(y < height); | |
39 *pixel = image[y * width + x]; | |
40 } | |
41 | |
42 T* getRow(int y) const { | |
43 return &image[y * width]; | |
44 } | |
45 | |
46 void writePixel(int x, int y, const T& pixel) { | |
47 SkASSERT(x >= 0); | |
48 SkASSERT(y >= 0); | |
49 SkASSERT(x < width); | |
50 SkASSERT(y < height); | |
51 image[y * width + x] = pixel; | |
52 } | |
53 }; | |
54 | |
55 typedef Image2D<float> ImageL; | |
56 typedef Image2D<RGB> ImageRGB; | |
57 typedef Image2D<LAB> ImageLAB; | |
58 | |
59 template<class T> | |
60 struct ImageArray | |
61 { | |
62 int slices; | |
63 Image2D<T>** image; | |
64 | |
65 ImageArray(int w, int h, int s) | |
66 : slices(s) { | |
67 SkASSERT(s > 0); | |
68 image = SkNEW_ARRAY(Image2D<T>*, s); | |
69 for (int sliceIndex = 0; sliceIndex < slices; sliceIndex++) { | |
70 image[sliceIndex] = SkNEW_ARGS(Image2D<T>, (w, h)); | |
71 } | |
72 } | |
73 | |
74 ~ImageArray() { | |
75 for (int sliceIndex = 0; sliceIndex < slices; sliceIndex++) { | |
76 SkDELETE(image[sliceIndex]); | |
77 } | |
78 SkDELETE_ARRAY(image); | |
79 } | |
80 | |
81 Image2D<T>* getLayer(int z) const { | |
82 SkASSERT(z >= 0); | |
83 SkASSERT(z < slices); | |
84 return image[z]; | |
85 } | |
86 }; | |
87 | |
88 typedef ImageArray<float> ImageL3D; | |
89 | |
90 | |
91 #define MAT_ROW_MULT(rc,gc,bc) r*rc + g*gc + b*bc | |
92 | |
93 | |
94 void adobergb_to_cielab(float r, float g, float b, LAB* lab) { | |
95 // Conversion of Adobe RGB to XYZ taken from from "Adobe RGB (1998) ColorIma
ge Encoding" | |
96 // URL:http://www.adobe.com/digitalimag/pdfs/AdobeRGB1998.pdf | |
97 // Section: 4.3.5.3 | |
98 // See Also: http://en.wikipedia.org/wiki/Adobe_rgb | |
99 float x = MAT_ROW_MULT(0.57667f, 0.18556f, 0.18823f); | |
100 float y = MAT_ROW_MULT(0.29734f, 0.62736f, 0.07529f); | |
101 float z = MAT_ROW_MULT(0.02703f, 0.07069f, 0.99134f); | |
102 | |
103 // The following is the white point in XYZ, so it's simply the row wise addi
tion of the above | |
104 // matrix. | |
105 const float xw = 0.5767f + 0.185556f + 0.188212f; | |
106 const float yw = 0.297361f + 0.627355f + 0.0752847f; | |
107 const float zw = 0.0270328f + 0.0706879f + 0.991248f; | |
108 | |
109 // This is the XYZ color point relative to the white point | |
110 float f[3] = { x / xw, y / yw, z / zw }; | |
111 | |
112 // Conversion from XYZ to LAB taken from | |
113 // http://en.wikipedia.org/wiki/CIELAB#Forward_transformation | |
114 for (int i = 0; i < 3; i++) { | |
115 if (f[i] >= 0.008856f) { | |
116 f[i] = SkPMetricUtil::get_cube_root(f[i]); | |
117 } else { | |
118 f[i] = 7.787f * f[i] + 4.0f / 29.0f; | |
119 } | |
120 } | |
121 lab->l = 116.0f * f[1] - 16.0f; | |
122 lab->a = 500.0f * (f[0] - f[1]); | |
123 lab->b = 200.0f * (f[1] - f[2]); | |
124 } | |
125 | |
126 /// Converts a 8888 bitmap to LAB color space and puts it into the output | |
127 static void bitmap_to_cielab(const SkBitmap* bitmap, ImageLAB* outImageLAB) { | |
128 SkASSERT(bitmap->config() == SkBitmap::kARGB_8888_Config); | |
129 | |
130 int width = bitmap->width(); | |
131 int height = bitmap->height(); | |
132 SkASSERT(outImageLAB->width == width); | |
133 SkASSERT(outImageLAB->height == height); | |
134 | |
135 bitmap->lockPixels(); | |
136 RGB rgb; | |
137 LAB lab; | |
138 for (int y = 0; y < height; y++) { | |
139 unsigned char* row = (unsigned char*)bitmap->getAddr(0, y); | |
140 for (int x = 0; x < width; x++) { | |
141 // Perform gamma correction which is assumed to be 2.2 | |
142 rgb.r = SkPMetricUtil::get_gamma(row[x * 4 + 2]); | |
143 rgb.g = SkPMetricUtil::get_gamma(row[x * 4 + 1]); | |
144 rgb.b = SkPMetricUtil::get_gamma(row[x * 4 + 0]); | |
145 adobergb_to_cielab(rgb.r, rgb.g, rgb.b, &lab); | |
146 outImageLAB->writePixel(x, y, lab); | |
147 } | |
148 } | |
149 bitmap->unlockPixels(); | |
150 } | |
151 | |
152 // From Barten SPIE 1989 | |
153 static float contrast_sensitivity(float cyclesPerDegree, float luminance) { | |
154 float a = 440.0f * powf(1.0f + 0.7f / luminance, -0.2f); | |
155 float b = 0.3f * powf(1 + 100.0 / luminance, 0.15f); | |
156 return a * | |
157 cyclesPerDegree * | |
158 expf(-b * cyclesPerDegree) * | |
159 sqrtf(1.0f + 0.06f * expf(b * cyclesPerDegree)); | |
160 } | |
161 | |
162 #if 0 | |
163 // We're keeping these around for reference and in case the lookup tables are no
longer desired. | |
164 // They are no longer called by any code in this file. | |
165 | |
166 // From Daly 1993 | |
167 static float visual_mask(float contrast) { | |
168 float x = powf(392.498f * contrast, 0.7f); | |
169 x = powf(0.0153f * x, 4.0f); | |
170 return powf(1.0f + x, 0.25f); | |
171 } | |
172 | |
173 // From Ward Larson Siggraph 1997 | |
174 static float threshold_vs_intensity(float adaptationLuminance) { | |
175 float logLum = log10f(adaptationLuminance); | |
176 float x; | |
177 if (logLum < -3.94f) { | |
178 x = -2.86f; | |
179 } else if (logLum < -1.44f) { | |
180 x = powf(0.405f * logLum + 1.6f, 2.18) - 2.86f; | |
181 } else if (logLum < -0.0184f) { | |
182 x = logLum - 0.395f; | |
183 } else if (logLum < 1.9f) { | |
184 x = powf(0.249f * logLum + 0.65f, 2.7f) - 0.72f; | |
185 } else { | |
186 x = logLum - 1.255f; | |
187 } | |
188 return powf(10.0f, x); | |
189 } | |
190 | |
191 #endif | |
192 | |
193 /// Simply takes the L channel from the input and puts it into the output | |
194 static void lab_to_l(const ImageLAB* imageLAB, ImageL* outImageL) { | |
195 for (int y = 0; y < imageLAB->height; y++) { | |
196 for (int x = 0; x < imageLAB->width; x++) { | |
197 LAB lab; | |
198 imageLAB->readPixel(x, y, &lab); | |
199 outImageL->writePixel(x, y, lab.l); | |
200 } | |
201 } | |
202 } | |
203 | |
204 /// Convolves an image with the given filter in one direction and saves it to th
e output image | |
205 static void convolve(const ImageL* imageL, bool vertical, ImageL* outImageL) { | |
206 SkASSERT(imageL->width == outImageL->width); | |
207 SkASSERT(imageL->height == outImageL->height); | |
208 | |
209 const float matrix[] = { 0.05f, 0.25f, 0.4f, 0.25f, 0.05f }; | |
210 const int matrixCount = sizeof(matrix) / sizeof(float); | |
211 const int radius = matrixCount / 2; | |
212 | |
213 // Keep track of what rows are being operated on for quick access. | |
214 float* rowPtrs[matrixCount]; // Because matrixCount is constant, this won't
create a VLA | |
215 for (int y = radius; y < matrixCount; y++) { | |
216 rowPtrs[y] = imageL->getRow(y - radius); | |
217 } | |
218 float* writeRow = outImageL->getRow(0); | |
219 | |
220 for (int y = 0; y < imageL->height; y++) { | |
221 for (int x = 0; x < imageL->width; x++) { | |
222 float lSum = 0.0f; | |
223 for (int xx = -radius; xx <= radius; xx++) { | |
224 int nx = x; | |
225 int ny = y; | |
226 | |
227 // We mirror at edges so that edge pixels that the filter weight
ing still makes | |
228 // sense. | |
229 if (vertical) { | |
230 ny += xx; | |
231 if (ny < 0) { | |
232 ny = -ny; | |
233 } | |
234 if (ny >= imageL->height) { | |
235 ny = imageL->height + (imageL->height - ny - 1); | |
236 } | |
237 } else { | |
238 nx += xx; | |
239 if (nx < 0) { | |
240 nx = -nx; | |
241 } | |
242 if (nx >= imageL->width) { | |
243 nx = imageL->width + (imageL->width - nx - 1); | |
244 } | |
245 } | |
246 | |
247 float weight = matrix[xx + radius]; | |
248 lSum += rowPtrs[ny - y + radius][nx] * weight; | |
249 } | |
250 writeRow[x] = lSum; | |
251 } | |
252 // As we move down, scroll the row pointers down with us | |
253 for (int y = 0; y < matrixCount - 1; y++) | |
254 { | |
255 rowPtrs[y] = rowPtrs[y + 1]; | |
256 } | |
257 rowPtrs[matrixCount - 1] += imageL->width; | |
258 writeRow += imageL->width; | |
259 } | |
260 } | |
261 | |
262 float pmetric(const ImageLAB* baselineLAB, const ImageLAB* testLAB, SkTDArray<Sk
IPoint>* poi) { | |
263 int width = baselineLAB->width; | |
264 int height = baselineLAB->height; | |
265 int maxLevels = (int)log2(width < height ? width : height); | |
266 | |
267 const float fov = M_PI / 180.0f * 45.0f; | |
268 float contrastSensitivityMax = contrast_sensitivity(3.248f, 100.0f); | |
269 float pixelsPerDegree = width / (2.0f * tanf(fov * 0.5f) * 180.0f / M_PI); | |
270 | |
271 ImageL3D baselineL(width, height, maxLevels); | |
272 ImageL3D testL(width, height, maxLevels); | |
273 ImageL scratchImageL(width, height); | |
274 float* cyclesPerDegree = SkNEW_ARRAY(float, maxLevels); | |
275 float* thresholdFactorFrequency = SkNEW_ARRAY(float, maxLevels - 2); | |
276 float* contrast = SkNEW_ARRAY(float, maxLevels - 2); | |
277 | |
278 lab_to_l(baselineLAB, baselineL.getLayer(0)); | |
279 lab_to_l(testLAB, testL.getLayer(0)); | |
280 | |
281 // Compute cpd - Cycles per degree on the pyramid | |
282 cyclesPerDegree[0] = 0.5f * pixelsPerDegree; | |
283 for (int levelIndex = 1; levelIndex < maxLevels; levelIndex++) { | |
284 cyclesPerDegree[levelIndex] = cyclesPerDegree[levelIndex - 1] * 0.5f; | |
285 } | |
286 | |
287 // Contrast sensitivity is based on image dimensions. Therefore it cannot be
statically | |
288 // generated. | |
289 float* contrastSensitivityTable = SkNEW_ARRAY(float, maxLevels * 1000); | |
290 for (int levelIndex = 0; levelIndex < maxLevels; levelIndex++) { | |
291 for (int csLum = 0; csLum < 1000; csLum++) { | |
292 contrastSensitivityTable[levelIndex * 1000 + csLum] = | |
293 contrast_sensitivity(cyclesPerDegree[levelIndex], (float)csLum / 10.0
f + 1e-5f); | |
294 } | |
295 } | |
296 | |
297 // Compute G - The convolved lum for the baseline | |
298 for (int levelIndex = 1; levelIndex < maxLevels; levelIndex++) { | |
299 convolve(baselineL.getLayer(levelIndex - 1), false, &scratchImageL); | |
300 convolve(&scratchImageL, true, baselineL.getLayer(levelIndex)); | |
301 } | |
302 for (int levelIndex = 1; levelIndex < maxLevels; levelIndex++) { | |
303 convolve(testL.getLayer(levelIndex - 1), false, &scratchImageL); | |
304 convolve(&scratchImageL, true, testL.getLayer(levelIndex)); | |
305 } | |
306 | |
307 // Compute F_freq - The elevation f | |
308 for (int levelIndex = 0; levelIndex < maxLevels - 2; levelIndex++) { | |
309 float cpd = cyclesPerDegree[levelIndex]; | |
310 thresholdFactorFrequency[levelIndex] = contrastSensitivityMax / | |
311 contrast_sensitivity(cpd, 100.0f)
; | |
312 } | |
313 | |
314 int failures = 0; | |
315 // Calculate F | |
316 for (int y = 0; y < height; y++) { | |
317 for (int x = 0; x < width; x++) { | |
318 float lBaseline; | |
319 float lTest; | |
320 baselineL.getLayer(0)->readPixel(x, y, &lBaseline); | |
321 testL.getLayer(0)->readPixel(x, y, &lTest); | |
322 | |
323 float avgLBaseline; | |
324 float avgLTest; | |
325 baselineL.getLayer(maxLevels - 1)->readPixel(x, y, &avgLBaseline); | |
326 testL.getLayer(maxLevels - 1)->readPixel(x, y, &avgLTest); | |
327 | |
328 float lAdapt = 0.5f * (avgLBaseline + avgLTest); | |
329 if (lAdapt < 1e-5) { | |
330 lAdapt = 1e-5; | |
331 } | |
332 | |
333 float contrastSum = 0.0f; | |
334 for (int levelIndex = 0; levelIndex < maxLevels - 2; levelIndex++) { | |
335 float baselineL0, baselineL1, baselineL2; | |
336 float testL0, testL1, testL2; | |
337 baselineL.getLayer(levelIndex + 0)->readPixel(x, y, &baselineL0)
; | |
338 testL. getLayer(levelIndex + 0)->readPixel(x, y, &testL0); | |
339 baselineL.getLayer(levelIndex + 1)->readPixel(x, y, &baselineL1)
; | |
340 testL. getLayer(levelIndex + 1)->readPixel(x, y, &testL1); | |
341 baselineL.getLayer(levelIndex + 2)->readPixel(x, y, &baselineL2)
; | |
342 testL. getLayer(levelIndex + 2)->readPixel(x, y, &testL2); | |
343 | |
344 float baselineContrast1 = fabsf(baselineL0 - baselineL1); | |
345 float testContrast1 = fabsf(testL0 - testL1); | |
346 float numerator = (baselineContrast1 > testContrast1) ? | |
347 baselineContrast1 : testContrast1; | |
348 | |
349 float baselineContrast2 = fabsf(baselineL2); | |
350 float testContrast2 = fabsf(testL2); | |
351 float denominator = (baselineContrast2 > testContrast2) ? | |
352 baselineContrast2 : testContrast2; | |
353 | |
354 // Avoid divides by close to zero | |
355 if (denominator < 1e-5) { | |
356 denominator = 1e-5; | |
357 } | |
358 contrast[levelIndex] = numerator / denominator; | |
359 contrastSum += contrast[levelIndex]; | |
360 } | |
361 | |
362 if (contrastSum < 1e-5) { | |
363 contrastSum = 1e-5; | |
364 } | |
365 | |
366 float F = 0.0f; | |
367 for (int levelIndex = 0; levelIndex < maxLevels - 2; levelIndex++) { | |
368 float contrastSensitivity = contrastSensitivityTable[levelIndex
* 1000 + | |
369 (int)(lAdap
t * 10.0)]; | |
370 float mask = SkPMetricUtil::get_visual_mask(contrast[levelIndex]
* | |
371 contrastSensitivity)
; | |
372 | |
373 F += contrast[levelIndex] + | |
374 thresholdFactorFrequency[levelIndex] * mask / contrastSum; | |
375 } | |
376 | |
377 if (F < 1.0f) { | |
378 F = 1.0f; | |
379 } | |
380 | |
381 if (F > 10.0f) { | |
382 F = 10.0f; | |
383 } | |
384 | |
385 | |
386 bool isFailure = false; | |
387 if (fabsf(lBaseline - lTest) > F * SkPMetricUtil::get_threshold_vs_i
ntensity(lAdapt)) { | |
388 isFailure = true; | |
389 } else { | |
390 LAB baselineColor; | |
391 LAB testColor; | |
392 baselineLAB->readPixel(x, y, &baselineColor); | |
393 testLAB->readPixel(x, y, &testColor); | |
394 float contrastA = baselineColor.a - testColor.a; | |
395 float contrastB = baselineColor.b - testColor.b; | |
396 float colorScale = 1.0f; | |
397 if (lAdapt < 10.0f) { | |
398 colorScale = lAdapt / 10.0f; | |
399 } | |
400 colorScale *= colorScale; | |
401 | |
402 if ((contrastA * contrastA + contrastB * contrastB) * colorScale
> F) | |
403 { | |
404 isFailure = true; | |
405 } | |
406 } | |
407 | |
408 if (isFailure) { | |
409 failures++; | |
410 poi->push()->set(x, y); | |
411 } | |
412 } | |
413 } | |
414 | |
415 SkDELETE_ARRAY(cyclesPerDegree); | |
416 SkDELETE_ARRAY(contrast); | |
417 SkDELETE_ARRAY(thresholdFactorFrequency); | |
418 SkDELETE_ARRAY(contrastSensitivityTable); | |
419 return 1.0 - (double)failures / (width * height); | |
420 } | |
421 | |
422 const char* SkPMetric::getName() { | |
423 return "perceptual"; | |
424 } | |
425 | |
426 int SkPMetric::queueDiff(SkBitmap* baseline, SkBitmap* test) { | |
427 double startTime = get_seconds(); | |
428 int diffID = fQueuedDiffs.count(); | |
429 QueuedDiff& diff = fQueuedDiffs.push_back(); | |
430 diff.result = 0.0; | |
431 | |
432 // Ensure the images are comparable | |
433 if (baseline->width() != test->width() || baseline->height() != test->height
() || | |
434 baseline->width() <= 0 || baseline->height() <= 0) { | |
435 diff.finished = true; | |
436 return diffID; | |
437 } | |
438 | |
439 ImageLAB baselineLAB(baseline->width(), baseline->height()); | |
440 ImageLAB testLAB(baseline->width(), baseline->height()); | |
441 | |
442 bitmap_to_cielab(baseline, &baselineLAB); | |
443 bitmap_to_cielab(test, &testLAB); | |
444 | |
445 diff.result = pmetric(&baselineLAB, &testLAB, &diff.poi); | |
446 | |
447 SkDebugf("Time: %f\n", (get_seconds() - startTime)); | |
448 | |
449 return diffID; | |
450 } | |
451 | |
452 | |
453 void SkPMetric::deleteDiff(int id) { | |
454 | |
455 } | |
456 | |
457 bool SkPMetric::isFinished(int id) { | |
458 return fQueuedDiffs[id].finished; | |
459 } | |
460 | |
461 double SkPMetric::getResult(int id) { | |
462 return fQueuedDiffs[id].result; | |
463 } | |
464 | |
465 int SkPMetric::getPointsOfInterestCount(int id) { | |
466 return fQueuedDiffs[id].poi.count(); | |
467 } | |
468 | |
469 SkIPoint* SkPMetric::getPointsOfInterest(int id) { | |
470 return fQueuedDiffs[id].poi.begin(); | |
471 } | |
OLD | NEW |