-// CLASS IDENTICATION
-// BackprojectIntDiff2
-//
-// PURPOSE
-// Integer version of BackprojectDiff2
-
-void
-BackprojectIntDiff2::BackprojectView (const double* const filteredProj, const double view_angle)
-{
- double theta = view_angle;
-
- CubicInterpolator* pCubicInterp = NULL;
- if (interpType == Backprojector::INTERP_CUBIC)
- pCubicInterp = new CubicInterpolator (filteredProj, nDet);
-
- static const kint32 scale = 1 << 16;
- static const double dScale = scale;
- static const kint32 halfScale = scale / 2;
-
- const kint32 det_dx = nearest<kint32> (xInc * cos (theta) / detInc * scale);
- const kint32 det_dy = nearest<kint32> (yInc * sin (theta) / detInc * scale);
-
- // calculate L for first point in image (0, 0)
- kint32 detPosColStart = nearest<kint32> (start_r * cos (theta - start_phi) / detInc * 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 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 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]);
- } else if (interpType = Backprojector::INTERP_CUBIC) {
- double p = iDetCenter + (static_cast<double>(curDetPos) / scale); // position along detector
- if (p >= 0 && p < nDet)
- *pImCol++ += pCubicInterp->interpolate (p);
- }
- } // end for y
- } // end for x
-
- if (interpType == Backprojector::INTERP_CUBIC)
- delete pCubicInterp;
-}