From 3cf63d335ebe392a8c77f0c90395c18150647aeb Mon Sep 17 00:00:00 2001 From: Willem Jan Palenstijn Date: Fri, 29 Mar 2019 21:21:29 +0100 Subject: Adjust adjoint to line integral scaling --- src/CudaReconstructionAlgorithm2D.cpp | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) (limited to 'src') diff --git a/src/CudaReconstructionAlgorithm2D.cpp b/src/CudaReconstructionAlgorithm2D.cpp index 1e81390..939a026 100644 --- a/src/CudaReconstructionAlgorithm2D.cpp +++ b/src/CudaReconstructionAlgorithm2D.cpp @@ -309,8 +309,7 @@ void CCudaReconstructionAlgorithm2D::run(int _iNrIterations) m_bAlgoInit = true; } - float fPixelSize = volgeom.getPixelLengthX(); - float fSinogramScale = 1.0f/(fPixelSize*fPixelSize); + float fSinogramScale = 1.0f; ok = m_pAlgo->copyDataToGPU(m_pSinogram->getDataConst(), m_pSinogram->getGeometry()->getDetectorCount(), fSinogramScale, m_pReconstruction->getDataConst(), volgeom.getGridColCount(), -- cgit v1.2.3 From 48f4e7b165404a0375db300b9fe59da92edf1cce Mon Sep 17 00:00:00 2001 From: Willem Jan Palenstijn Date: Sat, 30 Mar 2019 20:49:42 +0100 Subject: Adjust FBP to line integral scaling --- src/CudaFilteredBackProjectionAlgorithm.cpp | 7 +++++++ src/CudaReconstructionAlgorithm2D.cpp | 4 +--- src/FilteredBackProjectionAlgorithm.cpp | 8 ++++++-- 3 files changed, 14 insertions(+), 5 deletions(-) (limited to 'src') diff --git a/src/CudaFilteredBackProjectionAlgorithm.cpp b/src/CudaFilteredBackProjectionAlgorithm.cpp index 88e235b..69140fc 100644 --- a/src/CudaFilteredBackProjectionAlgorithm.cpp +++ b/src/CudaFilteredBackProjectionAlgorithm.cpp @@ -151,6 +151,13 @@ void CCudaFilteredBackProjectionAlgorithm::initCUDAAlgorithm() if (!ok) { ASTRA_ERROR("CCudaFilteredBackProjectionAlgorithm: Failed to set short-scan mode"); } + + const CVolumeGeometry2D& volGeom = *m_pProjector->getVolumeGeometry(); + float fPixelArea = volGeom.getPixelArea(); + ok &= pFBP->setReconstructionScale(1.0f/fPixelArea); + if (!ok) { + ASTRA_ERROR("CCudaFilteredBackProjectionAlgorithm: Failed to set reconstruction scale"); + } } diff --git a/src/CudaReconstructionAlgorithm2D.cpp b/src/CudaReconstructionAlgorithm2D.cpp index 939a026..6730cea 100644 --- a/src/CudaReconstructionAlgorithm2D.cpp +++ b/src/CudaReconstructionAlgorithm2D.cpp @@ -309,9 +309,7 @@ void CCudaReconstructionAlgorithm2D::run(int _iNrIterations) m_bAlgoInit = true; } - float fSinogramScale = 1.0f; - - ok = m_pAlgo->copyDataToGPU(m_pSinogram->getDataConst(), m_pSinogram->getGeometry()->getDetectorCount(), fSinogramScale, + ok = m_pAlgo->copyDataToGPU(m_pSinogram->getDataConst(), m_pSinogram->getGeometry()->getDetectorCount(), m_pReconstruction->getDataConst(), volgeom.getGridColCount(), m_bUseReconstructionMask ? m_pReconstructionMask->getDataConst() : 0, volgeom.getGridColCount(), m_bUseSinogramMask ? m_pSinogramMask->getDataConst() : 0, m_pSinogram->getGeometry()->getDetectorCount()); diff --git a/src/FilteredBackProjectionAlgorithm.cpp b/src/FilteredBackProjectionAlgorithm.cpp index 423dc6c..95bef3c 100644 --- a/src/FilteredBackProjectionAlgorithm.cpp +++ b/src/FilteredBackProjectionAlgorithm.cpp @@ -264,8 +264,12 @@ void CFilteredBackProjectionAlgorithm::run(int _iNrIterations) DefaultBPPolicy(m_pReconstruction, &filteredSinogram)); // Scale data - int iAngleCount = m_pProjector->getProjectionGeometry()->getProjectionAngleCount(); - (*m_pReconstruction) *= (PI/2)/iAngleCount; + const CVolumeGeometry2D& volGeom = *m_pProjector->getVolumeGeometry(); + const CProjectionGeometry2D& projGeom = *m_pProjector->getProjectionGeometry(); + + int iAngleCount = projGeom.getProjectionAngleCount(); + float fPixelArea = volGeom.getPixelArea(); + (*m_pReconstruction) *= PI/(2*iAngleCount*fPixelArea); m_pReconstruction->updateStatistics(); } -- cgit v1.2.3 From 4fba091ff811ac5e6035cf95efdaae6681e6cc46 Mon Sep 17 00:00:00 2001 From: Willem Jan Palenstijn Date: Mon, 1 Apr 2019 00:09:28 +0200 Subject: Add error check for non-parallel FBP --- src/FilteredBackProjectionAlgorithm.cpp | 5 +++++ 1 file changed, 5 insertions(+) (limited to 'src') diff --git a/src/FilteredBackProjectionAlgorithm.cpp b/src/FilteredBackProjectionAlgorithm.cpp index 95bef3c..6b4093d 100644 --- a/src/FilteredBackProjectionAlgorithm.cpp +++ b/src/FilteredBackProjectionAlgorithm.cpp @@ -167,6 +167,11 @@ bool CFilteredBackProjectionAlgorithm::initialize(const Config& _cfg) m_filterConfig = getFilterConfigForAlgorithm(_cfg, this); + const CParallelProjectionGeometry2D* parprojgeom = dynamic_cast(m_pSinogram->getGeometry()); + if (!parprojgeom) { + ASTRA_ERROR("FBP currently only supports parallel projection geometries."); + return false; + } // TODO: check that the angles are linearly spaced between 0 and pi -- cgit v1.2.3 From 975537c8c68b115399af522cb1a0a1e731eca576 Mon Sep 17 00:00:00 2001 From: Willem Jan Palenstijn Date: Tue, 2 Apr 2019 15:56:57 +0200 Subject: Fix fan-beam FBP scaling --- src/GeometryUtil2D.cpp | 9 ++++++--- 1 file changed, 6 insertions(+), 3 deletions(-) (limited to 'src') diff --git a/src/GeometryUtil2D.cpp b/src/GeometryUtil2D.cpp index e09a3bc..806572f 100644 --- a/src/GeometryUtil2D.cpp +++ b/src/GeometryUtil2D.cpp @@ -28,6 +28,7 @@ along with the ASTRA Toolbox. If not, see . #include "astra/GeometryUtil2D.h" #include +#include namespace astra { @@ -158,14 +159,16 @@ bool getFanParameters(const SFanProjection &proj, unsigned int iProjDets, float // project origin on detector line ( == project source on detector line) double t = (- proj.fDetSX) * proj.fDetUX + (- proj.fDetSY) * proj.fDetUY; + t /= (proj.fDetUX * proj.fDetUX + proj.fDetUY * proj.fDetUY); fOffset = (float)t - 0.5*iProjDets; - // TODO: CHECKME fOriginDetector = sqrt((proj.fDetSX + t * proj.fDetUX)*(proj.fDetSX + t * proj.fDetUX) + (proj.fDetSY + t * proj.fDetUY)*(proj.fDetSY + t * proj.fDetUY)); - //float fAngle = atan2(proj.fDetSX + t * proj.fDetUX - proj.fSrcX, proj.fDetSY + t * proj.fDetUY); // TODO: Fix order + sign - fAngle = atan2(proj.fDetUY, proj.fDetUX); // TODO: Check order + sign + fAngle = atan2(proj.fDetUY, proj.fDetUX); + + //fprintf(stderr, "getFanParams: s = (%f,%f) d = (%f,%f) u = (%f,%f)\n", proj.fSrcX, proj.fSrcY, proj.fDetSX, proj.fDetSY, proj.fDetUX, proj.fDetUY); + //fprintf(stderr, "getFanParams: fOS = %f, fOD = %f, detsize = %f, offset = %f (t = %f), angle = %f\n", fOriginSource, fOriginDetector, fDetSize, fOffset, t, fAngle); return true; } -- cgit v1.2.3 From 4b576ee6ad461b162dd64538f54df47b65b37972 Mon Sep 17 00:00:00 2001 From: Willem Jan Palenstijn Date: Wed, 3 Apr 2019 22:31:55 +0200 Subject: Remove obsolete DensityWeighting option --- src/CompositeGeometryManager.cpp | 4 +--- src/CudaProjector3D.cpp | 8 -------- 2 files changed, 1 insertion(+), 11 deletions(-) (limited to 'src') diff --git a/src/CompositeGeometryManager.cpp b/src/CompositeGeometryManager.cpp index 1319a87..822f746 100644 --- a/src/CompositeGeometryManager.cpp +++ b/src/CompositeGeometryManager.cpp @@ -1462,12 +1462,10 @@ static bool doJob(const CCompositeGeometryManager::TJobSet::const_iterator& iter Cuda3DProjectionKernel projKernel = ker3d_default; int detectorSuperSampling = 1; int voxelSuperSampling = 1; - bool densityWeighting = false; if (projector) { projKernel = projector->getProjectionKernel(); detectorSuperSampling = projector->getDetectorSuperSampling(); voxelSuperSampling = projector->getVoxelSuperSampling(); - densityWeighting = projector->getDensityWeighting(); } size_t inx, iny, inz; @@ -1513,7 +1511,7 @@ static bool doJob(const CCompositeGeometryManager::TJobSet::const_iterator& iter ASTRA_DEBUG("CCompositeGeometryManager::doJobs: doing BP"); - ok = astraCUDA3d::BP(((CCompositeGeometryManager::CProjectionPart*)j.pInput.get())->pGeom, srcMem->hnd, ((CCompositeGeometryManager::CVolumePart*)j.pOutput.get())->pGeom, dstMem->hnd, voxelSuperSampling, densityWeighting); + ok = astraCUDA3d::BP(((CCompositeGeometryManager::CProjectionPart*)j.pInput.get())->pGeom, srcMem->hnd, ((CCompositeGeometryManager::CVolumePart*)j.pOutput.get())->pGeom, dstMem->hnd, voxelSuperSampling); if (!ok) ASTRA_ERROR("Error performing sub-BP"); ASTRA_DEBUG("CCompositeGeometryManager::doJobs: BP done"); } diff --git a/src/CudaProjector3D.cpp b/src/CudaProjector3D.cpp index 3ea7043..e5c55cc 100644 --- a/src/CudaProjector3D.cpp +++ b/src/CudaProjector3D.cpp @@ -67,7 +67,6 @@ void CCudaProjector3D::_clear() m_iVoxelSuperSampling = 1; m_iDetectorSuperSampling = 1; m_iGPUIndex = -1; - m_bDensityWeighting = false; } //---------------------------------------------------------------------------------------- @@ -132,13 +131,6 @@ bool CCudaProjector3D::initialize(const Config& _cfg) m_iDetectorSuperSampling = (int)_cfg.self.getOptionNumerical("DetectorSuperSampling", 1); CC.markOptionParsed("DetectorSuperSampling"); - if (dynamic_cast(m_pProjectionGeometry) || - dynamic_cast(m_pProjectionGeometry)) - { - m_bDensityWeighting = _cfg.self.getOptionBool("DensityWeighting", false); - CC.markOptionParsed("DensityWeighting"); - } - m_iGPUIndex = (int)_cfg.self.getOptionNumerical("GPUindex", -1); m_iGPUIndex = (int)_cfg.self.getOptionNumerical("GPUIndex", m_iGPUIndex); CC.markOptionParsed("GPUIndex"); -- cgit v1.2.3 From 744605ed826196aa6e0e3e3b5e945e50d830ed3a Mon Sep 17 00:00:00 2001 From: Willem Jan Palenstijn Date: Wed, 3 Apr 2019 22:32:11 +0200 Subject: Add feature flags for changed scaling behaviour --- src/Features.cpp | 6 ++++++ 1 file changed, 6 insertions(+) (limited to 'src') diff --git a/src/Features.cpp b/src/Features.cpp index 9114131..09a3499 100644 --- a/src/Features.cpp +++ b/src/Features.cpp @@ -34,6 +34,12 @@ _AstraExport bool hasFeature(const std::string &flag) { if (flag == "cuda") { return cudaEnabled(); } + if (flag == "projectors_scaled_as_line_integrals") { + return true; + } + if (flag == "fan_cone_BP_density_weighting_by_default") { + return true; + } return false; } -- cgit v1.2.3 From 9d7e93d57183d5a1e1efd113067bdcc8d388bc50 Mon Sep 17 00:00:00 2001 From: Willem Jan Palenstijn Date: Thu, 11 Apr 2019 16:52:19 +0200 Subject: Fix crash in FBP_CUDA when called without projector --- src/CudaFilteredBackProjectionAlgorithm.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) (limited to 'src') diff --git a/src/CudaFilteredBackProjectionAlgorithm.cpp b/src/CudaFilteredBackProjectionAlgorithm.cpp index 69140fc..c1d3dc8 100644 --- a/src/CudaFilteredBackProjectionAlgorithm.cpp +++ b/src/CudaFilteredBackProjectionAlgorithm.cpp @@ -152,7 +152,7 @@ void CCudaFilteredBackProjectionAlgorithm::initCUDAAlgorithm() ASTRA_ERROR("CCudaFilteredBackProjectionAlgorithm: Failed to set short-scan mode"); } - const CVolumeGeometry2D& volGeom = *m_pProjector->getVolumeGeometry(); + const CVolumeGeometry2D& volGeom = *m_pReconstruction->getGeometry(); float fPixelArea = volGeom.getPixelArea(); ok &= pFBP->setReconstructionScale(1.0f/fPixelArea); if (!ok) { -- cgit v1.2.3