Sarah Parker | 4dc0f1b | 2016-08-09 17:40:53 -0700 | [diff] [blame] | 1 | /* |
Wan-Teh Chang | 0289a9f | 2022-12-14 09:33:45 -0800 | [diff] [blame] | 2 | * Copyright (c) 2016, Alliance for Open Media. All rights reserved |
Sarah Parker | 4dc0f1b | 2016-08-09 17:40:53 -0700 | [diff] [blame] | 3 | * |
Wan-Teh Chang | 0289a9f | 2022-12-14 09:33:45 -0800 | [diff] [blame] | 4 | * This source code is subject to the terms of the BSD 2 Clause License and |
| 5 | * the Alliance for Open Media Patent License 1.0. If the BSD 2 Clause License |
| 6 | * was not distributed with this source code in the LICENSE file, you can |
| 7 | * obtain it at www.aomedia.org/license/software. If the Alliance for Open |
| 8 | * Media Patent License 1.0 was not distributed with this source code in the |
| 9 | * PATENTS file, you can obtain it at www.aomedia.org/license/patent. |
Sarah Parker | 4dc0f1b | 2016-08-09 17:40:53 -0700 | [diff] [blame] | 10 | */ |
Rachel Barker | 53e3a69 | 2022-11-30 15:45:11 +0000 | [diff] [blame] | 11 | |
Sarah Parker | 4dc0f1b | 2016-08-09 17:40:53 -0700 | [diff] [blame] | 12 | #include <memory.h> |
| 13 | #include <math.h> |
| 14 | #include <time.h> |
| 15 | #include <stdio.h> |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 16 | #include <stdbool.h> |
Rachel Barker | a02796a | 2023-03-02 20:06:47 +0000 | [diff] [blame] | 17 | #include <string.h> |
Sarah Parker | 4dc0f1b | 2016-08-09 17:40:53 -0700 | [diff] [blame] | 18 | #include <assert.h> |
| 19 | |
Rachel Barker | 53e3a69 | 2022-11-30 15:45:11 +0000 | [diff] [blame] | 20 | #include "aom_dsp/flow_estimation/ransac.h" |
Yaowu Xu | 9a0d1b7 | 2021-07-13 09:56:50 -0700 | [diff] [blame] | 21 | #include "aom_dsp/mathutils.h" |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 22 | #include "aom_mem/aom_mem.h" |
Rachel Barker | 53e3a69 | 2022-11-30 15:45:11 +0000 | [diff] [blame] | 23 | |
| 24 | // TODO(rachelbarker): Remove dependence on code in av1/encoder/ |
Alex Converse | 9d068c1 | 2017-08-03 11:48:19 -0700 | [diff] [blame] | 25 | #include "av1/encoder/random.h" |
Sarah Parker | 4dc0f1b | 2016-08-09 17:40:53 -0700 | [diff] [blame] | 26 | |
Sarah Parker | 4dc0f1b | 2016-08-09 17:40:53 -0700 | [diff] [blame] | 27 | #define MAX_MINPTS 4 |
Sarah Parker | 4dc0f1b | 2016-08-09 17:40:53 -0700 | [diff] [blame] | 28 | #define MINPTS_MULTIPLIER 5 |
| 29 | |
Debargha Mukherjee | 3b16e00 | 2018-09-19 13:44:05 -0700 | [diff] [blame] | 30 | #define INLIER_THRESHOLD 1.25 |
Rachel Barker | 3604269 | 2023-01-11 18:09:43 +0000 | [diff] [blame] | 31 | #define INLIER_THRESHOLD_SQUARED (INLIER_THRESHOLD * INLIER_THRESHOLD) |
Rachel Barker | 41f1bb0 | 2023-02-13 16:22:38 +0000 | [diff] [blame] | 32 | #define NUM_TRIALS 20 |
Debargha Mukherjee | 5dfa930 | 2017-02-10 05:00:08 -0800 | [diff] [blame] | 33 | |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 34 | // Flag to enable functions for finding TRANSLATION type models. |
| 35 | // |
| 36 | // These modes are not considered currently due to a spec bug (see comments |
| 37 | // in gm_get_motion_vector() in av1/common/mv.h). Thus we don't need to compile |
| 38 | // the corresponding search functions, but it is nice to keep the source around |
| 39 | // but disabled, for completeness. |
| 40 | #define ALLOW_TRANSLATION_MODELS 0 |
| 41 | |
Sarah Parker | 4dc0f1b | 2016-08-09 17:40:53 -0700 | [diff] [blame] | 42 | //////////////////////////////////////////////////////////////////////////////// |
| 43 | // ransac |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 44 | typedef bool (*IsDegenerateFunc)(double *p); |
| 45 | typedef bool (*FindTransformationFunc)(int points, const double *points1, |
| 46 | const double *points2, double *params); |
| 47 | typedef void (*ProjectPointsFunc)(const double *mat, const double *points, |
| 48 | double *proj, int n, int stride_points, |
| 49 | int stride_proj); |
Sarah Parker | 97fa6da | 2016-09-23 11:17:27 -0700 | [diff] [blame] | 50 | |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 51 | // vtable-like structure which stores all of the information needed by RANSAC |
| 52 | // for a particular model type |
| 53 | typedef struct { |
| 54 | IsDegenerateFunc is_degenerate; |
| 55 | FindTransformationFunc find_transformation; |
| 56 | ProjectPointsFunc project_points; |
| 57 | int minpts; |
| 58 | } RansacModelInfo; |
| 59 | |
| 60 | #if ALLOW_TRANSLATION_MODELS |
| 61 | static void project_points_translation(const double *mat, const double *points, |
| 62 | double *proj, int n, int stride_points, |
| 63 | int stride_proj) { |
Sarah Parker | 97fa6da | 2016-09-23 11:17:27 -0700 | [diff] [blame] | 64 | int i; |
| 65 | for (i = 0; i < n; ++i) { |
| 66 | const double x = *(points++), y = *(points++); |
Debargha Mukherjee | 8db4c77 | 2016-11-07 12:54:21 -0800 | [diff] [blame] | 67 | *(proj++) = x + mat[0]; |
| 68 | *(proj++) = y + mat[1]; |
Sarah Parker | 97fa6da | 2016-09-23 11:17:27 -0700 | [diff] [blame] | 69 | points += stride_points - 2; |
| 70 | proj += stride_proj - 2; |
| 71 | } |
| 72 | } |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 73 | #endif // ALLOW_TRANSLATION_MODELS |
Sarah Parker | 97fa6da | 2016-09-23 11:17:27 -0700 | [diff] [blame] | 74 | |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 75 | static void project_points_affine(const double *mat, const double *points, |
| 76 | double *proj, int n, int stride_points, |
| 77 | int stride_proj) { |
Sarah Parker | 97fa6da | 2016-09-23 11:17:27 -0700 | [diff] [blame] | 78 | int i; |
| 79 | for (i = 0; i < n; ++i) { |
| 80 | const double x = *(points++), y = *(points++); |
Debargha Mukherjee | 8db4c77 | 2016-11-07 12:54:21 -0800 | [diff] [blame] | 81 | *(proj++) = mat[2] * x + mat[3] * y + mat[0]; |
| 82 | *(proj++) = mat[4] * x + mat[5] * y + mat[1]; |
Sarah Parker | 97fa6da | 2016-09-23 11:17:27 -0700 | [diff] [blame] | 83 | points += stride_points - 2; |
| 84 | proj += stride_proj - 2; |
| 85 | } |
| 86 | } |
| 87 | |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 88 | #if ALLOW_TRANSLATION_MODELS |
| 89 | static bool find_translation(int np, const double *pts1, const double *pts2, |
| 90 | double *params) { |
| 91 | double sumx = 0; |
| 92 | double sumy = 0; |
Debargha Mukherjee | afe7c5f | 2017-06-22 13:47:23 -0700 | [diff] [blame] | 93 | |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 94 | for (int i = 0; i < np; ++i) { |
| 95 | double dx = *(pts2++); |
| 96 | double dy = *(pts2++); |
| 97 | double sx = *(pts1++); |
| 98 | double sy = *(pts1++); |
Debargha Mukherjee | e6eb3b5 | 2017-02-26 08:50:56 -0800 | [diff] [blame] | 99 | |
| 100 | sumx += dx - sx; |
| 101 | sumy += dy - sy; |
| 102 | } |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 103 | |
| 104 | params[0] = sumx / np; |
| 105 | params[1] = sumy / np; |
| 106 | params[2] = 1; |
| 107 | params[3] = 0; |
| 108 | params[4] = 0; |
| 109 | params[5] = 1; |
| 110 | return true; |
| 111 | } |
| 112 | #endif // ALLOW_TRANSLATION_MODELS |
| 113 | |
| 114 | static bool find_rotzoom(int np, const double *pts1, const double *pts2, |
| 115 | double *params) { |
| 116 | const int n = 4; // Size of least-squares problem |
| 117 | double mat[4 * 4]; // Accumulator for A'A |
| 118 | double y[4]; // Accumulator for A'b |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 119 | double a[4]; // Single row of A |
| 120 | double b; // Single element of b |
| 121 | |
| 122 | least_squares_init(mat, y, n); |
| 123 | for (int i = 0; i < np; ++i) { |
| 124 | double dx = *(pts2++); |
| 125 | double dy = *(pts2++); |
| 126 | double sx = *(pts1++); |
| 127 | double sy = *(pts1++); |
| 128 | |
Rachel Barker | 24ec16a | 2023-02-09 13:47:35 +0000 | [diff] [blame] | 129 | a[0] = 1; |
| 130 | a[1] = 0; |
| 131 | a[2] = sx; |
| 132 | a[3] = sy; |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 133 | b = dx; |
| 134 | least_squares_accumulate(mat, y, a, b, n); |
| 135 | |
Rachel Barker | 24ec16a | 2023-02-09 13:47:35 +0000 | [diff] [blame] | 136 | a[0] = 0; |
| 137 | a[1] = 1; |
| 138 | a[2] = sy; |
| 139 | a[3] = -sx; |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 140 | b = dy; |
| 141 | least_squares_accumulate(mat, y, a, b, n); |
| 142 | } |
| 143 | |
Rachel Barker | 24ec16a | 2023-02-09 13:47:35 +0000 | [diff] [blame] | 144 | // Fill in params[0] .. params[3] with output model |
| 145 | if (!least_squares_solve(mat, y, params, n)) { |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 146 | return false; |
| 147 | } |
| 148 | |
Rachel Barker | 24ec16a | 2023-02-09 13:47:35 +0000 | [diff] [blame] | 149 | // Fill in remaining parameters |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 150 | params[4] = -params[3]; |
| 151 | params[5] = params[2]; |
| 152 | |
| 153 | return true; |
Debargha Mukherjee | e6eb3b5 | 2017-02-26 08:50:56 -0800 | [diff] [blame] | 154 | } |
| 155 | |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 156 | static bool find_affine(int np, const double *pts1, const double *pts2, |
| 157 | double *params) { |
Rachel Barker | 24ec16a | 2023-02-09 13:47:35 +0000 | [diff] [blame] | 158 | // Note: The least squares problem for affine models is 6-dimensional, |
| 159 | // but it splits into two independent 3-dimensional subproblems. |
| 160 | // Solving these two subproblems separately and recombining at the end |
| 161 | // results in less total computation than solving the 6-dimensional |
| 162 | // problem directly. |
| 163 | // |
| 164 | // The two subproblems correspond to all the parameters which contribute |
| 165 | // to the x output of the model, and all the parameters which contribute |
| 166 | // to the y output, respectively. |
Debargha Mukherjee | e6eb3b5 | 2017-02-26 08:50:56 -0800 | [diff] [blame] | 167 | |
Rachel Barker | 24ec16a | 2023-02-09 13:47:35 +0000 | [diff] [blame] | 168 | const int n = 3; // Size of each least-squares problem |
| 169 | double mat[2][3 * 3]; // Accumulator for A'A |
| 170 | double y[2][3]; // Accumulator for A'b |
| 171 | double x[2][3]; // Output vector |
| 172 | double a[2][3]; // Single row of A |
| 173 | double b[2]; // Single element of b |
| 174 | |
| 175 | least_squares_init(mat[0], y[0], n); |
| 176 | least_squares_init(mat[1], y[1], n); |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 177 | for (int i = 0; i < np; ++i) { |
| 178 | double dx = *(pts2++); |
| 179 | double dy = *(pts2++); |
| 180 | double sx = *(pts1++); |
| 181 | double sy = *(pts1++); |
Debargha Mukherjee | e6eb3b5 | 2017-02-26 08:50:56 -0800 | [diff] [blame] | 182 | |
Rachel Barker | 24ec16a | 2023-02-09 13:47:35 +0000 | [diff] [blame] | 183 | a[0][0] = 1; |
| 184 | a[0][1] = sx; |
| 185 | a[0][2] = sy; |
| 186 | b[0] = dx; |
| 187 | least_squares_accumulate(mat[0], y[0], a[0], b[0], n); |
Debargha Mukherjee | e6eb3b5 | 2017-02-26 08:50:56 -0800 | [diff] [blame] | 188 | |
Rachel Barker | 24ec16a | 2023-02-09 13:47:35 +0000 | [diff] [blame] | 189 | a[1][0] = 1; |
| 190 | a[1][1] = sx; |
| 191 | a[1][2] = sy; |
| 192 | b[1] = dy; |
| 193 | least_squares_accumulate(mat[1], y[1], a[1], b[1], n); |
Debargha Mukherjee | e6eb3b5 | 2017-02-26 08:50:56 -0800 | [diff] [blame] | 194 | } |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 195 | |
Rachel Barker | 24ec16a | 2023-02-09 13:47:35 +0000 | [diff] [blame] | 196 | if (!least_squares_solve(mat[0], y[0], x[0], n)) { |
| 197 | return false; |
| 198 | } |
| 199 | if (!least_squares_solve(mat[1], y[1], x[1], n)) { |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 200 | return false; |
Debargha Mukherjee | e6eb3b5 | 2017-02-26 08:50:56 -0800 | [diff] [blame] | 201 | } |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 202 | |
| 203 | // Rearrange least squares result to form output model |
Rachel Barker | 24ec16a | 2023-02-09 13:47:35 +0000 | [diff] [blame] | 204 | params[0] = x[0][0]; |
| 205 | params[1] = x[1][0]; |
| 206 | params[2] = x[0][1]; |
| 207 | params[3] = x[0][2]; |
| 208 | params[4] = x[1][1]; |
| 209 | params[5] = x[1][2]; |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 210 | |
| 211 | return true; |
Debargha Mukherjee | e6eb3b5 | 2017-02-26 08:50:56 -0800 | [diff] [blame] | 212 | } |
| 213 | |
emilkeyder@google.com | f347763 | 2017-03-01 16:29:14 -0500 | [diff] [blame] | 214 | typedef struct { |
| 215 | int num_inliers; |
Rachel Barker | 3604269 | 2023-01-11 18:09:43 +0000 | [diff] [blame] | 216 | double sse; // Sum of squared errors of inliers |
emilkeyder@google.com | f347763 | 2017-03-01 16:29:14 -0500 | [diff] [blame] | 217 | int *inlier_indices; |
| 218 | } RANSAC_MOTION; |
| 219 | |
| 220 | // Return -1 if 'a' is a better motion, 1 if 'b' is better, 0 otherwise. |
| 221 | static int compare_motions(const void *arg_a, const void *arg_b) { |
| 222 | const RANSAC_MOTION *motion_a = (RANSAC_MOTION *)arg_a; |
| 223 | const RANSAC_MOTION *motion_b = (RANSAC_MOTION *)arg_b; |
| 224 | |
| 225 | if (motion_a->num_inliers > motion_b->num_inliers) return -1; |
| 226 | if (motion_a->num_inliers < motion_b->num_inliers) return 1; |
Rachel Barker | 3604269 | 2023-01-11 18:09:43 +0000 | [diff] [blame] | 227 | if (motion_a->sse < motion_b->sse) return -1; |
| 228 | if (motion_a->sse > motion_b->sse) return 1; |
emilkeyder@google.com | f347763 | 2017-03-01 16:29:14 -0500 | [diff] [blame] | 229 | return 0; |
| 230 | } |
| 231 | |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 232 | static bool is_better_motion(const RANSAC_MOTION *motion_a, |
| 233 | const RANSAC_MOTION *motion_b) { |
emilkeyder@google.com | f347763 | 2017-03-01 16:29:14 -0500 | [diff] [blame] | 234 | return compare_motions(motion_a, motion_b) < 0; |
| 235 | } |
| 236 | |
| 237 | static void copy_points_at_indices(double *dest, const double *src, |
| 238 | const int *indices, int num_points) { |
| 239 | for (int i = 0; i < num_points; ++i) { |
| 240 | const int index = indices[i]; |
| 241 | dest[i * 2] = src[index * 2]; |
| 242 | dest[i * 2 + 1] = src[index * 2 + 1]; |
| 243 | } |
| 244 | } |
| 245 | |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 246 | // Returns true on success, false on error |
| 247 | static bool ransac_internal(const Correspondence *matched_points, int npoints, |
| 248 | MotionModel *motion_models, int num_desired_motions, |
| 249 | const RansacModelInfo *model_info) { |
Rachel Barker | a02796a | 2023-03-02 20:06:47 +0000 | [diff] [blame] | 250 | assert(npoints >= 0); |
emilkeyder@google.com | f347763 | 2017-03-01 16:29:14 -0500 | [diff] [blame] | 251 | int i = 0; |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 252 | int minpts = model_info->minpts; |
| 253 | bool ret_val = true; |
emilkeyder@google.com | f347763 | 2017-03-01 16:29:14 -0500 | [diff] [blame] | 254 | |
Sarah Parker | efa6582 | 2016-10-11 12:29:07 -0700 | [diff] [blame] | 255 | unsigned int seed = (unsigned int)npoints; |
Sarah Parker | 4dc0f1b | 2016-08-09 17:40:53 -0700 | [diff] [blame] | 256 | |
Sarah Parker | efa6582 | 2016-10-11 12:29:07 -0700 | [diff] [blame] | 257 | int indices[MAX_MINPTS] = { 0 }; |
Sarah Parker | 4dc0f1b | 2016-08-09 17:40:53 -0700 | [diff] [blame] | 258 | |
emilkeyder@google.com | f347763 | 2017-03-01 16:29:14 -0500 | [diff] [blame] | 259 | double *points1, *points2; |
| 260 | double *corners1, *corners2; |
Rachel Barker | 8e99722 | 2023-03-02 15:12:41 +0000 | [diff] [blame] | 261 | double *projected_corners; |
Sarah Parker | 4dc0f1b | 2016-08-09 17:40:53 -0700 | [diff] [blame] | 262 | |
emilkeyder@google.com | f347763 | 2017-03-01 16:29:14 -0500 | [diff] [blame] | 263 | // Store information for the num_desired_motions best transformations found |
| 264 | // and the worst motion among them, as well as the motion currently under |
| 265 | // consideration. |
| 266 | RANSAC_MOTION *motions, *worst_kept_motion = NULL; |
| 267 | RANSAC_MOTION current_motion; |
| 268 | |
| 269 | // Store the parameters and the indices of the inlier points for the motion |
| 270 | // currently under consideration. |
| 271 | double params_this_motion[MAX_PARAMDIM]; |
| 272 | |
Sarah Parker | fa75ae0 | 2016-10-26 12:48:01 -0700 | [diff] [blame] | 273 | if (npoints < minpts * MINPTS_MULTIPLIER || npoints == 0) { |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 274 | return false; |
Sarah Parker | 4dc0f1b | 2016-08-09 17:40:53 -0700 | [diff] [blame] | 275 | } |
| 276 | |
Rachel Barker | a02796a | 2023-03-02 20:06:47 +0000 | [diff] [blame] | 277 | int min_inliers = AOMMAX((int)(MIN_INLIER_PROB * npoints), minpts); |
| 278 | |
emilkeyder@google.com | f347763 | 2017-03-01 16:29:14 -0500 | [diff] [blame] | 279 | points1 = (double *)aom_malloc(sizeof(*points1) * npoints * 2); |
| 280 | points2 = (double *)aom_malloc(sizeof(*points2) * npoints * 2); |
Sarah Parker | f9a961c | 2016-09-06 11:25:04 -0700 | [diff] [blame] | 281 | corners1 = (double *)aom_malloc(sizeof(*corners1) * npoints * 2); |
Sarah Parker | f9a961c | 2016-09-06 11:25:04 -0700 | [diff] [blame] | 282 | corners2 = (double *)aom_malloc(sizeof(*corners2) * npoints * 2); |
Rachel Barker | 8e99722 | 2023-03-02 15:12:41 +0000 | [diff] [blame] | 283 | projected_corners = |
| 284 | (double *)aom_malloc(sizeof(*projected_corners) * npoints * 2); |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 285 | motions = |
Rachel Barker | dfd8046 | 2023-03-02 20:49:41 +0000 | [diff] [blame] | 286 | (RANSAC_MOTION *)aom_calloc(num_desired_motions, sizeof(RANSAC_MOTION)); |
James Zern | 9a0dd4b | 2022-05-17 18:20:49 -0700 | [diff] [blame] | 287 | |
Rachel Barker | 79fc28c | 2023-03-02 19:41:08 +0000 | [diff] [blame] | 288 | // Allocate one large buffer which will be carved up to store the inlier |
| 289 | // indices for the current motion plus the num_desired_motions many |
| 290 | // output models |
| 291 | // This allows us to keep the allocation/deallocation logic simple, without |
| 292 | // having to (for example) check that `motions` is non-null before allocating |
| 293 | // the inlier arrays |
| 294 | int *inlier_buffer = (int *)aom_malloc(sizeof(*inlier_buffer) * npoints * |
| 295 | (num_desired_motions + 1)); |
James Zern | 9a0dd4b | 2022-05-17 18:20:49 -0700 | [diff] [blame] | 296 | |
Rachel Barker | 8e99722 | 2023-03-02 15:12:41 +0000 | [diff] [blame] | 297 | if (!(points1 && points2 && corners1 && corners2 && projected_corners && |
Rachel Barker | 79fc28c | 2023-03-02 19:41:08 +0000 | [diff] [blame] | 298 | motions && inlier_buffer)) { |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 299 | ret_val = false; |
| 300 | goto finish_ransac; |
| 301 | } |
| 302 | |
Rachel Barker | 79fc28c | 2023-03-02 19:41:08 +0000 | [diff] [blame] | 303 | // Once all our allocations are known-good, we can fill in our structures |
| 304 | worst_kept_motion = motions; |
| 305 | |
| 306 | for (i = 0; i < num_desired_motions; ++i) { |
| 307 | motions[i].inlier_indices = inlier_buffer + i * npoints; |
Rachel Barker | 79fc28c | 2023-03-02 19:41:08 +0000 | [diff] [blame] | 308 | } |
Rachel Barker | dfd8046 | 2023-03-02 20:49:41 +0000 | [diff] [blame] | 309 | memset(¤t_motion, 0, sizeof(current_motion)); |
Rachel Barker | 79fc28c | 2023-03-02 19:41:08 +0000 | [diff] [blame] | 310 | current_motion.inlier_indices = inlier_buffer + num_desired_motions * npoints; |
Rachel Barker | 79fc28c | 2023-03-02 19:41:08 +0000 | [diff] [blame] | 311 | |
emilkeyder@google.com | f347763 | 2017-03-01 16:29:14 -0500 | [diff] [blame] | 312 | for (i = 0; i < npoints; ++i) { |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 313 | corners1[2 * i + 0] = matched_points[i].x; |
| 314 | corners1[2 * i + 1] = matched_points[i].y; |
| 315 | corners2[2 * i + 0] = matched_points[i].rx; |
| 316 | corners2[2 * i + 1] = matched_points[i].ry; |
Sarah Parker | 4dc0f1b | 2016-08-09 17:40:53 -0700 | [diff] [blame] | 317 | } |
Sarah Parker | 4dc0f1b | 2016-08-09 17:40:53 -0700 | [diff] [blame] | 318 | |
Rachel Barker | 41f1bb0 | 2023-02-13 16:22:38 +0000 | [diff] [blame] | 319 | for (int trial_count = 0; trial_count < NUM_TRIALS; trial_count++) { |
Rachel Barker | febe1a0 | 2023-03-07 15:10:17 +0000 | [diff] [blame] | 320 | lcg_pick(npoints, minpts, indices, &seed); |
emilkeyder@google.com | f347763 | 2017-03-01 16:29:14 -0500 | [diff] [blame] | 321 | |
Rachel Barker | 41f1bb0 | 2023-02-13 16:22:38 +0000 | [diff] [blame] | 322 | copy_points_at_indices(points1, corners1, indices, minpts); |
| 323 | copy_points_at_indices(points2, corners2, indices, minpts); |
emilkeyder@google.com | f347763 | 2017-03-01 16:29:14 -0500 | [diff] [blame] | 324 | |
Rachel Barker | 41f1bb0 | 2023-02-13 16:22:38 +0000 | [diff] [blame] | 325 | if (model_info->is_degenerate(points1)) { |
| 326 | continue; |
Sarah Parker | 4dc0f1b | 2016-08-09 17:40:53 -0700 | [diff] [blame] | 327 | } |
| 328 | |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 329 | if (!model_info->find_transformation(minpts, points1, points2, |
| 330 | params_this_motion)) { |
Sarah Parker | 4dc0f1b | 2016-08-09 17:40:53 -0700 | [diff] [blame] | 331 | continue; |
| 332 | } |
| 333 | |
Rachel Barker | 8e99722 | 2023-03-02 15:12:41 +0000 | [diff] [blame] | 334 | model_info->project_points(params_this_motion, corners1, projected_corners, |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 335 | npoints, 2, 2); |
Sarah Parker | 4dc0f1b | 2016-08-09 17:40:53 -0700 | [diff] [blame] | 336 | |
Rachel Barker | 41f1bb0 | 2023-02-13 16:22:38 +0000 | [diff] [blame] | 337 | current_motion.num_inliers = 0; |
Rachel Barker | 3604269 | 2023-01-11 18:09:43 +0000 | [diff] [blame] | 338 | double sse = 0.0; |
Sarah Parker | 4dc0f1b | 2016-08-09 17:40:53 -0700 | [diff] [blame] | 339 | for (i = 0; i < npoints; ++i) { |
Rachel Barker | 8e99722 | 2023-03-02 15:12:41 +0000 | [diff] [blame] | 340 | double dx = projected_corners[i * 2] - corners2[i * 2]; |
| 341 | double dy = projected_corners[i * 2 + 1] - corners2[i * 2 + 1]; |
Rachel Barker | 3604269 | 2023-01-11 18:09:43 +0000 | [diff] [blame] | 342 | double squared_error = dx * dx + dy * dy; |
Sarah Parker | 4dc0f1b | 2016-08-09 17:40:53 -0700 | [diff] [blame] | 343 | |
Rachel Barker | 3604269 | 2023-01-11 18:09:43 +0000 | [diff] [blame] | 344 | if (squared_error < INLIER_THRESHOLD_SQUARED) { |
emilkeyder@google.com | f347763 | 2017-03-01 16:29:14 -0500 | [diff] [blame] | 345 | current_motion.inlier_indices[current_motion.num_inliers++] = i; |
Rachel Barker | 3604269 | 2023-01-11 18:09:43 +0000 | [diff] [blame] | 346 | sse += squared_error; |
Sarah Parker | 4dc0f1b | 2016-08-09 17:40:53 -0700 | [diff] [blame] | 347 | } |
| 348 | } |
| 349 | |
Rachel Barker | a02796a | 2023-03-02 20:06:47 +0000 | [diff] [blame] | 350 | if (current_motion.num_inliers < min_inliers) { |
Rachel Barker | 41f1bb0 | 2023-02-13 16:22:38 +0000 | [diff] [blame] | 351 | // Reject models with too few inliers |
| 352 | continue; |
| 353 | } |
| 354 | |
| 355 | current_motion.sse = sse; |
| 356 | if (is_better_motion(¤t_motion, worst_kept_motion)) { |
| 357 | // This motion is better than the worst currently kept motion. Remember |
| 358 | // the inlier points and sse. The parameters for each kept motion |
| 359 | // will be recomputed later using only the inliers. |
| 360 | worst_kept_motion->num_inliers = current_motion.num_inliers; |
| 361 | worst_kept_motion->sse = current_motion.sse; |
Rachel Barker | dfd8046 | 2023-03-02 20:49:41 +0000 | [diff] [blame] | 362 | |
| 363 | // Rather than copying the (potentially many) inlier indices from |
| 364 | // current_motion.inlier_indices to worst_kept_motion->inlier_indices, |
| 365 | // we can swap the underlying pointers. |
| 366 | // |
| 367 | // This is okay because the next time current_motion.inlier_indices |
| 368 | // is used will be in the next trial, where we ignore its previous |
| 369 | // contents anyway. And both arrays will be deallocated together at the |
| 370 | // end of this function, so there are no lifetime issues. |
| 371 | int *tmp = worst_kept_motion->inlier_indices; |
| 372 | worst_kept_motion->inlier_indices = current_motion.inlier_indices; |
| 373 | current_motion.inlier_indices = tmp; |
| 374 | |
Rachel Barker | 41f1bb0 | 2023-02-13 16:22:38 +0000 | [diff] [blame] | 375 | // Determine the new worst kept motion and its num_inliers and sse. |
| 376 | for (i = 0; i < num_desired_motions; ++i) { |
| 377 | if (is_better_motion(worst_kept_motion, &motions[i])) { |
| 378 | worst_kept_motion = &motions[i]; |
emilkeyder@google.com | f347763 | 2017-03-01 16:29:14 -0500 | [diff] [blame] | 379 | } |
Sarah Parker | 4dc0f1b | 2016-08-09 17:40:53 -0700 | [diff] [blame] | 380 | } |
| 381 | } |
Sarah Parker | 4dc0f1b | 2016-08-09 17:40:53 -0700 | [diff] [blame] | 382 | } |
emilkeyder@google.com | f347763 | 2017-03-01 16:29:14 -0500 | [diff] [blame] | 383 | |
| 384 | // Sort the motions, best first. |
| 385 | qsort(motions, num_desired_motions, sizeof(RANSAC_MOTION), compare_motions); |
| 386 | |
| 387 | // Recompute the motions using only the inliers. |
| 388 | for (i = 0; i < num_desired_motions; ++i) { |
Rachel Barker | 41f1bb0 | 2023-02-13 16:22:38 +0000 | [diff] [blame] | 389 | int num_inliers = motions[i].num_inliers; |
| 390 | if (num_inliers > 0) { |
| 391 | assert(num_inliers >= minpts); |
| 392 | |
Debargha Mukherjee | afe7c5f | 2017-06-22 13:47:23 -0700 | [diff] [blame] | 393 | copy_points_at_indices(points1, corners1, motions[i].inlier_indices, |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 394 | num_inliers); |
Debargha Mukherjee | afe7c5f | 2017-06-22 13:47:23 -0700 | [diff] [blame] | 395 | copy_points_at_indices(points2, corners2, motions[i].inlier_indices, |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 396 | num_inliers); |
emilkeyder@google.com | f347763 | 2017-03-01 16:29:14 -0500 | [diff] [blame] | 397 | |
Rachel Barker | a02796a | 2023-03-02 20:06:47 +0000 | [diff] [blame] | 398 | if (!model_info->find_transformation(num_inliers, points1, points2, |
| 399 | motion_models[i].params)) { |
| 400 | // In the unlikely event that this model fitting fails, |
| 401 | // we don't have a good fallback. So just clear the output |
| 402 | // model and move on |
| 403 | memcpy(motion_models[i].params, kIdentityParams, |
Rachel Barker | 0d9df42 | 2023-03-07 20:51:55 +0000 | [diff] [blame] | 404 | MAX_PARAMDIM * sizeof(*(motion_models[i].params))); |
Rachel Barker | a02796a | 2023-03-02 20:06:47 +0000 | [diff] [blame] | 405 | motion_models[i].num_inliers = 0; |
| 406 | continue; |
| 407 | } |
Sarah Parker | 0c74c85 | 2019-05-09 15:15:32 -0700 | [diff] [blame] | 408 | |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 409 | // Populate inliers array |
| 410 | for (int j = 0; j < num_inliers; j++) { |
| 411 | int index = motions[i].inlier_indices[j]; |
| 412 | const Correspondence *corr = &matched_points[index]; |
| 413 | motion_models[i].inliers[2 * j + 0] = (int)rint(corr->x); |
| 414 | motion_models[i].inliers[2 * j + 1] = (int)rint(corr->y); |
| 415 | } |
Rachel Barker | a02796a | 2023-03-02 20:06:47 +0000 | [diff] [blame] | 416 | motion_models[i].num_inliers = num_inliers; |
| 417 | } else { |
| 418 | memcpy(motion_models[i].params, kIdentityParams, |
Rachel Barker | 0d9df42 | 2023-03-07 20:51:55 +0000 | [diff] [blame] | 419 | MAX_PARAMDIM * sizeof(*(motion_models[i].params))); |
Rachel Barker | a02796a | 2023-03-02 20:06:47 +0000 | [diff] [blame] | 420 | motion_models[i].num_inliers = 0; |
Debargha Mukherjee | afe7c5f | 2017-06-22 13:47:23 -0700 | [diff] [blame] | 421 | } |
emilkeyder@google.com | f347763 | 2017-03-01 16:29:14 -0500 | [diff] [blame] | 422 | } |
| 423 | |
Sarah Parker | 4dc0f1b | 2016-08-09 17:40:53 -0700 | [diff] [blame] | 424 | finish_ransac: |
Rachel Barker | 79fc28c | 2023-03-02 19:41:08 +0000 | [diff] [blame] | 425 | aom_free(inlier_buffer); |
| 426 | aom_free(motions); |
Rachel Barker | 8e99722 | 2023-03-02 15:12:41 +0000 | [diff] [blame] | 427 | aom_free(projected_corners); |
Rachel Barker | 79fc28c | 2023-03-02 19:41:08 +0000 | [diff] [blame] | 428 | aom_free(corners2); |
| 429 | aom_free(corners1); |
| 430 | aom_free(points2); |
| 431 | aom_free(points1); |
emilkeyder@google.com | f347763 | 2017-03-01 16:29:14 -0500 | [diff] [blame] | 432 | |
Sarah Parker | 4dc0f1b | 2016-08-09 17:40:53 -0700 | [diff] [blame] | 433 | return ret_val; |
| 434 | } |
| 435 | |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 436 | static bool is_collinear3(double *p1, double *p2, double *p3) { |
Sarah Parker | 4dc0f1b | 2016-08-09 17:40:53 -0700 | [diff] [blame] | 437 | static const double collinear_eps = 1e-3; |
| 438 | const double v = |
| 439 | (p2[0] - p1[0]) * (p3[1] - p1[1]) - (p2[1] - p1[1]) * (p3[0] - p1[0]); |
| 440 | return fabs(v) < collinear_eps; |
| 441 | } |
| 442 | |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 443 | #if ALLOW_TRANSLATION_MODELS |
| 444 | static bool is_degenerate_translation(double *p) { |
Sarah Parker | 4dc0f1b | 2016-08-09 17:40:53 -0700 | [diff] [blame] | 445 | return (p[0] - p[2]) * (p[0] - p[2]) + (p[1] - p[3]) * (p[1] - p[3]) <= 2; |
| 446 | } |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 447 | #endif // ALLOW_TRANSLATION_MODELS |
Sarah Parker | 4dc0f1b | 2016-08-09 17:40:53 -0700 | [diff] [blame] | 448 | |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 449 | static bool is_degenerate_affine(double *p) { |
Sarah Parker | 4dc0f1b | 2016-08-09 17:40:53 -0700 | [diff] [blame] | 450 | return is_collinear3(p, p + 2, p + 4); |
| 451 | } |
| 452 | |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 453 | static const RansacModelInfo ransac_model_info[TRANS_TYPES] = { |
| 454 | // IDENTITY |
| 455 | { NULL, NULL, NULL, 0 }, |
| 456 | // TRANSLATION |
| 457 | #if ALLOW_TRANSLATION_MODELS |
| 458 | { is_degenerate_translation, find_translation, project_points_translation, |
| 459 | 3 }, |
| 460 | #else |
| 461 | { NULL, NULL, NULL, 0 }, |
| 462 | #endif |
| 463 | // ROTZOOM |
| 464 | { is_degenerate_affine, find_rotzoom, project_points_affine, 3 }, |
| 465 | // AFFINE |
| 466 | { is_degenerate_affine, find_affine, project_points_affine, 3 }, |
| 467 | }; |
Sarah Parker | 4dc0f1b | 2016-08-09 17:40:53 -0700 | [diff] [blame] | 468 | |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 469 | // Returns true on success, false on error |
| 470 | bool ransac(const Correspondence *matched_points, int npoints, |
| 471 | TransformationType type, MotionModel *motion_models, |
| 472 | int num_desired_motions) { |
| 473 | #if ALLOW_TRANSLATION_MODELS |
| 474 | assert(type > IDENTITY && type < TRANS_TYPES); |
| 475 | #else |
| 476 | assert(type > TRANSLATION && type < TRANS_TYPES); |
| 477 | #endif // ALLOW_TRANSLATION_MODELS |
Sarah Parker | f9a961c | 2016-09-06 11:25:04 -0700 | [diff] [blame] | 478 | |
Rachel Barker | 738c4f4 | 2023-01-11 17:36:18 +0000 | [diff] [blame] | 479 | return ransac_internal(matched_points, npoints, motion_models, |
| 480 | num_desired_motions, &ransac_model_info[type]); |
Yaowu Xu | 0e4509c | 2019-05-02 15:59:17 -0700 | [diff] [blame] | 481 | } |