DCTTest.cpp

#include <stdio.h> #include <math.h> #include "DCTTest.h" #define PI 3.14159265359 void runTests() {    printf("DCT Test\n\nAAN constants:\n");    printf("AANDCT4 %f\nAANDCT6 %f\nAANDCT2PLUS6 %f\nAANDCT2MINUS6NEGATED %f\nAANDCTSCALE0 %f\nAANDCTSCALE1 %f\nAANDCTSCALE2 %f\nAANDCTSCALE3 %f\nAANDCTSCALE4 %f\nAANDCTSCALE5 %f\nAANDCTSCALE6 %f\nAANDCTSCALE7 %f\nAANIDCTSCALE0 %f\nAANIDCTSCALE1 %f\nAANIDCTSCALE2 %f\nAANIDCTSCALE3 %f\nAANIDCTSCALE4 %f\nAANIDCTSCALE5 %f\nAANIDCTSCALE6 %f\nAANIDCTSCALE7 %f\n\n",          cos(4.0 * PI / 16.0),          cos(6.0 * PI / 16.0),          cos(2.0 * PI / 16.0) + cos(6.0 * PI / 16.0),          - cos(2.0 * PI / 16.0) + cos(6.0 * PI / 16.0),          1.0 / (2 * sqrt(2) * cos(0.0)),          1.0 / (cos(1.0 * PI / 16.0) * 4),          1.0 / (cos(2.0 * PI / 16.0) * 4),          1.0 / (cos(3.0 * PI / 16.0) * 4),          1.0 / (cos(4.0 * PI / 16.0) * 4),          1.0 / (cos(5.0 * PI / 16.0) * 4),          1.0 / (cos(6.0 * PI / 16.0) * 4),          1.0 / (cos(7.0 * PI / 16.0) * 4),          4.0 * sqrt(2),          4.0 / (cos(1.0 * PI / 16)),          4.0 / (cos(2.0 * PI / 16)),          4.0 / (cos(3.0 * PI / 16)),          4.0 / (cos(4.0 * PI / 16)),          4.0 / (cos(5.0 * PI / 16)),          4.0 / (cos(6.0 * PI / 16)),          4.0 / (cos(7.0 * PI / 16)));        unsigned char flat128Input[8][8];    signed short flatTestInput[8][8];    unsigned char twoComponents[8][8];    for (int y = 0; y < 8; y++)    {       for (int x = 0; x < 8; x++)       {          flatTestInput[y][x] = x == 0 && y == 0 ? 1264 : 0;          flat128Input[y][x] = 128;       }       twoComponents[y][0] = y + 105;       twoComponents[y][1] = y + 102;       twoComponents[y][2] = y + 97;       twoComponents[y][3] = y + 91;       twoComponents[y][4] = y + 84;       twoComponents[y][5] = y + 78;       twoComponents[y][6] = y + 73;       twoComponents[y][7] = y + 70;    }        printf("Reference DCT using double\n");    double values[8][8];    for (int y = 0; y < 8; y++)    {       for (int x = 0; x < 8; x++)       {          values[y][x] = refDCTEvaluateElement(x, y, twoComponents);          printf("%4f ", values[y][x]);       }              printf("\n");    }    printf("\n");        printf("Reference IDCT using double\n");    for (int y = 0; y < 8; y++)    {       for (int x = 0; x < 8; x++)       {          double val = refIDCTEvaluateElement(x, y, values);          printf("%4f ", val);       }              printf("\n");    }    printf("\n");        printf("AAN DCT using float\n");    signed short int output[8][8];    aanDCTEvaluateAll(twoComponents, output);    for (int y = 0; y < 8; y++)    {       for (int x = 0; x < 8; x++)       {          printf("%d ", output[y][x]);       }              printf("\n");    }    printf("\n");        printf("AAN IDCT using float\n");    unsigned char output2[8][8];    aanIDCTEvaluateAll(output, output2);    for (int y = 0; y < 8; y++)    {       for (int x = 0; x < 8; x++)       {          printf("%d ", output2[y][x]);       }       printf("\n");    }    printf("\n");        printf("AAN IDCT using float on flat test input\n");    unsigned char output3[8][8];    aanIDCTEvaluateAll(flatTestInput, output3);    for (int y = 0; y < 8; y++)    {       for (int x = 0; x < 8; x++)       {          printf("%d ", output3[y][x]);       }       printf("\n");    }    printf("\n");    } double refDCTEvaluateElement(unsigned char u, unsigned char v, unsigned char input[8][8]) {    double output = 0.0;        for (int x = 0; x < 8; x++)    {       for (int y = 0; y < 8; y++)       {          output += input[y][x] * cos((1 + (x << 1)) * u * PI / 16.0) * cos((1 + (y << 1)) * v * PI / 16.0);       }    }        if (v == 0)       output *= 1.0 / sqrt(2);        if (u == 0)       output *= 1.0 / sqrt(2);        output *= 0.25;        return output; } double refIDCTEvaluateElement(unsigned char x, unsigned char y, double input[8][8]) {    double output = 0.0;        for (int u = 0; u < 8; u++)    {       for (int v = 0; v < 8; v++)       {          output += input[v][u] * cos((1 + (x << 1)) * u * PI / 16.0) * cos((1 + (y << 1)) * v * PI / 16.0) *             (u == 0 ? 1.0 / sqrt(2) : 1) * (v == 0 ? 1.0 / sqrt(2) : 1);       }    }        output *= 0.25;        return output; } #define AANDCT4 0.707107 #define AANDCT6 0.382683 #define AANDCT2PLUS6 1.306563 #define AANDCT2MINUS6NEGATED -0.541196 #define AANDCTSCALE0 0.353553 #define AANDCTSCALE1 0.254898 #define AANDCTSCALE2 0.270598 #define AANDCTSCALE3 0.300672 #define AANDCTSCALE4 0.353553 #define AANDCTSCALE5 0.449988 #define AANDCTSCALE6 0.653281 #define AANDCTSCALE7 1.281458 // One-dimensional AAN DCT, including post-scaling. #define AAN_DCT_1D {\   const float s07 = SOURCE(0) + SOURCE(7);\   const float s16 = SOURCE(1) + SOURCE(6);\   const float s25 = SOURCE(2) + SOURCE(5);\   const float s34 = SOURCE(3) + SOURCE(4);\   const float d07 = SOURCE(0) - SOURCE(7);\   const float d16 = SOURCE(1) - SOURCE(6);\   const float d25 = SOURCE(2) - SOURCE(5);\   const float d34 = SOURCE(3) - SOURCE(4);\   const float tmp0 = s07 + s34;\   const float tmp1 = s16 + s25;\   const float tmp2 = s07 - s34;\   const float tmp3 = s16 - s25;\   const float tmp4 = d16 + d25;\   const float tmp6 = d07 + d16;\   const float tmp7 = -d25 - d34;\   const float m5 = AANDCT4 * (tmp3 + tmp2);\   const float m6 = AANDCT4 * tmp4;\   const float m7 = AANDCT6 * (tmp6 + tmp7);\   const float p5 = d07 + m6;\   const float p6 = d07 - m6;\   const float p7 = (AANDCT2PLUS6 * tmp6) - m7;\   const float p8 = (AANDCT2MINUS6NEGATED * tmp7) - m7;\   DESTINATION(0, AANDCTSCALE0 * (tmp0 + tmp1));\   DESTINATION(1, AANDCTSCALE1 * (p5 + p7));\   DESTINATION(2, AANDCTSCALE2 * (tmp2 + m5));\   DESTINATION(3, AANDCTSCALE3 * (p6 - p8));\   DESTINATION(4, AANDCTSCALE4 * (tmp0 - tmp1));\   DESTINATION(5, AANDCTSCALE5 * (p6 + p8));\   DESTINATION(6, AANDCTSCALE6 * (tmp2 - m5));\   DESTINATION(7, AANDCTSCALE7 * (p5 - p7));\}; #define AANIDCTSCALE0 5.656854 #define AANIDCTSCALE1 4.078365 #define AANIDCTSCALE2 4.329569 #define AANIDCTSCALE3 4.810759 #define AANIDCTSCALE4 5.656854 #define AANIDCTSCALE5 7.199810 #define AANIDCTSCALE6 10.452504 #define AANIDCTSCALE7 20.503324 #define AAN_IDCT_1D {\   const float F0 = SOURCE(0) * AANIDCTSCALE0;\   const float F1 = SOURCE(1) * AANIDCTSCALE1;\   const float F2 = SOURCE(2) * AANIDCTSCALE2;\   const float F3 = SOURCE(3) * AANIDCTSCALE3;\   const float F4 = SOURCE(4) * AANIDCTSCALE4;\   const float F5 = SOURCE(5) * AANIDCTSCALE5;\   const float F6 = SOURCE(6) * AANIDCTSCALE6;\   const float F7 = SOURCE(7) * AANIDCTSCALE7;\   const float tmp0 = F0 + F4;\   const float tmp1 = F0 - F4;\   const float tmp3 = AANDCT4 * (F2 - F6);\   const float tmp2 = F6 + F2 + tmp3;\   const float s07 = tmp2 + tmp0;\   const float s16 = tmp1 + tmp3;\   const float s25 = tmp1 - tmp3;\   const float s34 = tmp0 - tmp2;\   const float p6 = F5 + F3;\   const float p8 = F5 - F3;\   const float p5 = F1 + F7;\   const float p7 = F1 - F7;\   const float m6 = p5 - p6;\   const float m7 = (-p7 - p8) * AANDCT6;\   const float tmp7 = p8 * AANDCT2MINUS6NEGATED + m7;\   const float tmp4 = m6 * AANDCT4;\   const float tmp6 = p7 * AANDCT2PLUS6 + m7;\   const float d07 = tmp6 + p6 + p5;\   const float d16 = tmp6 + tmp4;\   const float d25 = tmp4 - tmp7;\   const float d34 = tmp7;\   DESTINATION(0, (s07 + d07));\   DESTINATION(1, (s16 + d16));\   DESTINATION(2, (s25 + d25));\   DESTINATION(3, (s34 - d34));\   DESTINATION(4, (s34 + d34));\   DESTINATION(5, (s25 - d25));\   DESTINATION(6, (s16 - d16));\   DESTINATION(7, (s07 - d07));\}; void aanDCTEvaluateAll(unsigned char input[8][8], signed short int output[8][8]) {    float tmp[8][8]; #define SOURCE(x) input[x][i] #define DESTINATION(x,y) tmp[x][i] = (y)    for (int i = 0; i < 8; i++)       AAN_DCT_1D #undef SOURCE #undef DESTINATION        #define SOURCE(x) tmp[i][x] #define DESTINATION(x,y) output[x][i] = lrintf(y);    for (int i = 0; i < 8; i++)       AAN_DCT_1D #undef SOURCE #undef DESTINATION } void aanIDCTEvaluateAll(signed short int input[8][8], unsigned char output[8][8]) {    float tmp[8][8]; #define SOURCE(x) input[x][i] #define DESTINATION(x,y) tmp[x][i] = (y) / 16.0f    for (int i = 0; i < 8; i++)       AAN_IDCT_1D #undef SOURCE #undef DESTINATION        #define SOURCE(x) tmp[i][x] #define DESTINATION(x,y) output[x][i] = lrintf(y) >> 4;    for (int i = 0; i < 8; i++)       AAN_IDCT_1D #undef SOURCE #undef DESTINATION }