API Reference

This part of the documentation describes each function, class, and method in detail.

Class structure

The chart below is a class diagram. It shows the inheritance and composition relationships between the classes. For example, Ptycho is an Operator which is composed of Convolution and Propagation (which are also Operators).

_images/api-class-diagram.svg

Modules

operators

Defines reference implementations for all operators.

All of the solvers, rely on core operators including forward and adjoint operators. These core operators may have multiple implementations based on different backends e.g. CUDA, OpenCL, NumPy. The classes in this module prescribe an interface and reference implementation upon which specific solvers are based. In this way, multiple solvers (e.g. ePIE, gradient descent, SIRT) implemented in Python can share the same core operators and can be upgraded to better operators in the future.

All operator methods should take NumPy arrays as inputs. This is a design decision which was made because clients of the operators library should not need to be concerned about memory locality which is necessary complexity when implementating operators for specialized hardware (such as GPUs).

class tike.operators.Operator(array_module=None, asnumpy=None, **kwargs)[source]

Bases: abc.ABC

A base class for Operators.

An Operator is a context manager which provides the basic functions (forward and adjoint) required solve an inverse problem.

array_module

A module which provides a NumPy interface that will be compatible with the arrays passed to the operator.

Type:Python module
asnumpy

This function converts the arrays of the type processed by this operator into NumPy arrays.

Type:function
adj(**kwargs)[source]

Perform the adjoint operator.

fwd(**kwargs)[source]

Perform the forward operator.

class tike.operators.Convolution(probe_shape, nscan, nz, n, ntheta, **kwargs)[source]

Bases: tike.operators.operator.Operator

2D Convolution operator with linear interpolation.

adj(nearplane, scan)[source]

Combine probe shaped patches into a psi shaped grid by addition.

fwd(psi, scan)[source]

Extract probe shaped patches from the psi at each scan position.

The patches within the bounds of psi are linearly interpolated, and indices outside the bounds of psi are not allowed.

class tike.operators.Propagation(nwaves, detector_shape, probe_shape, model='gaussian', **kwargs)[source]

Bases: tike.operators.operator.Operator

A base class for Fourier-based free-space propagation.

adj(farplane, **kwargs)[source]

Adjoint Fourier-based free-space propagation operator.

fwd(nearplane, **kwargs)[source]

Forward Fourier-based free-space propagation operator.

class tike.operators.Ptycho(detector_shape, probe_shape, nscan, nz, n, ntheta=1, model='gaussian', nmodes=1, propagation=<class 'tike.operators.propagation.Propagation'>, diffraction=<class 'tike.operators.convolution.Convolution'>, **kwargs)[source]

Bases: tike.operators.operator.Operator

A base class for ptychography solvers.

This class is a context manager which provides the basic operators required to implement a ptychography solver. Specific implementations of this class can either inherit from this class or just provide the same interface.

Solver implementations should inherit from PtychoBacked which is an alias for whichever Ptycho implementation is selected at import time.

nscan

The number of scan positions at each angular view. (Assumed to be the same for all views.)

Type:int
probe_shape

The pixel width and height of the (square) probe illumination.

Type:int
detector_shape

The pixel width and height of the (square) detector grid.

Type:int
nz, n

The pixel width and height of the reconstructed grid.

Type:int
ntheta

The number of angular partitions of the data.

Type:int
ptheta

The number of angular partitions to process together simultaneously.

Type:int
Parameters:
  • psi ((ntheta, nz, n) complex64) – The complex wavefront modulation of the object.
  • probe ((ntheta, probe_shape, probe_shape) complex64) – The complex illumination function.
  • data, farplane ((ntheta, nscan, detector_shape, detector_shape) complex64) – data is the square of the absolute value of farplane. data is the intensity of the farplane.
  • scan ((ntheta, nscan, 2) float32) – Coordinates of the minimum corner of the probe grid for each measurement in the coordinate system of psi. Vertical coordinates first, horizontal coordinates second.
adj(farplane, probe, scan, **kwargs)[source]

Perform the adjoint operator.

adj_probe(farplane, scan, psi, **kwargs)[source]
fwd(probe, scan, psi, **kwargs)[source]

Perform the forward operator.

class tike.operators.Tomo(angles, ntheta, nz, n, centers)[source]

Bases: tike.operators.operator.Operator

A base class for tomography solvers.

This class is a context manager which provides the basic operators required to implement a tomography solver. Specific implementations of this class can either inherit from this class or just provide the same interface.

