warped_motion.c 62.9 KB
Newer Older
1
/*
2
 * Copyright (c) 2016, Alliance for Open Media. All rights reserved
3
 *
4 5 6 7 8 9
 * This source code is subject to the terms of the BSD 2 Clause License and
 * the Alliance for Open Media Patent License 1.0. If the BSD 2 Clause License
 * was not distributed with this source code in the LICENSE file, you can
 * obtain it at www.aomedia.org/license/software. If the Alliance for Open
 * Media Patent License 1.0 was not distributed with this source code in the
 * PATENTS file, you can obtain it at www.aomedia.org/license/patent.
10 11 12 13 14 15 16 17
 */

#include <stdio.h>
#include <stdlib.h>
#include <memory.h>
#include <math.h>
#include <assert.h>

18
#include "av1/common/warped_motion.h"
19

20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88
/* clang-format off */
static const int error_measure_lut[512] = {
  // power 0.6
  255, 254, 254, 253, 253, 252, 251, 251,
  250, 250, 249, 248, 248, 247, 247, 246,
  245, 245, 244, 243, 243, 242, 242, 241,
  240, 240, 239, 238, 238, 237, 237, 236,
  235, 235, 234, 233, 233, 232, 231, 231,
  230, 230, 229, 228, 228, 227, 226, 226,
  225, 224, 224, 223, 222, 222, 221, 220,
  220, 219, 218, 218, 217, 216, 216, 215,
  214, 214, 213, 212, 212, 211, 210, 210,
  209, 208, 208, 207, 206, 206, 205, 204,
  203, 203, 202, 201, 201, 200, 199, 199,
  198, 197, 196, 196, 195, 194, 194, 193,
  192, 191, 191, 190, 189, 188, 188, 187,
  186, 185, 185, 184, 183, 182, 182, 181,
  180, 179, 179, 178, 177, 176, 176, 175,
  174, 173, 173, 172, 171, 170, 169, 169,
  168, 167, 166, 165, 165, 164, 163, 162,
  161, 161, 160, 159, 158, 157, 156, 156,
  155, 154, 153, 152, 151, 151, 150, 149,
  148, 147, 146, 145, 145, 144, 143, 142,
  141, 140, 139, 138, 137, 137, 136, 135,
  134, 133, 132, 131, 130, 129, 128, 127,
  126, 125, 124, 123, 122, 121, 120, 119,
  118, 117, 116, 115, 114, 113, 112, 111,
  110, 109, 108, 107, 106, 105, 104, 103,
  102, 100,  99,  98,  97,  96,  95,  94,
  92, 91, 90, 89, 88, 86, 85, 84,
  83, 81, 80, 79, 77, 76, 75, 73,
  72, 71, 69, 68, 66, 65, 63, 62,
  60, 59, 57, 55, 54, 52, 50, 48,
  47, 45, 43, 41, 39, 37, 34, 32,
  29, 27, 24, 21, 18, 14,  9,  0,
  9, 14, 18, 21, 24, 27, 29, 32,
  34, 37, 39, 41, 43, 45, 47, 48,
  50, 52, 54, 55, 57, 59, 60, 62,
  63, 65, 66, 68, 69, 71, 72, 73,
  75, 76, 77, 79, 80, 81, 83, 84,
  85, 86, 88, 89, 90, 91, 92, 94,
  95,  96,  97,  98,  99, 100, 102, 103,
  104, 105, 106, 107, 108, 109, 110, 111,
  112, 113, 114, 115, 116, 117, 118, 119,
  120, 121, 122, 123, 124, 125, 126, 127,
  128, 129, 130, 131, 132, 133, 134, 135,
  136, 137, 137, 138, 139, 140, 141, 142,
  143, 144, 145, 145, 146, 147, 148, 149,
  150, 151, 151, 152, 153, 154, 155, 156,
  156, 157, 158, 159, 160, 161, 161, 162,
  163, 164, 165, 165, 166, 167, 168, 169,
  169, 170, 171, 172, 173, 173, 174, 175,
  176, 176, 177, 178, 179, 179, 180, 181,
  182, 182, 183, 184, 185, 185, 186, 187,
  188, 188, 189, 190, 191, 191, 192, 193,
  194, 194, 195, 196, 196, 197, 198, 199,
  199, 200, 201, 201, 202, 203, 203, 204,
  205, 206, 206, 207, 208, 208, 209, 210,
  210, 211, 212, 212, 213, 214, 214, 215,
  216, 216, 217, 218, 218, 219, 220, 220,
  221, 222, 222, 223, 224, 224, 225, 226,
  226, 227, 228, 228, 229, 230, 230, 231,
  231, 232, 233, 233, 234, 235, 235, 236,
  237, 237, 238, 238, 239, 240, 240, 241,
  242, 242, 243, 243, 244, 245, 245, 246,
  247, 247, 248, 248, 249, 250, 250, 251,
  251, 252, 253, 253, 254, 254, 255, 255,
};
/* clang-format on */
89

90
static ProjectPointsFunc get_project_points_type(TransformationType type) {
91
  switch (type) {
92 93 94 95
    case HOMOGRAPHY: return project_points_homography;
    case AFFINE: return project_points_affine;
    case ROTZOOM: return project_points_rotzoom;
    case TRANSLATION: return project_points_translation;
clang-format's avatar
clang-format committed
96
    default: assert(0); return NULL;
97 98 99
  }
}

100
void project_points_translation(int32_t *mat, int *points, int *proj,
101 102 103
                                const int n, const int stride_points,
                                const int stride_proj, const int subsampling_x,
                                const int subsampling_y) {
104 105 106 107
  int i;
  for (i = 0; i < n; ++i) {
    const int x = *(points++), y = *(points++);
    if (subsampling_x)
108
      *(proj++) = ROUND_POWER_OF_TWO_SIGNED(
109
          ((x * (1 << (WARPEDMODEL_PREC_BITS + 1))) + mat[0]),
110
          WARPEDDIFF_PREC_BITS + 1);
111
    else
112
      *(proj++) = ROUND_POWER_OF_TWO_SIGNED(
113
          ((x * (1 << WARPEDMODEL_PREC_BITS)) + mat[0]), WARPEDDIFF_PREC_BITS);
114
    if (subsampling_y)
115
      *(proj++) = ROUND_POWER_OF_TWO_SIGNED(
116
          ((y * (1 << (WARPEDMODEL_PREC_BITS + 1))) + mat[1]),
117
          WARPEDDIFF_PREC_BITS + 1);
118
    else
119
      *(proj++) = ROUND_POWER_OF_TWO_SIGNED(
120
          ((y * (1 << WARPEDMODEL_PREC_BITS))) + mat[1], WARPEDDIFF_PREC_BITS);
121 122 123 124 125
    points += stride_points - 2;
    proj += stride_proj - 2;
  }
}

126
void project_points_rotzoom(int32_t *mat, int *points, int *proj, const int n,
127 128
                            const int stride_points, const int stride_proj,
                            const int subsampling_x, const int subsampling_y) {
129 130 131 132
  int i;
  for (i = 0; i < n; ++i) {
    const int x = *(points++), y = *(points++);
    if (subsampling_x)
133
      *(proj++) = ROUND_POWER_OF_TWO_SIGNED(
134 135
          mat[2] * 2 * x + mat[3] * 2 * y + mat[0] +
              (mat[2] + mat[3] - (1 << WARPEDMODEL_PREC_BITS)) / 2,
136 137
          WARPEDDIFF_PREC_BITS + 1);
    else
138
      *(proj++) = ROUND_POWER_OF_TWO_SIGNED(mat[2] * x + mat[3] * y + mat[0],
139
                                            WARPEDDIFF_PREC_BITS);
140
    if (subsampling_y)
141
      *(proj++) = ROUND_POWER_OF_TWO_SIGNED(
142 143
          -mat[3] * 2 * x + mat[2] * 2 * y + mat[1] +
              (-mat[3] + mat[2] - (1 << WARPEDMODEL_PREC_BITS)) / 2,
144 145
          WARPEDDIFF_PREC_BITS + 1);
    else
146
      *(proj++) = ROUND_POWER_OF_TWO_SIGNED(-mat[3] * x + mat[2] * y + mat[1],
147
                                            WARPEDDIFF_PREC_BITS);
148 149 150 151 152
    points += stride_points - 2;
    proj += stride_proj - 2;
  }
}

