Commit ace5a174 authored by Pierre Paleo's avatar Pierre Paleo
Browse files

FBP OK with different filters

parent 19e4178a
......@@ -89,7 +89,7 @@ def get_shape_dtype(arr):
elif isinstance(arr, cuda.Array):
return cuarray_shape_dtype(arr)
else:
raise ValueError("Unknown array type")
raise ValueError("Unknown array type %s" % str(type(arr)))
def copy_array(dst, src, check=False, src_dtype=None):
......
......@@ -56,16 +56,16 @@ class Backprojector(CudaProcessing):
- "axis_correction": None
"""
super().__init__(**cuda_options)
self.configure_extra_options(extra_options=extra_options)
self.init_geometry(sino_shape, slice_shape, angles, rot_center)
self._configure_extra_options(extra_options=extra_options)
self._init_geometry(sino_shape, slice_shape, angles, rot_center)
self._init_filter(filter_name)
self.allocate_memory()
self.compute_angles()
self.compile_kernels()
self.bind_texture()
self._allocate_memory()
self._compute_angles()
self._compile_kernels()
self._bind_textures()
def configure_extra_options(self, extra_options={}):
def _configure_extra_options(self, extra_options={}):
self._axis_array = None
self.extra_options = {
"padding_mode": "zeros",
......@@ -75,7 +75,7 @@ class Backprojector(CudaProcessing):
self._axis_array = self.extra_options["axis_correction"]
def init_geometry(self, slice_shape, angles, rot_center):
def _init_geometry(self, sino_shape, slice_shape, angles, rot_center):
self.sino_shape = sino_shape
if len(sino_shape) == 2:
n_angles, dwidth = sino_shape
......@@ -98,6 +98,8 @@ class Backprojector(CudaProcessing):
self.n_sinos = n_sinos
self.n_angles = n_angles
self.dwidth = dwidth
if angles is None:
angles = n_angles
if np.isscalar(angles):
angles = np.linspace(0, np.pi, angles, False)
else:
......@@ -108,15 +110,16 @@ class Backprojector(CudaProcessing):
self.axis_pos = self.rot_center
def allocate_memory(self):
self._d_slice = garray.zeros(self.slice_shape, "f")
self._d_sino = garray.zeros(self.sino_shape, "f")
def _allocate_memory(self):
self._d_sino_cua = cuda.np_to_array(np.zeros(self.sino_shape, "f"), "C")
# 1D textures are not supported in pycuda
self.h_msin = np.zeros((1, self.n_angles), "f")
self.h_cos = np.zeros((1, self.n_angles), "f")
self._d_sino = garray.zeros(self.sino_shape, "f")
# ~ self._d_slice = garray.zeros(self.slice_shape, "f")
self._init_arrays_to_none(["_d_slice"])
def compute_angles(self):
def _compute_angles(self):
self.h_cos[0] = np.cos(self.angles).astype("f")
self.h_msin[0] = (-np.sin(self.angles)).astype("f")
self._d_msin = garray.to_gpu(self.h_msin[0])
......@@ -131,7 +134,7 @@ class Backprojector(CudaProcessing):
ctx=self.ctx,
)
def compile_kernels(self):
def _compile_kernels(self):
# Configure backprojector
fname = get_cuda_srcfile("backproj.cu")
self.gpu_projector = CudaKernel("backproj", filename=fname)
......@@ -147,12 +150,13 @@ class Backprojector(CudaProcessing):
)
def bind_texture(self):
def _bind_textures(self):
self.texref_proj.set_array(self._d_sino_cua)
def set_output(self, output, check=False):
def _set_output(self, output, check=False):
if output is None:
self._allocate_array("_d_slice", self.slice_shape, dtype=np.float32)
return self._d_slice
if check:
assert output.dtype == np.float32
......@@ -163,15 +167,14 @@ class Backprojector(CudaProcessing):
return output
# TODO support output=<numpy array>
def backproj(self, sino, output=None, do_checks=True):
copy_array(self._d_sino_cua, sino, check=do_checks)
d_slice = self.set_output(output, check=do_checks)
d_slice = self._set_output(output, check=do_checks)
self.gpu_projector(
d_slice,
self.n_angles,
self.dwidth_x,
self.dwidth,
self.axis_pos,
self.n_x,
self.n_y,
......
......@@ -172,7 +172,6 @@ class SinoFilter(CudaProcessing):
grid=self._pad_grid,
block=self._pad_block
)
else: # zeros
self.d_sino_padded.fill(0)
if self.ndim == 2:
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment