Skip to content
Snippets Groups Projects
Commit be937958 authored by Stephan Seitz's avatar Stephan Seitz
Browse files

Cherry-pick changes from interpolation-refactoring

parent 7cc8270b
1 merge request!135CI: Replace minimal-ubuntu job with ubuntu
Pipeline #21296 failed with stage
in 6 minutes and 42 seconds
...@@ -11,6 +11,8 @@ import itertools ...@@ -11,6 +11,8 @@ import itertools
from os.path import dirname, join from os.path import dirname, join
import numpy as np import numpy as np
import pycuda.autoinit # NOQA
import pycuda.gpuarray as gpuarray
import pytest import pytest
import sympy import sympy
...@@ -108,93 +110,79 @@ def test_rotate_interpolation(address_mode): ...@@ -108,93 +110,79 @@ def test_rotate_interpolation(address_mode):
pyconrad.imshow(out, "out " + address_mode) pyconrad.imshow(out, "out " + address_mode)
@pytest.mark.parametrize('address_mode', ['border', 'wrap', 'clamp', 'mirror']) @pytest.mark.parametrize('dtype', (np.int32, np.float32, np.float64))
def test_rotate_interpolation_gpu(address_mode): @pytest.mark.parametrize('address_mode', ('border', 'wrap', 'clamp', 'mirror'))
pytest.importorskip('pycuda') @pytest.mark.parametrize('use_textures', ('use_textures', False))
import pycuda.autoinit # NOQA def test_rotate_interpolation_gpu(dtype, address_mode, use_textures):
import pycuda.gpuarray as gpuarray
rotation_angle = sympy.pi / 5 rotation_angle = sympy.pi / 5
scale = 1 scale = 1
previous_result = None if dtype == np.int32:
for dtype in [np.int32, np.float32, np.float64]: lenna_gpu = gpuarray.to_gpu(
if dtype == np.int32: np.ascontiguousarray(lenna * 255, dtype))
lenna_gpu = gpuarray.to_gpu( else:
np.ascontiguousarray(lenna * 255, dtype)) lenna_gpu = gpuarray.to_gpu(
else: np.ascontiguousarray(lenna, dtype))
lenna_gpu = gpuarray.to_gpu( x_f, y_f = pystencils.fields('x,y: %s [2d]' % type_map[dtype], ghost_layers=0)
np.ascontiguousarray(lenna, dtype))
for use_textures in [True, False]: transformed = scale * \
x_f, y_f = pystencils.fields('x,y: %s [2d]' % type_map[dtype], ghost_layers=0) sympy.rot_axis3(rotation_angle)[:2, :2] * sympy.Matrix((x_, y_)) - sympy.Matrix([2, 2])
assignments = pystencils.AssignmentCollection({
transformed = scale * sympy.rot_axis3(rotation_angle)[:2, :2] * \ y_f.center(): LinearInterpolator(x_f, address_mode=address_mode).at(transformed)
sympy.Matrix((x_, y_)) - sympy.Matrix([2, 2]) })
assignments = pystencils.AssignmentCollection({ print(assignments)
y_f.center(): LinearInterpolator(x_f, address_mode=address_mode).at(transformed) ast = pystencils.create_kernel(assignments, target='gpu', use_textures_for_interpolation=use_textures)
}) print(ast)
print(assignments) print(pystencils.show_code(ast))
ast = pystencils.create_kernel(assignments, target='gpu', use_textures_for_interpolation=use_textures) kernel = ast.compile()
print(ast)
print(pystencils.show_code(ast)) out = gpuarray.zeros_like(lenna_gpu)
kernel = ast.compile() kernel(x=lenna_gpu, y=out)
pyconrad.imshow(out,
out = gpuarray.zeros_like(lenna_gpu) f"out {address_mode} texture:{use_textures} {type_map[dtype]}")
kernel(x=lenna_gpu, y=out) skimage.io.imsave(f"/tmp/out {address_mode} texture:{use_textures} {type_map[dtype]}.tif",
pyconrad.imshow(out, np.ascontiguousarray(out.get(), np.float32))
f"out {address_mode} texture:{use_textures} {type_map[dtype]}")
skimage.io.imsave(f"/tmp/out {address_mode} texture:{use_textures} {type_map[dtype]}.tif",
np.ascontiguousarray(out.get(), np.float32))
if previous_result is not None:
try:
assert np.allclose(previous_result[4:-4, 4:-4], out.get()[4:-4, 4:-4], rtol=100, atol=1e-3)
except AssertionError: # NOQA
print("Max error: %f" % np.max(previous_result - out.get()))
# pyconrad.imshow(previous_result - out.get(), "Difference image")
# raise e
previous_result = out.get()
@pytest.mark.parametrize('address_mode', ['border', 'wrap', 'clamp', 'mirror']) @pytest.mark.parametrize('address_mode', ['border', 'wrap', 'clamp', 'mirror'])
def test_shift_interpolation_gpu(address_mode): @pytest.mark.parametrize('dtype', [np.float64, np.float32, np.int32])
pytest.importorskip('pycuda') @pytest.mark.parametrize('use_textures', ('use_textures', False,))
import pycuda.autoinit # NOQA def test_shift_interpolation_gpu(address_mode, dtype, use_textures):
import pycuda.gpuarray as gpuarray
rotation_angle = 0 # sympy.pi / 5 rotation_angle = 0 # sympy.pi / 5
scale = 1 scale = 1
# shift = - sympy.Matrix([1.5, 1.5]) # shift = - sympy.Matrix([1.5, 1.5])
shift = sympy.Matrix((0.0, 0.0)) shift = sympy.Matrix((0.0, 0.0))
for dtype in [np.float64, np.float32, np.int32]: if dtype == np.int32:
if dtype == np.int32: lenna_gpu = gpuarray.to_gpu(
lenna_gpu = gpuarray.to_gpu( np.ascontiguousarray(lenna * 255, dtype))
np.ascontiguousarray(lenna * 255, dtype)) else:
else: lenna_gpu = gpuarray.to_gpu(
lenna_gpu = gpuarray.to_gpu( np.ascontiguousarray(lenna, dtype))
np.ascontiguousarray(lenna, dtype))
for use_textures in [True, False]: x_f, y_f = pystencils.fields('x,y: %s [2d]' % type_map[dtype], ghost_layers=0)
x_f, y_f = pystencils.fields('x,y: %s [2d]' % type_map[dtype], ghost_layers=0)
if use_textures:
if use_textures: transformed = scale * sympy.rot_axis3(rotation_angle)[:2, :2] * sympy.Matrix((x_, y_)) + shift
transformed = scale * sympy.rot_axis3(rotation_angle)[:2, :2] * sympy.Matrix((x_, y_)) + shift else:
else: transformed = scale * sympy.rot_axis3(rotation_angle)[:2, :2] * sympy.Matrix((x_, y_)) + shift
transformed = scale * sympy.rot_axis3(rotation_angle)[:2, :2] * sympy.Matrix((x_, y_)) + shift assignments = pystencils.AssignmentCollection({
assignments = pystencils.AssignmentCollection({ y_f.center(): LinearInterpolator(x_f, address_mode=address_mode).at(transformed)
y_f.center(): LinearInterpolator(x_f, address_mode=address_mode).at(transformed) })
}) # print(assignments)
# print(assignments) ast = pystencils.create_kernel(assignments, target='gpu', use_textures_for_interpolation=use_textures)
ast = pystencils.create_kernel(assignments, target='gpu', use_textures_for_interpolation=use_textures) # print(ast)
# print(ast) print(pystencils.show_code(ast))
print(pystencils.show_code(ast)) kernel = ast.compile()
kernel = ast.compile()
out = gpuarray.zeros_like(lenna_gpu)
out = gpuarray.zeros_like(lenna_gpu) kernel(x=lenna_gpu, y=out)
kernel(x=lenna_gpu, y=out) pyconrad.imshow(out,
pyconrad.imshow(out, f"out {address_mode} texture:{use_textures} {type_map[dtype]}")
f"out {address_mode} texture:{use_textures} {type_map[dtype]}") skimage.io.imsave(f"/tmp/out {address_mode} texture:{use_textures} {type_map[dtype]}.tif",
skimage.io.imsave(f"/tmp/out {address_mode} texture:{use_textures} {type_map[dtype]}.tif", np.ascontiguousarray(out.get(), np.float32))
np.ascontiguousarray(out.get(), np.float32))
@pytest.mark.parametrize('address_mode', ['border', 'clamp']) @pytest.mark.parametrize('address_mode', ['border', 'clamp'])
...@@ -222,7 +210,7 @@ def test_rotate_interpolation_size_change(address_mode): ...@@ -222,7 +210,7 @@ def test_rotate_interpolation_size_change(address_mode):
@pytest.mark.parametrize('address_mode, target', @pytest.mark.parametrize('address_mode, target',
itertools.product(['border', 'wrap', 'clamp', 'mirror'], ['cpu', 'gpu'])) itertools.product(['border', 'wrap', 'clamp', 'mirror'], ['cpu']))
def test_field_interpolated(address_mode, target): def test_field_interpolated(address_mode, target):
x_f, y_f = pystencils.fields('x,y: float64 [2d]') x_f, y_f = pystencils.fields('x,y: float64 [2d]')
...@@ -230,7 +218,7 @@ def test_field_interpolated(address_mode, target): ...@@ -230,7 +218,7 @@ def test_field_interpolated(address_mode, target):
y_f.center(): x_f.interpolated_access([0.5 * x_ + 2.7, 0.25 * y_ + 7.2], address_mode=address_mode) y_f.center(): x_f.interpolated_access([0.5 * x_ + 2.7, 0.25 * y_ + 7.2], address_mode=address_mode)
}) })
print(assignments) print(assignments)
ast = pystencils.create_kernel(assignments) ast = pystencils.create_kernel(assignments, target=target)
print(ast) print(ast)
print(pystencils.show_code(ast)) print(pystencils.show_code(ast))
kernel = ast.compile() kernel = ast.compile()
...@@ -238,5 +226,11 @@ def test_field_interpolated(address_mode, target): ...@@ -238,5 +226,11 @@ def test_field_interpolated(address_mode, target):
out = np.zeros_like(lenna) out = np.zeros_like(lenna)
kernel(x=lenna, y=out) kernel(x=lenna, y=out)
pyconrad.imshow(out, "out " + address_mode) pyconrad.imshow(out, "out " + address_mode)
kernel(x=lenna, y=out)
pyconrad.imshow(out, "out " + address_mode)
def test_spatial_derivative():
x, y = pystencils.fields('x, y: float32[2d]')
tx, ty = pystencils.fields('t_x, t_y: float32[2d]')
diff = sympy.diff(x.interpolated_access((tx.center, ty.center)), tx.center)
print("diff: " + str(diff))
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