Skip to content

Commit

Permalink
MAINT: update doctests and fix some minor things
Browse files Browse the repository at this point in the history
  • Loading branch information
Holger Kohr committed May 26, 2017
1 parent 9872b75 commit fb50f78
Show file tree
Hide file tree
Showing 15 changed files with 202 additions and 83 deletions.
1 change: 0 additions & 1 deletion examples/tomo/checks/check_axes_cone2d_bp.py
Original file line number Diff line number Diff line change
Expand Up @@ -39,7 +39,6 @@
grid = odl.RectGrid(np.linspace(0, 2 * np.pi, 360, endpoint=False))
angle_partition = odl.uniform_partition_fromgrid(grid)

# Make detector large enough to cover the object
# Make detector large enough to cover the object
src_radius = 500
det_radius = 1000
Expand Down
5 changes: 2 additions & 3 deletions examples/tomo/checks/check_axes_cone2d_fp.py
Original file line number Diff line number Diff line change
Expand Up @@ -6,8 +6,8 @@
Both pairs of plots of ODL projections and NumPy axis sums should look
similar in the sense that they should show the same features in the
right arrangement (not flipped, rotated, etc.), but differing in
magnification.
right arrangement (not flipped, rotated, etc.), up to differences in
magnification and other cone effects.
"""

import matplotlib.pyplot as plt
Expand Down Expand Up @@ -42,7 +42,6 @@
grid = odl.RectGrid([0, np.pi / 2, np.pi, 3 * np.pi / 2])
angle_partition = odl.uniform_partition_fromgrid(grid)

# Make detector large enough to cover the object
# Make detector large enough to cover the object
src_radius = 500
det_radius = 1000
Expand Down
6 changes: 3 additions & 3 deletions examples/tomo/checks/check_axes_cone3d_bp.py
Original file line number Diff line number Diff line change
Expand Up @@ -48,7 +48,7 @@
axis=[0, 0, 1])

# Create projections and back-projection
ray_trafo = odl.tomo.RayTransform(reco_space, geometry)
ray_trafo = odl.tomo.RayTransform(reco_space, geometry, impl='astra_cuda')
proj_data = ray_trafo(phantom)
backproj = ray_trafo.adjoint(proj_data)
backproj.show('Backprojection, axis = [0, 0, 1], middle z slice',
Expand All @@ -65,7 +65,7 @@
axis=[0, 1, 0])

# Create projections and back-projection
ray_trafo = odl.tomo.RayTransform(reco_space, geometry)
ray_trafo = odl.tomo.RayTransform(reco_space, geometry, impl='astra_cuda')
proj_data = ray_trafo(phantom)
backproj = ray_trafo.adjoint(proj_data)
backproj.show('Backprojection, axis = [0, 1, 0], middle y slice',
Expand All @@ -82,7 +82,7 @@
axis=[1, 0, 0])

# Create projections and back-projection
ray_trafo = odl.tomo.RayTransform(reco_space, geometry)
ray_trafo = odl.tomo.RayTransform(reco_space, geometry, impl='astra_cuda')
proj_data = ray_trafo(phantom)
backproj = ray_trafo.adjoint(proj_data)
backproj.show('Backprojection, axis = [1, 0, 0], almost max x slice',
Expand Down
12 changes: 6 additions & 6 deletions examples/tomo/checks/check_axes_cone3d_fp.py
Original file line number Diff line number Diff line change
Expand Up @@ -4,10 +4,10 @@
back-ends, a check is needed to confirm that the translation steps are
done correctly.
All pairs of plots of ODL projections and NumPy axis sums should look
Both pairs of plots of ODL projections and NumPy axis sums should look
similar in the sense that they should show the same features in the
right arrangement (not flipped, rotated, etc.), but differing in
magnification.
right arrangement (not flipped, rotated, etc.), up to differences in
magnification and other cone effects.
"""

import matplotlib.pyplot as plt
Expand Down Expand Up @@ -61,7 +61,7 @@
assert np.allclose(geometry.src_to_det_init, [0, 1, 0])