Solver implementations should inherit from TomoBackend which is an alias for whichever TomoCore implementation is selected at import time.

ntheta

The number of projections.

Type:int
n, nz

The pixel width and height of the projection.

Type:int
Parameters:
  • obj ((nz, n, n) complex64) – The complex object to be transformed or recovered.
  • tomo ((ntheta, nz, n) complex64) – The radon transform of obj.
  • angles ((ntheta, ) float32) – The radian angles at which the radon transform is sampled.
  • centers ((nz, ) float32) – The center of rotation in obj pixels for each slice along z.
adj(tomo, **kwargs)[source]

Perform the adjoint Radon transform (R*).

See help(TomoCore) for more information.

array_module = None
asnumpy = None
fwd(obj, **kwargs)[source]

Perform the forward Radon transform (R).

See help(TomoCore) for more information.

run(tomo, obj, **kwargs)[source]

Implement a specific tomography solving algorithm.

See help(TomoCore) for more information.

opt

Generic implementations of optimization routines.

Generic implementations of optimization algorithm such as conjugate gradient and line search that can be reused between domain specific modules. In, the future, this module may be replaced by Operator Discretization Library (ODL) solvers library.

tike.opt.conjugate_gradient(array_module, x, cost_function, grad, num_iter=1)[source]

Use conjugate gradient to estimate x.

Parameters:
  • array_module (module) – The Python module that will provide array operations.
  • x (array_like) – The object to be recovered.
  • cost_function (func(x) -> float) – The function being minimized to recover x.
  • grad (func(x) -> array_like) – The gradient of cost_function.
  • num_iter (int) – The number of steps to take.
tike.opt.direction_dy(grad0, grad1, dir)[source]

Return the Dai-Yuan search direction.

Parameters:
  • grad0 (array_like) – The gradient from the previous step.
  • grad1 (array_like) – The gradient from this step.
  • dir (array_like) – The previous search direction.

Return a new step_length using a backtracking line search.

Parameters:
  • f (function(x)) – The function being optimized.
  • x (vector) – The current position.
  • d (vector) – The search direction.
Returns:

  • step_length (float) – The optimal step length along d.
  • cost (float) – The new value of the cost function after stepping along d.

References

https://en.wikipedia.org/wiki/Backtracking_line_search

ptycho

Provides ptychography solvers.

The reference implementation uses NumPy’s FFT library. Select a non-default backend by setting the TIKE_PTYCHO_BACKEND environment variable.

Coordinate Systems

v, h are the horizontal and vertical directions perpendicular to the probe direction where positive directions are to the right and up.

Functions

Each function in this module should have the following interface:

Parameters:
  • data ((T, P, V, H) numpy.array float32) – An array of detector intensities for each of the P positions at T viewing angles. The grid of each detector is H pixels wide (the horizontal direction) and V pixels tall (the vertical direction).
  • probe ((T, P, M, V, H) numpy.array complex64) – The illuminations of the probes.
  • psi ((T, V, H) numpy.array complex64) – The object transmission function.
  • scan ((T, P, 2) numpy.array float32) – The scanning positions with vertical coordinate listed before horizontal coordinates.
  • kwargs (dict) – Keyword arguments specific to this function. **kwargs should always be included so that extra parameters are ignored instead of raising an error.
tike.ptycho.ptycho.gaussian(size, rin=0.8, rout=1.0)[source]

Return a complex gaussian probe distribution.

Illumination probe represented on a 2D regular grid.

A finite-extent circular shaped probe is represented as a complex wave. The intensity of the probe is maximum at the center and damps to zero at the borders of the frame.

Parameters:
  • size (int) – The side length of the distribution
  • rin (float [0, 1) < rout) – The inner radius of the distribution where the dampening of the intensity will start.
  • rout (float (0, 1] > rin) – The outer radius of the distribution where the intensity will reach zero.
tike.ptycho.ptycho.reconstruct(data, probe, scan, psi, algorithm, num_iter=1, rtol=0.001, **kwargs)[source]

Reconstruct the psi and probe using the given algorithm.

Parameters:
  • algorithm (string) – The name of one algorithms to use for reconstructing.
  • rtol (float) – Terminate early if the relative decrease of the cost function is less than this amount.
tike.ptycho.ptycho.simulate(detector_shape, probe, scan, psi, nmode=1, **kwargs)[source]

Propagate the wavefront to the detector.

Return real-valued intensities measured by the detector.

