x264的哈达玛变换（hadamard）的实现和优化

x264的哈达玛变换（hadamard）的实现和优化hadamard 变换理论
很多⽹页都有介绍，我就不拷贝了，给两个链接。

下⾯的是harvey mudd college 的⼀个“计算机图像处理分析”课件中哈达玛变换的⼀个章节
（JASON GARRETT-GLASER x264的主开发就是在这个学校上过学阿。

很棒的⼀个⼯程类⼤
学）
/e161/lectures/wht/index.html
⼤概提纲：
1.介绍了hadamard 矩阵的定义
2.快速hadamard变换算法（hadamard orderd）
3. Sequency Ordered hadamard 矩阵的定义 --h264使⽤了这个定义
4. 快速hadamard变换算法（Sequency orderd）
/xkfz007/articles/2616143.html
这个是中⽂的⼀个博客“X264中SATD实现分析”
优点是-中⽂，哈哈。

步骤也⽐较详细，缺点是关于h264使⽤的sequency ordered那部分交代的
不清楚。

JM86中hadamard的实现
函数 SATD ⾥
if (use_hadamard) //++ 如果采⽤了Hadamard变换，则先对残差块进⾏Hadamard变换，然后将变换后的16个残差值取绝对值相加作为代价
{
/*===== hadamard transform =====*/
m[ 0] = d[ 0] + d[12];
m[ 4] = d[ 4] + d[ 8];
m[ 8] = d[ 4] - d[ 8];
m[12] = d[ 0] - d[12];
m[ 1] = d[ 1] + d[13];
m[ 5] = d[ 5] + d[ 9];
m[ 9] = d[ 5] - d[ 9];
m[13] = d[ 1] - d[13];
m[ 2] = d[ 2] + d[14];
m[ 6] = d[ 6] + d[10];
m[10] = d[ 6] - d[10];
m[14] = d[ 2] - d[14];
m[ 3] = d[ 3] + d[15];
m[ 7] = d[ 7] + d[11];
m[11] = d[ 7] - d[11];
m[15] = d[ 3] - d[15];
d[ 0] = m[ 0] + m[ 4];
d[ 8] = m[ 0] - m[ 4];
d[ 4] = m[ 8] + m[12];
d[12] = m[12] - m[ 8];
d[ 1] = m[ 1] + m[ 5];
d[ 9] = m[ 1] - m[ 5];
d[ 5] = m[ 9] + m[13];
d[13] = m[13] - m[ 9];
d[ 2] = m[ 2] + m[ 6];
d[10] = m[ 2] - m[ 6];
d[ 6] = m[10] + m[14];
d[14] = m[14] - m[10];
d[ 3] = m[ 3] + m[ 7];
d[11] = m[ 3] - m[ 7];
d[ 7] = m[11] + m[15];
d[15] = m[15] - m[11];
m[ 0] = d[ 0] + d[ 3];
m[ 1] = d[ 1] + d[ 2];
m[ 2] = d[ 1] - d[ 2];
m[ 3] = d[ 0] - d[ 3];
m[ 4] = d[ 4] + d[ 7];
m[ 5] = d[ 5] + d[ 6];
m[ 6] = d[ 5] - d[ 6];
m[ 7] = d[ 4] - d[ 7];
m[ 8] = d[ 8] + d[11];
m[ 9] = d[ 9] + d[10];
m[10] = d[ 9] - d[10];
m[11] = d[ 8] - d[11];
m[12] = d[12] + d[15];
m[13] = d[13] + d[14];
m[14] = d[13] - d[14];
m[15] = d[12] - d[15];
d[ 0] = m[ 0] + m[ 1];
d[ 1] = m[ 0] - m[ 1];
d[ 2] = m[ 2] + m[ 3];
d[ 3] = m[ 3] - m[ 2];
d[ 4] = m[ 4] + m[ 5];
d[ 5] = m[ 4] - m[ 5];
d[ 6] = m[ 6] + m[ 7];
d[ 7] = m[ 7] - m[ 6];
d[ 8] = m[ 8] + m[ 9];
d[ 9] = m[ 8] - m[ 9];
d[10] = m[10] + m[11];
d[11] = m[11] - m[10];
d[12] = m[12] + m[13];
d[13] = m[12] - m[13];
d[14] = m[14] + m[15];
d[15] = m[15] - m[14];
/*===== sum up =====*/
for (dd=diff[k=0]; k<16; dd=diff[++k])
{
satd += (dd < 0 ? -dd : dd);
}
satd >>= 1;
}
这个算法过程很清楚，参见上⾯的 "快速hadamard变换算法（hadamard orderd）"
可能有点搞不清楚的地⽅就是，⼤部分资料讲的是⼀维的哈达吗的变换，⼆维的都只是交代⼀下先变换⼆维图像⾏，再在变换后的基础上在变换列。

