专栏首页Visual Codex[AV1] Constrained Directional Enhancement Filter

[AV1] Constrained Directional Enhancement Filter

CDEF原理

在AV1视频标准中,CDEF为非线性空间滤波器,该滤波器以8x8为基本单位,通过沿着物体的方向进行滤波,从而消除和减弱振铃(ringing artifact)效应,从而提升重建图像的质量。

在Encoder端,CDEF滤波一共分三部分组成,分别为

  • CDEF direction 处理
  • CDEF Block 处理
  • CDEF Filter 处理

CDEF direction 处理

首先,CDEF一共定义了八个方向,如下图所示:

沿着方向,一个8x8的块由很多行(由k表示)的像素组成,每个方向的块的行数不一定相等,同一个块的每一行的像素个数也不一定相等,其值分别为

Direction

The number of lines (k)

d = 0

15

d = 1

11

d = 2

8

d = 3

15

d = 4

11

d = 5

11

d = 6

8

d = 7

11

d = 8

11

(上述数据在代码中以partial_sum_hv,partial_sum_diag,partial_sum_alt的数组形式出现)

CDEF direction 处理主要由以下几个步骤组成

第一步:计算每个方向上像素与平均像素的差值。在AV1协议中,若像素的比特深度超过8比特,则要将像素比特先缩放到8比特,然后再减去128(这样做的目的据说是为了使计算的整数只需要一个32位的有符号整型就可以了),在代码里的 cdef_find_dir_c 函数里可以找到如下代码

const int px = (img[x] >> bitdepth_min_8) - 128;

对于每一个方向的每一行,求其像素值与行像素均值的差的平方和。

下面的公式是每一行的均值。(d代表八个方向其中之一,k代表这个方向中的第k行,P代表k行的所有像素集合,x是k行的一个像素,N是k行像素的个数)

则其像素值与行像素均值的差的平方和SSE即为

把上面两个式子结合化简后可以得到

第二步就是为了求上式的极值,上式中,其中k行的x的平方和是固定值,所以想得到E的最小值,只需要求右边的第二项的最大值即可,即

第三步就是取其右边第二项的最大值(也为Cost)对应的方向,即我们需要求的方向。所以在计算的时候,Cost最大的direction即为最佳方向

    int best_dir = 0;
    unsigned best_cost = cost[0];
    for (int n = 1; n < 8; n++) {
        if (cost[n] > best_cost) {
            best_cost = cost[n];
            best_dir = n;
        }
    }

这里面有一个问题,就是每个方向的每行像素个数其实是不一样的,为了使其归一化,计算的时候会对每个方向加上权重。

CDEF direction 处理

至此为止,最佳的方向已经选定,但是这个方向并不需要传输到decoder去,因为这个计算过程与original frame并没有关系,在decoder可以依靠同样的方式计算出来最佳的方向。

下面就是说确定方向后,怎么滤波的问题了。

CDEF 整体流程图

核心代码分析

函数cdef_find_dir_c

这个函数是在decoder端,返回值是一个int型的变量:best_dir,其为当前8x8块的方向,之后调用滤波函数时,filter tap的选择要靠这个参数来抉择。

