Discrete Cosine Transformation DCT C

Question

Im trying to implement direct and inverse discrete cosine transform (DCT) to C. The code is to transfer one input block of pixels to the transformation matrix through the dct () function, and then return to the original pixel values ​​through the idct () function. See the attached code. My conclusion from idct is the sequential values ​​244, 116, 244, 116, etc. From looking at idct values, it doesn't look like my program works. Can someone help me and give me an idea of ​​what results I should expect after each function? Obviously, after idct, I should get prety close to the original input matrix.

thanks

# include <stdio.h> # define PI 3.14 void dct(float [][]); // Function prototypes void idct(float [][]); // Function prototypes void dct(float inMatrix[8][8]){ double dct, Cu, sum, Cv; int i, j, u, h = 0, v; FILE * fp = fopen("mydata.csv", "w"); float dctMatrix[8][8], greyLevel; for (u = 0; u < 8; ++u) { for (v = 0; v < 8; ++v) { if (u == 0) { Cu = 1.0 / sqrt(2.0); } else { Cu = 1.0; } if (v == 0) { Cv = 1.0 / sqrt(2.0); } else { Cu = (1.0); } sum = 0.0; for (i = 0; i < 8; i++) { for (j = 0; j < 8; j++) { // Level around 0 greyLevel = inMatrix[i][j]; dct = greyLevel * cos((2 * i + 1) * u * PI / 16.0) * cos((2 * j + 1) * v * PI / 16.0); sum += dct; } } dctMatrix[u][v] = 0.25 * Cu * Cv * sum; fprintf(fp, "\n %f", dctMatrix[u][v]); } fprintf(fp, "\n"); } idct(dctMatrix); } void idct(float dctMatrix[8][8]){ double idct, Cu, sum, Cv; int i, j, u, v; float idctMatrix[8][8], greyLevel; FILE * fp = fopen("mydata.csv", "a"); fprintf(fp, "\n Inverse DCT"); for (i = 0; i < 8; ++i) { for (j = 0; j < 8; ++j) { sum = 0.0; for (u = 0; u < 8; u++) { for (v = 0; v < 8; v++) { if (u == 0) { Cu = 1.0 / sqrt(2.0); } else { Cu = 1.0; } if (v == 0) { Cv = 1.0 / sqrt(2.0); } else { Cu = (1.0); } // Level around 0 greyLevel = dctMatrix[u][v]; idct = (greyLevel * cos((2 * i + 1) * u * PI / 16.0) * cos((2 * j + 1) * v * PI / 16.0)); sum += idct; } } idctMatrix[i][j] = 0.25 * Cu * Cv * sum; fprintf(fp, "\n %f", idctMatrix[i][j]); } fprintf(fp, "\n"); } } int main() { float testBlockA[8][8] = { {255, 255, 255, 255, 255, 255, 255, 255}, {255, 255, 255, 255, 255, 255, 255, 255}, {255, 255, 255, 255, 255, 255, 255, 255}, {255, 255, 255, 255, 255, 255, 255, 255}, {255, 255, 255, 255, 255, 255, 255, 255}, {255, 255, 255, 255, 255, 255, 255, 255}, {255, 255, 255, 255, 255, 255, 255, 255}, {255, 255, 255, 255, 255, 255, 255, 255} }, testBlockB[8][8] = {{255, 0, 255, 0, 255, 0, 255, 0}, {0, 255, 0, 255, 0, 255, 0, 255}, {255, 0, 255, 0, 255, 0, 255, 0}, {0, 255, 0, 255, 0, 255, 0, 255}, {255, 0, 255, 0, 255, 0, 255, 0}, {0, 255, 0, 255, 0, 255, 0, 255}, {255, 0, 255, 0, 255, 0, 255, 0}, {0, 255, 0, 255, 0, 255, 0, 255} }; dct(testBlockB); } 

Answer 1

