图像的正交变换---离散余弦变换

/*************************************************************************

*

* 函数名称：

*   DCT()

*

* 参数:

*   double * f				- 指向时域值的指针

*   double * F				- 指向频域值的指针

*   r						－2的幂数

*

* 返回值:

*   无。

*

* 说明:

*   该函数用来实现一维快速离散余弦变换

*

************************************************************************/

void DCT(double*f,double*F,int r);

/*************************************************************************

*

* 函数名称：

*   IDCT()

*

* 参数:

*   double * F				- 指向频域值的指针

*   double * f				- 指向时域值的指针

*   r						－2的幂数

*

* 返回值:

*   无。

*

* 说明:

*   该函数实现一维快速离散余弦逆变换

*

************************************************************************/

void IDCT(double*F,double*f,int r);

/*************************************************************************

*

* 函数名称：

*   FreqDCT()

*

* 参数:

*   double* f			- 输入的时域序列

*   double* F			- 输出的频域序列

*	 LONG width			- 图象宽度

*	 LONG height		- 图象高度

*

* 返回值:

*   BOOL               - 成功返回TRUE，否则返回FALSE。

*

* 说明:

*   该函数进行二维快速离散余弦变换。

*

************************************************************************/

BOOL FreqDCT(double*f,double*F, LONG width, LONG height);

/*************************************************************************

*

* 函数名称：

*   IFreqDCT()

*

* 参数:

*   double* f			- 输入的时域序列

*   double* F			- 输出的频域序列

*	 LONG width			- 图象宽度

*	 LONG height		- 图象高度

*

* 返回值:

*   BOOL               - 成功返回TRUE，否则返回FALSE。

*

* 说明:

*   该函数进行二维快速离散余弦逆变换。

*

************************************************************************/

BOOL IFreqDCT(double*f,double*F, LONG lWidth, LONG lHeight);

/*************************************************************************

*

* 函数名称：

*   DIBDCT()

*

* 参数:

*   BYTE* bmp,LONG width,LONG height

*

* 返回值:

*   BOOL               - 成功返回TRUE，否则返回FALSE。

*

* 说明:

*   该函数用来对图像进行离散余弦变换。

*

************************************************************************/

BOOLBmpDCT(BYTE* bmp,LONG width,LONG height);

voidMyProcess::DCT(double*f,double*F,int r)

{

// 循环变量

int		i;

// 中间变量

double	dTemp;

// 计算离散余弦变换点数

LONG N =1<<r;

// 申请并分配内存

complex<double>*XIn;

complex<double>*XOut;

XIn=new complex<double>[N *2];

XOut=new complex<double>[N *2];

// 赋初值为0

memset(XIn,0,sizeof(complex<double>)* N *2);

memset(XOut,0,sizeof(complex<double>)* N *2);

// 将时域点写入数组X

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

XIn[i]= complex<double>(*(f + i),0);

// 调用快速付立叶变换

FFT(XIn,XOut, r +1);

// 调整系数

dTemp =1/ sqrt(N);

// 求F[0]

F[0]=XOut[0].real()* dTemp;

dTemp *= sqrt(2);

// 求F[u]

for(i =1; i < N; i++)

*(F + i)=(XOut[i].real()* cos(i * PI /(N *2))+XOut[i].imag()* sin(i * PI /(N *2)))* dTemp;

// 释放内存

delete[]XIn;

delete[]XOut;

}

voidMyProcess::IDCT(double*F,double*f,int r)

{

// 循环变量

int		i;

// 中间变量

double	dTemp, d0;

// 计算离散余弦变换点数

LONG N =1<<r;

// 分配内存

complex<double>*XIn;

complex<double>*XOut;

XIn=new complex<double>[N *2];

XOut=new complex<double>[N *2];

// 赋初值为0

memset(XIn,0,sizeof(complex<double>)* N *2);

memset(XOut,0,sizeof(complex<double>)* N *2);

// 将频域变换后点写入数组X

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

XIn[i]= complex<double>(F[i]* cos(i * PI /(N *2)), F[i]* sin(i * PI /(N *2)));

// 调用快速付立叶反变换

IFFT(XIn,XOut, r +1);

// 调整系数

dTemp = sqrt(2.0/ N);

d0 =(sqrt(1.0/ N)- dTemp)* F[0];

// 计算f(x)

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

f[i]= d0 +XOut[i].real()* dTemp *2* N;

// 释放内存

delete[]XIn;

delete[]XOut;

}

BOOL MyProcess::FreqDCT(double*f,double*F, LONG width, LONG height)

