r636: Optimized Rebinning, Added Reconstruct with Rebinning option

diff --git a/ChangeLog b/ChangeLog
index d723fa9bd649f8065c03cadd5db58ca417935fbf..a9677d80ea8eb54a96f08e25d65238ea071b094d 100644 (file)
--- a/ChangeLog
+++ b/ChangeLog
@@ -21,6 +21,9 @@
        * ctsim: Added imagefile export to DICOM files.
 
        * ctsim: Added clipboard cut/copy/paste for image files.
+
+       * ctsim: Added Reconstruction with Rebinning for faster
+       divergent beam reconstructions.
        
        * ctsim: Added plot t-theta sampling to projection file menu.
        

diff --git a/include/projections.h b/include/projections.h
index b86f448a9edd091f7c45248c5d4b355c1c6b9a6b..c88fca709552c694e375ced81711486e91732df6 100644 (file)
--- a/include/projections.h
+++ b/include/projections.h
@@ -9,7 +9,7 @@
 **  This is part of the CTSim program
 **  Copyright (c) 1983-2001 Kevin Rosenberg
 **
-**  $Id: projections.h,v 1.32 2001/03/13 04:44:25 kevin Exp $
+**  $Id: projections.h,v 1.33 2001/03/13 08:24:41 kevin Exp $
 **
 **
 **  This program is free software; you can redistribute it and/or modify
@@ -65,7 +65,7 @@ public:
     THETA_RANGE_FOLD_TO_PI,
   };
 
-  ParallelRaysums (Projections* pProjections, int iThetaRange);
+  ParallelRaysums (const Projections* pProjections, int iThetaRange);
   ~ParallelRaysums ();
 
   typedef std::vector<ParallelRaysumCoordinate*> CoordinateContainer;
@@ -84,6 +84,7 @@ private:
   CoordinateContainer m_vecpCoordinates;
   CoordinateContainer m_vecpSortedByT;
   CoordinateContainer m_vecpSortedByTheta;
+  ParallelRaysumCoordinate* m_pCoordinates;
   int m_iNumCoordinates;
   int m_iNumView;
   int m_iNumDet;
@@ -128,7 +129,7 @@ class Projections
   bool detarrayRead (fnetorderstream& fs, DetectorArray& darray, const int view_num);
   bool detarrayWrite (fnetorderstream& fs, const DetectorArray& darray, const int view_num);
 
-  Projections* interpolateToParallel();
+  Projections* interpolateToParallel() const;
 
   bool convertPolar (ImageFile& rIF, int iInterpolation);
   bool convertFFTPolar (ImageFile& rIF, int iInterpolation, int iZeropad);

diff --git a/include/reconstruct.h b/include/reconstruct.h
index 61df5257df5f38d542fb19b2b39f418cfacacdd5..473fa66d6701c16b918a604e8faef51f6975cea6 100644 (file)
--- a/include/reconstruct.h
+++ b/include/reconstruct.h
@@ -8,7 +8,7 @@
 **  This is part of the CTSim program
 **  Copyright (c) 1983-2001 Kevin Rosenberg
 **
-**  $Id: reconstruct.h,v 1.8 2001/03/11 15:27:30 kevin Exp $
+**  $Id: reconstruct.h,v 1.9 2001/03/13 08:24:41 kevin Exp $
 **
 **  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
