1 // Copyright (c) 2013 The Chromium Authors. All rights reserved. 2 // Use of this source code is governed by a BSD-style license that can be 3 // found in the LICENSE file. 4 5 #include <functional> 6 #include <numeric> 7 #include <vector> 8 9 #include "base/basictypes.h" 10 #include "base/file_util.h" 11 #include "base/files/file_path.h" 12 #include "base/logging.h" 13 #include "base/time/time.h" 14 #include "skia/ext/convolver.h" 15 #include "skia/ext/recursive_gaussian_convolution.h" 16 #include "testing/gtest/include/gtest/gtest.h" 17 #include "third_party/skia/include/core/SkPoint.h" 18 #include "third_party/skia/include/core/SkRect.h" 19 20 namespace { 21 22 int ComputeRowStride(int width, int channel_count, int stride_slack) { 23 return width * channel_count + stride_slack; 24 } 25 26 SkIPoint MakeImpulseImage(std::vector<unsigned char>* image, 27 int width, 28 int height, 29 int channel_index, 30 int channel_count, 31 int stride_slack) { 32 const int src_row_stride = ComputeRowStride( 33 width, channel_count, stride_slack); 34 const int src_byte_count = src_row_stride * height; 35 const int signal_x = width / 2; 36 const int signal_y = height / 2; 37 38 image->resize(src_byte_count, 0); 39 const int non_zero_pixel_index = 40 signal_y * src_row_stride + signal_x * channel_count + channel_index; 41 (*image)[non_zero_pixel_index] = 255; 42 return SkIPoint::Make(signal_x, signal_y); 43 } 44 45 SkIRect MakeBoxImage(std::vector<unsigned char>* image, 46 int width, 47 int height, 48 int channel_index, 49 int channel_count, 50 int stride_slack, 51 int box_width, 52 int box_height, 53 unsigned char value) { 54 const int src_row_stride = ComputeRowStride( 55 width, channel_count, stride_slack); 56 const int src_byte_count = src_row_stride * height; 57 const SkIRect box = SkIRect::MakeXYWH((width - box_width) / 2, 58 (height - box_height) / 2, 59 box_width, box_height); 60 61 image->resize(src_byte_count, 0); 62 for (int y = box.top(); y < box.bottom(); ++y) { 63 for (int x = box.left(); x < box.right(); ++x) 64 (*image)[y * src_row_stride + x * channel_count + channel_index] = value; 65 } 66 67 return box; 68 } 69 70 int ComputeBoxSum(const std::vector<unsigned char>& image, 71 const SkIRect& box, 72 int image_width) { 73 // Compute the sum of all pixels in the box. Assume byte stride 1 and row 74 // stride same as image_width. 75 int sum = 0; 76 for (int y = box.top(); y < box.bottom(); ++y) { 77 for (int x = box.left(); x < box.right(); ++x) 78 sum += image[y * image_width + x]; 79 } 80 81 return sum; 82 } 83 84 } // namespace 85 86 namespace skia { 87 88 TEST(RecursiveGaussian, SmoothingMethodComparison) { 89 static const int kImgWidth = 512; 90 static const int kImgHeight = 220; 91 static const int kChannelIndex = 3; 92 static const int kChannelCount = 3; 93 static const int kStrideSlack = 22; 94 95 std::vector<unsigned char> input; 96 SkISize image_size = SkISize::Make(kImgWidth, kImgHeight); 97 MakeImpulseImage( 98 &input, kImgWidth, kImgHeight, kChannelIndex, kChannelCount, 99 kStrideSlack); 100 101 // Destination will be a single channel image with stide matching width. 102 const int dest_row_stride = kImgWidth; 103 const int dest_byte_count = dest_row_stride * kImgHeight; 104 std::vector<unsigned char> intermediate(dest_byte_count); 105 std::vector<unsigned char> intermediate2(dest_byte_count); 106 std::vector<unsigned char> control(dest_byte_count); 107 std::vector<unsigned char> output(dest_byte_count); 108 109 const int src_row_stride = ComputeRowStride( 110 kImgWidth, kChannelCount, kStrideSlack); 111 112 const float kernel_sigma = 2.5f; 113 ConvolutionFilter1D filter; 114 SetUpGaussianConvolutionKernel(&filter, kernel_sigma, false); 115 // Process the control image. 116 SingleChannelConvolveX1D(&input[0], src_row_stride, 117 kChannelIndex, kChannelCount, 118 filter, image_size, 119 &intermediate[0], dest_row_stride, 0, 1, false); 120 SingleChannelConvolveY1D(&intermediate[0], dest_row_stride, 0, 1, 121 filter, image_size, 122 &control[0], dest_row_stride, 0, 1, false); 123 124 // Now try the same using the other method. 125 RecursiveFilter recursive_filter(kernel_sigma, RecursiveFilter::FUNCTION); 126 SingleChannelRecursiveGaussianY(&input[0], src_row_stride, 127 kChannelIndex, kChannelCount, 128 recursive_filter, image_size, 129 &intermediate2[0], dest_row_stride, 130 0, 1, false); 131 SingleChannelRecursiveGaussianX(&intermediate2[0], dest_row_stride, 0, 1, 132 recursive_filter, image_size, 133 &output[0], dest_row_stride, 0, 1, false); 134 135 // We cannot expect the results to be really the same. In particular, 136 // the standard implementation is computed in completely fixed-point, while 137 // recursive is done in floating point and squeezed back into char*. On top 138 // of that, its characteristics are a bit different (consult the paper). 139 EXPECT_NEAR(std::accumulate(intermediate.begin(), intermediate.end(), 0), 140 std::accumulate(intermediate2.begin(), intermediate2.end(), 0), 141 50); 142 int difference_count = 0; 143 std::vector<unsigned char>::const_iterator i1, i2; 144 for (i1 = control.begin(), i2 = output.begin(); 145 i1 != control.end(); ++i1, ++i2) { 146 if ((*i1 != 0) != (*i2 != 0)) 147 difference_count++; 148 } 149 150 EXPECT_LE(difference_count, 44); // 44 is 2 * PI * r (r == 7, spot size). 151 } 152 153 TEST(RecursiveGaussian, SmoothingImpulse) { 154 static const int kImgWidth = 200; 155 static const int kImgHeight = 300; 156 static const int kChannelIndex = 3; 157 static const int kChannelCount = 3; 158 static const int kStrideSlack = 22; 159 160 std::vector<unsigned char> input; 161 SkISize image_size = SkISize::Make(kImgWidth, kImgHeight); 162 const SkIPoint centre_point = MakeImpulseImage( 163 &input, kImgWidth, kImgHeight, kChannelIndex, kChannelCount, 164 kStrideSlack); 165 166 // Destination will be a single channel image with stide matching width. 167 const int dest_row_stride = kImgWidth; 168 const int dest_byte_count = dest_row_stride * kImgHeight; 169 std::vector<unsigned char> intermediate(dest_byte_count); 170 std::vector<unsigned char> output(dest_byte_count); 171 172 const int src_row_stride = ComputeRowStride( 173 kImgWidth, kChannelCount, kStrideSlack); 174 175 const float kernel_sigma = 5.0f; 176 RecursiveFilter recursive_filter(kernel_sigma, RecursiveFilter::FUNCTION); 177 SingleChannelRecursiveGaussianY(&input[0], src_row_stride, 178 kChannelIndex, kChannelCount, 179 recursive_filter, image_size, 180 &intermediate[0], dest_row_stride, 181 0, 1, false); 182 SingleChannelRecursiveGaussianX(&intermediate[0], dest_row_stride, 0, 1, 183 recursive_filter, image_size, 184 &output[0], dest_row_stride, 0, 1, false); 185 186 // Check we got the expected impulse response. 187 const int cx = centre_point.x(); 188 const int cy = centre_point.y(); 189 unsigned char value_x = output[dest_row_stride * cy + cx]; 190 unsigned char value_y = value_x; 191 EXPECT_GT(value_x, 0); 192 for (int offset = 0; 193 offset < std::min(kImgWidth, kImgHeight) && (value_y > 0 || value_x > 0); 194 ++offset) { 195 // Symmetricity and monotonicity along X. 196 EXPECT_EQ(output[dest_row_stride * cy + cx - offset], 197 output[dest_row_stride * cy + cx + offset]); 198 EXPECT_LE(output[dest_row_stride * cy + cx - offset], value_x); 199 value_x = output[dest_row_stride * cy + cx - offset]; 200 201 // Symmetricity and monotonicity along Y. 202 EXPECT_EQ(output[dest_row_stride * (cy - offset) + cx], 203 output[dest_row_stride * (cy + offset) + cx]); 204 EXPECT_LE(output[dest_row_stride * (cy - offset) + cx], value_y); 205 value_y = output[dest_row_stride * (cy - offset) + cx]; 206 207 // Symmetricity along X/Y (not really assured, but should be close). 208 EXPECT_NEAR(value_x, value_y, 1); 209 } 210 211 // Smooth the inverse now. 212 std::vector<unsigned char> output2(dest_byte_count); 213 std::transform(input.begin(), input.end(), input.begin(), 214 std::bind1st(std::minus<unsigned char>(), 255U)); 215 SingleChannelRecursiveGaussianY(&input[0], src_row_stride, 216 kChannelIndex, kChannelCount, 217 recursive_filter, image_size, 218 &intermediate[0], dest_row_stride, 219 0, 1, false); 220 SingleChannelRecursiveGaussianX(&intermediate[0], dest_row_stride, 0, 1, 221 recursive_filter, image_size, 222 &output2[0], dest_row_stride, 0, 1, false); 223 // The image should be the reverse of output, but permitting for rounding 224 // we will only claim that wherever output is 0, output2 should be 255. 225 // There still can be differences at the edges of the object. 226 std::vector<unsigned char>::const_iterator i1, i2; 227 int difference_count = 0; 228 for (i1 = output.begin(), i2 = output2.begin(); 229 i1 != output.end(); ++i1, ++i2) { 230 // The line below checks (*i1 == 0 <==> *i2 == 255). 231 if ((*i1 != 0 && *i2 == 255) && ! (*i1 == 0 && *i2 != 255)) 232 ++difference_count; 233 } 234 EXPECT_LE(difference_count, 8); 235 } 236 237 TEST(RecursiveGaussian, FirstDerivative) { 238 static const int kImgWidth = 512; 239 static const int kImgHeight = 1024; 240 static const int kChannelIndex = 2; 241 static const int kChannelCount = 4; 242 static const int kStrideSlack = 22; 243 static const int kBoxSize = 400; 244 245 std::vector<unsigned char> input; 246 const SkISize image_size = SkISize::Make(kImgWidth, kImgHeight); 247 const SkIRect box = MakeBoxImage( 248 &input, kImgWidth, kImgHeight, kChannelIndex, kChannelCount, 249 kStrideSlack, kBoxSize, kBoxSize, 200); 250 251 // Destination will be a single channel image with stide matching width. 252 const int dest_row_stride = kImgWidth; 253 const int dest_byte_count = dest_row_stride * kImgHeight; 254 std::vector<unsigned char> output_x(dest_byte_count); 255 std::vector<unsigned char> output_y(dest_byte_count); 256 std::vector<unsigned char> output(dest_byte_count); 257 258 const int src_row_stride = ComputeRowStride( 259 kImgWidth, kChannelCount, kStrideSlack); 260 261 const float kernel_sigma = 3.0f; 262 const int spread = 4 * kernel_sigma; 263 RecursiveFilter recursive_filter(kernel_sigma, 264 RecursiveFilter::FIRST_DERIVATIVE); 265 SingleChannelRecursiveGaussianX(&input[0], src_row_stride, 266 kChannelIndex, kChannelCount, 267 recursive_filter, image_size, 268 &output_x[0], dest_row_stride, 269 0, 1, true); 270 SingleChannelRecursiveGaussianY(&input[0], src_row_stride, 271 kChannelIndex, kChannelCount, 272 recursive_filter, image_size, 273 &output_y[0], dest_row_stride, 274 0, 1, true); 275 276 // In test code we can assume adding the two up should do fine. 277 std::vector<unsigned char>::const_iterator ix, iy; 278 std::vector<unsigned char>::iterator target; 279 for (target = output.begin(), ix = output_x.begin(), iy = output_y.begin(); 280 target < output.end(); ++target, ++ix, ++iy) { 281 *target = *ix + *iy; 282 } 283 284 SkIRect inflated_rect(box); 285 inflated_rect.outset(spread, spread); 286 SkIRect deflated_rect(box); 287 deflated_rect.inset(spread, spread); 288 289 int image_total = ComputeBoxSum(output, 290 SkIRect::MakeWH(kImgWidth, kImgHeight), 291 kImgWidth); 292 int box_inflated = ComputeBoxSum(output, inflated_rect, kImgWidth); 293 int box_deflated = ComputeBoxSum(output, deflated_rect, kImgWidth); 294 EXPECT_EQ(box_deflated, 0); 295 EXPECT_EQ(image_total, box_inflated); 296 297 // Try inverted image. Behaviour should be very similar (modulo rounding). 298 std::transform(input.begin(), input.end(), input.begin(), 299 std::bind1st(std::minus<unsigned char>(), 255U)); 300 SingleChannelRecursiveGaussianX(&input[0], src_row_stride, 301 kChannelIndex, kChannelCount, 302 recursive_filter, image_size, 303 &output_x[0], dest_row_stride, 304 0, 1, true); 305 SingleChannelRecursiveGaussianY(&input[0], src_row_stride, 306 kChannelIndex, kChannelCount, 307 recursive_filter, image_size, 308 &output_y[0], dest_row_stride, 309 0, 1, true); 310 311 for (target = output.begin(), ix = output_x.begin(), iy = output_y.begin(); 312 target < output.end(); ++target, ++ix, ++iy) { 313 *target = *ix + *iy; 314 } 315 316 image_total = ComputeBoxSum(output, 317 SkIRect::MakeWH(kImgWidth, kImgHeight), 318 kImgWidth); 319 box_inflated = ComputeBoxSum(output, inflated_rect, kImgWidth); 320 box_deflated = ComputeBoxSum(output, deflated_rect, kImgWidth); 321 322 EXPECT_EQ(box_deflated, 0); 323 EXPECT_EQ(image_total, box_inflated); 324 } 325 326 TEST(RecursiveGaussian, SecondDerivative) { 327 static const int kImgWidth = 700; 328 static const int kImgHeight = 500; 329 static const int kChannelIndex = 0; 330 static const int kChannelCount = 2; 331 static const int kStrideSlack = 42; 332 static const int kBoxSize = 200; 333 334 std::vector<unsigned char> input; 335 SkISize image_size = SkISize::Make(kImgWidth, kImgHeight); 336 const SkIRect box = MakeBoxImage( 337 &input, kImgWidth, kImgHeight, kChannelIndex, kChannelCount, 338 kStrideSlack, kBoxSize, kBoxSize, 200); 339 340 // Destination will be a single channel image with stide matching width. 341 const int dest_row_stride = kImgWidth; 342 const int dest_byte_count = dest_row_stride * kImgHeight; 343 std::vector<unsigned char> output_x(dest_byte_count); 344 std::vector<unsigned char> output_y(dest_byte_count); 345 std::vector<unsigned char> output(dest_byte_count); 346 347 const int src_row_stride = ComputeRowStride( 348 kImgWidth, kChannelCount, kStrideSlack); 349 350 const float kernel_sigma = 5.0f; 351 const int spread = 8 * kernel_sigma; 352 RecursiveFilter recursive_filter(kernel_sigma, 353 RecursiveFilter::SECOND_DERIVATIVE); 354 SingleChannelRecursiveGaussianX(&input[0], src_row_stride, 355 kChannelIndex, kChannelCount, 356 recursive_filter, image_size, 357 &output_x[0], dest_row_stride, 358 0, 1, true); 359 SingleChannelRecursiveGaussianY(&input[0], src_row_stride, 360 kChannelIndex, kChannelCount, 361 recursive_filter, image_size, 362 &output_y[0], dest_row_stride, 363 0, 1, true); 364 365 // In test code we can assume adding the two up should do fine. 366 std::vector<unsigned char>::const_iterator ix, iy; 367 std::vector<unsigned char>::iterator target; 368 for (target = output.begin(),ix = output_x.begin(), iy = output_y.begin(); 369 target < output.end(); ++target, ++ix, ++iy) { 370 *target = *ix + *iy; 371 } 372 373 int image_total = ComputeBoxSum(output, 374 SkIRect::MakeWH(kImgWidth, kImgHeight), 375 kImgWidth); 376 int box_inflated = ComputeBoxSum(output, 377 SkIRect::MakeLTRB(box.left() - spread, 378 box.top() - spread, 379 box.right() + spread, 380 box.bottom() + spread), 381 kImgWidth); 382 int box_deflated = ComputeBoxSum(output, 383 SkIRect::MakeLTRB(box.left() + spread, 384 box.top() + spread, 385 box.right() - spread, 386 box.bottom() - spread), 387 kImgWidth); 388 // Since second derivative is not really used and implemented mostly 389 // for the sake of completeness, we do not verify the detail (that dip 390 // in the middle). But it is there. 391 EXPECT_EQ(box_deflated, 0); 392 EXPECT_EQ(image_total, box_inflated); 393 } 394 395 } // namespace skia 396