{

// 循环变量

LONG	i;

LONG	j;

LONG    k;

// 进行离散余弦变换的宽度和高度（2的整数次方）

LONG w =1;

LONG h =1;

int wp =0;

int hp =0;

// 计算进行离散余弦变换的宽度和高度（2的整数次方）

while(w < width/3)

{

w *=2;

wp++;

}

while(h < height)

{

h *=2;

hp++;

}

// 分配内存

double*TempIn=newdouble[h];

double*TempOut=newdouble[h];

// 对y方向进行离散余弦变换

for(i =0; i < w *3; i++)

{

// 抽取数据

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

TempIn[j]= f[j * w *3+ i];

// 一维快速离散余弦变换

DCT(TempIn,TempOut, hp);

// 保存变换结果

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

f[j * w *3+ i]=TempOut[j];

}

// 释放内存

deleteTempIn;

deleteTempOut;

// 分配内存

TempIn=newdouble[w];

TempOut=newdouble[w];

// 对x方向进行快速离散余弦变换

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

{

for(k =0; k <3; k++)

{

// 抽取数据

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

TempIn[j]= f[i * w *3+ j *3+ k];

// 一维快速离散余弦变换

DCT(TempIn,TempOut, wp);

// 保存变换结果

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

F[i * w *3+ j *3+ k]=TempOut[j];

}

}

// 释放内存

deleteTempIn;

deleteTempOut;

return TRUE;

}

BOOL MyProcess::IFreqDCT(double*f,double*F, LONG width, LONG height)

{

// 循环变量

LONG	i;

LONG	j;

LONG    k;

// 进行离散余弦变换的宽度和高度（2的整数次方）

LONG w =1;

LONG h =1;

int wp =0;

int hp =0;

// 计算进行付立叶变换的宽度和高度（2的整数次方）

while(w < width/3)

{

w *=2;

wp++;

}

while(h < height)

{

h *=2;

hp++;

}

// 分配内存

double*TempIn=newdouble[w];

double*TempOut=newdouble[w];

// 对x方向进行快速付立叶变换

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

{

for(k =0; k <3; k++)

{

// 抽取数据

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

TempIn[j]= F[i * w *3+ j *3+ k];

// 一维快速傅立叶变换

IDCT(TempIn,TempOut, wp);

// 保存变换结果

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

F[i * w *3+ j *3+ k]=TempOut[j];

}

}

// 释放内存

deleteTempIn;

deleteTempOut;

TempIn=newdouble[h];

TempOut=newdouble[h];

// 对y方向进行快速付立叶变换

for(i =0; i < w *3; i++)

{

// 抽取数据

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

TempIn[j]= F[j * w *3+ i];

// 一维快速傅立叶变换

IDCT(TempIn,TempOut, hp);

// 保存变换结果

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

F[j * w *3+ i]=TempOut[j];

}

// 释放内存

deleteTempIn;

deleteTempOut;

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

{

for(j =0; j < w *3; j++)

{

if(i < height && j < width)

*(f + i * width + j)= F[i * w *3+ j];

}

}

return TRUE;

}

BOOL MyProcess::BmpDCT(BYTE* bmp,LONG width,LONG height)

{

// 进行离散余弦变换的宽度和高度（2的整数次方）

LONG	i;

LONG	j;

// 进行离散余弦变换的宽度和高度（2的整数次方）

LONG w =1;

LONG h =1;

int wp =0;

int hp =0;

// 计算进行离散余弦变换的宽度和高度（2的整数次方）

while(w < width/3)

{

w *=2;

wp++;

}

while(h < height)

{

h *=2;

hp++;

}

// 分配内存

double*f =newdouble[w * h *3];

double*F =newdouble[w * h *3];

// 向时域赋值并补零

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

{

for(j =0; j < w *3; j++)

{

if(i < height && j < width)

f[i * w *3+ j]= bmp[width * i + j];

else

f[w * i *3+ j]=0.0f;

}

}

// 进行频谱分析

if(FreqDCT(f, F,width, height)== FALSE)

return FALSE;

// 更新所有象素

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

{

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

{

// 判断是否超过255

if(fabs(F[i * w *3+ j])>255)

{

// 对于超过的，直接设置为255

bmp[width *(height -1- i)+ j]=255;

}

else

{

// 如果没有超过，则按实际计算结果赋值

bmp[width *(height -1- i)+ j]= fabs(F[i * w *3+ j]);

}

}

}

// 释放内存

delete[]f;

delete[] F;

// 返回

return TRUE;

}