153
void project_points_affine(int32_t *mat, int *points, int *proj, const int n,
154 155
                           const int stride_points, const int stride_proj,
                           const int subsampling_x, const int subsampling_y) {
156 157 158 159
  int i;
  for (i = 0; i < n; ++i) {
    const int x = *(points++), y = *(points++);
    if (subsampling_x)
160
      *(proj++) = ROUND_POWER_OF_TWO_SIGNED(
161 162
          mat[2] * 2 * x + mat[3] * 2 * y + mat[0] +
              (mat[2] + mat[3] - (1 << WARPEDMODEL_PREC_BITS)) / 2,
163 164
          WARPEDDIFF_PREC_BITS + 1);
    else
165
      *(proj++) = ROUND_POWER_OF_TWO_SIGNED(mat[2] * x + mat[3] * y + mat[0],
166
                                            WARPEDDIFF_PREC_BITS);
167
    if (subsampling_y)
168
      *(proj++) = ROUND_POWER_OF_TWO_SIGNED(
169
          mat[4] * 2 * x + mat[5] * 2 * y + mat[1] +
Sarah Parker's avatar
Sarah Parker committed
170
              (mat[4] + mat[5] - (1 << WARPEDMODEL_PREC_BITS)) / 2,
171 172
          WARPEDDIFF_PREC_BITS + 1);
    else
173
      *(proj++) = ROUND_POWER_OF_TWO_SIGNED(mat[4] * x + mat[5] * y + mat[1],
174
                                            WARPEDDIFF_PREC_BITS);
175 176 177 178 179
    points += stride_points - 2;
    proj += stride_proj - 2;
  }
}

180
void project_points_homography(int32_t *mat, int *points, int *proj,
181 182 183
                               const int n, const int stride_points,
                               const int stride_proj, const int subsampling_x,
                               const int subsampling_y) {
184
  int i;
185 186
  int64_t x, y, Z;
  int64_t xp, yp;
187 188 189 190 191
  for (i = 0; i < n; ++i) {
    x = *(points++), y = *(points++);
    x = (subsampling_x ? 4 * x + 1 : 2 * x);
    y = (subsampling_y ? 4 * y + 1 : 2 * y);

192 193
    Z = (mat[6] * x + mat[7] * y + (1 << (WARPEDMODEL_ROW3HOMO_PREC_BITS + 1)));
    xp = (mat[2] * x + mat[3] * y + 2 * mat[0]) *
194 195
         (1 << (WARPEDPIXEL_PREC_BITS + WARPEDMODEL_ROW3HOMO_PREC_BITS -
                WARPEDMODEL_PREC_BITS));
196
    yp = (mat[4] * x + mat[5] * y + 2 * mat[1]) *
197 198
         (1 << (WARPEDPIXEL_PREC_BITS + WARPEDMODEL_ROW3HOMO_PREC_BITS -
                WARPEDMODEL_PREC_BITS));
199 200 201 202

    xp = xp > 0 ? (xp + Z / 2) / Z : (xp - Z / 2) / Z;
    yp = yp > 0 ? (yp + Z / 2) / Z : (yp - Z / 2) / Z;

clang-format's avatar
clang-format committed
203 204
    if (subsampling_x) xp = (xp - (1 << (WARPEDPIXEL_PREC_BITS - 1))) / 2;
    if (subsampling_y) yp = (yp - (1 << (WARPEDPIXEL_PREC_BITS - 1))) / 2;
205 206 207 208 209 210 211 212
    *(proj++) = xp;
    *(proj++) = yp;

    points += stride_points - 2;
    proj += stride_proj - 2;
  }
}

Yue Chen's avatar
Yue Chen committed
213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235
// 'points' are at original scale, output 'proj's are scaled up by
// 1 << WARPEDPIXEL_PREC_BITS
void project_points(WarpedMotionParams *wm_params, int *points, int *proj,
                    const int n, const int stride_points, const int stride_proj,
                    const int subsampling_x, const int subsampling_y) {
  switch (wm_params->wmtype) {
    case AFFINE:
      project_points_affine(wm_params->wmmat, points, proj, n, stride_points,
                            stride_proj, subsampling_x, subsampling_y);
      break;
    case ROTZOOM:
      project_points_rotzoom(wm_params->wmmat, points, proj, n, stride_points,
                             stride_proj, subsampling_x, subsampling_y);
      break;
    case HOMOGRAPHY:
      project_points_homography(wm_params->wmmat, points, proj, n,
                                stride_points, stride_proj, subsampling_x,
                                subsampling_y);
      break;
    default: assert(0 && "Invalid warped motion type!"); return;
  }
}

236
static const int16_t
clang-format's avatar
clang-format committed
237 238 239 240 241 242 243 244 245 246 247 248 249 250 251 252 253 254 255 256 257 258 259 260 261 262 263 264 265 266 267 268 269 270
    filter_ntap[WARPEDPIXEL_PREC_SHIFTS][WARPEDPIXEL_FILTER_TAPS] = {
      { 0, 0, 128, 0, 0, 0 },      { 0, -1, 128, 2, -1, 0 },
      { 1, -3, 127, 4, -1, 0 },    { 1, -4, 126, 6, -2, 1 },
      { 1, -5, 126, 8, -3, 1 },    { 1, -6, 125, 11, -4, 1 },
      { 1, -7, 124, 13, -4, 1 },   { 2, -8, 123, 15, -5, 1 },
      { 2, -9, 122, 18, -6, 1 },   { 2, -10, 121, 20, -6, 1 },
      { 2, -11, 120, 22, -7, 2 },  { 2, -12, 119, 25, -8, 2 },
      { 3, -13, 117, 27, -8, 2 },  { 3, -13, 116, 29, -9, 2 },
      { 3, -14, 114, 32, -10, 3 }, { 3, -15, 113, 35, -10, 2 },
      { 3, -15, 111, 37, -11, 3 }, { 3, -16, 109, 40, -11, 3 },
      { 3, -16, 108, 42, -12, 3 }, { 4, -17, 106, 45, -13, 3 },
      { 4, -17, 104, 47, -13, 3 }, { 4, -17, 102, 50, -14, 3 },
      { 4, -17, 100, 52, -14, 3 }, { 4, -18, 98, 55, -15, 4 },
      { 4, -18, 96, 58, -15, 3 },  { 4, -18, 94, 60, -16, 4 },
      { 4, -18, 91, 63, -16, 4 },  { 4, -18, 89, 65, -16, 4 },
      { 4, -18, 87, 68, -17, 4 },  { 4, -18, 85, 70, -17, 4 },
      { 4, -18, 82, 73, -17, 4 },  { 4, -18, 80, 75, -17, 4 },
      { 4, -18, 78, 78, -18, 4 },  { 4, -17, 75, 80, -18, 4 },
      { 4, -17, 73, 82, -18, 4 },  { 4, -17, 70, 85, -18, 4 },
      { 4, -17, 68, 87, -18, 4 },  { 4, -16, 65, 89, -18, 4 },
      { 4, -16, 63, 91, -18, 4 },  { 4, -16, 60, 94, -18, 4 },
      { 3, -15, 58, 96, -18, 4 },  { 4, -15, 55, 98, -18, 4 },
      { 3, -14, 52, 100, -17, 4 }, { 3, -14, 50, 102, -17, 4 },
      { 3, -13, 47, 104, -17, 4 }, { 3, -13, 45, 106, -17, 4 },
      { 3, -12, 42, 108, -16, 3 }, { 3, -11, 40, 109, -16, 3 },
      { 3, -11, 37, 111, -15, 3 }, { 2, -10, 35, 113, -15, 3 },
      { 3, -10, 32, 114, -14, 3 }, { 2, -9, 29, 116, -13, 3 },
      { 2, -8, 27, 117, -13, 3 },  { 2, -8, 25, 119, -12, 2 },
      { 2, -7, 22, 120, -11, 2 },  { 1, -6, 20, 121, -10, 2 },
      { 1, -6, 18, 122, -9, 2 },   { 1, -5, 15, 123, -8, 2 },
      { 1, -4, 13, 124, -7, 1 },   { 1, -4, 11, 125, -6, 1 },
      { 1, -3, 8, 126, -5, 1 },    { 1, -2, 6, 126, -4, 1 },
      { 0, -1, 4, 127, -3, 1 },    { 0, -1, 2, 128, -1, 0 },
    };
