diff options
Diffstat (limited to 'src')
-rw-r--r-- | src/CudaFDKAlgorithm3D.cpp | 4 | ||||
-rw-r--r-- | src/CudaFilteredBackProjectionAlgorithm.cpp | 2 | ||||
-rw-r--r-- | src/Filters.cpp | 66 |
3 files changed, 67 insertions, 5 deletions
diff --git a/src/CudaFDKAlgorithm3D.cpp b/src/CudaFDKAlgorithm3D.cpp index 4d8fc5a..24ed04f 100644 --- a/src/CudaFDKAlgorithm3D.cpp +++ b/src/CudaFDKAlgorithm3D.cpp @@ -35,9 +35,9 @@ along with the ASTRA Toolbox. If not, see <http://www.gnu.org/licenses/>. #include "astra/CompositeGeometryManager.h" #include "astra/Logging.h" +#include "astra/Filters.h" #include "astra/cuda/3d/astra3d.h" -#include "astra/cuda/2d/fft.h" #include "astra/cuda/3d/util3d.h" using namespace std; @@ -156,7 +156,7 @@ bool CCudaFDKAlgorithm3D::initialize(const Config& _cfg) const CProjectionGeometry3D* projgeom = m_pSinogram->getGeometry(); const CProjectionGeometry2D* filtgeom = pFilterData->getGeometry(); int iPaddedDetCount = calcNextPowerOfTwo(2 * projgeom->getDetectorColCount()); - int iHalfFFTSize = astraCUDA::calcFFTFourierSize(iPaddedDetCount); + int iHalfFFTSize = calcFFTFourierSize(iPaddedDetCount); if(filtgeom->getDetectorCount()!=iHalfFFTSize || filtgeom->getProjectionAngleCount()!=projgeom->getProjectionCount()){ ASTRA_ERROR("Filter size does not match required size (%i angles, %i detectors)",projgeom->getProjectionCount(),iHalfFFTSize); return false; diff --git a/src/CudaFilteredBackProjectionAlgorithm.cpp b/src/CudaFilteredBackProjectionAlgorithm.cpp index 90b831e..88e235b 100644 --- a/src/CudaFilteredBackProjectionAlgorithm.cpp +++ b/src/CudaFilteredBackProjectionAlgorithm.cpp @@ -176,6 +176,8 @@ bool CCudaFilteredBackProjectionAlgorithm::check() // check pixel supersampling ASTRA_CONFIG_CHECK(m_iPixelSuperSampling >= 0, "FBP_CUDA", "PixelSuperSampling must be a non-negative integer."); + ASTRA_CONFIG_CHECK(checkCustomFilterSize(m_filterConfig, *m_pSinogram->getGeometry()), "FBP_CUDA", "Filter size mismatch"); + // success m_bIsInitialized = true; diff --git a/src/Filters.cpp b/src/Filters.cpp index bbd7cfd..3ee3749 100644 --- a/src/Filters.cpp +++ b/src/Filters.cpp @@ -504,14 +504,15 @@ SFilterConfig getFilterConfigForAlgorithm(const Config& _cfg, CAlgorithm *_alg) int id = node.getContentInt(); const CFloat32ProjectionData2D * pFilterData = dynamic_cast<CFloat32ProjectionData2D*>(CData2DManager::getSingleton().get(id)); c.m_iCustomFilterWidth = pFilterData->getGeometry()->getDetectorCount(); - int iFilterProjectionCount = pFilterData->getGeometry()->getProjectionAngleCount(); + c.m_iCustomFilterHeight = pFilterData->getGeometry()->getProjectionAngleCount(); - c.m_pfCustomFilter = new float[c.m_iCustomFilterWidth * iFilterProjectionCount]; - memcpy(c.m_pfCustomFilter, pFilterData->getDataConst(), sizeof(float) * c.m_iCustomFilterWidth * iFilterProjectionCount); + c.m_pfCustomFilter = new float[c.m_iCustomFilterWidth * c.m_iCustomFilterHeight]; + memcpy(c.m_pfCustomFilter, pFilterData->getDataConst(), sizeof(float) * c.m_iCustomFilterWidth * c.m_iCustomFilterHeight); } else { c.m_iCustomFilterWidth = 0; + c.m_iCustomFilterHeight = 0; c.m_pfCustomFilter = NULL; } CC.markNodeParsed("FilterSinogramId"); // TODO: Only for some types! @@ -545,4 +546,63 @@ SFilterConfig getFilterConfigForAlgorithm(const Config& _cfg, CAlgorithm *_alg) return c; } +static int calcNextPowerOfTwo(int n) +{ + int x = 1; + while (x < n) + x *= 2; + + return x; +} + +// Because the input is real, the Fourier transform is symmetric. +// CUFFT only outputs the first half (ignoring the redundant second half), +// and expects the same as input for the IFFT. +int calcFFTFourierSize(int _iFFTRealSize) +{ + int iFFTFourierSize = _iFFTRealSize / 2 + 1; + + return iFFTFourierSize; +} + +bool checkCustomFilterSize(const SFilterConfig &_cfg, const CProjectionGeometry2D &_geom) { + int iExpectedWidth = -1, iExpectedHeight = 1; + + switch (_cfg.m_eType) { + case FILTER_ERROR: + ASTRA_ERROR("checkCustomFilterSize: internal error; FILTER_ERROR passed"); + return false; + case FILTER_NONE: + return true; + case FILTER_SINOGRAM: + iExpectedHeight = _geom.getProjectionAngleCount(); + // fallthrough + case FILTER_PROJECTION: + { + int iPaddedDetCount = calcNextPowerOfTwo(2 * _geom.getDetectorCount()); + iExpectedWidth = calcFFTFourierSize(iPaddedDetCount); + } + if (_cfg.m_iCustomFilterWidth != iExpectedWidth || + _cfg.m_iCustomFilterHeight != iExpectedHeight) + { + ASTRA_ERROR("filter size mismatch: %dx%d (received) is not %dx%d (expected)", _cfg.m_iCustomFilterHeight, _cfg.m_iCustomFilterWidth, iExpectedHeight, iExpectedWidth); + return false; + } + return true; + case FILTER_RSINOGRAM: + iExpectedHeight = _geom.getProjectionAngleCount(); + // fallthrough + case FILTER_RPROJECTION: + if (_cfg.m_iCustomFilterHeight != iExpectedHeight) + { + ASTRA_ERROR("filter size mismatch: %dx%d (received) is not %dxX (expected)", _cfg.m_iCustomFilterHeight, _cfg.m_iCustomFilterWidth, iExpectedHeight); + return false; + } + return true; + default: + // Non-custom filters; nothing to check. + return true; + } +} + } |