C - DCT 구현해보기
Updated:
제가 C언어를 공부하면서 배웠던 문제를 공유하겠습니다.
이번 post는 DCT를 구현해보는 저장하는 문제입니다.
입력의 lena파일은 512*512 size의 raw 파일입니다.
파일 입,출력 포스팅에 있는 주석들은 생략하였습니다.
영상에 대해 DCT를 수행해보고 그 결과로 Inverse DCT를 수행해서 영상을 복원해봅시다.
8x8 block과 전체 영상에 대해 DCT를 수행했습니다. 양자화도 진행했습니다.
DCT의 수식적인 부분은 코드에서 확인할 수 있으니 설명하지 않겠습니다!
DCT를 사용하는 힌트는 아래와 같습니다 :)
- 우리가 접하는 모든 영상은 고주파 성분보다 저주파 성분을 훨씬 많이 가지고 있음
- 주파수 성분으로 영상을 분할한 뒤 파일을 압축한다면??
#pragma warning(disable:4996)
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#define PI 3.141592654
#define DIMENTION 8 // 4, 8, 16, 32, 64
#define DIMENTION1 512
#define HEI 512
#define WID 512
// 더블 함수 포인터
unsigned char **UCalloc(int width, int height)
{
int i, j;
unsigned char **ptr;
if ((ptr = (unsigned char**)malloc(height * sizeof(unsigned char*))) == NULL)
{
printf("\nMemory allocation failure\n");
exit(1);
}
for (i = 0; i<height; i++)
{
if ((ptr[i] = (unsigned char*)malloc(width * sizeof(unsigned char))) == NULL)
{
printf("\nMemory allocation failure\n");
exit(1);
}
}
for (i = 0; i<height; i++)
for (j = 0; j<width; j++)
ptr[i][j] = 0;
printf("\nMEMORY ALLOCATION(UNSIGNED CHAR) OK!\n");
return ptr;
} // 더블 함수 포인터
// 닫는 free 함수 묶은 것
void UCfree(unsigned char **ptr, int height)
{
int i;
for (i = 0; i<height; i++)
free(ptr[i]);
free(ptr);
printf("\nMEMORY FREE(UNSIGNED CHAR) OK!\n");
}
// filename을 읽어서 source(2차원배열)로 저장, ([height][width]) 후 닫음
void Readfile(char *filename, unsigned char **source, int width, int height)
{
int i, j;
FILE *readf;
if ((readf = fopen(filename, "rb")) == NULL)
{
fprintf(stderr, "ERROR : File cannot open : %s \n", filename);
exit(-1);
}
for (i = 0; i<height; i++)
for (j = 0; j<width; j++)
source[i][j] = (unsigned char)getc(readf);
printf("\n%sFILE READING OK!\n", filename);
fclose(readf);
}
// result(2차원 배열)를 filename 에 저장 후 닫음
void Writefile(char *filename, unsigned char **result, int width, int height)
{
int i, j;
FILE *writef;
if ((writef = fopen(filename, "wb")) == NULL)
{
fprintf(stderr, "ERROR : File cannot open : %s \n", filename);
exit(-1);
}
for (i = 0; i<height; i++)
for (j = 0; j<width; j++)
putc((unsigned char)result[i][j], writef);
printf("\n%sFILE WRITING OK!\n", filename);
fclose(writef);
}
double trans[HEI][WID];
double transcp[HEI][WID];
double transqua[HEI][WID];
double trans_quant[HEI][WID];
double trans_w[HEI][WID];
double trans_wh[HEI][WID];
double trans_quantoidct[HEI][WID];
int main(void)
{
int i, j, k, l, m, n;
double cu, cv, result;
//double MSE, SUM, PSNR = 0;
unsigned char Qaunt[8][8] = {
{ 16,11,10,16,24,40,51,61 },
{ 12,12,14,19,26,58,60,55 },
{ 14,13,16,24,40,57,69,56 },
{ 14,17,22,29,51,87,80,62 },
{ 18,22,37,56,68,109,103,77 },
{ 24,35,55,64,81,104,113,92 },
{ 49,64,78,87,103,121,120,101 },
{ 72,82,95,98,112,100,103,99 }
};
unsigned char **image = UCalloc(512, 512);
unsigned char **image_1 = UCalloc(512, 512);
unsigned char **dct = UCalloc(512, 512);
unsigned char **tran_quant = UCalloc(512, 512);
unsigned char **dct_w = UCalloc(512, 512);
unsigned char **image1 = UCalloc(512, 512);
unsigned char **image2 = UCalloc(512,512);
unsigned char **image3 = UCalloc(512, 512);
Readfile("Lena_512.raw", image, 512, 512);
Readfile("Lena_512.raw", image_1, 512, 512);
//Forward DCT
for (i = 0; i<(512 / DIMENTION); i++)
{
for (j = 0; j<(512 / DIMENTION); j++)
{
for (k = 0; k<DIMENTION; k++)
{
for (l = 0; l<DIMENTION; l++)
{
result = 0.0;
for (m = 0; m<DIMENTION; m++)
{
for (n = 0; n<DIMENTION; n++)
{
result += (double)cos((((2 * (double)m + 1)*(double)k*PI) / (DIMENTION * 2)))*cos((((2 * (double)n + 1)*(double)l*PI) / (DIMENTION * 2)))*(double)image[(i*DIMENTION) + m][(j*DIMENTION) + n];
}
}
if (k == 0)
cu = 1 / sqrt(2);
else
cu = 1;
if (l == 0)
cv = 1 / sqrt(2);
else
cv = 1;
trans[(i*DIMENTION) + k][(j*DIMENTION) + l] = result*((cu*cv) / 4);
}
}
}
}
// DCT된 상태 출력
for (i = 0; i<512; i++)
{
for (j = 0; j<512; j++)
{
if (trans[i][j]<0)
transcp[i][j] = 0;
else
transcp[i][j] = trans[i][j];
if (transcp[i][j]>255)
transcp[i][j] = 255;
dct[i][j] = (unsigned char)floor(transcp[i][j] + 0.5); // 반올림
}
}
//printf("%f\n",PSNR);
Writefile("8X8_FDCT_Lena.raw", dct, 512, 512);
//////////////////////////////////////////////// 양자화 시킨 DCT
for (i = 0; i < 512; i += 8)
{
for (j = 0; j < 512; j += 8)
{
for (k = 0; k < DIMENTION; k++)
{
for (l = 0; l < DIMENTION; l++)
{
transqua[i + k][j + l] =(trans[i + k][j + l] / (double)Qaunt[k][l]);
}
}
}
}
for (i = 0; i<512; i++)
{
for (j = 0; j<512; j++)
{
if (transqua[i][j]<0)
trans_quant[i][j] = 0;
else
trans_quant[i][j] = transqua[i][j];
if (trans_quant[i][j]>255)
trans_quant[i][j] = 255;
tran_quant[i][j] = (unsigned char)floor(trans_quant[i][j] + 0.5); // 반올림
}
}
Writefile("Quant_FDCT_Lena.raw", tran_quant, 512, 512);
UCfree(dct, 512);
UCfree(tran_quant, 512);
/////////////////////////////////////////////전체에 걸친 DCT
for (i = 0; i<(512 / DIMENTION1); i++)
{
for (j = 0; j<(512 / DIMENTION1); j++)
{
for (k = 0; k<DIMENTION1; k++)
{
for (l = 0; l<DIMENTION1; l++)
{
result = 0.0;
for (m = 0; m<DIMENTION1; m++)
{
for (n = 0; n<DIMENTION1; n++)
{
result += (double)cos((((2 * (double)m + 1)*(double)k*PI) / (DIMENTION1 * 2)))*cos((((2 * (double)n + 1)*(double)l*PI) / (DIMENTION1 * 2)))*(double)image_1[(i*DIMENTION1) + m][(j*DIMENTION1) + n];
}
}
if (k == 0)
cu = 1 / sqrt(2);
else
cu = 1;
if (l == 0)
cv = 1 / sqrt(2);
else
cv = 1;
trans_w[(i*DIMENTION1) + k][(j*DIMENTION1) + l] = result*((cu*cv) / 256);
}
}
}
}
for (i = 0; i<512; i++)
{
for (j = 0; j<512; j++)
{
if (trans_w[i][j]<0)
trans_wh[i][j] = 0;
else
trans_wh[i][j] = trans_w[i][j];
if (trans_wh[i][j]>255)
trans_wh[i][j] = 255;
dct_w[i][j] = (unsigned char)floor(trans_wh[i][j] + 0.5); // 반올림
}
}
Writefile("FDCT_Lena_w.raw", dct_w, 512, 512);
UCfree(dct_w, 512);
//Inverse DCT
for (i = 0; i<(512 / DIMENTION); i++)
{
for (j = 0; j<(512 / DIMENTION); j++)
{
for (n = 0; n<DIMENTION; n++)// x
{
for (m = 0; m<DIMENTION; m++) // y
{
result = 0.0;
for (l = 0; l<DIMENTION; l++) // v
{
for (k = 0; k<DIMENTION; k++) // u
{
if (k == 0)
cu = 1 / sqrt(2);
else
cu = 1;
if (l == 0)
cv = 1 / sqrt(2);
else
cv = 1;
result += ((cu*cv) / 4)*cos((((2 * (double)n + 1)*(double)k*PI) / (DIMENTION * 2)))*cos((((2 * (double)m + 1)*(double)l*PI) / (DIMENTION * 2)))*trans[(i*DIMENTION) + k][(j*DIMENTION) + l];
}
}
image1[(i*DIMENTION) + n][(j*DIMENTION) + m] = (unsigned char)floor(result + 0.5);
}
}
}
}
Writefile("8X8_IDCT_Lena.raw", image1, 512, 512);
UCfree(image1, 512);
//양자화가 포함된 Inverse DCT
for (i = 0; i < 512; i += 8)
{
for (j = 0; j < 512; j += 8)
{
for (k = 0; k < DIMENTION; k++)
{
for (l = 0; l < DIMENTION; l++)
{
trans_quantoidct[i + k][j + l] = transqua[i + k][j + l] * (double)Qaunt[k][l];
}
}
}
}
for (i = 0; i<(512 / DIMENTION); i++)
{
for (j = 0; j<(512 / DIMENTION); j++)
{
for (n = 0; n<DIMENTION; n++)// x
{
for (m = 0; m<DIMENTION; m++) // y
{
result = 0.0;
for (l = 0; l<DIMENTION; l++) // v
{
for (k = 0; k<DIMENTION; k++) // u
{
if (k == 0)
cu = 1 / sqrt(2);
else
cu = 1;
if (l == 0)
cv = 1 / sqrt(2);
else
cv = 1;
result += ((cu*cv) / 4)*cos((((2 * (double)n + 1)*(double)k*PI) / (DIMENTION * 2)))*cos((((2 * (double)m + 1)*(double)l*PI) / (DIMENTION * 2)))*trans_quantoidct[(i*DIMENTION) + k][(j*DIMENTION) + l];
}
}
image2[(i*DIMENTION) + n][(j*DIMENTION) + m] = (unsigned char)floor(result + 0.5);
}
}
}
}
Writefile("Quant_IDCT_Lena.raw", image2, 512, 512);
UCfree(image2, 512);
//
// 전체에 걸친 IDCT
for (i = 0; i<(512 / DIMENTION1); i++)
{
for (j = 0; j<(512 / DIMENTION1); j++)
{
for (n = 0; n<DIMENTION1; n++)// x
{
for (m = 0; m<DIMENTION1; m++) // y
{
result = 0.0;
for (l = 0; l<DIMENTION1; l++) // v
{
for (k = 0; k<DIMENTION1; k++) // u
{
if (k == 0)
cu = 1 / sqrt(2);
else
cu = 1;
if (l == 0)
cv = 1 / sqrt(2);
else
cv = 1;
result += ((cu*cv) / 256)*cos((((2 * (double)n + 1)*(double)k*PI) / (DIMENTION1 * 2)))*cos((((2 * (double)m + 1)*(double)l*PI) / (DIMENTION1 * 2)))*trans_w[(i*DIMENTION1) + k][(j*DIMENTION1) + l];
}
}
image3[(i*DIMENTION1) + n][(j*DIMENTION1) + m] = (unsigned char)floor(result + 0.5);
}
}
}
}
Writefile("IDCT_Lena_w.raw", image3, 512, 512);
UCfree(image3, 512);
// 원본과 복구한 그림 비교
k = 1;
for (i = 0; i<512; i++)
{
for (j = 0; j<512; j++)
{
if (image[i][j] != dct[i][j])
{
printf("두 파일의 %d행 %d열 픽셀이 틀립니다.\n", i + 1, j + 1);
k = 0;
}
}
}
// 비교 후 무결성 여부 발표
if (k == 0)
printf("원복과 복구한 이미지가 같지 않습니다.\n\n");
else
printf("원복과 복구한 이미지가 동일합니다.\n\n");
return 0;
}
Leave a comment