Last active
November 3, 2015 17:43
-
-
Save domitry/495850648ce3d480f294 to your computer and use it in GitHub Desktop.
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
#include <stdio.h> | |
#include <stdlib.h> | |
#include <math.h> | |
#include <GL/glut.h> | |
#include "image.h" | |
#define PI 3.1415926 | |
#define P(x,y) x+y*WD | |
void dft(double *A, double *re, double *im){ | |
double x_re[WD*HT], x_im[WD*HT]; | |
for(int v=0;v<HT;v++) | |
for(int u=0;u<WD;u++){ | |
x_re[P(u,v)]=0.0; | |
x_im[P(u,v)]=0.0; | |
for(int x=0;x<WD;x++){ | |
x_re[P(u,v)]+=A[P(x,v)]*cos((-2.0*PI*u*x)/WD); | |
x_im[P(u,v)]+=A[P(x,v)]*sin((-2.0*PI*u*x)/WD); | |
} | |
} | |
for(int u=0;u<WD;u++) | |
for(int v=0;v<HT;v++){ | |
re[P(u,v)]=0.0; | |
im[P(u,v)]=0.0; | |
for(int y=0;y<HT;y++){ | |
re[P(u,v)]+=x_re[P(u,y)]*cos((-2.0*PI*v*y)/HT)-x_im[P(u,y)]*sin((-2.0*PI*v*y)/HT); | |
im[P(u,v)]+=x_re[P(u,y)]*sin((-2.0*PI*v*y)/HT)+x_im[P(u,y)]*cos((-2.0*PI*v*y)/HT); | |
} | |
} | |
} | |
void idft(double *re, double *im, double *A){ | |
double x_re[WD*HT], x_im[WD*HT]; | |
for(int v=0;v<HT;v++){ | |
for(int u=0;u<WD;u++){ | |
x_re[P(u,v)]=0.0; | |
x_im[P(u,v)]=0.0; | |
for(int x=0;x<WD;x++){ | |
x_re[P(u,v)]+=re[P(x,v)]*cos((2.0*PI*u*x)/WD)-im[P(x,v)]*sin((2.0*PI*u*x)/WD); | |
x_im[P(u,v)]+=re[P(x,v)]*sin((2.0*PI*u*x)/WD)+im[P(x,v)]*cos((2.0*PI*u*x)/WD); | |
} | |
} | |
} | |
for(int u=0;u<WD;u++){ | |
for(int v=0;v<WD;v++){ | |
A[P(u,v)]=0.0; | |
for(int y=0;y<HT;y++){ | |
A[P(u,v)]+=x_re[P(u,y)]*cos((2.0*PI*v*y)/HT)-x_im[P(u,y)]*sin((2.0*PI*v*y)/HT); | |
} | |
} | |
} | |
for(int x=0;x<WD;x++) | |
for(int y=0;y<HT;y++) | |
A[P(x,y)]/=(WD*HT); | |
} | |
double max_double(double *buff, int size){ | |
double ret=-1e30; | |
int i; | |
for(i=0;i<size;i++)ret=(ret < buff[i] ? buff[i] : ret); | |
return ret; | |
} | |
double min_double(double *buff, int size){ | |
double ret=1e30; | |
int i; | |
for(i=0;i<size;i++)ret=(ret > buff[i] ? buff[i] : ret); | |
return ret; | |
} | |
void flip(double *arr){ | |
int y; | |
for(y=0; y<HT/2; y++){ | |
int x, src=y, dst=(HT-y-1); | |
double tmp[WD]; | |
for(x=0;x<WD;x++)tmp[x]=arr[src*WD+x]; | |
for(x=0;x<WD;x++)arr[src*WD+x]=arr[dst*WD+x]; | |
for(x=0;x<WD;x++)arr[dst*WD+x]=tmp[x]; | |
} | |
} | |
int main(int argc, char **argv){ | |
double re[WD*HT], im[WD*HT], tempbuf[WD*HT], ABS[WD*HT]; | |
read_image("sample5.ras"); printf("read\n"); | |
for(int i=0;i<WD*HT;i++)tempbuf[i]=buf[i]; | |
flip(tempbuf); | |
dft(tempbuf,re,im); | |
idft(re,im,tempbuf); | |
for(int i=0;i<WD*HT;i++)buf[i]=(unsigned char)tempbuf[i]; | |
/*// debug from here // | |
for(i=0; i<WD*HT; i++){ | |
ABS[i] = log(sqrt(pow(re[i], 2) + pow(im[i], 2))); | |
} | |
for(i=0;i<WD*HT;i++){ | |
double max_ = max_double(ABS, WD*HT); | |
double min_ = min_double(ABS, WD*HT); | |
buf[i]=(unsigned char)((511*(ABS[i]-min_))/(max_-min_)); | |
} | |
// debug upto here //*/ | |
write_image("output.ras"); printf("write\n"); | |
glutInit(&argc, argv); | |
glutInitDisplayMode(GLUT_DOUBLE|GLUT_RGB|GLUT_DEPTH); | |
glutInitWindowSize(WD, HT); | |
glutInitWindowPosition(100, 100); | |
glutCreateWindow(argv[0]); | |
init(); | |
glutDisplayFunc(display); | |
glutReshapeFunc(reshape); | |
glutKeyboardFunc(keyboard); | |
glutMainLoop(); | |
return 0; | |
} |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment