dct变换源代码

#define EXTRA_NAME "@dcttransform."

#include "loadbmp.h"

#define Point(x,y) lpPoints[(x)+(y)*nWidth]

void FFT(COMPLEX * TD, COMPLEX * FD, int power)

{

int count;

int i,j,k,bfsize,p;

double angle;

COMPLEX *W,*X1,*X2,*X;

count=1<<power;

W=(COMPLEX *)malloc(sizeof(COMPLEX)*count/2);

X1=(COMPLEX *)malloc(sizeof(COMPLEX)*count);

X2=(COMPLEX *)malloc(sizeof(COMPLEX)*count);

for(i=0;i<count/2;i++)

{

angle=-i*pi*2/count;

W[i].re=cos(angle);

W[i].im=sin(angle);

}

memcpy(X1,TD,sizeof(COMPLEX)*count);

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

{

for(j=0;j<1<<k;j++)

{

bfsize=1<<(power-k);

for(i=0;i<bfsize/2;i++)

{

p=j*bfsize;

X2[i+p]=Add(X1[i+p],X1[i+p+bfsize/2]);

X2[i+p+bfsize/2]=Mul(Sub(X1[i+p],X1[i+p+bfsize/2]),W[i*(1<<k)]);

}

}

X=X1;

X1=X2;

X2=X;

}

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

{

p=0;

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

{

if (j&(1<<i)) p+=1<<(power-i-1);

}

FD[j]=X1[p];

}

free(W);

free(X1);

free(X2);

}

void DCT_Pre(double *f, double *F, int power)

{

int i,count;

COMPLEX *X;

double s;

count=1<<power;

X=(COMPLEX *)malloc(sizeof(COMPLEX)*count*2);

memset(X,0,sizeof(COMPLEX)*count*2);

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

{

X[i].re=f[i];

}

FFT(X,X,power+1);

s=1/sqrt(count);

F[0]=X[0].re*s;

s*=sqrt(2);

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

{

F[i]=(X[i].re*cos(i*pi/(count*2))+X[i].im*sin(i*pi/(count*2)))*s;

}

free(X);

}

void Dct()

{

int w=1,h=1,wp=0,hp=0;

while(w*2<=nWidth)

{

w*=2;

wp++;

}

while(h*2<=nHeight)

{

h*=2;

hp++;

}

int x,y;

BYTE *lpPoints=new BYTE[nWidth*nHeight];

GetPoints(lpPoints);

double *f=new double[w*h];

double *W=new double[w*h];

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

{

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

{

f[x+y*w]=Point(x,y);

}

}

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

{

DCT_Pre(&f[w*y],&W[w*y],wp);

}

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

{

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

{

f[x*h+y]=W[x+w*y];

}

}

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

{

DCT_Pre(&f[x*h],&W[x*h],hp);

}

double a;

memset(lpPoints,0,nWidth*nHeight);

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

{

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

{

a=fabs(W[x*h+y]);

if (a>255) a=255;

Point(x,nHeight-y-1)=(BYTE)(a);

}

}

delete f;

delete W;

PutPoints(lpPoints);

delete lpPoints;

}

void main(int argc, char *argv[])

{

if(argc==2)

FileName=argv[1];

else

return;

OpenFile();

Dct();

SaveAs();

}