X-Git-Url: http://git.kpe.io/?a=blobdiff_plain;f=libctsim%2Fbackprojectors.cpp;h=bfd594e8697ef0075b98ce6bce5b325261b4b082;hb=b7cbad1230585741079fad4881a8d67be2592100;hp=c7a8f27f47f7bd34b1523af9d852760ff967cf35;hpb=d3fa225aa232e132cc198672c4fc148f96a1ab8c;p=ctsim.git

diff --git a/libctsim/backprojectors.cpp b/libctsim/backprojectors.cpp
index c7a8f27..bfd594e 100644
--- a/libctsim/backprojectors.cpp
+++ b/libctsim/backprojectors.cpp
@@ -8,7 +8,7 @@
 **  This is part of the CTSim program
 **  Copyright (c) 1983-2001 Kevin Rosenberg
 **
-**  $Id: backprojectors.cpp,v 1.26 2001/02/11 04:56:37 kevin Exp $
+**  $Id: backprojectors.cpp,v 1.30 2001/03/07 19:02:44 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
@@ -32,7 +32,7 @@ const int Backprojector::BPROJ_TABLE = 1;
 const int Backprojector::BPROJ_DIFF = 2;
 const int Backprojector::BPROJ_IDIFF = 3;
 
-const char* Backprojector::s_aszBackprojectName[] = 
+const char* const Backprojector::s_aszBackprojectName[] = 
 {
   {"trig"},
   {"table"},
@@ -40,7 +40,7 @@ const char* Backprojector::s_aszBackprojectName[] =
   {"idiff"},
 };
 
