logo标签列表

离散余弦变换（DCT）的实现C语言代码

cdct
9
我正在尝试在C语言中实现正向和反向离散余弦变换（DCT）。该代码通过dct()函数将单个输入像素块转换为变换矩阵，然后通过idct()函数将其转换回原始像素值。请参见附加的代码。我的idct输出是连续的244、116、244、116等数值。从idct值的外观来看，我的程序似乎没有起作用。有人能帮我解决问题并给我一个关于每个函数后应该得到什么结果的想法吗？显然，在idct之后，我应该得到非常接近原始输入矩阵的结果。

谢谢

 # 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);
}
你好Paul，根据Stack Overflow的指南，插入代码时请“缩进四个空格”，或者在列表中使用“预先格式化的块，缩进八个空格”..我按照这两种方法进行了操作，但输出结果却不符合预期..也许你能为我提供一些指导？ - user915071
只需选中代码并单击代码格式化按钮 {} - 这将自动缩进4个空格以获得所需的格式。我最初已经为您完成了此操作，但不知何故，您设法撤消了它。 - Paul R
也许开发人员会关心编辑说明，因为那不同于所给出的说明。我不得不在输入问题一分钟后进行编辑，但与此同时你已经格式化了代码。 - user915071
3个回答
9
在if语句中的Cv常量赋值至少有两个错别字。
    if (v == 0) {
        Cv = 1.0 / sqrt(2.0);
    } else {
        Cu = (1.0); // << this should be Cv = 1.0
    }

虽然没有仔细检查，但使用关于余弦变换的德语维基百科，以下代码可以正常工作... 我不想花时间去研究你是如何定义转换常数的。 我猜你需要确保使用正确的常数和反函数：

#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;
}

编辑 离散余弦变换（DCT）基于维基百科中的DCT-II公式的叉积。 离散余弦逆变换（IDCT）基于DCT-III公式的叉积，并且每个维度的归一化因子为2/N（因为这是文本中提到的DCT-II的反函数）。 编辑 我非常确定你版本中逆dct中的因子应该是sqrt(2)而不是1/sqrt(2)。

嗨Bort，我尝试在DEV-cpp中编译它，但是我得到了以下错误：“在函数float ** calloc_mat（int，int）'中：错误：从void *'到float **'的无效转换错误：从void *'到`float **'的无效转换”...也就是说，它在上面的代码的第13和14行抱怨。有什么想法吗？顺便问一下，你为什么要使用指向指针的指针？ - user915071
我使用指针，因为它们允许我重复使用适用于大于8x8和小于8x8的图像的代码。我假设您可以通过对calloc返回的void指针进行转换来修复错误，尽管在现代编译器中这不再是必需的（我使用的是gcc 4.5）。将强制转换添加为float **m = (float**) calloc(dimX, sizeof(float*));和float *p = (float *) calloc(dimX*dimY, sizeof(float));。 - Bort
它能正常工作。使用以下命令进行编译和运行： gcc dct.c -o dct -lm; ./dct - Juan Carlos Kuri Pinto
2

您不需要

#include <math.h>

这可能意味着编译器对数学函数做出了一些不正确的假设，例如它们都返回int。请注意，您调用的所有函数都需要在某个地方声明，C语言中并没有像sin()或printf()这样的“内置”函数（当然，对于后者，您正确地包含了stdin.h）。
另外，一旦您包含了<math.h>，您可以使用M_PI。
嗨，unwind，我已经进行了调整，但是我得到了相同的输出...根据GNU C： "Function: double cos (double x)此函数返回x的余弦值，其中x以弧度给出。返回值在-1到1的范围内。" 我应该将dct方程中的所有int值转换为double吗？另外，该函数返回弧度值，所以我应该转换回float吗？ - user915071
1
除了之前关于Cv常量（在dct()和idct()函数中）的错别字的答案外，您在反DCT公式（第2个）中使用不正确。 您每次在循环中都必须乘以Cv和Cu。因此，idct()的正确代码应为：
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");      
    }    
 }

在这种情况下，输出值接近于255、0、255、0等。
网页内容由stack overflow 提供, 点击上面的
可以查看英文原文，
原文链接