我当时糊涂了半天（要吃核桃了），那现在看到这个实现，应该很清楚了，实际上就是算两次，
第⼀次以把4x4横向划分为4个1x4⽮量，做⼀维变换，总共做四次，这样以⾏做变换就完成了，
第⼆次，在得到的4x4矩阵上再纵向划分为4个1x4⽮量，做⼀维变换。

总共做四次。

这样就做完了⼆维的变换了。

x264未优化c语⾔版本的hadamard变换实现
在版本 commit 5dc0aae2f900064d1f58579929a2285ab289a436 ，也就是最开始的版本
在函数 pixel_satd_wxh 中
for( d = 0; d < 4; d++ )
{
int s01, s23;
int d01, d23;
s01 = diff[d][0] + diff[d][1]; s23 = diff[d][2] + diff[d][3];
d01 = diff[d][0] - diff[d][1]; d23 = diff[d][2] - diff[d][3];
tmp[d][0] = s01 + s23;
tmp[d][1] = s01 - s23;
tmp[d][2] = d01 - d23;
tmp[d][3] = d01 + d23;
}
for( d = 0; d < 4; d++ )
{
int s01, s23;
int d01, d23;
s01 = tmp[0][d] + tmp[1][d]; s23 = tmp[2][d] + tmp[3][d];
d01 = tmp[0][d] - tmp[1][d]; d23 = tmp[2][d] - tmp[3][d];
i_satd += abs( s01 + s23 ) + abs( s01 - s23 ) + abs( d01 - d23 ) + abs( d01 + d23 );
}
可以看到算法有些和JM86的不⼀样，因为他⽤了哈达吗变换的Sequency orderd。

说实话，在课件⾥的”快速hadamard变换算法（Sequency orderd）“ 说明，我没看很懂阿。

有⾼⼿能讲解下吗？但是这段代码⼤致不难理解，也是典型的碟形算法。

只是次序和JM86不⼀样估计是因为采⽤了Sequency ordered的原因。

之所以采⽤这种ordered，我理解是因为Sequency ordered 是按照hadamard变换矩阵的符号变换次序重拍列了。

变化最少的在嘴上⾯，最多的在最下⾯⼀列。

这样应该能量更加集中在了左上⾓，更⽅便后⾯的压缩。

x264优化c语⾔版本的hadamard变换实现参见patch “
1.6x faster satd_c (and sa8d and hadamard_ac) with pseudo-simd.”
static NOINLINE int x264_pixel_satd_4x4( pixel *pix1, intptr_t i_pix1, pixel *pix2, intptr_t i_pix2 )
{
sum2_t tmp[4][2];
sum2_t a0, a1, a2, a3, b0, b1;
sum2_t sum = 0;
for( int i = 0; i < 4; i++, pix1 += i_pix1, pix2 += i_pix2 )
{
a0 = pix1[0] - pix2[0];
a1 = pix1[1] - pix2[1];
b0 = (a0+a1) + ((a0-a1)<<BITS_PER_SUM);
a2 = pix1[2] - pix2[2];
a3 = pix1[3] - pix2[3];
b1 = (a2+a3) + ((a2-a3)<<BITS_PER_SUM);
tmp[i][0] = b0 + b1;
tmp[i][1] = b0 - b1;
}
for( int i = 0; i < 2; i++ )
{
HADAMARD4( a0, a1, a2, a3, tmp[0][i], tmp[1][i], tmp[2][i], tmp[3][i] );
a0 = abs2(a0) + abs2(a1) + abs2(a2) + abs2(a3);
sum += ((sum_t)a0) + (a0>>BITS_PER_SUM);
}
return sum >> 1;
}
仔细观察可以发现，快速的哈达吗变换算法是对称的。

所以这个优化的主要思想是第⼀趟⾏变换的时候把数据重排列，⾼位和低位都放置好第⼀次变换后的数据，
(类似，tmp[i][0] = x|x ; tmp[i][1] = x|x 这样的布局) 这样在后⾯的列变换的时候就可以并⾏了处理两个数据
正常的4x4哈达吗快速变换的算法的复杂度是，8（⼀维哈达吗变换的加法次数）*4（⾏） + 8 * 4（列变换）=64
⽽新的算法是 4*8 + 2*8 = 48次。

加速⽐1.3倍。