X-Git-Url: http://git.kpe.io/?a=blobdiff_plain;f=libctsim%2Fprojections.cpp;h=8674647146a8c4f64f160c4323ef7be974b304d9;hb=82a12d49ecd2c09c301e6513b549c358f715e764;hp=64730029714c60dac9fb70dcc9630ff81448e60d;hpb=4f15a69a3150f180bd93fcbe1c87dbaca92a77ad;p=ctsim.git

diff --git a/libctsim/projections.cpp b/libctsim/projections.cpp
index 6473002..8674647 100644
--- a/libctsim/projections.cpp
+++ b/libctsim/projections.cpp
@@ -8,7 +8,7 @@
 **  This is part of the CTSim program
 **  Copyright (c) 1983-2001 Kevin Rosenberg
 **
-**  $Id: projections.cpp,v 1.82 2003/03/23 18:37:42 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
@@ -37,15 +37,15 @@ const int Projections::POLAR_INTERP_BICUBIC = 2;
 
 const char* const Projections::s_aszInterpName[] = 
 {
-  {"nearest"},
-  {"bilinear"},
+  "nearest",
+  "bilinear",
 //  {"bicubic"},
 };
 
 const char* const Projections::s_aszInterpTitle[] = 
 {
-  {"Nearest"},
-  {"Bilinear"},
+  "Nearest",
+  "Bilinear",
 //  {"Bicubic"},
 };
 
@@ -1019,12 +1019,12 @@ Projections::convertFFTPolar (ImageFile& rIF, int iInterpolationID, int iZeropad
 
   int iInterpDet = static_cast<int>(static_cast<double>(sqrt(nx*nx+ny*ny)));
   int iNumInterpDetWithZeros = ProcessSignal::addZeropadFactor (iInterpDet, iZeropad);
-  double dProjScale = iInterpDet / (pProj->viewDiameter() * 0.5);
+  double dProjScale = iInterpDet / (pProj->viewDiameter() * 0.05);
   double dZeropadRatio = static_cast<double>(iNumInterpDetWithZeros) / static_cast<double>(iInterpDet);
 
-  fftw_plan plan = fftw_create_plan (iNumInterpDetWithZeros, FFTW_FORWARD, FFTW_IN_PLACE | FFTW_ESTIMATE | FFTW_USE_WISDOM);
+  fftw_complex* pcIn = static_cast<fftw_complex*> (fftw_malloc (sizeof(fftw_complex) * iNumInterpDetWithZeros));
+  fftw_plan plan = fftw_plan_dft_1d (iNumInterpDetWithZeros, pcIn, pcIn, FFTW_FORWARD, FFTW_ESTIMATE);
 
-  fftw_complex* pcIn = new fftw_complex [iNumInterpDetWithZeros];
   std::complex<double>** ppcDetValue = new std::complex<double>* [pProj->m_nView];
   //double dInterpScale = (pProj->m_nDet-1) / static_cast<double>(iInterpDet-1);
   double dInterpScale = pProj->m_nDet / static_cast<double>(iInterpDet);
@@ -1040,28 +1040,30 @@ Projections::convertFFTPolar (ImageFile& rIF, int iInterpolationID, int iZeropad
     LinearInterpolator<DetectorValue> projInterp (detval, pProj->m_nDet);
     for (int iDet = 0; iDet < iInterpDet; iDet++) {
       double dInterpPos = (m_nDet / 2.) + (iDet - dMidPoint) * dInterpScale;
-      pcIn[iDet].re = projInterp.interpolate (dInterpPos) * dProjScale;
-      pcIn[iDet].im = 0;
+      pcIn[iDet][0] = projInterp.interpolate (dInterpPos) * dProjScale;
+      pcIn[iDet][1] = 0;
     }
 
     Fourier::shuffleFourierToNaturalOrder (pcIn, iInterpDet);
     if (iZerosAdded > 0) {
-       for (int iDet1 = iInterpDet -1; iDet1 >= iMidPoint; iDet1--)
-        pcIn[iDet1+iZerosAdded] = pcIn[iDet1];
+      for (int iDet1 = iInterpDet -1; iDet1 >= iMidPoint; iDet1--) {
+        pcIn[iDet1+iZerosAdded][0] = pcIn[iDet1][0];
+        pcIn[iDet1+iZerosAdded][1] = pcIn[iDet1][1];
+      }
       for (int iDet2 = iMidPoint; iDet2 < iMidPoint + iZerosAdded; iDet2++) 
-        pcIn[iDet2].re = pcIn[iDet2].im = 0;
+        pcIn[iDet2][0] = pcIn[iDet2][1] = 0;
     }
 
-    fftw_one (plan, pcIn, NULL);
+    fftw_execute (plan);
 
     ppcDetValue[iView] = new std::complex<double> [iNumInterpDetWithZeros];
     for (int iD = 0; iD < iNumInterpDetWithZeros; iD++) {
-      ppcDetValue[iView][iD] = std::complex<double> (pcIn[iD].re * dFFTScale, pcIn[iD].im * dFFTScale); 
+      ppcDetValue[iView][iD] = std::complex<double> (pcIn[iD][0] * dFFTScale, pcIn[iD][1] * dFFTScale); 
     }
 
     Fourier::shuffleFourierToNaturalOrder (ppcDetValue[iView], iNumInterpDetWithZeros);
   }
-  delete [] pcIn;
+  fftw_free(pcIn) ;
 
   fftw_destroy_plan (plan);