If statements have at least two typos in the Cv constant:

  if (v == 0) { Cv = 1.0 / sqrt(2.0); } else { Cu = (1.0); // << this should be Cv = 1.0 } 

Did not check too correctly. Using the German cosplay transform wikipedia , the following code works ... I didn’t want to waste time figuring out how you determined which conversion constant. I think you need to make sure that you are using the correct constants and inverse functions:

 #include <stdio.h> #include <math.h> #include <stdlib.h> void dct(float **DCTMatrix, float **Matrix, int N, int M); void write_mat(FILE *fp, float **testRes, int N, int M); void idct(float **Matrix, float **DCTMatrix, int N, int M); float **calloc_mat(int dimX, int dimY); void free_mat(float **p); float **calloc_mat(int dimX, int dimY){ float **m = calloc(dimX, sizeof(float*)); float *p = calloc(dimX*dimY, sizeof(float)); int i; for(i=0; i <dimX;i++){ m[i] = &p[i*dimY]; } return m; } void free_mat(float **m){ free(m[0]); free(m); } void write_mat(FILE *fp, float **m, int N, int M){ int i, j; for(i =0; i< N; i++){ fprintf(fp, "%f", m[i][0]); for(j = 1; j < M; j++){ fprintf(fp, "\t%f", m[i][j]); } fprintf(fp, "\n"); } fprintf(fp, "\n"); } void dct(float **DCTMatrix, float **Matrix, int N, int M){ int i, j, u, v; for (u = 0; u < N; ++u) { for (v = 0; v < M; ++v) { DCTMatrix[u][v] = 0; for (i = 0; i < N; i++) { for (j = 0; j < M; j++) { DCTMatrix[u][v] += Matrix[i][j] * cos(M_PI/((float)N)*(i+1./2.)*u)*cos(M_PI/((float)M)*(j+1./2.)*v); } } } } } void idct(float **Matrix, float **DCTMatrix, int N, int M){ int i, j, u, v; for (u = 0; u < N; ++u) { for (v = 0; v < M; ++v) { Matrix[u][v] = 1/4.*DCTMatrix[0][0]; for(i = 1; i < N; i++){ Matrix[u][v] += 1/2.*DCTMatrix[i][0]; } for(j = 1; j < M; j++){ Matrix[u][v] += 1/2.*DCTMatrix[0][j]; } for (i = 1; i < N; i++) { for (j = 1; j < M; j++) { Matrix[u][v] += DCTMatrix[i][j] * cos(M_PI/((float)N)*(u+1./2.)*i)*cos(M_PI/((float)M)*(v+1./2.)*j); } } Matrix[u][v] *= 2./((float)N)*2./((float)M); } } } int main() { float testBlockA[8][8] = { {255, 255, 255, 255, 255, 255, 255, 255}, {255, 255, 255, 255, 255, 255, 255, 255}, {255, 255, 255, 255, 255, 255, 255, 255}, {255, 255, 255, 255, 255, 255, 255, 255}, {255, 255, 255, 255, 255, 255, 255, 255}, {255, 255, 255, 255, 255, 255, 255, 255}, {255, 255, 255, 255, 255, 255, 255, 255}, {255, 255, 255, 255, 255, 255, 255, 255} }, testBlockB[8][8] = {{255, 0, 255, 0, 255, 0, 255, 0}, {0, 255, 0, 255, 0, 255, 0, 255}, {255, 0, 255, 0, 255, 0, 255, 0}, {0, 255, 0, 255, 0, 255, 0, 255}, {255, 0, 255, 0, 255, 0, 255, 0}, {0, 255, 0, 255, 0, 255, 0, 255}, {255, 0, 255, 0, 255, 0, 255, 0}, {0, 255, 0, 255, 0, 255, 0, 255} }; FILE * fp = fopen("mydata.csv", "w"); int dimX = 8, dimY = 8; int i, j; float **testBlock = calloc_mat(dimX, dimY); float **testDCT = calloc_mat(dimX, dimY); float **testiDCT = calloc_mat(dimX, dimY); for(i = 0; i<dimX; i++){ for(j = 0; j<dimY; j++){ testBlock[i][j] = testBlockB[i][j]; } } dct(testDCT, testBlock, dimX, dimY); write_mat(fp, testDCT, dimX, dimY); idct(testiDCT, testDCT, dimX, dimY); write_mat(fp, testiDCT, dimX, dimY); fclose(fp); free_mat(testBlock); free_mat(testDCT); free_mat(testiDCT); return 0; } 

Edit dct is based on the cross-product of the DCT-II formula on the wiki. idct is based on the cross product of the DCT-III formula with a normalization factor of 2 / N per size (since this is the opposite of DCT-II, as indicated in the text). Edit I am sure that the coefficient in the inverse dct should be sqrt (2), not 1 / sqrt (2) in your version.

Answer 2

You do not

 #include <math.h> 

This may mean that the compiler accepts things about mathematical functions that are not true, for example, that they all return an int . Note that all functions that you call must be declared somewhere, C no longer has a "built-in" sin() than a built-in printf() (for the latter, you include stdin.h , of course).

Alternatively, you can use M_PI after including <math.h> .

Answer 3

In addition to the previous answer to typos in the constant Cv (as in the dct () and idct () functions, you used the reverse DCT (2nd) formula incorrectly. You had to multiply by Cv and Cu every time in the loop. So, the correct idct code () it should be:

 void idct(float dctMatrix[8][8]){ double idct, Cu, sum, Cv; int i, j, u, v; float idctMatrix[8][8], greyLevel; FILE * fp = fopen("mydata.csv", "a"); fprintf(fp, "\n Inverse DCT"); for (i = 0; i < 8; ++i) { for (j = 0; j < 8; ++j) { sum = 0.0; for (u = 0; u < 8; u++) { for (v = 0; v < 8; v++) { if (u == 0) { Cu = 1.0 / sqrt(2.0); } else { Cu = 1.0; } if (v == 0) { Cv = 1.0 / sqrt(2.0); } else { Cv = (1.0); //mistake was here - the same is in dct() } greyLevel = dctMatrix[u][v]; // Multiply by Cv and Cu here! idct = (greyLevel * Cu * Cv * cos((2 * i + 1) * u * PI / 16.0) * cos((2 * j + 1) * v * PI / 16.0)); sum += idct; } } // not "* Cv * Cu" here! idctMatrix[i][j] = 0.25 * sum; fprintf(fp, "\n %f", idctMatrix[i][j]); } fprintf(fp, "\n"); } } 

In this case, the output values ​​are close to 255, 0, 255, 0, etc.