271 272 273 274 275 276 277 278 279 280 281

static int32_t do_ntap_filter(int32_t *p, int x) {
  int i;
  int32_t sum = 0;
  for (i = 0; i < WARPEDPIXEL_FILTER_TAPS; ++i) {
    sum += p[i - WARPEDPIXEL_FILTER_TAPS / 2 + 1] * filter_ntap[x][i];
  }
  return sum;
}

static int32_t do_cubic_filter(int32_t *p, int x) {
clang-format's avatar
clang-format committed
282
  if (x == 0) {
283
    return p[0] * (1 << WARPEDPIXEL_FILTER_BITS);
284
  } else if (x == (1 << WARPEDPIXEL_PREC_BITS)) {
285
    return p[1] * (1 << WARPEDPIXEL_FILTER_BITS);
286
  } else {
287
    const int64_t v1 = (int64_t)x * x * x * (3 * (p[0] - p[1]) + p[2] - p[-1]);
288 289 290
    const int64_t v2 = x * x * (2 * p[-1] - 5 * p[0] + 4 * p[1] - p[2]);
    const int64_t v3 = x * (p[1] - p[-1]);
    const int64_t v4 = 2 * p[0];
291
    return (int32_t)ROUND_POWER_OF_TWO_SIGNED(
292 293 294
        (v4 * (1 << (3 * WARPEDPIXEL_PREC_BITS))) +
            (v3 * (1 << (2 * WARPEDPIXEL_PREC_BITS))) +
            (v2 * (1 << WARPEDPIXEL_PREC_BITS)) + v1,
295 296 297 298 299 300 301 302 303 304 305 306 307 308 309
        3 * WARPEDPIXEL_PREC_BITS + 1 - WARPEDPIXEL_FILTER_BITS);
  }
}

static INLINE void get_subcolumn(int taps, uint8_t *ref, int32_t *col,
                                 int stride, int x, int y_start) {
  int i;
  for (i = 0; i < taps; ++i) {
    col[i] = ref[(i + y_start) * stride + x];
  }
}

static uint8_t bi_ntap_filter(uint8_t *ref, int x, int y, int stride) {
  int32_t val, arr[WARPEDPIXEL_FILTER_TAPS];
  int k;
clang-format's avatar
clang-format committed
310 311
  int i = (int)x >> WARPEDPIXEL_PREC_BITS;
  int j = (int)y >> WARPEDPIXEL_PREC_BITS;
312 313 314 315 316 317
  for (k = 0; k < WARPEDPIXEL_FILTER_TAPS; ++k) {
    int32_t arr_temp[WARPEDPIXEL_FILTER_TAPS];
    get_subcolumn(WARPEDPIXEL_FILTER_TAPS, ref, arr_temp, stride,
                  i + k + 1 - WARPEDPIXEL_FILTER_TAPS / 2,
                  j + 1 - WARPEDPIXEL_FILTER_TAPS / 2);
    arr[k] = do_ntap_filter(arr_temp + WARPEDPIXEL_FILTER_TAPS / 2 - 1,
318
                            y - (j * (1 << WARPEDPIXEL_PREC_BITS)));
319 320
  }
  val = do_ntap_filter(arr + WARPEDPIXEL_FILTER_TAPS / 2 - 1,
321
                       x - (i * (1 << WARPEDPIXEL_PREC_BITS)));
322 323 324 325 326 327 328
  val = ROUND_POWER_OF_TWO_SIGNED(val, WARPEDPIXEL_FILTER_BITS * 2);
  return (uint8_t)clip_pixel(val);
}

static uint8_t bi_cubic_filter(uint8_t *ref, int x, int y, int stride) {
  int32_t val, arr[4];
  int k;
clang-format's avatar
clang-format committed
329 330
  int i = (int)x >> WARPEDPIXEL_PREC_BITS;
  int j = (int)y >> WARPEDPIXEL_PREC_BITS;
331 332
  for (k = 0; k < 4; ++k) {
    int32_t arr_temp[4];
clang-format's avatar
clang-format committed
333
    get_subcolumn(4, ref, arr_temp, stride, i + k - 1, j - 1);
334 335
    arr[k] =
        do_cubic_filter(arr_temp + 1, y - (j * (1 << WARPEDPIXEL_PREC_BITS)));
336
  }
337
  val = do_cubic_filter(arr + 1, x - (i * (1 << WARPEDPIXEL_PREC_BITS)));
338 339 340 341 342 343 344
  val = ROUND_POWER_OF_TWO_SIGNED(val, WARPEDPIXEL_FILTER_BITS * 2);
  return (uint8_t)clip_pixel(val);
}

static uint8_t bi_linear_filter(uint8_t *ref, int x, int y, int stride) {
  const int ix = x >> WARPEDPIXEL_PREC_BITS;
  const int iy = y >> WARPEDPIXEL_PREC_BITS;
345 346
  const int sx = x - (ix * (1 << WARPEDPIXEL_PREC_BITS));
  const int sy = y - (iy * (1 << WARPEDPIXEL_PREC_BITS));
347 348 349
  int32_t val;
  val = ROUND_POWER_OF_TWO_SIGNED(
      ref[iy * stride + ix] * (WARPEDPIXEL_PREC_SHIFTS - sy) *
clang-format's avatar
clang-format committed
350 351 352 353
              (WARPEDPIXEL_PREC_SHIFTS - sx) +
          ref[iy * stride + ix + 1] * (WARPEDPIXEL_PREC_SHIFTS - sy) * sx +
          ref[(iy + 1) * stride + ix] * sy * (WARPEDPIXEL_PREC_SHIFTS - sx) +
          ref[(iy + 1) * stride + ix + 1] * sy * sx,
354 355 356 357
      WARPEDPIXEL_PREC_BITS * 2);
  return (uint8_t)clip_pixel(val);
}

clang-format's avatar
clang-format committed
358 359
static uint8_t warp_interpolate(uint8_t *ref, int x, int y, int width,
                                int height, int stride) {
360 361
  int ix = x >> WARPEDPIXEL_PREC_BITS;
  int iy = y >> WARPEDPIXEL_PREC_BITS;
362 363
  int sx = x - (ix * (1 << WARPEDPIXEL_PREC_BITS));
  int sy = y - (iy * (1 << WARPEDPIXEL_PREC_BITS));
364 365
  int32_t v;

clang-format's avatar
clang-format committed
366 367
  if (ix < 0 && iy < 0)
    return ref[0];
Yue Chen's avatar
Yue Chen committed
368
  else if (ix < 0 && iy >= height - 1)
369
    return ref[(height - 1) * stride];
Yue Chen's avatar
Yue Chen committed
370
  else if (ix >= width - 1 && iy < 0)
371
    return ref[width - 1];
Yue Chen's avatar
Yue Chen committed
372
  else if (ix >= width - 1 && iy >= height - 1)
373 374 375 376
    return ref[(height - 1) * stride + (width - 1)];
  else if (ix < 0) {
    v = ROUND_POWER_OF_TWO_SIGNED(
        ref[iy * stride] * (WARPEDPIXEL_PREC_SHIFTS - sy) +
clang-format's avatar
clang-format committed
377
            ref[(iy + 1) * stride] * sy,
378 379 380 381
        WARPEDPIXEL_PREC_BITS);
    return clip_pixel(v);
  } else if (iy < 0) {
    v = ROUND_POWER_OF_TWO_SIGNED(
clang-format's avatar
clang-format committed
382
        ref[ix] * (WARPEDPIXEL_PREC_SHIFTS - sx) + ref[ix + 1] * sx,
383 384
        WARPEDPIXEL_PREC_BITS);
    return clip_pixel(v);
Yue Chen's avatar
Yue Chen committed
385
  } else if (ix >= width - 1) {
386 387
    v = ROUND_POWER_OF_TWO_SIGNED(
        ref[iy * stride + width - 1] * (WARPEDPIXEL_PREC_SHIFTS - sy) +
clang-format's avatar
clang-format committed
388
            ref[(iy + 1) * stride + width - 1] * sy,
389 390
        WARPEDPIXEL_PREC_BITS);
    return clip_pixel(v);
Yue Chen's avatar
Yue Chen committed
391
  } else if (iy >= height - 1) {
392 393
    v = ROUND_POWER_OF_TWO_SIGNED(
        ref[(height - 1) * stride + ix] * (WARPEDPIXEL_PREC_SHIFTS - sx) +
clang-format's avatar
clang-format committed
394
            ref[(height - 1) * stride + ix + 1] * sx,
395 396 397 398 399 400 401
        WARPEDPIXEL_PREC_BITS);
    return clip_pixel(v);
  } else if (ix >= WARPEDPIXEL_FILTER_TAPS / 2 - 1 &&
             iy >= WARPEDPIXEL_FILTER_TAPS / 2 - 1 &&
             ix < width - WARPEDPIXEL_FILTER_TAPS / 2 &&
             iy < height - WARPEDPIXEL_FILTER_TAPS / 2) {
    return bi_ntap_filter(ref, x, y, stride);
clang-format's avatar
clang-format committed
402
  } else if (ix >= 1 && iy >= 1 && ix < width - 2 && iy < height - 2) {
403 404 405 406 407 408
    return bi_cubic_filter(ref, x, y, stride);
  } else {
    return bi_linear_filter(ref, x, y, stride);
  }
}

