in libvmaf/src/feature/integer_adm.c [1266:1645]
static float adm_cm(AdmBuffer *buf, int w, int h, int src_stride, int csf_a_stride,
double adm_norm_view_dist, int adm_ref_display_height)
{
const adm_dwt_band_t *src = &buf->decouple_r;
const adm_dwt_band_t *csf_f = &buf->csf_f;
const adm_dwt_band_t *csf_a = &buf->csf_a;
// for ADM: scales goes from 0 to 3 but in noise floor paper, it goes from
// 1 to 4 (from finest scale to coarsest scale).
// 0 is scale zero passed to dwt_quant_step
const float factor1 = dwt_quant_step(&dwt_7_9_YCbCr_threshold[0], 0, 1, adm_norm_view_dist, adm_ref_display_height);
const float factor2 = dwt_quant_step(&dwt_7_9_YCbCr_threshold[0], 0, 2, adm_norm_view_dist, adm_ref_display_height);
const float rfactor1[3] = { 1.0f / factor1, 1.0f / factor1, 1.0f / factor2 };
/**
* rfactor is converted to fixed-point for scale0 and stored in i_rfactor
* multiplied by 2^21 for rfactor[0,1] and by 2^23 for rfactor[2].
* For adm_norm_view_dist 3.0 and adm_ref_display_height 1080,
* i_rfactor is around { 36453,36453,49417 }
*/
uint16_t i_rfactor[3];
if (fabs(adm_norm_view_dist * adm_ref_display_height - DEFAULT_ADM_NORM_VIEW_DIST * DEFAULT_ADM_REF_DISPLAY_HEIGHT) < 1.0e-8) {
i_rfactor[0] = 36453;
i_rfactor[1] = 36453;
i_rfactor[2] = 49417;
}
else {
const double pow2_21 = pow(2, 21);
const double pow2_23 = pow(2, 23);
i_rfactor[0] = (uint16_t) (rfactor1[0] * pow2_21);
i_rfactor[1] = (uint16_t) (rfactor1[1] * pow2_21);
i_rfactor[2] = (uint16_t) (rfactor1[2] * pow2_23);
}
const int32_t shift_xhsq = 29;
const int32_t shift_xvsq = 29;
const int32_t shift_xdsq = 30;
const int32_t add_shift_xhsq = 268435456;
const int32_t add_shift_xvsq = 268435456;
const int32_t add_shift_xdsq = 536870912;
const uint32_t shift_xhcub = (uint32_t)ceil(log2(w) - 4);
const uint32_t add_shift_xhcub = (uint32_t)pow(2, (shift_xhcub - 1));
const uint32_t shift_xvcub = (uint32_t)ceil(log2(w) - 4);
const uint32_t add_shift_xvcub = (uint32_t)pow(2, (shift_xvcub - 1));
const uint32_t shift_xdcub = (uint32_t)ceil(log2(w) - 3);
const uint32_t add_shift_xdcub = (uint32_t)pow(2, (shift_xdcub - 1));
const uint32_t shift_inner_accum = (uint32_t)ceil(log2(h));
const uint32_t add_shift_inner_accum = (uint32_t)pow(2, (shift_inner_accum - 1));
const int32_t shift_xhsub = 10;
const int32_t shift_xvsub = 10;
const int32_t shift_xdsub = 12;
int16_t *angles[3] = { csf_a->band_h, csf_a->band_v, csf_a->band_d };
int16_t *flt_angles[3] = { csf_f->band_h, csf_f->band_v, csf_f->band_d };
/* The computation of the scales is not required for the regions which lie
* outside the frame borders
*/
int left = w * ADM_BORDER_FACTOR - 0.5;
int top = h * ADM_BORDER_FACTOR - 0.5;
int right = w - left;
int bottom = h - top;
const int start_col = (left > 1) ? left : 1;
const int end_col = (right < (w - 1)) ? right : (w - 1);
const int start_row = (top > 1) ? top : 1;
const int end_row = (bottom < (h - 1)) ? bottom : (h - 1);
int i, j;
int64_t val;
int32_t xh, xv, xd, thr;
int32_t xh_sq, xv_sq, xd_sq;
int64_t accum_h = 0, accum_v = 0, accum_d = 0;
int64_t accum_inner_h = 0, accum_inner_v = 0, accum_inner_d = 0;
/* i=0,j=0 */
if ((top <= 0) && (left <= 0))
{
xh = (int32_t)src->band_h[0] * i_rfactor[0];
xv = (int32_t)src->band_v[0] * i_rfactor[1];
xd = (int32_t)src->band_d[0] * i_rfactor[2];
ADM_CM_THRESH_S_0_0(angles, flt_angles, csf_a_stride, &thr, w, h, 0, 0);
//thr is shifted to make it's Q format equivalent to xh,xv,xd
/**
* max value of xh_sq and xv_sq is 1301381973 and that of xd_sq is 1195806729
*
* max(val before shift for h and v) is 9.995357299 * —10^17.
* 9.995357299 * —10^17 * 2^4 is close to 2^64.
* Hence shift is done based on width subtracting 4
*
* max(val before shift for d) is 1.355006643 * —10^18
* 1.355006643 * —10^18 * 2^3 is close to 2^64
* Hence shift is done based on width subtracting 3
*/
ADM_CM_ACCUM_ROUND(xh, thr, shift_xhsub, xh_sq, add_shift_xhsq, shift_xhsq, val,
add_shift_xhcub, shift_xhcub, accum_inner_h);
ADM_CM_ACCUM_ROUND(xv, thr, shift_xvsub, xv_sq, add_shift_xvsq, shift_xvsq, val,
add_shift_xvcub, shift_xvcub, accum_inner_v);
ADM_CM_ACCUM_ROUND(xd, thr, shift_xdsub, xd_sq, add_shift_xdsq, shift_xdsq, val,
add_shift_xdcub, shift_xdcub, accum_inner_d);
}
/* i=0, j */
if (top <= 0) {
for (j = start_col; j < end_col; ++j) {
xh = src->band_h[j] * i_rfactor[0];
xv = src->band_v[j] * i_rfactor[1];
xd = src->band_d[j] * i_rfactor[2];
ADM_CM_THRESH_S_0_J(angles, flt_angles, csf_a_stride, &thr, w, h, 0, j);
ADM_CM_ACCUM_ROUND(xh, thr, shift_xhsub, xh_sq, add_shift_xhsq, shift_xhsq, val,
add_shift_xhcub, shift_xhcub, accum_inner_h);
ADM_CM_ACCUM_ROUND(xv, thr, shift_xvsub, xv_sq, add_shift_xvsq, shift_xvsq, val,
add_shift_xvcub, shift_xvcub, accum_inner_v);
ADM_CM_ACCUM_ROUND(xd, thr, shift_xdsub, xd_sq, add_shift_xdsq, shift_xdsq, val,
add_shift_xdcub, shift_xdcub, accum_inner_d);
}
}
/* i=0,j=w-1 */
if ((top <= 0) && (right > (w - 1)))
{
xh = src->band_h[w - 1] * i_rfactor[0];
xv = src->band_v[w - 1] * i_rfactor[1];
xd = src->band_d[w - 1] * i_rfactor[2];
ADM_CM_THRESH_S_0_W_M_1(angles, flt_angles, csf_a_stride, &thr, w, h, 0, (w - 1));
ADM_CM_ACCUM_ROUND(xh, thr, shift_xhsub, xh_sq, add_shift_xhsq, shift_xhsq, val,
add_shift_xhcub, shift_xhcub, accum_inner_h);
ADM_CM_ACCUM_ROUND(xv, thr, shift_xvsub, xv_sq, add_shift_xvsq, shift_xvsq, val,
add_shift_xvcub, shift_xvcub, accum_inner_v);
ADM_CM_ACCUM_ROUND(xd, thr, shift_xdsub, xd_sq, add_shift_xdsq, shift_xdsq, val,
add_shift_xdcub, shift_xdcub, accum_inner_d);
}
//Shift is done based on height
accum_h += (accum_inner_h + add_shift_inner_accum) >> shift_inner_accum;
accum_v += (accum_inner_v + add_shift_inner_accum) >> shift_inner_accum;
accum_d += (accum_inner_d + add_shift_inner_accum) >> shift_inner_accum;
if ((left > 0) && (right <= (w - 1))) /* Completely within frame */
{
for (i = start_row; i < end_row; ++i) {
accum_inner_h = 0;
accum_inner_v = 0;
accum_inner_d = 0;
for (j = start_col; j < end_col; ++j) {
xh = src->band_h[i * src_stride + j] * i_rfactor[0];
xv = src->band_v[i * src_stride + j] * i_rfactor[1];
xd = src->band_d[i * src_stride + j] * i_rfactor[2];
ADM_CM_THRESH_S_I_J(angles, flt_angles, csf_a_stride, &thr, w, h, i, j);
ADM_CM_ACCUM_ROUND(xh, thr, shift_xhsub, xh_sq, add_shift_xhsq, shift_xhsq, val,
add_shift_xhcub, shift_xhcub, accum_inner_h);
ADM_CM_ACCUM_ROUND(xv, thr, shift_xvsub, xv_sq, add_shift_xvsq, shift_xvsq, val,
add_shift_xvcub, shift_xvcub, accum_inner_v);
ADM_CM_ACCUM_ROUND(xd, thr, shift_xdsub, xd_sq, add_shift_xdsq, shift_xdsq, val,
add_shift_xdcub, shift_xdcub, accum_inner_d);
}
accum_h += (accum_inner_h + add_shift_inner_accum) >> shift_inner_accum;
accum_v += (accum_inner_v + add_shift_inner_accum) >> shift_inner_accum;
accum_d += (accum_inner_d + add_shift_inner_accum) >> shift_inner_accum;
}
}
else if ((left <= 0) && (right <= (w - 1))) /* Right border within frame, left outside */
{
for (i = start_row; i < end_row; ++i) {
accum_inner_h = 0;
accum_inner_v = 0;
accum_inner_d = 0;
/* j = 0 */
xh = src->band_h[i * src_stride] * i_rfactor[0];
xv = src->band_v[i * src_stride] * i_rfactor[1];
xd = src->band_d[i * src_stride] * i_rfactor[2];
ADM_CM_THRESH_S_I_0(angles, flt_angles, csf_a_stride, &thr, w, h, i, 0);
ADM_CM_ACCUM_ROUND(xh, thr, shift_xhsub, xh_sq, add_shift_xhsq, shift_xhsq, val,
add_shift_xhcub, shift_xhcub, accum_inner_h);
ADM_CM_ACCUM_ROUND(xv, thr, shift_xvsub, xv_sq, add_shift_xvsq, shift_xvsq, val,
add_shift_xvcub, shift_xvcub, accum_inner_v);
ADM_CM_ACCUM_ROUND(xd, thr, shift_xdsub, xd_sq, add_shift_xdsq, shift_xdsq, val,
add_shift_xdcub, shift_xdcub, accum_inner_d);
/* j within frame */
for (j = start_col; j < end_col; ++j) {
xh = src->band_h[i * src_stride + j] * i_rfactor[0];
xv = src->band_v[i * src_stride + j] * i_rfactor[1];
xd = src->band_d[i * src_stride + j] * i_rfactor[2];
ADM_CM_THRESH_S_I_J(angles, flt_angles, csf_a_stride, &thr, w, h, i, j);
ADM_CM_ACCUM_ROUND(xh, thr, shift_xhsub, xh_sq, add_shift_xhsq, shift_xhsq, val,
add_shift_xhcub, shift_xhcub, accum_inner_h);
ADM_CM_ACCUM_ROUND(xv, thr, shift_xvsub, xv_sq, add_shift_xvsq, shift_xvsq, val,
add_shift_xvcub, shift_xvcub, accum_inner_v);
ADM_CM_ACCUM_ROUND(xd, thr, shift_xdsub, xd_sq, add_shift_xdsq, shift_xdsq, val,
add_shift_xdcub, shift_xdcub, accum_inner_d);
}
accum_h += (accum_inner_h + add_shift_inner_accum) >> shift_inner_accum;
accum_v += (accum_inner_v + add_shift_inner_accum) >> shift_inner_accum;
accum_d += (accum_inner_d + add_shift_inner_accum) >> shift_inner_accum;
}
}
else if ((left > 0) && (right > (w - 1))) /* Left border within frame, right outside */
{
for (i = start_row; i < end_row; ++i) {
accum_inner_h = 0;
accum_inner_v = 0;
accum_inner_d = 0;
/* j within frame */
for (j = start_col; j < end_col; ++j) {
xh = src->band_h[i * src_stride + j] * i_rfactor[0];
xv = src->band_v[i * src_stride + j] * i_rfactor[1];
xd = src->band_d[i * src_stride + j] * i_rfactor[2];
ADM_CM_THRESH_S_I_J(angles, flt_angles, csf_a_stride, &thr, w, h, i, j);
ADM_CM_ACCUM_ROUND(xh, thr, shift_xhsub, xh_sq, add_shift_xhsq, shift_xhsq, val,
add_shift_xhcub, shift_xhcub, accum_inner_h);
ADM_CM_ACCUM_ROUND(xv, thr, shift_xvsub, xv_sq, add_shift_xvsq, shift_xvsq, val,
add_shift_xvcub, shift_xvcub, accum_inner_v);
ADM_CM_ACCUM_ROUND(xd, thr, shift_xdsub, xd_sq, add_shift_xdsq, shift_xdsq, val,
add_shift_xdcub, shift_xdcub, accum_inner_d);
}
/* j = w-1 */
xh = src->band_h[i * src_stride + w - 1] * i_rfactor[0];
xv = src->band_v[i * src_stride + w - 1] * i_rfactor[1];
xd = src->band_d[i * src_stride + w - 1] * i_rfactor[2];
ADM_CM_THRESH_S_I_W_M_1(angles, flt_angles, csf_a_stride, &thr, w, h, i, (w - 1));
ADM_CM_ACCUM_ROUND(xh, thr, shift_xhsub, xh_sq, add_shift_xhsq, shift_xhsq, val,
add_shift_xhcub, shift_xhcub, accum_inner_h);
ADM_CM_ACCUM_ROUND(xv, thr, shift_xvsub, xv_sq, add_shift_xvsq, shift_xvsq, val,
add_shift_xvcub, shift_xvcub, accum_inner_v);
ADM_CM_ACCUM_ROUND(xd, thr, shift_xdsub, xd_sq, add_shift_xdsq, shift_xdsq, val,
add_shift_xdcub, shift_xdcub, accum_inner_d);
accum_h += (accum_inner_h + add_shift_inner_accum) >> shift_inner_accum;
accum_v += (accum_inner_v + add_shift_inner_accum) >> shift_inner_accum;
accum_d += (accum_inner_d + add_shift_inner_accum) >> shift_inner_accum;
}
}
else /* Both borders outside frame */
{
for (i = start_row; i < end_row; ++i) {
accum_inner_h = 0;
accum_inner_v = 0;
accum_inner_d = 0;
/* j = 0 */
xh = src->band_h[i * src_stride] * i_rfactor[0];
xv = src->band_v[i * src_stride] * i_rfactor[1];
xd = src->band_d[i * src_stride] * i_rfactor[2];
ADM_CM_THRESH_S_I_0(angles, flt_angles, csf_a_stride, &thr, w, h, i, 0);
ADM_CM_ACCUM_ROUND(xh, thr, shift_xhsub, xh_sq, add_shift_xhsq, shift_xhsq, val,
add_shift_xhcub, shift_xhcub, accum_inner_h);
ADM_CM_ACCUM_ROUND(xv, thr, shift_xvsub, xv_sq, add_shift_xvsq, shift_xvsq, val,
add_shift_xvcub, shift_xvcub, accum_inner_v);
ADM_CM_ACCUM_ROUND(xd, thr, shift_xdsub, xd_sq, add_shift_xdsq, shift_xdsq, val,
add_shift_xdcub, shift_xdcub, accum_inner_d);
/* j within frame */
for (j = start_col; j < end_col; ++j) {
xh = src->band_h[i * src_stride + j] * i_rfactor[0];
xv = src->band_v[i * src_stride + j] * i_rfactor[1];
xd = src->band_d[i * src_stride + j] * i_rfactor[2];
ADM_CM_THRESH_S_I_J(angles, flt_angles, csf_a_stride, &thr, w, h, i, j);
ADM_CM_ACCUM_ROUND(xh, thr, shift_xhsub, xh_sq, add_shift_xhsq, shift_xhsq, val,
add_shift_xhcub, shift_xhcub, accum_inner_h);
ADM_CM_ACCUM_ROUND(xv, thr, shift_xvsub, xv_sq, add_shift_xvsq, shift_xvsq, val,
add_shift_xvcub, shift_xvcub, accum_inner_v);
ADM_CM_ACCUM_ROUND(xd, thr, shift_xdsub, xd_sq, add_shift_xdsq, shift_xdsq, val,
add_shift_xdcub, shift_xdcub, accum_inner_d);
}
/* j = w-1 */
xh = src->band_h[i * src_stride + w - 1] * i_rfactor[0];
xv = src->band_v[i * src_stride + w - 1] * i_rfactor[1];
xd = src->band_d[i * src_stride + w - 1] * i_rfactor[2];
ADM_CM_THRESH_S_I_W_M_1(angles, flt_angles, csf_a_stride, &thr, w, h, i, (w - 1));
ADM_CM_ACCUM_ROUND(xh, thr, shift_xhsub, xh_sq, add_shift_xhsq, shift_xhsq, val,
add_shift_xhcub, shift_xhcub, accum_inner_h);
ADM_CM_ACCUM_ROUND(xv, thr, shift_xvsub, xv_sq, add_shift_xvsq, shift_xvsq, val,
add_shift_xvcub, shift_xvcub, accum_inner_v);
ADM_CM_ACCUM_ROUND(xd, thr, shift_xdsub, xd_sq, add_shift_xdsq, shift_xdsq, val,
add_shift_xdcub, shift_xdcub, accum_inner_d);
accum_h += (accum_inner_h + add_shift_inner_accum) >> shift_inner_accum;
accum_v += (accum_inner_v + add_shift_inner_accum) >> shift_inner_accum;
accum_d += (accum_inner_d + add_shift_inner_accum) >> shift_inner_accum;
}
}
accum_inner_h = 0;
accum_inner_v = 0;
accum_inner_d = 0;
/* i=h-1,j=0 */
if ((bottom > (h - 1)) && (left <= 0))
{
xh = src->band_h[(h - 1) * src_stride] * i_rfactor[0];
xv = src->band_v[(h - 1) * src_stride] * i_rfactor[1];
xd = src->band_d[(h - 1) * src_stride] * i_rfactor[2];
ADM_CM_THRESH_S_H_M_1_0(angles, flt_angles, csf_a_stride, &thr, w, h, (h - 1), 0);
ADM_CM_ACCUM_ROUND(xh, thr, shift_xhsub, xh_sq, add_shift_xhsq, shift_xhsq, val,
add_shift_xhcub, shift_xhcub, accum_inner_h);
ADM_CM_ACCUM_ROUND(xv, thr, shift_xvsub, xv_sq, add_shift_xvsq, shift_xvsq, val,
add_shift_xvcub, shift_xvcub, accum_inner_v);
ADM_CM_ACCUM_ROUND(xd, thr, shift_xdsub, xd_sq, add_shift_xdsq, shift_xdsq, val,
add_shift_xdcub, shift_xdcub, accum_inner_d);
}
/* i=h-1,j */
if (bottom > (h - 1)) {
for (j = start_col; j < end_col; ++j) {
xh = src->band_h[(h - 1) * src_stride + j] * i_rfactor[0];
xv = src->band_v[(h - 1) * src_stride + j] * i_rfactor[1];
xd = src->band_d[(h - 1) * src_stride + j] * i_rfactor[2];
ADM_CM_THRESH_S_H_M_1_J(angles, flt_angles, csf_a_stride, &thr, w, h, (h - 1), j);
ADM_CM_ACCUM_ROUND(xh, thr, shift_xhsub, xh_sq, add_shift_xhsq, shift_xhsq, val,
add_shift_xhcub, shift_xhcub, accum_inner_h);
ADM_CM_ACCUM_ROUND(xv, thr, shift_xvsub, xv_sq, add_shift_xvsq, shift_xvsq, val,
add_shift_xvcub, shift_xvcub, accum_inner_v);
ADM_CM_ACCUM_ROUND(xd, thr, shift_xdsub, xd_sq, add_shift_xdsq, shift_xdsq, val,
add_shift_xdcub, shift_xdcub, accum_inner_d);
}
}
/* i-h-1,j=w-1 */
if ((bottom > (h - 1)) && (right > (w - 1)))
{
xh = src->band_h[(h - 1) * src_stride + w - 1] * i_rfactor[0];
xv = src->band_v[(h - 1) * src_stride + w - 1] * i_rfactor[1];
xd = src->band_d[(h - 1) * src_stride + w - 1] * i_rfactor[2];
ADM_CM_THRESH_S_H_M_1_W_M_1(angles, flt_angles, csf_a_stride, &thr, w, h,
(h - 1), (w - 1));
ADM_CM_ACCUM_ROUND(xh, thr, shift_xhsub, xh_sq, add_shift_xhsq, shift_xhsq, val,
add_shift_xhcub, shift_xhcub, accum_inner_h);
ADM_CM_ACCUM_ROUND(xv, thr, shift_xvsub, xv_sq, add_shift_xvsq, shift_xvsq, val,
add_shift_xvcub, shift_xvcub, accum_inner_v);
ADM_CM_ACCUM_ROUND(xd, thr, shift_xdsub, xd_sq, add_shift_xdsq, shift_xdsq, val,
add_shift_xdcub, shift_xdcub, accum_inner_d);
}
accum_h += (accum_inner_h + add_shift_inner_accum) >> shift_inner_accum;
accum_v += (accum_inner_v + add_shift_inner_accum) >> shift_inner_accum;
accum_d += (accum_inner_d + add_shift_inner_accum) >> shift_inner_accum;
/**
* For h and v total shifts pending from last stage is 6 rfactor[0,1] has 21 shifts
* => after cubing (6+21)*3=81 after squaring shifted by 29
* hence pending is 52-shift's done based on width and height
*
* For d total shifts pending from last stage is 6 rfactor[2] has 23 shifts
* => after cubing (6+23)*3=87 after squaring shifted by 30
* hence pending is 57-shift's done based on width and height
*/
float f_accum_h = (float)(accum_h / pow(2, (52 - shift_xhcub - shift_inner_accum)));
float f_accum_v = (float)(accum_v / pow(2, (52 - shift_xvcub - shift_inner_accum)));
float f_accum_d = (float)(accum_d / pow(2, (57 - shift_xdcub - shift_inner_accum)));
float num_scale_h = powf(f_accum_h, 1.0f / 3.0f) + powf((bottom - top) *
(right - left) / 32.0f, 1.0f / 3.0f);
float num_scale_v = powf(f_accum_v, 1.0f / 3.0f) + powf((bottom - top) *
(right - left) / 32.0f, 1.0f / 3.0f);
float num_scale_d = powf(f_accum_d, 1.0f / 3.0f) + powf((bottom - top) *
(right - left) / 32.0f, 1.0f / 3.0f);
return (num_scale_h + num_scale_v + num_scale_d);
}