psnr.c (16195B)
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_config.h" 16 #include "config/aom_dsp_rtcd.h" 17 18 #include "aom_dsp/psnr.h" 19 #include "aom_scale/yv12config.h" 20 21 #if CONFIG_INTERNAL_STATS 22 #define STATIC 23 #else 24 #define STATIC static 25 #endif // CONFIG_INTERNAL_STATS 26 27 STATIC double aom_sse_to_psnr(double samples, double peak, double sse) { 28 if (sse > 0.0) { 29 const double psnr = 10.0 * log10(samples * peak * peak / sse); 30 return psnr > MAX_PSNR ? MAX_PSNR : psnr; 31 } else { 32 return MAX_PSNR; 33 } 34 } 35 36 #undef STATIC 37 38 static int64_t encoder_sse(const uint8_t *a, int a_stride, const uint8_t *b, 39 int b_stride, int w, int h) { 40 int i, j; 41 int64_t sse = 0; 42 43 for (i = 0; i < h; i++) { 44 for (j = 0; j < w; j++) { 45 const int diff = a[j] - b[j]; 46 sse += diff * diff; 47 } 48 49 a += a_stride; 50 b += b_stride; 51 } 52 return sse; 53 } 54 55 #if CONFIG_AV1_HIGHBITDEPTH 56 static int64_t encoder_highbd_sse(const uint8_t *a8, int a_stride, 57 const uint8_t *b8, int b_stride, int w, 58 int h) { 59 const uint16_t *a = CONVERT_TO_SHORTPTR(a8); 60 const uint16_t *b = CONVERT_TO_SHORTPTR(b8); 61 int64_t sse = 0; 62 for (int i = 0; i < h; ++i) { 63 for (int j = 0; j < w; ++j) { 64 const int diff = a[j] - b[j]; 65 sse += diff * diff; 66 } 67 a += a_stride; 68 b += b_stride; 69 } 70 return sse; 71 } 72 73 #endif // CONFIG_AV1_HIGHBITDEPTH 74 75 static int64_t get_sse(const uint8_t *a, int a_stride, const uint8_t *b, 76 int b_stride, int width, int height) { 77 const int dw = width % 16; 78 const int dh = height % 16; 79 int64_t total_sse = 0; 80 int x, y; 81 82 if (dw > 0) { 83 total_sse += encoder_sse(&a[width - dw], a_stride, &b[width - dw], b_stride, 84 dw, height); 85 } 86 87 if (dh > 0) { 88 total_sse += 89 encoder_sse(&a[(height - dh) * a_stride], a_stride, 90 &b[(height - dh) * b_stride], b_stride, width - dw, dh); 91 } 92 93 for (y = 0; y < height / 16; ++y) { 94 const uint8_t *pa = a; 95 const uint8_t *pb = b; 96 for (x = 0; x < width / 16; ++x) { 97 total_sse += aom_sse(pa, a_stride, pb, b_stride, 16, 16); 98 99 pa += 16; 100 pb += 16; 101 } 102 103 a += 16 * a_stride; 104 b += 16 * b_stride; 105 } 106 107 return total_sse; 108 } 109 110 #if CONFIG_AV1_HIGHBITDEPTH 111 static int64_t highbd_get_sse_shift(const uint8_t *a8, int a_stride, 112 const uint8_t *b8, int b_stride, int width, 113 int height, unsigned int input_shift) { 114 const uint16_t *a = CONVERT_TO_SHORTPTR(a8); 115 const uint16_t *b = CONVERT_TO_SHORTPTR(b8); 116 int64_t total_sse = 0; 117 int x, y; 118 for (y = 0; y < height; ++y) { 119 for (x = 0; x < width; ++x) { 120 int64_t diff; 121 diff = (a[x] >> input_shift) - (b[x] >> input_shift); 122 total_sse += diff * diff; 123 } 124 a += a_stride; 125 b += b_stride; 126 } 127 return total_sse; 128 } 129 130 static int64_t highbd_get_sse(const uint8_t *a, int a_stride, const uint8_t *b, 131 int b_stride, int width, int height) { 132 int64_t total_sse = 0; 133 int x, y; 134 const int dw = width % 16; 135 const int dh = height % 16; 136 137 if (dw > 0) { 138 total_sse += encoder_highbd_sse(&a[width - dw], a_stride, &b[width - dw], 139 b_stride, dw, height); 140 } 141 if (dh > 0) { 142 total_sse += encoder_highbd_sse(&a[(height - dh) * a_stride], a_stride, 143 &b[(height - dh) * b_stride], b_stride, 144 width - dw, dh); 145 } 146 147 for (y = 0; y < height / 16; ++y) { 148 const uint8_t *pa = a; 149 const uint8_t *pb = b; 150 for (x = 0; x < width / 16; ++x) { 151 total_sse += aom_highbd_sse(pa, a_stride, pb, b_stride, 16, 16); 152 pa += 16; 153 pb += 16; 154 } 155 a += 16 * a_stride; 156 b += 16 * b_stride; 157 } 158 return total_sse; 159 } 160 #endif // CONFIG_AV1_HIGHBITDEPTH 161 162 uint64_t aom_get_y_var(const YV12_BUFFER_CONFIG *a, int hstart, int width, 163 int vstart, int height) { 164 return aom_var_2d_u8(a->y_buffer + vstart * a->y_stride + hstart, a->y_stride, 165 width, height) / 166 (width * height); 167 } 168 169 uint64_t aom_get_u_var(const YV12_BUFFER_CONFIG *a, int hstart, int width, 170 int vstart, int height) { 171 return aom_var_2d_u8(a->u_buffer + vstart * a->uv_stride + hstart, 172 a->uv_stride, width, height) / 173 (width * height); 174 } 175 176 uint64_t aom_get_v_var(const YV12_BUFFER_CONFIG *a, int hstart, int width, 177 int vstart, int height) { 178 return aom_var_2d_u8(a->v_buffer + vstart * a->uv_stride + hstart, 179 a->uv_stride, width, height) / 180 (width * height); 181 } 182 183 int64_t aom_get_y_sse_part(const YV12_BUFFER_CONFIG *a, 184 const YV12_BUFFER_CONFIG *b, int hstart, int width, 185 int vstart, int height) { 186 return get_sse(a->y_buffer + vstart * a->y_stride + hstart, a->y_stride, 187 b->y_buffer + vstart * b->y_stride + hstart, b->y_stride, 188 width, height); 189 } 190 191 int64_t aom_get_y_sse(const YV12_BUFFER_CONFIG *a, 192 const YV12_BUFFER_CONFIG *b) { 193 assert(a->y_crop_width == b->y_crop_width); 194 assert(a->y_crop_height == b->y_crop_height); 195 196 return get_sse(a->y_buffer, a->y_stride, b->y_buffer, b->y_stride, 197 a->y_crop_width, a->y_crop_height); 198 } 199 200 int64_t aom_get_u_sse_part(const YV12_BUFFER_CONFIG *a, 201 const YV12_BUFFER_CONFIG *b, int hstart, int width, 202 int vstart, int height) { 203 return get_sse(a->u_buffer + vstart * a->uv_stride + hstart, a->uv_stride, 204 b->u_buffer + vstart * b->uv_stride + hstart, b->uv_stride, 205 width, height); 206 } 207 208 int64_t aom_get_u_sse(const YV12_BUFFER_CONFIG *a, 209 const YV12_BUFFER_CONFIG *b) { 210 assert(a->uv_crop_width == b->uv_crop_width); 211 assert(a->uv_crop_height == b->uv_crop_height); 212 213 return get_sse(a->u_buffer, a->uv_stride, b->u_buffer, b->uv_stride, 214 a->uv_crop_width, a->uv_crop_height); 215 } 216 217 int64_t aom_get_v_sse_part(const YV12_BUFFER_CONFIG *a, 218 const YV12_BUFFER_CONFIG *b, int hstart, int width, 219 int vstart, int height) { 220 return get_sse(a->v_buffer + vstart * a->uv_stride + hstart, a->uv_stride, 221 b->v_buffer + vstart * b->uv_stride + hstart, b->uv_stride, 222 width, height); 223 } 224 225 int64_t aom_get_v_sse(const YV12_BUFFER_CONFIG *a, 226 const YV12_BUFFER_CONFIG *b) { 227 assert(a->uv_crop_width == b->uv_crop_width); 228 assert(a->uv_crop_height == b->uv_crop_height); 229 230 return get_sse(a->v_buffer, a->uv_stride, b->v_buffer, b->uv_stride, 231 a->uv_crop_width, a->uv_crop_height); 232 } 233 234 #if CONFIG_AV1_HIGHBITDEPTH 235 uint64_t aom_highbd_get_y_var(const YV12_BUFFER_CONFIG *a, int hstart, 236 int width, int vstart, int height) { 237 return aom_var_2d_u16(a->y_buffer + vstart * a->y_stride + hstart, 238 a->y_stride, width, height) / 239 (width * height); 240 } 241 242 uint64_t aom_highbd_get_u_var(const YV12_BUFFER_CONFIG *a, int hstart, 243 int width, int vstart, int height) { 244 return aom_var_2d_u16(a->u_buffer + vstart * a->uv_stride + hstart, 245 a->uv_stride, width, height) / 246 (width * height); 247 } 248 249 uint64_t aom_highbd_get_v_var(const YV12_BUFFER_CONFIG *a, int hstart, 250 int width, int vstart, int height) { 251 return aom_var_2d_u16(a->v_buffer + vstart * a->uv_stride + hstart, 252 a->uv_stride, width, height) / 253 (width * height); 254 } 255 256 int64_t aom_highbd_get_y_sse_part(const YV12_BUFFER_CONFIG *a, 257 const YV12_BUFFER_CONFIG *b, int hstart, 258 int width, int vstart, int height) { 259 return highbd_get_sse( 260 a->y_buffer + vstart * a->y_stride + hstart, a->y_stride, 261 b->y_buffer + vstart * b->y_stride + hstart, b->y_stride, width, height); 262 } 263 264 int64_t aom_highbd_get_y_sse(const YV12_BUFFER_CONFIG *a, 265 const YV12_BUFFER_CONFIG *b) { 266 assert(a->y_crop_width == b->y_crop_width); 267 assert(a->y_crop_height == b->y_crop_height); 268 assert((a->flags & YV12_FLAG_HIGHBITDEPTH) != 0); 269 assert((b->flags & YV12_FLAG_HIGHBITDEPTH) != 0); 270 271 return highbd_get_sse(a->y_buffer, a->y_stride, b->y_buffer, b->y_stride, 272 a->y_crop_width, a->y_crop_height); 273 } 274 275 int64_t aom_highbd_get_u_sse_part(const YV12_BUFFER_CONFIG *a, 276 const YV12_BUFFER_CONFIG *b, int hstart, 277 int width, int vstart, int height) { 278 return highbd_get_sse(a->u_buffer + vstart * a->uv_stride + hstart, 279 a->uv_stride, 280 b->u_buffer + vstart * b->uv_stride + hstart, 281 b->uv_stride, width, height); 282 } 283 284 int64_t aom_highbd_get_u_sse(const YV12_BUFFER_CONFIG *a, 285 const YV12_BUFFER_CONFIG *b) { 286 assert(a->uv_crop_width == b->uv_crop_width); 287 assert(a->uv_crop_height == b->uv_crop_height); 288 assert((a->flags & YV12_FLAG_HIGHBITDEPTH) != 0); 289 assert((b->flags & YV12_FLAG_HIGHBITDEPTH) != 0); 290 291 return highbd_get_sse(a->u_buffer, a->uv_stride, b->u_buffer, b->uv_stride, 292 a->uv_crop_width, a->uv_crop_height); 293 } 294 295 int64_t aom_highbd_get_v_sse_part(const YV12_BUFFER_CONFIG *a, 296 const YV12_BUFFER_CONFIG *b, int hstart, 297 int width, int vstart, int height) { 298 return highbd_get_sse(a->v_buffer + vstart * a->uv_stride + hstart, 299 a->uv_stride, 300 b->v_buffer + vstart * b->uv_stride + hstart, 301 b->uv_stride, width, height); 302 } 303 304 int64_t aom_highbd_get_v_sse(const YV12_BUFFER_CONFIG *a, 305 const YV12_BUFFER_CONFIG *b) { 306 assert(a->uv_crop_width == b->uv_crop_width); 307 assert(a->uv_crop_height == b->uv_crop_height); 308 assert((a->flags & YV12_FLAG_HIGHBITDEPTH) != 0); 309 assert((b->flags & YV12_FLAG_HIGHBITDEPTH) != 0); 310 311 return highbd_get_sse(a->v_buffer, a->uv_stride, b->v_buffer, b->uv_stride, 312 a->uv_crop_width, a->uv_crop_height); 313 } 314 #endif // CONFIG_AV1_HIGHBITDEPTH 315 316 int64_t aom_get_sse_plane(const YV12_BUFFER_CONFIG *a, 317 const YV12_BUFFER_CONFIG *b, int plane, int highbd) { 318 #if CONFIG_AV1_HIGHBITDEPTH 319 if (highbd) { 320 switch (plane) { 321 case 0: return aom_highbd_get_y_sse(a, b); 322 case 1: return aom_highbd_get_u_sse(a, b); 323 case 2: return aom_highbd_get_v_sse(a, b); 324 default: assert(plane >= 0 && plane <= 2); return 0; 325 } 326 } else { 327 switch (plane) { 328 case 0: return aom_get_y_sse(a, b); 329 case 1: return aom_get_u_sse(a, b); 330 case 2: return aom_get_v_sse(a, b); 331 default: assert(plane >= 0 && plane <= 2); return 0; 332 } 333 } 334 #else 335 (void)highbd; 336 switch (plane) { 337 case 0: return aom_get_y_sse(a, b); 338 case 1: return aom_get_u_sse(a, b); 339 case 2: return aom_get_v_sse(a, b); 340 default: assert(plane >= 0 && plane <= 2); return 0; 341 } 342 #endif 343 } 344 345 #if CONFIG_AV1_HIGHBITDEPTH 346 void aom_calc_highbd_psnr(const YV12_BUFFER_CONFIG *a, 347 const YV12_BUFFER_CONFIG *b, PSNR_STATS *psnr, 348 uint32_t bit_depth, uint32_t in_bit_depth) { 349 assert(a->y_crop_width == b->y_crop_width); 350 assert(a->y_crop_height == b->y_crop_height); 351 assert(a->uv_crop_width == b->uv_crop_width); 352 assert(a->uv_crop_height == b->uv_crop_height); 353 const int widths[3] = { a->y_crop_width, a->uv_crop_width, a->uv_crop_width }; 354 const int heights[3] = { a->y_crop_height, a->uv_crop_height, 355 a->uv_crop_height }; 356 const int a_strides[3] = { a->y_stride, a->uv_stride, a->uv_stride }; 357 const int b_strides[3] = { b->y_stride, b->uv_stride, b->uv_stride }; 358 int i; 359 uint64_t total_sse = 0; 360 uint32_t total_samples = 0; 361 #if CONFIG_LIBVMAF_PSNR_PEAK 362 double peak = (double)(255 << (in_bit_depth - 8)); 363 #else 364 double peak = (double)((1 << in_bit_depth) - 1); 365 #endif // CONFIG_LIBVMAF_PSNR_PEAK 366 const unsigned int input_shift = bit_depth - in_bit_depth; 367 368 for (i = 0; i < 3; ++i) { 369 const int w = widths[i]; 370 const int h = heights[i]; 371 const uint32_t samples = w * h; 372 uint64_t sse; 373 if (a->flags & YV12_FLAG_HIGHBITDEPTH) { 374 if (input_shift) { 375 sse = highbd_get_sse_shift(a->buffers[i], a_strides[i], b->buffers[i], 376 b_strides[i], w, h, input_shift); 377 } else { 378 sse = highbd_get_sse(a->buffers[i], a_strides[i], b->buffers[i], 379 b_strides[i], w, h); 380 } 381 } else { 382 sse = get_sse(a->buffers[i], a_strides[i], b->buffers[i], b_strides[i], w, 383 h); 384 } 385 psnr->sse[1 + i] = sse; 386 psnr->samples[1 + i] = samples; 387 psnr->psnr[1 + i] = aom_sse_to_psnr(samples, peak, (double)sse); 388 389 total_sse += sse; 390 total_samples += samples; 391 } 392 393 psnr->sse[0] = total_sse; 394 psnr->samples[0] = total_samples; 395 psnr->psnr[0] = 396 aom_sse_to_psnr((double)total_samples, peak, (double)total_sse); 397 398 // Compute PSNR based on stream bit depth 399 if ((a->flags & YV12_FLAG_HIGHBITDEPTH) && (in_bit_depth < bit_depth)) { 400 #if CONFIG_LIBVMAF_PSNR_PEAK 401 peak = (double)(255 << (bit_depth - 8)); 402 #else 403 peak = (double)((1 << bit_depth) - 1); 404 #endif // CONFIG_LIBVMAF_PSNR_PEAK 405 total_sse = 0; 406 total_samples = 0; 407 for (i = 0; i < 3; ++i) { 408 const int w = widths[i]; 409 const int h = heights[i]; 410 const uint32_t samples = w * h; 411 uint64_t sse; 412 sse = highbd_get_sse(a->buffers[i], a_strides[i], b->buffers[i], 413 b_strides[i], w, h); 414 psnr->sse_hbd[1 + i] = sse; 415 psnr->samples_hbd[1 + i] = samples; 416 psnr->psnr_hbd[1 + i] = aom_sse_to_psnr(samples, peak, (double)sse); 417 total_sse += sse; 418 total_samples += samples; 419 } 420 421 psnr->sse_hbd[0] = total_sse; 422 psnr->samples_hbd[0] = total_samples; 423 psnr->psnr_hbd[0] = 424 aom_sse_to_psnr((double)total_samples, peak, (double)total_sse); 425 } 426 } 427 #endif 428 429 void aom_calc_psnr(const YV12_BUFFER_CONFIG *a, const YV12_BUFFER_CONFIG *b, 430 PSNR_STATS *psnr) { 431 assert(a->y_crop_width == b->y_crop_width); 432 assert(a->y_crop_height == b->y_crop_height); 433 assert(a->uv_crop_width == b->uv_crop_width); 434 assert(a->uv_crop_height == b->uv_crop_height); 435 static const double peak = 255.0; 436 const int widths[3] = { a->y_crop_width, a->uv_crop_width, a->uv_crop_width }; 437 const int heights[3] = { a->y_crop_height, a->uv_crop_height, 438 a->uv_crop_height }; 439 const int a_strides[3] = { a->y_stride, a->uv_stride, a->uv_stride }; 440 const int b_strides[3] = { b->y_stride, b->uv_stride, b->uv_stride }; 441 int i; 442 uint64_t total_sse = 0; 443 uint32_t total_samples = 0; 444 445 for (i = 0; i < 3; ++i) { 446 const int w = widths[i]; 447 const int h = heights[i]; 448 const uint32_t samples = w * h; 449 const uint64_t sse = 450 get_sse(a->buffers[i], a_strides[i], b->buffers[i], b_strides[i], w, h); 451 psnr->sse[1 + i] = sse; 452 psnr->samples[1 + i] = samples; 453 psnr->psnr[1 + i] = aom_sse_to_psnr(samples, peak, (double)sse); 454 455 total_sse += sse; 456 total_samples += samples; 457 } 458 459 psnr->sse[0] = total_sse; 460 psnr->samples[0] = total_samples; 461 psnr->psnr[0] = 462 aom_sse_to_psnr((double)total_samples, peak, (double)total_sse); 463 }