409 410 411 412 413 414 415 416 417 418 419 420 421 422 423 424 425 426 427 428 429 430 431 432 433 434 435 436 437 438 439 440 441 442 443 444 445 446 447 448 449 450 451 452 453 454 455 456 457 458 459 460 461 462 463 464 465 466 467 468 469 470 471 472 473 474 475 476 477 478 479 480 481 482 483 484 485 486 487 488 489 490 491 492 493 494 495 496 497 498 499 500 501 502 503 504 505 506 507 508 509 510 511 512 513 514 515 516 517 518 519
// For warping, we really use a 6-tap filter, but we do blocks of 8 pixels
// at a time. The zoom/rotation/shear in the model are applied to the
// "fractional" position of each pixel, which therefore varies within
// [-1, 2) * WARPEDPIXEL_PREC_SHIFTS.
// We need an extra 2 taps to fit this in, for a total of 8 taps.
/* clang-format off */
static const int16_t warped_filter[WARPEDPIXEL_PREC_SHIFTS*3][8] = {
  // [-1, 0)
  { 0,   0, 128,   0,   0, 0, 0, 0 }, { 0, - 1, 128,   2, - 1, 0, 0, 0 },
  { 1, - 3, 127,   4, - 1, 0, 0, 0 }, { 1, - 4, 126,   6, - 2, 1, 0, 0 },
  { 1, - 5, 126,   8, - 3, 1, 0, 0 }, { 1, - 6, 125,  11, - 4, 1, 0, 0 },
  { 1, - 7, 124,  13, - 4, 1, 0, 0 }, { 2, - 8, 123,  15, - 5, 1, 0, 0 },
  { 2, - 9, 122,  18, - 6, 1, 0, 0 }, { 2, -10, 121,  20, - 6, 1, 0, 0 },
  { 2, -11, 120,  22, - 7, 2, 0, 0 }, { 2, -12, 119,  25, - 8, 2, 0, 0 },
  { 3, -13, 117,  27, - 8, 2, 0, 0 }, { 3, -13, 116,  29, - 9, 2, 0, 0 },
  { 3, -14, 114,  32, -10, 3, 0, 0 }, { 3, -15, 113,  35, -10, 2, 0, 0 },
  { 3, -15, 111,  37, -11, 3, 0, 0 }, { 3, -16, 109,  40, -11, 3, 0, 0 },
  { 3, -16, 108,  42, -12, 3, 0, 0 }, { 4, -17, 106,  45, -13, 3, 0, 0 },
  { 4, -17, 104,  47, -13, 3, 0, 0 }, { 4, -17, 102,  50, -14, 3, 0, 0 },
  { 4, -17, 100,  52, -14, 3, 0, 0 }, { 4, -18,  98,  55, -15, 4, 0, 0 },
  { 4, -18,  96,  58, -15, 3, 0, 0 }, { 4, -18,  94,  60, -16, 4, 0, 0 },
  { 4, -18,  91,  63, -16, 4, 0, 0 }, { 4, -18,  89,  65, -16, 4, 0, 0 },
  { 4, -18,  87,  68, -17, 4, 0, 0 }, { 4, -18,  85,  70, -17, 4, 0, 0 },
  { 4, -18,  82,  73, -17, 4, 0, 0 }, { 4, -18,  80,  75, -17, 4, 0, 0 },
  { 4, -18,  78,  78, -18, 4, 0, 0 }, { 4, -17,  75,  80, -18, 4, 0, 0 },
  { 4, -17,  73,  82, -18, 4, 0, 0 }, { 4, -17,  70,  85, -18, 4, 0, 0 },
  { 4, -17,  68,  87, -18, 4, 0, 0 }, { 4, -16,  65,  89, -18, 4, 0, 0 },
  { 4, -16,  63,  91, -18, 4, 0, 0 }, { 4, -16,  60,  94, -18, 4, 0, 0 },
  { 3, -15,  58,  96, -18, 4, 0, 0 }, { 4, -15,  55,  98, -18, 4, 0, 0 },
  { 3, -14,  52, 100, -17, 4, 0, 0 }, { 3, -14,  50, 102, -17, 4, 0, 0 },
  { 3, -13,  47, 104, -17, 4, 0, 0 }, { 3, -13,  45, 106, -17, 4, 0, 0 },
  { 3, -12,  42, 108, -16, 3, 0, 0 }, { 3, -11,  40, 109, -16, 3, 0, 0 },
  { 3, -11,  37, 111, -15, 3, 0, 0 }, { 2, -10,  35, 113, -15, 3, 0, 0 },
  { 3, -10,  32, 114, -14, 3, 0, 0 }, { 2, - 9,  29, 116, -13, 3, 0, 0 },
  { 2, - 8,  27, 117, -13, 3, 0, 0 }, { 2, - 8,  25, 119, -12, 2, 0, 0 },
  { 2, - 7,  22, 120, -11, 2, 0, 0 }, { 1, - 6,  20, 121, -10, 2, 0, 0 },
  { 1, - 6,  18, 122, - 9, 2, 0, 0 }, { 1, - 5,  15, 123, - 8, 2, 0, 0 },
  { 1, - 4,  13, 124, - 7, 1, 0, 0 }, { 1, - 4,  11, 125, - 6, 1, 0, 0 },
  { 1, - 3,   8, 126, - 5, 1, 0, 0 }, { 1, - 2,   6, 126, - 4, 1, 0, 0 },
  { 0, - 1,   4, 127, - 3, 1, 0, 0 }, { 0, - 1,   2, 128, - 1, 0, 0, 0 },

  // [0, 1)
  { 0, 0,   0, 128,   0,   0, 0, 0 }, { 0, 0, - 1, 128,   2, - 1, 0, 0 },
  { 0, 1, - 3, 127,   4, - 1, 0, 0 }, { 0, 1, - 4, 126,   6, - 2, 1, 0 },
  { 0, 1, - 5, 126,   8, - 3, 1, 0 }, { 0, 1, - 6, 125,  11, - 4, 1, 0 },
  { 0, 1, - 7, 124,  13, - 4, 1, 0 }, { 0, 2, - 8, 123,  15, - 5, 1, 0 },
  { 0, 2, - 9, 122,  18, - 6, 1, 0 }, { 0, 2, -10, 121,  20, - 6, 1, 0 },
  { 0, 2, -11, 120,  22, - 7, 2, 0 }, { 0, 2, -12, 119,  25, - 8, 2, 0 },
  { 0, 3, -13, 117,  27, - 8, 2, 0 }, { 0, 3, -13, 116,  29, - 9, 2, 0 },
  { 0, 3, -14, 114,  32, -10, 3, 0 }, { 0, 3, -15, 113,  35, -10, 2, 0 },
  { 0, 3, -15, 111,  37, -11, 3, 0 }, { 0, 3, -16, 109,  40, -11, 3, 0 },
  { 0, 3, -16, 108,  42, -12, 3, 0 }, { 0, 4, -17, 106,  45, -13, 3, 0 },
  { 0, 4, -17, 104,  47, -13, 3, 0 }, { 0, 4, -17, 102,  50, -14, 3, 0 },
  { 0, 4, -17, 100,  52, -14, 3, 0 }, { 0, 4, -18,  98,  55, -15, 4, 0 },
  { 0, 4, -18,  96,  58, -15, 3, 0 }, { 0, 4, -18,  94,  60, -16, 4, 0 },
  { 0, 4, -18,  91,  63, -16, 4, 0 }, { 0, 4, -18,  89,  65, -16, 4, 0 },
  { 0, 4, -18,  87,  68, -17, 4, 0 }, { 0, 4, -18,  85,  70, -17, 4, 0 },
  { 0, 4, -18,  82,  73, -17, 4, 0 }, { 0, 4, -18,  80,  75, -17, 4, 0 },
  { 0, 4, -18,  78,  78, -18, 4, 0 }, { 0, 4, -17,  75,  80, -18, 4, 0 },
  { 0, 4, -17,  73,  82, -18, 4, 0 }, { 0, 4, -17,  70,  85, -18, 4, 0 },
  { 0, 4, -17,  68,  87, -18, 4, 0 }, { 0, 4, -16,  65,  89, -18, 4, 0 },
  { 0, 4, -16,  63,  91, -18, 4, 0 }, { 0, 4, -16,  60,  94, -18, 4, 0 },
  { 0, 3, -15,  58,  96, -18, 4, 0 }, { 0, 4, -15,  55,  98, -18, 4, 0 },
  { 0, 3, -14,  52, 100, -17, 4, 0 }, { 0, 3, -14,  50, 102, -17, 4, 0 },
  { 0, 3, -13,  47, 104, -17, 4, 0 }, { 0, 3, -13,  45, 106, -17, 4, 0 },
  { 0, 3, -12,  42, 108, -16, 3, 0 }, { 0, 3, -11,  40, 109, -16, 3, 0 },
  { 0, 3, -11,  37, 111, -15, 3, 0 }, { 0, 2, -10,  35, 113, -15, 3, 0 },
  { 0, 3, -10,  32, 114, -14, 3, 0 }, { 0, 2, - 9,  29, 116, -13, 3, 0 },
  { 0, 2, - 8,  27, 117, -13, 3, 0 }, { 0, 2, - 8,  25, 119, -12, 2, 0 },
  { 0, 2, - 7,  22, 120, -11, 2, 0 }, { 0, 1, - 6,  20, 121, -10, 2, 0 },
  { 0, 1, - 6,  18, 122, - 9, 2, 0 }, { 0, 1, - 5,  15, 123, - 8, 2, 0 },
  { 0, 1, - 4,  13, 124, - 7, 1, 0 }, { 0, 1, - 4,  11, 125, - 6, 1, 0 },
  { 0, 1, - 3,   8, 126, - 5, 1, 0 }, { 0, 1, - 2,   6, 126, - 4, 1, 0 },
  { 0, 0, - 1,   4, 127, - 3, 1, 0 }, { 0, 0, - 1,   2, 128, - 1, 0, 0 },

  // [1, 2)
  { 0, 0, 0,   0, 128,   0,   0, 0 }, { 0, 0, 0, - 1, 128,   2, - 1, 0 },
  { 0, 0, 1, - 3, 127,   4, - 1, 0 }, { 0, 0, 1, - 4, 126,   6, - 2, 1 },
  { 0, 0, 1, - 5, 126,   8, - 3, 1 }, { 0, 0, 1, - 6, 125,  11, - 4, 1 },
  { 0, 0, 1, - 7, 124,  13, - 4, 1 }, { 0, 0, 2, - 8, 123,  15, - 5, 1 },
  { 0, 0, 2, - 9, 122,  18, - 6, 1 }, { 0, 0, 2, -10, 121,  20, - 6, 1 },
  { 0, 0, 2, -11, 120,  22, - 7, 2 }, { 0, 0, 2, -12, 119,  25, - 8, 2 },
  { 0, 0, 3, -13, 117,  27, - 8, 2 }, { 0, 0, 3, -13, 116,  29, - 9, 2 },
  { 0, 0, 3, -14, 114,  32, -10, 3 }, { 0, 0, 3, -15, 113,  35, -10, 2 },
  { 0, 0, 3, -15, 111,  37, -11, 3 }, { 0, 0, 3, -16, 109,  40, -11, 3 },
  { 0, 0, 3, -16, 108,  42, -12, 3 }, { 0, 0, 4, -17, 106,  45, -13, 3 },
  { 0, 0, 4, -17, 104,  47, -13, 3 }, { 0, 0, 4, -17, 102,  50, -14, 3 },
  { 0, 0, 4, -17, 100,  52, -14, 3 }, { 0, 0, 4, -18,  98,  55, -15, 4 },
  { 0, 0, 4, -18,  96,  58, -15, 3 }, { 0, 0, 4, -18,  94,  60, -16, 4 },
  { 0, 0, 4, -18,  91,  63, -16, 4 }, { 0, 0, 4, -18,  89,  65, -16, 4 },
  { 0, 0, 4, -18,  87,  68, -17, 4 }, { 0, 0, 4, -18,  85,  70, -17, 4 },
  { 0, 0, 4, -18,  82,  73, -17, 4 }, { 0, 0, 4, -18,  80,  75, -17, 4 },
  { 0, 0, 4, -18,  78,  78, -18, 4 }, { 0, 0, 4, -17,  75,  80, -18, 4 },
  { 0, 0, 4, -17,  73,  82, -18, 4 }, { 0, 0, 4, -17,  70,  85, -18, 4 },
  { 0, 0, 4, -17,  68,  87, -18, 4 }, { 0, 0, 4, -16,  65,  89, -18, 4 },
  { 0, 0, 4, -16,  63,  91, -18, 4 }, { 0, 0, 4, -16,  60,  94, -18, 4 },
  { 0, 0, 3, -15,  58,  96, -18, 4 }, { 0, 0, 4, -15,  55,  98, -18, 4 },
  { 0, 0, 3, -14,  52, 100, -17, 4 }, { 0, 0, 3, -14,  50, 102, -17, 4 },
  { 0, 0, 3, -13,  47, 104, -17, 4 }, { 0, 0, 3, -13,  45, 106, -17, 4 },
  { 0, 0, 3, -12,  42, 108, -16, 3 }, { 0, 0, 3, -11,  40, 109, -16, 3 },
  { 0, 0, 3, -11,  37, 111, -15, 3 }, { 0, 0, 2, -10,  35, 113, -15, 3 },
  { 0, 0, 3, -10,  32, 114, -14, 3 }, { 0, 0, 2, - 9,  29, 116, -13, 3 },
  { 0, 0, 2, - 8,  27, 117, -13, 3 }, { 0, 0, 2, - 8,  25, 119, -12, 2 },
  { 0, 0, 2, - 7,  22, 120, -11, 2 }, { 0, 0, 1, - 6,  20, 121, -10, 2 },
  { 0, 0, 1, - 6,  18, 122, - 9, 2 }, { 0, 0, 1, - 5,  15, 123, - 8, 2 },
  { 0, 0, 1, - 4,  13, 124, - 7, 1 }, { 0, 0, 1, - 4,  11, 125, - 6, 1 },
  { 0, 0, 1, - 3,   8, 126, - 5, 1 }, { 0, 0, 1, - 2,   6, 126, - 4, 1 },
  { 0, 0, 0, - 1,   4, 127, - 3, 1 }, { 0, 0, 0, - 1,   2, 128, - 1, 0 },
};
/* clang-format on */

