Add flexible volgeom3d support to astraCudaFP

3 files changed, 39 insertions, 273 deletions

diff --git a/cuda/3d/astra3d.cu b/cuda/3d/astra3d.cu
index 97bebf4..b2375f3 100644
--- a/cuda/3d/astra3d.cu
+++ b/cuda/3d/astra3d.cu
@@ -1103,179 +1103,31 @@ float AstraCGLS3d::computeDiffNorm()
 
 
 
-bool astraCudaConeFP(const float* pfVolume, float* pfProjections,
-                     unsigned int iVolX,
-                     unsigned int iVolY,
-                     unsigned int iVolZ,
-                     unsigned int iProjAngles,
-                     unsigned int iProjU,
-                     unsigned int iProjV,
-                     float fOriginSourceDistance,
-                     float fOriginDetectorDistance,
-                     float fDetUSize,
-                     float fDetVSize,
-                     const float *pfAngles,
-                     int iGPUIndex, int iDetectorSuperSampling)
-{
-	if (iVolX == 0 || iVolY == 0 || iVolZ == 0)
-		return false;
-	if (iProjAngles == 0 || iProjU == 0 || iProjV == 0 || pfAngles == 0)
-		return false;
-
-	SConeProjection* p = genConeProjections(iProjAngles,
-                                            iProjU, iProjV,
-                                            fOriginSourceDistance,
-                                            fOriginDetectorDistance,
-                                            fDetUSize, fDetVSize,
-                                            pfAngles);
-
-	bool ok;
-	ok = astraCudaConeFP(pfVolume, pfProjections, iVolX, iVolY, iVolZ,
-	                     iProjAngles, iProjU, iProjV, p, iGPUIndex, iDetectorSuperSampling);
-
-	delete[] p;
-
-	return ok;
-}
-
-bool astraCudaConeFP(const float* pfVolume, float* pfProjections,
-                     unsigned int iVolX,
-                     unsigned int iVolY,
-                     unsigned int iVolZ,
-                     unsigned int iProjAngles,
-                     unsigned int iProjU,
-                     unsigned int iProjV,
-                     const SConeProjection *pfAngles,
-                     int iGPUIndex, int iDetectorSuperSampling)
+bool astraCudaFP(const float* pfVolume, float* pfProjections,
+                 const CVolumeGeometry3D* pVolGeom,
+                 const CProjectionGeometry3D* pProjGeom,
+                 int iGPUIndex, int iDetectorSuperSampling,
+                 Cuda3DProjectionKernel projKernel)
 {
 	SDimensions3D dims;
 
-	dims.iVolX = iVolX;
-	dims.iVolY = iVolY;
-	dims.iVolZ = iVolZ;
-	if (iVolX == 0 || iVolY == 0 || iVolZ == 0)
-		return false;
-
-	dims.iProjAngles = iProjAngles;
-	dims.iProjU = iProjU;
-	dims.iProjV = iProjV;
-
-	if (iProjAngles == 0 || iProjU == 0 || iProjV == 0 || pfAngles == 0)
+	bool ok = convertAstraGeometry_dims(pVolGeom, pProjGeom, dims);
+	if (!ok)
 		return false;
 
 	dims.iRaysPerDetDim = iDetectorSuperSampling;
-
 	if (iDetectorSuperSampling == 0)
 		return false;
 
-	if (iGPUIndex != -1) {
-		cudaSetDevice(iGPUIndex);
-		cudaError_t err = cudaGetLastError();
-
-		// Ignore errors caused by calling cudaSetDevice multiple times
-		if (err != cudaSuccess && err != cudaErrorSetOnActiveProcess)
-			return false;
-	}
-
-	cudaPitchedPtr D_volumeData = allocateVolumeData(dims);
-	bool ok = D_volumeData.ptr;
-	if (!ok)
-		return false;
-
-	cudaPitchedPtr D_projData = allocateProjectionData(dims);
-	ok = D_projData.ptr;
-	if (!ok) {
-		cudaFree(D_volumeData.ptr);
-		return false;
-	}
-
-	ok &= copyVolumeToDevice(pfVolume, D_volumeData, dims, dims.iVolX);
-
-	ok &= zeroProjectionData(D_projData, dims);
-
-	if (!ok) {
-		cudaFree(D_volumeData.ptr);
-		cudaFree(D_projData.ptr);
-		return false;
-	}
-
-	ok &= ConeFP(D_volumeData, D_projData, dims, pfAngles, 1.0f);
-
-	ok &= copyProjectionsFromDevice(pfProjections, D_projData,
-	                                dims, dims.iProjU);
-
-
-	cudaFree(D_volumeData.ptr);
-	cudaFree(D_projData.ptr);
-
-	return ok;
-
-}
-
-bool astraCudaPar3DFP(const float* pfVolume, float* pfProjections,
-                      unsigned int iVolX,
-                      unsigned int iVolY,
-                      unsigned int iVolZ,
-                      unsigned int iProjAngles,
-                      unsigned int iProjU,
-                      unsigned int iProjV,
-                      float fDetUSize,
-                      float fDetVSize,
-                      const float *pfAngles,
-                      int iGPUIndex, int iDetectorSuperSampling,
-                      Cuda3DProjectionKernel projKernel)
-{
-	if (iVolX == 0 || iVolY == 0 || iVolZ == 0)
-		return false;
-	if (iProjAngles == 0 || iProjU == 0 || iProjV == 0 || pfAngles == 0)
-		return false;
-
-	SPar3DProjection* p = genPar3DProjections(iProjAngles,
-                                             iProjU, iProjV,
-                                             fDetUSize, fDetVSize,
-                                             pfAngles);
-
-	bool ok;
-	ok = astraCudaPar3DFP(pfVolume, pfProjections, iVolX, iVolY, iVolZ,
-	                      iProjAngles, iProjU, iProjV, p, iGPUIndex, iDetectorSuperSampling,
-	                      projKernel);
-
-	delete[] p;
-
-	return ok;
-}
-
-
-bool astraCudaPar3DFP(const float* pfVolume, float* pfProjections,
-                      unsigned int iVolX,
-                      unsigned int iVolY,
-                      unsigned int iVolZ,
-                      unsigned int iProjAngles,
-                      unsigned int iProjU,
-                      unsigned int iProjV,
-                      const SPar3DProjection *pfAngles,
-                      int iGPUIndex, int iDetectorSuperSampling,
-                      Cuda3DProjectionKernel projKernel)
-{
-	SDimensions3D dims;
-
-	dims.iVolX = iVolX;
-	dims.iVolY = iVolY;
-	dims.iVolZ = iVolZ;
-	if (iVolX == 0 || iVolY == 0 || iVolZ == 0)
-		return false;
-
-	dims.iProjAngles = iProjAngles;
-	dims.iProjU = iProjU;
-	dims.iProjV = iProjV;
+	SPar3DProjection* pParProjs;
+	SConeProjection* pConeProjs;
 
-	if (iProjAngles == 0 || iProjU == 0 || iProjV == 0 || pfAngles == 0)
-		return false;
+	float outputScale;
 
-	dims.iRaysPerDetDim = iDetectorSuperSampling;
+	ok = convertAstraGeometry(pVolGeom, pProjGeom,
+	                          pParProjs, pConeProjs,
+	                          outputScale);
 
-	if (iDetectorSuperSampling == 0)
-		return false;
 
 	if (iGPUIndex != -1) {
 		cudaSetDevice(iGPUIndex);
@@ -1288,7 +1140,7 @@ bool astraCudaPar3DFP(const float* pfVolume, float* pfProjections,
 
 
 	cudaPitchedPtr D_volumeData = allocateVolumeData(dims);
-	bool ok = D_volumeData.ptr;
+	ok = D_volumeData.ptr;
 	if (!ok)
 		return false;
 
@@ -1309,15 +1161,25 @@ bool astraCudaPar3DFP(const float* pfVolume, float* pfProjections,
 		return false;
 	}
 
-	switch (projKernel) {
-	case ker3d_default:
-		ok &= Par3DFP(D_volumeData, D_projData, dims, pfAngles, 1.0f);
-		break;
-	case ker3d_sum_square_weights:
-		ok &= Par3DFP_SumSqW(D_volumeData, D_projData, dims, pfAngles, 1.0f);
-		break;
-	default:
-		assert(false);
+	if (pParProjs) {
+		switch (projKernel) {
+		case ker3d_default:
+			ok &= Par3DFP(D_volumeData, D_projData, dims, pParProjs, 1.0f);
+			break;
+		case ker3d_sum_square_weights:
+			ok &= Par3DFP_SumSqW(D_volumeData, D_projData, dims, pParProjs, 1.0f);
+			break;
+		default:
+			assert(false);
+		}
+	} else {
+		switch (projKernel) {
+		case ker3d_default:
+			ok &= ConeFP(D_volumeData, D_projData, dims, pConeProjs, 1.0f);
+			break;
+		default:
+			assert(false);
+		}
 	}
 
 	ok &= copyProjectionsFromDevice(pfProjections, D_projData,
diff --git a/cuda/3d/astra3d.h b/cuda/3d/astra3d.h
index 5464d2f..6c3fcfb 100644
--- a/cuda/3d/astra3d.h
+++ b/cuda/3d/astra3d.h
@@ -282,52 +282,9 @@ protected:
 };
 
 
-
-_AstraExport bool astraCudaConeFP(const float* pfVolume, float* pfProjections,
-                     unsigned int iVolX,
-                     unsigned int iVolY,
-                     unsigned int iVolZ,
-                     unsigned int iProjAngles,
-                     unsigned int iProjU,
-                     unsigned int iProjV,
-                     float fOriginSourceDistance,
-                     float fOriginDetectorDistance,
-                     float fDetUSize,
-                     float fDetVSize,
-                     const float *pfAngles,
-                     int iGPUIndex, int iDetectorSuperSampling);
-
-_AstraExport bool astraCudaConeFP(const float* pfVolume, float* pfProjections,
-                     unsigned int iVolX,
-                     unsigned int iVolY,
-                     unsigned int iVolZ,
-                     unsigned int iProjAngles,
-                     unsigned int iProjU,
-                     unsigned int iProjV,
-                     const SConeProjection *pfAngles,
-                     int iGPUIndex, int iDetectorSuperSampling);
-
-_AstraExport bool astraCudaPar3DFP(const float* pfVolume, float* pfProjections,
-                      unsigned int iVolX,
-                      unsigned int iVolY,
-                      unsigned int iVolZ,
-                      unsigned int iProjAngles,
-                      unsigned int iProjU,
-                      unsigned int iProjV,
-                      float fDetUSize,
-                      float fDetVSize,
-                      const float *pfAngles,
-                      int iGPUIndex, int iDetectorSuperSampling,
-                      Cuda3DProjectionKernel projKernel);
-
-_AstraExport bool astraCudaPar3DFP(const float* pfVolume, float* pfProjections,
-                      unsigned int iVolX,
-                      unsigned int iVolY,
-                      unsigned int iVolZ,
-                      unsigned int iProjAngles,
-                      unsigned int iProjU,
-                      unsigned int iProjV,
-                      const SPar3DProjection *pfAngles,
+_AstraExport bool astraCudaFP(const float* pfVolume, float* pfProjections,
+                      const CVolumeGeometry3D* pVolGeom,
+                      const CProjectionGeometry3D* pProjGeom,
                       int iGPUIndex, int iDetectorSuperSampling,
                       Cuda3DProjectionKernel projKernel);
 
diff --git a/src/CudaForwardProjectionAlgorithm3D.cpp b/src/CudaForwardProjectionAlgorithm3D.cpp
index bb122e0..914ee2f 100644
--- a/src/CudaForwardProjectionAlgorithm3D.cpp
+++ b/src/CudaForwardProjectionAlgorithm3D.cpp
@@ -239,10 +239,6 @@ void CCudaForwardProjectionAlgorithm3D::run(int)
 	assert(m_bIsInitialized);
 
 	const CProjectionGeometry3D* projgeom = m_pProjections->getGeometry();
-	const CConeProjectionGeometry3D* conegeom = dynamic_cast<const CConeProjectionGeometry3D*>(projgeom);
-	const CParallelProjectionGeometry3D* par3dgeom = dynamic_cast<const CParallelProjectionGeometry3D*>(projgeom);
-	const CConeVecProjectionGeometry3D* conevecgeom = dynamic_cast<const CConeVecProjectionGeometry3D*>(projgeom);
-	const CParallelVecProjectionGeometry3D* parvec3dgeom = dynamic_cast<const CParallelVecProjectionGeometry3D*>(projgeom);
 	const CVolumeGeometry3D& volgeom = *m_pVolume->getGeometry();
 
 	Cuda3DProjectionKernel projKernel = ker3d_default;
@@ -270,58 +266,9 @@ void CCudaForwardProjectionAlgorithm3D::run(int)
 	}
 #endif
 
-	if (conegeom) {
-		astraCudaConeFP(m_pVolume->getDataConst(), m_pProjections->getData(),
-		                volgeom.getGridColCount(),
-		                volgeom.getGridRowCount(),
-		                volgeom.getGridSliceCount(),
-		                conegeom->getProjectionCount(),
-		                conegeom->getDetectorColCount(),
-		                conegeom->getDetectorRowCount(),
-		                conegeom->getOriginSourceDistance(),
-		                conegeom->getOriginDetectorDistance(),
-		                conegeom->getDetectorSpacingX(),
-		                conegeom->getDetectorSpacingY(),
-		                conegeom->getProjectionAngles(),
-		                m_iGPUIndex, m_iDetectorSuperSampling);
-	} else if (par3dgeom) {
-		astraCudaPar3DFP(m_pVolume->getDataConst(), m_pProjections->getData(),
-		                 volgeom.getGridColCount(),
-		                 volgeom.getGridRowCount(),
-		                 volgeom.getGridSliceCount(),
-		                 par3dgeom->getProjectionCount(),
-		                 par3dgeom->getDetectorColCount(),
-		                 par3dgeom->getDetectorRowCount(),
-		                 par3dgeom->getDetectorSpacingX(),
-		                 par3dgeom->getDetectorSpacingY(),
-		                 par3dgeom->getProjectionAngles(),
-		                 m_iGPUIndex, m_iDetectorSuperSampling,
-		                 projKernel);
-	} else if (parvec3dgeom) {
-		astraCudaPar3DFP(m_pVolume->getDataConst(), m_pProjections->getData(),
-		                 volgeom.getGridColCount(),
-		                 volgeom.getGridRowCount(),
-		                 volgeom.getGridSliceCount(),
-		                 parvec3dgeom->getProjectionCount(),
-		                 parvec3dgeom->getDetectorColCount(),
-		                 parvec3dgeom->getDetectorRowCount(),
-		                 parvec3dgeom->getProjectionVectors(),
-		                 m_iGPUIndex, m_iDetectorSuperSampling,
-		                 projKernel);
-	} else if (conevecgeom) {
-		astraCudaConeFP(m_pVolume->getDataConst(), m_pProjections->getData(),
-		                volgeom.getGridColCount(),
-		                volgeom.getGridRowCount(),
-		                volgeom.getGridSliceCount(),
-		                conevecgeom->getProjectionCount(),
-		                conevecgeom->getDetectorColCount(),
-		                conevecgeom->getDetectorRowCount(),
-		                conevecgeom->getProjectionVectors(),
-		                m_iGPUIndex, m_iDetectorSuperSampling);
-	} else {
-		ASTRA_ASSERT(false);
-	}
-
+	astraCudaFP(m_pVolume->getDataConst(), m_pProjections->getData(),
+	            &volgeom, projgeom,
+	            m_iGPUIndex, m_iDetectorSuperSampling, projKernel);
 }