aom_dsp/psnr.c

   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 <assert.h>
  13 #include <math.h>
  14
  15 #include "config/aom_dsp_rtcd.h"
  16
  17 #include "aom_dsp/psnr.h"
  18 #include "aom_scale/yv12config.h"
  19
  20 double aom_sse_to_psnr(double samples, double peak, double sse) {
  21   if (sse > 0.0) {
  22     const double psnr = 10.0 * log10(samples * peak * peak / sse);
  23     return psnr > MAX_PSNR ? MAX_PSNR : psnr;
  24   } else {
  25     return MAX_PSNR;
  26   }
  27 }
  28
  29 /* TODO(yaowu): The block_variance calls the unoptimized versions of variance()
  30  * and highbd_8_variance(). It should not.
  31  */
  32 static void encoder_variance(const uint8_t *a, int a_stride, const uint8_t *b,
  33                              int b_stride, int w, int h, unsigned int *sse,
  34                              int *sum) {
  35   int i, j;
  36
  37   *sum = 0;
  38   *sse = 0;
  39
  40   for (i = 0; i < h; i++) {
  41     for (j = 0; j < w; j++) {
  42       const int diff = a[j] - b[j];
  43       *sum += diff;
  44       *sse += diff * diff;
  45     }
  46
  47     a += a_stride;
  48     b += b_stride;
  49   }
  50 }
  51
  52 static void encoder_highbd_variance64(const uint8_t *a8, int a_stride,
  53                                       const uint8_t *b8, int b_stride, int w,
  54                                       int h, uint64_t *sse, int64_t *sum) {
  55   const uint16_t *a = CONVERT_TO_SHORTPTR(a8);
  56   const uint16_t *b = CONVERT_TO_SHORTPTR(b8);
  57   int64_t tsum = 0;
  58   uint64_t tsse = 0;
  59   for (int i = 0; i < h; ++i) {
  60     int32_t lsum = 0;
  61     for (int j = 0; j < w; ++j) {
  62       const int diff = a[j] - b[j];
  63       lsum += diff;
  64       tsse += (uint32_t)(diff * diff);
  65     }
  66     tsum += lsum;
  67     a += a_stride;
  68     b += b_stride;
  69   }
  70   *sum = tsum;
  71   *sse = tsse;
  72 }
  73
  74 static void encoder_highbd_8_variance(const uint8_t *a8, int a_stride,
  75                                       const uint8_t *b8, int b_stride, int w,
  76                                       int h, unsigned int *sse, int *sum) {
  77   uint64_t sse_long = 0;
  78   int64_t sum_long = 0;
  79   encoder_highbd_variance64(a8, a_stride, b8, b_stride, w, h, &sse_long,
  80                             &sum_long);
  81   *sse = (unsigned int)sse_long;
  82   *sum = (int)sum_long;
  83 }
  84
  85 static int64_t get_sse(const uint8_t *a, int a_stride, const uint8_t *b,
  86                        int b_stride, int width, int height) {
  87   const int dw = width % 16;
  88   const int dh = height % 16;
  89   int64_t total_sse = 0;
  90   unsigned int sse = 0;
  91   int sum = 0;
  92   int x, y;
  93
  94   if (dw > 0) {
  95     encoder_variance(&a[width - dw], a_stride, &b[width - dw], b_stride, dw,
  96                      height, &sse, &sum);
  97     total_sse += sse;
  98   }
  99
 100   if (dh > 0) {
 101     encoder_variance(&a[(height - dh) * a_stride], a_stride,
 102                      &b[(height - dh) * b_stride], b_stride, width - dw, dh,
 103                      &sse, &sum);
 104     total_sse += sse;
 105   }
 106
 107   for (y = 0; y < height / 16; ++y) {
 108     const uint8_t *pa = a;
 109     const uint8_t *pb = b;
 110     for (x = 0; x < width / 16; ++x) {
 111       aom_mse16x16(pa, a_stride, pb, b_stride, &sse);
 112       total_sse += sse;
 113
 114       pa += 16;
 115       pb += 16;
 116     }
 117
 118     a += 16 * a_stride;
 119     b += 16 * b_stride;
 120   }
 121
 122   return total_sse;
 123 }
 124
 125 static int64_t highbd_get_sse_shift(const uint8_t *a8, int a_stride,
 126                                     const uint8_t *b8, int b_stride, int width,
 127                                     int height, unsigned int input_shift) {
 128   const uint16_t *a = CONVERT_TO_SHORTPTR(a8);
 129   const uint16_t *b = CONVERT_TO_SHORTPTR(b8);
 130   int64_t total_sse = 0;
 131   int x, y;
 132   for (y = 0; y < height; ++y) {
 133     for (x = 0; x < width; ++x) {
 134       int64_t diff;
 135       diff = (a[x] >> input_shift) - (b[x] >> input_shift);
 136       total_sse += diff * diff;
 137     }
 138     a += a_stride;
 139     b += b_stride;
 140   }
 141   return total_sse;
 142 }
 143
 144 static int64_t highbd_get_sse(const uint8_t *a, int a_stride, const uint8_t *b,
 145                               int b_stride, int width, int height) {
 146   int64_t total_sse = 0;
 147   int x, y;
 148   const int dw = width % 16;
 149   const int dh = height % 16;
 150   unsigned int sse = 0;
 151   int sum = 0;
 152   if (dw > 0) {
 153     encoder_highbd_8_variance(&a[width - dw], a_stride, &b[width - dw],
 154                               b_stride, dw, height, &sse, &sum);
 155     total_sse += sse;
 156   }
 157   if (dh > 0) {
 158     encoder_highbd_8_variance(&a[(height - dh) * a_stride], a_stride,
 159                               &b[(height - dh) * b_stride], b_stride,
 160                               width - dw, dh, &sse, &sum);
 161     total_sse += sse;
 162   }
 163   for (y = 0; y < height / 16; ++y) {
 164     const uint8_t *pa = a;
 165     const uint8_t *pb = b;
 166     for (x = 0; x < width / 16; ++x) {
 167       aom_highbd_8_mse16x16(pa, a_stride, pb, b_stride, &sse);
 168       total_sse += sse;
 169       pa += 16;
 170       pb += 16;
 171     }
 172     a += 16 * a_stride;
 173     b += 16 * b_stride;
 174   }
 175   return total_sse;
 176 }
 177
 178 int64_t aom_get_y_sse_part(const YV12_BUFFER_CONFIG *a,
 179                            const YV12_BUFFER_CONFIG *b, int hstart, int width,
 180                            int vstart, int height) {
 181   return get_sse(a->y_buffer + vstart * a->y_stride + hstart, a->y_stride,
 182                  b->y_buffer + vstart * b->y_stride + hstart, b->y_stride,
 183                  width, height);
 184 }
 185
 186 int64_t aom_get_y_sse(const YV12_BUFFER_CONFIG *a,
 187                       const YV12_BUFFER_CONFIG *b) {
 188   assert(a->y_crop_width == b->y_crop_width);
 189   assert(a->y_crop_height == b->y_crop_height);
 190
 191   return get_sse(a->y_buffer, a->y_stride, b->y_buffer, b->y_stride,
 192                  a->y_crop_width, a->y_crop_height);
 193 }
 194
 195 int64_t aom_get_u_sse_part(const YV12_BUFFER_CONFIG *a,
 196                            const YV12_BUFFER_CONFIG *b, int hstart, int width,
 197                            int vstart, int height) {
 198   return get_sse(a->u_buffer + vstart * a->uv_stride + hstart, a->uv_stride,
 199                  b->u_buffer + vstart * b->uv_stride + hstart, b->uv_stride,
 200                  width, height);
 201 }
 202
 203 int64_t aom_get_u_sse(const YV12_BUFFER_CONFIG *a,
 204                       const YV12_BUFFER_CONFIG *b) {
 205   assert(a->uv_crop_width == b->uv_crop_width);
 206   assert(a->uv_crop_height == b->uv_crop_height);
 207
 208   return get_sse(a->u_buffer, a->uv_stride, b->u_buffer, b->uv_stride,
 209                  a->uv_crop_width, a->uv_crop_height);
 210 }
 211
 212 int64_t aom_get_v_sse_part(const YV12_BUFFER_CONFIG *a,
 213                            const YV12_BUFFER_CONFIG *b, int hstart, int width,
 214                            int vstart, int height) {
 215   return get_sse(a->v_buffer + vstart * a->uv_stride + hstart, a->uv_stride,
 216                  b->v_buffer + vstart * b->uv_stride + hstart, b->uv_stride,
 217                  width, height);
 218 }
 219
 220 int64_t aom_get_v_sse(const YV12_BUFFER_CONFIG *a,
 221                       const YV12_BUFFER_CONFIG *b) {
 222   assert(a->uv_crop_width == b->uv_crop_width);
 223   assert(a->uv_crop_height == b->uv_crop_height);
 224
 225   return get_sse(a->v_buffer, a->uv_stride, b->v_buffer, b->uv_stride,
 226                  a->uv_crop_width, a->uv_crop_height);
 227 }
 228
 229 int64_t aom_highbd_get_y_sse_part(const YV12_BUFFER_CONFIG *a,
 230                                   const YV12_BUFFER_CONFIG *b, int hstart,
 231                                   int width, int vstart, int height) {
 232   return highbd_get_sse(
 233       a->y_buffer + vstart * a->y_stride + hstart, a->y_stride,
 234       b->y_buffer + vstart * b->y_stride + hstart, b->y_stride, width, height);
 235 }
 236
 237 int64_t aom_highbd_get_y_sse(const YV12_BUFFER_CONFIG *a,
 238                              const YV12_BUFFER_CONFIG *b) {
 239   assert(a->y_crop_width == b->y_crop_width);
 240   assert(a->y_crop_height == b->y_crop_height);
 241   assert((a->flags & YV12_FLAG_HIGHBITDEPTH) != 0);
 242   assert((b->flags & YV12_FLAG_HIGHBITDEPTH) != 0);
 243
 244   return highbd_get_sse(a->y_buffer, a->y_stride, b->y_buffer, b->y_stride,
 245                         a->y_crop_width, a->y_crop_height);
 246 }
 247
 248 int64_t aom_highbd_get_u_sse_part(const YV12_BUFFER_CONFIG *a,
 249                                   const YV12_BUFFER_CONFIG *b, int hstart,
 250                                   int width, int vstart, int height) {
 251   return highbd_get_sse(a->u_buffer + vstart * a->uv_stride + hstart,
 252                         a->uv_stride,
 253                         b->u_buffer + vstart * b->uv_stride + hstart,
 254                         b->uv_stride, width, height);
 255 }
 256
 257 int64_t aom_highbd_get_u_sse(const YV12_BUFFER_CONFIG *a,
 258                              const YV12_BUFFER_CONFIG *b) {
 259   assert(a->uv_crop_width == b->uv_crop_width);
 260   assert(a->uv_crop_height == b->uv_crop_height);
 261   assert((a->flags & YV12_FLAG_HIGHBITDEPTH) != 0);
 262   assert((b->flags & YV12_FLAG_HIGHBITDEPTH) != 0);
 263
 264   return highbd_get_sse(a->u_buffer, a->uv_stride, b->u_buffer, b->uv_stride,
 265                         a->uv_crop_width, a->uv_crop_height);
 266 }
 267
 268 int64_t aom_highbd_get_v_sse_part(const YV12_BUFFER_CONFIG *a,
 269                                   const YV12_BUFFER_CONFIG *b, int hstart,
 270                                   int width, int vstart, int height) {
 271   return highbd_get_sse(a->v_buffer + vstart * a->uv_stride + hstart,
 272                         a->uv_stride,
 273                         b->v_buffer + vstart * b->uv_stride + hstart,
 274                         b->uv_stride, width, height);
 275 }
 276
 277 int64_t aom_highbd_get_v_sse(const YV12_BUFFER_CONFIG *a,
 278                              const YV12_BUFFER_CONFIG *b) {
 279   assert(a->uv_crop_width == b->uv_crop_width);
 280   assert(a->uv_crop_height == b->uv_crop_height);
 281   assert((a->flags & YV12_FLAG_HIGHBITDEPTH) != 0);
 282   assert((b->flags & YV12_FLAG_HIGHBITDEPTH) != 0);
 283
 284   return highbd_get_sse(a->v_buffer, a->uv_stride, b->v_buffer, b->uv_stride,
 285                         a->uv_crop_width, a->uv_crop_height);
 286 }
 287
 288 int64_t aom_get_sse_plane(const YV12_BUFFER_CONFIG *a,
 289                           const YV12_BUFFER_CONFIG *b, int plane, int highbd) {
 290   if (highbd) {
 291     switch (plane) {
 292       case 0: return aom_highbd_get_y_sse(a, b);
 293       case 1: return aom_highbd_get_u_sse(a, b);
 294       case 2: return aom_highbd_get_v_sse(a, b);
 295       default: assert(plane >= 0 && plane <= 2); return 0;
 296     }
 297   }
 298   (void)highbd;
 299   switch (plane) {
 300     case 0: return aom_get_y_sse(a, b);
 301     case 1: return aom_get_u_sse(a, b);
 302     case 2: return aom_get_v_sse(a, b);
 303     default: assert(plane >= 0 && plane <= 2); return 0;
 304   }
 305 }
 306
 307 void aom_calc_highbd_psnr(const YV12_BUFFER_CONFIG *a,
 308                           const YV12_BUFFER_CONFIG *b, PSNR_STATS *psnr,
 309                           uint32_t bit_depth, uint32_t in_bit_depth) {
 310   const int widths[3] = { a->y_crop_width, a->uv_crop_width, a->uv_crop_width };
 311   const int heights[3] = { a->y_crop_height, a->uv_crop_height,
 312                            a->uv_crop_height };
 313   const int a_strides[3] = { a->y_stride, a->uv_stride, a->uv_stride };
 314   const int b_strides[3] = { b->y_stride, b->uv_stride, b->uv_stride };
 315   int i;
 316   uint64_t total_sse = 0;
 317   uint32_t total_samples = 0;
 318   const double peak = (double)((1 << in_bit_depth) - 1);
 319   const unsigned int input_shift = bit_depth - in_bit_depth;
 320
 321   for (i = 0; i < 3; ++i) {
 322     const int w = widths[i];
 323     const int h = heights[i];
 324     const uint32_t samples = w * h;
 325     uint64_t sse;
 326     if (a->flags & YV12_FLAG_HIGHBITDEPTH) {
 327       if (input_shift) {
 328         sse = highbd_get_sse_shift(a->buffers[i], a_strides[i], b->buffers[i],
 329                                    b_strides[i], w, h, input_shift);
 330       } else {
 331         sse = highbd_get_sse(a->buffers[i], a_strides[i], b->buffers[i],
 332                              b_strides[i], w, h);
 333       }
 334     } else {
 335       sse = get_sse(a->buffers[i], a_strides[i], b->buffers[i], b_strides[i], w,
 336                     h);
 337     }
 338     psnr->sse[1 + i] = sse;
 339     psnr->samples[1 + i] = samples;
 340     psnr->psnr[1 + i] = aom_sse_to_psnr(samples, peak, (double)sse);
 341
 342     total_sse += sse;
 343     total_samples += samples;
 344   }
 345
 346   psnr->sse[0] = total_sse;
 347   psnr->samples[0] = total_samples;
 348   psnr->psnr[0] =
 349       aom_sse_to_psnr((double)total_samples, peak, (double)total_sse);
 350 }
 351
 352 void aom_calc_psnr(const YV12_BUFFER_CONFIG *a, const YV12_BUFFER_CONFIG *b,
 353                    PSNR_STATS *psnr) {
 354   static const double peak = 255.0;
 355   const int widths[3] = { a->y_crop_width, a->uv_crop_width, a->uv_crop_width };
 356   const int heights[3] = { a->y_crop_height, a->uv_crop_height,
 357                            a->uv_crop_height };
 358   const int a_strides[3] = { a->y_stride, a->uv_stride, a->uv_stride };
 359   const int b_strides[3] = { b->y_stride, b->uv_stride, b->uv_stride };
 360   int i;
 361   uint64_t total_sse = 0;
 362   uint32_t total_samples = 0;
 363
 364   for (i = 0; i < 3; ++i) {
 365     const int w = widths[i];
 366     const int h = heights[i];
 367     const uint32_t samples = w * h;
 368     const uint64_t sse =
 369         get_sse(a->buffers[i], a_strides[i], b->buffers[i], b_strides[i], w, h);
 370     psnr->sse[1 + i] = sse;
 371     psnr->samples[1 + i] = samples;
 372     psnr->psnr[1 + i] = aom_sse_to_psnr(samples, peak, (double)sse);
 373
 374     total_sse += sse;
 375     total_samples += samples;
 376   }
 377
 378   psnr->sse[0] = total_sse;
 379   psnr->samples[0] = total_samples;
 380   psnr->psnr[0] =
 381       aom_sse_to_psnr((double)total_samples, peak, (double)total_sse);
 382 }