/ocr/ocrservice/jni/opticalflow/feature_detector.cpp
C++ | 382 lines | 216 code | 75 blank | 91 comment | 30 complexity | 1b03698e37a3ba59df9d247f52e7002f MD5 | raw file
1/* 2 * Copyright 2011, Google Inc. 3 * 4 * Licensed under the Apache License, Version 2.0 (the "License"); 5 * you may not use this file except in compliance with the License. 6 * You may obtain a copy of the License at 7 * 8 * http://www.apache.org/licenses/LICENSE-2.0 9 * 10 * Unless required by applicable law or agreed to in writing, software 11 * distributed under the License is distributed on an "AS IS" BASIS, 12 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. 13 * See the License for the specific language governing permissions and 14 * limitations under the License. 15 */ 16 17// Author: Andrew Harp 18// 19// Various feature detecting functions. 20 21#include <float.h> 22 23#include "utils.h" 24#include "time_log.h" 25 26#include "image.h" 27#include "feature_detector.h" 28 29// Threshold for pixels to be considered different. 30#define FAST_DIFF_AMOUNT 10 31 32// How far from edge of frame to stop looking for FAST features. 33#define FAST_BORDER_BUFFER 20 34 35// Minimum enforced distance between detected features. 36// Default 37#define MIN_FEATURE_DIST_NORMAL 24 38 39// Regions selected as "interesting" (aka have annotations) can have denser 40// feature coverage. 41#define MIN_FEATURE_DIST_INTEREST 12 42 43// How many FAST qualifying pixels must be connected to a pixel for it to be 44// considered a candidate feature for Harris filtering. 45#define MIN_NUM_CONNECTED 8 46 47// Size of the window to integrate over for Harris filtering. 48// Compare to WINDOW_SIZE in optical_flow.h. 49#define HARRIS_WINDOW_SIZE 2 50 51// Arbitrary parameter for how picky Harris filter is, the higher the more 52// discriminating. 53#define SENSITIVITY 0.2f 54 55namespace flow { 56 57void scoreFeatures(const Image<int32>& I_x, const Image<int32>& I_y, 58 const int32 num_candidates, 59 Point2D* const candidate_features) { 60 // Score all the features 61 for (int32 i = 0; i < num_candidates; ++i) { 62 Point2D* const feature = candidate_features + i; 63 feature->score = harrisFilter(I_x, I_y, feature->x, feature->y); 64 } 65} 66 67// Quicksorts detected features by score and then selects them such that 68// they are separated by a minimum distance. 69int32 sortAndSelect(const int32 num_candidates, const int32 max_features, 70 const Image<bool>& interest_map, 71 Point2D* const candidate_features, 72 Point2D* const final_features, 73 Image<uint8>* const best_feature_map) { 74 qsort(candidate_features, num_candidates); 75 76#ifdef SANITY_CHECKS 77 // Verify that the array got sorted. 78 float32 last_score = -FLT_MAX; 79 for (int32 i = 0; i < num_candidates; ++i) { 80 const float32 curr_score = (candidate_features + i)->score; 81 82 // Scores should be monotonically increasing. 83 CHECK(last_score <= curr_score, 84 "Quicksort failure! %d: %.5f > %d: %.5f", 85 i - 1, last_score, i, curr_score); 86 87 last_score = curr_score; 88 } 89#endif 90 91 best_feature_map->clear(false); 92 93 int32 num_features = 0; 94 95 for (int32 i = num_candidates - 1; i >= 0; --i) { 96 const Point2D& candidate = candidate_features[i]; 97 98 // Since features are sorted, the first 0 or less value means we can stop 99 // looking. 100 if (candidate.score <= 0.0f) { 101 break; 102 } 103 104 // Lookup whether this feature is in an interest region. If so, other 105 // features may appear closer to it than normal. 106 const int32 distance = interest_map.getPixel(candidate.x, candidate.y) ? 107 MIN_FEATURE_DIST_INTEREST : MIN_FEATURE_DIST_NORMAL; 108 109 if (markImage(candidate.x, candidate.y, distance, best_feature_map)) { 110 final_features[num_features] = candidate; 111 num_features++; 112 113 if (num_features >= max_features) { 114 break; 115 } 116 } 117 } 118 119 return num_features; 120} 121 122// Walks along the given circle checking for pixels above or below the center. 123// Returns a score, or 0 if the feature did not pass the criteria. 124// 125// Parameters: 126// circle_perimeter: the circumference in pixels of the circle. 127// threshold: the minimum number of contiguous pixels that must be above or 128// below the center value. 129// center_ptr: the location of the center pixel in memory 130// offsets: the relative offsets from the center pixel of the edge pixels. 131inline int32 testCircle(const int32 circle_perimeter, const int32 threshold, 132 const uint8* const center_ptr, 133 const int32* offsets) { 134 // Get the actual value of the center pixel for easier reference later on. 135 const int32 center_value = static_cast<int32>(*center_ptr); 136 137 // Number of total pixels to check. Have to wrap around some in case 138 // the contiguous section is split by the array edges. 139 const int32 num_total = circle_perimeter + threshold - 1; 140 141 int32 num_above = 0; 142 int32 above_diff = 0; 143 144 int32 num_below = 0; 145 int32 below_diff = 0; 146 147 // Used to tell when this is definitely not going to meet the threshold so we 148 // can early abort. 149 int32 minimum_by_now = threshold - num_total + 1; 150 151 // Go through every pixel along the perimeter of the circle, and then around 152 // again a little bit. 153 for (int32 i = 0; i < num_total; ++i) { 154 // This should be faster than mod. 155 const int32 perim_index = i < circle_perimeter ? i : i - circle_perimeter; 156 157 // This gets the value of the current pixel along the perimeter by using 158 // a precomputed offset. 159 const int32 curr_value = 160 static_cast<int32>(center_ptr[offsets[perim_index]]); 161 162 const int32 difference = curr_value - center_value; 163 164 if (difference > FAST_DIFF_AMOUNT) { 165 above_diff += difference; 166 ++num_above; 167 168 num_below = 0; 169 below_diff = 0; 170 171 if (num_above >= threshold) { 172 return above_diff; 173 } 174 } else if (difference < -FAST_DIFF_AMOUNT) { 175 below_diff += difference; 176 ++num_below; 177 178 num_above = 0; 179 above_diff = 0; 180 181 if (num_below >= threshold) { 182 return below_diff; 183 } 184 } else { 185 num_above = 0; 186 num_below = 0; 187 above_diff = 0; 188 below_diff = 0; 189 } 190 191 // See if there's any chance of making the threshold. 192 if (max(num_above, num_below) < minimum_by_now) { 193 // Didn't pass. 194 return 0; 195 } 196 ++minimum_by_now; 197 } 198 199 // Didn't pass. 200 return 0; 201} 202 203// Creates features in a regular grid, regardless of image contents. 204int32 seedFeatures(const Image<uint8>& frame, 205 const int32 num_x, const int32 num_y, 206 const float32 left, const float32 top, 207 const float32 right, const float32 bottom, 208 const int32 type, Point2D* const features) { 209 int32 num_features = 0; 210 211 const float32 step_x = ((right - left) / (num_x - 1)); 212 const float32 step_y = ((bottom - top) / (num_y - 1)); 213 214 for (int32 x = 0; x < num_x; ++x) { 215 for (int32 y = 0; y < num_y; ++y) { 216 const int32 x_pos = x * step_x + left; 217 const int32 y_pos = y * step_y + top; 218 if (inRange(x_pos, 0, frame.width_less_one_) && 219 inRange(y_pos, 0, frame.height_less_one_)) { 220 Point2D* const feature = features + num_features; 221 feature->x = x_pos; 222 feature->y = y_pos; 223 feature->type = type; 224 225 ++num_features; 226 } 227 } 228 } 229 230 return num_features; 231} 232 233 234// Returns how likely a point in the image is to be a corner. 235float32 harrisFilter(const Image<int32>& I_x, const Image<int32>& I_y, 236 const int32 x, const int32 y) { 237 // Image gradient matrix. 238 float32 G[] = { 0, 0, 0, 0 }; 239 calculateG(HARRIS_WINDOW_SIZE, x, y, I_x, I_y, G); 240 241 const float32 g_sum = G[0] + G[1] + G[2] + G[3]; 242 243 const float32 a = G[0] / g_sum; 244 const float32 b = G[1] / g_sum; 245 const float32 c = G[2] / g_sum; 246 const float32 d = G[3] / g_sum; 247 248 const float32 det = a * d - b * c; 249 const float32 trace = a + d; 250 251 const float32 inner = square(trace) - 4 * det; 252 253 if (inner >= 0.0f) { 254 const float32 square_root_inner = sqrtf(inner); 255 const float32 eig1 = (trace + square_root_inner) / 2.0f; 256 const float32 eig2 = (trace - square_root_inner) / 2.0f; 257 return eig1 * eig2 - SENSITIVITY * square(eig1 + eig2); 258 } 259 260 // Way negative. 261 return -100.0f; 262} 263 264// FAST feature detector. 265int32 findFastFeatures(const Image<uint8>& frame, const int32 max_num_features, 266 Point2D* const features, 267 Image<uint8>* const best_feature_map) { 268 /* 269 // Reference for a circle of diameter 7. 270 const int32 circle[] = {0, 0, 1, 1, 1, 0, 0, 271 0, 1, 0, 0, 0, 1, 0, 272 1, 0, 0, 0, 0, 0, 1, 273 1, 0, 0, 0, 0, 0, 1, 274 1, 0, 0, 0, 0, 0, 1, 275 0, 1, 0, 0, 0, 1, 0, 276 0, 0, 1, 1, 1, 0, 0}; 277 const int32 circle_offset[] = 278 {2, 3, 4, 8, 12, 14, 20, 21, 27, 28, 34, 36, 40, 44, 45, 46}; 279 */ 280 281 // Quick test of compass directions. Any length 16 circle with a break of up 282 // to 4 pixels will have at least 3 of these 4 pixels active. 283 static const int32 short_circle_perimeter = 4; 284 static const int32 short_threshold = 3; 285 static const int32 short_circle_x[] = { -3, 0, +3, 0 }; 286 static const int32 short_circle_y[] = { 0, -3, 0, +3 }; 287 288 // Precompute image offsets. 289 int32 short_offsets[short_circle_perimeter]; 290 for (int i = 0; i < short_circle_perimeter; ++i) { 291 short_offsets[i] = short_circle_x[i] + short_circle_y[i] * frame.getWidth(); 292 } 293 294 // Large circle values. 295 static const int32 full_circle_perimeter = 16; 296 static const int32 full_threshold = 12; 297 static const int32 full_circle_x[] = 298 { -1, 0, +1, +2, +3, +3, +3, +2, +1, +0, -1, -2, -3, -3, -3, -2 }; 299 static const int32 full_circle_y[] = 300 { -3, -3, -3, -2, -1, 0, +1, +2, +3, +3, +3, +2, +1, +0, -1, -2 }; 301 302 // Precompute image offsets. 303 int32 full_offsets[full_circle_perimeter]; 304 for (int i = 0; i < full_circle_perimeter; ++i) { 305 full_offsets[i] = full_circle_x[i] + full_circle_y[i] * frame.getWidth(); 306 } 307 308 const int frame_width = frame.getWidth(); 309 310 const int end_y = frame.getHeight() - FAST_BORDER_BUFFER; 311 const int end_x = frame.getWidth() - FAST_BORDER_BUFFER; 312 313 best_feature_map->clear(0); 314 315 // Loop through once to find FAST feature clumps. 316 for (int32 img_y = FAST_BORDER_BUFFER; img_y < end_y; ++img_y) { 317 const uint8* curr_pixel_ptr = 318 frame.getPixelPtrConst(FAST_BORDER_BUFFER, img_y); 319 320 for (int32 img_x = FAST_BORDER_BUFFER; img_x < end_x; ++img_x) { 321 // Only insert it if it meets the quick minimum requirements test. 322 if (testCircle(short_circle_perimeter, short_threshold, 323 curr_pixel_ptr, short_offsets) != 0) { 324 325 // Longer test for actual feature score.. 326 const int32 fast_score = testCircle(full_circle_perimeter, 327 full_threshold, 328 curr_pixel_ptr, 329 full_offsets); 330 331 // Non-zero score means the feature was found. 332 if (fast_score != 0) { 333 uint8* const center_ptr = best_feature_map->getPixelPtr(img_x, img_y); 334 335 // Increase the feature count on this pixel and the pixels in all 336 // 4 cardinal directions. 337 *center_ptr += 5; 338 *(center_ptr - 1) += 1; 339 *(center_ptr + 1) += 1; 340 *(center_ptr - frame_width) += 1; 341 *(center_ptr + frame_width) += 1; 342 } 343 } 344 345 ++curr_pixel_ptr; 346 } // x 347 } // y 348 349 timeLog("Found FAST features"); 350 351 int32 num_features = 0; 352 // Loop through again and Harris filter pixels in the center of clumps. 353 // We can shrink the window by 1 pixel on every side. 354 for (int32 img_y = FAST_BORDER_BUFFER + 1; img_y < end_y - 1; ++img_y) { 355 const int32 start_x = FAST_BORDER_BUFFER + 1; 356 357 const uint8* curr_pixel_ptr = 358 best_feature_map->getPixelPtrConst(start_x, img_y); 359 360 for (int32 img_x = start_x; img_x < end_x - 1; ++img_x) { 361 if (*curr_pixel_ptr >= MIN_NUM_CONNECTED) { 362 Point2D* const feature = features + num_features; 363 feature->x = img_x; 364 feature->y = img_y; 365 feature->score = 0; 366 feature->type = FEATURE_FAST; 367 368 ++num_features; 369 if (num_features >= max_num_features) { 370 return num_features; 371 } 372 } 373 374 ++curr_pixel_ptr; 375 } // x 376 } // y 377 378 timeLog("Filtered FAST features"); 379 return num_features; 380} 381 382} // namespace flow