static int cdef_find_dir_c(const pixel *img, const ptrdiff_t stride, unsigned *const var HIGHBD_DECL_SUFFIX)
{
    const int bitdepth_min_8 = bitdepth_from_max(bitdepth_max) - 8;
    int partial_sum_hv[2][8] = { { 0 } };
    int partial_sum_diag[2][15] = { { 0 } };
    int partial_sum_alt[4][11] = { { 0 } };

    for (int y = 0; y < 8; y++) {
        for (int x = 0; x < 8; x++) {
            const int px = (img[x] >> bitdepth_min_8) - 128;

            partial_sum_diag[0][     y       +  x      ] += px;
            partial_sum_alt [0][     y       + (x >> 1)] += px;
            partial_sum_hv  [0][     y                 ] += px;
            partial_sum_alt [1][3 +  y       - (x >> 1)] += px;
            partial_sum_diag[1][7 +  y       -  x      ] += px;
            partial_sum_alt [2][3 - (y >> 1) +  x      ] += px;
            partial_sum_hv  [1][                x      ] += px;
            partial_sum_alt [3][    (y >> 1) +  x      ] += px;
        }
        img += PXSTRIDE(stride);
    }

    unsigned cost[8] = { 0 };
    for (int n = 0; n < 8; n++) {
        cost[2] += partial_sum_hv[0][n] * partial_sum_hv[0][n];
        cost[6] += partial_sum_hv[1][n] * partial_sum_hv[1][n];
    }
    cost[2] *= 105; // 这个105即是上面理论部分说的权重 weight
    cost[6] *= 105;

    static const uint16_t div_table[7] = { 840, 420, 280, 210, 168, 140, 120 };
    for (int n = 0; n < 7; n++) {
        const int d = div_table[n];
        cost[0] += (partial_sum_diag[0][n]      * partial_sum_diag[0][n] +
                    partial_sum_diag[0][14 - n] * partial_sum_diag[0][14 - n]) * d;
        cost[4] += (partial_sum_diag[1][n]      * partial_sum_diag[1][n] +
                    partial_sum_diag[1][14 - n] * partial_sum_diag[1][14 - n]) * d;
    }
    cost[0] += partial_sum_diag[0][7] * partial_sum_diag[0][7] * 105;
    cost[4] += partial_sum_diag[1][7] * partial_sum_diag[1][7] * 105;

    for (int n = 0; n < 4; n++) {
        unsigned *const cost_ptr = &cost[n * 2 + 1];
        for (int m = 0; m < 5; m++)
            *cost_ptr += partial_sum_alt[n][3 + m] * partial_sum_alt[n][3 + m];
        *cost_ptr *= 105;
        for (int m = 0; m < 3; m++) {
            const int d = div_table[2 * m + 1];
            *cost_ptr += (partial_sum_alt[n][m]      * partial_sum_alt[n][m] +
                          partial_sum_alt[n][10 - m] * partial_sum_alt[n][10 - m]) * d;
        }
    }

    int best_dir = 0;
    unsigned best_cost = cost[0];
    for (int n = 1; n < 8; n++) {
        if (cost[n] > best_cost) {
            best_cost = cost[n];
            best_dir = n;
        }
    }

    *var = (best_cost - (cost[best_dir ^ 4])) >> 10;
    return best_dir;
}

函数cdef_filter_block_c

下面这个函数:cdef_filter_block_c 是调用cdef后filter每一个块的函数(包括4x4, 4x8, 8x8)

static NOINLINE void cdef_filter_block_c(pixel *dst, const ptrdiff_t dst_stride, const pixel (*left)[2], const pixel *const top,
                    const int pri_strength, const int sec_strength,
                    const int dir, const int damping, const int w, int h, const enum CdefEdgeFlags edges HIGHBD_DECL_SUFFIX)
