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 <stdio.h>
13 : #include <stdlib.h>
14 : #include <math.h>
15 : #include <assert.h>
16 : #include "aom_dsp_rtcd.h"
17 : #include "EbWarpedMotion.h"
18 :
19 : #define WARP_ERROR_BLOCK 32
20 :
21 : /* clang-format off */
22 : const int error_measure_lut[512] = {
23 : // pow 0.7
24 : 16384, 16339, 16294, 16249, 16204, 16158, 16113, 16068,
25 : 16022, 15977, 15932, 15886, 15840, 15795, 15749, 15703,
26 : 15657, 15612, 15566, 15520, 15474, 15427, 15381, 15335,
27 : 15289, 15242, 15196, 15149, 15103, 15056, 15010, 14963,
28 : 14916, 14869, 14822, 14775, 14728, 14681, 14634, 14587,
29 : 14539, 14492, 14445, 14397, 14350, 14302, 14254, 14206,
30 : 14159, 14111, 14063, 14015, 13967, 13918, 13870, 13822,
31 : 13773, 13725, 13676, 13628, 13579, 13530, 13481, 13432,
32 : 13383, 13334, 13285, 13236, 13187, 13137, 13088, 13038,
33 : 12988, 12939, 12889, 12839, 12789, 12739, 12689, 12639,
34 : 12588, 12538, 12487, 12437, 12386, 12335, 12285, 12234,
35 : 12183, 12132, 12080, 12029, 11978, 11926, 11875, 11823,
36 : 11771, 11719, 11667, 11615, 11563, 11511, 11458, 11406,
37 : 11353, 11301, 11248, 11195, 11142, 11089, 11036, 10982,
38 : 10929, 10875, 10822, 10768, 10714, 10660, 10606, 10552,
39 : 10497, 10443, 10388, 10333, 10279, 10224, 10168, 10113,
40 : 10058, 10002, 9947, 9891, 9835, 9779, 9723, 9666,
41 : 9610, 9553, 9497, 9440, 9383, 9326, 9268, 9211,
42 : 9153, 9095, 9037, 8979, 8921, 8862, 8804, 8745,
43 : 8686, 8627, 8568, 8508, 8449, 8389, 8329, 8269,
44 : 8208, 8148, 8087, 8026, 7965, 7903, 7842, 7780,
45 : 7718, 7656, 7593, 7531, 7468, 7405, 7341, 7278,
46 : 7214, 7150, 7086, 7021, 6956, 6891, 6826, 6760,
47 : 6695, 6628, 6562, 6495, 6428, 6361, 6293, 6225,
48 : 6157, 6089, 6020, 5950, 5881, 5811, 5741, 5670,
49 : 5599, 5527, 5456, 5383, 5311, 5237, 5164, 5090,
50 : 5015, 4941, 4865, 4789, 4713, 4636, 4558, 4480,
51 : 4401, 4322, 4242, 4162, 4080, 3998, 3916, 3832,
52 : 3748, 3663, 3577, 3490, 3402, 3314, 3224, 3133,
53 : 3041, 2948, 2854, 2758, 2661, 2562, 2461, 2359,
54 : 2255, 2148, 2040, 1929, 1815, 1698, 1577, 1452,
55 : 1323, 1187, 1045, 894, 731, 550, 339, 0,
56 : 339, 550, 731, 894, 1045, 1187, 1323, 1452,
57 : 1577, 1698, 1815, 1929, 2040, 2148, 2255, 2359,
58 : 2461, 2562, 2661, 2758, 2854, 2948, 3041, 3133,
59 : 3224, 3314, 3402, 3490, 3577, 3663, 3748, 3832,
60 : 3916, 3998, 4080, 4162, 4242, 4322, 4401, 4480,
61 : 4558, 4636, 4713, 4789, 4865, 4941, 5015, 5090,
62 : 5164, 5237, 5311, 5383, 5456, 5527, 5599, 5670,
63 : 5741, 5811, 5881, 5950, 6020, 6089, 6157, 6225,
64 : 6293, 6361, 6428, 6495, 6562, 6628, 6695, 6760,
65 : 6826, 6891, 6956, 7021, 7086, 7150, 7214, 7278,
66 : 7341, 7405, 7468, 7531, 7593, 7656, 7718, 7780,
67 : 7842, 7903, 7965, 8026, 8087, 8148, 8208, 8269,
68 : 8329, 8389, 8449, 8508, 8568, 8627, 8686, 8745,
69 : 8804, 8862, 8921, 8979, 9037, 9095, 9153, 9211,
70 : 9268, 9326, 9383, 9440, 9497, 9553, 9610, 9666,
71 : 9723, 9779, 9835, 9891, 9947, 10002, 10058, 10113,
72 : 10168, 10224, 10279, 10333, 10388, 10443, 10497, 10552,
73 : 10606, 10660, 10714, 10768, 10822, 10875, 10929, 10982,
74 : 11036, 11089, 11142, 11195, 11248, 11301, 11353, 11406,
75 : 11458, 11511, 11563, 11615, 11667, 11719, 11771, 11823,
76 : 11875, 11926, 11978, 12029, 12080, 12132, 12183, 12234,
77 : 12285, 12335, 12386, 12437, 12487, 12538, 12588, 12639,
78 : 12689, 12739, 12789, 12839, 12889, 12939, 12988, 13038,
79 : 13088, 13137, 13187, 13236, 13285, 13334, 13383, 13432,
80 : 13481, 13530, 13579, 13628, 13676, 13725, 13773, 13822,
81 : 13870, 13918, 13967, 14015, 14063, 14111, 14159, 14206,
82 : 14254, 14302, 14350, 14397, 14445, 14492, 14539, 14587,
83 : 14634, 14681, 14728, 14775, 14822, 14869, 14916, 14963,
84 : 15010, 15056, 15103, 15149, 15196, 15242, 15289, 15335,
85 : 15381, 15427, 15474, 15520, 15566, 15612, 15657, 15703,
86 : 15749, 15795, 15840, 15886, 15932, 15977, 16022, 16068,
87 : 16113, 16158, 16204, 16249, 16294, 16339, 16384, 16384,
88 : };
89 : /* clang-format on */
90 :
91 : // For warping, we really use a 6-tap filter, but we do blocks of 8 pixels
92 : // at a time. The zoom/rotation/shear in the model are applied to the
93 : // "fractional" position of each pixel, which therefore varies within
94 : // [-1, 2) * WARPEDPIXEL_PREC_SHIFTS.
95 : // We need an extra 2 taps to fit this in, for a total of 8 taps.
96 : /* clang-format off */
97 : EB_ALIGN(16) const int16_t eb_warped_filter[WARPEDPIXEL_PREC_SHIFTS * 3 + 1][8] = {
98 : #if WARPEDPIXEL_PREC_BITS == 6
99 : // [-1, 0)
100 : { 0, 0, 127, 1, 0, 0, 0, 0 }, { 0, - 1, 127, 2, 0, 0, 0, 0 },
101 : { 1, - 3, 127, 4, - 1, 0, 0, 0 }, { 1, - 4, 126, 6, - 2, 1, 0, 0 },
102 : { 1, - 5, 126, 8, - 3, 1, 0, 0 }, { 1, - 6, 125, 11, - 4, 1, 0, 0 },
103 : { 1, - 7, 124, 13, - 4, 1, 0, 0 }, { 2, - 8, 123, 15, - 5, 1, 0, 0 },
104 : { 2, - 9, 122, 18, - 6, 1, 0, 0 }, { 2, -10, 121, 20, - 6, 1, 0, 0 },
105 : { 2, -11, 120, 22, - 7, 2, 0, 0 }, { 2, -12, 119, 25, - 8, 2, 0, 0 },
106 : { 3, -13, 117, 27, - 8, 2, 0, 0 }, { 3, -13, 116, 29, - 9, 2, 0, 0 },
107 : { 3, -14, 114, 32, -10, 3, 0, 0 }, { 3, -15, 113, 35, -10, 2, 0, 0 },
108 : { 3, -15, 111, 37, -11, 3, 0, 0 }, { 3, -16, 109, 40, -11, 3, 0, 0 },
109 : { 3, -16, 108, 42, -12, 3, 0, 0 }, { 4, -17, 106, 45, -13, 3, 0, 0 },
110 : { 4, -17, 104, 47, -13, 3, 0, 0 }, { 4, -17, 102, 50, -14, 3, 0, 0 },
111 : { 4, -17, 100, 52, -14, 3, 0, 0 }, { 4, -18, 98, 55, -15, 4, 0, 0 },
112 : { 4, -18, 96, 58, -15, 3, 0, 0 }, { 4, -18, 94, 60, -16, 4, 0, 0 },
113 : { 4, -18, 91, 63, -16, 4, 0, 0 }, { 4, -18, 89, 65, -16, 4, 0, 0 },
114 : { 4, -18, 87, 68, -17, 4, 0, 0 }, { 4, -18, 85, 70, -17, 4, 0, 0 },
115 : { 4, -18, 82, 73, -17, 4, 0, 0 }, { 4, -18, 80, 75, -17, 4, 0, 0 },
116 : { 4, -18, 78, 78, -18, 4, 0, 0 }, { 4, -17, 75, 80, -18, 4, 0, 0 },
117 : { 4, -17, 73, 82, -18, 4, 0, 0 }, { 4, -17, 70, 85, -18, 4, 0, 0 },
118 : { 4, -17, 68, 87, -18, 4, 0, 0 }, { 4, -16, 65, 89, -18, 4, 0, 0 },
119 : { 4, -16, 63, 91, -18, 4, 0, 0 }, { 4, -16, 60, 94, -18, 4, 0, 0 },
120 : { 3, -15, 58, 96, -18, 4, 0, 0 }, { 4, -15, 55, 98, -18, 4, 0, 0 },
121 : { 3, -14, 52, 100, -17, 4, 0, 0 }, { 3, -14, 50, 102, -17, 4, 0, 0 },
122 : { 3, -13, 47, 104, -17, 4, 0, 0 }, { 3, -13, 45, 106, -17, 4, 0, 0 },
123 : { 3, -12, 42, 108, -16, 3, 0, 0 }, { 3, -11, 40, 109, -16, 3, 0, 0 },
124 : { 3, -11, 37, 111, -15, 3, 0, 0 }, { 2, -10, 35, 113, -15, 3, 0, 0 },
125 : { 3, -10, 32, 114, -14, 3, 0, 0 }, { 2, - 9, 29, 116, -13, 3, 0, 0 },
126 : { 2, - 8, 27, 117, -13, 3, 0, 0 }, { 2, - 8, 25, 119, -12, 2, 0, 0 },
127 : { 2, - 7, 22, 120, -11, 2, 0, 0 }, { 1, - 6, 20, 121, -10, 2, 0, 0 },
128 : { 1, - 6, 18, 122, - 9, 2, 0, 0 }, { 1, - 5, 15, 123, - 8, 2, 0, 0 },
129 : { 1, - 4, 13, 124, - 7, 1, 0, 0 }, { 1, - 4, 11, 125, - 6, 1, 0, 0 },
130 : { 1, - 3, 8, 126, - 5, 1, 0, 0 }, { 1, - 2, 6, 126, - 4, 1, 0, 0 },
131 : { 0, - 1, 4, 127, - 3, 1, 0, 0 }, { 0, 0, 2, 127, - 1, 0, 0, 0 },
132 : // [0, 1)
133 : { 0, 0, 0, 127, 1, 0, 0, 0}, { 0, 0, -1, 127, 2, 0, 0, 0},
134 : { 0, 1, -3, 127, 4, -2, 1, 0}, { 0, 1, -5, 127, 6, -2, 1, 0},
135 : { 0, 2, -6, 126, 8, -3, 1, 0}, {-1, 2, -7, 126, 11, -4, 2, -1},
136 : {-1, 3, -8, 125, 13, -5, 2, -1}, {-1, 3, -10, 124, 16, -6, 3, -1},
137 : {-1, 4, -11, 123, 18, -7, 3, -1}, {-1, 4, -12, 122, 20, -7, 3, -1},
138 : {-1, 4, -13, 121, 23, -8, 3, -1}, {-2, 5, -14, 120, 25, -9, 4, -1},
139 : {-1, 5, -15, 119, 27, -10, 4, -1}, {-1, 5, -16, 118, 30, -11, 4, -1},
140 : {-2, 6, -17, 116, 33, -12, 5, -1}, {-2, 6, -17, 114, 35, -12, 5, -1},
141 : {-2, 6, -18, 113, 38, -13, 5, -1}, {-2, 7, -19, 111, 41, -14, 6, -2},
142 : {-2, 7, -19, 110, 43, -15, 6, -2}, {-2, 7, -20, 108, 46, -15, 6, -2},
143 : {-2, 7, -20, 106, 49, -16, 6, -2}, {-2, 7, -21, 104, 51, -16, 7, -2},
144 : {-2, 7, -21, 102, 54, -17, 7, -2}, {-2, 8, -21, 100, 56, -18, 7, -2},
145 : {-2, 8, -22, 98, 59, -18, 7, -2}, {-2, 8, -22, 96, 62, -19, 7, -2},
146 : {-2, 8, -22, 94, 64, -19, 7, -2}, {-2, 8, -22, 91, 67, -20, 8, -2},
147 : {-2, 8, -22, 89, 69, -20, 8, -2}, {-2, 8, -22, 87, 72, -21, 8, -2},
148 : {-2, 8, -21, 84, 74, -21, 8, -2}, {-2, 8, -22, 82, 77, -21, 8, -2},
149 : {-2, 8, -21, 79, 79, -21, 8, -2}, {-2, 8, -21, 77, 82, -22, 8, -2},
150 : {-2, 8, -21, 74, 84, -21, 8, -2}, {-2, 8, -21, 72, 87, -22, 8, -2},
151 : {-2, 8, -20, 69, 89, -22, 8, -2}, {-2, 8, -20, 67, 91, -22, 8, -2},
152 : {-2, 7, -19, 64, 94, -22, 8, -2}, {-2, 7, -19, 62, 96, -22, 8, -2},
153 : {-2, 7, -18, 59, 98, -22, 8, -2}, {-2, 7, -18, 56, 100, -21, 8, -2},
154 : {-2, 7, -17, 54, 102, -21, 7, -2}, {-2, 7, -16, 51, 104, -21, 7, -2},
155 : {-2, 6, -16, 49, 106, -20, 7, -2}, {-2, 6, -15, 46, 108, -20, 7, -2},
156 : {-2, 6, -15, 43, 110, -19, 7, -2}, {-2, 6, -14, 41, 111, -19, 7, -2},
157 : {-1, 5, -13, 38, 113, -18, 6, -2}, {-1, 5, -12, 35, 114, -17, 6, -2},
158 : {-1, 5, -12, 33, 116, -17, 6, -2}, {-1, 4, -11, 30, 118, -16, 5, -1},
159 : {-1, 4, -10, 27, 119, -15, 5, -1}, {-1, 4, -9, 25, 120, -14, 5, -2},
160 : {-1, 3, -8, 23, 121, -13, 4, -1}, {-1, 3, -7, 20, 122, -12, 4, -1},
161 : {-1, 3, -7, 18, 123, -11, 4, -1}, {-1, 3, -6, 16, 124, -10, 3, -1},
162 : {-1, 2, -5, 13, 125, -8, 3, -1}, {-1, 2, -4, 11, 126, -7, 2, -1},
163 : { 0, 1, -3, 8, 126, -6, 2, 0}, { 0, 1, -2, 6, 127, -5, 1, 0},
164 : { 0, 1, -2, 4, 127, -3, 1, 0}, { 0, 0, 0, 2, 127, -1, 0, 0},
165 : // [1, 2)
166 : { 0, 0, 0, 1, 127, 0, 0, 0 }, { 0, 0, 0, - 1, 127, 2, 0, 0 },
167 : { 0, 0, 1, - 3, 127, 4, - 1, 0 }, { 0, 0, 1, - 4, 126, 6, - 2, 1 },
168 : { 0, 0, 1, - 5, 126, 8, - 3, 1 }, { 0, 0, 1, - 6, 125, 11, - 4, 1 },
169 : { 0, 0, 1, - 7, 124, 13, - 4, 1 }, { 0, 0, 2, - 8, 123, 15, - 5, 1 },
170 : { 0, 0, 2, - 9, 122, 18, - 6, 1 }, { 0, 0, 2, -10, 121, 20, - 6, 1 },
171 : { 0, 0, 2, -11, 120, 22, - 7, 2 }, { 0, 0, 2, -12, 119, 25, - 8, 2 },
172 : { 0, 0, 3, -13, 117, 27, - 8, 2 }, { 0, 0, 3, -13, 116, 29, - 9, 2 },
173 : { 0, 0, 3, -14, 114, 32, -10, 3 }, { 0, 0, 3, -15, 113, 35, -10, 2 },
174 : { 0, 0, 3, -15, 111, 37, -11, 3 }, { 0, 0, 3, -16, 109, 40, -11, 3 },
175 : { 0, 0, 3, -16, 108, 42, -12, 3 }, { 0, 0, 4, -17, 106, 45, -13, 3 },
176 : { 0, 0, 4, -17, 104, 47, -13, 3 }, { 0, 0, 4, -17, 102, 50, -14, 3 },
177 : { 0, 0, 4, -17, 100, 52, -14, 3 }, { 0, 0, 4, -18, 98, 55, -15, 4 },
178 : { 0, 0, 4, -18, 96, 58, -15, 3 }, { 0, 0, 4, -18, 94, 60, -16, 4 },
179 : { 0, 0, 4, -18, 91, 63, -16, 4 }, { 0, 0, 4, -18, 89, 65, -16, 4 },
180 : { 0, 0, 4, -18, 87, 68, -17, 4 }, { 0, 0, 4, -18, 85, 70, -17, 4 },
181 : { 0, 0, 4, -18, 82, 73, -17, 4 }, { 0, 0, 4, -18, 80, 75, -17, 4 },
182 : { 0, 0, 4, -18, 78, 78, -18, 4 }, { 0, 0, 4, -17, 75, 80, -18, 4 },
183 : { 0, 0, 4, -17, 73, 82, -18, 4 }, { 0, 0, 4, -17, 70, 85, -18, 4 },
184 : { 0, 0, 4, -17, 68, 87, -18, 4 }, { 0, 0, 4, -16, 65, 89, -18, 4 },
185 : { 0, 0, 4, -16, 63, 91, -18, 4 }, { 0, 0, 4, -16, 60, 94, -18, 4 },
186 : { 0, 0, 3, -15, 58, 96, -18, 4 }, { 0, 0, 4, -15, 55, 98, -18, 4 },
187 : { 0, 0, 3, -14, 52, 100, -17, 4 }, { 0, 0, 3, -14, 50, 102, -17, 4 },
188 : { 0, 0, 3, -13, 47, 104, -17, 4 }, { 0, 0, 3, -13, 45, 106, -17, 4 },
189 : { 0, 0, 3, -12, 42, 108, -16, 3 }, { 0, 0, 3, -11, 40, 109, -16, 3 },
190 : { 0, 0, 3, -11, 37, 111, -15, 3 }, { 0, 0, 2, -10, 35, 113, -15, 3 },
191 : { 0, 0, 3, -10, 32, 114, -14, 3 }, { 0, 0, 2, - 9, 29, 116, -13, 3 },
192 : { 0, 0, 2, - 8, 27, 117, -13, 3 }, { 0, 0, 2, - 8, 25, 119, -12, 2 },
193 : { 0, 0, 2, - 7, 22, 120, -11, 2 }, { 0, 0, 1, - 6, 20, 121, -10, 2 },
194 : { 0, 0, 1, - 6, 18, 122, - 9, 2 }, { 0, 0, 1, - 5, 15, 123, - 8, 2 },
195 : { 0, 0, 1, - 4, 13, 124, - 7, 1 }, { 0, 0, 1, - 4, 11, 125, - 6, 1 },
196 : { 0, 0, 1, - 3, 8, 126, - 5, 1 }, { 0, 0, 1, - 2, 6, 126, - 4, 1 },
197 : { 0, 0, 0, - 1, 4, 127, - 3, 1 }, { 0, 0, 0, 0, 2, 127, - 1, 0 },
198 : // dummy (replicate row index 191)
199 : { 0, 0, 0, 0, 2, 127, - 1, 0 },
200 : #elif WARPEDPIXEL_PREC_BITS == 5
201 : // [-1, 0)
202 : {0, 0, 127, 1, 0, 0, 0, 0}, {1, -3, 127, 4, -1, 0, 0, 0},
203 : {1, -5, 126, 8, -3, 1, 0, 0}, {1, -7, 124, 13, -4, 1, 0, 0},
204 : {2, -9, 122, 18, -6, 1, 0, 0}, {2, -11, 120, 22, -7, 2, 0, 0},
205 : {3, -13, 117, 27, -8, 2, 0, 0}, {3, -14, 114, 32, -10, 3, 0, 0},
206 : {3, -15, 111, 37, -11, 3, 0, 0}, {3, -16, 108, 42, -12, 3, 0, 0},
207 : {4, -17, 104, 47, -13, 3, 0, 0}, {4, -17, 100, 52, -14, 3, 0, 0},
208 : {4, -18, 96, 58, -15, 3, 0, 0}, {4, -18, 91, 63, -16, 4, 0, 0},
209 : {4, -18, 87, 68, -17, 4, 0, 0}, {4, -18, 82, 73, -17, 4, 0, 0},
210 : {4, -18, 78, 78, -18, 4, 0, 0}, {4, -17, 73, 82, -18, 4, 0, 0},
211 : {4, -17, 68, 87, -18, 4, 0, 0}, {4, -16, 63, 91, -18, 4, 0, 0},
212 : {3, -15, 58, 96, -18, 4, 0, 0}, {3, -14, 52, 100, -17, 4, 0, 0},
213 : {3, -13, 47, 104, -17, 4, 0, 0}, {3, -12, 42, 108, -16, 3, 0, 0},
214 : {3, -11, 37, 111, -15, 3, 0, 0}, {3, -10, 32, 114, -14, 3, 0, 0},
215 : {2, -8, 27, 117, -13, 3, 0, 0}, {2, -7, 22, 120, -11, 2, 0, 0},
216 : {1, -6, 18, 122, -9, 2, 0, 0}, {1, -4, 13, 124, -7, 1, 0, 0},
217 : {1, -3, 8, 126, -5, 1, 0, 0}, {0, -1, 4, 127, -3, 1, 0, 0},
218 : // [0, 1)
219 : { 0, 0, 0, 127, 1, 0, 0, 0}, { 0, 1, -3, 127, 4, -2, 1, 0},
220 : { 0, 2, -6, 126, 8, -3, 1, 0}, {-1, 3, -8, 125, 13, -5, 2, -1},
221 : {-1, 4, -11, 123, 18, -7, 3, -1}, {-1, 4, -13, 121, 23, -8, 3, -1},
222 : {-1, 5, -15, 119, 27, -10, 4, -1}, {-2, 6, -17, 116, 33, -12, 5, -1},
223 : {-2, 6, -18, 113, 38, -13, 5, -1}, {-2, 7, -19, 110, 43, -15, 6, -2},
224 : {-2, 7, -20, 106, 49, -16, 6, -2}, {-2, 7, -21, 102, 54, -17, 7, -2},
225 : {-2, 8, -22, 98, 59, -18, 7, -2}, {-2, 8, -22, 94, 64, -19, 7, -2},
226 : {-2, 8, -22, 89, 69, -20, 8, -2}, {-2, 8, -21, 84, 74, -21, 8, -2},
227 : {-2, 8, -21, 79, 79, -21, 8, -2}, {-2, 8, -21, 74, 84, -21, 8, -2},
228 : {-2, 8, -20, 69, 89, -22, 8, -2}, {-2, 7, -19, 64, 94, -22, 8, -2},
229 : {-2, 7, -18, 59, 98, -22, 8, -2}, {-2, 7, -17, 54, 102, -21, 7, -2},
230 : {-2, 6, -16, 49, 106, -20, 7, -2}, {-2, 6, -15, 43, 110, -19, 7, -2},
231 : {-1, 5, -13, 38, 113, -18, 6, -2}, {-1, 5, -12, 33, 116, -17, 6, -2},
232 : {-1, 4, -10, 27, 119, -15, 5, -1}, {-1, 3, -8, 23, 121, -13, 4, -1},
233 : {-1, 3, -7, 18, 123, -11, 4, -1}, {-1, 2, -5, 13, 125, -8, 3, -1},
234 : { 0, 1, -3, 8, 126, -6, 2, 0}, { 0, 1, -2, 4, 127, -3, 1, 0},
235 : // [1, 2)
236 : {0, 0, 0, 1, 127, 0, 0, 0}, {0, 0, 1, -3, 127, 4, -1, 0},
237 : {0, 0, 1, -5, 126, 8, -3, 1}, {0, 0, 1, -7, 124, 13, -4, 1},
238 : {0, 0, 2, -9, 122, 18, -6, 1}, {0, 0, 2, -11, 120, 22, -7, 2},
239 : {0, 0, 3, -13, 117, 27, -8, 2}, {0, 0, 3, -14, 114, 32, -10, 3},
240 : {0, 0, 3, -15, 111, 37, -11, 3}, {0, 0, 3, -16, 108, 42, -12, 3},
241 : {0, 0, 4, -17, 104, 47, -13, 3}, {0, 0, 4, -17, 100, 52, -14, 3},
242 : {0, 0, 4, -18, 96, 58, -15, 3}, {0, 0, 4, -18, 91, 63, -16, 4},
243 : {0, 0, 4, -18, 87, 68, -17, 4}, {0, 0, 4, -18, 82, 73, -17, 4},
244 : {0, 0, 4, -18, 78, 78, -18, 4}, {0, 0, 4, -17, 73, 82, -18, 4},
245 : {0, 0, 4, -17, 68, 87, -18, 4}, {0, 0, 4, -16, 63, 91, -18, 4},
246 : {0, 0, 3, -15, 58, 96, -18, 4}, {0, 0, 3, -14, 52, 100, -17, 4},
247 : {0, 0, 3, -13, 47, 104, -17, 4}, {0, 0, 3, -12, 42, 108, -16, 3},
248 : {0, 0, 3, -11, 37, 111, -15, 3}, {0, 0, 3, -10, 32, 114, -14, 3},
249 : {0, 0, 2, -8, 27, 117, -13, 3}, {0, 0, 2, -7, 22, 120, -11, 2},
250 : {0, 0, 1, -6, 18, 122, -9, 2}, {0, 0, 1, -4, 13, 124, -7, 1},
251 : {0, 0, 1, -3, 8, 126, -5, 1}, {0, 0, 0, -1, 4, 127, -3, 1},
252 : // dummy (replicate row index 95)
253 : {0, 0, 0, -1, 4, 127, -3, 1},
254 : #endif // WARPEDPIXEL_PREC_BITS == 6
255 : };
256 :
257 : /* clang-format on */
258 :
259 : #define DIV_LUT_PREC_BITS 14
260 : #define DIV_LUT_BITS 8
261 : #define DIV_LUT_NUM (1 << DIV_LUT_BITS)
262 :
263 : static const uint16_t div_lut[DIV_LUT_NUM + 1] = {
264 : 16384, 16320, 16257, 16194, 16132, 16070, 16009, 15948, 15888, 15828, 15768,
265 : 15709, 15650, 15592, 15534, 15477, 15420, 15364, 15308, 15252, 15197, 15142,
266 : 15087, 15033, 14980, 14926, 14873, 14821, 14769, 14717, 14665, 14614, 14564,
267 : 14513, 14463, 14413, 14364, 14315, 14266, 14218, 14170, 14122, 14075, 14028,
268 : 13981, 13935, 13888, 13843, 13797, 13752, 13707, 13662, 13618, 13574, 13530,
269 : 13487, 13443, 13400, 13358, 13315, 13273, 13231, 13190, 13148, 13107, 13066,
270 : 13026, 12985, 12945, 12906, 12866, 12827, 12788, 12749, 12710, 12672, 12633,
271 : 12596, 12558, 12520, 12483, 12446, 12409, 12373, 12336, 12300, 12264, 12228,
272 : 12193, 12157, 12122, 12087, 12053, 12018, 11984, 11950, 11916, 11882, 11848,
273 : 11815, 11782, 11749, 11716, 11683, 11651, 11619, 11586, 11555, 11523, 11491,
274 : 11460, 11429, 11398, 11367, 11336, 11305, 11275, 11245, 11215, 11185, 11155,
275 : 11125, 11096, 11067, 11038, 11009, 10980, 10951, 10923, 10894, 10866, 10838,
276 : 10810, 10782, 10755, 10727, 10700, 10673, 10645, 10618, 10592, 10565, 10538,
277 : 10512, 10486, 10460, 10434, 10408, 10382, 10356, 10331, 10305, 10280, 10255,
278 : 10230, 10205, 10180, 10156, 10131, 10107, 10082, 10058, 10034, 10010, 9986,
279 : 9963, 9939, 9916, 9892, 9869, 9846, 9823, 9800, 9777, 9754, 9732,
280 : 9709, 9687, 9664, 9642, 9620, 9598, 9576, 9554, 9533, 9511, 9489,
281 : 9468, 9447, 9425, 9404, 9383, 9362, 9341, 9321, 9300, 9279, 9259,
282 : 9239, 9218, 9198, 9178, 9158, 9138, 9118, 9098, 9079, 9059, 9039,
283 : 9020, 9001, 8981, 8962, 8943, 8924, 8905, 8886, 8867, 8849, 8830,
284 : 8812, 8793, 8775, 8756, 8738, 8720, 8702, 8684, 8666, 8648, 8630,
285 : 8613, 8595, 8577, 8560, 8542, 8525, 8508, 8490, 8473, 8456, 8439,
286 : 8422, 8405, 8389, 8372, 8355, 8339, 8322, 8306, 8289, 8273, 8257,
287 : 8240, 8224, 8208, 8192,
288 : };
289 :
290 : // Decomposes a divisor D such that 1/D = y/2^shift, where y is returned
291 : // at precision of DIV_LUT_PREC_BITS along with the shift.
292 737783 : static int16_t resolve_divisor_64(uint64_t D, int16_t *shift) {
293 : int64_t f;
294 1464570 : *shift = (int16_t)((D >> 32) ? get_msb((unsigned int)(D >> 32)) + 32
295 726790 : : get_msb((unsigned int)D));
296 : // e is obtained from D after resetting the most significant 1 bit.
297 737782 : const int64_t e = D - ((uint64_t)1 << *shift);
298 : // Get the most significant DIV_LUT_BITS (8) bits of e into f
299 737782 : if (*shift > DIV_LUT_BITS)
300 737781 : f = ROUND_POWER_OF_TWO_64(e, *shift - DIV_LUT_BITS);
301 : else
302 1 : f = e << (DIV_LUT_BITS - *shift);
303 737782 : assert(f <= DIV_LUT_NUM);
304 737782 : *shift += DIV_LUT_PREC_BITS;
305 : // Use f as lookup into the precomputed table of multipliers
306 737782 : return div_lut[f];
307 : }
308 :
309 829496 : static int16_t resolve_divisor_32(uint32_t D, int16_t *shift) {
310 : int32_t f;
311 829496 : *shift = get_msb(D);
312 : // e is obtained from D after resetting the most significant 1 bit.
313 829495 : const int32_t e = D - ((uint32_t)1 << *shift);
314 : // Get the most significant DIV_LUT_BITS (8) bits of e into f
315 829495 : if (*shift > DIV_LUT_BITS)
316 829496 : f = ROUND_POWER_OF_TWO(e, *shift - DIV_LUT_BITS);
317 : else
318 0 : f = e << (DIV_LUT_BITS - *shift);
319 829495 : assert(f <= DIV_LUT_NUM);
320 829495 : *shift += DIV_LUT_PREC_BITS;
321 : // Use f as lookup into the precomputed table of multipliers
322 829495 : return div_lut[f];
323 : }
324 :
325 829495 : static int is_affine_valid(const EbWarpedMotionParams *const wm) {
326 829495 : const int32_t *mat = wm->wmmat;
327 829495 : return (mat[2] > 0);
328 : }
329 :
330 829487 : static int is_affine_shear_allowed(int16_t alpha, int16_t beta, int16_t gamma,
331 : int16_t delta) {
332 829487 : if ((4 * abs(alpha) + 7 * abs(beta) >= (1 << WARPEDMODEL_PREC_BITS)) ||
333 782617 : (4 * abs(gamma) + 4 * abs(delta) >= (1 << WARPEDMODEL_PREC_BITS)))
334 50470 : return 0;
335 : else
336 779017 : return 1;
337 : }
338 :
339 : // Returns 1 on success or 0 on an invalid affine set
340 829495 : int eb_get_shear_params(EbWarpedMotionParams *wm) {
341 829495 : const int32_t *mat = wm->wmmat;
342 829495 : if (!is_affine_valid(wm)) return 0;
343 829494 : wm->alpha =
344 829495 : clamp(mat[2] - (1 << WARPEDMODEL_PREC_BITS), INT16_MIN, INT16_MAX);
345 829494 : wm->beta = clamp(mat[3], INT16_MIN, INT16_MAX);
346 : int16_t shift;
347 829495 : int16_t y = resolve_divisor_32(abs(mat[2]), &shift) * (mat[2] < 0 ? -1 : 1);
348 829497 : int64_t v = ((int64_t)mat[4] * (1 << WARPEDMODEL_PREC_BITS)) * y;
349 829487 : wm->gamma =
350 829497 : clamp((int)ROUND_POWER_OF_TWO_SIGNED_64(v, shift), INT16_MIN, INT16_MAX);
351 829487 : v = ((int64_t)mat[3] * mat[4]) * y;
352 829487 : wm->delta = clamp(mat[5] - (int)ROUND_POWER_OF_TWO_SIGNED_64(v, shift) -
353 : (1 << WARPEDMODEL_PREC_BITS),
354 : INT16_MIN, INT16_MAX);
355 :
356 829482 : wm->alpha = ROUND_POWER_OF_TWO_SIGNED(wm->alpha, WARP_PARAM_REDUCE_BITS) *
357 : (1 << WARP_PARAM_REDUCE_BITS);
358 829482 : wm->beta = ROUND_POWER_OF_TWO_SIGNED(wm->beta, WARP_PARAM_REDUCE_BITS) *
359 : (1 << WARP_PARAM_REDUCE_BITS);
360 829482 : wm->gamma = ROUND_POWER_OF_TWO_SIGNED(wm->gamma, WARP_PARAM_REDUCE_BITS) *
361 : (1 << WARP_PARAM_REDUCE_BITS);
362 829482 : wm->delta = ROUND_POWER_OF_TWO_SIGNED(wm->delta, WARP_PARAM_REDUCE_BITS) *
363 : (1 << WARP_PARAM_REDUCE_BITS);
364 :
365 829482 : if (!is_affine_shear_allowed(wm->alpha, wm->beta, wm->gamma, wm->delta))
366 50471 : return 0;
367 :
368 779013 : return 1;
369 : }
370 :
371 0 : static INLINE int highbd_error_measure(int err, int bd) {
372 0 : const int b = bd - 8;
373 0 : const int bmask = (1 << b) - 1;
374 0 : const int v = (1 << b);
375 0 : err = abs(err);
376 0 : const int e1 = err >> b;
377 0 : const int e2 = err & bmask;
378 0 : return error_measure_lut[255 + e1] * (v - e2) +
379 0 : error_measure_lut[256 + e1] * e2;
380 : }
381 :
382 : /* Note: For an explanation of the warp algorithm, and some notes on bit widths
383 : for hardware implementations, see the comments above eb_av1_warp_affine_c
384 : */
385 0 : void eb_av1_highbd_warp_affine_c(const int32_t *mat, const uint16_t *ref,
386 : int width, int height, int stride, uint16_t *pred,
387 : int p_col, int p_row, int p_width, int p_height,
388 : int p_stride, int subsampling_x,
389 : int subsampling_y, int bd,
390 : ConvolveParams *conv_params, int16_t alpha,
391 : int16_t beta, int16_t gamma, int16_t delta) {
392 : int32_t tmp[15 * 8];
393 0 : const int reduce_bits_horiz =
394 0 : conv_params->round_0 +
395 0 : AOMMAX(bd + FILTER_BITS - conv_params->round_0 - 14, 0);
396 0 : const int reduce_bits_vert = conv_params->is_compound
397 : ? conv_params->round_1
398 0 : : 2 * FILTER_BITS - reduce_bits_horiz;
399 0 : const int max_bits_horiz = bd + FILTER_BITS + 1 - reduce_bits_horiz;
400 0 : const int offset_bits_horiz = bd + FILTER_BITS - 1;
401 0 : const int offset_bits_vert = bd + 2 * FILTER_BITS - reduce_bits_horiz;
402 0 : const int round_bits =
403 0 : 2 * FILTER_BITS - conv_params->round_0 - conv_params->round_1;
404 0 : const int offset_bits = bd + 2 * FILTER_BITS - conv_params->round_0;
405 : (void)max_bits_horiz;
406 0 : assert(IMPLIES(conv_params->is_compound, conv_params->dst != NULL));
407 :
408 0 : for (int i = p_row; i < p_row + p_height; i += 8) {
409 0 : for (int j = p_col; j < p_col + p_width; j += 8) {
410 : // Calculate the center of this 8x8 block,
411 : // project to luma coordinates (if in a subsampled chroma plane),
412 : // apply the affine transformation,
413 : // then convert back to the original coordinates (if necessary)
414 0 : const int32_t src_x = (j + 4) << subsampling_x;
415 0 : const int32_t src_y = (i + 4) << subsampling_y;
416 0 : const int32_t dst_x = mat[2] * src_x + mat[3] * src_y + mat[0];
417 0 : const int32_t dst_y = mat[4] * src_x + mat[5] * src_y + mat[1];
418 0 : const int32_t x4 = dst_x >> subsampling_x;
419 0 : const int32_t y4 = dst_y >> subsampling_y;
420 :
421 0 : const int32_t ix4 = x4 >> WARPEDMODEL_PREC_BITS;
422 0 : int32_t sx4 = x4 & ((1 << WARPEDMODEL_PREC_BITS) - 1);
423 0 : const int32_t iy4 = y4 >> WARPEDMODEL_PREC_BITS;
424 0 : int32_t sy4 = y4 & ((1 << WARPEDMODEL_PREC_BITS) - 1);
425 :
426 0 : sx4 += alpha * (-4) + beta * (-4);
427 0 : sy4 += gamma * (-4) + delta * (-4);
428 :
429 0 : sx4 &= ~((1 << WARP_PARAM_REDUCE_BITS) - 1);
430 0 : sy4 &= ~((1 << WARP_PARAM_REDUCE_BITS) - 1);
431 :
432 : // Horizontal filter
433 0 : for (int k = -7; k < 8; ++k) {
434 0 : const int iy = clamp(iy4 + k, 0, height - 1);
435 :
436 0 : int sx = sx4 + beta * (k + 4);
437 0 : for (int l = -4; l < 4; ++l) {
438 0 : int ix = ix4 + l - 3;
439 0 : const int offs = ROUND_POWER_OF_TWO(sx, WARPEDDIFF_PREC_BITS) +
440 : WARPEDPIXEL_PREC_SHIFTS;
441 0 : assert(offs >= 0 && offs <= WARPEDPIXEL_PREC_SHIFTS * 3);
442 0 : const int16_t *coeffs = eb_warped_filter[offs];
443 :
444 0 : int32_t sum = 1 << offset_bits_horiz;
445 0 : for (int m = 0; m < 8; ++m) {
446 0 : const int sample_x = clamp(ix + m, 0, width - 1);
447 0 : sum += ref[iy * stride + sample_x] * coeffs[m];
448 : }
449 0 : sum = ROUND_POWER_OF_TWO(sum, reduce_bits_horiz);
450 0 : assert(0 <= sum && sum < (1 << max_bits_horiz));
451 0 : tmp[(k + 7) * 8 + (l + 4)] = sum;
452 0 : sx += alpha;
453 : }
454 : }
455 :
456 : // Vertical filter
457 0 : for (int k = -4; k < AOMMIN(4, p_row + p_height - i - 4); ++k) {
458 0 : int sy = sy4 + delta * (k + 4);
459 0 : for (int l = -4; l < AOMMIN(4, p_col + p_width - j - 4); ++l) {
460 0 : const int offs = ROUND_POWER_OF_TWO(sy, WARPEDDIFF_PREC_BITS) +
461 : WARPEDPIXEL_PREC_SHIFTS;
462 0 : assert(offs >= 0 && offs <= WARPEDPIXEL_PREC_SHIFTS * 3);
463 0 : const int16_t *coeffs = eb_warped_filter[offs];
464 :
465 0 : int32_t sum = 1 << offset_bits_vert;
466 0 : for (int m = 0; m < 8; ++m)
467 0 : sum += tmp[(k + m + 4) * 8 + (l + 4)] * coeffs[m];
468 0 : if (conv_params->is_compound) {
469 0 : ConvBufType *p =
470 : &conv_params
471 0 : ->dst[(i - p_row + k + 4) * conv_params->dst_stride +
472 0 : (j - p_col + l + 4)];
473 0 : sum = ROUND_POWER_OF_TWO(sum, reduce_bits_vert);
474 0 : if (conv_params->do_average) {
475 0 : uint16_t *dst16 =
476 0 : &pred[(i - p_row + k + 4) * p_stride + (j - p_col + l + 4)];
477 0 : int32_t tmp32 = *p;
478 0 : if (conv_params->use_jnt_comp_avg) {
479 0 : tmp32 = tmp32 * conv_params->fwd_offset +
480 0 : sum * conv_params->bck_offset;
481 0 : tmp32 = tmp32 >> DIST_PRECISION_BITS;
482 : } else {
483 0 : tmp32 += sum;
484 0 : tmp32 = tmp32 >> 1;
485 : }
486 0 : tmp32 = tmp32 - (1 << (offset_bits - conv_params->round_1)) -
487 0 : (1 << (offset_bits - conv_params->round_1 - 1));
488 0 : *dst16 =
489 0 : clip_pixel_highbd(ROUND_POWER_OF_TWO(tmp32, round_bits), bd);
490 : } else
491 0 : *p = sum;
492 : } else {
493 0 : uint16_t *p =
494 0 : &pred[(i - p_row + k + 4) * p_stride + (j - p_col + l + 4)];
495 0 : sum = ROUND_POWER_OF_TWO(sum, reduce_bits_vert);
496 0 : assert(0 <= sum && sum < (1 << (bd + 2)));
497 0 : *p = clip_pixel_highbd(sum - (1 << (bd - 1)) - (1 << bd), bd);
498 : }
499 0 : sy += gamma;
500 : }
501 : }
502 : }
503 : }
504 0 : }
505 :
506 0 : static void highbd_warp_plane(EbWarpedMotionParams *wm, const uint8_t *const ref8,
507 : int width, int height, int stride,
508 : const uint8_t *const pred8, int p_col, int p_row,
509 : int p_width, int p_height, int p_stride,
510 : int subsampling_x, int subsampling_y, int bd,
511 : ConvolveParams *conv_params) {
512 0 : assert(wm->wmtype <= AFFINE);
513 0 : if (wm->wmtype == ROTZOOM) {
514 0 : wm->wmmat[5] = wm->wmmat[2];
515 0 : wm->wmmat[4] = -wm->wmmat[3];
516 : }
517 0 : const int32_t *const mat = wm->wmmat;
518 0 : const int16_t alpha = wm->alpha;
519 0 : const int16_t beta = wm->beta;
520 0 : const int16_t gamma = wm->gamma;
521 0 : const int16_t delta = wm->delta;
522 :
523 0 : const uint16_t *const ref = (uint16_t *)ref8;
524 0 : uint16_t *pred = (uint16_t *)pred8;
525 0 : eb_av1_highbd_warp_affine_c(mat, ref, width, height, stride, pred, p_col, p_row,
526 : p_width, p_height, p_stride, subsampling_x,
527 : subsampling_y, bd, conv_params, alpha, beta, gamma,
528 : delta);
529 0 : }
530 :
531 0 : static int64_t highbd_frame_error(const uint16_t *const ref, int stride,
532 : const uint16_t *const dst, int p_width,
533 : int p_height, int p_stride, int bd) {
534 0 : int64_t sum_error = 0;
535 0 : for (int i = 0; i < p_height; ++i) {
536 0 : for (int j = 0; j < p_width; ++j) {
537 0 : sum_error +=
538 0 : highbd_error_measure(dst[j + i * p_stride] - ref[j + i * stride], bd);
539 : }
540 : }
541 0 : return sum_error;
542 : }
543 :
544 0 : static int64_t highbd_warp_error(
545 : EbWarpedMotionParams *wm, const uint8_t *const ref8, int width, int height,
546 : int stride, const uint8_t *const dst8, int p_col, int p_row, int p_width,
547 : int p_height, int p_stride, int subsampling_x, int subsampling_y, int bd,
548 : int64_t best_error) {
549 0 : int64_t gm_sumerr = 0;
550 0 : const int error_bsize_w = AOMMIN(p_width, WARP_ERROR_BLOCK);
551 0 : const int error_bsize_h = AOMMIN(p_height, WARP_ERROR_BLOCK);
552 : uint16_t tmp[WARP_ERROR_BLOCK * WARP_ERROR_BLOCK];
553 :
554 0 : ConvolveParams conv_params = get_conv_params(0, 0, 0, bd);
555 0 : conv_params.use_jnt_comp_avg = 0;
556 0 : for (int i = p_row; i < p_row + p_height; i += WARP_ERROR_BLOCK) {
557 0 : for (int j = p_col; j < p_col + p_width; j += WARP_ERROR_BLOCK) {
558 : // avoid warping extra 8x8 blocks in the padded region of the frame
559 : // when p_width and p_height are not multiples of WARP_ERROR_BLOCK
560 0 : const int warp_w = AOMMIN(error_bsize_w, p_col + p_width - j);
561 0 : const int warp_h = AOMMIN(error_bsize_h, p_row + p_height - i);
562 0 : highbd_warp_plane(wm, ref8, width, height, stride,
563 0 : CONVERT_TO_BYTEPTR(tmp), j, i, warp_w, warp_h,
564 : WARP_ERROR_BLOCK, subsampling_x, subsampling_y, bd,
565 : &conv_params);
566 :
567 0 : gm_sumerr += highbd_frame_error(
568 0 : tmp, WARP_ERROR_BLOCK, CONVERT_TO_SHORTPTR(dst8) + j + i * p_stride,
569 : warp_w, warp_h, p_stride, bd);
570 0 : if (gm_sumerr > best_error) return gm_sumerr;
571 : }
572 : }
573 0 : return gm_sumerr;
574 : }
575 :
576 : /* The warp filter for ROTZOOM and AFFINE models works as follows:
577 : * Split the input into 8x8 blocks
578 : * For each block, project the point (4, 4) within the block, to get the
579 : overall block position. Split into integer and fractional coordinates,
580 : maintaining full WARPEDMODEL precision
581 : * Filter horizontally: Generate 15 rows of 8 pixels each. Each pixel gets a
582 : variable horizontal offset. This means that, while the rows of the
583 : intermediate buffer align with the rows of the *reference* image, the
584 : columns align with the columns of the *destination* image.
585 : * Filter vertically: Generate the output block (up to 8x8 pixels, but if the
586 : destination is too small we crop the output at this stage). Each pixel has
587 : a variable vertical offset, so that the resulting rows are aligned with
588 : the rows of the destination image.
589 :
590 : To accomplish these alignments, we factor the warp matrix as a
591 : product of two shear / asymmetric zoom matrices:
592 : / a b \ = / 1 0 \ * / 1+alpha beta \
593 : \ c d / \ gamma 1+delta / \ 0 1 /
594 : where a, b, c, d are wmmat[2], wmmat[3], wmmat[4], wmmat[5] respectively.
595 : The horizontal shear (with alpha and beta) is applied first,
596 : then the vertical shear (with gamma and delta) is applied second.
597 :
598 : The only limitation is that, to fit this in a fixed 8-tap filter size,
599 : the fractional pixel offsets must be at most +-1. Since the horizontal filter
600 : generates 15 rows of 8 columns, and the initial point we project is at (4, 4)
601 : within the block, the parameters must satisfy
602 : 4 * |alpha| + 7 * |beta| <= 1 and 4 * |gamma| + 4 * |delta| <= 1
603 : for this filter to be applicable.
604 :
605 : Note: This function assumes that the caller has done all of the relevant
606 : checks, ie. that we have a ROTZOOM or AFFINE model, that wm[4] and wm[5]
607 : are set appropriately (if using a ROTZOOM model), and that alpha, beta,
608 : gamma, delta are all in range.
609 :
610 : TODO(david.barker): Maybe support scaled references?
611 : */
612 : /* A note on hardware implementation:
613 : The warp filter is intended to be implementable using the same hardware as
614 : the high-precision convolve filters from the loop-restoration and
615 : convolve-round experiments.
616 :
617 : For a single filter stage, considering all of the coefficient sets for the
618 : warp filter and the regular convolution filter, an input in the range
619 : [0, 2^k - 1] is mapped into the range [-56 * (2^k - 1), 184 * (2^k - 1)]
620 : before rounding.
621 :
622 : Allowing for some changes to the filter coefficient sets, call the range
623 : [-64 * 2^k, 192 * 2^k]. Then, if we initialize the accumulator to 64 * 2^k,
624 : we can replace this by the range [0, 256 * 2^k], which can be stored in an
625 : unsigned value with 8 + k bits.
626 :
627 : This allows the derivation of the appropriate bit widths and offsets for
628 : the various intermediate values: If
629 :
630 : F := FILTER_BITS = 7 (or else the above ranges need adjusting)
631 : So a *single* filter stage maps a k-bit input to a (k + F + 1)-bit
632 : intermediate value.
633 : H := ROUND0_BITS
634 : V := VERSHEAR_REDUCE_PREC_BITS
635 : (and note that we must have H + V = 2*F for the output to have the same
636 : scale as the input)
637 :
638 : then we end up with the following offsets and ranges:
639 : Horizontal filter: Apply an offset of 1 << (bd + F - 1), sum fits into a
640 : uint{bd + F + 1}
641 : After rounding: The values stored in 'tmp' fit into a uint{bd + F + 1 - H}.
642 : Vertical filter: Apply an offset of 1 << (bd + 2*F - H), sum fits into a
643 : uint{bd + 2*F + 2 - H}
644 : After rounding: The final value, before undoing the offset, fits into a
645 : uint{bd + 2}.
646 :
647 : Then we need to undo the offsets before clamping to a pixel. Note that,
648 : if we do this at the end, the amount to subtract is actually independent
649 : of H and V:
650 :
651 : offset to subtract = (1 << ((bd + F - 1) - H + F - V)) +
652 : (1 << ((bd + 2*F - H) - V))
653 : == (1 << (bd - 1)) + (1 << bd)
654 :
655 : This allows us to entirely avoid clamping in both the warp filter and
656 : the convolve-round experiment. As of the time of writing, the Wiener filter
657 : from loop-restoration can encode a central coefficient up to 216, which
658 : leads to a maximum value of about 282 * 2^k after applying the offset.
659 : So in that case we still need to clamp.
660 : */
661 0 : void eb_av1_warp_affine_c(const int32_t *mat, const uint8_t *ref, int width,
662 : int height, int stride, uint8_t *pred, int p_col,
663 : int p_row, int p_width, int p_height, int p_stride,
664 : int subsampling_x, int subsampling_y,
665 : ConvolveParams *conv_params, int16_t alpha, int16_t beta,
666 : int16_t gamma, int16_t delta) {
667 : int32_t tmp[15 * 8];
668 0 : const int bd = 8;
669 0 : const int reduce_bits_horiz = conv_params->round_0;
670 0 : const int reduce_bits_vert = conv_params->is_compound
671 : ? conv_params->round_1
672 0 : : 2 * FILTER_BITS - reduce_bits_horiz;
673 0 : const int max_bits_horiz = bd + FILTER_BITS + 1 - reduce_bits_horiz;
674 0 : const int offset_bits_horiz = bd + FILTER_BITS - 1;
675 0 : const int offset_bits_vert = bd + 2 * FILTER_BITS - reduce_bits_horiz;
676 0 : const int round_bits =
677 0 : 2 * FILTER_BITS - conv_params->round_0 - conv_params->round_1;
678 0 : const int offset_bits = bd + 2 * FILTER_BITS - conv_params->round_0;
679 : (void)max_bits_horiz;
680 0 : assert(IMPLIES(conv_params->is_compound, conv_params->dst != NULL));
681 0 : assert(IMPLIES(conv_params->do_average, conv_params->is_compound));
682 :
683 0 : for (int i = p_row; i < p_row + p_height; i += 8) {
684 0 : for (int j = p_col; j < p_col + p_width; j += 8) {
685 : // Calculate the center of this 8x8 block,
686 : // project to luma coordinates (if in a subsampled chroma plane),
687 : // apply the affine transformation,
688 : // then convert back to the original coordinates (if necessary)
689 0 : const int32_t src_x = (j + 4) << subsampling_x;
690 0 : const int32_t src_y = (i + 4) << subsampling_y;
691 0 : const int32_t dst_x = mat[2] * src_x + mat[3] * src_y + mat[0];
692 0 : const int32_t dst_y = mat[4] * src_x + mat[5] * src_y + mat[1];
693 0 : const int32_t x4 = dst_x >> subsampling_x;
694 0 : const int32_t y4 = dst_y >> subsampling_y;
695 :
696 0 : int32_t ix4 = x4 >> WARPEDMODEL_PREC_BITS;
697 0 : int32_t sx4 = x4 & ((1 << WARPEDMODEL_PREC_BITS) - 1);
698 0 : int32_t iy4 = y4 >> WARPEDMODEL_PREC_BITS;
699 0 : int32_t sy4 = y4 & ((1 << WARPEDMODEL_PREC_BITS) - 1);
700 :
701 0 : sx4 += alpha * (-4) + beta * (-4);
702 0 : sy4 += gamma * (-4) + delta * (-4);
703 :
704 0 : sx4 &= ~((1 << WARP_PARAM_REDUCE_BITS) - 1);
705 0 : sy4 &= ~((1 << WARP_PARAM_REDUCE_BITS) - 1);
706 :
707 : // Horizontal filter
708 0 : for (int k = -7; k < 8; ++k) {
709 : // Clamp to top/bottom edge of the frame
710 0 : const int iy = clamp(iy4 + k, 0, height - 1);
711 :
712 0 : int sx = sx4 + beta * (k + 4);
713 :
714 0 : for (int l = -4; l < 4; ++l) {
715 0 : int ix = ix4 + l - 3;
716 : // At this point, sx = sx4 + alpha * l + beta * k
717 0 : const int offs = ROUND_POWER_OF_TWO(sx, WARPEDDIFF_PREC_BITS) +
718 : WARPEDPIXEL_PREC_SHIFTS;
719 0 : assert(offs >= 0 && offs <= WARPEDPIXEL_PREC_SHIFTS * 3);
720 0 : const int16_t *coeffs = eb_warped_filter[offs];
721 :
722 0 : int32_t sum = 1 << offset_bits_horiz;
723 0 : for (int m = 0; m < 8; ++m) {
724 : // Clamp to left/right edge of the frame
725 0 : const int sample_x = clamp(ix + m, 0, width - 1);
726 :
727 0 : sum += ref[iy * stride + sample_x] * coeffs[m];
728 : }
729 0 : sum = ROUND_POWER_OF_TWO(sum, reduce_bits_horiz);
730 0 : assert(0 <= sum && sum < (1 << max_bits_horiz));
731 0 : tmp[(k + 7) * 8 + (l + 4)] = sum;
732 0 : sx += alpha;
733 : }
734 : }
735 :
736 : // Vertical filter
737 0 : for (int k = -4; k < AOMMIN(4, p_row + p_height - i - 4); ++k) {
738 0 : int sy = sy4 + delta * (k + 4);
739 0 : for (int l = -4; l < AOMMIN(4, p_col + p_width - j - 4); ++l) {
740 : // At this point, sy = sy4 + gamma * l + delta * k
741 0 : const int offs = ROUND_POWER_OF_TWO(sy, WARPEDDIFF_PREC_BITS) +
742 : WARPEDPIXEL_PREC_SHIFTS;
743 0 : assert(offs >= 0 && offs <= WARPEDPIXEL_PREC_SHIFTS * 3);
744 0 : const int16_t *coeffs = eb_warped_filter[offs];
745 :
746 0 : int32_t sum = 1 << offset_bits_vert;
747 0 : for (int m = 0; m < 8; ++m)
748 0 : sum += tmp[(k + m + 4) * 8 + (l + 4)] * coeffs[m];
749 0 : if (conv_params->is_compound) {
750 0 : ConvBufType *p =
751 : &conv_params
752 0 : ->dst[(i - p_row + k + 4) * conv_params->dst_stride +
753 0 : (j - p_col + l + 4)];
754 0 : sum = ROUND_POWER_OF_TWO(sum, reduce_bits_vert);
755 0 : if (conv_params->do_average) {
756 0 : uint8_t *dst8 =
757 0 : &pred[(i - p_row + k + 4) * p_stride + (j - p_col + l + 4)];
758 0 : int32_t tmp32 = *p;
759 0 : if (conv_params->use_jnt_comp_avg) {
760 0 : tmp32 = tmp32 * conv_params->fwd_offset +
761 0 : sum * conv_params->bck_offset;
762 0 : tmp32 = tmp32 >> DIST_PRECISION_BITS;
763 : } else {
764 0 : tmp32 += sum;
765 0 : tmp32 = tmp32 >> 1;
766 : }
767 0 : tmp32 = tmp32 - (1 << (offset_bits - conv_params->round_1)) -
768 0 : (1 << (offset_bits - conv_params->round_1 - 1));
769 0 : *dst8 = clip_pixel(ROUND_POWER_OF_TWO(tmp32, round_bits));
770 : } else
771 0 : *p = sum;
772 : } else {
773 0 : uint8_t *p =
774 0 : &pred[(i - p_row + k + 4) * p_stride + (j - p_col + l + 4)];
775 0 : sum = ROUND_POWER_OF_TWO(sum, reduce_bits_vert);
776 0 : assert(0 <= sum && sum < (1 << (bd + 2)));
777 0 : *p = clip_pixel(sum - (1 << (bd - 1)) - (1 << bd));
778 : }
779 0 : sy += gamma;
780 : }
781 : }
782 : }
783 : }
784 0 : }
785 :
786 26608800 : static void warp_plane(EbWarpedMotionParams *wm, const uint8_t *const ref,
787 : int width, int height, int stride, uint8_t *pred,
788 : int p_col, int p_row, int p_width, int p_height,
789 : int p_stride, int subsampling_x, int subsampling_y,
790 : ConvolveParams *conv_params) {
791 26608800 : assert(wm->wmtype <= AFFINE);
792 26608800 : if (wm->wmtype == ROTZOOM) {
793 25138600 : wm->wmmat[5] = wm->wmmat[2];
794 25138600 : wm->wmmat[4] = -wm->wmmat[3];
795 : }
796 26608800 : const int32_t *const mat = wm->wmmat;
797 26608800 : const int16_t alpha = wm->alpha;
798 26608800 : const int16_t beta = wm->beta;
799 26608800 : const int16_t gamma = wm->gamma;
800 26608800 : const int16_t delta = wm->delta;
801 26608800 : eb_av1_warp_affine(mat, ref, width, height, stride, pred, p_col, p_row, p_width,
802 : p_height, p_stride, subsampling_x, subsampling_y, conv_params,
803 : alpha, beta, gamma, delta);
804 26608400 : }
805 :
806 0 : int64_t eb_av1_calc_frame_error_c(const uint8_t *const ref, int stride,
807 : const uint8_t *const dst, int p_width, int p_height, int p_stride) {
808 0 : int64_t sum_error = 0;
809 0 : for (int i = 0; i < p_height; ++i) {
810 0 : for (int j = 0; j < p_width; ++j) {
811 0 : sum_error +=
812 0 : (int64_t)error_measure(dst[j + i * p_stride] - ref[j + i * stride]);
813 : }
814 : }
815 0 : return sum_error;
816 : }
817 :
818 91233 : static int64_t warp_error(EbWarpedMotionParams *wm, const uint8_t *const ref,
819 : int width, int height, int stride,
820 : const uint8_t *const dst, int p_col, int p_row,
821 : int p_width, int p_height, int p_stride,
822 : int subsampling_x, int subsampling_y,
823 : int64_t best_error) {
824 91233 : int64_t gm_sumerr = 0;
825 : int warp_w, warp_h;
826 91233 : int error_bsize_w = AOMMIN(p_width, WARP_ERROR_BLOCK);
827 91233 : int error_bsize_h = AOMMIN(p_height, WARP_ERROR_BLOCK);
828 : uint8_t tmp[WARP_ERROR_BLOCK * WARP_ERROR_BLOCK];
829 91233 : ConvolveParams conv_params = get_conv_params(0, 0, 0, 8);
830 91233 : conv_params.use_jnt_comp_avg = 0;
831 :
832 997091 : for (int i = p_row; i < p_row + p_height; i += WARP_ERROR_BLOCK) {
833 19931400 : for (int j = p_col; j < p_col + p_width; j += WARP_ERROR_BLOCK) {
834 : // avoid warping extra 8x8 blocks in the padded region of the frame
835 : // when p_width and p_height are not multiples of WARP_ERROR_BLOCK
836 19025600 : warp_w = AOMMIN(error_bsize_w, p_col + p_width - j);
837 19025600 : warp_h = AOMMIN(error_bsize_h, p_row + p_height - i);
838 19025600 : warp_plane(wm, ref, width, height, stride, tmp, j, i, warp_w, warp_h,
839 : WARP_ERROR_BLOCK, subsampling_x, subsampling_y, &conv_params);
840 :
841 19024300 : gm_sumerr += eb_av1_calc_frame_error(tmp, WARP_ERROR_BLOCK, dst + j + i * p_stride,
842 : warp_w, warp_h, p_stride);
843 19025500 : if (gm_sumerr > best_error) return gm_sumerr;
844 : }
845 : }
846 2214 : return gm_sumerr;
847 : }
848 :
849 113 : int64_t eb_av1_frame_error(int use_hbd, int bd, const uint8_t *ref, int stride,
850 : uint8_t *dst, int p_width, int p_height, int p_stride) {
851 113 : if (use_hbd) {
852 0 : return highbd_frame_error(CONVERT_TO_SHORTPTR(ref), stride,
853 0 : CONVERT_TO_SHORTPTR(dst), p_width, p_height,
854 : p_stride, bd);
855 : }
856 113 : return eb_av1_calc_frame_error(ref, stride, dst, p_width, p_height, p_stride);
857 : }
858 :
859 91233 : int64_t eb_av1_warp_error(EbWarpedMotionParams *wm, int use_hbd, int bd,
860 : const uint8_t *ref, int width, int height, int stride,
861 : uint8_t *dst, int p_col, int p_row, int p_width,
862 : int p_height, int p_stride, int subsampling_x,
863 : int subsampling_y, int64_t best_error) {
864 91233 : if (wm->wmtype <= AFFINE)
865 91233 : if (!eb_get_shear_params(wm)) return 1;
866 91233 : if (use_hbd)
867 0 : return highbd_warp_error(wm, ref, width, height, stride, dst, p_col, p_row,
868 : p_width, p_height, p_stride, subsampling_x,
869 : subsampling_y, bd, best_error);
870 91233 : return warp_error(wm, ref, width, height, stride, dst, p_col, p_row, p_width,
871 : p_height, p_stride, subsampling_x, subsampling_y,
872 : best_error);
873 : }
874 :
875 7584170 : void eb_av1_warp_plane(EbWarpedMotionParams *wm, int use_hbd, int bd,
876 : const uint8_t *ref, int width, int height, int stride,
877 : uint8_t *pred, int p_col, int p_row, int p_width,
878 : int p_height, int p_stride, int subsampling_x,
879 : int subsampling_y, ConvolveParams *conv_params) {
880 7584170 : if (use_hbd)
881 0 : highbd_warp_plane(wm, ref, width, height, stride, pred, p_col, p_row,
882 : p_width, p_height, p_stride, subsampling_x, subsampling_y,
883 : bd, conv_params);
884 : else
885 7584170 : warp_plane(wm, ref, width, height, stride, pred, p_col, p_row, p_width,
886 : p_height, p_stride, subsampling_x, subsampling_y, conv_params);
887 7585150 : }
888 :
889 : #define LS_MV_MAX 256 // max mv in 1/8-pel
890 : // Use LS_STEP = 8 so that 2 less bits needed for A, Bx, By.
891 : #define LS_STEP 8
892 :
893 : // Assuming LS_MV_MAX is < MAX_SB_SIZE * 8,
894 : // the precision needed is:
895 : // (MAX_SB_SIZE_LOG2 + 3) [for sx * sx magnitude] +
896 : // (MAX_SB_SIZE_LOG2 + 4) [for sx * dx magnitude] +
897 : // 1 [for sign] +
898 : // LEAST_SQUARES_SAMPLES_MAX_BITS
899 : // [for adding up to LEAST_SQUARES_SAMPLES_MAX samples]
900 : // The value is 23
901 : #define LS_MAT_RANGE_BITS \
902 : ((MAX_SB_SIZE_LOG2 + 4) * 2 + LEAST_SQUARES_SAMPLES_MAX_BITS)
903 :
904 : // Bit-depth reduction from the full-range
905 : #define LS_MAT_DOWN_BITS 2
906 :
907 : // bits range of A, Bx and By after downshifting
908 : #define LS_MAT_BITS (LS_MAT_RANGE_BITS - LS_MAT_DOWN_BITS)
909 : #define LS_MAT_MIN (-(1 << (LS_MAT_BITS - 1)))
910 : #define LS_MAT_MAX ((1 << (LS_MAT_BITS - 1)) - 1)
911 :
912 : // By setting LS_STEP = 8, the least 2 bits of every elements in A, Bx, By are
913 : // 0. So, we can reduce LS_MAT_RANGE_BITS(2) bits here.
914 : #define LS_SQUARE(a) \
915 : (((a) * (a)*4 + (a)*4 * LS_STEP + LS_STEP * LS_STEP * 2) >> \
916 : (2 + LS_MAT_DOWN_BITS))
917 : #define LS_PRODUCT1(a, b) \
918 : (((a) * (b)*4 + ((a) + (b)) * 2 * LS_STEP + LS_STEP * LS_STEP) >> \
919 : (2 + LS_MAT_DOWN_BITS))
920 : #define LS_PRODUCT2(a, b) \
921 : (((a) * (b)*4 + ((a) + (b)) * 2 * LS_STEP + LS_STEP * LS_STEP * 2) >> \
922 : (2 + LS_MAT_DOWN_BITS))
923 :
924 : #define USE_LIMITED_PREC_MULT 0
925 :
926 : #if USE_LIMITED_PREC_MULT
927 :
928 : #define MUL_PREC_BITS 16
929 : static uint16_t resolve_multiplier_64(uint64_t D, int16_t *shift) {
930 : int msb = 0;
931 : uint16_t mult = 0;
932 : *shift = 0;
933 : if (D != 0) {
934 : msb = (int16_t)((D >> 32) ? get_msb((unsigned int)(D >> 32)) + 32
935 : : get_msb((unsigned int)D));
936 : if (msb >= MUL_PREC_BITS) {
937 : mult = (uint16_t)ROUND_POWER_OF_TWO_64(D, msb + 1 - MUL_PREC_BITS);
938 : *shift = msb + 1 - MUL_PREC_BITS;
939 : } else {
940 : mult = (uint16_t)D;
941 : *shift = 0;
942 : }
943 : }
944 : return mult;
945 : }
946 :
947 : static int32_t get_mult_shift_ndiag(int64_t Px, int16_t iDet, int shift) {
948 : int32_t ret;
949 : int16_t mshift;
950 : uint16_t Mul = resolve_multiplier_64(llabs(Px), &mshift);
951 : int32_t v = (int32_t)Mul * (int32_t)iDet * (Px < 0 ? -1 : 1);
952 : shift -= mshift;
953 : if (shift > 0) {
954 : return (int32_t)clamp(ROUND_POWER_OF_TWO_SIGNED(v, shift),
955 : -WARPEDMODEL_NONDIAGAFFINE_CLAMP + 1,
956 : WARPEDMODEL_NONDIAGAFFINE_CLAMP - 1);
957 : } else {
958 : return (int32_t)clamp(v * (1 << (-shift)),
959 : -WARPEDMODEL_NONDIAGAFFINE_CLAMP + 1,
960 : WARPEDMODEL_NONDIAGAFFINE_CLAMP - 1);
961 : }
962 : return ret;
963 : }
964 :
965 : static int32_t get_mult_shift_diag(int64_t Px, int16_t iDet, int shift) {
966 : int16_t mshift;
967 : uint16_t Mul = resolve_multiplier_64(llabs(Px), &mshift);
968 : int32_t v = (int32_t)Mul * (int32_t)iDet * (Px < 0 ? -1 : 1);
969 : shift -= mshift;
970 : if (shift > 0) {
971 : return (int32_t)clamp(
972 : ROUND_POWER_OF_TWO_SIGNED(v, shift),
973 : (1 << WARPEDMODEL_PREC_BITS) - WARPEDMODEL_NONDIAGAFFINE_CLAMP + 1,
974 : (1 << WARPEDMODEL_PREC_BITS) + WARPEDMODEL_NONDIAGAFFINE_CLAMP - 1);
975 : } else {
976 : return (int32_t)clamp(
977 : v * (1 << (-shift)),
978 : (1 << WARPEDMODEL_PREC_BITS) - WARPEDMODEL_NONDIAGAFFINE_CLAMP + 1,
979 : (1 << WARPEDMODEL_PREC_BITS) + WARPEDMODEL_NONDIAGAFFINE_CLAMP - 1);
980 : }
981 : }
982 :
983 : #else
984 :
985 1475510 : static int32_t get_mult_shift_ndiag(int64_t Px, int16_t iDet, int shift) {
986 1475510 : int64_t v = Px * (int64_t)iDet;
987 1475510 : return (int32_t)clamp64(ROUND_POWER_OF_TWO_SIGNED_64(v, shift),
988 : -WARPEDMODEL_NONDIAGAFFINE_CLAMP + 1,
989 : WARPEDMODEL_NONDIAGAFFINE_CLAMP - 1);
990 : }
991 :
992 1475490 : static int32_t get_mult_shift_diag(int64_t Px, int16_t iDet, int shift) {
993 1475490 : int64_t v = Px * (int64_t)iDet;
994 2950980 : return (int32_t)clamp64(
995 1475490 : ROUND_POWER_OF_TWO_SIGNED_64(v, shift),
996 : (1 << WARPEDMODEL_PREC_BITS) - WARPEDMODEL_NONDIAGAFFINE_CLAMP + 1,
997 : (1 << WARPEDMODEL_PREC_BITS) + WARPEDMODEL_NONDIAGAFFINE_CLAMP - 1);
998 : }
999 : #endif // USE_LIMITED_PREC_MULT
1000 :
1001 773046 : static int find_affine_int(int np, const int *pts1, const int *pts2,
1002 : BlockSize bsize, int mvy, int mvx,
1003 : EbWarpedMotionParams *wm, int mi_row, int mi_col) {
1004 773046 : int32_t A[2][2] = { { 0, 0 }, { 0, 0 } };
1005 773046 : int32_t Bx[2] = { 0, 0 };
1006 773046 : int32_t By[2] = { 0, 0 };
1007 : int i;
1008 :
1009 773046 : const int bw = block_size_wide[bsize];
1010 773046 : const int bh = block_size_high[bsize];
1011 773046 : const int rsuy = (AOMMAX(bh, MI_SIZE) / 2 - 1);
1012 773046 : const int rsux = (AOMMAX(bw, MI_SIZE) / 2 - 1);
1013 773046 : const int suy = rsuy * 8;
1014 773046 : const int sux = rsux * 8;
1015 773046 : const int duy = suy + mvy;
1016 773046 : const int dux = sux + mvx;
1017 773046 : const int isuy = (mi_row * MI_SIZE + rsuy);
1018 773046 : const int isux = (mi_col * MI_SIZE + rsux);
1019 :
1020 : // Assume the center pixel of the block has exactly the same motion vector
1021 : // as transmitted for the block. First shift the origin of the source
1022 : // points to the block center, and the origin of the destination points to
1023 : // the block center added to the motion vector transmitted.
1024 : // Let (xi, yi) denote the source points and (xi', yi') denote destination
1025 : // points after origin shfifting, for i = 0, 1, 2, .... n-1.
1026 : // Then if P = [x0, y0,
1027 : // x1, y1
1028 : // x2, y1,
1029 : // ....
1030 : // ]
1031 : // q = [x0', x1', x2', ... ]'
1032 : // r = [y0', y1', y2', ... ]'
1033 : // the least squares problems that need to be solved are:
1034 : // [h1, h2]' = inv(P'P)P'q and
1035 : // [h3, h4]' = inv(P'P)P'r
1036 : // where the affine transformation is given by:
1037 : // x' = h1.x + h2.y
1038 : // y' = h3.x + h4.y
1039 : //
1040 : // The loop below computes: A = P'P, Bx = P'q, By = P'r
1041 : // We need to just compute inv(A).Bx and inv(A).By for the solutions.
1042 : // Contribution from neighbor block
1043 2269110 : for (i = 0; i < np; i++) {
1044 1496060 : const int dx = pts2[i * 2] - dux;
1045 1496060 : const int dy = pts2[i * 2 + 1] - duy;
1046 1496060 : const int sx = pts1[i * 2] - sux;
1047 1496060 : const int sy = pts1[i * 2 + 1] - suy;
1048 : // (TODO)yunqing: This comparison wouldn't be necessary if the sample
1049 : // selection is done in find_samples(). Also, global offset can be removed
1050 : // while collecting samples.
1051 1496060 : if (abs(sx - dx) < LS_MV_MAX && abs(sy - dy) < LS_MV_MAX) {
1052 1460760 : A[0][0] += LS_SQUARE(sx);
1053 1460760 : A[0][1] += LS_PRODUCT1(sx, sy);
1054 1460760 : A[1][1] += LS_SQUARE(sy);
1055 1460760 : Bx[0] += LS_PRODUCT2(sx, dx);
1056 1460760 : Bx[1] += LS_PRODUCT1(sy, dx);
1057 1460760 : By[0] += LS_PRODUCT1(sx, dy);
1058 1460760 : By[1] += LS_PRODUCT2(sy, dy);
1059 : }
1060 : }
1061 :
1062 : // Just for debugging, and can be removed later.
1063 773046 : assert(A[0][0] >= LS_MAT_MIN && A[0][0] <= LS_MAT_MAX);
1064 773065 : assert(A[0][1] >= LS_MAT_MIN && A[0][1] <= LS_MAT_MAX);
1065 773071 : assert(A[1][1] >= LS_MAT_MIN && A[1][1] <= LS_MAT_MAX);
1066 773072 : assert(Bx[0] >= LS_MAT_MIN && Bx[0] <= LS_MAT_MAX);
1067 773073 : assert(Bx[1] >= LS_MAT_MIN && Bx[1] <= LS_MAT_MAX);
1068 773076 : assert(By[0] >= LS_MAT_MIN && By[0] <= LS_MAT_MAX);
1069 773077 : assert(By[1] >= LS_MAT_MIN && By[1] <= LS_MAT_MAX);
1070 :
1071 : int64_t Det;
1072 : int16_t iDet, shift;
1073 :
1074 : // Compute Determinant of A
1075 773077 : Det = (int64_t)A[0][0] * A[1][1] - (int64_t)A[0][1] * A[0][1];
1076 773077 : if (Det == 0) return 1;
1077 737786 : iDet = resolve_divisor_64(llabs(Det), &shift) * (Det < 0 ? -1 : 1);
1078 737780 : shift -= WARPEDMODEL_PREC_BITS;
1079 737780 : if (shift < 0) {
1080 0 : iDet <<= (-shift);
1081 0 : shift = 0;
1082 : }
1083 :
1084 : int64_t Px[2], Py[2];
1085 :
1086 : // These divided by the Det, are the least squares solutions
1087 737780 : Px[0] = (int64_t)A[1][1] * Bx[0] - (int64_t)A[0][1] * Bx[1];
1088 737780 : Px[1] = -(int64_t)A[0][1] * Bx[0] + (int64_t)A[0][0] * Bx[1];
1089 737780 : Py[0] = (int64_t)A[1][1] * By[0] - (int64_t)A[0][1] * By[1];
1090 737780 : Py[1] = -(int64_t)A[0][1] * By[0] + (int64_t)A[0][0] * By[1];
1091 :
1092 737780 : wm->wmmat[2] = get_mult_shift_diag(Px[0], iDet, shift);
1093 737773 : wm->wmmat[3] = get_mult_shift_ndiag(Px[1], iDet, shift);
1094 737763 : wm->wmmat[4] = get_mult_shift_ndiag(Py[0], iDet, shift);
1095 737746 : wm->wmmat[5] = get_mult_shift_diag(Py[1], iDet, shift);
1096 :
1097 : // Note: In the vx, vy expressions below, the max value of each of the
1098 : // 2nd and 3rd terms are (2^16 - 1) * (2^13 - 1). That leaves enough room
1099 : // for the first term so that the overall sum in the worst case fits
1100 : // within 32 bits overall.
1101 737746 : int32_t vx = mvx * (1 << (WARPEDMODEL_PREC_BITS - 3)) -
1102 737746 : (isux * (wm->wmmat[2] - (1 << WARPEDMODEL_PREC_BITS)) +
1103 737746 : isuy * wm->wmmat[3]);
1104 737746 : int32_t vy = mvy * (1 << (WARPEDMODEL_PREC_BITS - 3)) -
1105 737746 : (isux * wm->wmmat[4] +
1106 737746 : isuy * (wm->wmmat[5] - (1 << WARPEDMODEL_PREC_BITS)));
1107 737747 : wm->wmmat[0] =
1108 737746 : clamp(vx, -WARPEDMODEL_TRANS_CLAMP, WARPEDMODEL_TRANS_CLAMP - 1);
1109 737737 : wm->wmmat[1] =
1110 737747 : clamp(vy, -WARPEDMODEL_TRANS_CLAMP, WARPEDMODEL_TRANS_CLAMP - 1);
1111 :
1112 737737 : wm->wmmat[6] = wm->wmmat[7] = 0;
1113 737737 : return 0;
1114 : }
1115 :
1116 773044 : EbBool eb_find_projection(
1117 : int np,
1118 : int *pts1,
1119 : int *pts2,
1120 : BlockSize bsize,
1121 : int mvy,
1122 : int mvx,
1123 : EbWarpedMotionParams *wm_params,
1124 : int mi_row,
1125 : int mi_col)
1126 : {
1127 773044 : if (find_affine_int(np, pts1, pts2, bsize, mvy, mvx, wm_params,
1128 : mi_row, mi_col))
1129 : {
1130 35291 : return 1;
1131 : }
1132 :
1133 : // check compatibility with the fast warp filter
1134 737737 : if (!eb_get_shear_params(wm_params))
1135 50471 : return 1;
1136 :
1137 687261 : return 0;
1138 : }
|
