P
Ich habe mir für FFTW einen Wrapper gebastelt der inzwischen auch sehr schön FFTs rechnet. Allerdings habe ich damit zwei Probleme.
void CFFTW::FFT2D(CComplex* InArray, CComplex* OutArray)
{
if (m_Init != true)
return;
for (int ix=0; ix < m_size.x; ix++)
{
for (int iy=0; iy < m_size.y; iy++)
{
in[ArrPos(ix,iy)][0] = InArray[ArrPos(ix,iy)].re;
in[ArrPos(ix,iy)][1] = InArray[ArrPos(ix,iy)].im;
}
}
fftw_execute(plan);
for (int ix=0; ix < m_size.x; ix++)
{
for (int iy=0; iy < m_size.y; iy++)
{
OutArray[ArrPos(ix,iy)].re = out[ArrPos(ix,iy)][0];
OutArray[ArrPos(ix,iy)].im = out[ArrPos(ix,iy)][1];
}
}
Cycle2D(OutArray, m_size.x);
}
Zum einen habe ich eine Klasse CComplex die nicht gleich std::Complex ist. Daher kopiere ich im Moment die Arrays mit den for Schleifen um.
1. Kann man das effizienter machen
2. Wenn ich CComplex von std::Complex ableiten würde könnte ich dann das Array direkt verwenden obwohl es den Typ CComplex hat ?
Zum anderen wird das Ergebnis der FFT verschoben dargestellt, also statt
1 2
3 4
sehe ich
4 3
2 1
das ändere ich mit Cycle2D - auch hier wäre die Frage ob man das schneller bewerkstelligen könnte, bzw warum die FFT das so ausgibt.
// *********************************
// Cycle changes 012 345 to 345 012
// *********************************
void CFFTW::Cycle(CComplex *A1,long N)
{
Cycle(A1,N,N/2);
}
void CFFTW::Cycle(CComplex *A1,long N,long del)
{
long i;
CComplex *D;
D=new CComplex [N];
if (del>=N) del-=N;
if (del<=-N) del+=N;
if(del>0 && del<N)
{
for (i=N-del;i<N;i++) D[i-N+del]=A1[i];
for (i=N-1;i>=del;i--) A1[i]=A1[i-del];
for (i=0;i<del;i++) A1[i]=D[i];
}
else if(del<0 && del > -N)
{
del=-del;
for (i=0;i<del;i++) D[i]=A1[i];
for (i=0;i<N-del;i++) A1[i]=A1[i+del];
for (i=N-del;i<N;i++) A1[i]=D[i-N+del];
}
delete [] D;
}
// **************************************
// 1 2 4 3
// 3 4 -> 2 1
// **************************************
void CFFTW::Cycle2D(CComplex *A,long N)
{
CComplex swap;
long ix,iy;
for (ix=0;ix<N;ix++) Cycle(A+ix*N,N); // Spalten shiften
for (iy=0;iy<N;iy++) // Zeilen einzeln shiften
for (ix=0;ix<N/2;ix++)
{
swap=A[(ix+N/2)*N+iy];
A[(ix+N/2)*N+iy]=A[ix*N+iy];
A[ix*N+iy]=swap;
}
}
Matthias