// 参数中
// 1. dst 为滤波结果,完成CDEF后的结果将存入到此指针指向的内存中
// 2. dst_stride 为目标内存的 stride
// 3. left 为左侧的像素值(一列)
// 4. top 为上方的像素值(一排)
// 5. pri_strength 是第一级的滤波强度
// 6. sec_strength 是第二级的滤波强度
// 7. dir 由之前的 find_dir 函数计算出的滤波方向
// 8. damping 指滤波的threshold,详见 Spec 7.15.3的计算步骤
// 9. w, h 为滤波块的长和宽,一般是8x8
// 10. edges 指边缘情况,ABOVE, LEFT, RIGHT, BOTTOM 四个方向是否存在,由相应的宏定义
{
    // 1. 为滤波当前块做准备
    const ptrdiff_t tmp_stride = 12;
    assert((w == 4 || w == 8) && (h == 4 || h == 8));
    int16_t tmp_buf[144]; // 12*12 is the maximum value of tmp_stride * (h + 4)
    int16_t *tmp = tmp_buf + 2 * tmp_stride + 2;
	padding(tmp, tmp_stride, dst, dst_stride, left, top, w, h, edges);
	
    // 2. 滤波开始
    // 2.1 进行一级滤波
	if (pri_strength) {
    	const int bitdepth_min_8 = bitdepth_from_max(bitdepth_max) - 8;
    	const int pri_tap = 4 - ((pri_strength >> bitdepth_min_8) & 1);
    	const int pri_shift = imax(0, damping - ulog2(pri_strength));
		// 2.1.1 进行一级滤波的同时也进行二级滤波
        if (sec_strength) {
        	const int sec_shift = imax(0, damping - ulog2(sec_strength));
        	do {
            	for (int x = 0; x < w; x++) {
                	const int px = dst[x];
                	int sum = 0;
                	int max = px, min = px;
                	int pri_tap_k = pri_tap;
                	for (int k = 0; k < 2; k++) {
                    	const int off1 = dav1d_cdef_directions[dir + 2][k]; // 这个off1是中心元素离primary tap的滤波像素取值点的距离。
	                    const int p0 = tmp[x + off1]; // 获取到primary tap的上方像素位置
    	                const int p1 = tmp[x - off1]; // 获取到primary tap的下方像素位置
        	            sum += pri_tap_k * constrain(p0 - px, pri_strength, pri_shift); // 滤波计算(上)
            	        sum += pri_tap_k * constrain(p1 - px, pri_strength, pri_shift); // 滤波计算(下)
                	    // if pri_tap_k == 4 then it becomes 2 else it remains 3
	                    pri_tap_k = (pri_tap_k & 3) | 2;
    	                min = umin(p0, min);
        	            max = imax(p0, max);
            	        min = umin(p1, min);
                	    max = imax(p1, max);
                    	const int off2 = dav1d_cdef_directions[dir + 4][k]; // secondary tap的位置与中心像素点的位置差(方向1)
	                    const int off3 = dav1d_cdef_directions[dir + 0][k]; // secondary tap的位置与中心像素点的位置差(方向2)
    	                const int s0 = tmp[x + off2]; // 获取到secondary tap的像素位置(方向1上)
        	            const int s1 = tmp[x - off2]; // 获取到secondary tap的像素位置(方向1下)
            	        const int s2 = tmp[x + off3]; // 获取到secondary tap的像素位置(方向2上)
                	    const int s3 = tmp[x - off3]; // 获取到secondary tap的像素位置(方向2下)
	                    // sec_tap starts at 2 and becomes 1
    	                const int sec_tap = 2 - k;
        	            sum += sec_tap * constrain(s0 - px, sec_strength, sec_shift); // 滤波计算
            	        sum += sec_tap * constrain(s1 - px, sec_strength, sec_shift); // 滤波计算
                	    sum += sec_tap * constrain(s2 - px, sec_strength, sec_shift); // 滤波计算
                    	sum += sec_tap * constrain(s3 - px, sec_strength, sec_shift); // 滤波计算
                        min = umin(s0, min);
    	                max = imax(s0, max);
        	            min = umin(s1, min);
            	        max = imax(s1, max);
                	    min = umin(s2, min);
                    	max = imax(s2, max);
	                    min = umin(s3, min);
    	                max = imax(s3, max);
        	        }
                    // 最后的结果必须不超过所有参与运算的像素值的最大值,不小于所有参与运算的像素值的最小值,每一步的min,max的计算就是在比较像素值然后                     // 最后得出所有参与运算的像素的值的上下界
            	    dst[x] = iclip(px + ((sum - (sum < 0) + 8) >> 4), min, max);
            	}
            	dst += PXSTRIDE(dst_stride);
            	tmp += tmp_stride;
	        } while (--h);
        // 2.1.2 仅进行一级滤波
    	} else { 
        	do {
            	for (int x = 0; x < w; x++) {
                	const int px = dst[x];
                	int sum = 0;
	                int pri_tap_k = pri_tap;
    	            for (int k = 0; k < 2; k++) {
        	            const int off = dav1d_cdef_directions[dir + 2][k]; // dir
            	        const int p0 = tmp[x + off];
                	    const int p1 = tmp[x - off];
                    	sum += pri_tap_k * constrain(p0 - px, pri_strength, pri_shift);
	                    sum += pri_tap_k * constrain(p1 - px, pri_strength, pri_shift);
    	                pri_tap_k = (pri_tap_k & 3) | 2;
        	        }
            	    dst[x] = px + ((sum - (sum < 0) + 8) >> 4);
	            }	
    	        dst += PXSTRIDE(dst_stride);
        	    tmp += tmp_stride;
	        } while (--h);
	    }
    // 2.2 仅进行二级滤波
	} else {
    	assert(sec_strength);
	    const int sec_shift = imax(0, damping - ulog2(sec_strength));
    	do {
        	for (int x = 0; x < w; x++) {
            	const int px = dst[x];
	            int sum = 0;
    	        for (int k = 0; k < 2; k++) {
        	        const int off1 = dav1d_cdef_directions[dir + 4][k]; // dir + 2
            	    const int off2 = dav1d_cdef_directions[dir + 0][k]; // dir - 2
                	const int s0 = tmp[x + off1];
	                const int s1 = tmp[x - off1];
    	            const int s2 = tmp[x + off2];
        	        const int s3 = tmp[x - off2];
            	    const int sec_tap = 2 - k;
	                sum += sec_tap * constrain(s0 - px, sec_strength, sec_shift);
    	            sum += sec_tap * constrain(s1 - px, sec_strength, sec_shift);
        	        sum += sec_tap * constrain(s2 - px, sec_strength, sec_shift);
            	    sum += sec_tap * constrain(s3 - px, sec_strength, sec_shift);
	            }
    	        dst[x] = px + ((sum - (sum < 0) + 8) >> 4);
        	}
	        dst += PXSTRIDE(dst_stride);
    	    tmp += tmp_stride;
	    } while (--h);
	}
}

