summaryrefslogtreecommitdiff
path: root/silx/opencl/test
diff options
context:
space:
mode:
Diffstat (limited to 'silx/opencl/test')
-rw-r--r--silx/opencl/test/__init__.py4
-rw-r--r--silx/opencl/test/test_backprojection.py21
-rw-r--r--silx/opencl/test/test_convolution.py263
-rw-r--r--silx/opencl/test/test_sparse.py192
4 files changed, 479 insertions, 1 deletions
diff --git a/silx/opencl/test/__init__.py b/silx/opencl/test/__init__.py
index f3121d5..2e90e66 100644
--- a/silx/opencl/test/__init__.py
+++ b/silx/opencl/test/__init__.py
@@ -38,6 +38,8 @@ from ..codec import test as test_codec
from . import test_image
from . import test_kahan
from . import test_stats
+from . import test_convolution
+from . import test_sparse
def suite():
@@ -52,6 +54,8 @@ def suite():
test_suite.addTests(test_image.suite())
test_suite.addTests(test_kahan.suite())
test_suite.addTests(test_stats.suite())
+ test_suite.addTests(test_convolution.suite())
+ test_suite.addTests(test_sparse.suite())
# Allow to remove sift from the project
test_base_dir = os.path.dirname(__file__)
sift_dir = os.path.join(test_base_dir, "..", "sift")
diff --git a/silx/opencl/test/test_backprojection.py b/silx/opencl/test/test_backprojection.py
index 7ab416d..b2f2070 100644
--- a/silx/opencl/test/test_backprojection.py
+++ b/silx/opencl/test/test_backprojection.py
@@ -45,7 +45,7 @@ except ImportError:
from ..common import ocl
if ocl:
from .. import backprojection
- from ..sinofilter import compute_fourier_filter
+ from ...image.tomography import compute_fourier_filter
from silx.test.utils import utilstest
logger = logging.getLogger(__name__)
@@ -199,11 +199,30 @@ class TestFBP(unittest.TestCase):
"Something wrong with FBP(filter=%s)" % filter_name
)
+ @unittest.skipUnless(ocl and mako, "pyopencl is missing")
+ def test_fbp_oddsize(self):
+ # Generate a 513-sinogram.
+ # The padded width will be nextpow(513*2).
+ # silx [0.10, 0.10.1] will give 1029, which makes R2C transform fail.
+ sino = np.pad(self.sino, ((0, 0), (1, 0)), mode='edge')
+ B = backprojection.Backprojection(sino.shape, axis_position=self.fbp.axis_pos+1)
+ res = B(sino)
+ # Compare with self.reference_rec. Tolerance is high as backprojector
+ # is not fully shift-invariant.
+ errmax = np.max(np.abs(clip_circle(res[1:, 1:] - self.reference_rec)))
+ self.assertLess(
+ errmax, 1.e-1,
+ "Something wrong with FBP on odd-sized sinogram"
+ )
+
+
+
def suite():
testSuite = unittest.TestSuite()
testSuite.addTest(TestFBP("test_fbp"))
testSuite.addTest(TestFBP("test_fbp_filters"))
+ testSuite.addTest(TestFBP("test_fbp_oddsize"))
return testSuite
diff --git a/silx/opencl/test/test_convolution.py b/silx/opencl/test/test_convolution.py
new file mode 100644
index 0000000..8acd385
--- /dev/null
+++ b/silx/opencl/test/test_convolution.py
@@ -0,0 +1,263 @@
+#!/usr/bin/env python
+# coding: utf-8
+# /*##########################################################################
+#
+# Copyright (c) 2019 European Synchrotron Radiation Facility
+#
+# Permission is hereby granted, free of charge, to any person obtaining a copy
+# of this software and associated documentation files (the "Software"), to deal
+# in the Software without restriction, including without limitation the rights
+# to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
+# copies of the Software, and to permit persons to whom the Software is
+# furnished to do so, subject to the following conditions:
+#
+# The above copyright notice and this permission notice shall be included in
+# all copies or substantial portions of the Software.
+#
+# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
+# IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+# FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
+# AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+# LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
+# OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
+# THE SOFTWARE.
+#
+# ###########################################################################*/
+
+"""
+Test of the Convolution class.
+"""
+
+from __future__ import division, print_function
+
+__authors__ = ["Pierre Paleo"]
+__contact__ = "pierre.paleo@esrf.fr"
+__license__ = "MIT"
+__copyright__ = "2019 European Synchrotron Radiation Facility, Grenoble, France"
+__date__ = "15/02/2019"
+
+import logging
+from itertools import product
+import numpy as np
+from silx.utils.testutils import parameterize
+try:
+ from scipy.ndimage import convolve, convolve1d
+ from scipy.misc import ascent
+ scipy_convolve = convolve
+ scipy_convolve1d = convolve1d
+except ImportError:
+ scipy_convolve = None
+import unittest
+from ..common import ocl
+if ocl:
+ import pyopencl as cl
+ import pyopencl.array as parray
+ from ..convolution import Convolution, gaussian_kernel
+logger = logging.getLogger(__name__)
+
+
+@unittest.skipUnless(ocl and scipy_convolve, "PyOpenCl/scipy is missing")
+class TestConvolution(unittest.TestCase):
+
+ @classmethod
+ def setUpClass(cls):
+ super(TestConvolution, cls).setUpClass()
+ cls.image = np.ascontiguousarray(ascent()[:, :511], dtype="f")
+ cls.data1d = cls.image[0]
+ cls.data2d = cls.image
+ cls.data3d = np.tile(cls.image[224:-224, 224:-224], (62, 1, 1))
+ cls.kernel1d = gaussian_kernel(1.)
+ cls.kernel2d = np.outer(cls.kernel1d, cls.kernel1d)
+ cls.kernel3d = np.multiply.outer(cls.kernel2d, cls.kernel1d)
+ cls.ctx = ocl.create_context()
+ cls.tol = {
+ "1D": 1e-4,
+ "2D": 1e-3,
+ "3D": 1e-3,
+ }
+
+ @classmethod
+ def tearDownClass(cls):
+ cls.data1d = cls.data2d = cls.data3d = cls.image = None
+ cls.kernel1d = cls.kernel2d = cls.kernel3d = None
+
+ @staticmethod
+ def compare(arr1, arr2):
+ return np.max(np.abs(arr1 - arr2))
+
+ @staticmethod
+ def print_err(conv):
+ errmsg = str(
+ """
+ Something wrong with %s
+ mode=%s, texture=%s
+ """
+ % (conv.use_case_desc, conv.mode, conv.use_textures)
+ )
+ return errmsg
+
+ def __init__(self, methodName='runTest', param=None):
+ unittest.TestCase.__init__(self, methodName)
+ self.param = param
+ self.mode = param["boundary_handling"]
+ logger.debug(
+ """
+ Testing convolution with boundary_handling=%s,
+ use_textures=%s, input_device=%s, output_device=%s
+ """
+ % (
+ self.mode, param["use_textures"],
+ param["input_on_device"], param["output_on_device"]
+ )
+ )
+
+ def instantiate_convol(self, shape, kernel, axes=None):
+ if (self.mode == "constant") and (
+ not(self.param["use_textures"])
+ or (self.ctx.devices[0].type == cl._cl.device_type.CPU)
+ ):
+ self.skipTest("mode=constant not implemented without textures")
+ C = Convolution(
+ shape, kernel,
+ mode=self.mode,
+ ctx=self.ctx,
+ axes=axes,
+ extra_options={"dont_use_textures": not(self.param["use_textures"])}
+ )
+ return C
+
+ def get_data_and_kernel(self, test_name):
+ dims = {
+ "test_1D": (1, 1),
+ "test_separable_2D": (2, 1),
+ "test_separable_3D": (3, 1),
+ "test_nonseparable_2D": (2, 2),
+ "test_nonseparable_3D": (3, 3),
+ }
+ dim_data = {
+ 1: self.data1d,
+ 2: self.data2d,
+ 3: self.data3d
+ }
+ dim_kernel = {
+ 1: self.kernel1d,
+ 2: self.kernel2d,
+ 3: self.kernel3d,
+ }
+ dd, kd = dims[test_name]
+ return dim_data[dd], dim_kernel[kd]
+
+ def get_reference_function(self, test_name):
+ ref_func = {
+ "test_1D":
+ lambda x, y : scipy_convolve1d(x, y, mode=self.mode),
+ "test_separable_2D":
+ lambda x, y : scipy_convolve1d(
+ scipy_convolve1d(x, y, mode=self.mode, axis=1),
+ y, mode=self.mode, axis=0
+ ),
+ "test_separable_3D":
+ lambda x, y: scipy_convolve1d(
+ scipy_convolve1d(
+ scipy_convolve1d(x, y, mode=self.mode, axis=2),
+ y, mode=self.mode, axis=1),
+ y, mode=self.mode, axis=0
+ ),
+ "test_nonseparable_2D":
+ lambda x, y: scipy_convolve(x, y, mode=self.mode),
+ "test_nonseparable_3D":
+ lambda x, y : scipy_convolve(x, y, mode=self.mode),
+ }
+ return ref_func[test_name]
+
+ def template_test(self, test_name):
+ data, kernel = self.get_data_and_kernel(test_name)
+ conv = self.instantiate_convol(data.shape, kernel)
+ if self.param["input_on_device"]:
+ data_ref = parray.to_device(conv.queue, data)
+ else:
+ data_ref = data
+ if self.param["output_on_device"]:
+ d_res = parray.zeros_like(conv.data_out)
+ res = d_res
+ else:
+ res = None
+ res = conv(data_ref, output=res)
+ if self.param["output_on_device"]:
+ res = res.get()
+ ref_func = self.get_reference_function(test_name)
+ ref = ref_func(data, kernel)
+ metric = self.compare(res, ref)
+ logger.info("%s: max error = %.2e" % (test_name, metric))
+ tol = self.tol[str("%dD" % kernel.ndim)]
+ self.assertLess(metric, tol, self.print_err(conv))
+
+ def test_1D(self):
+ self.template_test("test_1D")
+
+ def test_separable_2D(self):
+ self.template_test("test_separable_2D")
+
+ def test_separable_3D(self):
+ self.template_test("test_separable_3D")
+
+ def test_nonseparable_2D(self):
+ self.template_test("test_nonseparable_2D")
+
+ def test_nonseparable_3D(self):
+ self.template_test("test_nonseparable_3D")
+
+ def test_batched_2D(self):
+ """
+ Test batched (nonseparable) 2D convolution on 3D data.
+ In this test: batch along "z" (axis 0)
+ """
+ data = self.data3d
+ kernel = self.kernel2d
+ conv = self.instantiate_convol(data.shape, kernel, axes=(0,))
+ res = conv(data) # 3D
+ ref = scipy_convolve(data[0], kernel, mode=self.mode) # 2D
+
+ std = np.std(res, axis=0)
+ std_max = np.max(np.abs(std))
+ self.assertLess(std_max, self.tol["2D"], self.print_err(conv))
+ metric = self.compare(res[0], ref)
+ logger.info("test_nonseparable_3D: max error = %.2e" % metric)
+ self.assertLess(metric, self.tol["2D"], self.print_err(conv))
+
+
+def test_convolution():
+ boundary_handling_ = ["reflect", "nearest", "wrap", "constant"]
+ use_textures_ = [True, False]
+ input_on_device_ = [True, False]
+ output_on_device_ = [True, False]
+ testSuite = unittest.TestSuite()
+
+ param_vals = list(product(
+ boundary_handling_,
+ use_textures_,
+ input_on_device_,
+ output_on_device_
+ ))
+ for boundary_handling, use_textures, input_dev, output_dev in param_vals:
+ testcase = parameterize(
+ TestConvolution,
+ param={
+ "boundary_handling": boundary_handling,
+ "input_on_device": input_dev,
+ "output_on_device": output_dev,
+ "use_textures": use_textures,
+ }
+ )
+ testSuite.addTest(testcase)
+ return testSuite
+
+
+
+def suite():
+ testSuite = test_convolution()
+ return testSuite
+
+
+if __name__ == '__main__':
+ unittest.main(defaultTest="suite")
diff --git a/silx/opencl/test/test_sparse.py b/silx/opencl/test/test_sparse.py
new file mode 100644
index 0000000..56f1ba4
--- /dev/null
+++ b/silx/opencl/test/test_sparse.py
@@ -0,0 +1,192 @@
+#!/usr/bin/env python
+# coding: utf-8
+# /*##########################################################################
+#
+# Copyright (c) 2018-2019 European Synchrotron Radiation Facility
+#
+# Permission is hereby granted, free of charge, to any person obtaining a copy
+# of this software and associated documentation files (the "Software"), to deal
+# in the Software without restriction, including without limitation the rights
+# to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
+# copies of the Software, and to permit persons to whom the Software is
+# furnished to do so, subject to the following conditions:
+#
+# The above copyright notice and this permission notice shall be included in
+# all copies or substantial portions of the Software.
+#
+# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
+# IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+# FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
+# AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+# LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
+# OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
+# THE SOFTWARE.
+#
+# ###########################################################################*/
+"""Test of the sparse module"""
+
+import numpy as np
+import unittest
+import logging
+from itertools import product
+from ..common import ocl
+if ocl:
+ import pyopencl.array as parray
+ from silx.opencl.sparse import CSR
+try:
+ import scipy.sparse as sp
+except ImportError:
+ sp = None
+logger = logging.getLogger(__name__)
+
+
+
+def generate_sparse_random_data(
+ shape=(1000,),
+ data_min=0, data_max=100,
+ density=0.1,
+ use_only_integers=True,
+ dtype="f"):
+ """
+ Generate random sparse data where.
+
+ Parameters
+ ------------
+ shape: tuple
+ Output data shape.
+ data_min: int or float
+ Minimum value of data
+ data_max: int or float
+ Maximum value of data
+ density: float
+ Density of non-zero elements in the output data.
+ Low value of density mean low number of non-zero elements.
+ use_only_integers: bool
+ If set to True, the output data items will be primarily integers,
+ possibly casted to float if dtype is a floating-point type.
+ This can be used for ease of debugging.
+ dtype: str or numpy.dtype
+ Output data type
+ """
+ mask = np.random.binomial(1, density, size=shape)
+ if use_only_integers:
+ d = np.random.randint(data_min, high=data_max, size=shape)
+ else:
+ d = data_min + (data_max - data_min) * np.random.rand(*shape)
+ return (d * mask).astype(dtype)
+
+
+
+@unittest.skipUnless(ocl and sp, "PyOpenCl/scipy is missing")
+class TestCSR(unittest.TestCase):
+ """Test CSR format"""
+
+ def setUp(self):
+ self.array = generate_sparse_random_data(shape=(512, 511))
+ # Compute reference sparsification
+ a_s = sp.csr_matrix(self.array)
+ self.ref_data = a_s.data
+ self.ref_indices = a_s.indices
+ self.ref_indptr = a_s.indptr
+ self.ref_nnz = a_s.nnz
+ # Test possible configurations
+ input_on_device = [False, True]
+ output_on_device = [False, True]
+ self._test_configs = list(product(input_on_device, output_on_device))
+
+
+ def test_sparsification(self):
+ for input_on_device, output_on_device in self._test_configs:
+ self._test_sparsification(input_on_device, output_on_device)
+
+
+ def _test_sparsification(self, input_on_device, output_on_device):
+ current_config = "input on device: %s, output on device: %s" % (
+ str(input_on_device), str(output_on_device)
+ )
+ # Sparsify on device
+ csr = CSR(self.array.shape)
+ if input_on_device:
+ # The array has to be flattened
+ arr = parray.to_device(csr.queue, self.array.ravel())
+ else:
+ arr = self.array
+ if output_on_device:
+ d_data = parray.zeros_like(csr.data)
+ d_indices = parray.zeros_like(csr.indices)
+ d_indptr = parray.zeros_like(csr.indptr)
+ output = (d_data, d_indices, d_indptr)
+ else:
+ output = None
+ data, indices, indptr = csr.sparsify(arr, output=output)
+ if output_on_device:
+ data = data.get()
+ indices = indices.get()
+ indptr = indptr.get()
+ # Compare
+ nnz = self.ref_nnz
+ self.assertTrue(
+ np.allclose(data[:nnz], self.ref_data),
+ "something wrong with sparsified data (%s)"
+ % current_config
+ )
+ self.assertTrue(
+ np.allclose(indices[:nnz], self.ref_indices),
+ "something wrong with sparsified indices (%s)"
+ % current_config
+ )
+ self.assertTrue(
+ np.allclose(indptr, self.ref_indptr),
+ "something wrong with sparsified indices pointers (indptr) (%s)"
+ % current_config
+ )
+
+
+ def test_desparsification(self):
+ for input_on_device, output_on_device in self._test_configs:
+ self._test_desparsification(input_on_device, output_on_device)
+
+
+ def _test_desparsification(self, input_on_device, output_on_device):
+ current_config = "input on device: %s, output on device: %s" % (
+ str(input_on_device), str(output_on_device)
+ )
+ # De-sparsify on device
+ csr = CSR(self.array.shape, max_nnz=self.ref_nnz)
+ if input_on_device:
+ data = parray.to_device(csr.queue, self.ref_data)
+ indices = parray.to_device(csr.queue, self.ref_indices)
+ indptr = parray.to_device(csr.queue, self.ref_indptr)
+ else:
+ data = self.ref_data
+ indices = self.ref_indices
+ indptr = self.ref_indptr
+ if output_on_device:
+ d_arr = parray.zeros_like(csr.array)
+ output = d_arr
+ else:
+ output = None
+ arr = csr.densify(data, indices, indptr, output=output)
+ if output_on_device:
+ arr = arr.get()
+ # Compare
+ self.assertTrue(
+ np.allclose(arr.reshape(self.array.shape), self.array),
+ "something wrong with densified data (%s)"
+ % current_config
+ )
+
+
+
+def suite():
+ suite = unittest.TestSuite()
+ suite.addTest(
+ unittest.defaultTestLoader.loadTestsFromTestCase(TestCSR)
+ )
+ return suite
+
+
+if __name__ == '__main__':
+ unittest.main(defaultTest="suite")
+
+
generated by cgit v1.2.3 (git 2.25.1) at 2024-05-23 18:40:31 +0000