Linear Operators (LinOp)

\[\]

This section contains linear operator classes which all derive from the abstract class LinOp.

LinOpBroadcast

class LinOp.LinOpBroadcast

Bases: LinOp

LinOpSum: Broacasting linear operator which broadcast the elements of a variable along given directions. $$\mathrm{H} : \mathrm{x} \mapsto \mathrm{y_{k,l}} = \mathrm{x}_{k} $$

Parameters:

sz – size of $\mathrm{y}$ the output of the LinOpBroadcast.
index – dimensions along which vector will be broadcasted

LinOpBroadcast is the Adjoint of the LinOpSum operator

All attributes of parent class LinOp are inherited.

Example S=LinOpBroadcast(sz,index)

See also LinOp, Map

Method Summary

apply_(this, x): Reimplemented from parent class LinOp. $$\mathrm{H}^* : \mathrm{x} \mapsto \mathrm{y_{k,l}} = \mathrm{x}_{k} \; \forall l$$

applyAdjoint_(this, x): Reimplemented from parent class LinOp. $$\mathrm{H} : \mathrm{x} \mapsto \mathrm{y_k} = \sum_l \mathrm{x}_{k,l} $$

applyHtH_(this, x): Reimplemented from parent class LinOp.

makeAdjoint_(this): Reimplemented from parent class LinOp.

makeHtH_(this): Reimplemented from parent class LinOp.

LinOpBroadcastMatrix

class LinOp.LinOpBroadcastMatrix

Bases: LinOp

LinOpBroadCastMatrix which applies a matrix along a given dimension

Build a LinOp from a given matrix

Parameters:

M – of size M x N
sz – input size of the LinOpBroadCastMatrix
index – dimension along which the matrix is applied (sz(index) must be equal to N)

All attributes of parent class LinOp are inherited.

If M is 2x3, sz=[2 3 2] and index=2, then the matrix M is applied to each x(i,:,j) leading to an output of size [2 2 2].

Example H=LinOpBroadcastMatrix(M,sz,index)

See also LinOp, Map

Method Summary

apply_(this, x): Reimplemented from parent class LinOp.

applyAdjoint_(this, x): Reimplemented from parent class LinOp.

applyInverse_(this, x): Reimplemented from parent class LinOp.

applyHtH_(this, x): Reimplemented from parent class LinOp.

applyHHt_(this, x): Reimplemented from parent class LinOp.

makeAdjoint_(this): Reimplemented from parent class LinOp.

makeHtH_(this): Reimplemented from parent class LinOp.

makeHHt_(this): Reimplemented from parent class LinOp.

makeInversion_(this): Reimplemented from parent class LinOp.

plus_(this, G): Reimplemented from parent class LinOp.

makeComposition_(this, G): Reimplemented from parent class LinOp.

LinOpConv

class LinOp.LinOpConv

Bases: LinOp

LinOpConv: Convolution operator

Parameters:

mtf – Fourier transform of Point Spread Function
isReal – if true (default) the result of the convolution should be real
index – dimensions along which the convolution is performed (the MTF/PSF must have a comformable size
'MTF' – keyword to provide MTF (default)
'PSF' – keyword to provide PSF instead of MTF
'Centered' – if the PSF is centered in the image
'Pad' – if the PSF must be padded to the size SZ with the value padvalue (default 0)
'useRFT' – keyword to use the real-to-half-complex fourier transformation (works with a given ‘PSF’, not ‘MTF’)

All attributes of parent class LinOp are inherited.

Example H=LinOpConv(mtf,isReal,index)

Example H=LinOpConv(‘MTF’,mtf,isReal,index)

Example H=LinOpConv(‘PSF’, psf,isReal,index)

Example H=LinOpConv(‘PSF’, psf,isReal,index,’useRFT’)

Example H=LinOpConv(‘PSF’, psf,isReal,index,’Centered’)

Example H=LinOpConv(‘PSF’, psf,isReal,index,’Pad’,sz,padvalue)

Example H=LinOpConv(‘PSF’, psf,isReal,index,’Centered’,’Pad’,sz,padvalue)

Note the order of keywords do not matters. For example the two following lines are equivalent

H=LinOpConv(‘PSF’, psf,isReal,index,’Pad’,sz,padvalue,’useRFT’)

H=LinOpConv(‘PSF’, psf,isReal,index,’useRFT’,’Pad’,sz,padvalue)

See also LinOp, Map

Method Summary

apply_(this, x): Reimplemented from parent class LinOp.

applyAdjoint_(this, x): Reimplemented from parent class LinOp.

applyHtH_(this, x): Reimplemented from parent class LinOp.

applyHHt_(this, x): Reimplemented from parent class LinOp.

applyInverse_(this, x): Reimplemented from parent class LinOp.

applyAdjointInverse_(this, x): Reimplemented from parent class LinOp.

plus_(this, G): Reimplemented from parent class LinOp.

makeAdjoint_(this): Reimplemented from parent class LinOp.

makeHtH_(this): Reimplemented from parent class LinOp.

makeHHt_(this): Reimplemented from parent class LinOp.

makeInversion_(this): Reimplemented from parent class LinOp

makeComposition_(this, H): Reimplemented from parent class LinOp

LinOpDiag

class LinOp.LinOpDiag

Bases: LinOp

LinOpDiag: Diagonal operator $$ \mathrm{Hx}= \mathrm{\mathbf{diag}(w)x}$$ where $\mathrm{w} \in \mathbb{R}^N$ or $\mathbb{C}^N$ is a vector containing the diagonal elements of $\mathrm{H}$.

Parameters:

diag – elements of the diagonal (Non-singleton dimensions of diag and sz must be consistent)
sz – size (if not the size of the given diag).

All attributes of parent class LinOp are inherited.

If size(diag)=[2 1] and sz=[2 2], then LinOpDiag multiplies each column of a given 2x2 matrix with the diag vector (see bsxfun).

Example D=LinOpDiag(sz,diag)

See also LinOp, Map

Method Summary

apply_(this, x): Reimplemented from parent class LinOp.

applyAdjoint_(this, x): Reimplemented from parent class LinOp.

applyHtH_(this, x): Reimplemented from parent class LinOp.

applyHHt_(this, x): Reimplemented from parent class LinOp.

applyInverse_(this, x): Reimplemented from parent class LinOp.

applyAdjointInverse_(this, x): Reimplemented from parent class LinOp.

plus_(this, G): Reimplemented from parent class LinOp.

makeAdjoint_(this): Reimplemented from parent class LinOp.

makeHHt_(this): Reimplemented from parent class LinOp.

makeHtH_(this): Reimplemented from parent class LinOp.

mpower_(this, p): Reimplemented from parent class LinOp

makeInversion_(this): Reimplemented from parent class LinOp.

makeComposition_(this, G): Reimplemented from parent class LinOp.

LinOpDCT

class LinOp.LinOpDCT

Bases: LinOp

LinOpDCT : Discrete cosine transform

Parameters:: sz – sizein of the operator.

All attributes of parent class LinOp are inherited.

Note Currently only implemented in 2D

Example DCT=LinOpDCT(sz)

LinOpGrad

class LinOp.LinOpGrad

Bases: LinOp

LinOpGrad: Gradient linear operator (Finite differences)

Parameters:

sz – sizein of the gradient operator
index – dimension along which the gradient is computed (all by default)
bc – boundary condition: ‘circular’ (default), ‘zeros’, ‘mirror’
res – vector containing the resolution along each dimension (default all 1)
useRFT – use RFT when defining the LinOpConv associated to $\mathrm{H^TH}$

All attributes of parent class LinOp are inherited.

Note When circular boundary conditions are selected, the method: makeHtH (or equivalently the composition H'*H) returns a convolution linear operator LinOp

Example G = LinOpGrad(sz,index,bc,res)

LinOpHess

class LinOp.LinOpHess

Bases: LinOp

LinOpHess: Hessian linear operator (finite differences)

Parameters:

sz – sizein of the gradient operator
bc – boundary condition: ‘circular’ (default), ‘mirror’
useRFT – use RFT when defining the LinOpConv associated to $\mathrm{H^TH}$
index – dimension along which the hessian is computed (all by default)

The output size is [sz,lidx(lidx+1)/2] where lidx is the length of the index vector. The last dimension of length lidx(lidx+1)/2 represents the coefficients of upper-right part of the symmetric hessian matrix ordered from top to bottom and left to right. For instance:

in 3D with index=[1 2 3], [d^2F/dxx;d^2F/dxy;d^2F/dxz;d^2F/dyy;d^2F/dyz;d^2F/dzz]
in 3D with index=[2 3], [d^2F/dyy;d^2F/dyz;d^2F/dzz]

All attributes of parent class LinOp are inherited.

Note When circular boundary conditions are selected, the method makeHtH (or equivalently the composition H'*H) returns a convolution linear operator LinOp

Example H = LinOpHess(sz,bc,index)

LinOpMatrix

class LinOp.LinOpMatrix

Bases: LinOp

LinOpMatrix : Matrix operator

Build a LinOp from a given matrix

Parameters:

M – matrix
index –

All attributes of parent class LinOp are inherited.

Example H=LinOpMatrix(M,index)

See also LinOp, Map

Method Summary

apply_(this, x): Reimplemented from parent class LinOp.

applyAdjoint_(this, x): Reimplemented from parent class LinOp.

LinOpSDFT

class LinOp.LinOpSDFT

Bases: LinOp

LinOpSDFT : Sliced Discrete Fourier operator

Parameters:

sz – sizein of the operator.
unitary – boolean true when normalized DFT (default false)
pad – padding size (see the doc of fftn function).
index – index along wich dimension are computed the FFT (default all)

All attributes of parent class LinOp are inherited.

Example SDFT=LinOpSDFT(sz,index,unitary, pad)

See also LinOp, Map, fftn, ifftn, Sfft, iSfft

Method Summary

apply_(this, x): Reimplemented from parent class LinOp.

applyAdjoint_(this, x): Reimplemented from parent class LinOp.

applyInverse_(this, x): Reimplemented from parent class LinOp.

applyHtH_(this, x): Reimplemented from parent class LinOp.

applyHHt_(this, x): Reimplemented from parent class LinOp.

applyAdjointInverse_(this, x): Reimplemented from parent class LinOp.

makeHHt_(this): Reimplemented from parent class LinOp.

makeHtH_(this): Reimplemented from parent class LinOp.

makeInversion_(this): Reimplemented from parent class LinOp.

LinOpShape

class LinOp.LinOpShape

Bases: LinOp

LinOpShape: reshaping operator

Reshape an array of size sizein in a array of size sizeout

Parameters:

sizein – input size
sizeout – output size

Example R=LinOpShape(sizein,sizeout)

See also LinOp, Map

Method Summary

apply_(this, x): Reimplemented from parent class LinOp.

applyAdjoint_(this, x): Reimplemented from parent class LinOp.

applyHHt_(~, x): Reimplemented from parent class LinOp.

applyHtH_(~, x): Reimplemented from parent class LinOp.

applyInverse_(this, x): Reimplemented from parent class LinOp.

applyAdjointInverse_(this, x): Reimplemented from parent class LinOp.

makeHHt_(this): Reimplemented from parent class LinOp.

makeHtH_(this): Reimplemented from parent class LinOp.

LinOpSum

class LinOp.LinOpSum

Bases: LinOp

LinOpSum: Summation linear operator which sums the elements of a variable along given directions. $$\mathrm{H} : \mathrm{x} \mapsto \mathrm{y_k} = \sum_l \mathrm{x}_{k,l} $$

Parameters:

sz – size of $\mathrm{x}$ on which the LinOpSum applies.
index – dimensions along which the sum will be performed (inner sum over l)

All attributes of parent class LinOp are inherited.

Example S=LinOpSum(sz,index)

See also LinOp, Map

Method Summary

apply_(this, x): Reimplemented from parent class LinOp.

applyAdjoint_(this, x): Reimplemented from parent class LinOp. $$\mathrm{H}^* : \mathrm{x} \mapsto \mathrm{y_{k,l}} = \mathrm{x}_{k} \; \forall l$$

makeAdjoint_(this): Reimplemented from parent class LinOp.

applyHHt_(this, x): Reimplemented from parent class LinOp.

makeHHt_(this): Reimplemented from parent class LinOp.

makeComposition_(this, H): Reimplemented from parent class LinOp

LinOpXRay

class LinOp.LinOpXRay

Bases: LinOp

LinOpXray: 2D, parallel x-ray operator based on MATLAB’s radon function.

Parameters:

x – structure defining the reconstruction domain (always 2D)
x.size – 1x2 vector giving dimensions of the domain in pixels, must be of the form N*ones(1,2).
x.step – 1x2 vector giving the spacing of the pixels (in an apropirate unit of length, e.g. milimeteres).
thetas – 1x(number of projections) vector giving the projection directions in radians. Specifically, [cos(theta) sin(theta)] points in the direction of the positive axis in the projection coordinate system.

Note We take rows as $$x$$ and columns as $$y$$ in the reconstruction domain. Note The center of the reconstruction (and rotation axis) is assumed to be at floor((x.size+1)/2).

Note The center of each projection is at ceil(sizeout/2).

Note Be warned that the adjoint is only accurate to around 20 dB for for images of size 500x500 px. The error is smaller for larger images.

Method Summary

apply_(this, x): apply the operator

applyAdjoint_(this, c): apply the adjoint

applyHtH_(this, x): apply the HtH matrix

SelectorLinOps

LinOpSelector

class LinOp.SelectorLinOps.LinOpSelector

Bases: LinOp

Selector linear operator which extracts specific entries of a vector $$\mathrm{H} : \mathrm{x} \mapsto \mathrm{x}_{\mathrm{sel}} $$ where $\mathrm{sel}$ is a subset of selected indexes.

Parameters:: sel – boolean with true at selected positions

All attributes of parent class LinOp are inherited.

Note The output of the apply() method will be a colunm vector whatever the given parameter sel. For specific selectors (i.e., downsampling or patch see LinOpDownsample and LinOpSelectorPatch)

Example S=LinOpSelector(sel)

Method Summary

apply_(this, x): Reimplemented from parent class LinOpSelector.

applyAdjoint_(this, x): Reimplemented from parent class LinOpSelector.

applyHtH_(this, x): Reimplemented from parent class LinOpSelector.

applyHHt_(~, x): Reimplemented from parent class LinOp.

makeHHt_(this): Reimplemented from parent class LinOp.

makeHtH_(this): Reimplemented from parent class LinOp.

LinOpDownsample

class LinOp.SelectorLinOps.LinOpDownsample

Bases: LinOp.SelectorLinOps.LinOpSelector

Downsampling linear operator $$\mathrm{H} : \mathrm{x} \mapsto \mathrm{y}_k =\mathrm{x}_{(k-1)d+1} $$ where $d $ is the downsampling factor.

Parameters:

sz – size of $\mathrm{x}$ on which the LinOpDownsample applies.
df – array containing the downsampling factor in each direction (must divide the corresponding entry of sz)
first – array containing the index of the first element in each direction

All attributes of the parent class LinOpSelector are inherited.

Example D=LinOpDownsample(sz,df,first)

Method Summary

makeComposition_(this, H): Reimplemented from parent class LinOp

LinOpSelectorPatch

class LinOp.SelectorLinOps.LinOpSelectorPatch

Bases: LinOp.SelectorLinOps.LinOpSelector

Patch Selector linear operator which extracts a patch from the given vector $$\mathrm{H} : \mathrm{x} \mapsto \mathrm{x}_{[i_{min}:i_{max}]}$$ where $ i_{min} $ and $ i_{max} $ are indexes corresponding to the first and last elements of the patch.

Parameters:

sz – size of $\mathrm{x}$ on which the LinOpSelectorPatch applies.
idxmin – array containing the first kept index in each direction
idxmax – array containing the last kept index in each direction
squeezeflag – true to squeeze singleton dimensions (default: false)

All attributes of parent class LinOpSelector are inherited.

Example S=LinOpSelectorPatch(sz,idxmin,idxmax)

See also LinOp, LinOpSelector, LinOpDownsampling.