Skip to content

Commit

Permalink
ENH: implement tensor-valued DiscreteLp
Browse files Browse the repository at this point in the history
Changes in detail:
- Add dtype with shape to DiscreteLp (mostly __repr__, factory
  functions and some downstream methods). As a consequence,
  `shape_[in,out]` and `ndim_[in,out]` are added for the
  different types of axes, as well as `scalar_dtype`.
- Add `PerAxisWeighting` and make it the default for
  `DiscreteLp` type spaces. Reason: this way the `tspace`
  knows how to deal with removed axes etc. This is important
  for a smooth experience with indexing and reductions over
  axes.
  Helpers for slicing weightings help structure this task.
- Implement `__getitem__` for `TensorSpace` and `DiscreteLp`,
  including (hopefully) reasonable propagation of weights.
  The new `simulate_slicing` utility function simplifies
  this task.
- Allow indexing with ODL tensors of boolean or integer dtype.
- Implement correct weighting for backprojections with
  non-uniform angles, using per-axis weighting and a new
  helper `adjoint_weightings` to apply the weightings in an
  efficient way.
  The correct weighting from the range of `RayTransform`
  is determined by the new `proj_space_weighting` helper.
- Change the space `_*_impl` methods to always expect and
  return Numpy arrays, and adapt the calling code.
- Change behavior of `norm` and `dist` to ignoring weights
  for `exponent=inf`, in accordance with math.
- Improve speed of `all_equal` for comparison of arrays.
- Account for `None` entries in indices in the
  `normalized_index_expression` helper, thus allowing
  creation of new axes.
- Remove `dicsr_sequence_space`, it was largely unused and
  just a maintenance burden. Use a regular `uniform-discr`
  from zero to `shape` instead.
- Remove `Weighting.equiv()` mehtods, never used and hard
  to maintain (n^2 possibilities).
- Remove the (largely useless) `_weighting` helper to create
  weighting instances since it would have been ambiguous
  with sequences of scalars (array or per axis?). Also remove
  the `npy_weighted_*` functions, they were useless, too.
- Remove some dead code from tomo/util.
- A bunch of minor fixes, as usual.

Closes: odlgroup#908, odlgroup#907, odlgroup#1113, odlgroup#965, odlgroup#286, odlgroup#267, odlgroup#1001
  • Loading branch information
Holger Kohr committed Dec 11, 2017
1 parent 9ba1c8d commit 53d02c4
Show file tree
Hide file tree
Showing 22 changed files with 1,783 additions and 817 deletions.
2 changes: 1 addition & 1 deletion examples/tomo/ray_trafo_parallel_2d.py
Original file line number Diff line number Diff line change
Expand Up @@ -28,7 +28,7 @@
# projection data (or any element in the projection space).
backproj = ray_trafo.adjoint(proj_data)

# Shows a slice of the phantom, projections, and reconstruction
# Show a slice of the phantom, projections, and reconstruction
phantom.show(title='Phantom')
proj_data.show(title='Projection data (sinogram)')
backproj.show(title='Back-projected data', force_show=True)
9 changes: 5 additions & 4 deletions odl/discr/discr_mappings.py
Original file line number Diff line number Diff line change
Expand Up @@ -79,10 +79,11 @@ def __init__(self, map_type, fspace, partition, tspace, linear=False):
'of the function set {}'
''.format(partition, fspace.domain, fspace))

if tspace.shape != partition.shape:
raise ValueError('`tspace.shape` not equal to `partition.shape`: '
'{} != {}'
''.format(tspace.shape, partition.shape))
if tspace.shape != fspace.out_shape + partition.shape:
raise ValueError(
'`tspace.shape` not equal to '
'`fspace.out_shape + partition.shape`: {} != {}'
''.format(tspace.shape, partition.shape))

domain = fspace if map_type == 'sampling' else tspace
range = tspace if map_type == 'sampling' else fspace
Expand Down
13 changes: 10 additions & 3 deletions odl/discr/discretization.py
Original file line number Diff line number Diff line change
Expand Up @@ -375,6 +375,13 @@ def __eq__(self, other):
def __getitem__(self, indices):
"""Return ``self[indices]``.
.. note::
This is a basic implementation that does little more than
propagating the indexing to the `tensor` attribute. In
particular, the result will **not** be wrapped as a
`DiscretizedSpaceElement` again. Subclasses should take
care of that by overriding ``__getitem__``.
Parameters
----------
indices : int or `slice`
Expand Down Expand Up @@ -417,9 +424,9 @@ def sampling(self, ufunc, **kwargs):
Parameters
----------
ufunc : ``self.space.fspace`` element
The continuous function that should be samplingicted.
The continuous function that should be sampled.
kwargs :
Additional arugments for the sampling operator implementation
Additional arguments for the sampling operator implementation.
Examples
--------
Expand All @@ -445,7 +452,7 @@ def interpolation(self):
Returns
-------
interpolation_op : `FunctionSpaceMapping`
Operatior representing a continuous interpolation of this
Operator representing a continuous interpolation of this
element.
Examples
Expand Down
Loading

0 comments on commit 53d02c4

Please sign in to comment.