summaryrefslogtreecommitdiffstats
diff options
context:
space:
mode:
authorWillem Jan Palenstijn <Willem.Jan.Palenstijn@cwi.nl>2015-03-09 15:43:56 +0100
committerWillem Jan Palenstijn <Willem.Jan.Palenstijn@cwi.nl>2015-03-09 15:43:56 +0100
commit5304d08cd1ab7b8d778c367912934376eb92370f (patch)
treef0af13b4d9320fdd3af12f36a14e66feffebb90e
parent248639b4fee8659a4106dcc44d721149a1885018 (diff)
downloadastra-5304d08cd1ab7b8d778c367912934376eb92370f.tar.gz
astra-5304d08cd1ab7b8d778c367912934376eb92370f.tar.bz2
astra-5304d08cd1ab7b8d778c367912934376eb92370f.tar.xz
astra-5304d08cd1ab7b8d778c367912934376eb92370f.zip
Allow non-centered volume geometry in SIRT3D and CGLS3D
-rw-r--r--cuda/3d/astra3d.cu284
-rw-r--r--cuda/3d/astra3d.h90
-rw-r--r--src/CudaCglsAlgorithm3D.cpp39
-rw-r--r--src/CudaSirtAlgorithm3D.cpp38
4 files changed, 90 insertions, 361 deletions
diff --git a/cuda/3d/astra3d.cu b/cuda/3d/astra3d.cu
index f672d6c..426f3a0 100644
--- a/cuda/3d/astra3d.cu
+++ b/cuda/3d/astra3d.cu
@@ -182,6 +182,20 @@ static bool convertAstraGeometry_internal(const CVolumeGeometry3D* pVolGeom,
}
+void convertAstraGeometry_dims(const CVolumeGeometry3D* pVolGeom,
+ const CProjectionGeometry3D* pProjGeom,
+ SDimensions3D& dims)
+{
+ dims.iVolX = pVolGeom->getGridColCount();
+ dims.iVolY = pVolGeom->getGridRowCount();
+ dims.iVolZ = pVolGeom->getGridSliceCount();
+ dims.iProjAngles = pProjGeom->getProjectionCount();
+ dims.iProjU = pProjGeom->getDetectorColCount(),
+ dims.iProjV = pProjGeom->getDetectorRowCount(),
+ dims.iRaysPerDetDim = 1;
+ dims.iRaysPerVoxelDim = 1;
+}
+
bool convertAstraGeometry(const CVolumeGeometry3D* pVolGeom,
const CParallelProjectionGeometry3D* pProjGeom,
@@ -370,127 +384,55 @@ AstraSIRT3d::~AstraSIRT3d()
pData = 0;
}
-bool AstraSIRT3d::setReconstructionGeometry(unsigned int iVolX,
- unsigned int iVolY,
- unsigned int iVolZ/*,
- float fPixelSize = 1.0f*/)
+bool AstraSIRT3d::setGeometry(const CVolumeGeometry3D* pVolGeom,
+ const CProjectionGeometry3D* pProjGeom)
{
if (pData->initialized)
return false;
- pData->dims.iVolX = iVolX;
- pData->dims.iVolY = iVolY;
- pData->dims.iVolZ = iVolZ;
+ convertAstraGeometry_dims(pVolGeom, pProjGeom, pData->dims);
- return (iVolX > 0 && iVolY > 0 && iVolZ > 0);
-}
-
-
-bool AstraSIRT3d::setPar3DGeometry(unsigned int iProjAngles,
- unsigned int iProjU,
- unsigned int iProjV,
- const SPar3DProjection* projs)
-{
- if (pData->initialized)
+ if (pData->dims.iVolX <= 0 || pData->dims.iVolX <= 0 || pData->dims.iVolX <= 0)
return false;
-
- pData->dims.iProjAngles = iProjAngles;
- pData->dims.iProjU = iProjU;
- pData->dims.iProjV = iProjV;
-
- if (iProjAngles == 0 || iProjU == 0 || iProjV == 0 || projs == 0)
+ if (pData->dims.iProjAngles <= 0 || pData->dims.iProjU <= 0 || pData->dims.iProjV <= 0)
return false;
- pData->parprojs = new SPar3DProjection[iProjAngles];
- memcpy(pData->parprojs, projs, iProjAngles * sizeof(projs[0]));
-
- pData->projType = PROJ_PARALLEL;
-
- return true;
-}
-
-bool AstraSIRT3d::setPar3DGeometry(unsigned int iProjAngles,
- unsigned int iProjU,
- unsigned int iProjV,
- float fDetUSize,
- float fDetVSize,
- const float *pfAngles)
-{
- if (pData->initialized)
- return false;
-
- if (iProjAngles == 0 || iProjU == 0 || iProjV == 0 || pfAngles == 0)
- return false;
-
- SPar3DProjection* p = genPar3DProjections(iProjAngles,
- iProjU, iProjV,
- fDetUSize, fDetVSize,
- pfAngles);
- pData->dims.iProjAngles = iProjAngles;
- pData->dims.iProjU = iProjU;
- pData->dims.iProjV = iProjV;
-
- pData->parprojs = p;
- pData->projType = PROJ_PARALLEL;
-
- return true;
-}
-
-
+ const CConeProjectionGeometry3D* conegeom = dynamic_cast<const CConeProjectionGeometry3D*>(pProjGeom);
+ const CParallelProjectionGeometry3D* par3dgeom = dynamic_cast<const CParallelProjectionGeometry3D*>(pProjGeom);
+ const CParallelVecProjectionGeometry3D* parvec3dgeom = dynamic_cast<const CParallelVecProjectionGeometry3D*>(pProjGeom);
+ const CConeVecProjectionGeometry3D* conevec3dgeom = dynamic_cast<const CConeVecProjectionGeometry3D*>(pProjGeom);
-bool AstraSIRT3d::setConeGeometry(unsigned int iProjAngles,
- unsigned int iProjU,
- unsigned int iProjV,
- const SConeProjection* projs)
-{
- if (pData->initialized)
- return false;
+ float outputScale;
+ bool ok;
- pData->dims.iProjAngles = iProjAngles;
- pData->dims.iProjU = iProjU;
- pData->dims.iProjV = iProjV;
+ pData->projs = 0;
+ pData->parprojs = 0;
+
+ if (conegeom) {
+ ok = convertAstraGeometry(pVolGeom, conegeom, pData->projs, outputScale);
+ pData->projType = PROJ_PARALLEL;
+ } else if (conevec3dgeom) {
+ ok = convertAstraGeometry(pVolGeom, conevec3dgeom, pData->projs, outputScale);
+ pData->projType = PROJ_PARALLEL;
+ } else if (par3dgeom) {
+ ok = convertAstraGeometry(pVolGeom, par3dgeom, pData->parprojs, outputScale);
+ pData->projType = PROJ_CONE;
+ } else if (parvec3dgeom) {
+ ok = convertAstraGeometry(pVolGeom, parvec3dgeom, pData->parprojs, outputScale);
+ pData->projType = PROJ_CONE;
+ } else {
+ ok = false;
+ }
- if (iProjAngles == 0 || iProjU == 0 || iProjV == 0 || projs == 0)
+ if (!ok)
return false;
- pData->projs = new SConeProjection[iProjAngles];
- memcpy(pData->projs, projs, iProjAngles * sizeof(projs[0]));
- pData->projType = PROJ_CONE;
+ // TODO: Handle outputScale
return true;
}
-bool AstraSIRT3d::setConeGeometry(unsigned int iProjAngles,
- unsigned int iProjU,
- unsigned int iProjV,
- float fOriginSourceDistance,
- float fOriginDetectorDistance,
- float fDetUSize,
- float fDetVSize,
- const float *pfAngles)
-{
- if (pData->initialized)
- return false;
-
- if (iProjAngles == 0 || iProjU == 0 || iProjV == 0 || pfAngles == 0)
- return false;
-
- SConeProjection* p = genConeProjections(iProjAngles,
- iProjU, iProjV,
- fOriginSourceDistance,
- fOriginDetectorDistance,
- fDetUSize, fDetVSize,
- pfAngles);
- pData->dims.iProjAngles = iProjAngles;
- pData->dims.iProjU = iProjU;
- pData->dims.iProjV = iProjV;
-
- pData->projs = p;
- pData->projType = PROJ_CONE;
-
- return true;
-}
bool AstraSIRT3d::enableSuperSampling(unsigned int iVoxelSuperSampling,
unsigned int iDetectorSuperSampling)
@@ -837,125 +779,51 @@ AstraCGLS3d::~AstraCGLS3d()
pData = 0;
}
-bool AstraCGLS3d::setReconstructionGeometry(unsigned int iVolX,
- unsigned int iVolY,
- unsigned int iVolZ/*,
- float fPixelSize = 1.0f*/)
+bool AstraCGLS3d::setGeometry(const CVolumeGeometry3D* pVolGeom,
+ const CProjectionGeometry3D* pProjGeom)
{
if (pData->initialized)
return false;
- pData->dims.iVolX = iVolX;
- pData->dims.iVolY = iVolY;
- pData->dims.iVolZ = iVolZ;
+ convertAstraGeometry_dims(pVolGeom, pProjGeom, pData->dims);
- return (iVolX > 0 && iVolY > 0 && iVolZ > 0);
-}
-
-
-bool AstraCGLS3d::setPar3DGeometry(unsigned int iProjAngles,
- unsigned int iProjU,
- unsigned int iProjV,
- const SPar3DProjection* projs)
-{
- if (pData->initialized)
+ if (pData->dims.iVolX <= 0 || pData->dims.iVolX <= 0 || pData->dims.iVolX <= 0)
return false;
-
- pData->dims.iProjAngles = iProjAngles;
- pData->dims.iProjU = iProjU;
- pData->dims.iProjV = iProjV;
-
- if (iProjAngles == 0 || iProjU == 0 || iProjV == 0 || projs == 0)
+ if (pData->dims.iProjAngles <= 0 || pData->dims.iProjU <= 0 || pData->dims.iProjV <= 0)
return false;
- pData->parprojs = new SPar3DProjection[iProjAngles];
- memcpy(pData->parprojs, projs, iProjAngles * sizeof(projs[0]));
+ const CConeProjectionGeometry3D* conegeom = dynamic_cast<const CConeProjectionGeometry3D*>(pProjGeom);
+ const CParallelProjectionGeometry3D* par3dgeom = dynamic_cast<const CParallelProjectionGeometry3D*>(pProjGeom);
+ const CParallelVecProjectionGeometry3D* parvec3dgeom = dynamic_cast<const CParallelVecProjectionGeometry3D*>(pProjGeom);
+ const CConeVecProjectionGeometry3D* conevec3dgeom = dynamic_cast<const CConeVecProjectionGeometry3D*>(pProjGeom);
- pData->projType = PROJ_PARALLEL;
-
- return true;
-}
-
-bool AstraCGLS3d::setPar3DGeometry(unsigned int iProjAngles,
- unsigned int iProjU,
- unsigned int iProjV,
- float fDetUSize,
- float fDetVSize,
- const float *pfAngles)
-{
- if (pData->initialized)
- return false;
-
- if (iProjAngles == 0 || iProjU == 0 || iProjV == 0 || pfAngles == 0)
- return false;
-
- SPar3DProjection* p = genPar3DProjections(iProjAngles,
- iProjU, iProjV,
- fDetUSize, fDetVSize,
- pfAngles);
- pData->dims.iProjAngles = iProjAngles;
- pData->dims.iProjU = iProjU;
- pData->dims.iProjV = iProjV;
-
- pData->parprojs = p;
- pData->projType = PROJ_PARALLEL;
-
- return true;
-}
-
-
-
-bool AstraCGLS3d::setConeGeometry(unsigned int iProjAngles,
- unsigned int iProjU,
- unsigned int iProjV,
- const SConeProjection* projs)
-{
- if (pData->initialized)
- return false;
-
- pData->dims.iProjAngles = iProjAngles;
- pData->dims.iProjU = iProjU;
- pData->dims.iProjV = iProjV;
-
- if (iProjAngles == 0 || iProjU == 0 || iProjV == 0 || projs == 0)
- return false;
-
- pData->projs = new SConeProjection[iProjAngles];
- memcpy(pData->projs, projs, iProjAngles * sizeof(projs[0]));
-
- pData->projType = PROJ_CONE;
-
- return true;
-}
+ float outputScale;
+ bool ok;
-bool AstraCGLS3d::setConeGeometry(unsigned int iProjAngles,
- unsigned int iProjU,
- unsigned int iProjV,
- float fOriginSourceDistance,
- float fOriginDetectorDistance,
- float fDetUSize,
- float fDetVSize,
- const float *pfAngles)
-{
- if (pData->initialized)
- return false;
+ pData->projs = 0;
+ pData->parprojs = 0;
+
+ if (conegeom) {
+ ok = convertAstraGeometry(pVolGeom, conegeom, pData->projs, outputScale);
+ pData->projType = PROJ_PARALLEL;
+ } else if (conevec3dgeom) {
+ ok = convertAstraGeometry(pVolGeom, conevec3dgeom, pData->projs, outputScale);
+ pData->projType = PROJ_PARALLEL;
+ } else if (par3dgeom) {
+ ok = convertAstraGeometry(pVolGeom, par3dgeom, pData->parprojs, outputScale);
+ pData->projType = PROJ_CONE;
+ } else if (parvec3dgeom) {
+ ok = convertAstraGeometry(pVolGeom, parvec3dgeom, pData->parprojs, outputScale);
+ pData->projType = PROJ_CONE;
+ } else {
+ ok = false;
+ }
- if (iProjAngles == 0 || iProjU == 0 || iProjV == 0 || pfAngles == 0)
+ if (!ok)
return false;
- SConeProjection* p = genConeProjections(iProjAngles,
- iProjU, iProjV,
- fOriginSourceDistance,
- fOriginDetectorDistance,
- fDetUSize, fDetVSize,
- pfAngles);
-
- pData->dims.iProjAngles = iProjAngles;
- pData->dims.iProjU = iProjU;
- pData->dims.iProjV = iProjV;
- pData->projs = p;
- pData->projType = PROJ_CONE;
+ // TODO: Handle outputScale
return true;
}
diff --git a/cuda/3d/astra3d.h b/cuda/3d/astra3d.h
index 47e252e..cab5479 100644
--- a/cuda/3d/astra3d.h
+++ b/cuda/3d/astra3d.h
@@ -42,7 +42,12 @@ enum Cuda3DProjectionKernel {
ker3d_sum_square_weights
};
-
+class CProjectionGeometry3D;
+class CParallelProjectionGeometry3D;
+class CParallelVecProjectionGeometry3D;
+class CConeProjectionGeometry3D;
+class CConeVecProjectionGeometry3D;
+class CVolumeGeometry3D;
class AstraSIRT3d_internal;
@@ -52,37 +57,9 @@ public:
AstraSIRT3d();
~AstraSIRT3d();
- // Set the number of pixels in the reconstruction rectangle,
- // and the length of the edge of a pixel.
- // Volume pixels are assumed to be square.
- // This must be called before setting the projection geometry.
- bool setReconstructionGeometry(unsigned int iVolX,
- unsigned int iVolY,
- unsigned int iVolZ/*,
- float fPixelSize = 1.0f*/);
-
- bool setConeGeometry(unsigned int iProjAngles,
- unsigned int iProjU,
- unsigned int iProjV,
- const SConeProjection* projs);
- bool setConeGeometry(unsigned int iProjAngles,
- unsigned int iProjU,
- unsigned int iProjV,
- float fOriginSourceDistance,
- float fOriginDetectorDistance,
- float fSourceZ,
- float fDetSize,
- const float *pfAngles);
- bool setPar3DGeometry(unsigned int iProjAngles,
- unsigned int iProjU,
- unsigned int iProjV,
- const SPar3DProjection* projs);
- bool setPar3DGeometry(unsigned int iProjAngles,
- unsigned int iProjU,
- unsigned int iProjV,
- float fSourceZ,
- float fDetSize,
- const float *pfAngles);
+ // Set the volume and projection geometry
+ bool setGeometry(const CVolumeGeometry3D* pVolGeom,
+ const CProjectionGeometry3D* pProjGeom);
// Enable supersampling.
//
@@ -197,37 +174,9 @@ public:
AstraCGLS3d();
~AstraCGLS3d();
- // Set the number of pixels in the reconstruction rectangle,
- // and the length of the edge of a pixel.
- // Volume pixels are assumed to be square.
- // This must be called before setting the projection geometry.
- bool setReconstructionGeometry(unsigned int iVolX,
- unsigned int iVolY,
- unsigned int iVolZ/*,
- float fPixelSize = 1.0f*/);
-
- bool setConeGeometry(unsigned int iProjAngles,
- unsigned int iProjU,
- unsigned int iProjV,
- const SConeProjection* projs);
- bool setConeGeometry(unsigned int iProjAngles,
- unsigned int iProjU,
- unsigned int iProjV,
- float fOriginSourceDistance,
- float fOriginDetectorDistance,
- float fSourceZ,
- float fDetSize,
- const float *pfAngles);
- bool setPar3DGeometry(unsigned int iProjAngles,
- unsigned int iProjU,
- unsigned int iProjV,
- const SPar3DProjection* projs);
- bool setPar3DGeometry(unsigned int iProjAngles,
- unsigned int iProjU,
- unsigned int iProjV,
- float fSourceZ,
- float fDetSize,
- const float *pfAngles);
+ // Set the volume and projection geometry
+ bool setGeometry(const CVolumeGeometry3D* pVolGeom,
+ const CProjectionGeometry3D* pProjGeom);
// Enable supersampling.
//
@@ -466,21 +415,6 @@ _AstraExport bool astraCudaFDK(float* pfVolume, const float* pfProjections,
bool bShortScan,
int iGPUIndex, int iVoxelSuperSampling);
-_AstraExport bool convertAstraGeometry(const CVolumeGeometry3D* pVolGeom,
- const CParallelProjectionGeometry3D* pProjGeom,
- SPar3DProjection*& pProjs, float& fOutputScale);
-
-_AstraExport bool convertAstraGeometry(const CVolumeGeometry3D* pVolGeom,
- const CParallelVecProjectionGeometry3D* pProjGeom,
- SPar3DProjection*& pProjs, float& fOutputScale);
-
-_AstraExport bool convertAstraGeometry(const CVolumeGeometry3D* pVolGeom,
- const CConeProjectionGeometry3D* pProjGeom,
- SConeProjection*& pProjs, float& fOutputScale);
-
-_AstraExport bool convertAstraGeometry(const CVolumeGeometry3D* pVolGeom,
- const CConeVecProjectionGeometry3D* pProjGeom,
- SConeProjection*& pProjs, float& fOutputScale);
}
diff --git a/src/CudaCglsAlgorithm3D.cpp b/src/CudaCglsAlgorithm3D.cpp
index a5500d6..3677458 100644
--- a/src/CudaCglsAlgorithm3D.cpp
+++ b/src/CudaCglsAlgorithm3D.cpp
@@ -171,9 +171,6 @@ void CCudaCglsAlgorithm3D::run(int _iNrIterations)
ASTRA_ASSERT(m_bIsInitialized);
const CProjectionGeometry3D* projgeom = m_pSinogram->getGeometry();
- const CConeProjectionGeometry3D* conegeom = dynamic_cast<const CConeProjectionGeometry3D*>(projgeom);
- const CParallelVecProjectionGeometry3D* parvec3dgeom = dynamic_cast<const CParallelVecProjectionGeometry3D*>(projgeom);
- const CConeVecProjectionGeometry3D* conevec3dgeom = dynamic_cast<const CConeVecProjectionGeometry3D*>(projgeom);
const CVolumeGeometry3D& volgeom = *m_pReconstruction->getGeometry();
bool ok = true;
@@ -182,41 +179,7 @@ void CCudaCglsAlgorithm3D::run(int _iNrIterations)
ok &= m_pCgls->setGPUIndex(m_iGPUIndex);
- ok &= m_pCgls->setReconstructionGeometry(volgeom.getGridColCount(),
- volgeom.getGridRowCount(),
- volgeom.getGridSliceCount());
-/*
- unsigned int iProjAngles,
- unsigned int iProjU,
- unsigned int iProjV,
- float fOriginSourceDistance,
- float fOriginDetectorDistance,
- float fDetUSize,
- float fDetVSize,
- const float *pfAngles)
-*/
- if (conegeom) {
- ok &= m_pCgls->setConeGeometry(conegeom->getProjectionCount(),
- conegeom->getDetectorColCount(),
- conegeom->getDetectorRowCount(),
- conegeom->getOriginSourceDistance(),
- conegeom->getOriginDetectorDistance(),
- conegeom->getDetectorSpacingX(),
- conegeom->getDetectorSpacingY(),
- conegeom->getProjectionAngles());
- } else if (parvec3dgeom) {
- ok &= m_pCgls->setPar3DGeometry(parvec3dgeom->getProjectionCount(),
- parvec3dgeom->getDetectorColCount(),
- parvec3dgeom->getDetectorRowCount(),
- parvec3dgeom->getProjectionVectors());
- } else if (conevec3dgeom) {
- ok &= m_pCgls->setConeGeometry(conevec3dgeom->getProjectionCount(),
- conevec3dgeom->getDetectorColCount(),
- conevec3dgeom->getDetectorRowCount(),
- conevec3dgeom->getProjectionVectors());
- } else {
- ASTRA_ASSERT(false);
- }
+ ok &= m_pCgls->setGeometry(&volgeom, projgeom);
ok &= m_pCgls->enableSuperSampling(m_iVoxelSuperSampling, m_iDetectorSuperSampling);
diff --git a/src/CudaSirtAlgorithm3D.cpp b/src/CudaSirtAlgorithm3D.cpp
index da83c7e..d67778f 100644
--- a/src/CudaSirtAlgorithm3D.cpp
+++ b/src/CudaSirtAlgorithm3D.cpp
@@ -172,10 +172,6 @@ void CCudaSirtAlgorithm3D::run(int _iNrIterations)
ASTRA_ASSERT(m_bIsInitialized);
const CProjectionGeometry3D* projgeom = m_pSinogram->getGeometry();
- const CConeProjectionGeometry3D* conegeom = dynamic_cast<const CConeProjectionGeometry3D*>(projgeom);
- const CParallelProjectionGeometry3D* par3dgeom = dynamic_cast<const CParallelProjectionGeometry3D*>(projgeom);
- const CParallelVecProjectionGeometry3D* parvec3dgeom = dynamic_cast<const CParallelVecProjectionGeometry3D*>(projgeom);
- const CConeVecProjectionGeometry3D* conevec3dgeom = dynamic_cast<const CConeVecProjectionGeometry3D*>(projgeom);
const CVolumeGeometry3D& volgeom = *m_pReconstruction->getGeometry();
bool ok = true;
@@ -184,39 +180,7 @@ void CCudaSirtAlgorithm3D::run(int _iNrIterations)
ok &= m_pSirt->setGPUIndex(m_iGPUIndex);
- ok &= m_pSirt->setReconstructionGeometry(volgeom.getGridColCount(),
- volgeom.getGridRowCount(),
- volgeom.getGridSliceCount());
-
- if (conegeom) {
- ok &= m_pSirt->setConeGeometry(conegeom->getProjectionCount(),
- conegeom->getDetectorColCount(),
- conegeom->getDetectorRowCount(),
- conegeom->getOriginSourceDistance(),
- conegeom->getOriginDetectorDistance(),
- conegeom->getDetectorSpacingX(),
- conegeom->getDetectorSpacingY(),
- conegeom->getProjectionAngles());
- } else if (par3dgeom) {
- ok &= m_pSirt->setPar3DGeometry(par3dgeom->getProjectionCount(),
- par3dgeom->getDetectorColCount(),
- par3dgeom->getDetectorRowCount(),
- par3dgeom->getDetectorSpacingX(),
- par3dgeom->getDetectorSpacingY(),
- par3dgeom->getProjectionAngles());
- } else if (parvec3dgeom) {
- ok &= m_pSirt->setPar3DGeometry(parvec3dgeom->getProjectionCount(),
- parvec3dgeom->getDetectorColCount(),
- parvec3dgeom->getDetectorRowCount(),
- parvec3dgeom->getProjectionVectors());
- } else if (conevec3dgeom) {
- ok &= m_pSirt->setConeGeometry(conevec3dgeom->getProjectionCount(),
- conevec3dgeom->getDetectorColCount(),
- conevec3dgeom->getDetectorRowCount(),
- conevec3dgeom->getProjectionVectors());
- } else {
- ASTRA_ASSERT(false);
- }
+ ok &= m_pSirt->setGeometry(&volgeom, projgeom);
ok &= m_pSirt->enableSuperSampling(m_iVoxelSuperSampling, m_iDetectorSuperSampling);