# Create projections
ray_trafo = odl.tomo.RayTransform(reco_space, geometry)
ray_trafo = odl.tomo.RayTransform(reco_space, geometry, impl='astra_cuda')
proj_data = ray_trafo(phantom)

# Axes in this image are (x, z). This corresponds to
Expand Down Expand Up @@ -107,7 +107,7 @@
assert np.allclose(geometry.src_to_det_init, [0, 0, -1])

# Create projections
ray_trafo = odl.tomo.RayTransform(reco_space, geometry)
ray_trafo = odl.tomo.RayTransform(reco_space, geometry, impl='astra_cuda')
proj_data = ray_trafo(phantom)

# Axes in this image are (x, y). This corresponds to:
Expand Down Expand Up @@ -153,7 +153,7 @@
assert np.allclose(geometry.src_to_det_init, [0, 1, 0])

# Create projections
ray_trafo = odl.tomo.RayTransform(reco_space, geometry)
ray_trafo = odl.tomo.RayTransform(reco_space, geometry, impl='astra_cuda')
proj_data = ray_trafo(phantom)

# Axes in this image are (y, x). This corresponds to
Expand Down
2 changes: 1 addition & 1 deletion examples/tomo/checks/check_axes_parallel2d_fp.py
Original file line number Diff line number Diff line change
Expand Up @@ -64,7 +64,7 @@
assert np.allclose(geometry.det_pos_init, [0, 1])

# Create projections
ray_trafo = odl.tomo.RayTransform(reco_space, geometry, impl='astra_cuda')
ray_trafo = odl.tomo.RayTransform(reco_space, geometry, impl='astra_cpu')
proj_data = ray_trafo(phantom)

# Axis in this image is x. This corresponds to 0 degrees.
Expand Down
6 changes: 3 additions & 3 deletions examples/tomo/checks/check_axes_parallel3d_bp.py
Original file line number Diff line number Diff line change
Expand Up @@ -41,7 +41,7 @@
axis=[0, 0, 1])

# Create projections and back-projection
ray_trafo = odl.tomo.RayTransform(reco_space, geometry)
ray_trafo = odl.tomo.RayTransform(reco_space, geometry, impl='astra_cuda')
proj_data = ray_trafo(phantom)
backproj = ray_trafo.adjoint(proj_data)
backproj.show('Backprojection, axis = [0, 0, 1], middle z slice',
Expand All @@ -57,7 +57,7 @@
axis=[0, 1, 0])

# Create projections and back-projection
ray_trafo = odl.tomo.RayTransform(reco_space, geometry)
ray_trafo = odl.tomo.RayTransform(reco_space, geometry, impl='astra_cuda')
proj_data = ray_trafo(phantom)
backproj = ray_trafo.adjoint(proj_data)
backproj.show('Backprojection, axis = [0, 1, 0], middle y slice',
Expand All @@ -73,7 +73,7 @@
axis=[1, 0, 0])

