C实现 离散余弦转换DCT 源代码

算法特点：

1：将调用FFT-1D和FFT-2D实现DCT-1D和DCT-2D

2:初始化函数将任意长度数组进行调整使其满足FFT处理的2的整数次幂。

3:由于核心算法采用的是蝶形运算的FFT-1D，所以运算速度特别快。

4：FFT算法的源代码参照： C实现 快速傅立叶变换FFT FFT-1D FFT-2D 源代码

源代码：

void InitTDAndFD(double *&TD,double *&FD,unsigned char *&data,int lWidth,int lHeight)//初始化TD和FD（DCT操作）
{
	int w=1;
	int h=1;
	int wp=0;//存储w的2的幂数
	int hp=0;//存储h的2的幂数
 
	//计算刚好大于或等于lWidth，lHeight的2的整数次幂，和相应的幂数///
	while(w<lWidth)
	{
		w=w*2;
		wp++;
	}
	while(h<lHeight)
	{
		h=h*2;
		hp++;
	}	
	TD=new double[w*h];
	FD=new double[w*h];
	memset(TD,0,sizeof(double)*w*h);
	memset(FD,0,sizeof(double)*w*h);

	for(int i=0;i<h;i++)
	{
		if(i<lHeight)
		{
			for(int j=0;j<w;j++)
			{
				if(j< lWidth)
				{
					TD[i*w+j]=(double)data[i*lWidth+j];//将char数据，准换为实数为data数据，虚数为0的复数
				}
				else
				{
					TD[i*w+j]=0.0;//对于超出原数据的数据进行补零操作
				}

			}
		}
		else
		{
			for(int j=0;j<w;j++)
			{
				TD[i*w+j]=0.0;//对于超出原数据的数据进行补零操作
		
			}
		}
	}

}

void DCT_1D(double *&TD,double *&FD,int Len)//DCT变换前后的数据都是实数
{
	double dtmp=0;
	int l=1;
	int lp=0; 

	while(l<Len)
	{
		l=l*2;
		lp++;
	}
	int N=1<<lp;
	complex<double> *xIn;//复数数组的定义
	complex<double> *xOut;

	xIn =new complex<double>[N*2];//为数组分配内存
	xOut=new complex<double>[N*2];

	memset(xIn,0,sizeof(complex<double>)*N*2);//将分配的内存块中的数据的初始值设为0
	memset(xOut,0,sizeof(complex<double>)*N*2);
	
	for(int i=0;i<N;i++)
	xIn[i]=complex<double>(TD[i],0);//将x由实数转化为复数数据类型

	FFT_1D(xIn,xOut,2*N);//调用FFT_1D算出xOut后在对其进行适当的调整，其中点数fft的2倍（r+1）

	dtmp=1/sqrt((double)N);//调整系数

	FD[0]=xOut[0].real()*dtmp;//FD[0]的值进行调节时乘以1/sqrt(N)

	dtmp=dtmp*sqrt(2.0);//FD[n]n>0时乘以sqrt(2/N)

	for(int i=0;i<N;i++)
	{
		FD[i]=(xOut[i].real()*cos(i*PI/(N*2))+xOut[i].imag()*sin(i*PI/(N*2)))*dtmp;//应用到正余弦公式简单推导即可得DCT值
	}
	delete xIn;
	delete xOut;

}

void DCT_2D(double *&TD, double *&FD, int lWidth, int lHeight)//二维DCT变换
{
	
	int w=1;
	int h=1;
	int wp=0;//存储w的2的幂数
	int hp=0;//存储h的2的幂数

//计算刚好大于或等于lWidth，lHeight的2的整数次幂，和相应的幂数///
	while(w<lWidth)
	{
		w=w*2;
		wp++;
	}
	while(h<lHeight)
	{
		h=h*2;
		hp++;
	}
对y方向进行离散余弦变换/

	double* TmpIn = new double[h];
	double* TmpOut = new double[h];//为临时数据分配内存

	// 对y方向进行离散余弦变换
	for (int i = 0; i < w; i++){
		// 抽取数据
		for (int j = 0; j < h; j++)
			TmpIn[j] = TD[j * w + i];
		
		// 一维快速离散余弦变换求纵向的DCT
		DCT_1D(TmpIn, TmpOut,lHeight);

		// 保存变换结果
		for (int j = 0; j < h; j++)
			TD[j * w + i] = TmpOut[j];
	}
	
	// 释放内存
	delete TmpIn;
	delete TmpOut;
///对x方向进行快速离散余弦变换//
	// 分配内存
	TmpIn = new double[w];
	TmpOut = new double[w];
	
	// 对x方向进行快速离散余弦变换
	for (int i = 0; i < h; i++){
		for (int k = 0; k < 1; k++){
			// 抽取数据
			for (int j = 0; j < w; j++)
				TmpIn[j] = TD[i * w + j + k];

			// 一维快速离散余弦变换
			DCT_1D(TmpIn, TmpOut, lWidth);

			// 保存变换结果
			for (int j = 0; j < w; j++)
				FD[i * w + j + k] = TmpOut[j];
		}
	}

	// 释放内存
	delete TmpIn;
	delete TmpOut;
}