- for (int iy = 0; iy < ny; iy++, curDetPos += det_dy) {
- if (interpType == Backprojector::INTERP_NEAREST) {
- int iDetPos = (curDetPos + halfScale) >> scaleShift;
-
- if (iDetPos < 0 || iDetPos >= nDet) // check for impossible: index outside of raysum pos
- errorIndexOutsideDetector (ix, iy, theta, curDetPos, iDetPos);
- else
- *pImCol++ += filteredProj[iDetPos];
- } else if (interpType == Backprojector::INTERP_LINEAR) {
- kint32 detPosFloor = curDetPos / scale;
- kint32 detPosRemainder = curDetPos % scale;
- if (detPosRemainder < 0) {
- detPosFloor--;
- detPosRemainder += scale;
- }
- int iDetPos = iDetCenter + detPosFloor;
- double frac = detPosRemainder / dScale;
- if (iDetPos < 0 || iDetPos >= nDet - 1)
- errorIndexOutsideDetector (ix, iy, theta, curDetPos, iDetPos);
- else
- *pImCol++ += ((1.-frac) * filteredProj[iDetPos] + frac * filteredProj[iDetPos+1]);
- }
- } // end for y
- } // end for x
+ if (interpType == Backprojector::INTERP_NEAREST) {
+ for (int iy = 0; iy < ny; iy++, curDetPos += det_dy) {
+ int iDetPos = (curDetPos + halfScale) >> 16;
+ assert(iDetPos >= 0 && iDetPos < nDet);
+ *pImCol++ += filteredProj[iDetPos];
+ } // end for iy
+ } else if (interpType == Backprojector::INTERP_LINEAR) {
+ for (int iy = 0; iy < ny; iy++, curDetPos += det_dy) {
+ kint32 iDetPos = curDetPos >> scaleShift;
+ kint32 detRemainder = curDetPos & scaleBitmask;
+#if 0
+ double frac = detRemainder * (1. / 65536.);
+ assert(iDetPos >= 0 && iDetPos < nDet - 1);
+ *pImCol++ += ((1.-frac) * filteredProj[iDetPos]) + (frac * filteredProj[iDetPos+1]);
+#else
+ assert(iDetPos >= 0 && iDetPos < nDet - 1);
+ const double* const detPointer = &filteredProj[iDetPos];
+ *pImCol++ += (((scale-detRemainder) * *detPointer) + (detRemainder * *(detPointer+1))) * dInvScale;
+#endif
+ } // end for iy
+ } //end linear
+ } // end for ix