Redo the forward 4x4 dct

The new fdct lowers the round trip sum squared error for a
4x4 block ~0.12. or ~0.008/pixel. For reference, the old
matrix multiply version has average round trip error 1.46
for a 4x4 block.

Thanks to "derf" for his suggestions and references.

Change-Id: I5559d1e81d333b319404ab16b336b739f87afc79
This commit is contained in:
Yaowu Xu
2010-06-16 12:52:18 -07:00
parent a5906668a3
commit d0dd01b8ce
14 changed files with 118 additions and 562 deletions

View File

@@ -11,163 +11,54 @@
#include <math.h>
static const short dct_matrix2[4][4] =
{
{ 23170, 30274, 23170, 12540 },
{ 23170, 12540, -23170, -30274 },
{ 23170, -12540, -23170, 30274 },
{ 23170, -30274, 23170, -12540 }
};
static const short dct_matrix1[4][4] =
{
{ 23170, 23170, 23170, 23170 },
{ 30274, 12540, -12540, -30274 },
{ 23170, -23170, -23170, 23170 },
{ 12540, -30274, 30274, -12540 }
};
#define _1STSTAGESHIFT 14
#define _1STSTAGEROUNDING (1<<( _1STSTAGESHIFT-1))
#define _2NDSTAGESHIFT 16
#define _2NDSTAGEROUNDING (1<<( _2NDSTAGESHIFT-1))
// using matrix multiply
void vp8_short_fdct4x4_c(short *input, short *output, int pitch)
{
int i, j, k;
short temp[4][4];
int sumtemp;
pitch >>= 1;
for (i = 0; i < 4; i++)
{
for (j = 0; j < 4; j++)
{
sumtemp = 0;
for (k = 0; k < 4; k++)
{
sumtemp += input[i*pitch+k] * dct_matrix2[k][j];
}
temp[i][j] = (short)((sumtemp + _1STSTAGEROUNDING) >> _1STSTAGESHIFT);
}
}
for (i = 0; i < 4; i++)
{
for (j = 0; j < 4; j++)
{
sumtemp = 0;
for (k = 0; k < 4; k++)
{
sumtemp += dct_matrix1[i][ k] * temp[k][ j];
}
output[i*4+j] = (short)((sumtemp + _2NDSTAGEROUNDING) >> _2NDSTAGESHIFT);
}
}
}
void vp8_short_fdct8x4_c(short *input, short *output, int pitch)
{
vp8_short_fdct4x4_c(input, output, pitch);
vp8_short_fdct4x4_c(input + 4, output + 16, pitch);
}
static const signed short x_c1 = 60547;
static const signed short x_c2 = 46341;
static const signed short x_c3 = 25080;
void vp8_fast_fdct4x4_c(short *input, short *output, int pitch)
{
int i;
int a1, b1, c1, d1;
int a2, b2, c2, d2;
short *ip = input;
short *op = output;
int temp1, temp2;
for (i = 0; i < 4; i++)
{
a1 = (ip[0] + ip[3]) * 2;
b1 = (ip[1] + ip[2]) * 2;
c1 = (ip[1] - ip[2]) * 2;
d1 = (ip[0] - ip[3]) * 2;
a1 = ((ip[0] + ip[3])<<3);
b1 = ((ip[1] + ip[2])<<3);
c1 = ((ip[1] - ip[2])<<3);
d1 = ((ip[0] - ip[3])<<3);
temp1 = a1 + b1;
temp2 = a1 - b1;
op[0] = a1 + b1;
op[2] = a1 - b1;
op[0] = ((temp1 * x_c2) >> 16) + temp1;
op[2] = ((temp2 * x_c2) >> 16) + temp2;
temp1 = (c1 * x_c3) >> 16;
temp2 = ((d1 * x_c1) >> 16) + d1;
op[1] = temp1 + temp2;
temp1 = (d1 * x_c3) >> 16;
temp2 = ((c1 * x_c1) >> 16) + c1;
op[3] = temp1 - temp2;
op[1] = (c1 * 2217 + d1 * 5352 + 14500)>>12;
op[3] = (d1 * 2217 - c1 * 5352 + 7500)>>12;
ip += pitch / 2;
op += 4;
}
}
ip = output;
op = output;
for (i = 0; i < 4; i++)
{
a1 = ip[0] + ip[12];
b1 = ip[4] + ip[8];
c1 = ip[4] - ip[8];
d1 = ip[0] - ip[12];
op[0] = ( a1 + b1 + 7)>>4;
op[8] = ( a1 - b1 + 7)>>4;
temp1 = a1 + b1;
temp2 = a1 - b1;
a2 = ((temp1 * x_c2) >> 16) + temp1;
c2 = ((temp2 * x_c2) >> 16) + temp2;
temp1 = (c1 * x_c3) >> 16;
temp2 = ((d1 * x_c1) >> 16) + d1;
b2 = temp1 + temp2;
temp1 = (d1 * x_c3) >> 16;
temp2 = ((c1 * x_c1) >> 16) + c1;
d2 = temp1 - temp2;
op[0] = (a2 + 1) >> 1;
op[4] = (b2 + 1) >> 1;
op[8] = (c2 + 1) >> 1;
op[12] = (d2 + 1) >> 1;
op[4] =((c1 * 2217 + d1 * 5352 + 12000)>>16) + (d1!=0);
op[12] = (d1 * 2217 - c1 * 5352 + 51000)>>16;
ip++;
op++;
}
}
void vp8_fast_fdct8x4_c(short *input, short *output, int pitch)
void vp8_short_fdct8x4_c(short *input, short *output, int pitch)
{
vp8_fast_fdct4x4_c(input, output, pitch);
vp8_fast_fdct4x4_c(input + 4, output + 16, pitch);
vp8_short_fdct4x4_c(input, output, pitch);
vp8_short_fdct4x4_c(input + 4, output + 16, pitch);
}
void vp8_short_walsh4x4_c(short *input, short *output, int pitch)