为滤波当前块做准备

以8x8的块为例,这个做准备的过程主要是构建一个12x12的区域。

为什么是12x12的大小呢?

因为块的大小为8x8,依据下图,针对于每一个pixel的滤波都需要一个5x5的块(其中中部的深灰色块为当前块),那么在滤波左边缘(left edge)右边缘(right edge)上边缘(above edge)和下边缘(bottom edge)的时候,5x5的这个mask会有一部分延申到8x8块的外面,而延申的距离最大为2,所以8x8的块,上下左右均延申2的话,即得到12x12的块。

12x12的块滤波的时候的情况如下

padding过程

下面是padding过程的代码,主要实现的功能是把上面8x8的块的外面、12x12块的里面的空白区域进行

static void padding(int16_t *tmp, const ptrdiff_t tmp_stride,
                    const pixel *src, const ptrdiff_t src_stride,
                    const pixel (*left)[2], const pixel *top,
                    const int w, const int h,
                    const enum CdefEdgeFlags edges)
{
    // 1. 根据当前块的左右上下是否是8x8块的边界来决定是否要padding,然后移动x_start和y_start
    // fill extended input buffer
    int x_start = -2, x_end = w + 2, y_start = -2, y_end = h + 2;
    if (!(edges & CDEF_HAVE_TOP)) {
        fill(tmp - 2 - 2 * tmp_stride, tmp_stride, w + 4, 2);
        y_start = 0;
    }
    if (!(edges & CDEF_HAVE_BOTTOM)) {
        fill(tmp + h * tmp_stride - 2, tmp_stride, w + 4, 2);
        y_end -= 2;
    }
    if (!(edges & CDEF_HAVE_LEFT)) {
        fill(tmp + y_start * tmp_stride - 2, tmp_stride, 2, y_end - y_start);
        x_start = 0;
    }
    if (!(edges & CDEF_HAVE_RIGHT)) {
        fill(tmp + y_start * tmp_stride + w, tmp_stride, 2, y_end - y_start);
        x_end -= 2;
    }
    
    // 2. 对于没用padding的部分,拷贝原本的像素值到12x12块的对应区域
	for (int y = y_start; y < 0; y++) {
    	for (int x = x_start; x < x_end; x++)
        	tmp[x + y * tmp_stride] = top[x];
	    top += PXSTRIDE(src_stride);
	}
	for (int y = 0; y < h; y++)
    	for (int x = x_start; x < 0; x++)
        	tmp[x + y * tmp_stride] = left[y][2 + x];
	for (int y = 0; y < y_end; y++) {
    	for (int x = (y < h) ? 0 : x_start; x < x_end; x++)
        	tmp[x] = src[x];
	    src += PXSTRIDE(src_stride);
    	tmp += tmp_stride;
	}
}

