+ *pImCol++ += ((1.-frac) * filteredProj[iDetPos] + frac * filteredProj[iDetPos+1]);
+ }
+ } // end for y
+ } // end for x
+}
+
+// CLASS IDENTICATION
+// BackprojectIntDiff3
+//
+// PURPOSE
+// Highly optimized version of BackprojectIntDiff2
+
+void
+BackprojectIntDiff3::BackprojectView (const double* const filteredProj, const double view_angle)
+{
+ double theta = - view_angle; // add half PI to view angle to get perpendicular theta angle
+ static const int scaleShift = 16;
+ static const kint32 scale = (1 << scaleShift);
+ static const double dScale = scale;
+ static const kint32 halfScale = scale / 2;
+
+ const kint32 det_dx = nearest<kint32> (xInc * sin (theta) / detInc * scale);
+ const kint32 det_dy = nearest<kint32> (yInc * cos (theta) / detInc * scale);
+
+ // calculate L for first point in image (0, 0)
+ kint32 detPosColStart = nearest<kint32> ((start_r * cos (theta - start_phi) / detInc + iDetCenter) * scale);
+
+ for (int ix = 0; ix < nx; ix++, detPosColStart += det_dx) {
+ kint32 curDetPos = detPosColStart;
+ ImageFileColumn pImCol = v[ix];
+
+ 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]);