@@ -45,13 +45,14 @@ struct ReconstructionROI {
 class Reconstructor 
 {
  private:
-    const Projections& m_rProj;
+    const Projections& m_rOriginalProj;
+    const Projections* m_pProj;
     ImageFile& m_rImagefile;
     ProcessSignal* m_pProcessSignal;
     Backprojector* m_pBackprojector;
     int m_nFilteredProjections;
     int m_iTrace;
-
+    const bool m_bRebinToParallel;
     bool m_bFail;
     std::string m_strFailMessage;
 
@@ -61,7 +62,7 @@ class Reconstructor
     Reconstructor (const Projections& rProj, ImageFile& rIF, const char* const filterName, double filt_param, 
       const char* const filterMethodName, const int zeropad, const char* filterGenerationName, 
       const char* const interpName, int interpFactor, const char* const backprojectName, const int trace, 
-      ReconstructionROI* pROI = NULL, SGP* pSGP = NULL);
+      ReconstructionROI* pROI = NULL, bool bRebinToParallel = false, SGP* pSGP = NULL);
 
     ~Reconstructor ();
 

diff --git a/libctsim/projections.cpp b/libctsim/projections.cpp
index 60f9467e06cbb2426ed2db355b9c163248ffdbd7..3eab99c4fe9a8a8725bc2108492eb86302e50a64 100644 (file)
--- 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.62 2001/03/13 04:44:25 kevin Exp $
+**  $Id: projections.cpp,v 1.63 2001/03/13 08:24:41 kevin Exp $
 **
 **  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
@@ -895,59 +895,53 @@ Projections::initFromSomatomAR_STAR (int iNViews, int iNDets, unsigned char* pDa
   m_dFocalLength = 510;
   m_dSourceDetectorLength = 890;
   m_detInc = convertDegreesToRadians (3.06976 / 60);
+  m_dFanBeamAngle = (iNDets + 1) * m_detInc;
   m_detStart = -(m_dFanBeamAngle / 2);
   m_rotInc = TWOPI / static_cast<double>(iNViews);
   m_rotStart = HALFPI;
-  m_dFanBeamAngle = (iNDets + 1) * m_detInc;
   m_dViewDiameter = sin (m_dFanBeamAngle / 2) * m_dFocalLength * 2;
 
   if (! ((iNViews == 750 && lDataLength == 1560000L) || (iNViews == 950 && lDataLength == 1976000L) 
                 || (iNViews == 1500 && lDataLength == 3120000)))
     return false;
 
-  int iCenter = (iNDets - 1) / 2; // change from (Nm+1)/2 because of 0 vs. 1 indexing
+  double dCenter = (iNDets - 1.) / 2.; // change from (Nm+1)/2 because of 0 vs. 1 indexing
   double* pdCosScale = new double [iNDets];
   for (int i = 0; i < iNDets; i++)
-    pdCosScale[i] = cos ((i - iCenter) * m_detInc);
+    pdCosScale[i] = 1. / cos ((i - dCenter) * m_detInc);
 
   long lDataPos = 0;
   for (int iv = 0; iv < iNViews; iv++) {
     unsigned char* pArgBase = pData + lDataPos;
-    unsigned char* p = pArgBase+0;
-    SwapBytes4IfLittleEndian (p);
+    unsigned char* p = pArgBase+0; SwapBytes4IfLittleEndian (p);
     long lProjNumber = *reinterpret_cast<long*>(p);
 
-    p = pArgBase+20;
-    SwapBytes4IfLittleEndian (p);
+    p = pArgBase+20;  SwapBytes4IfLittleEndian (p);
     long lEscale = *reinterpret_cast<long*>(p);
 
-    p = pArgBase+28;
-    SwapBytes4IfLittleEndian (p);
+    p = pArgBase+28;  SwapBytes4IfLittleEndian (p);
     long lTime = *reinterpret_cast<long*>(p);
 
-    p = pArgBase + 4;
-    SwapBytes4IfLittleEndian (p);
+    p = pArgBase + 4; SwapBytes4IfLittleEndian (p);
     double dAlpha = *reinterpret_cast<float*>(p) + HALFPI;
 
-    p = pArgBase+12;
-    SwapBytes4IfLittleEndian (p);
+    p = pArgBase+12; SwapBytes4IfLittleEndian (p);
     double dAlign = *reinterpret_cast<float*>(p);
 
-    p = pArgBase + 16;
-    SwapBytes4IfLittleEndian (p);
+    p = pArgBase + 16; SwapBytes4IfLittleEndian (p);
     double dMaxValue = *reinterpret_cast<float*>(p);
 
     DetectorArray& detArray = getDetectorArray (iv);
     detArray.setViewAngle (dAlpha);
     DetectorValue* detval = detArray.detValues();
 
-    double dViewScale = 2294.4871 * ::pow (2.0, -lEscale);
+    double dViewScale = 1. / (2294.4871 * ::pow (2.0, -lEscale));
     lDataPos += 32;
     for (int id = 0; id < iNDets; id++) {
-      int iV = pData[lDataPos+1] + 256 * pData[lDataPos];
+      int iV = pData[lDataPos+1] + (pData[lDataPos] << 8);
       if (iV > 32767)   // two's complement signed conversion
         iV = iV - 65536;
-      detval[id] = iV / (dViewScale * pdCosScale[id]);
+      detval[id] = iV * dViewScale * pdCosScale[id];
       lDataPos += 2;
     }
   }
@@ -957,10 +951,10 @@ Projections::initFromSomatomAR_STAR (int iNViews, int iNDets, unsigned char* pDa
 }
 
 Projections*
-Projections::interpolateToParallel ()
+Projections::interpolateToParallel () const
 {
   if (m_geometry == Scanner::GEOMETRY_PARALLEL)
-    return this;
+    return const_cast<Projections*>(this);
 
   int nDet = m_nDet;
   int nView = m_nView;
@@ -1045,9 +1039,9 @@ Projections::interpolateToParallel ()
 //
 ///////////////////////////////////////////////////////////////////////////////
 
-ParallelRaysums::ParallelRaysums (Projections* pProjections, int iThetaRange)
+ParallelRaysums::ParallelRaysums (const Projections* pProjections, int iThetaRange)
 : m_iNumCoordinates(0), m_iNumView(pProjections->nView()), m_iNumDet(pProjections->nDet()),
-  m_iThetaRange (iThetaRange)
+  m_iThetaRange (iThetaRange), m_pCoordinates(NULL)
 {
   int iGeometry = pProjections->geometry();
   double dDetInc = pProjections->detInc();
@@ -1055,9 +1049,10 @@ ParallelRaysums::ParallelRaysums (Projections* pProjections, int iThetaRange)
   double dFocalLength = pProjections->focalLength();
 
   m_iNumCoordinates =  m_iNumView * m_iNumDet;
+  m_pCoordinates = new ParallelRaysumCoordinate [m_iNumCoordinates];
   m_vecpCoordinates.reserve (m_iNumCoordinates);
   for (int i = 0; i < m_iNumCoordinates; i++)
-    m_vecpCoordinates[i] = new ParallelRaysumCoordinate;
+    m_vecpCoordinates[i] = m_pCoordinates + i;
 
   int iCoordinate = 0;
   for (int iV = 0; iV < m_iNumView; iV++) {
@@ -1096,8 +1091,7 @@ ParallelRaysums::ParallelRaysums (Projections* pProjections, int iThetaRange)
 
 ParallelRaysums::~ParallelRaysums()
 {
-  for (int i = 0; i < m_iNumCoordinates; i++)
-    delete m_vecpCoordinates[i];
+  delete m_pCoordinates;
 }
 
 ParallelRaysums::CoordinateContainer&

diff --git a/libctsim/reconstruct.cpp b/libctsim/reconstruct.cpp
index c3b6e8b65eb6777cace756dfe7c95544d52a9c70..8caed206c2df46997d38eb77b3c3ef674b644db3 100644 (file)
--- a/libctsim/reconstruct.cpp
+++ b/libctsim/reconstruct.cpp
@@ -8,7 +8,7 @@
 **  This is part of the CTSim program
 **  Copyright (c) 1983-2001 Kevin Rosenberg
 **
-**  $Id: reconstruct.cpp,v 1.16 2001/03/11 15:27:30 kevin Exp $
+**  $Id: reconstruct.cpp,v 1.17 2001/03/13 08:24:41 kevin Exp $
 **
 **  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
@@ -53,11 +53,13 @@ Reconstructor::Reconstructor (const Projections& rProj, ImageFile& rIF, const ch
                               double filt_param, const char* const filterMethodName, const int zeropad, 
                               const char* filterGenerationName, const char* const interpName, 
                               int interpFactor, const char* const backprojectName, const int iTrace, 
-                              ReconstructionROI* pROI, SGP* pSGP)
-  : m_rProj(rProj), m_rImagefile(rIF), m_pProcessSignal(0), m_pBackprojector(0), m_iTrace(iTrace), 
-    m_bFail(false), m_adPlotXAxis(0)
+                              ReconstructionROI* pROI, bool bRebinToParallel, SGP* pSGP)
+  : m_rOriginalProj(rProj), 
+    m_pProj(bRebinToParallel ? m_rOriginalProj.interpolateToParallel() : &m_rOriginalProj),
+    m_rImagefile(rIF), m_pProcessSignal(0), m_pBackprojector(0), 
+    m_iTrace(iTrace), m_bRebinToParallel(bRebinToParallel), m_bFail(false), m_adPlotXAxis(0)
 {
-  m_nFilteredProjections = m_rProj.nDet() * interpFactor;
+  m_nFilteredProjections = m_pProj->nDet() * interpFactor;
 
 #ifdef HAVE_BSPLINE_INTERP
   int spline_order = 0, zoom_factor = 0;
@@ -69,10 +71,10 @@ Reconstructor::Reconstructor (const Projections& rProj, ImageFile& rIF, const ch
   }
 #endif
 
-  double filterBW = 1. / m_rProj.detInc();
-  m_pProcessSignal = new ProcessSignal (filterName, filterMethodName, filterBW, m_rProj.detInc(), 
-    m_rProj.nDet(), filt_param, "spatial", filterGenerationName, zeropad, interpFactor, iTrace, 
-    m_rProj.geometry(), m_rProj.focalLength(), m_rProj.sourceDetectorLength(), pSGP);
+  double filterBW = 1. / m_pProj->detInc();
+  m_pProcessSignal = new ProcessSignal (filterName, filterMethodName, filterBW, m_pProj->detInc(), 
+    m_pProj->nDet(), filt_param, "spatial", filterGenerationName, zeropad, interpFactor, iTrace, 
+    m_pProj->geometry(), m_pProj->focalLength(), m_pProj->sourceDetectorLength(), pSGP);
 
   if (m_pProcessSignal->fail()) {
     m_bFail = true;
@@ -82,7 +84,7 @@ Reconstructor::Reconstructor (const Projections& rProj, ImageFile& rIF, const ch
     return;
   }
 
-  m_pBackprojector = new Backprojector (m_rProj, m_rImagefile, backprojectName, interpName, interpFactor, pROI);
+  m_pBackprojector = new Backprojector (*m_pProj, m_rImagefile, backprojectName, interpName, interpFactor, pROI);
   if (m_pBackprojector->fail()) {
     m_bFail = true;
     m_strFailMessage = "Error creating backprojector: ";
@@ -93,17 +95,20 @@ Reconstructor::Reconstructor (const Projections& rProj, ImageFile& rIF, const ch
   }
 
 #ifdef HAVE_SGP
-  m_adPlotXAxis = new double [m_rProj.nDet()];
-  double x = - ((m_rProj.nDet() - 1) / 2) * m_rProj.detInc();
-  double xInc = m_rProj.detInc();
+  m_adPlotXAxis = new double [m_pProj->nDet()];
+  double x = - ((m_pProj->nDet() - 1) / 2) * m_pProj->detInc();
+  double xInc = m_pProj->detInc();
 
-  for (int i = 0; i < m_rProj.nDet(); i++, x += xInc)
+  for (int i = 0; i < m_pProj->nDet(); i++, x += xInc)
     m_adPlotXAxis[i] = x;
 #endif
 }
 
 Reconstructor::~Reconstructor ()
 {
+  if (m_bRebinToParallel)
+    delete m_pProj;
+
   delete m_pBackprojector;
   delete m_pProcessSignal;
   delete m_adPlotXAxis;
@@ -139,7 +144,7 @@ Reconstructor::plotFilter (SGP* pSGP)
 void
 Reconstructor::reconstructAllViews ()
 {
-  reconstructView (0, m_rProj.nView());
+  reconstructView (0, m_pProj->nView());
   postProcessing();
 }
 
@@ -156,24 +161,24 @@ Reconstructor::reconstructView (int iStartView, int iViewCount, SGP* pSGP, bool
   double* adFilteredProj = new double [m_nFilteredProjections];   // filtered projections
 
   if (iViewCount <= 0)
-    iViewCount = m_rProj.nView() - iStartView;
+    iViewCount = m_pProj->nView() - iStartView;
       
   for (int iView = iStartView; iView < (iStartView + iViewCount); iView++)  {
     if (m_iTrace == Trace::TRACE_CONSOLE) 
-               std::cout <<"Reconstructing view " << iView << " (last = " << m_rProj.nView() - 1 << ")\n";
+               std::cout <<"Reconstructing view " << iView << " (last = " << m_pProj->nView() - 1 << ")\n";
       
-    const DetectorArray& rDetArray = m_rProj.getDetectorArray (iView);
+    const DetectorArray& rDetArray = m_pProj->getDetectorArray (iView);
     const DetectorValue* detval = rDetArray.detValues();
 
     m_pProcessSignal->filterSignal (detval, adFilteredProj);
 
 #ifdef HAVE_BSPLINE_INTERP
     if (interp_type == I_BSPLINE) 
-       bspline (m_rProj.nDet(), zoom_factor, spline_order, adFilteredProj, adFilteredProj);
+       bspline (m_pProj->nDet(), zoom_factor, spline_order, adFilteredProj, adFilteredProj);
     
 #ifdef HAVE_SGP
     if (trace >= Trace::TRACE_PLOT && interp_type == I_BSPLINE && pSGP) {
-       bspline (m_rProj.nDet(), zoom_factor, spline_order, adFilteredProj, adFilteredProj);
+       bspline (m_pProj->nDet(), zoom_factor, spline_order, adFilteredProj, adFilteredProj);
       ezplot_1d (adFilteredProj, m_nFilteredProjections);
     }
 #endif
@@ -201,14 +206,14 @@ Reconstructor::reconstructView (int iStartView, int iViewCount, SGP* pSGP, bool
       ezplotProj.ezset ("box.");
       ezplotProj.ezset ("grid.");
 #if 0  // workaround c++ optimizer bug, now disabled by using /O1 in code
-      double* pdDetval = new double [m_rProj.nDet()];
-      for (unsigned int id = 0; id < m_rProj.nDet(); id++) {
+      double* pdDetval = new double [m_pProj->nDet()];
+      for (unsigned int id = 0; id < m_pProj->nDet(); id++) {
         pdDetval[id] = detval[id];
       }
-      ezplotProj.addCurve (m_adPlotXAxis, pdDetval, m_rProj.nDet());
+      ezplotProj.addCurve (m_adPlotXAxis, pdDetval, m_pProj->nDet());
       delete pdDetval;
 #else
-      ezplotProj.addCurve (m_adPlotXAxis, detval, m_rProj.nDet());
+      ezplotProj.addCurve (m_adPlotXAxis, detval, m_pProj->nDet());
 #endif
       pSGP->setTextPointSize (12);
       ezplotProj.plot (pSGP);

diff --git a/msvc/ctsim/ctsim.plg b/msvc/ctsim/ctsim.plg
index ba431eb85ec37e7b8d19808ed6aaf65994979c7d..5c51373d3d61b0874af767712bdb4b0b43786867 100644 (file)
--- a/msvc/ctsim/ctsim.plg
+++ b/msvc/ctsim/ctsim.plg
@@ -6,19 +6,13 @@
 --------------------Configuration: ctsim - Win32 Debug--------------------
 </h3>
 <h3>Command Lines</h3>
-Creating temporary file "C:\DOCUME~1\kevin\LOCALS~1\Temp\RSPD2.tmp" with contents
+Creating temporary file "C:\DOCUME~1\kevin\LOCALS~1\Temp\RSP15B.tmp" with contents
 [
 /nologo /G6 /MTd /W3 /Gm /Gi /GR /GX /Zi /Od /Gy /I "\wx2.2.5\include" /I "..\..\..\fftw-2.1.3\fftw" /I "\wx2.2.5\src\png" /I "\wx2.2.5\src\zlib" /I "..\..\include" /I "..\..\getopt" /I "..\..\..\fftw-2.1.3\rfftw" /I "\dicom\ctn\include" /D VERSION=\"3.0.0beta1\" /D "_DEBUG" /D "__WXMSW__" /D "HAVE_SGP" /D "HAVE_PNG" /D "HAVE_WXWINDOWS" /D "WIN32" /D "_WINDOWS" /D "_MBCS" /D "HAVE_STRING_H" /D "HAVE_FFTW" /D "HAVE_RFFTW" /D "HAVE_GETOPT_H" /D "MSVC" /D "__WIN95__" /D "__WIN32__" /D WINVER=0x0400 /D "STRICT" /D CTSIMVERSION=\"3.1.0\" /D "HAVE_CTN_DICOM" /D CTSIMVERSION=\"3.0.0alpha5\" /FR"Debug/" /Fp"Debug/ctsim.pch" /YX"ctsim.h" /Fo"Debug/" /Fd"Debug/" /FD /GZ /c 
-"C:\ctsim\src\docs.cpp"
+"C:\ctsim\src\views.cpp"
 ]
-Creating command line "cl.exe @C:\DOCUME~1\kevin\LOCALS~1\Temp\RSPD2.tmp" 
-Creating temporary file "C:\DOCUME~1\kevin\LOCALS~1\Temp\RSPD3.tmp" with contents
-[
-/nologo /G6 /MTd /W3 /Gm /Gi /GR /GX /Zi /Od /Gy /I "\wx2.2.5\include" /I "..\..\..\fftw-2.1.3\fftw" /I "\wx2.2.5\src\png" /I "\wx2.2.5\src\zlib" /I "..\..\include" /I "..\..\getopt" /I "..\..\..\fftw-2.1.3\rfftw" /I "\dicom\ctn\include" /D VERSION=\"3.0.0beta1\" /D "_DEBUG" /D "__WXMSW__" /D "HAVE_SGP" /D "HAVE_PNG" /D "HAVE_WXWINDOWS" /D "WIN32" /D "_WINDOWS" /D "_MBCS" /D "HAVE_STRING_H" /D "HAVE_FFTW" /D "HAVE_RFFTW" /D "HAVE_GETOPT_H" /D "MSVC" /D "__WIN95__" /D "__WIN32__" /D WINVER=0x0400 /D "STRICT" /D CTSIMVERSION=\"3.1.0\" /D "HAVE_CTN_DICOM" /FR"Debug/" /Fp"Debug/ctsim.pch" /YX"ctsim.h" /Fo"Debug/" /Fd"Debug/" /FD /GZ /c 
-"C:\ctsim\src\graph3dview.cpp"
-]
-Creating command line "cl.exe @C:\DOCUME~1\kevin\LOCALS~1\Temp\RSPD3.tmp" 
-Creating temporary file "C:\DOCUME~1\kevin\LOCALS~1\Temp\RSPD4.tmp" with contents
+Creating command line "cl.exe @C:\DOCUME~1\kevin\LOCALS~1\Temp\RSP15B.tmp" 
+Creating temporary file "C:\DOCUME~1\kevin\LOCALS~1\Temp\RSP15C.tmp" with contents
 [
 winmm.lib rpcrt4.lib ws2_32.lib ../libctsim/Debug/libctsim.lib libcmtd.lib ..\..\..\fftw-2.1.3\Win32\FFTW2st\Debug\FFTW2st.lib ..\..\..\fftw-2.1.3\Win32\RFFTW2st\Debug\RFFTW2st.lib wxd.lib xpmd.lib tiffd.lib zlibd.lib pngd.lib comctl32.lib kernel32.lib user32.lib gdi32.lib winspool.lib comdlg32.lib advapi32.lib shell32.lib ole32.lib oleaut32.lib uuid.lib odbc32.lib odbccp32.lib opengl32.lib glu32.lib htmlhelp.lib ctn_lib.lib /nologo /subsystem:windows /incremental:yes /pdb:"Debug/ctsim.pdb" /debug /machine:I386 /out:"Debug/ctsim.exe" /pdbtype:sept /libpath:"\wx2.2.5\lib" /libpath:"\dicom\ctn\winctn\ctn_lib\Debug" 
 .\Debug\backgroundmgr.obj
@@ -45,12 +39,10 @@ winmm.lib rpcrt4.lib ws2_32.lib ../libctsim/Debug/libctsim.lib libcmtd.lib ..\..
 \wx2.2.5\lib\zlibd.lib
 \wx2.2.5\lib\tiffd.lib
 ]
-Creating command line "link.exe @C:\DOCUME~1\kevin\LOCALS~1\Temp\RSPD4.tmp"
+Creating command line "link.exe @C:\DOCUME~1\kevin\LOCALS~1\Temp\RSP15C.tmp"
 <h3>Output Window</h3>
 Compiling...
-docs.cpp
-Compiling...
-graph3dview.cpp
+views.cpp
 Linking...
 
 

diff --git a/src/ctsim.h b/src/ctsim.h
index ef100ad474a3ba13b1fdc8404b5ed2f90fa21389..b7d20a9d31a7e6c5654085bd33ba6cad47fd953c 100644 (file)
--- a/src/ctsim.h
+++ b/src/ctsim.h
@@ -9,7 +9,7 @@
 **  This is part of the CTSim program
 **  Copyright (c) 1983-2001 Kevin Rosenberg
 **
-**  $Id: ctsim.h,v 1.58 2001/03/11 17:55:29 kevin Exp $
+**  $Id: ctsim.h,v 1.59 2001/03/13 08:24:41 kevin Exp $
 **
 **  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
@@ -306,6 +306,7 @@ enum {
     
     PJMENU_FILE_PROPERTIES,
     PJMENU_RECONSTRUCT_FBP,
+    PJMENU_RECONSTRUCT_FBP_REBIN,
     PJMENU_RECONSTRUCT_FOURIER,
     PJMENU_CONVERT_POLAR,
     PJMENU_CONVERT_FFT_POLAR,

diff --git a/src/threadrecon.cpp b/src/threadrecon.cpp
index e8deca41e6809f065135d49d6eaa946d88d7e144..2fde6279e7c0deca3f2b64856eea5b29e5841746 100644 (file)
--- a/src/threadrecon.cpp
+++ b/src/threadrecon.cpp
@@ -9,7 +9,7 @@
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2001 Kevin Rosenberg
 **
-**  $Id: threadrecon.cpp,v 1.24 2001/03/11 15:27:30 kevin Exp $
+**  $Id: threadrecon.cpp,v 1.25 2001/03/13 08:24:41 kevin Exp $
 **
 **  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
@@ -54,18 +54,23 @@
 ReconstructorSupervisorThread::ReconstructorSupervisorThread (ProjectionFileView* pProjView, int iNX, int iNY, 
    const char* pszFilterName, double dFilterParam, const char* pszFilterMethod, int iZeropad, 
    const char* pszFilterGenerationName, const char* pszInterpName, int iInterpParam, 
-   const char* pszBackprojectName, const char* const pszLabel, ReconstructionROI* pROI)
+   const char* pszBackprojectName, const char* const pszLabel, ReconstructionROI* pROI, bool bRebinToParallel)
 :   SupervisorThread(), m_pProjView(pProjView), m_iNX(iNX), m_iNY(iNY), m_strFilterName(pszFilterName), m_dFilterParam(dFilterParam), 
   m_strFilterMethod(pszFilterMethod), m_iZeropad(iZeropad), m_strFilterGenerationName(pszFilterGenerationName), 
   m_strInterpName(pszInterpName), m_iInterpParam(iInterpParam), m_strBackprojectName(pszBackprojectName), 
-  m_strLabel(pszLabel), m_reconROI(*pROI)
+  m_strLabel(pszLabel), m_reconROI(*pROI), m_bRebinToParallel(bRebinToParallel)
 {
 }
 
 wxThread::ExitCode
 ReconstructorSupervisorThread::Entry()
 {
-  ReconstructorSupervisor reconSupervisor (this, m_pProjView, m_iNX, m_iNY, 
+  Projections* pProj = &m_pProjView->GetDocument()->getProjections();
+
+  if (m_bRebinToParallel)
+    pProj = pProj->interpolateToParallel();
+
+  ReconstructorSupervisor reconSupervisor (this, pProj, m_pProjView, m_iNX, m_iNY, 
    m_strFilterName.c_str(), m_dFilterParam, m_strFilterMethod.c_str(), m_iZeropad, m_strFilterGenerationName.c_str(), 
    m_strInterpName.c_str(), m_iInterpParam, m_strBackprojectName.c_str(), m_strLabel.c_str(), &m_reconROI);
 
@@ -86,6 +91,9 @@ ReconstructorSupervisorThread::Entry()
          reconSupervisor.onDone();
   reconSupervisor.deleteWorkers();
 
+  if (m_bRebinToParallel)
+    delete pProj;
+
   return static_cast<wxThread::ExitCode>(0);
 }
 
@@ -101,12 +109,13 @@ ReconstructorSupervisorThread::OnExit()
 //
 /////////////////////////////////////////////////////////////////////
 
-ReconstructorSupervisor::ReconstructorSupervisor (SupervisorThread* pThread, ProjectionFileView* pProjView, 
-  int iImageNX, int iImageNY, const char* pszFilterName, double dFilterParam, const char* pszFilterMethod, 
-  int iZeropad, const char* pszFilterGenerationName, const char* pszInterpName, int iInterpParam,
-  const char* pszBackprojectName, const char* const pszLabel, ReconstructionROI* pROI)
+ReconstructorSupervisor::ReconstructorSupervisor (SupervisorThread* pThread, Projections* pProj,
+  ProjectionFileView* pProjView, int iImageNX, int iImageNY, const char* pszFilterName, double dFilterParam, 
+  const char* pszFilterMethod, int iZeropad, const char* pszFilterGenerationName, 
+  const char* pszInterpName, int iInterpParam, const char* pszBackprojectName, const char* const pszLabel, 
+  ReconstructionROI* pROI)
     : BackgroundSupervisor (pThread, pProjView->GetFrame(), pProjView->GetDocument(), "Reconstructing", pProjView->GetDocument()->getProjections().nView()),
-      m_pProjView(pProjView), m_pProjDoc(pProjView->GetDocument()), 
+      m_pProj(pProj), m_pProjView(pProjView), m_pProjDoc(pProjView->GetDocument()), 
       m_iImageNX(iImageNX), m_iImageNY(iImageNY), 
       m_pszFilterName(pszFilterName), m_dFilterParam(dFilterParam), m_pszFilterMethod(pszFilterMethod),
       m_iZeropad(iZeropad), m_pszFilterGenerationName(pszFilterGenerationName), m_pszInterpName(pszInterpName),
@@ -132,8 +141,8 @@ BackgroundWorkerThread*
 ReconstructorSupervisor::createWorker (int iThread, int iStartUnit, int iNumUnits)
 {
    ReconstructorWorker* pThread = new ReconstructorWorker (this, iThread, iStartUnit, iNumUnits);
-   pThread->SetParameters (m_pProjView, m_vecpChildImageFile[iThread], m_pszFilterName, m_dFilterParam, 
-     m_pszFilterMethod, m_iZeropad, m_pszFilterGenerationName, m_pszInterpName,
+   pThread->SetParameters (m_pProj, m_pProjView, m_vecpChildImageFile[iThread], m_pszFilterName, 
+     m_dFilterParam, m_pszFilterMethod, m_iZeropad, m_pszFilterGenerationName, m_pszInterpName,
      m_iInterpParam, m_pszBackprojectName, m_pReconROI);
 
    return pThread;
@@ -146,7 +155,7 @@ ReconstructorSupervisor::onDone()
   wxCriticalSectionLocker critsect (doneSection);
 
   ImageFile* pImageFile = getImageFile();
-  pImageFile->labelAdd (m_pProjView->GetDocument()->getProjections().getLabel());
+  pImageFile->labelAdd (m_pProj->getLabel());
   pImageFile->labelAdd (m_pszLabel, getTimerEnd());
 
   wxCommandEvent eventLog (wxEVT_COMMAND_MENU_SELECTED, MAINMENU_LOG_EVENT );
@@ -189,11 +198,12 @@ ReconstructorSupervisor::getImageFile()
 /////////////////////////////////////////////////////////////////////
 
 void
-ReconstructorWorker::SetParameters (ProjectionFileView* pProjView, ImageFile* pImageFile, 
+ReconstructorWorker::SetParameters (const Projections* pProj, ProjectionFileView* pProjView, ImageFile* pImageFile, 
  const char* pszFilterName, double dFilterParam, const char* pszFilterMethod, int iZeropad,
  const char* pszFilterGenerationName, const char* pszInterpName, int iInterpParam, 
  const char* pszBackprojectName, ReconstructionROI* pROI)
 {
+   m_pProj = pProj;
    m_pProjView = pProjView;
    m_pImageFile = pImageFile;
    m_pszFilterName = pszFilterName;
@@ -210,10 +220,9 @@ ReconstructorWorker::SetParameters (ProjectionFileView* pProjView, ImageFile* pI
 wxThread::ExitCode
 ReconstructorWorker::Entry ()
 {
-  Reconstructor* pReconstructor = new Reconstructor (m_pProjView->GetDocument()->getProjections(), 
-    *m_pImageFile, m_pszFilterName, m_dFilterParam, m_pszFilterMethod, m_iZeropad, 
-    m_pszFilterGenerationName, m_pszInterpName, m_iInterpParam, m_pszBackprojectName, Trace::TRACE_NONE,
-    m_pReconROI);
+  Reconstructor* pReconstructor = new Reconstructor (*m_pProj, *m_pImageFile, m_pszFilterName, 
+    m_dFilterParam, m_pszFilterMethod, m_iZeropad, m_pszFilterGenerationName, m_pszInterpName, 
+    m_iInterpParam, m_pszBackprojectName, Trace::TRACE_NONE, m_pReconROI, false);
 
   bool bFail = pReconstructor->fail();
   std::string failMsg;

diff --git a/src/threadrecon.h b/src/threadrecon.h
index 0b56a8c69cf33649fa29a4d7fd72ff767468c3e9..d0987caef9e5f8f77b81d7535aacf835f6f684d1 100644 (file)
--- a/src/threadrecon.h
+++ b/src/threadrecon.h
@@ -9,7 +9,7 @@
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2001 Kevin Rosenberg
 **
-**  $Id: threadrecon.h,v 1.13 2001/03/11 15:27:30 kevin Exp $
+**  $Id: threadrecon.h,v 1.14 2001/03/13 08:24:41 kevin Exp $
 **
 **  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
@@ -61,12 +61,13 @@ private:
   const std::string m_strBackprojectName;
   const std::string m_strLabel;
   ReconstructionROI m_reconROI;
+  const bool m_bRebinToParallel;
 
 public:
   ReconstructorSupervisorThread(ProjectionFileView* pProjView, int iNX, int iNY, const char* pszFilterName, 
    double dFilterParam, const char* pszFilterMethod, int iZeropad, const char* pszFilterGenerationName, 
    const char* pszInterpName, int iInterpParam, const char* pszBackprojectName, const char* const pszLabel,
-   ReconstructionROI* pROI);
+   ReconstructionROI* pROI, bool bRebinToParallel);
 
   virtual wxThread::ExitCode Entry();
 
@@ -79,6 +80,7 @@ class ReconstructorSupervisor : public BackgroundSupervisor {
 private:
 
   std::vector<ImageFile*> m_vecpChildImageFile;
+  const Projections* m_pProj;
   ProjectionFileView* m_pProjView;
   ProjectionFileDocument* m_pProjDoc;
     
@@ -97,10 +99,10 @@ private:
   ReconstructionROI* m_pReconROI;
 
 public:
-   ReconstructorSupervisor (SupervisorThread* pMyThread, ProjectionFileView* pProjView, int iNX, int iNY, const char* pszFilterName, 
-   double dFilterParam, const char* pszFilterMethod, int iZeropad, const char* pszFilterGenerationName, 
-   const char* pszInterpName, int iInterpParam, const char* pszBackprojectName, const char* const pszLabel,
-   ReconstructionROI* pReconROI);
+   ReconstructorSupervisor (SupervisorThread* pMyThread, Projections* pProj, ProjectionFileView* pProjView, 
+   int iNX, int iNY, const char* pszFilterName, double dFilterParam, const char* pszFilterMethod, int iZeropad, 
+   const char* pszFilterGenerationName, const char* pszInterpName, int iInterpParam, 
+   const char* pszBackprojectName, const char* const pszLabel, ReconstructionROI* pReconROI);
 
    virtual BackgroundWorkerThread* createWorker (int iThread, int iStartUnit, int iNumUnits);
 
@@ -113,10 +115,9 @@ public:
 };
 
 
-
-
 class ReconstructorWorker : public BackgroundWorkerThread {
 private:
+  const Projections* m_pProj;
   ProjectionFileView* m_pProjView;
   ImageFile* m_pImageFile;
   const char* m_pszFilterName;
@@ -134,7 +135,7 @@ public:
     : BackgroundWorkerThread (pSupervisor, iThread, iStartView, iNumViews)
   {}
   
-  void SetParameters (ProjectionFileView* pProjFile, ImageFile* pImageFile, 
+  void SetParameters (const Projections* pProj, ProjectionFileView* pProjFile, ImageFile* pImageFile, 
    const char* const pszFilterName, double dFilterParam, const char* const pszFilterMethod, 
    int iZeropad, const char* const pszFilterGenerationName, const char* const pszInterpName, int iInterpParam,
    const char* pszBackprojectName, ReconstructionROI* pROI);

diff --git a/src/views.cpp b/src/views.cpp
index dfce19050dd228534721978c7f8a1abca95674d3..8b40f35744b2fd2039288dfa165e1e3d2c4401a7 100644 (file)
--- a/src/views.cpp
+++ b/src/views.cpp
@@ -9,7 +9,7 @@
 **  This is part of the CTSim program
 **  Copyright (c) 1983-2001 Kevin Rosenberg
 **
-**  $Id: views.cpp,v 1.136 2001/03/13 04:44:25 kevin Exp $
+**  $Id: views.cpp,v 1.137 2001/03/13 08:24:41 kevin Exp $
 **
 **  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
@@ -2378,6 +2378,7 @@ IMPLEMENT_DYNAMIC_CLASS(ProjectionFileView, wxView)
 BEGIN_EVENT_TABLE(ProjectionFileView, wxView)
 EVT_MENU(PJMENU_FILE_PROPERTIES, ProjectionFileView::OnProperties)
 EVT_MENU(PJMENU_RECONSTRUCT_FBP, ProjectionFileView::OnReconstructFBP)
+EVT_MENU(PJMENU_RECONSTRUCT_FBP_REBIN, ProjectionFileView::OnReconstructFBPRebin)
 EVT_MENU(PJMENU_RECONSTRUCT_FOURIER, ProjectionFileView::OnReconstructFourier)
 EVT_MENU(PJMENU_CONVERT_POLAR, ProjectionFileView::OnConvertPolar)
 EVT_MENU(PJMENU_CONVERT_FFT_POLAR, ProjectionFileView::OnConvertFFTPolar)
@@ -2599,10 +2600,23 @@ ProjectionFileView::OnReconstructFourier (wxCommandEvent& event)
   wxMessageBox ("Fourier Reconstruction is not yet supported", "Unimplemented function");
 }
 
+void
+ProjectionFileView::OnReconstructFBPRebin (wxCommandEvent& event)
+{
+  Projections& rProj = GetDocument()->getProjections();
+  doReconstructFBP (rProj, true);
+}
+
 void
 ProjectionFileView::OnReconstructFBP (wxCommandEvent& event)
 {
-  const Projections& rProj = GetDocument()->getProjections();
+  Projections& rProj = GetDocument()->getProjections();
+  doReconstructFBP (rProj, false);
+}
+
+void
+ProjectionFileView::doReconstructFBP (const Projections& rProj, bool bRebinToParallel)
+{
   ReconstructionROI defaultROI;
   defaultROI.m_dXMin = -rProj.phmLen() / 2;
   defaultROI.m_dXMax = defaultROI.m_dXMin + rProj.phmLen();
@@ -2641,14 +2655,17 @@ ProjectionFileView::OnReconstructFBP (wxCommandEvent& event)
   
   std::ostringstream os;
   os << "Reconstruct " << rProj.getFilename() << ": xSize=" << m_iDefaultNX << ", ySize=" << m_iDefaultNY << ", Filter=" << optFilterName.c_str() << ", FilterParam=" << m_dDefaultFilterParam << ", FilterMethod=" << optFilterMethodName.c_str() << ", FilterGeneration=" << optFilterGenerationName.c_str() << ", Zeropad=" << m_iDefaultZeropad << ", Interpolation=" << optInterpName.c_str() << ", InterpolationParam=" << m_iDefaultInterpParam << ", Backprojection=" << optBackprojectName.c_str();
-  
+  if (bRebinToParallel)
+    os << "; Interpolate to Parallel";
+
   Timer timerRecon;
   ImageFile* pImageFile = NULL;
   if (m_iDefaultTrace > Trace::TRACE_CONSOLE) {
     pImageFile = new ImageFile (m_iDefaultNX, m_iDefaultNY);
     Reconstructor* pReconstructor = new Reconstructor (rProj, *pImageFile, optFilterName.c_str(), 
       m_dDefaultFilterParam, optFilterMethodName.c_str(), m_iDefaultZeropad, optFilterGenerationName.c_str(), 
-      optInterpName.c_str(), m_iDefaultInterpParam, optBackprojectName.c_str(), m_iDefaultTrace, &defaultROI);
+      optInterpName.c_str(), m_iDefaultInterpParam, optBackprojectName.c_str(), m_iDefaultTrace, 
+      &defaultROI, bRebinToParallel);
     
     ReconstructDialog* pDlgReconstruct = new ReconstructDialog (*pReconstructor, rProj, *pImageFile, m_iDefaultTrace, getFrameForChild());
     for (int iView = 0; iView < rProj.nView(); iView++) {
@@ -2671,10 +2688,10 @@ ProjectionFileView::OnReconstructFBP (wxCommandEvent& event)
   } else {
 #if HAVE_WXTHREADS
     if (theApp->getUseBackgroundTasks()) {
-      ReconstructorSupervisorThread* pReconstructor = new ReconstructorSupervisorThread (this, 
-        m_iDefaultNX, m_iDefaultNY, optFilterName.c_str(), 
-        m_dDefaultFilterParam, optFilterMethodName.c_str(), m_iDefaultZeropad, optFilterGenerationName.c_str(), 
-        optInterpName.c_str(), m_iDefaultInterpParam, optBackprojectName.c_str(), os.str().c_str(), &defaultROI);
+      ReconstructorSupervisorThread* pReconstructor = new ReconstructorSupervisorThread (this, m_iDefaultNX, 
+        m_iDefaultNY, optFilterName.c_str(), m_dDefaultFilterParam, optFilterMethodName.c_str(), 
+        m_iDefaultZeropad, optFilterGenerationName.c_str(), optInterpName.c_str(), m_iDefaultInterpParam, 
+        optBackprojectName.c_str(), os.str().c_str(), &defaultROI, bRebinToParallel);
       if (pReconstructor->Create() != wxTHREAD_NO_ERROR) {
         sys_error (ERR_SEVERE, "Error creating reconstructor thread");
         delete pReconstructor;
@@ -2687,11 +2704,12 @@ ProjectionFileView::OnReconstructFBP (wxCommandEvent& event)
 #endif
     {
       pImageFile = new ImageFile (m_iDefaultNX, m_iDefaultNY);
+      wxProgressDialog dlgProgress (wxString("Reconstruction"), wxString("Reconstruction Progress"), rProj.nView() + 1, getFrameForChild(), wxPD_CAN_ABORT );
       Reconstructor* pReconstructor = new Reconstructor (rProj, *pImageFile, optFilterName.c_str(), 
         m_dDefaultFilterParam, optFilterMethodName.c_str(), m_iDefaultZeropad, optFilterGenerationName.c_str(), 
-        optInterpName.c_str(), m_iDefaultInterpParam, optBackprojectName.c_str(), m_iDefaultTrace, &defaultROI);
+        optInterpName.c_str(), m_iDefaultInterpParam, optBackprojectName.c_str(), m_iDefaultTrace, 
+        &defaultROI, bRebinToParallel);
       
-      wxProgressDialog dlgProgress (wxString("Reconstruction"), wxString("Reconstruction Progress"), rProj.nView() + 1, getFrameForChild(), wxPD_CAN_ABORT );
       for (int iView = 0; iView < rProj.nView(); iView++) {
         pReconstructor->reconstructView (iView, 1);
         if (! dlgProgress.Update (iView + 1)) {
@@ -2800,6 +2818,7 @@ ProjectionFileView::CreateChildFrame(wxDocument *doc, wxView *view)
   
   wxMenu *reconstruct_menu = new wxMenu;
   reconstruct_menu->Append (PJMENU_RECONSTRUCT_FBP, "&Filtered Backprojection...\tCtrl-R", "Reconstruct image using filtered backprojection");
+  reconstruct_menu->Append (PJMENU_RECONSTRUCT_FBP_REBIN, "Filtered &Backprojection (Rebin to Parallel)...\tCtrl-B", "Reconstruct image using filtered backprojection");
   reconstruct_menu->Append (PJMENU_RECONSTRUCT_FOURIER, "&Fourier...\tCtrl-E", "Reconstruct image using inverse Fourier");
   reconstruct_menu->Enable (PJMENU_RECONSTRUCT_FOURIER, false);
   
@@ -2821,14 +2840,15 @@ ProjectionFileView::CreateChildFrame(wxDocument *doc, wxView *view)
   subframe->SetMenuBar(menu_bar);  
   subframe->Centre(wxBOTH);
   
-  wxAcceleratorEntry accelEntries[6];
+  wxAcceleratorEntry accelEntries[7];
   accelEntries[0].Set (wxACCEL_CTRL, static_cast<int>('L'), PJMENU_CONVERT_POLAR);
   accelEntries[1].Set (wxACCEL_CTRL, static_cast<int>('M'), PJMENU_CONVERT_FFT_POLAR);
   accelEntries[2].Set (wxACCEL_CTRL, static_cast<int>('R'), PJMENU_RECONSTRUCT_FBP);
-  accelEntries[3].Set (wxACCEL_CTRL, static_cast<int>('E'), PJMENU_RECONSTRUCT_FOURIER);
-  accelEntries[4].Set (wxACCEL_CTRL, static_cast<int>('I'), PJMENU_FILE_PROPERTIES);
-  accelEntries[5].Set (wxACCEL_CTRL, static_cast<int>('T'), PJMENU_PLOT_TTHETA_SAMPLING);
-  wxAcceleratorTable accelTable (6, accelEntries);
+  accelEntries[3].Set (wxACCEL_CTRL, static_cast<int>('B'), PJMENU_RECONSTRUCT_FBP_REBIN);
+  accelEntries[4].Set (wxACCEL_CTRL, static_cast<int>('E'), PJMENU_RECONSTRUCT_FOURIER);
+  accelEntries[5].Set (wxACCEL_CTRL, static_cast<int>('I'), PJMENU_FILE_PROPERTIES);
+  accelEntries[6].Set (wxACCEL_CTRL, static_cast<int>('T'), PJMENU_PLOT_TTHETA_SAMPLING);
+  wxAcceleratorTable accelTable (7, accelEntries);
   subframe->SetAcceleratorTable (accelTable);
   
   return subframe;

diff --git a/src/views.h b/src/views.h
index 9aca9648cc4565f8a58cb96e4fe0777eeb284563..69c2a1424ee4f3dfe4b67481fc3ed3fecbeddae5 100644 (file)
--- a/src/views.h
+++ b/src/views.h
@@ -9,7 +9,7 @@
 **  This is part of the CTSim program
 **  Copyright (c) 1983-2001 Kevin Rosenberg
 **
-**  $Id: views.h,v 1.50 2001/03/11 17:55:29 kevin Exp $
+**  $Id: views.h,v 1.51 2001/03/13 08:24:41 kevin Exp $
 **
 **  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
@@ -241,6 +241,7 @@ public:
   bool OnClose (bool deleteWindow = true);
   void OnProperties (wxCommandEvent& event);
   void OnReconstructFBP (wxCommandEvent& event);
+  void OnReconstructFBPRebin (wxCommandEvent& event);
   void OnReconstructFourier (wxCommandEvent& event);
   void OnConvertPolar (wxCommandEvent& event);
   void OnConvertFFTPolar (wxCommandEvent& event);
@@ -248,6 +249,8 @@ public:
   void OnConvertParallel (wxCommandEvent& event);
   void OnArtifactReduction (wxCommandEvent& event);
 
+  void doReconstructFBP (const Projections& rProj, bool bRebinToParallel);
+
 #if CTSIM_MDI
   wxDocMDIChildFrame* getFrame() { return m_pFrame; }
 #else