test/av1_txfm_test.cc

   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 "test/av1_txfm_test.h"
  14
  15 namespace libaom_test {
  16
  17 int get_txfm1d_size(TX_SIZE tx_size) { return tx_size_wide[tx_size]; }
  18
  19 void get_txfm1d_type(TX_TYPE txfm2d_type, TYPE_TXFM *type0, TYPE_TXFM *type1) {
  20   switch (txfm2d_type) {
  21     case DCT_DCT:
  22       *type0 = TYPE_DCT;
  23       *type1 = TYPE_DCT;
  24       break;
  25     case ADST_DCT:
  26       *type0 = TYPE_ADST;
  27       *type1 = TYPE_DCT;
  28       break;
  29     case DCT_ADST:
  30       *type0 = TYPE_DCT;
  31       *type1 = TYPE_ADST;
  32       break;
  33     case ADST_ADST:
  34       *type0 = TYPE_ADST;
  35       *type1 = TYPE_ADST;
  36       break;
  37     case FLIPADST_DCT:
  38       *type0 = TYPE_ADST;
  39       *type1 = TYPE_DCT;
  40       break;
  41     case DCT_FLIPADST:
  42       *type0 = TYPE_DCT;
  43       *type1 = TYPE_ADST;
  44       break;
  45     case FLIPADST_FLIPADST:
  46       *type0 = TYPE_ADST;
  47       *type1 = TYPE_ADST;
  48       break;
  49     case ADST_FLIPADST:
  50       *type0 = TYPE_ADST;
  51       *type1 = TYPE_ADST;
  52       break;
  53     case FLIPADST_ADST:
  54       *type0 = TYPE_ADST;
  55       *type1 = TYPE_ADST;
  56       break;
  57     case IDTX:
  58       *type0 = TYPE_IDTX;
  59       *type1 = TYPE_IDTX;
  60       break;
  61     case H_DCT:
  62       *type0 = TYPE_IDTX;
  63       *type1 = TYPE_DCT;
  64       break;
  65     case V_DCT:
  66       *type0 = TYPE_DCT;
  67       *type1 = TYPE_IDTX;
  68       break;
  69     case H_ADST:
  70       *type0 = TYPE_IDTX;
  71       *type1 = TYPE_ADST;
  72       break;
  73     case V_ADST:
  74       *type0 = TYPE_ADST;
  75       *type1 = TYPE_IDTX;
  76       break;
  77     case H_FLIPADST:
  78       *type0 = TYPE_IDTX;
  79       *type1 = TYPE_ADST;
  80       break;
  81     case V_FLIPADST:
  82       *type0 = TYPE_ADST;
  83       *type1 = TYPE_IDTX;
  84       break;
  85     default:
  86       *type0 = TYPE_DCT;
  87       *type1 = TYPE_DCT;
  88       assert(0);
  89       break;
  90   }
  91 }
  92
  93 double Sqrt2 = pow(2, 0.5);
  94 double invSqrt2 = 1 / pow(2, 0.5);
  95
  96 double dct_matrix(double n, double k, int size) {
  97   return cos(M_PI * (2 * n + 1) * k / (2 * size));
  98 }
  99
 100 void reference_dct_1d(const double *in, double *out, int size) {
 101   for (int k = 0; k < size; ++k) {
 102     out[k] = 0;
 103     for (int n = 0; n < size; ++n) {
 104       out[k] += in[n] * dct_matrix(n, k, size);
 105     }
 106     if (k == 0) out[k] = out[k] * invSqrt2;
 107   }
 108 }
 109
 110 void reference_idct_1d(const double *in, double *out, int size) {
 111   for (int k = 0; k < size; ++k) {
 112     out[k] = 0;
 113     for (int n = 0; n < size; ++n) {
 114       if (n == 0)
 115         out[k] += invSqrt2 * in[n] * dct_matrix(k, n, size);
 116       else
 117         out[k] += in[n] * dct_matrix(k, n, size);
 118     }
 119   }
 120 }
 121
 122 // TODO(any): Copied from the old 'fadst4' (same as the new 'av1_fadst4_new'
 123 // function). Should be replaced by a proper reference function that takes
 124 // 'double' input & output.
 125 static void fadst4_new(const tran_low_t *input, tran_low_t *output) {
 126   tran_high_t x0, x1, x2, x3;
 127   tran_high_t s0, s1, s2, s3, s4, s5, s6, s7;
 128
 129   x0 = input[0];
 130   x1 = input[1];
 131   x2 = input[2];
 132   x3 = input[3];
 133
 134   if (!(x0 | x1 | x2 | x3)) {
 135     output[0] = output[1] = output[2] = output[3] = 0;
 136     return;
 137   }
 138
 139   s0 = sinpi_1_9 * x0;
 140   s1 = sinpi_4_9 * x0;
 141   s2 = sinpi_2_9 * x1;
 142   s3 = sinpi_1_9 * x1;
 143   s4 = sinpi_3_9 * x2;
 144   s5 = sinpi_4_9 * x3;
 145   s6 = sinpi_2_9 * x3;
 146   s7 = x0 + x1 - x3;
 147
 148   x0 = s0 + s2 + s5;
 149   x1 = sinpi_3_9 * s7;
 150   x2 = s1 - s3 + s6;
 151   x3 = s4;
 152
 153   s0 = x0 + x3;
 154   s1 = x1;
 155   s2 = x2 - x3;
 156   s3 = x2 - x0 + x3;
 157
 158   // 1-D transform scaling factor is sqrt(2).
 159   output[0] = (tran_low_t)fdct_round_shift(s0);
 160   output[1] = (tran_low_t)fdct_round_shift(s1);
 161   output[2] = (tran_low_t)fdct_round_shift(s2);
 162   output[3] = (tran_low_t)fdct_round_shift(s3);
 163 }
 164
 165 void reference_adst_1d(const double *in, double *out, int size) {
 166   if (size == 4) {  // Special case.
 167     tran_low_t int_input[4];
 168     for (int i = 0; i < 4; ++i) {
 169       int_input[i] = static_cast<tran_low_t>(round(in[i]));
 170     }
 171     tran_low_t int_output[4];
 172     fadst4_new(int_input, int_output);
 173     for (int i = 0; i < 4; ++i) {
 174       out[i] = int_output[i];
 175     }
 176     return;
 177   }
 178
 179   for (int k = 0; k < size; ++k) {
 180     out[k] = 0;
 181     for (int n = 0; n < size; ++n) {
 182       out[k] += in[n] * sin(M_PI * (2 * n + 1) * (2 * k + 1) / (4 * size));
 183     }
 184   }
 185 }
 186
 187 void reference_idtx_1d(const double *in, double *out, int size) {
 188   double scale = 0;
 189   if (size == 4)
 190     scale = Sqrt2;
 191   else if (size == 8)
 192     scale = 2;
 193   else if (size == 16)
 194     scale = 2 * Sqrt2;
 195   else if (size == 32)
 196     scale = 4;
 197   else if (size == 64)
 198     scale = 4 * Sqrt2;
 199   for (int k = 0; k < size; ++k) {
 200     out[k] = in[k] * scale;
 201   }
 202 }
 203
 204 void reference_hybrid_1d(double *in, double *out, int size, int type) {
 205   if (type == TYPE_DCT)
 206     reference_dct_1d(in, out, size);
 207   else if (type == TYPE_ADST)
 208     reference_adst_1d(in, out, size);
 209   else
 210     reference_idtx_1d(in, out, size);
 211 }
 212
 213 double get_amplification_factor(TX_TYPE tx_type, TX_SIZE tx_size) {
 214   TXFM_2D_FLIP_CFG fwd_txfm_flip_cfg;
 215   av1_get_fwd_txfm_cfg(tx_type, tx_size, &fwd_txfm_flip_cfg);
 216   const int tx_width = tx_size_wide[fwd_txfm_flip_cfg.tx_size];
 217   const int tx_height = tx_size_high[fwd_txfm_flip_cfg.tx_size];
 218   const int8_t *shift = fwd_txfm_flip_cfg.shift;
 219   const int amplify_bit = shift[0] + shift[1] + shift[2];
 220   double amplify_factor =
 221       amplify_bit >= 0 ? (1 << amplify_bit) : (1.0 / (1 << -amplify_bit));
 222
 223   // For rectangular transforms, we need to multiply by an extra factor.
 224   const int rect_type = get_rect_tx_log_ratio(tx_width, tx_height);
 225   if (abs(rect_type) == 1) {
 226     amplify_factor *= pow(2, 0.5);
 227   }
 228   return amplify_factor;
 229 }
 230
 231 void reference_hybrid_2d(double *in, double *out, TX_TYPE tx_type,
 232                          TX_SIZE tx_size) {
 233   // Get transform type and size of each dimension.
 234   TYPE_TXFM type0;
 235   TYPE_TXFM type1;
 236   get_txfm1d_type(tx_type, &type0, &type1);
 237   const int tx_width = tx_size_wide[tx_size];
 238   const int tx_height = tx_size_high[tx_size];
 239
 240   double *const temp_in = new double[AOMMAX(tx_width, tx_height)];
 241   double *const temp_out = new double[AOMMAX(tx_width, tx_height)];
 242   double *const out_interm = new double[tx_width * tx_height];
 243   const int stride = tx_width;
 244
 245   // Transform columns.
 246   for (int c = 0; c < tx_width; ++c) {
 247     for (int r = 0; r < tx_height; ++r) {
 248       temp_in[r] = in[r * stride + c];
 249     }
 250     reference_hybrid_1d(temp_in, temp_out, tx_height, type0);
 251     for (int r = 0; r < tx_height; ++r) {
 252       out_interm[r * stride + c] = temp_out[r];
 253     }
 254   }
 255
 256   // Transform rows.
 257   for (int r = 0; r < tx_height; ++r) {
 258     reference_hybrid_1d(out_interm + r * stride, out + r * stride, tx_width,
 259                         type1);
 260   }
 261
 262   delete[] temp_in;
 263   delete[] temp_out;
 264   delete[] out_interm;
 265
 266   // These transforms use an approximate 2D DCT transform, by only keeping the
 267   // top-left quarter of the coefficients, and repacking them in the first
 268   // quarter indices.
 269   // TODO(urvang): Refactor this code.
 270   if (tx_width == 64 && tx_height == 64) {  // tx_size == TX_64X64
 271     // Zero out top-right 32x32 area.
 272     for (int row = 0; row < 32; ++row) {
 273       memset(out + row * 64 + 32, 0, 32 * sizeof(*out));
 274     }
 275     // Zero out the bottom 64x32 area.
 276     memset(out + 32 * 64, 0, 32 * 64 * sizeof(*out));
 277     // Re-pack non-zero coeffs in the first 32x32 indices.
 278     for (int row = 1; row < 32; ++row) {
 279       memcpy(out + row * 32, out + row * 64, 32 * sizeof(*out));
 280     }
 281   } else if (tx_width == 32 && tx_height == 64) {  // tx_size == TX_32X64
 282     // Zero out the bottom 32x32 area.
 283     memset(out + 32 * 32, 0, 32 * 32 * sizeof(*out));
 284     // Note: no repacking needed here.
 285   } else if (tx_width == 64 && tx_height == 32) {  // tx_size == TX_64X32
 286     // Zero out right 32x32 area.
 287     for (int row = 0; row < 32; ++row) {
 288       memset(out + row * 64 + 32, 0, 32 * sizeof(*out));
 289     }
 290     // Re-pack non-zero coeffs in the first 32x32 indices.
 291     for (int row = 1; row < 32; ++row) {
 292       memcpy(out + row * 32, out + row * 64, 32 * sizeof(*out));
 293     }
 294   } else if (tx_width == 16 && tx_height == 64) {  // tx_size == TX_16X64
 295     // Zero out the bottom 16x32 area.
 296     memset(out + 16 * 32, 0, 16 * 32 * sizeof(*out));
 297     // Note: no repacking needed here.
 298   } else if (tx_width == 64 && tx_height == 16) {  // tx_size == TX_64X16
 299     // Zero out right 32x16 area.
 300     for (int row = 0; row < 16; ++row) {
 301       memset(out + row * 64 + 32, 0, 32 * sizeof(*out));
 302     }
 303     // Re-pack non-zero coeffs in the first 32x16 indices.
 304     for (int row = 1; row < 16; ++row) {
 305       memcpy(out + row * 32, out + row * 64, 32 * sizeof(*out));
 306     }
 307   }
 308
 309   // Apply appropriate scale.
 310   const double amplify_factor = get_amplification_factor(tx_type, tx_size);
 311   for (int c = 0; c < tx_width; ++c) {
 312     for (int r = 0; r < tx_height; ++r) {
 313       out[r * stride + c] *= amplify_factor;
 314     }
 315   }
 316 }
 317
 318 template <typename Type>
 319 void fliplr(Type *dest, int width, int height, int stride) {
 320   for (int r = 0; r < height; ++r) {
 321     for (int c = 0; c < width / 2; ++c) {
 322       const Type tmp = dest[r * stride + c];
 323       dest[r * stride + c] = dest[r * stride + width - 1 - c];
 324       dest[r * stride + width - 1 - c] = tmp;
 325     }
 326   }
 327 }
 328
 329 template <typename Type>
 330 void flipud(Type *dest, int width, int height, int stride) {
 331   for (int c = 0; c < width; ++c) {
 332     for (int r = 0; r < height / 2; ++r) {
 333       const Type tmp = dest[r * stride + c];
 334       dest[r * stride + c] = dest[(height - 1 - r) * stride + c];
 335       dest[(height - 1 - r) * stride + c] = tmp;
 336     }
 337   }
 338 }
 339
 340 template <typename Type>
 341 void fliplrud(Type *dest, int width, int height, int stride) {
 342   for (int r = 0; r < height / 2; ++r) {
 343     for (int c = 0; c < width; ++c) {
 344       const Type tmp = dest[r * stride + c];
 345       dest[r * stride + c] = dest[(height - 1 - r) * stride + width - 1 - c];
 346       dest[(height - 1 - r) * stride + width - 1 - c] = tmp;
 347     }
 348   }
 349 }
 350
 351 template void fliplr<double>(double *dest, int width, int height, int stride);
 352 template void flipud<double>(double *dest, int width, int height, int stride);
 353 template void fliplrud<double>(double *dest, int width, int height, int stride);
 354
 355 int bd_arr[BD_NUM] = { 8, 10, 12 };
 356
 357 int8_t low_range_arr[BD_NUM] = { 18, 32, 32 };
 358 int8_t high_range_arr[BD_NUM] = { 32, 32, 32 };
 359
 360 void txfm_stage_range_check(const int8_t *stage_range, int stage_num,
 361                             int8_t cos_bit, int low_range, int high_range) {
 362   for (int i = 0; i < stage_num; ++i) {
 363     EXPECT_LE(stage_range[i], low_range);
 364     ASSERT_LE(stage_range[i] + cos_bit, high_range) << "stage = " << i;
 365   }
 366   for (int i = 0; i < stage_num - 1; ++i) {
 367     // make sure there is no overflow while doing half_btf()
 368     ASSERT_LE(stage_range[i + 1] + cos_bit, high_range) << "stage = " << i;
 369   }
 370 }
 371 }  // namespace libaom_test