# Create projections and back-projection
ray_trafo = odl.tomo.RayTransform(reco_space, geometry)
ray_trafo = odl.tomo.RayTransform(reco_space, geometry, impl='astra_cuda')
proj_data = ray_trafo(phantom)
backproj = ray_trafo.adjoint(proj_data)
backproj.show('Backprojection, axis = [1, 0, 0], almost max x slice',
Expand Down
6 changes: 3 additions & 3 deletions examples/tomo/checks/check_axes_parallel3d_fp.py
Original file line number Diff line number Diff line change
Expand Up @@ -56,7 +56,7 @@
assert np.allclose(geometry.det_pos_init, [0, 1, 0])

# Create projections
ray_trafo = odl.tomo.RayTransform(reco_space, geometry)
ray_trafo = odl.tomo.RayTransform(reco_space, geometry, impl='astra_cuda')
proj_data = ray_trafo(phantom)

# Axes in this image are (x, z). This corresponds to
Expand Down Expand Up @@ -101,7 +101,7 @@
assert np.allclose(geometry.det_pos_init, [0, 0, -1])

# Create projections
ray_trafo = odl.tomo.RayTransform(reco_space, geometry)
ray_trafo = odl.tomo.RayTransform(reco_space, geometry, impl='astra_cuda')
proj_data = ray_trafo(phantom)

# Axes in this image are (x, y). This corresponds to:
Expand Down Expand Up @@ -146,7 +146,7 @@
assert np.allclose(geometry.det_pos_init, [0, 1, 0])

# Create projections
ray_trafo = odl.tomo.RayTransform(reco_space, geometry)
ray_trafo = odl.tomo.RayTransform(reco_space, geometry, impl='astra_cuda')
proj_data = ray_trafo(phantom)

# Axes in this image are (y, x). This corresponds to
Expand Down
2 changes: 1 addition & 1 deletion odl/tomo/backends/astra_cuda.py
Original file line number Diff line number Diff line change
Expand Up @@ -275,7 +275,7 @@ def call_backward(self, proj_data, out=None):

# Copy result to CPU memory
if self.geometry.ndim == 2:
# Rotate 90 degrees counter-clockwise from coordinate system
# Rotate 90 degrees clockwise from coordinate system
# (rows, cols) to (x, y)
out[:] = np.rot90(self.out_array, -1)
elif self.geometry.ndim == 3:
Expand Down
6 changes: 3 additions & 3 deletions odl/tomo/backends/astra_setup.py
Original file line number Diff line number Diff line change
Expand Up @@ -237,7 +237,7 @@ def astra_conebeam_3d_geom_to_vec(geometry):
References
----------
.. _astra projection geometry documentation:
.. _ASTRA projection geometry documentation:
http://www.astra-toolbox.com/docs/geom3d.html#projection-geometries
"""
angles = geometry.angles
Expand Down Expand Up @@ -297,7 +297,7 @@ def astra_conebeam_2d_geom_to_vec(geometry):
References
----------
.. _astra projection geometry documentation:
.. _ASTRA projection geometry documentation:
http://www.astra-toolbox.com/docs/geom3d.html#projection-geometries
"""
angles = geometry.angles
Expand Down Expand Up @@ -350,7 +350,7 @@ def astra_parallel_3d_geom_to_vec(geometry):
References
----------
.. _astra projection geometry documentation:
.. _ASTRA projection geometry documentation:
http://www.astra-toolbox.com/docs/geom3d.html#projection-geometries
"""
angles = geometry.angles
Expand Down
57 changes: 49 additions & 8 deletions odl/tomo/geometry/conebeam.py
Original file line number Diff line number Diff line change
Expand Up @@ -43,6 +43,9 @@ class HelicalConeFlatGeometry(DivergentBeamGeometry, AxisOrientedGeometry):
the initial source-to-detector vector is ``(0, 1, 0)``, and the
initial detector axes are ``[(1, 0, 0), (0, 0, 1)]``.
For details, check `the online docs
<https://odlgroup.github.io/odl/guide/geometry_guide.html>`_.
See Also
--------
CircularConeFlatGeometry : Case with zero pitch
Expand Down Expand Up @@ -118,7 +121,6 @@ def __init__(self, apart, dpart, src_radius, det_radius, pitch,
Initialization with default parameters and some (arbitrary)
choices for pitch and radii:
>>> e_x, e_y, e_z = np.eye(3) # standard unit vectors
>>> apart = odl.uniform_partition(0, 4 * np.pi, 10)
>>> dpart = odl.uniform_partition([-1, -1], [1, 1], (20, 20))
>>> geom = HelicalConeFlatGeometry(
Expand All @@ -128,8 +130,12 @@ def __init__(self, apart, dpart, src_radius, det_radius, pitch,
>>> geom.det_refpoint(0)
array([ 0., 10., 0.])
>>> np.allclose(geom.src_position(2 * np.pi),
... geom.src_position(0) + (0, 0, 2)) # z shift due to pitch
... geom.src_position(0) + (0, 0, 2)) # z shift by pitch
True
Checking the default orientation:
>>> e_x, e_y, e_z = np.eye(3) # standard unit vectors
>>> np.allclose(geom.axis, e_z)
True
>>> np.allclose(geom.src_to_det_init, e_y)
Expand Down Expand Up @@ -235,10 +241,8 @@ def __init__(self, apart, dpart, src_radius, det_radius, pitch,
if det_axes_init is None:
vecs_to_transform.extend(default_det_axes_init)

transformed_vecs, matrix = transform_system(
transformed_vecs = transform_system(
axis, default_axis, vecs_to_transform, matrix=init_matrix)
# Store matrix computed from the vectors for reference
self._rot_matrix_from_vecs = matrix if init_matrix is None else None

axis = transformed_vecs[0]
if src_to_det_init is None:
Expand Down Expand Up @@ -341,8 +345,9 @@ def det_refpoint(self, angle):
For an angle ``phi``, the detector position is given by::
ref(phi) = det_rad * rot_matrix(phi) * src_to_det_init +
(pitch_offset + pitch * phi) * axis
det_ref(phi) = translation +
rot_matrix(phi) * (det_rad * src_to_det_init) +
(pitch_offset + pitch * phi) * axis
where ``src_to_det_init`` is the initial unit vector pointing
from source to detector.
Expand All @@ -361,6 +366,21 @@ def det_refpoint(self, angle):
See Also
--------
rotation_matrix
Examples
--------
With default arguments, the detector starts at ``det_rad * e_y``
and rotates to ``det_rad * (-e_x) + pitch/4 * e_z`` at
90 degrees:
>>> apart = odl.uniform_partition(0, 4 * np.pi, 10)
>>> dpart = odl.uniform_partition([-1, -1], [1, 1], (20, 20))
>>> geom = HelicalConeFlatGeometry(
... apart, dpart, src_radius=5, det_radius=10, pitch=2)
>>> np.allclose(geom.det_refpoint(0), [0, 10, 0])
True
>>> np.allclose(geom.det_refpoint(np.pi / 2), [-10, 0, 0.5])
True
"""
angle = float(angle)
if angle not in self.motion_params:
Expand Down Expand Up @@ -405,6 +425,21 @@ def src_position(self, angle):
See Also
--------
rotation_matrix
Examples
--------
With default arguments, the source starts at ``src_rad * (-e_y)``
and rotates to ``src_rad * e_x + pitch/4 * e_z`` at
90 degrees:
>>> apart = odl.uniform_partition(0, 4 * np.pi, 10)
>>> dpart = odl.uniform_partition([-1, -1], [1, 1], (20, 20))
>>> geom = HelicalConeFlatGeometry(
... apart, dpart, src_radius=5, det_radius=10, pitch=2)
>>> np.allclose(geom.src_position(0), [0, -5, 0])
True
>>> np.allclose(geom.src_position(np.pi / 2), [5, 0, 0.5])
True
"""
angle = float(angle)
if angle not in self.motion_params:
Expand Down Expand Up @@ -486,6 +521,9 @@ class CircularConeFlatGeometry(HelicalConeFlatGeometry):
the initial source-to-detector vector is ``(0, 1, 0)``, and the
initial detector axes are ``[(1, 0, 0), (0, 0, 1)]``.
For details, check `the online docs
<https://odlgroup.github.io/odl/guide/geometry_guide.html>`_.
See Also
--------
HelicalConeFlatGeometry : General case with motion in z direction
Expand Down Expand Up @@ -552,7 +590,6 @@ def __init__(self, apart, dpart, src_radius, det_radius, axis=(0, 0, 1),
Initialization with default parameters and some (arbitrary)
choices for the radii:
>>> e_x, e_y, e_z = np.eye(3) # standard unit vectors
>>> apart = odl.uniform_partition(0, 4 * np.pi, 10)
>>> dpart = odl.uniform_partition([-1, -1], [1, 1], (20, 20))
>>> geom = CircularConeFlatGeometry(
Expand All @@ -563,6 +600,10 @@ def __init__(self, apart, dpart, src_radius, det_radius, axis=(0, 0, 1),
array([ 0., 10., 0.])
>>> np.allclose(geom.src_position(2 * np.pi), geom.src_position(0))
True
Checking the default orientation:
>>> e_x, e_y, e_z = np.eye(3) # standard unit vectors
>>> np.allclose(geom.axis, e_z)
True
>>> np.allclose(geom.src_to_det_init, e_y)
Expand Down
Loading

0 comments on commit fb50f78

Please sign in to comment.