diff options
author | Willem Jan Palenstijn <wjp@usecode.org> | 2019-03-29 19:08:11 +0100 |
---|---|---|
committer | Willem Jan Palenstijn <Willem.Jan.Palenstijn@cwi.nl> | 2019-09-25 14:10:08 +0200 |
commit | 9c869d834ddc681299df9999787eb92bd6010dac (patch) | |
tree | c21e3efd97ca910388db92062b44dc9dfdd76955 /include/astra | |
parent | 15f53527e2f38e8eeb8bece79376b5e4686f3dd4 (diff) | |
download | astra-9c869d834ddc681299df9999787eb92bd6010dac.tar.gz astra-9c869d834ddc681299df9999787eb92bd6010dac.tar.bz2 astra-9c869d834ddc681299df9999787eb92bd6010dac.tar.xz astra-9c869d834ddc681299df9999787eb92bd6010dac.zip |
Adjust strip kernels to line integral scaling
Diffstat (limited to 'include/astra')
-rw-r--r-- | include/astra/FanFlatBeamStripKernelProjector2D.inl | 12 | ||||
-rw-r--r-- | include/astra/ParallelBeamStripKernelProjector2D.inl | 8 |
2 files changed, 14 insertions, 6 deletions
diff --git a/include/astra/FanFlatBeamStripKernelProjector2D.inl b/include/astra/FanFlatBeamStripKernelProjector2D.inl index 889d0a3..f5a688c 100644 --- a/include/astra/FanFlatBeamStripKernelProjector2D.inl +++ b/include/astra/FanFlatBeamStripKernelProjector2D.inl @@ -109,8 +109,12 @@ void CFanFlatBeamStripKernelProjector2D::projectBlock_internal(int _iProjFrom, i // POLICY: RAY PRIOR if (!p.rayPrior(iRayIndex)) continue; - float32 dist_srcDetPixSquared = projgeom->getSourceDetectorDistance() * projgeom->getSourceDetectorDistance() + - (iDetector + 0.5f - m_pProjectionGeometry->getDetectorCount()*0.5f) * (iDetector + 0.5f - m_pProjectionGeometry->getDetectorCount()*0.5f) * DW * DW; + float32 dist_srcDetPixSquared = projgeom->getSourceDetectorDistance() * projgeom->getSourceDetectorDistance() + (iDetector + 0.5f - m_pProjectionGeometry->getDetectorCount()*0.5f) * (iDetector + 0.5f - m_pProjectionGeometry->getDetectorCount()*0.5f) * DW * DW; + dist_srcDetPixSquared = dist_srcDetPixSquared * dist_srcDetPixSquared / (projgeom->getSourceDetectorDistance() * projgeom->getSourceDetectorDistance() * DW * DW); + + + + //float32 InvRayWidthSquared = (projgeom->getSourceDetectorDistance() * projgeom->getSourceDetectorDistance()) / dist_srcDetPixSquared; float32 sin_theta_left, cos_theta_left; float32 sin_theta_right, cos_theta_right; @@ -257,8 +261,8 @@ void CFanFlatBeamStripKernelProjector2D::projectBlock_internal(int _iProjFrom, i // POLICY: RAY PRIOR if (!p.rayPrior(iRayIndex)) continue; - float32 dist_srcDetPixSquared = projgeom->getSourceDetectorDistance() * projgeom->getSourceDetectorDistance() + - (iDetector + 0.5f - m_pProjectionGeometry->getDetectorCount()*0.5f) * (iDetector + 0.5f - m_pProjectionGeometry->getDetectorCount()*0.5f) * DW * DW; + float32 dist_srcDetPixSquared = projgeom->getSourceDetectorDistance() * projgeom->getSourceDetectorDistance() + (iDetector + 0.5f - m_pProjectionGeometry->getDetectorCount()*0.5f) * (iDetector + 0.5f - m_pProjectionGeometry->getDetectorCount()*0.5f) * DW * DW; + dist_srcDetPixSquared = dist_srcDetPixSquared * dist_srcDetPixSquared / (projgeom->getSourceDetectorDistance() * projgeom->getSourceDetectorDistance() * DW * DW); // get theta_l = alpha_left + theta and theta_r = alpha_right + theta float32 sin_theta_left, cos_theta_left; diff --git a/include/astra/ParallelBeamStripKernelProjector2D.inl b/include/astra/ParallelBeamStripKernelProjector2D.inl index 0d775b3..bed02ab 100644 --- a/include/astra/ParallelBeamStripKernelProjector2D.inl +++ b/include/astra/ParallelBeamStripKernelProjector2D.inl @@ -142,6 +142,10 @@ void CParallelBeamStripKernelProjector2D::projectBlock_internal(int _iProjFrom, const SParProjection * proj = &pVecProjectionGeometry->getProjectionVectors()[iAngle]; + const float32 rayWidth = fabs(proj->fDetUX * proj->fRayY - proj->fDetUY * proj->fRayX) / + sqrt(proj->fRayX * proj->fRayX + proj->fRayY * proj->fRayY); + const float32 relPixelArea = pixelArea / rayWidth; + bool vertical = fabs(proj->fRayX) < fabs(proj->fRayY); if (vertical) { RxOverRy = proj->fRayX/proj->fRayY; @@ -219,7 +223,7 @@ void CParallelBeamStripKernelProjector2D::projectBlock_internal(int _iProjFrom, else if (-S < offsetL) res -= 0.5f + offsetL; else if (-T < offsetL) res -= 0.5f*(offsetL+T)*(offsetL+T)*invTminS; - p.addWeight(iRayIndex, iVolumeIndex, pixelArea*res); + p.addWeight(iRayIndex, iVolumeIndex, relPixelArea*res); p.pixelPosterior(iVolumeIndex); } } @@ -272,7 +276,7 @@ void CParallelBeamStripKernelProjector2D::projectBlock_internal(int _iProjFrom, else if (-S < offsetL) res -= 0.5f + offsetL; else if (-T < offsetL) res -= 0.5f*(offsetL+T)*(offsetL+T)*invTminS; - p.addWeight(iRayIndex, iVolumeIndex, pixelArea*res); + p.addWeight(iRayIndex, iVolumeIndex, relPixelArea*res); p.pixelPosterior(iVolumeIndex); } } |