Line data Source code
1 : /*
2 : * Copyright (c) 2016, Alliance for Open Media. All rights reserved
3 : *
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.
10 : */
11 :
12 : #include <stdlib.h>
13 : #include <memory.h>
14 : #include <math.h>
15 :
16 : #include "corner_match.h"
17 :
18 : #define SEARCH_SZ 9
19 : #define SEARCH_SZ_BY2 ((SEARCH_SZ - 1) / 2)
20 :
21 : #define THRESHOLD_NCC 0.75
22 :
23 : /* Compute var(im) * MATCH_SZ_SQ over a MATCH_SZ by MATCH_SZ window of im,
24 : centered at (x, y).
25 : */
26 369020 : static double compute_variance(unsigned char *im, int stride, int x, int y) {
27 369020 : int sum = 0;
28 369020 : int sumsq = 0;
29 : int var;
30 : int i, j;
31 5165450 : for (i = 0; i < MATCH_SZ; ++i)
32 67147000 : for (j = 0; j < MATCH_SZ; ++j) {
33 62350600 : sum += im[(i + y - MATCH_SZ_BY2) * stride + (j + x - MATCH_SZ_BY2)];
34 62350600 : sumsq += im[(i + y - MATCH_SZ_BY2) * stride + (j + x - MATCH_SZ_BY2)] *
35 62350600 : im[(i + y - MATCH_SZ_BY2) * stride + (j + x - MATCH_SZ_BY2)];
36 : }
37 369020 : var = sumsq * MATCH_SZ_SQ - sum * sum;
38 369020 : return (double)var;
39 : }
40 :
41 : /* Compute corr(im1, im2) * MATCH_SZ * stddev(im1), where the
42 : correlation/standard deviation are taken over MATCH_SZ by MATCH_SZ windows
43 : of each image, centered at (x1, y1) and (x2, y2) respectively.
44 : */
45 104459000 : double av1_compute_cross_correlation_c(unsigned char *im1, int stride1, int x1,
46 : int y1, unsigned char *im2, int stride2,
47 : int x2, int y2) {
48 : int v1, v2;
49 104459000 : int sum1 = 0;
50 104459000 : int sum2 = 0;
51 104459000 : int sumsq2 = 0;
52 104459000 : int cross = 0;
53 : int var2, cov;
54 : int i, j;
55 1396870000 : for (i = 0; i < MATCH_SZ; ++i)
56 17670400000 : for (j = 0; j < MATCH_SZ; ++j) {
57 16378000000 : v1 = im1[(i + y1 - MATCH_SZ_BY2) * stride1 + (j + x1 - MATCH_SZ_BY2)];
58 16378000000 : v2 = im2[(i + y2 - MATCH_SZ_BY2) * stride2 + (j + x2 - MATCH_SZ_BY2)];
59 16378000000 : sum1 += v1;
60 16378000000 : sum2 += v2;
61 16378000000 : sumsq2 += v2 * v2;
62 16378000000 : cross += v1 * v2;
63 : }
64 104459000 : var2 = sumsq2 * MATCH_SZ_SQ - sum2 * sum2;
65 104459000 : cov = cross * MATCH_SZ_SQ - sum1 * sum2;
66 104459000 : return cov / sqrt((double)var2);
67 : }
68 :
69 1276410000 : static int is_eligible_point(int pointx, int pointy, int width, int height) {
70 1264830000 : return (pointx >= MATCH_SZ_BY2 && pointy >= MATCH_SZ_BY2 &&
71 2541240000 : pointx + MATCH_SZ_BY2 < width && pointy + MATCH_SZ_BY2 < height);
72 : }
73 :
74 1234900000 : static int is_eligible_distance(int point1x, int point1y, int point2x,
75 : int point2y, int width, int height) {
76 1234900000 : const int thresh = (width < height ? height : width) >> 4;
77 1234900000 : return ((point1x - point2x) * (point1x - point2x) +
78 1234900000 : (point1y - point2y) * (point1y - point2y)) <= thresh * thresh;
79 : }
80 :
81 113 : static void improve_correspondence(unsigned char *frm, unsigned char *ref,
82 : int width, int height, int frm_stride,
83 : int ref_stride,
84 : Correspondence *correspondences,
85 : int num_correspondences) {
86 : int i;
87 283034 : for (i = 0; i < num_correspondences; ++i) {
88 282909 : int x, y, best_x = 0, best_y = 0;
89 282909 : double best_match_ncc = 0.0;
90 2829090 : for (y = -SEARCH_SZ_BY2; y <= SEARCH_SZ_BY2; ++y) {
91 25460600 : for (x = -SEARCH_SZ_BY2; x <= SEARCH_SZ_BY2; ++x) {
92 : double match_ncc;
93 22914500 : if (!is_eligible_point(correspondences[i].rx + x,
94 22914400 : correspondences[i].ry + y, width, height))
95 156780 : continue;
96 22757700 : if (!is_eligible_distance(correspondences[i].x, correspondences[i].y,
97 22757700 : correspondences[i].rx + x,
98 22757700 : correspondences[i].ry + y, width, height))
99 31502 : continue;
100 22726200 : match_ncc = av1_compute_cross_correlation_c(
101 22726200 : frm, frm_stride, correspondences[i].x, correspondences[i].y, ref,
102 22726200 : ref_stride, correspondences[i].rx + x, correspondences[i].ry + y);
103 22726200 : if (match_ncc > best_match_ncc) {
104 2963740 : best_match_ncc = match_ncc;
105 2963740 : best_y = y;
106 2963740 : best_x = x;
107 : }
108 : }
109 : }
110 282921 : correspondences[i].rx += best_x;
111 282921 : correspondences[i].ry += best_y;
112 : }
113 282798 : for (i = 0; i < num_correspondences; ++i) {
114 282908 : int x, y, best_x = 0, best_y = 0;
115 282908 : double best_match_ncc = 0.0;
116 2828810 : for (y = -SEARCH_SZ_BY2; y <= SEARCH_SZ_BY2; ++y)
117 25458400 : for (x = -SEARCH_SZ_BY2; x <= SEARCH_SZ_BY2; ++x) {
118 : double match_ncc;
119 22913000 : if (!is_eligible_point(correspondences[i].x + x,
120 22912500 : correspondences[i].y + y, width, height))
121 160898 : continue;
122 22752400 : if (!is_eligible_distance(
123 22752100 : correspondences[i].x + x, correspondences[i].y + y,
124 22752100 : correspondences[i].rx, correspondences[i].ry, width, height))
125 31163 : continue;
126 22721200 : match_ncc = av1_compute_cross_correlation_c(
127 22721200 : ref, ref_stride, correspondences[i].rx, correspondences[i].ry, frm,
128 22721200 : frm_stride, correspondences[i].x + x, correspondences[i].y + y);
129 22720200 : if (match_ncc > best_match_ncc) {
130 2973280 : best_match_ncc = match_ncc;
131 2973280 : best_y = y;
132 2973280 : best_x = x;
133 : }
134 : }
135 282673 : correspondences[i].x += best_x;
136 282673 : correspondences[i].y += best_y;
137 : }
138 0 : }
139 :
140 113 : int av1_determine_correspondence(unsigned char *frm, int *frm_corners,
141 : int num_frm_corners, unsigned char *ref,
142 : int *ref_corners, int num_ref_corners,
143 : int width, int height, int frm_stride,
144 : int ref_stride, int *correspondence_pts) {
145 : // TODO(sarahparker) Improve this to include 2-way match
146 : int i, j;
147 113 : Correspondence *correspondences = (Correspondence *)correspondence_pts;
148 113 : int num_correspondences = 0;
149 381816 : for (i = 0; i < num_frm_corners; ++i) {
150 381704 : double best_match_ncc = 0.0;
151 : double template_norm;
152 381704 : int best_match_j = -1;
153 381704 : if (!is_eligible_point(frm_corners[2 * i], frm_corners[2 * i + 1], width,
154 : height))
155 12681 : continue;
156 1234350000 : for (j = 0; j < num_ref_corners; ++j) {
157 : double match_ncc;
158 1232250000 : if (!is_eligible_point(ref_corners[2 * j], ref_corners[2 * j + 1], width,
159 : height))
160 41499000 : continue;
161 1191310000 : if (!is_eligible_distance(frm_corners[2 * i], frm_corners[2 * i + 1],
162 1190720000 : ref_corners[2 * j], ref_corners[2 * j + 1],
163 : width, height))
164 1133270000 : continue;
165 58040100 : match_ncc = av1_compute_cross_correlation_c(
166 58040100 : frm, frm_stride, frm_corners[2 * i], frm_corners[2 * i + 1], ref,
167 58040100 : ref_stride, ref_corners[2 * j], ref_corners[2 * j + 1]);
168 59207600 : if (match_ncc > best_match_ncc) {
169 2020630 : best_match_ncc = match_ncc;
170 2020630 : best_match_j = j;
171 : }
172 : }
173 : // Note: We want to test if the best correlation is >= THRESHOLD_NCC,
174 : // but need to account for the normalization in
175 : // av1_compute_cross_correlation.
176 2095460 : template_norm = compute_variance(frm, frm_stride, frm_corners[2 * i],
177 2095460 : frm_corners[2 * i + 1]);
178 369022 : if (best_match_ncc > THRESHOLD_NCC * sqrt(template_norm)) {
179 282909 : correspondences[num_correspondences].x = frm_corners[2 * i];
180 282909 : correspondences[num_correspondences].y = frm_corners[2 * i + 1];
181 282909 : correspondences[num_correspondences].rx = ref_corners[2 * best_match_j];
182 282909 : correspondences[num_correspondences].ry =
183 282909 : ref_corners[2 * best_match_j + 1];
184 :
185 : /*printf("corresp: %d %d - %d %d\n",
186 : correspondences[num_correspondences].x,
187 : correspondences[num_correspondences].y,
188 : correspondences[num_correspondences].rx,
189 : correspondences[num_correspondences].ry);*/
190 :
191 282909 : num_correspondences++;
192 : }
193 : }
194 112 : improve_correspondence(frm, ref, width, height, frm_stride, ref_stride,
195 : correspondences, num_correspondences);
196 113 : return num_correspondences;
197 : }
|