X-Git-Url: http://git.kpe.io/?p=ctsim.git;a=blobdiff_plain;f=libctsim%2Ffourier.cpp;h=149b4de234b5fc2c4738295ce04df5c51e792dd6;hp=6fec9f5d5d95cb3756ab7196a8c3ba77f85ce960;hb=1a050c98763fbbc0662731b0b76953acede6f5d7;hpb=3ea498d51ce4597e9649cd21f155b51175ea0bea

diff --git a/libctsim/fourier.cpp b/libctsim/fourier.cpp
index 6fec9f5..149b4de 100644
--- a/libctsim/fourier.cpp
+++ b/libctsim/fourier.cpp
@@ -9,7 +9,7 @@
 **  This is part of the CTSim program
 **  Copyright (c) 1983-2001 Kevin Rosenberg
 **
-**  $Id: fourier.cpp,v 1.6 2001/03/21 21:45:31 kevin Exp $
+**  $Id$
 **
 **  This program is free software; you can redistribute it and/or modify
 **  it under the terms of the GNU General Public License (version 2) as
@@ -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
+