tike.ptycho.solvers.combined(operator, data, probe, scan, psi, recover_psi=True, recover_probe=True, **kwargs)[source]

Solve the ptychography problem using a combined approach.

See also

tike.ptycho.divided

tike.ptycho.solvers.divided(self, data, probe, scan, psi, recover_psi=True, recover_probe=False, recover_positions=False, nmodes=1, **kwargs)[source]

Solve near- and farfield- ptychography problems separately.

References

Michal Odstrcil, Andreas Menzel, and Manuel Guizar-Sicaros. Iteraive least-squares solver for generalized maximum-likelihood ptychography. Optics Express. 2018.

See also

tike.ptycho.combined

tomo

Provides tomography solvers.

The reference implementation uses NumPy’s FFT library. Select a non-default backend by setting the TIKE_TOMO_BACKEND environment variable.

Coordinate Systems

theta, v, h. v, h are the horizontal vertical directions perpendicular to the probe direction where positive directions are to the right and up. theta is the rotation angle around the vertical reconstruction space axis, z. z is parallel to v, and uses the right hand rule to determine reconstruction space coordinates z, x, y. theta is measured from the x axis, so when theta = 0, h is parallel to y.

Functions

Each public function in this module should have the following interface:

Parameters:
  • obj ((Z, X, Y, P) numpy.array float32) – An array of material properties. The first three dimensions Z, X, Y are spatial dimensions. The fourth dimension, P, holds properties at each grid position: refractive indices, attenuation coefficents, etc.
  • integrals ((M, V, H, P) numpy.array float32) – Integrals across the obj for each of the probe rays and P parameters.
  • theta, v, h ((M, ) numpy.array float32) – The min corner (theta, v, h) of the probe for each measurement.
  • kwargs – Keyword arguments specific to this function. **kwargs should always be included so that extra parameters are ignored instead of raising an error.
tike.tomo.tomo.reconstruct(obj, integrals, theta, algorithm=None, num_iter=1, **kwargs)[source]

Reconstruct the obj using the given algorithm.

tike.tomo.tomo.simulate(obj, theta, **kwargs)[source]

Compute line integrals over an obj.

Provides Solver implementations for a variety of algorithms.

class tike.tomo.solvers.ConjugateGradientTomoSolver(angles, ntheta, nz, n, centers)[source]

Bases: tike.operators.tomo.Tomo

Solve the ptychography problem using gradient descent.

adj(tomo, **kwargs)

Perform the adjoint Radon transform (R*).

See help(TomoCore) for more information.

array_module = None
asnumpy = None
fwd(obj, **kwargs)

Perform the forward Radon transform (R).

See help(TomoCore) for more information.

run(tomo, obj, theta, num_iter, rho=1.0, tau=0.0, reg=0j, K=(1+0j), **kwargs)[source]

Use conjugate gradient to estimate obj.

Parameters:
  • tomo (array-like float32) – Line integrals through the object.
  • obj (array-like float32) – The object to be recovered.
  • num_iter (int) – Number of steps to take.
  • rho, tau (float32) – Weights for data and variation components of the cost function
  • reg (complex64) – The regularizer for total variation

view

Define functions for plotting and viewing data of various types.

tike.view.plot_complex(Z, rmin=None, rmax=None, imin=None, imax=None)[source]

Plot real and imaginary parts of a 2D image size by side.

Takes parameters rmin, rmax, imin, imax to scale the ranges of the real and imaginary plots.

tike.view.plot_phase(Z, amin=None, amax=None)[source]

Plot the amplitude and phase of a 2D image side by side.

Takes parameters amin, amax to scale the range of the amplitude. The phase is scaled to the range -pi to pi.

tike.view.trajectory(x, y, connect=True, frame=None, pause=True, dt=1e-12)[source]

Plot a 2D trajectory.

tike.view.plot_footprint(theta, v, h)[source]

Plot 2D projections of the trajectory for each pair of axes.

tike.view.plot_trajectories(theta, v, h, t)[source]

Plot each trajectory as a function of time in the current figure.

Plots two subplots in the current figure. The top one shows horizonal and vertical position as a function of time and the bottom shows angular position as a function of time.

Returns:ax1, ax1b (axes) – Handles to the two axes
tike.view.plot_sino_coverage(theta, v, h, dwell=None, bins=[16, 8, 4], probe_grid=[[1]], probe_shape=(0, 0))[source]

Plot projections of minimum coverage in the sinogram space.