离散余弦变换的C++实现

期末会有一个大作业是搞这个东西，自己还不是完全明白，先转载过来，有空再看吧。。。。

最近网上有朋友因为要进行图像处理，而问及离散余弦正逆变换的问题，我在网上找了一些资料，网下也找了一些参考书，发现居然有很多错误，尤其是逆变换(IDCT)的公式都是错的！怪不得网友做完DCT后，再IDCT，却得不到原来的初始数据了(四舍五入的误差除外)，实在误人甚深。写此博文，以正其谬。

 

说明：

下面实现的代码，严格遵循离散余弦正逆变换的公式，仅供理解离散余弦逆变换的原理而写。并未考虑代码优化，也未考虑快速算法，快速算法见快速离散余弦变换代码实现(FDCT) 和快速逆离散余弦变换代码实现(FIDCT)。

离散余弦变换(DCT)公式：

 

 

离散余弦逆变换(IDCT)公式：

 

 

注：严格的离散余弦变换公式中，其中一个N可以写为M，且M和N可以不相等。这里由于为图像处理服务，所以都写成了N。

 

代码实现：

#include <iostream>

#include <math.h>

using namespace std;

 

#define NUM 8

#define PI 3.1415926

 

short round(double a)

{

         if (a >= 0)

         {

                   return (short)(a + 0.5);

         }

         else

         {

                   return (short)(a - 0.5);

         }

}

 

// DCT - Discrete Cosine Transform

void DCT(short data[NUM][NUM])

{

 

         short output[NUM][NUM];

         double ALPHA, BETA;

 

         short u = 0;

         short v = 0;

         short i = 0;

         short j = 0;

 

         for(u = 0; u < NUM; u++)

         {

                   for(v = 0; v < NUM; v++)

                   {

                            if(u == 0)

                            {

                                     ALPHA = sqrt(1.0 / NUM);

                            }

                            else

                            {

                                     ALPHA = sqrt(2.0 / NUM);

                            }

 

                            if(v == 0)

                            {

                                     BETA = sqrt(1.0 / NUM);

                            }

                            else

                            {

                                     BETA = sqrt(2.0 / NUM);

                            }

 

                            double tmp = 0.0;

                            for(i = 0; i < NUM; i++)

                            {

                                     for(j = 0; j < NUM; j++)

                                     {

                                               tmp += data[i][j] * cos((2*i+1)*u*PI/(2.0 * NUM)) * cos((2*j+1)*v*PI/(2.0 * NUM));

                                     }

                            }

                            output[u][v] = round(ALPHA * BETA * tmp);

                   }

         }

 

         memset(data, 0, NUM * NUM * sizeof(short));

         memcpy(data, output, NUM * NUM * sizeof(short));

}

 

// Inverse DCT

void IDCT(short data[NUM][NUM])

{

         short output[NUM][NUM];

         double ALPHA, BETA;

 

         short u = 0;

         short v = 0;

         short i = 0;

         short j = 0;

 

         for(i = 0; i < NUM; i++)

         {

                   for(short j = 0; j < NUM; j++)

                   {

                            double tmp = 0.0;

 

                            for(short u = 0; u < NUM; u++)

                            {

                                     for(v = 0; v < NUM; v++)

                                     {

                                               if(u == 0)

                                               {

                                                        ALPHA = sqrt(1.0 / NUM);

                                               }

                                               else

                                               {

                                                        ALPHA = sqrt(2.0 / NUM);

                                               }

                                               if(v == 0)

                                               {

                                                        BETA = sqrt(1.0 / NUM);

                                               }

                                               else

                                               {

                                                        BETA = sqrt(2.0 / NUM);

                                               }

                                               tmp += ALPHA * BETA * data[u][v] * cos((2*i+1)*u*PI/(2.0 * NUM)) * cos((2*j+1)*v*PI/(2.0 * NUM));

                                     }

                            }

                            output[i][j] = round(tmp);

                   }

         }

         memset(data, 0, NUM * NUM * sizeof(short));

         memcpy(data, output, NUM * NUM * sizeof(short));

}

 

int main(void)

{

         short i = 0;

         short j = 0;

         short u = 0;

         short v = 0;

        

         // 8 x 8 的图像数据

         short input[NUM][NUM] =

         {

                   {89, 101, 114, 125, 126, 115, 105, 96},

                   {97, 115, 131, 147, 149, 135, 123, 113},

                   {114, 134, 159, 178, 175, 164, 149, 137},

                   {121, 143, 177, 196, 201, 189, 165, 150},

                   {119, 141, 175, 201, 207, 186, 162, 144},

                   {107, 130, 165, 189, 192, 171, 144, 125},

                   {97, 119, 149, 171, 172, 145, 117, 96},

                   {88, 107, 136, 156, 155, 129, 97, 75}

         };

 

         DCT(input);

         cout << "The result of DCT:" << endl;

         for(u = 0; u < NUM; u++)

         {

                   for(v = 0; v < NUM; v++)

                   {

                            cout << input[u][v] << '/t';

                   }

                   cout << endl;

         }

 

         IDCT(input);

         cout << "The result of IDCT:" << endl;

         for(i = 0; i < NUM; i++)

         {

                   for(short j = 0; j < NUM; j++)

                   {

                            cout << input[i][j] << '/t';

                   }

                   cout << endl;

         }

 

         return 0;

}

运算结果如下：

The result of DCT:
1125    -32    -185     -7      2     -1     -2      2 
   -22    -16        45     -3     -2     0      -2     -2 
 -165     32        17     2        1    -1     -3       0 
      -7      -4          0     2        2    -1     -1       2 
     -2        0          0     3        0     0       2       1 
      3        1          1     -1      -2     0       2       0 
      0        0          2     -1      -1     2       1      -1 
      0        3          1     -1       2     1      -2       0 

The result of IDCT:
   89   101   114   125   126   115   105     96 
   97   115   131   147   149   135   123   113 
114   134   159   178   175   164   149   137 
121   143   177   196   201   189   165   150 
119   141   175   201   207   186   162   144 
107   130   165   189   192   171   144   125 
  96   119   150   171   172   145   116      96 
  88   107   136   156   155   129      97     75

转自http://blog.csdn.net/pathuang68/article/details/4193990