其中有一个fill函数,即把需要padding的部分填充INI16_MIN值。

static inline void fill(int16_t *tmp, const ptrdiff_t stride, const int w, const int h)
{
   /* Use a value that's a large positive number when interpreted as unsigned, and a large negative number when interpreted as signed. */
	for (int y = 0; y < h; y++) {
		for (int x = 0; x < w; x++)
		    tmp[x] = INT16_MIN;
			tmp += stride;
    }
}

实际滤波过程

这里重复放了一张之前放过的图,这张图表示的是一级滤波与二级滤波的滤波mask的形状与滤波时所取内部元素的分布图。

从图中可以知道,每一级的滤波都有8个tap,即8个方向,每一个8x8的块的方向均属于其中的一种,且primary tap和secondary tap是对应的,也就是说primary tap是d1,那secondary也是d1。而且primary tap与对应的secondary tap也是有关联的。

在primary tap中,把有值的方块连成的线的方向确定为当前块的方向的话。那么这个方向逆时针旋转45°即可得到secondary tap中的其中一个方向。

然后对这个方向取垂直中线方向即可得到secondary tap的第二个方向。

滤波计算

滤波计算主要是在函数 constrain 函数里进行。

static inline int apply_sign(const int v, const int s) 
{
    return s < 0 ? -v : v;
}

static inline int constrain(const int diff, const int threshold, const int shift)
    // 1. diff 是当前像素与滤波mask的像素之间的差,可能是正数,可能是负数
    // 2. threshold 是参数 damping
    // 3. shift 是像素差需要右移的位数
{
	const int adiff = abs(diff);
    return apply_sign(imin(adiff, imax(0, threshold - (adiff >> shift))), diff);
}

在constrain这个函数里,

本文参与腾讯云自媒体分享计划,欢迎正在阅读的你也加入,一起分享。

我来说两句

0 条评论
登录 后参与评论

相关文章

  • [AV1] AV1专栏介紹(正在更新中)

    这个专栏是我花了半年的时间阅读libav1,SVT-AV1以及dav1d的源码后摸索总结出来的AV1技能树,希望能帮助到你的AV1编解码器的学习。

    轻舞飞扬SR
  • ECCV 2020 的对抗相关论文(对抗生成、对抗攻击)

    本文汇总了ECCV 2020上部分对抗相关论文,后续公众号会随缘对一些paper做解读。感兴趣的同学,可先自行根据标题,搜索对应链接(有些paper可能未公布)...

    公众号机器学习与生成对抗网络
  • 视频编解码器 2020-比赛开始

    原文链接:https://blog.beamr.com/2020/05/28/video-codecs-in-2020-the-race-is-on/

    LiveVideoStack
  • 【ACM MM论文集】国际多媒体顶级会议ACM Multimedia 2017 Open Access Repository

    点击上方“专知”关注获取更多AI知识! 【导读】第25届ACM国际多媒体会议(ACM International Conference on Multimedi...

    WZEARW
  • 重磅纯干货 | 超级赞的语音识别/语音合成经典论文的路线图(1982-2018.5)

    网址:https://github.com/zzw922cn/awesome-speech-recognition-speech-synthesis-paper...

    用户7623498
  • 比较全面的3D数据处理建模等链接收集

    点云PCL博主
  • [计算机视觉论文速递] 2018-03-03

    通知:这篇推文很长,有32篇论文速递信息,涉及目标检测、图像分割、网络优化、人脸表情识别、SLAM和OCR等方向。 [1]《The 2018 DAVIS Cha...

    Amusi
  • 最全OCR相关资料整理

    最近看到一个非常赞的OCR相关资源,收集从2015.10.9到现在的一些OCR文献,github项目和博客资源等

    AI算法与图像处理
  • HDR关键技术:色调映射(三)

    HDR技术近年来发展迅猛,在未来将会成为图像与视频领域的主流。如何让HDR图像与视频也能够同时兼容地在现有的SDR显示设备上显示,是非常重要的技术挑战。色调映射...

    用户1324186

扫码关注云+社区

领取腾讯云代金券