diff options
-rw-r--r-- | cuda/3d/astra3d.cu | 4 | ||||
-rw-r--r-- | cuda/3d/astra3d.h | 2 | ||||
-rw-r--r-- | cuda/3d/fdk.cu | 12 | ||||
-rw-r--r-- | cuda/3d/fdk.h | 3 | ||||
-rw-r--r-- | include/astra/CudaFDKAlgorithm3D.h | 1 | ||||
-rw-r--r-- | src/CudaFDKAlgorithm3D.cpp | 33 |
6 files changed, 48 insertions, 7 deletions
diff --git a/cuda/3d/astra3d.cu b/cuda/3d/astra3d.cu index 5670873..35e3cd4 100644 --- a/cuda/3d/astra3d.cu +++ b/cuda/3d/astra3d.cu @@ -1311,7 +1311,7 @@ bool astraCudaFDK(float* pfVolume, const float* pfProjections, const CVolumeGeometry3D* pVolGeom, const CConeProjectionGeometry3D* pProjGeom, bool bShortScan, - int iGPUIndex, int iVoxelSuperSampling) + int iGPUIndex, int iVoxelSuperSampling, const float* filter) { SDimensions3D dims; @@ -1369,7 +1369,7 @@ bool astraCudaFDK(float* pfVolume, const float* pfProjections, // TODO: Offer interface for SrcZ, DetZ ok &= FDK(D_volumeData, D_projData, fOriginSourceDistance, fOriginDetectorDistance, 0, 0, fDetUSize, fDetVSize, - dims, pfAngles, bShortScan); + dims, pfAngles, bShortScan, filter); ok &= copyVolumeFromDevice(pfVolume, D_volumeData, dims, dims.iVolX); diff --git a/cuda/3d/astra3d.h b/cuda/3d/astra3d.h index 2137587..dde1347 100644 --- a/cuda/3d/astra3d.h +++ b/cuda/3d/astra3d.h @@ -314,7 +314,7 @@ _AstraExport bool astraCudaFDK(float* pfVolume, const float* pfProjections, const CVolumeGeometry3D* pVolGeom, const CConeProjectionGeometry3D* pProjGeom, bool bShortScan, - int iGPUIndex, int iVoxelSuperSampling); + int iGPUIndex, int iVoxelSuperSampling, const float* filter); } diff --git a/cuda/3d/fdk.cu b/cuda/3d/fdk.cu index 0e13be1..4899ad1 100644 --- a/cuda/3d/fdk.cu +++ b/cuda/3d/fdk.cu @@ -394,7 +394,8 @@ bool FDK(cudaPitchedPtr D_volumeData, cudaPitchedPtr D_projData, float fSrcOrigin, float fDetOrigin, float fSrcZ, float fDetZ, float fDetUSize, float fDetVSize, - const SDimensions3D& dims, const float* angles, bool bShortScan) + const SDimensions3D& dims, const float* angles, bool bShortScan, + const float* filter) { bool ok; // Generate filter @@ -412,7 +413,14 @@ bool FDK(cudaPitchedPtr D_volumeData, cufftComplex *pHostFilter = new cufftComplex[dims.iProjAngles * iHalfFFTSize]; memset(pHostFilter, 0, sizeof(cufftComplex) * dims.iProjAngles * iHalfFFTSize); - genFilter(FILTER_RAMLAK, 1.0f, dims.iProjAngles, pHostFilter, iPaddedDetCount, iHalfFFTSize); + if (filter==NULL){ + genFilter(FILTER_RAMLAK, 1.0f, dims.iProjAngles, pHostFilter, iPaddedDetCount, iHalfFFTSize); + }else{ + for(int i=0;i<dims.iProjAngles * iHalfFFTSize;i++){ + pHostFilter[i].x = filter[i]; + pHostFilter[i].y = 0; + } + } allocateComplexOnDevice(dims.iProjAngles, iHalfFFTSize, &D_filter); diff --git a/cuda/3d/fdk.h b/cuda/3d/fdk.h index de7fe53..7a9d318 100644 --- a/cuda/3d/fdk.h +++ b/cuda/3d/fdk.h @@ -43,7 +43,8 @@ bool FDK(cudaPitchedPtr D_volumeData, cudaPitchedPtr D_projData, float fSrcOrigin, float fDetOrigin, float fSrcZ, float fDetZ, float fDetUSize, float fDetVSize, - const SDimensions3D& dims, const float* angles, bool bShortScan); + const SDimensions3D& dims, const float* angles, bool bShortScan, + const float* filter); } diff --git a/include/astra/CudaFDKAlgorithm3D.h b/include/astra/CudaFDKAlgorithm3D.h index 63f07fd..477bf34 100644 --- a/include/astra/CudaFDKAlgorithm3D.h +++ b/include/astra/CudaFDKAlgorithm3D.h @@ -151,6 +151,7 @@ protected: int m_iGPUIndex; int m_iVoxelSuperSampling; + int m_iFilterDataId; bool m_bShortScan; void initializeFromProjector(); diff --git a/src/CudaFDKAlgorithm3D.cpp b/src/CudaFDKAlgorithm3D.cpp index b5ce545..bc15a3d 100644 --- a/src/CudaFDKAlgorithm3D.cpp +++ b/src/CudaFDKAlgorithm3D.cpp @@ -36,8 +36,11 @@ $Id$ #include "astra/Logging.h" #include "../cuda/3d/astra3d.h" +#include "../cuda/2d/fft.h" +#include "../cuda/3d/util3d.h" using namespace std; +using namespace astraCUDA3d; namespace astra { @@ -129,6 +132,28 @@ bool CCudaFDKAlgorithm3D::initialize(const Config& _cfg) CC.markOptionParsed("GPUIndex"); if (!_cfg.self.hasOption("GPUIndex")) CC.markOptionParsed("GPUindex"); + + // filter + if (_cfg.self.hasOption("FilterSinogramId")){ + m_iFilterDataId = (int)_cfg.self.getOptionInt("FilterSinogramId"); + const CFloat32ProjectionData2D * pFilterData = dynamic_cast<CFloat32ProjectionData2D*>(CData2DManager::getSingleton().get(m_iFilterDataId)); + if (!pFilterData){ + ASTRA_ERROR("Incorrect FilterSinogramId"); + return false; + } + const CProjectionGeometry3D* projgeom = m_pSinogram->getGeometry(); + const CProjectionGeometry2D* filtgeom = pFilterData->getGeometry(); + int iPaddedDetCount = calcNextPowerOfTwo(2 * projgeom->getDetectorColCount()); + int iHalfFFTSize = calcFFTFourSize(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; + } + }else + { + m_iFilterDataId = -1; + } + CC.markOptionParsed("FilterSinogramId"); @@ -196,10 +221,16 @@ void CCudaFDKAlgorithm3D::run(int _iNrIterations) bool ok = true; + + const float *filter = NULL; + if(m_iFilterDataId!=-1){ + const CFloat32ProjectionData2D * pFilterData = dynamic_cast<CFloat32ProjectionData2D*>(CData2DManager::getSingleton().get(m_iFilterDataId)); + filter = pFilterData->getDataConst(); + } ok = astraCudaFDK(pReconMem->getData(), pSinoMem->getDataConst(), &volgeom, conegeom, - m_bShortScan, m_iGPUIndex, m_iVoxelSuperSampling); + m_bShortScan, m_iGPUIndex, m_iVoxelSuperSampling, filter); ASTRA_ASSERT(ok); |