X-Git-Url: http://git.kpe.io/?a=blobdiff_plain;f=libctsim%2Ffourier.cpp;h=149b4de234b5fc2c4738295ce04df5c51e792dd6;hb=ee3b6b957bf6786bd0b4c15ce67c22381a8301fc;hp=83922d6c6af0530b89ab91100807ff952fea7f54;hpb=8a7697ce57b56cdc43698cd1241ad98d49f9b5ac;p=ctsim.git

diff --git a/libctsim/fourier.cpp b/libctsim/fourier.cpp
index 83922d6..149b4de 100644
--- a/libctsim/fourier.cpp
+++ b/libctsim/fourier.cpp
@@ -98,3 +98,84 @@ Fourier::shuffleNaturalToFourierOrder (ImageFile& im)
   }
   delete [] pRow;
 }
+
+#ifdef HAVE_FFTW
+void Fourier::shuffleNaturalToFourierOrder (fftw_complex* pVector, const int n)
+{
+  fftw_complex* pTemp = static_cast<fftw_complex*>(fftw_malloc(sizeof(fftw_complex) * n));
+  int i;
+
+  if (isOdd(n)) { // Odd
+    int iHalfN = (n - 1) / 2;
+
+    pTemp[0][0] = pVector[iHalfN][0];
+    pTemp[0][1] = pVector[iHalfN][1];
+    for (i = 0; i < iHalfN; i++) {
+      pTemp[i + 1][0] = pVector[i + 1 + iHalfN][0];
+      pTemp[i + 1][1] = pVector[i + 1 + iHalfN][1];
+    }
+    for (i = 0; i < iHalfN; i++) {
+      pTemp[i + iHalfN + 1][0] = pVector[i][0];
+      pTemp[i + iHalfN + 1][1] = pVector[i][1];
+    }
+  } else {     // Even
+    int iHalfN = n / 2;
+    pTemp[0][0] = pVector[iHalfN][0];
+    pTemp[0][1] = pVector[iHalfN][1];
+    for (i = 0; i < iHalfN - 1; i++) {
+      pTemp[i + 1][0] = pVector[i + iHalfN + 1][0];
+      pTemp[i + 1][1] = pVector[i + iHalfN + 1][1];
+    }
+    for (i = 0; i < iHalfN; i++) {
+      pTemp[i + iHalfN][0] = pVector[i][0];
+      pTemp[i + iHalfN][1] = pVector[i][1];
+    }
+  }
+
+  for (i = 0; i < n; i++) {
+    pVector[i][0] = pTemp[i][0];
+    pVector[i][1] = pTemp[i][1];
+  }
+  fftw_free(pTemp);
+}
+
+void Fourier::shuffleFourierToNaturalOrder (fftw_complex* pVector, const int n)
+{
+  fftw_complex* pTemp = static_cast<fftw_complex*>(fftw_malloc(sizeof(fftw_complex) * n));
+  int i;
+  if (isOdd(n)) { // Odd
+    int iHalfN = (n - 1) / 2;
+
+    pTemp[iHalfN][0] = pVector[0][0];
+    pTemp[iHalfN][1] = pVector[0][1];
+    for (i = 0; i < iHalfN; i++) {
+      pTemp[i + 1 + iHalfN][0] = pVector[i + 1][0];
+      pTemp[i + 1 + iHalfN][1] = pVector[i + 1][1];
+    }
+    for (i = 0; i < iHalfN; i++) {
+      pTemp[i][0] = pVector[i + iHalfN + 1][0];
+      pTemp[i][1] = pVector[i + iHalfN + 1][1];
+    }
+  } else {     // Even
+    int iHalfN = n / 2;
+    pTemp[iHalfN][0] = pVector[0][0];
+    pTemp[iHalfN][1] = pVector[0][1];
+    for (i = 0; i < iHalfN; i++) {
+      pTemp[i][0] = pVector[i + iHalfN][0];
+      pTemp[i][1] = pVector[i + iHalfN][1];
+    }
+    for (i = 0; i < iHalfN - 1; i++) {
+      pTemp[i + iHalfN + 1][0] = pVector[i+1][0];
+      pTemp[i + iHalfN + 1][1] = pVector[i+1][1];
+    }
+  }
+
+  for (i = 0; i < n; i++) {
+    pVector[i][0] = pTemp[i][0];
+    pVector[i][1] = pTemp[i][1];
+  }
+
+  fftw_free(pTemp);
+}
+#endif
+