c++ DFT 二维傅里叶变换

最近写了一个计算图像二维傅里叶变换的程序， c++， 采用索引表的方法计算。

程序如下

#include <opencv2/core/core.hpp>
#include <opencv2/highgui/highgui.hpp>
#include <opencv2/imgproc/imgproc.hpp>
using namespace cv;

#include <stdlib.h>
#include <stdio.h>

#define PI 3.14159265359

float *trigalesin;
float *trigalecos;

void calculate(int col, int row)
{
trigalesin  = new float[col*row];
trigalecos = new float[col*row];

for (int y = 0; y< row; y++)
{
for(int x = 0; x< col; x++)
{
float tmp = -2*PI*(1.0*x/col+1.0*y/row);

trigalesin[y*col+x] = sin(tmp);
trigalecos[y*col+x] = cos(tmp);
}
}
}

void getDFT(int *a, int arows, int acols, float **b, float **c)
{


*b = new float [arows*acols];
*c = new float [arows*acols];

for (int v = 0; v< arows; v++)
{
if (v*100%arows == 0)
{
printf("%d\n",v*100/arows );
}
for (int u = 0; u< acols; u++)
{
float tmpb = 0;
float tmpc = 0;



for (int y = 0; y< arows; y++)
{
for (int x = 0; x< acols; x++)
{
float t1 = a[y*acols+x];

if (t1 == 0)
{
continue;
}

// float t2 = -2*PI*(u*x/acols+v*y/arows);

int index1 = (u*x)%acols;
int index2 = (v*y)%arows;
int index = index2*acols+index1;


tmpb += t1*trigalesin[index];
tmpc += t1*trigalecos[index];
}
}

(*b)[v*acols+u] = tmpb;
(*c)[v*acols+u] = tmpc;
}
}
}

void getIDFT(float *a, float *b, float **c, int arows, int acols)
{
*c = new float [arows*acols];

for (int y = 0; y< arows; y++)
{
if (y*100%arows == 0)
{
printf("%d\n",y*100/arows );
}

for (int x = 0; x< acols; x++)
{
float tmpa = 0;
float tmpb = 0;

for (int v = 0; v< arows; v++)
{
for (int u = 0; u< acols; u++)
{
float tmpfa = a[v*acols+u];
float tmpfb = b[v*acols+u];

int index1 = (u*x)%acols;
int index2 = (v*y)%arows;
int index = index2*acols+index1;

float tmpea = -trigalesin[index];//sin(2*PI*(u*x/acols+v*y/arows));
float tmpeb = trigalecos[index];//cos(2*PI*(u*x/acols+v*y/arows));

tmpa += tmpfa*tmpea - tmpfb*tmpeb;
tmpb += tmpfa*tmpeb+tmpfb*tmpea;
}
}

(*c)[y*acols+x] = sqrt(tmpa*tmpa+tmpb*tmpb);
}
}
}

void getPadding(Mat a, Mat b, int **adata, int **bdata, int &rows, int &cols)
{
rows = pow(2, 1.0*(int)(log(b.rows*1.0)/log(2.0)+1));
cols = pow(2, 1.0*(int)(log(b.cols*1.0)/log(2.0)+1));

rows = b.rows;
cols = b.cols;

*adata = new int[rows*cols];
*bdata =  new int[rows*cols];



for (int i = 0; i< rows; i++)
{
for (int j = 0; j< cols; j++)
{
(*adata)[i*cols+j] = 0;
(*bdata)[i*cols+j] = 0;
}
}

for (int i = 0; i< a.rows; i++)
{
for (int j = 0; j< a.cols; j++)
{
(*adata)[i*cols+j] = a.at<uchar>(i, j);
}
}

for (int i = 0; i< b.rows; i++)
{
for (int j = 0; j< b.cols; j++)
{
(*bdata)[i*cols+j] = b.at<uchar>(i, j);
}
}
}

void getMul(float *ab, float *ac, float *bb, float *bc, float **cb, float **cc, int rows, int cols)
{
*cb = new float [rows*cols];
*cc = new float [rows*cols];

for (int i = 0; i< rows; i++)
{
for (int j = 0; j< cols; j++)
{
(*cb)[i*cols+j] = ab[i*cols+j]*bb[i*cols+j] - ac[i*cols+j]*bc[i*cols+j];
(*cc)[i*cols+j] = ab[i*cols+j]*bc[i*cols+j]+ac[i*cols+j]*bb[i*cols+j];
}
}
}

over.

另外， dft做出来了， 但是如何计算图像位置还不清楚。