Yaowu Xu's avatar
Yaowu Xu committed
520
#if CONFIG_AOM_HIGHBITDEPTH
521 522 523 524 525 526 527 528
static INLINE void highbd_get_subcolumn(int taps, uint16_t *ref, int32_t *col,
                                        int stride, int x, int y_start) {
  int i;
  for (i = 0; i < taps; ++i) {
    col[i] = ref[(i + y_start) * stride + x];
  }
}

clang-format's avatar
clang-format committed
529
static uint16_t highbd_bi_ntap_filter(uint16_t *ref, int x, int y, int stride,
530 531 532
                                      int bd) {
  int32_t val, arr[WARPEDPIXEL_FILTER_TAPS];
  int k;
clang-format's avatar
clang-format committed
533 534
  int i = (int)x >> WARPEDPIXEL_PREC_BITS;
  int j = (int)y >> WARPEDPIXEL_PREC_BITS;
535 536 537 538 539 540
  for (k = 0; k < WARPEDPIXEL_FILTER_TAPS; ++k) {
    int32_t arr_temp[WARPEDPIXEL_FILTER_TAPS];
    highbd_get_subcolumn(WARPEDPIXEL_FILTER_TAPS, ref, arr_temp, stride,
                         i + k + 1 - WARPEDPIXEL_FILTER_TAPS / 2,
                         j + 1 - WARPEDPIXEL_FILTER_TAPS / 2);
    arr[k] = do_ntap_filter(arr_temp + WARPEDPIXEL_FILTER_TAPS / 2 - 1,
541
                            y - (j * (1 << WARPEDPIXEL_PREC_BITS)));
542 543
  }
  val = do_ntap_filter(arr + WARPEDPIXEL_FILTER_TAPS / 2 - 1,
544
                       x - (i * (1 << WARPEDPIXEL_PREC_BITS)));
545
  val = ROUND_POWER_OF_TWO_SIGNED(val, WARPEDPIXEL_FILTER_BITS * 2);
546
  return (uint16_t)clip_pixel_highbd(val, bd);
547 548
}

