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 "config/av1_rtcd.h"
17
18 #include "av1/encoder/corner_match.h"
19
20 #define SEARCH_SZ 9
21 #define SEARCH_SZ_BY2 ((SEARCH_SZ - 1) / 2)
22
23 #define THRESHOLD_NCC 0.75
24
25 /* Compute var(im) * MATCH_SZ_SQ over a MATCH_SZ by MATCH_SZ window of im,
26 centered at (x, y).
27 */
compute_variance(unsigned char * im,int stride,int x,int y)28 static double compute_variance(unsigned char *im, int stride, int x, int y) {
29 int sum = 0;
30 int sumsq = 0;
31 int var;
32 int i, j;
33 for (i = 0; i < MATCH_SZ; ++i)
34 for (j = 0; j < MATCH_SZ; ++j) {
35 sum += im[(i + y - MATCH_SZ_BY2) * stride + (j + x - MATCH_SZ_BY2)];
36 sumsq += im[(i + y - MATCH_SZ_BY2) * stride + (j + x - MATCH_SZ_BY2)] *
37 im[(i + y - MATCH_SZ_BY2) * stride + (j + x - MATCH_SZ_BY2)];
38 }
39 var = sumsq * MATCH_SZ_SQ - sum * sum;
40 return (double)var;
41 }
42
43 /* Compute corr(im1, im2) * MATCH_SZ * stddev(im1), where the
44 correlation/standard deviation are taken over MATCH_SZ by MATCH_SZ windows
45 of each image, centered at (x1, y1) and (x2, y2) respectively.
46 */
compute_cross_correlation_c(unsigned char * im1,int stride1,int x1,int y1,unsigned char * im2,int stride2,int x2,int y2)47 double compute_cross_correlation_c(unsigned char *im1, int stride1, int x1,
48 int y1, unsigned char *im2, int stride2,
49 int x2, int y2) {
50 int v1, v2;
51 int sum1 = 0;
52 int sum2 = 0;
53 int sumsq2 = 0;
54 int cross = 0;
55 int var2, cov;
56 int i, j;
57 for (i = 0; i < MATCH_SZ; ++i)
58 for (j = 0; j < MATCH_SZ; ++j) {
59 v1 = im1[(i + y1 - MATCH_SZ_BY2) * stride1 + (j + x1 - MATCH_SZ_BY2)];
60 v2 = im2[(i + y2 - MATCH_SZ_BY2) * stride2 + (j + x2 - MATCH_SZ_BY2)];
61 sum1 += v1;
62 sum2 += v2;
63 sumsq2 += v2 * v2;
64 cross += v1 * v2;
65 }
66 var2 = sumsq2 * MATCH_SZ_SQ - sum2 * sum2;
67 cov = cross * MATCH_SZ_SQ - sum1 * sum2;
68 return cov / sqrt((double)var2);
69 }
70
is_eligible_point(int pointx,int pointy,int width,int height)71 static int is_eligible_point(int pointx, int pointy, int width, int height) {
72 return (pointx >= MATCH_SZ_BY2 && pointy >= MATCH_SZ_BY2 &&
73 pointx + MATCH_SZ_BY2 < width && pointy + MATCH_SZ_BY2 < height);
74 }
75
is_eligible_distance(int point1x,int point1y,int point2x,int point2y,int width,int height)76 static int is_eligible_distance(int point1x, int point1y, int point2x,
77 int point2y, int width, int height) {
78 const int thresh = (width < height ? height : width) >> 4;
79 return ((point1x - point2x) * (point1x - point2x) +
80 (point1y - point2y) * (point1y - point2y)) <= thresh * thresh;
81 }
82
improve_correspondence(unsigned char * frm,unsigned char * ref,int width,int height,int frm_stride,int ref_stride,Correspondence * correspondences,int num_correspondences)83 static void improve_correspondence(unsigned char *frm, unsigned char *ref,
84 int width, int height, int frm_stride,
85 int ref_stride,
86 Correspondence *correspondences,
87 int num_correspondences) {
88 int i;
89 for (i = 0; i < num_correspondences; ++i) {
90 int x, y, best_x = 0, best_y = 0;
91 double best_match_ncc = 0.0;
92 for (y = -SEARCH_SZ_BY2; y <= SEARCH_SZ_BY2; ++y) {
93 for (x = -SEARCH_SZ_BY2; x <= SEARCH_SZ_BY2; ++x) {
94 double match_ncc;
95 if (!is_eligible_point(correspondences[i].rx + x,
96 correspondences[i].ry + y, width, height))
97 continue;
98 if (!is_eligible_distance(correspondences[i].x, correspondences[i].y,
99 correspondences[i].rx + x,
100 correspondences[i].ry + y, width, height))
101 continue;
102 match_ncc = compute_cross_correlation(
103 frm, frm_stride, correspondences[i].x, correspondences[i].y, ref,
104 ref_stride, correspondences[i].rx + x, correspondences[i].ry + y);
105 if (match_ncc > best_match_ncc) {
106 best_match_ncc = match_ncc;
107 best_y = y;
108 best_x = x;
109 }
110 }
111 }
112 correspondences[i].rx += best_x;
113 correspondences[i].ry += best_y;
114 }
115 for (i = 0; i < num_correspondences; ++i) {
116 int x, y, best_x = 0, best_y = 0;
117 double best_match_ncc = 0.0;
118 for (y = -SEARCH_SZ_BY2; y <= SEARCH_SZ_BY2; ++y)
119 for (x = -SEARCH_SZ_BY2; x <= SEARCH_SZ_BY2; ++x) {
120 double match_ncc;
121 if (!is_eligible_point(correspondences[i].x + x,
122 correspondences[i].y + y, width, height))
123 continue;
124 if (!is_eligible_distance(
125 correspondences[i].x + x, correspondences[i].y + y,
126 correspondences[i].rx, correspondences[i].ry, width, height))
127 continue;
128 match_ncc = compute_cross_correlation(
129 ref, ref_stride, correspondences[i].rx, correspondences[i].ry, frm,
130 frm_stride, correspondences[i].x + x, correspondences[i].y + y);
131 if (match_ncc > best_match_ncc) {
132 best_match_ncc = match_ncc;
133 best_y = y;
134 best_x = x;
135 }
136 }
137 correspondences[i].x += best_x;
138 correspondences[i].y += best_y;
139 }
140 }
141
determine_correspondence(unsigned char * frm,int * frm_corners,int num_frm_corners,unsigned char * ref,int * ref_corners,int num_ref_corners,int width,int height,int frm_stride,int ref_stride,int * correspondence_pts)142 int determine_correspondence(unsigned char *frm, int *frm_corners,
143 int num_frm_corners, unsigned char *ref,
144 int *ref_corners, int num_ref_corners, int width,
145 int height, int frm_stride, int ref_stride,
146 int *correspondence_pts) {
147 // TODO(sarahparker) Improve this to include 2-way match
148 int i, j;
149 Correspondence *correspondences = (Correspondence *)correspondence_pts;
150 int num_correspondences = 0;
151 for (i = 0; i < num_frm_corners; ++i) {
152 double best_match_ncc = 0.0;
153 double template_norm;
154 int best_match_j = -1;
155 if (!is_eligible_point(frm_corners[2 * i], frm_corners[2 * i + 1], width,
156 height))
157 continue;
158 for (j = 0; j < num_ref_corners; ++j) {
159 double match_ncc;
160 if (!is_eligible_point(ref_corners[2 * j], ref_corners[2 * j + 1], width,
161 height))
162 continue;
163 if (!is_eligible_distance(frm_corners[2 * i], frm_corners[2 * i + 1],
164 ref_corners[2 * j], ref_corners[2 * j + 1],
165 width, height))
166 continue;
167 match_ncc = compute_cross_correlation(
168 frm, frm_stride, frm_corners[2 * i], frm_corners[2 * i + 1], ref,
169 ref_stride, ref_corners[2 * j], ref_corners[2 * j + 1]);
170 if (match_ncc > best_match_ncc) {
171 best_match_ncc = match_ncc;
172 best_match_j = j;
173 }
174 }
175 // Note: We want to test if the best correlation is >= THRESHOLD_NCC,
176 // but need to account for the normalization in compute_cross_correlation.
177 template_norm = compute_variance(frm, frm_stride, frm_corners[2 * i],
178 frm_corners[2 * i + 1]);
179 if (best_match_ncc > THRESHOLD_NCC * sqrt(template_norm)) {
180 correspondences[num_correspondences].x = frm_corners[2 * i];
181 correspondences[num_correspondences].y = frm_corners[2 * i + 1];
182 correspondences[num_correspondences].rx = ref_corners[2 * best_match_j];
183 correspondences[num_correspondences].ry =
184 ref_corners[2 * best_match_j + 1];
185 num_correspondences++;
186 }
187 }
188 improve_correspondence(frm, ref, width, height, frm_stride, ref_stride,
189 correspondences, num_correspondences);
190 return num_correspondences;
191 }
192