summaryrefslogtreecommitdiffstats
path: root/cuda
diff options
context:
space:
mode:
authorSuren A. Chilingaryan <csa@ipecompute4.ands.kit.edu>2022-07-26 23:29:12 +0200
committerSuren A. Chilingaryan <csa@ipecompute4.ands.kit.edu>2022-07-26 23:29:12 +0200
commit9888385e475a747250ab802ee137b76f2b87925d (patch)
tree9fd968221e25d71dec56cb285e0d5942b5751c38 /cuda
parent0bad868e66827f94ae93a909dc8ea330a8540a11 (diff)
downloadastra-9888385e475a747250ab802ee137b76f2b87925d.tar.gz
astra-9888385e475a747250ab802ee137b76f2b87925d.tar.bz2
astra-9888385e475a747250ab802ee137b76f2b87925d.tar.xz
astra-9888385e475a747250ab802ee137b76f2b87925d.zip
Half-precision back-/forward-projection for parallel geometryhalf
Diffstat (limited to 'cuda')
-rw-r--r--cuda/3d/par3d_bp.cu20
-rw-r--r--cuda/3d/par3d_fp.cu22
-rw-r--r--cuda/3d/rounding.h9
3 files changed, 49 insertions, 2 deletions
diff --git a/cuda/3d/par3d_bp.cu b/cuda/3d/par3d_bp.cu
index 27d95fe..7958ac9 100644
--- a/cuda/3d/par3d_bp.cu
+++ b/cuda/3d/par3d_bp.cu
@@ -54,6 +54,8 @@ __constant__ DevPar3DParams gC_C[g_MaxAngles];
__constant__ float gC_scale[g_MaxAngles];
+#include "rounding.h"
+
template<unsigned int ZSIZE>
__global__ void dev_par3D_BP(void* D_volData, unsigned int volPitch, cudaTextureObject_t tex, int startAngle, int angleOffset, const SDimensions3D dims, float fOutputScale)
{
@@ -100,10 +102,26 @@ __global__ void dev_par3D_BP(void* D_volData, unsigned int volPitch, cudaTexture
float fU = fCu.w + fX * fCu.x + fY * fCu.y + fZ * fCu.z;
float fV = fCv.w + fX * fCv.x + fY * fCv.y + fZ * fCv.z;
+// printf("%f %f\n", fU, fV);
for (int idx = 0; idx < ZSIZE; ++idx) {
+ float fVal;
+ textype h5 = texto(0.5f);
+ textype fU_ = texto(fU);
+ textype fUf_ = texto(floor(fU));
+ float fUf = floor(fU);
+
+ if ((fU - fUf) < 0.5f) {
+ textype fVal1 = texto(tex3D<float>(tex, fUf - 0.5f, fAngle, fV));
+ textype fVal2 = texto(tex3D<float>(tex, fUf + 0.5f, fAngle, fV));
+ fVal = texfrom(fVal1 + (fU_ + h5 - fUf_) * (fVal2 - fVal1));
+ } else {
+ textype fVal1 = texto(tex3D<float>(tex, fUf + 0.5f, fAngle, fV));
+ textype fVal2 = texto(tex3D<float>(tex, fUf + 1.5f, fAngle, fV));
+ fVal = texfrom(fVal1 + (fU_ - h5 - fUf_) * (fVal2 - fVal1));
+ }
- float fVal = tex3D<float>(tex, fU, fAngle, fV);
+// float fVal = tex3D<float>(tex, fU, fAngle, fV);
Z[idx] += fVal * fS;
fU += fCu.z;
diff --git a/cuda/3d/par3d_fp.cu b/cuda/3d/par3d_fp.cu
index fda6f93..075784b 100644
--- a/cuda/3d/par3d_fp.cu
+++ b/cuda/3d/par3d_fp.cu
@@ -146,6 +146,7 @@ bool transferConstants(const SPar3DProjection* angles, unsigned int iProjAngles)
// blockIdx: x = u/v detector
// y = angle block
+#include "rounding.h"
template<class COORD, class SCALE>
__global__ void par3D_FP_t(float* D_projData, unsigned int projPitch,
@@ -212,10 +213,28 @@ __global__ void par3D_FP_t(float* D_projData, unsigned int projPitch,
float f0 = startSlice + 0.5f;
float f1 = a1 * (startSlice - 0.5f*c.nSlices(dims) + 0.5f) + b1 + 0.5f*c.nDim1(dims) - 0.5f + 0.5f;
float f2 = a2 * (startSlice - 0.5f*c.nSlices(dims) + 0.5f) + b2 + 0.5f*c.nDim2(dims) - 0.5f + 0.5f;
+ //printf("%f, %f (%f), %f (%f)\n", f0, f1, a1, f2, a2); // Only f1 non linear
for (int s = startSlice; s < endSlice; ++s)
{
- fVal += c.tex(tex, f0, f1, f2);
+
+ textype h5 = texto(0.5f);
+ textype f1_ = texto(f1);
+ textype f1f_ = texto(floor(f1));
+ float f1f = floor(f1);
+
+ if ((f1 - f1f) < 0.5f) {
+ textype fVal1 = texto(c.tex(tex, f0, f1f - 0.5f, f2));
+ textype fVal2 = texto(c.tex(tex, f0, f1f + 0.5f, f2));
+ fVal += texfrom(fVal1 + (f1_ + h5 - f1f_) * (fVal2 - fVal1));
+// fVal += texfrom(__hfma(__hadd(h5,__hsub(f1_, f1f_)), __hsub(fVal2, fVal1), fVal1));
+ } else {
+ textype fVal1 = texto(c.tex(tex, f0, f1f + 0.5f, f2));
+ textype fVal2 = texto(c.tex(tex, f0, f1f + 1.5f, f2));
+ fVal += texfrom(fVal1 + (f1_ - h5 - f1f_) * (fVal2 - fVal1));
+ }
+
+// fVal += c.tex(tex, f0, f1, f2);
f0 += 1.0f;
f1 += a1;
f2 += a2;
@@ -308,6 +327,7 @@ __global__ void par3D_FP_SS_t(float* D_projData, unsigned int projPitch,
for (int s = startSlice; s < endSlice; ++s)
{
fVal += c.tex(tex, f0, f1, f2);
+
f0 += 1.0f;
f1 += a1;
f2 += a2;
diff --git a/cuda/3d/rounding.h b/cuda/3d/rounding.h
new file mode 100644
index 0000000..a263f84
--- /dev/null
+++ b/cuda/3d/rounding.h
@@ -0,0 +1,9 @@
+#include <cuda_fp16.h>
+
+#define texto __float2half
+#define texfrom __half2float
+#define textype half
+
+//#define texto
+//#define texfrom
+//#define textype float