clang-format's avatar
clang-format committed
549
static uint16_t highbd_bi_cubic_filter(uint16_t *ref, int x, int y, int stride,
550 551 552
                                       int bd) {
  int32_t val, arr[4];
  int k;
clang-format's avatar
clang-format committed
553 554
  int i = (int)x >> WARPEDPIXEL_PREC_BITS;
  int j = (int)y >> WARPEDPIXEL_PREC_BITS;
555 556
  for (k = 0; k < 4; ++k) {
    int32_t arr_temp[4];
clang-format's avatar
clang-format committed
557
    highbd_get_subcolumn(4, ref, arr_temp, stride, i + k - 1, j - 1);
558 559
    arr[k] =
        do_cubic_filter(arr_temp + 1, y - (j * (1 << WARPEDPIXEL_PREC_BITS)));
560
  }
561
  val = do_cubic_filter(arr + 1, x - (i * (1 << WARPEDPIXEL_PREC_BITS)));
562
  val = ROUND_POWER_OF_TWO_SIGNED(val, WARPEDPIXEL_FILTER_BITS * 2);
563
  return (uint16_t)clip_pixel_highbd(val, bd);
564 565
}

clang-format's avatar
clang-format committed
566
static uint16_t highbd_bi_linear_filter(uint16_t *ref, int x, int y, int stride,
567 568 569
                                        int bd) {
  const int ix = x >> WARPEDPIXEL_PREC_BITS;
  const int iy = y >> WARPEDPIXEL_PREC_BITS;
570 571
  const int sx = x - (ix * (1 << WARPEDPIXEL_PREC_BITS));
  const int sy = y - (iy * (1 << WARPEDPIXEL_PREC_BITS));
572 573 574
  int32_t val;
  val = ROUND_POWER_OF_TWO_SIGNED(
      ref[iy * stride + ix] * (WARPEDPIXEL_PREC_SHIFTS - sy) *
clang-format's avatar
clang-format committed
575 576 577 578
              (WARPEDPIXEL_PREC_SHIFTS - sx) +
          ref[iy * stride + ix + 1] * (WARPEDPIXEL_PREC_SHIFTS - sy) * sx +
          ref[(iy + 1) * stride + ix] * sy * (WARPEDPIXEL_PREC_SHIFTS - sx) +
          ref[(iy + 1) * stride + ix + 1] * sy * sx,
579
      WARPEDPIXEL_PREC_BITS * 2);
580
  return (uint16_t)clip_pixel_highbd(val, bd);
581 582
}

clang-format's avatar
clang-format committed
583 584
static uint16_t highbd_warp_interpolate(uint16_t *ref, int x, int y, int width,
                                        int height, int stride, int bd) {
585 586
  int ix = x >> WARPEDPIXEL_PREC_BITS;
  int iy = y >> WARPEDPIXEL_PREC_BITS;
587 588
  int sx = x - (ix * (1 << WARPEDPIXEL_PREC_BITS));
  int sy = y - (iy * (1 << WARPEDPIXEL_PREC_BITS));
589 590
  int32_t v;

clang-format's avatar
clang-format committed
591 592
  if (ix < 0 && iy < 0)
    return ref[0];
593 594 595 596 597 598 599 600 601
  else if (ix < 0 && iy > height - 1)
    return ref[(height - 1) * stride];
  else if (ix > width - 1 && iy < 0)
    return ref[width - 1];
  else if (ix > width - 1 && iy > height - 1)
    return ref[(height - 1) * stride + (width - 1)];
  else if (ix < 0) {
    v = ROUND_POWER_OF_TWO_SIGNED(
        ref[iy * stride] * (WARPEDPIXEL_PREC_SHIFTS - sy) +
clang-format's avatar
clang-format committed
602
            ref[(iy + 1) * stride] * sy,
603
        WARPEDPIXEL_PREC_BITS);
604
    return clip_pixel_highbd(v, bd);
605 606
  } else if (iy < 0) {
    v = ROUND_POWER_OF_TWO_SIGNED(
clang-format's avatar
clang-format committed
607
        ref[ix] * (WARPEDPIXEL_PREC_SHIFTS - sx) + ref[ix + 1] * sx,
608
        WARPEDPIXEL_PREC_BITS);
609
    return clip_pixel_highbd(v, bd);
610 611 612
  } else if (ix > width - 1) {
    v = ROUND_POWER_OF_TWO_SIGNED(
        ref[iy * stride + width - 1] * (WARPEDPIXEL_PREC_SHIFTS - sy) +
clang-format's avatar
clang-format committed
613
            ref[(iy + 1) * stride + width - 1] * sy,
614
        WARPEDPIXEL_PREC_BITS);
615
    return clip_pixel_highbd(v, bd);
616 617 618
  } else if (iy > height - 1) {
    v = ROUND_POWER_OF_TWO_SIGNED(
        ref[(height - 1) * stride + ix] * (WARPEDPIXEL_PREC_SHIFTS - sx) +
clang-format's avatar
clang-format committed
619
            ref[(height - 1) * stride + ix + 1] * sx,
620
        WARPEDPIXEL_PREC_BITS);
621
    return clip_pixel_highbd(v, bd);
622 623 624 625 626
  } else if (ix >= WARPEDPIXEL_FILTER_TAPS / 2 - 1 &&
             iy >= WARPEDPIXEL_FILTER_TAPS / 2 - 1 &&
             ix < width - WARPEDPIXEL_FILTER_TAPS / 2 &&
             iy < height - WARPEDPIXEL_FILTER_TAPS / 2) {
    return highbd_bi_ntap_filter(ref, x, y, stride, bd);
clang-format's avatar
clang-format committed
627
  } else if (ix >= 1 && iy >= 1 && ix < width - 2 && iy < height - 2) {
628 629 630 631 632 633
    return highbd_bi_cubic_filter(ref, x, y, stride, bd);
  } else {
    return highbd_bi_linear_filter(ref, x, y, stride, bd);
  }
}

634 635 636
static INLINE int highbd_error_measure(int err, int bd) {
  const int b = bd - 8;
  const int bmask = (1 << b) - 1;
637
  const int v = (1 << b);
638 639 640 641 642 643 644 645
  int e1, e2;
  err = abs(err);
  e1 = err >> b;
  e2 = err & bmask;
  return error_measure_lut[255 + e1] * (v - e2) +
         error_measure_lut[256 + e1] * e2;
}