-const char* Backprojector::s_aszBackprojectTitle[] = 
+const char* const Backprojector::s_aszBackprojectTitle[] = 
 {
   {"Direct Trigometric"},
   {"Trigometric Table"},
@@ -62,7 +62,7 @@ const int Backprojector::INTERP_2BSPLINE = 6;
 const int Backprojector::INTERP_3BSPLINE = 7;
 #endif
 
-const char* Backprojector::s_aszInterpName[] = 
+const char* const Backprojector::s_aszInterpName[] = 
 {
   {"nearest"},
   {"linear"},
@@ -78,7 +78,7 @@ const char* Backprojector::s_aszInterpName[] =
 #endif
 };
 
-const char* Backprojector::s_aszInterpTitle[] = 
+const char* const Backprojector::s_aszInterpTitle[] = 
 {
   {"Nearest"},
   {"Linear"},
@@ -113,6 +113,13 @@ Backprojector::BackprojectView (const double* const viewData, const double viewA
     m_pBackprojectImplem->BackprojectView (viewData, viewAngle);
 }
 
+void 
+Backprojector::PostProcessing()
+{
+  if (m_pBackprojectImplem != NULL)
+    m_pBackprojectImplem->PostProcessing();
+}
+
 Backprojector::~Backprojector ()
 {
   delete m_pBackprojectImplem;
@@ -254,7 +261,7 @@ Backprojector::convertInterpIDToTitle (const int interpID)
 //   Pure virtual base class for all backprojectors.
 
 Backproject::Backproject (const Projections& proj, ImageFile& im, int interpType, const int interpFactor)
-: proj(proj), im(im), interpType(interpType), m_interpFactor(interpFactor)
+: proj(proj), im(im), interpType(interpType), m_interpFactor(interpFactor), m_bPostProcessingDone(false)
 {
   detInc = proj.detInc();
   nDet = proj.nDet();
@@ -282,11 +289,18 @@ Backproject::Backproject (const Projections& proj, ImageFile& im, int interpType
   yInc = (yMax - yMin) / ny;
   
   m_dFocalLength = proj.focalLength();
+  m_dSourceDetectorLength = proj.sourceDetectorLength();
 }
 
 Backproject::~Backproject ()
 {}
 
+void
+Backproject::PostProcessing()
+{
+  m_bPostProcessingDone = true;
+}
+
 void
 Backproject::ScaleImageByRotIncrement ()
 {
@@ -351,7 +365,7 @@ BackprojectTrig::BackprojectView (const double* const filteredProj, const double
         double frac = p - pFloor;	// fraction distance from det
         if (iDetPos >= 0 && iDetPos < nDet - 1)
           v[ix][iy] += rotScale * ((1-frac) * filteredProj[iDetPos] + frac * filteredProj[iDetPos+1]);
-      } else if (interpType = Backprojector::INTERP_CUBIC) {
+      } else if (interpType == Backprojector::INTERP_CUBIC) {
         double p = iDetCenter + (L / detInc);	// position along detector
         if (p >= 0 && p < nDet)
           v[ix][iy] += rotScale * pCubicInterp->interpolate (p);
@@ -389,7 +403,15 @@ BackprojectTable::BackprojectTable (const Projections& proj, ImageFile& im, int
 
 BackprojectTable::~BackprojectTable ()
 {
-  ScaleImageByRotIncrement();
+}
+
+void
+BackprojectTable::PostProcessing()
+{
+  if (! m_bPostProcessingDone) {
+    ScaleImageByRotIncrement();
+    m_bPostProcessingDone = true;
+  }
 }
 
 void
@@ -419,7 +441,7 @@ BackprojectTable::BackprojectView (const double* const filteredProj, const doubl
         double frac = dPos - dPosFloor;	// fraction distance from det 
         if (iDetPos >= 0 && iDetPos < nDet - 1)
           pImCol[iy] += ((1-frac) * filteredProj[iDetPos] + frac * filteredProj[iDetPos+1]);
-      } else if (interpType = Backprojector::INTERP_CUBIC) {
+      } else if (interpType == Backprojector::INTERP_CUBIC) {
         double p = iDetCenter + (L / detInc);	// position along detector
         if (p >= 0 && p < nDet)
           pImCol[iy] += pCubicInterp->interpolate (p);
@@ -451,11 +473,18 @@ BackprojectDiff::BackprojectDiff (const Projections& proj, ImageFile& im, int in
   im.arrayDataClear();
 }
 
-BackprojectDiff::~BackprojectDiff()
+BackprojectDiff::~BackprojectDiff ()
 {
-  ScaleImageByRotIncrement();
 }
 
+void
+BackprojectDiff::PostProcessing()
+{
+  if (! m_bPostProcessingDone) {
+    ScaleImageByRotIncrement();
+    m_bPostProcessingDone = true;
+  }
+}
 
 void
 BackprojectDiff::BackprojectView (const double* const filteredProj, const double view_angle)
@@ -489,7 +518,7 @@ BackprojectDiff::BackprojectView (const double* const filteredProj, const double
         double frac = curDetPos - detPosFloor;	// fraction distance from det 
         if (iDetPos > 0 && iDetPos < nDet - 1)
           *pImCol++ += filteredProj[iDetPos] + (frac * (filteredProj[iDetPos+1] - filteredProj[iDetPos]));
-      } else if (interpType = Backprojector::INTERP_CUBIC) {
+      } else if (interpType == Backprojector::INTERP_CUBIC) {
         double p = iDetCenter + curDetPos;	// position along detector
         if (p >= 0 && p < nDet)
           *pImCol++  += pCubicInterp->interpolate (p);
@@ -560,7 +589,7 @@ BackprojectIntDiff::BackprojectView (const double* const filteredProj, const dou
         if (iDetPos >= 0 && iDetPos <= iLastDet)
           *pImCol++ += filteredProj[iDetPos] + (detRemainder * deltaFilteredProj[iDetPos]);
       }	// end for iy
-    } else if (interpType = Backprojector::INTERP_CUBIC) {
+    } else if (interpType == Backprojector::INTERP_CUBIC) {
       for (int iy = 0; iy < ny; iy++, curDetPos += det_dy) {
         *pImCol++ += pCubicInterp->interpolate (static_cast<double>(curDetPos) / 65536);
       }
@@ -636,10 +665,8 @@ BackprojectEquilinear::BackprojectView (const double* const filteredProj, const
       
       double dU = (m_dFocalLength + rsin_t) / m_dFocalLength;
       double dDetPos =  rcos_t / dU;
-      // double to scale for imaginary detector that passes through origin
-      // of phantom, see Kak-Slaney Figure 3.22. This assumes that the detector is also
-      // located focal-length away from the origin.
-      dDetPos *= 2; 
+      // Scale for imaginary detector that passes through origin of phantom, see Kak-Slaney Figure 3.22. 
+      dDetPos *= m_dSourceDetectorLength / m_dFocalLength; 
       double dPos = dDetPos / detInc;  // position along detector array 
 
       if (interpType == Backprojector::INTERP_NEAREST) {