#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<?W=(COMPLEX *)malloc(sizeof(COMPLEX)*count/2);
?X1=(COMPLEX *)malloc(sizeof(COMPLEX)*count);
?X2=(COMPLEX *)malloc(sizeof(COMPLEX)*count);
?for(i=0;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?{
??for(j=0;j<1<??{
???bfsize=1<<(power-k);
???for(i=0;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<???}
??}
??X=X1;
??X1=X2;
??X2=X;
?}
?for(j=0;j?{
??p=0;
??for(i=0;i??{
???if (j&(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<?X=(COMPLEX *)malloc(sizeof(COMPLEX)*count*2);
?memset(X,0,sizeof(COMPLEX)*count*2);
?for(i=0;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?{
??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?{
??for(x=0;x??{
???f[x+y*w]=Point(x,y);
??}
?}

?for(y=0;y?{
??DCT_Pre(&f[w*y],&W[w*y],wp);
?}
?for(y=0;y?{
??for(x=0;x??{
???f[x*h+y]=W[x+w*y];
??}
?}
?for(x=0;x?{
??DCT_Pre(&f[x*h],&W[x*h],hp);
?}
?double a;
?memset(lpPoints,0,nWidth*nHeight);

?for(y=0;y?{
??for(x=0;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();
}