646 647 648 649 650 651 652
static void highbd_warp_plane_old(WarpedMotionParams *wm, uint8_t *ref8,
                                  int width, int height, int stride,
                                  uint8_t *pred8, int p_col, int p_row,
                                  int p_width, int p_height, int p_stride,
                                  int subsampling_x, int subsampling_y,
                                  int x_scale, int y_scale, int bd,
                                  int ref_frm) {
653
  int i, j;
654
  ProjectPointsFunc projectpoints = get_project_points_type(wm->wmtype);
655
  uint16_t *pred = CONVERT_TO_SHORTPTR(pred8);
656
  uint16_t *ref = CONVERT_TO_SHORTPTR(ref8);
657
  if (projectpoints == NULL) return;
658 659 660 661 662
  for (i = p_row; i < p_row + p_height; ++i) {
    for (j = p_col; j < p_col + p_width; ++j) {
      int in[2], out[2];
      in[0] = j;
      in[1] = i;
663
      projectpoints(wm->wmmat, in, out, 1, 2, 2, subsampling_x, subsampling_y);
664 665
      out[0] = ROUND_POWER_OF_TWO_SIGNED(out[0] * x_scale, 4);
      out[1] = ROUND_POWER_OF_TWO_SIGNED(out[1] * y_scale, 4);
666 667 668 669 670 671 672 673 674
      if (ref_frm)
        pred[(j - p_col) + (i - p_row) * p_stride] = ROUND_POWER_OF_TWO(
            pred[(j - p_col) + (i - p_row) * p_stride] +
                highbd_warp_interpolate(ref, out[0], out[1], width, height,
                                        stride, bd),
            1);
      else
        pred[(j - p_col) + (i - p_row) * p_stride] = highbd_warp_interpolate(
            ref, out[0], out[1], width, height, stride, bd);
675 676 677 678
    }
  }
}

clang-format's avatar
clang-format committed
679 680 681 682 683
static void highbd_warp_plane(WarpedMotionParams *wm, uint8_t *ref8, int width,
                              int height, int stride, uint8_t *pred8, int p_col,
                              int p_row, int p_width, int p_height,
                              int p_stride, int subsampling_x,
                              int subsampling_y, int x_scale, int y_scale,
684
                              int bd, int ref_frm) {
685 686 687 688 689 690 691 692 693 694 695 696 697 698 699 700 701 702 703 704 705 706 707 708 709 710 711 712 713 714 715 716 717 718 719 720 721 722 723 724 725 726 727 728 729 730 731 732 733 734 735 736 737 738 739 740 741 742 743 744 745 746 747 748 749 750 751 752 753 754 755 756 757 758 759 760 761 762 763 764 765 766 767 768 769 770 771 772 773 774 775 776 777 778 779 780 781 782 783 784 785 786 787 788 789 790 791 792 793 794 795 796 797 798 799
  if (wm->wmtype == ROTZOOM) {
    wm->wmmat[5] = wm->wmmat[2];
    wm->wmmat[4] = -wm->wmmat[3];
  }
  if (wm->wmtype == ROTZOOM || wm->wmtype == AFFINE) {
    int32_t tmp[15 * 8];
    int i, j, k, l, m;
    int32_t *mat = wm->wmmat;

    int32_t alpha = mat[2] - (1 << WARPEDMODEL_PREC_BITS);
    int32_t beta = mat[3];
    int32_t gamma = ((int64_t)mat[4] << WARPEDMODEL_PREC_BITS) / mat[2];
    int32_t delta = mat[5] -
                    (((int64_t)mat[3] * mat[4] + (mat[2] / 2)) / mat[2]) -
                    (1 << WARPEDMODEL_PREC_BITS);
    uint16_t *ref = CONVERT_TO_SHORTPTR(ref8);
    uint16_t *pred = CONVERT_TO_SHORTPTR(pred8);

    if ((4 * abs(alpha) + 7 * abs(beta) > (1 << WARPEDMODEL_PREC_BITS)) ||
        (4 * abs(gamma) + 7 * abs(delta) > (1 << WARPEDMODEL_PREC_BITS))) {
      assert(0 && "Warped motion model is incompatible with new warp filter");
      highbd_warp_plane_old(wm, ref8, width, height, stride, pred8, p_col,
                            p_row, p_width, p_height, p_stride, subsampling_x,
                            subsampling_y, x_scale, y_scale, bd, ref_frm);
      return;
    }

    for (i = p_row; i < p_row + p_height; i += 8) {
      for (j = p_col; j < p_col + p_width; j += 8) {
        int32_t x4, y4, ix4, sx4, iy4, sy4;
        if (subsampling_x)
          x4 = ROUND_POWER_OF_TWO_SIGNED(
              mat[2] * 2 * (j + 4) + mat[3] * 2 * (i + 4) + mat[0] +
                  (mat[2] + mat[3] - (1 << WARPEDMODEL_PREC_BITS)) / 2,
              1);
        else
          x4 = mat[2] * (j + 4) + mat[3] * (i + 4) + mat[0];

        if (subsampling_y)
          y4 = ROUND_POWER_OF_TWO_SIGNED(
              mat[4] * 2 * (j + 4) + mat[5] * 2 * (i + 4) + mat[1] +
                  (mat[4] + mat[5] - (1 << WARPEDMODEL_PREC_BITS)) / 2,
              1);
        else
          y4 = mat[4] * (j + 4) + mat[5] * (i + 4) + mat[1];

        ix4 = x4 >> WARPEDMODEL_PREC_BITS;
        sx4 = x4 & ((1 << WARPEDMODEL_PREC_BITS) - 1);
        iy4 = y4 >> WARPEDMODEL_PREC_BITS;
        sy4 = y4 & ((1 << WARPEDMODEL_PREC_BITS) - 1);

        // Horizontal filter
        for (k = -7; k < 8; ++k) {
          int iy = iy4 + k;
          if (iy < 0)
            iy = 0;
          else if (iy > height - 1)
            iy = height - 1;

          for (l = -4; l < 4; ++l) {
            int ix = ix4 + l;
            int sx = ROUND_POWER_OF_TWO_SIGNED(sx4 + alpha * l + beta * k,
                                               WARPEDDIFF_PREC_BITS);
            const int16_t *coeffs = warped_filter[sx + WARPEDPIXEL_PREC_SHIFTS];
            int32_t sum = 0;
            for (m = 0; m < 8; ++m) {
              if (ix + m - 3 < 0)
                sum += ref[iy * stride] * coeffs[m];
              else if (ix + m - 3 > width - 1)
                sum += ref[iy * stride + width - 1] * coeffs[m];
              else
                sum += ref[iy * stride + ix + m - 3] * coeffs[m];
            }
            tmp[(k + 7) * 8 + (l + 4)] = sum;
          }
        }

        // Vertical filter
        for (k = -4; k < AOMMIN(4, p_row + p_height - i - 4); ++k) {
          for (l = -4; l < AOMMIN(4, p_col + p_width - j - 4); ++l) {
            uint16_t *p =
                &pred[(i - p_row + k + 4) * p_stride + (j - p_col + l + 4)];
            int sy = ROUND_POWER_OF_TWO_SIGNED(sy4 + gamma * l + delta * k,
                                               WARPEDDIFF_PREC_BITS);
            const int16_t *coeffs = warped_filter[sy + WARPEDPIXEL_PREC_SHIFTS];
            int32_t sum = 0;
            for (m = 0; m < 8; ++m) {
              sum += tmp[(k + m + 4) * 8 + (l + 4)] * coeffs[m];
            }
            sum = clip_pixel_highbd(
                ROUND_POWER_OF_TWO_SIGNED(sum, 2 * WARPEDPIXEL_FILTER_BITS),
                bd);
            if (ref_frm)
              *p = ROUND_POWER_OF_TWO_SIGNED(*p + sum, 1);
            else
              *p = sum;
          }
        }
      }
    }
  } else {
    highbd_warp_plane_old(wm, ref8, width, height, stride, pred8, p_col, p_row,
                          p_width, p_height, p_stride, subsampling_x,
                          subsampling_y, x_scale, y_scale, bd, ref_frm);
  }
}

