#include <assert.h>
#include <math.h>
#include "config/aom_dsp_rtcd.h"
#include "aom_dsp/ssim.h"
#include "aom_ports/mem.h"
void aom_ssim_parms_8x8_c(const uint8_t *s, int sp, const uint8_t *r, int rp,
uint32_t *sum_s, uint32_t *sum_r, uint32_t *sum_sq_s,
uint32_t *sum_sq_r, uint32_t *sum_sxr) { … }
static const int64_t cc1 = …;
static const int64_t cc2 = …;
static const int64_t cc1_10 = …;
static const int64_t cc2_10 = …;
static const int64_t cc1_12 = …;
static const int64_t cc2_12 = …;
static double similarity(uint32_t sum_s, uint32_t sum_r, uint32_t sum_sq_s,
uint32_t sum_sq_r, uint32_t sum_sxr, int count,
uint32_t bd) { … }
static double ssim_8x8(const uint8_t *s, int sp, const uint8_t *r, int rp) { … }
double aom_ssim2(const uint8_t *img1, const uint8_t *img2, int stride_img1,
int stride_img2, int width, int height) { … }
#if CONFIG_INTERNAL_STATS
void aom_lowbd_calc_ssim(const YV12_BUFFER_CONFIG *source,
const YV12_BUFFER_CONFIG *dest, double *weight,
double *fast_ssim) {
double abc[3];
for (int i = 0; i < 3; ++i) {
const int is_uv = i > 0;
abc[i] = aom_ssim2(source->buffers[i], dest->buffers[i],
source->strides[is_uv], dest->strides[is_uv],
source->crop_widths[is_uv], source->crop_heights[is_uv]);
}
*weight = 1;
*fast_ssim = abc[0] * .8 + .1 * (abc[1] + abc[2]);
}
static double ssimv_similarity(const Ssimv *sv, int64_t n) {
const int64_t c1 = (cc1 * n * n) >> 12;
const int64_t c2 = (cc2 * n * n) >> 12;
const double l = 1.0 * (2 * sv->sum_s * sv->sum_r + c1) /
(sv->sum_s * sv->sum_s + sv->sum_r * sv->sum_r + c1);
const double v = (2.0 * n * sv->sum_sxr - 2 * sv->sum_s * sv->sum_r + c2) /
(n * sv->sum_sq_s - sv->sum_s * sv->sum_s +
n * sv->sum_sq_r - sv->sum_r * sv->sum_r + c2);
return l * v;
}
static double ssimv_similarity2(const Ssimv *sv, int64_t n) {
const int64_t c1 = (cc1 * n * n) >> 12;
const int64_t c2 = (cc2 * n * n) >> 12;
const double mean_diff = (1.0 * sv->sum_s - sv->sum_r) / n;
const double l = (255 * 255 - mean_diff * mean_diff + c1) / (255 * 255 + c1);
const double v = (2.0 * n * sv->sum_sxr - 2 * sv->sum_s * sv->sum_r + c2) /
(n * sv->sum_sq_s - sv->sum_s * sv->sum_s +
n * sv->sum_sq_r - sv->sum_r * sv->sum_r + c2);
return l * v;
}
static void ssimv_parms(uint8_t *img1, int img1_pitch, uint8_t *img2,
int img2_pitch, Ssimv *sv) {
aom_ssim_parms_8x8(img1, img1_pitch, img2, img2_pitch, &sv->sum_s, &sv->sum_r,
&sv->sum_sq_s, &sv->sum_sq_r, &sv->sum_sxr);
}
double aom_get_ssim_metrics(uint8_t *img1, int img1_pitch, uint8_t *img2,
int img2_pitch, int width, int height, Ssimv *sv2,
Metrics *m, int do_inconsistency) {
double dssim_total = 0;
double ssim_total = 0;
double ssim2_total = 0;
double inconsistency_total = 0;
int i, j;
int c = 0;
double norm;
double old_ssim_total = 0;
for (i = 0; i < height;
i += 4, img1 += img1_pitch * 4, img2 += img2_pitch * 4) {
for (j = 0; j < width; j += 4, ++c) {
Ssimv sv = { 0, 0, 0, 0, 0, 0 };
double ssim;
double ssim2;
double dssim;
uint32_t var_new;
uint32_t var_old;
uint32_t mean_new;
uint32_t mean_old;
double ssim_new;
double ssim_old;
if (j + 8 <= width && i + 8 <= height) {
ssimv_parms(img1 + j, img1_pitch, img2 + j, img2_pitch, &sv);
}
ssim = ssimv_similarity(&sv, 64);
ssim2 = ssimv_similarity2(&sv, 64);
sv.ssim = ssim2;
dssim = 255 * 255 * (1 - ssim2) / 2;
var_new = sv.sum_sq_s - sv.sum_s * sv.sum_s / 64;
var_old = sv2[c].sum_sq_s - sv2[c].sum_s * sv2[c].sum_s / 64;
mean_new = sv.sum_s;
mean_old = sv2[c].sum_s;
ssim_new = sv.ssim;
ssim_old = sv2[c].ssim;
if (do_inconsistency) {
static const double kScaling = 4. * 4 * 255 * 255;
static const double c1 = 1, c2 = 1, c3 = 1;
const double variance_term =
(2.0 * var_old * var_new + c1) /
(1.0 * var_old * var_old + 1.0 * var_new * var_new + c1);
const double mean_term =
(2.0 * mean_old * mean_new + c2) /
(1.0 * mean_old * mean_old + 1.0 * mean_new * mean_new + c2);
double ssim_term =
pow((2.0 * ssim_old * ssim_new + c3) /
(ssim_old * ssim_old + ssim_new * ssim_new + c3),
5);
double this_inconsistency;
if (ssim_term > 1) ssim_term = 1;
this_inconsistency = (1 - ssim_term) * variance_term * mean_term;
this_inconsistency *= kScaling;
inconsistency_total += this_inconsistency;
}
sv2[c] = sv;
ssim_total += ssim;
ssim2_total += ssim2;
dssim_total += dssim;
old_ssim_total += ssim_old;
}
old_ssim_total += 0;
}
norm = 1. / (width / 4) / (height / 4);
ssim_total *= norm;
ssim2_total *= norm;
m->ssim2 = ssim2_total;
m->ssim = ssim_total;
if (old_ssim_total == 0) inconsistency_total = 0;
m->ssimc = inconsistency_total;
m->dssim = dssim_total;
return inconsistency_total;
}
#endif
#if CONFIG_AV1_HIGHBITDEPTH
void aom_highbd_ssim_parms_8x8_c(const uint16_t *s, int sp, const uint16_t *r,
int rp, uint32_t *sum_s, uint32_t *sum_r,
uint32_t *sum_sq_s, uint32_t *sum_sq_r,
uint32_t *sum_sxr) {
int i, j;
for (i = 0; i < 8; i++, s += sp, r += rp) {
for (j = 0; j < 8; j++) {
*sum_s += s[j];
*sum_r += r[j];
*sum_sq_s += s[j] * s[j];
*sum_sq_r += r[j] * r[j];
*sum_sxr += s[j] * r[j];
}
}
}
static double highbd_ssim_8x8(const uint16_t *s, int sp, const uint16_t *r,
int rp, uint32_t bd, uint32_t shift) {
uint32_t sum_s = 0, sum_r = 0, sum_sq_s = 0, sum_sq_r = 0, sum_sxr = 0;
aom_highbd_ssim_parms_8x8(s, sp, r, rp, &sum_s, &sum_r, &sum_sq_s, &sum_sq_r,
&sum_sxr);
return similarity(sum_s >> shift, sum_r >> shift, sum_sq_s >> (2 * shift),
sum_sq_r >> (2 * shift), sum_sxr >> (2 * shift), 64, bd);
}
double aom_highbd_ssim2(const uint8_t *img1, const uint8_t *img2,
int stride_img1, int stride_img2, int width, int height,
uint32_t bd, uint32_t shift) {
int i, j;
int samples = 0;
double ssim_total = 0;
for (i = 0; i <= height - 8;
i += 4, img1 += stride_img1 * 4, img2 += stride_img2 * 4) {
for (j = 0; j <= width - 8; j += 4) {
double v = highbd_ssim_8x8(CONVERT_TO_SHORTPTR(img1 + j), stride_img1,
CONVERT_TO_SHORTPTR(img2 + j), stride_img2, bd,
shift);
ssim_total += v;
samples++;
}
}
ssim_total /= samples;
return ssim_total;
}
#if CONFIG_INTERNAL_STATS
void aom_highbd_calc_ssim(const YV12_BUFFER_CONFIG *source,
const YV12_BUFFER_CONFIG *dest, double *weight,
uint32_t bd, uint32_t in_bd, double *fast_ssim) {
assert(bd >= in_bd);
uint32_t shift = bd - in_bd;
double abc[3];
for (int i = 0; i < 3; ++i) {
const int is_uv = i > 0;
abc[i] = aom_highbd_ssim2(source->buffers[i], dest->buffers[i],
source->strides[is_uv], dest->strides[is_uv],
source->crop_widths[is_uv],
source->crop_heights[is_uv], in_bd, shift);
}
weight[0] = 1;
fast_ssim[0] = abc[0] * .8 + .1 * (abc[1] + abc[2]);
if (bd > in_bd) {
shift = 0;
for (int i = 0; i < 3; ++i) {
const int is_uv = i > 0;
abc[i] = aom_highbd_ssim2(source->buffers[i], dest->buffers[i],
source->strides[is_uv], dest->strides[is_uv],
source->crop_widths[is_uv],
source->crop_heights[is_uv], bd, shift);
}
weight[1] = 1;
fast_ssim[1] = abc[0] * .8 + .1 * (abc[1] + abc[2]);
}
}
#endif
#endif
#if CONFIG_INTERNAL_STATS
void aom_calc_ssim(const YV12_BUFFER_CONFIG *orig,
const YV12_BUFFER_CONFIG *recon, const uint32_t bit_depth,
const uint32_t in_bit_depth, int is_hbd, double *weight,
double *frame_ssim2) {
#if CONFIG_AV1_HIGHBITDEPTH
if (is_hbd) {
aom_highbd_calc_ssim(orig, recon, weight, bit_depth, in_bit_depth,
frame_ssim2);
return;
}
#else
(void)bit_depth;
(void)in_bit_depth;
(void)is_hbd;
#endif
aom_lowbd_calc_ssim(orig, recon, weight, frame_ssim2);
}
#endif