r184: *** empty log message ***

[ctsim.git] / libctsim / backprojectors.cpp

diff --git a/libctsim/backprojectors.cpp b/libctsim/backprojectors.cpp
index 7a1c55bccf080b304695190e9d5b7df096512c85..71aa8d4824722d2940d874b0b87666007121954f 100644 (file)
--- a/libctsim/backprojectors.cpp
+++ b/libctsim/backprojectors.cpp
@@ -8,7 +8,7 @@
 **  This is part of the CTSim program
 **  Copyright (C) 1983-2000 Kevin Rosenberg
 **
-**  $Id: backprojectors.cpp,v 1.11 2000/07/23 01:49:03 kevin Exp $
+**  $Id: backprojectors.cpp,v 1.12 2000/08/25 15:59:13 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
@@ -286,18 +286,20 @@ Backproject::ScaleImageByRotIncrement ()
 
 void Backproject::errorIndexOutsideDetector (int ix, int iy, double theta, double r, double phi, double L, int iDetPos)
 {
-    printf ("r=%f, phi=%f\n", r, phi);
+    sys_error (ERR_WARNING, "r=%f, phi=%f", r, phi);
     errorIndexOutsideDetector (ix, iy, theta, L, iDetPos);
 }
 
 void Backproject::errorIndexOutsideDetector (int ix, int iy, double theta, double L, int iDetPos)
 {
-    printf ("ix=%d, iy=%d\n", ix, iy);
-    printf ("theta=%f, L=%f, detInc=%f\n", theta, L, detInc);
-    printf ("proj.ndet=%d, proj.detInc=%.4f, iDetCenter=%d\n", nDet, detInc, iDetCenter);
-    printf ("xMin=%15.8f, xMax=%15.8f, xInc=%15.8f\n", xMin, xMax, xInc);
-    printf ("yMin=%15.8f, yMax=%15.8f, yInc=%15.8f\n", yMin, yMax, yInc);
-    sys_error (ERR_WARNING, "iDetPos index outside bounds: %d [backprojector]", iDetPos);
+  ostringstream os;
+  os << "ix=" << ix << ", iy=" << iy << ", theta=" << theta << ", L=" << L << ", detinc=" << detInc << "\n";
+  os << "ndet=" << nDet << ", detInc=" << detInc << ", iDetCenter=" << iDetCenter << "\n";
+  os << "xMin=" << xMin << ", xMax=" << xMax << ", xInc=" << xInc << "\n";
+  os << "yMin=" << yMin << ", yMax=" << yMax << ", yInc=" << yInc << "\n";
+  os << "iDetPos index outside bounds: " << iDetPos << " [backprojector]";;
+
+  sys_error (ERR_WARNING, os.str().c_str());
 }
 
 
@@ -543,7 +545,7 @@ BackprojectIntDiff2::BackprojectView (const double* const filteredProj, const do
        int detPosNearest = (curDetPos >= 0 ? ((curDetPos + halfScale) / scale) : ((curDetPos - halfScale) / scale));
        int iDetPos = iDetCenter + detPosNearest;       // calc index in the filtered raysum vector 
 
-       if (iDetPos < 0 || iDetPos >= nDet)  // check for impossible: index outside of raysum pos 
+       if (iDetPos < 0 || iDetPos >= nDet)  // check for index outside of raysum pos 
            errorIndexOutsideDetector (ix, iy, theta, curDetPos, iDetPos);
        else
          *pImCol++ += filteredProj[iDetPos];
@@ -588,12 +590,14 @@ BackprojectIntDiff3::BackprojectView (const double* const filteredProj, const do
   kint32 detPosColStart = nearest<kint32> ((start_r * cos (theta - start_phi) / detInc + iDetCenter) * scale);
        
   // precalculate scaled difference for linear interpolation
-  double deltaFilteredProj [nDet - 1];
+  double deltaFilteredProj [nDet];
   if (interpType == Backprojector::INTERP_LINEAR) {
     for (int i = 0; i < nDet - 1; i++)
       deltaFilteredProj[i] = (filteredProj[i+1] - filteredProj[i]) * dInvScale;
   }
+  deltaFilteredProj[nDet - 1] = 0;  // last detector
 
+  int iLastDet = nDet - 1;
   for (int ix = 0; ix < nx; ix++, detPosColStart += det_dx) {
     kint32 curDetPos = detPosColStart;
     ImageFileColumn pImCol = v[ix];
@@ -601,21 +605,21 @@ BackprojectIntDiff3::BackprojectView (const double* const filteredProj, const do
     if (interpType == Backprojector::INTERP_NEAREST) {
       for (int iy = 0; iy < ny; iy++, curDetPos += det_dy) {
        const int iDetPos = (curDetPos + halfScale) >> 16;
-       assert(iDetPos >= 0 && iDetPos < nDet);
-       *pImCol++ += filteredProj[iDetPos];
+       if (iDetPos >= 0 && iDetPos <= iLastDet)
+         *pImCol++ += filteredProj[iDetPos];
       }        // end for iy
     } else if (interpType == Backprojector::INTERP_FREQ_PREINTERPOLATION) {
       for (int iy = 0; iy < ny; iy++, curDetPos += det_dy) {
        const int iDetPos = ((curDetPos + halfScale) >> 16) * m_interpFactor;
-       assert(iDetPos >= 0 && iDetPos < nDet);
+       if (iDetPos >= 0 && iDetPos <= iLastDet)
        *pImCol++ += filteredProj[iDetPos];
       }        // end for iy
     } else if (interpType == Backprojector::INTERP_LINEAR) {
       for (int iy = 0; iy < ny; iy++, curDetPos += det_dy) {
        const kint32 iDetPos = curDetPos >> scaleShift;
        const kint32 detRemainder = curDetPos & scaleBitmask;
-       assert(iDetPos >= 0 && iDetPos < nDet - 1);
-       *pImCol++ += filteredProj[iDetPos] + (detRemainder * deltaFilteredProj[iDetPos]);
+       if (iDetPos >= 0 && iDetPos <= iLastDet)
+         *pImCol++ += filteredProj[iDetPos] + (detRemainder * deltaFilteredProj[iDetPos]);
       }        // end for iy
     } //end linear
   } // end for ix