static double highbd_warp_erroradv(WarpedMotionParams *wm, uint8_t *ref8,
                                   int width, int height, int stride,
                                   uint8_t *dst8, int p_col, int p_row,
                                   int p_width, int p_height, int p_stride,
                                   int subsampling_x, int subsampling_y,
                                   int x_scale, int y_scale, int bd) {
  int gm_err = 0, no_gm_err = 0;
  int64_t gm_sumerr = 0, no_gm_sumerr = 0;
800
  int i, j;
801 802
  uint16_t *tmp = aom_malloc(p_width * p_height * sizeof(*tmp));
  uint16_t *dst = CONVERT_TO_SHORTPTR(dst8);
803
  uint16_t *ref = CONVERT_TO_SHORTPTR(ref8);
804 805 806 807 808 809 810 811 812 813
  highbd_warp_plane(wm, ref8, width, height, stride, CONVERT_TO_BYTEPTR(tmp),
                    p_col, p_row, p_width, p_height, p_width, subsampling_x,
                    subsampling_y, x_scale, y_scale, bd, 0);
  for (i = 0; i < p_height; ++i) {
    for (j = 0; j < p_width; ++j) {
      gm_err = dst[j + i * p_stride] - tmp[j + i * p_width];
      no_gm_err =
          dst[j + i * p_stride] - ref[(j + p_col) + (i + p_row) * stride];
      gm_sumerr += highbd_error_measure(gm_err, bd);
      no_gm_sumerr += highbd_error_measure(no_gm_err, bd);
814 815
    }
  }
816 817
  aom_free(tmp);
  return (double)gm_sumerr / no_gm_sumerr;
818
}
Yaowu Xu's avatar
Yaowu Xu committed
819
#endif  // CONFIG_AOM_HIGHBITDEPTH
820

821 822 823 824
static INLINE int error_measure(int err) {
  return error_measure_lut[255 + err];
}

825 826 827 828 829
static void warp_plane_old(WarpedMotionParams *wm, uint8_t *ref, int width,
                           int height, int stride, uint8_t *pred, int p_col,
                           int p_row, int p_width, int p_height, int p_stride,
                           int subsampling_x, int subsampling_y, int x_scale,
                           int y_scale, int ref_frm) {
830
  int i, j;
831
  ProjectPointsFunc projectpoints = get_project_points_type(wm->wmtype);
832
  if (projectpoints == NULL) return;
833 834 835 836 837
  for (i = p_row; i < p_row + p_height; ++i) {
    for (j = p_col; j < p_col + p_width; ++j) {
      int in[2], out[2];
      in[0] = j;
      in[1] = i;
838
      projectpoints(wm->wmmat, in, out, 1, 2, 2, subsampling_x, subsampling_y);
839 840
      out[0] = ROUND_POWER_OF_TWO_SIGNED(out[0] * x_scale, 4);
      out[1] = ROUND_POWER_OF_TWO_SIGNED(out[1] * y_scale, 4);
841 842 843 844 845 846 847 848
      if (ref_frm)
        pred[(j - p_col) + (i - p_row) * p_stride] = ROUND_POWER_OF_TWO(
            pred[(j - p_col) + (i - p_row) * p_stride] +
                warp_interpolate(ref, out[0], out[1], width, height, stride),
            1);
      else
        pred[(j - p_col) + (i - p_row) * p_stride] =
            warp_interpolate(ref, out[0], out[1], width, height, stride);
849 850 851 852
    }
  }
}

853 854 855 856 857 858 859 860 861 862 863 864 865 866 867 868 869 870 871 872 873 874 875 876 877 878 879 880 881 882
/* The warp filter for ROTZOOM and AFFINE models works as follows:
   * Split the input into 8x8 blocks
   * For each block, project the point (4, 4) within the block, to get the
     overall block position. Split into integer and fractional coordinates,
     maintaining full WARPEDMODEL precision
   * Filter horizontally: Generate 15 rows of 8 pixels each. Each pixel gets a
     variable horizontal offset. This means that, while the rows of the
     intermediate buffer align with the rows of the *reference* image, the
     columns align with the columns of the *destination* image.
   * Filter vertically: Generate the output block (up to 8x8 pixels, but if the
     destination is too small we crop the output at this stage). Each pixel has
     a variable vertical offset, so that the resulting rows are aligned with
     the rows of the destination image.

   To accomplish these alignments, we factor the warp matrix as a
   product of two shear / asymmetric zoom matrices:
   / a b \  = /   1       0    \ * / 1+alpha  beta \
   \ c d /    \ gamma  1+delta /   \    0      1   /
   where a, b, c, d are wmmat[2], wmmat[3], wmmat[4], wmmat[5] respectively.
   The second shear (with alpha and beta) is applied by the horizontal filter,
   then the first shear (with gamma and delta) is applied by the vertical
   filter.

   The only limitation is that, to fit this in a fixed 8-tap filter size,
   the fractional pixel offsets must be at most +-1. Since the horizontal filter
   generates 15 rows of 8 columns, and the initial point we project is at (4, 4)
   within the block, the parameters must satisfy
   4 * |alpha| + 7 * |beta| <= 1   and   4 * |gamma| + 7 * |delta| <= 1
   for this filter to be applicable.
*/
883 884 885 886
static void warp_plane(WarpedMotionParams *wm, uint8_t *ref, int width,
                       int height, int stride, uint8_t *pred, int p_col,
                       int p_row, int p_width, int p_height, int p_stride,
                       int subsampling_x, int subsampling_y, int x_scale,
887
                       int y_scale, int ref_frm) {
888 889 890 891 892 893 894 895 896 897 898 899 900 901 902 903 904 905 906 907 908 909 910 911 912 913 914 915 916 917 918 919 920 921 922 923 924 925 926 927 928 929 930 931 932 933 934 935 936 937 938 939 940 941 942 943 944 945 946 947 948 949 950 951 952 953 954 955 956 957 958 959 960 961 962 963 964 965 966 967 968 969 970 971 972 973 974 975 976 977 978 979 980 981 982 983 984 985 986 987 988 989 990 991 992 993 994 995 996 997 998
  if (wm->wmtype == ROTZOOM) {
    wm->wmmat[5] = wm->wmmat[2];
    wm->wmmat[4] = -wm->wmmat[3];
  }
  if (wm->wmtype == ROTZOOM || wm->wmtype == AFFINE) {
    int32_t tmp[15 * 8];
    int i, j, k, l, m;
    int32_t *mat = wm->wmmat;

    int32_t alpha = mat[2] - (1 << WARPEDMODEL_PREC_BITS);
    int32_t beta = mat[3];
    int32_t gamma = ((int64_t)mat[4] << WARPEDMODEL_PREC_BITS) / mat[2];
    int32_t delta = mat[5] -
                    (((int64_t)mat[3] * mat[4] + (mat[2] / 2)) / mat[2]) -
                    (1 << WARPEDMODEL_PREC_BITS);

    if ((4 * abs(alpha) + 7 * abs(beta) > (1 << WARPEDMODEL_PREC_BITS)) ||
        (4 * abs(gamma) + 7 * abs(delta) > (1 << WARPEDMODEL_PREC_BITS))) {
      assert(0 && "Warped motion model is incompatible with new warp filter");
      warp_plane_old(wm, ref, width, height, stride, pred, p_col, p_row,
                     p_width, p_height, p_stride, subsampling_x, subsampling_y,
                     x_scale, y_scale, ref_frm);
      return;
    }

    for (i = p_row; i < p_row + p_height; i += 8) {
      for (j = p_col; j < p_col + p_width; j += 8) {
        int32_t x4, y4, ix4, sx4, iy4, sy4;
        if (subsampling_x)
          x4 = ROUND_POWER_OF_TWO_SIGNED(
              mat[2] * 2 * (j + 4) + mat[3] * 2 * (i + 4) + mat[0] +
                  (mat[2] + mat[3] - (1 << WARPEDMODEL_PREC_BITS)) / 2,
              1);
        else
          x4 = mat[2] * (j + 4) + mat[3] * (i + 4) + mat[0];

        if (subsampling_y)
          y4 = ROUND_POWER_OF_TWO_SIGNED(
              mat[4] * 2 * (j + 4) + mat[5] * 2 * (i + 4) + mat[1] +
                  (mat[4] + mat[5] - (1 << WARPEDMODEL_PREC_BITS)) / 2,
              1);
        else
          y4 = mat[4]