summaryrefslogtreecommitdiffstats
diff options
context:
space:
mode:
-rw-r--r--cuda/3d/astra3d.cu2
-rw-r--r--cuda/3d/fdk.cu12
-rw-r--r--cuda/3d/fdk.h3
-rw-r--r--cuda/3d/mem3d.cu4
-rw-r--r--cuda/3d/mem3d.h2
-rw-r--r--include/astra/CompositeGeometryManager.h4
-rw-r--r--include/astra/CudaFDKAlgorithm3D.h1
-rw-r--r--src/CompositeGeometryManager.cpp6
-rw-r--r--src/CudaFDKAlgorithm3D.cpp37
9 files changed, 57 insertions, 14 deletions
diff --git a/cuda/3d/astra3d.cu b/cuda/3d/astra3d.cu
index 91f4530..3bd56ea 100644
--- a/cuda/3d/astra3d.cu
+++ b/cuda/3d/astra3d.cu
@@ -1289,6 +1289,4 @@ bool astraCudaBP_SIRTWeighted(float* pfVolume,
}
-
-
}
diff --git a/cuda/3d/fdk.cu b/cuda/3d/fdk.cu
index 19ef338..bbac0be 100644
--- a/cuda/3d/fdk.cu
+++ b/cuda/3d/fdk.cu
@@ -281,7 +281,8 @@ bool FDK_Filter(cudaPitchedPtr D_projData,
bool FDK(cudaPitchedPtr D_volumeData,
cudaPitchedPtr D_projData,
const SConeProjection* angles,
- const SDimensions3D& dims, SProjectorParams3D params, bool bShortScan)
+ const SDimensions3D& dims, SProjectorParams3D params, bool bShortScan,
+ const float* pfFilter)
{
bool ok;
// Generate filter
@@ -332,7 +333,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 (pfFilter == 0){
+ genFilter(FILTER_RAMLAK, 1.0f, dims.iProjAngles, pHostFilter, iPaddedDetCount, iHalfFFTSize);
+ } else {
+ for (int i = 0; i < dims.iProjAngles * iHalfFFTSize; i++) {
+ pHostFilter[i].x = pfFilter[i];
+ pHostFilter[i].y = 0;
+ }
+ }
allocateComplexOnDevice(dims.iProjAngles, iHalfFFTSize, &D_filter);
diff --git a/cuda/3d/fdk.h b/cuda/3d/fdk.h
index 8d2a128..fc0df20 100644
--- a/cuda/3d/fdk.h
+++ b/cuda/3d/fdk.h
@@ -42,7 +42,8 @@ bool FDK_PreWeight(cudaPitchedPtr D_projData,
bool FDK(cudaPitchedPtr D_volumeData,
cudaPitchedPtr D_projData,
const SConeProjection* angles,
- const SDimensions3D& dims, SProjectorParams3D params, bool bShortScan);
+ const SDimensions3D& dims, SProjectorParams3D params, bool bShortScan,
+ const float* filter);
}
diff --git a/cuda/3d/mem3d.cu b/cuda/3d/mem3d.cu
index 4fb30a2..cb15ab9 100644
--- a/cuda/3d/mem3d.cu
+++ b/cuda/3d/mem3d.cu
@@ -279,7 +279,7 @@ bool BP(const astra::CProjectionGeometry3D* pProjGeom, MemHandle3D projData, con
}
-bool FDK(const astra::CProjectionGeometry3D* pProjGeom, MemHandle3D projData, const astra::CVolumeGeometry3D* pVolGeom, MemHandle3D volData, bool bShortScan)
+bool FDK(const astra::CProjectionGeometry3D* pProjGeom, MemHandle3D projData, const astra::CVolumeGeometry3D* pVolGeom, MemHandle3D volData, bool bShortScan, const float *pfFilter)
{
SDimensions3D dims;
SProjectorParams3D params;
@@ -298,7 +298,7 @@ bool FDK(const astra::CProjectionGeometry3D* pProjGeom, MemHandle3D projData, co
if (!ok || !pConeProjs)
return false;
- ok &= FDK(volData.d->ptr, projData.d->ptr, pConeProjs, dims, params, bShortScan);
+ ok &= FDK(volData.d->ptr, projData.d->ptr, pConeProjs, dims, params, bShortScan, pfFilter);
return ok;
diff --git a/cuda/3d/mem3d.h b/cuda/3d/mem3d.h
index 2d0fb27..bb8b091 100644
--- a/cuda/3d/mem3d.h
+++ b/cuda/3d/mem3d.h
@@ -95,7 +95,7 @@ bool FP(const astra::CProjectionGeometry3D* pProjGeom, MemHandle3D projData, con
bool BP(const astra::CProjectionGeometry3D* pProjGeom, MemHandle3D projData, const astra::CVolumeGeometry3D* pVolGeom, MemHandle3D volData, int iVoxelSuperSampling);
-bool FDK(const astra::CProjectionGeometry3D* pProjGeom, MemHandle3D projData, const astra::CVolumeGeometry3D* pVolGeom, MemHandle3D volData, bool bShortScan);
+bool FDK(const astra::CProjectionGeometry3D* pProjGeom, MemHandle3D projData, const astra::CVolumeGeometry3D* pVolGeom, MemHandle3D volData, bool bShortScan, const float *pfFilter = 0);
}
diff --git a/include/astra/CompositeGeometryManager.h b/include/astra/CompositeGeometryManager.h
index 064370a..d30ffd1 100644
--- a/include/astra/CompositeGeometryManager.h
+++ b/include/astra/CompositeGeometryManager.h
@@ -57,6 +57,7 @@ struct SGPUParams {
struct SFDKSettings {
bool bShortScan;
+ const float *pfFilter;
};
@@ -161,7 +162,8 @@ public:
bool doBP(CProjector3D *pProjector, CFloat32VolumeData3DMemory *pVolData,
CFloat32ProjectionData3DMemory *pProjData);
bool doFDK(CProjector3D *pProjector, CFloat32VolumeData3DMemory *pVolData,
- CFloat32ProjectionData3DMemory *pProjData, bool bShortScan);
+ CFloat32ProjectionData3DMemory *pProjData, bool bShortScan,
+ const float *pfFilter = 0);
bool doFP(CProjector3D *pProjector, const std::vector<CFloat32VolumeData3DMemory *>& volData, const std::vector<CFloat32ProjectionData3DMemory *>& projData);
bool doBP(CProjector3D *pProjector, const std::vector<CFloat32VolumeData3DMemory *>& volData, const std::vector<CFloat32ProjectionData3DMemory *>& projData);
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/CompositeGeometryManager.cpp b/src/CompositeGeometryManager.cpp
index 7c4f8e6..5879aec 100644
--- a/src/CompositeGeometryManager.cpp
+++ b/src/CompositeGeometryManager.cpp
@@ -1042,7 +1042,8 @@ bool CCompositeGeometryManager::doBP(CProjector3D *pProjector, CFloat32VolumeDat
bool CCompositeGeometryManager::doFDK(CProjector3D *pProjector, CFloat32VolumeData3DMemory *pVolData,
- CFloat32ProjectionData3DMemory *pProjData, bool bShortScan)
+ CFloat32ProjectionData3DMemory *pProjData, bool bShortScan,
+ const float *pfFilter)
{
if (!dynamic_cast<CConeProjectionGeometry3D*>(pProjData->getGeometry())) {
ASTRA_ERROR("CCompositeGeometryManager::doFDK: cone geometry required");
@@ -1052,6 +1053,7 @@ bool CCompositeGeometryManager::doFDK(CProjector3D *pProjector, CFloat32VolumeDa
SJob job = createJobBP(pProjector, pVolData, pProjData);
job.eType = SJob::JOB_FDK;
job.FDKSettings.bShortScan = bShortScan;
+ job.FDKSettings.pfFilter = pfFilter;
TJobList L;
L.push_back(job);
@@ -1288,7 +1290,7 @@ static bool doJob(const CCompositeGeometryManager::TJobSet::const_iterator& iter
} else {
ASTRA_DEBUG("CCompositeGeometryManager::doJobs: doing FDK");
- ok = astraCUDA3d::FDK(((CCompositeGeometryManager::CProjectionPart*)j.pInput.get())->pGeom, inputMem, ((CCompositeGeometryManager::CVolumePart*)j.pOutput.get())->pGeom, outputMem, j.FDKSettings.bShortScan);
+ ok = astraCUDA3d::FDK(((CCompositeGeometryManager::CProjectionPart*)j.pInput.get())->pGeom, inputMem, ((CCompositeGeometryManager::CVolumePart*)j.pOutput.get())->pGeom, outputMem, j.FDKSettings.bShortScan, j.FDKSettings.pfFilter);
if (!ok) ASTRA_ERROR("Error performing sub-FDK");
ASTRA_DEBUG("CCompositeGeometryManager::doJobs: FDK done");
}
diff --git a/src/CudaFDKAlgorithm3D.cpp b/src/CudaFDKAlgorithm3D.cpp
index c7c8ed5..15abed4 100644
--- a/src/CudaFDKAlgorithm3D.cpp
+++ b/src/CudaFDKAlgorithm3D.cpp
@@ -37,8 +37,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 {
@@ -141,6 +144,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");
@@ -206,20 +231,26 @@ void CCudaFDKAlgorithm3D::run(int _iNrIterations)
CFloat32VolumeData3DMemory* pReconMem = dynamic_cast<CFloat32VolumeData3DMemory*>(m_pReconstruction);
ASTRA_ASSERT(pReconMem);
+ const float *filter = NULL;
+ if (m_iFilterDataId != -1) {
+ const CFloat32ProjectionData2D *pFilterData = dynamic_cast<CFloat32ProjectionData2D*>(CData2DManager::getSingleton().get(m_iFilterDataId));
+ if (pFilterData)
+ filter = pFilterData->getDataConst();
+ }
#if 0
bool ok = true;
-
+
ok = astraCudaFDK(pReconMem->getData(), pSinoMem->getDataConst(),
&volgeom, conegeom,
- m_bShortScan, m_iGPUIndex, m_iVoxelSuperSampling);
+ m_bShortScan, m_iGPUIndex, m_iVoxelSuperSampling, filter);
ASTRA_ASSERT(ok);
#endif
CCompositeGeometryManager cgm;
- cgm.doFDK(m_pProjector, pReconMem, pSinoMem, m_bShortScan);
+ cgm.doFDK(m_pProjector, pReconMem, pSinoMem, m_bShortScan, filter);