edu.stanford.rsl.konrad.volume3d
Class VolumeOperator

java.lang.Object
  extended by edu.stanford.rsl.konrad.volume3d.VolumeOperator
Direct Known Subclasses:
CUDAVolumeOperator, ParallelVolumeOperator

public class VolumeOperator
extends java.lang.Object

Class to model operation on volumes. While this class offers many simple operations which can be applied to volumes like add, and multiply, it also contains a suite of operations tailored to filtering volumes based on quadrature filters. While a quadrature filter is complex in spatial domain, it can often be described using only its real part in frequency domain. Many of the filters used here are based on this feature. The filter functions are modeled as real volumes. However, they are designed to be applied in frequency domain for image processing. All of these filters %preamble{\usepackage{amsmath}} $F_\text{orient}(\mathbf{u})$ are based on a radial formulation:

%preamble{\usepackage{amsmath}} \begin{align*}
$\mathbf{u}$ denotes the coordinate in frequency space and $\rho = |\mathbf{u}|$ is its L2-norm. %preamble{\usepackage{amsmath}} $D_k(\mathbf{u})$ controls the isotropy of the filter. For isotropic filters it is simply 1. For anisotropic filters the direction of anisotropy $\hat{\mathbf{n}}_k$ is required:
%preamble{\usepackage{amsmath}} \begin{align*}
The VolumeOperator contains also methods to obtain directionally well spaced $\hat{\mathbf{n}}_k$.

Note that any method which starts with create will allocate memory. A call of Volume3D.destory() is recommended after it's use. This will free the memory in device-dependent implementations.

Author:
akmaier

Nested Class Summary
static class VolumeOperator.FILTER_TYPE
          Sets the behavior of the filter.
 
Constructor Summary
VolumeOperator()
           
 
Method Summary
 int abs(Volume3D vol)
          Determines the absolute volume of the input volume.
If the input volume is real Math.abs() is called for each element.
If the input volume is complex the power spectrum is computed.
 int addScalar(Volume3D vol, float realPart, float imagPart)
          Adds a scalar to the Volume3D.
 int addVolume(Volume3D vol1, Volume3D vol2)
          Adds the second volume to the first volume.
 int addVolume(Volume3D vol1, Volume3D vol2, double weight)
          Adds the second volume %preamble{\usepackage{amsmath}} $y_j$ multiplied by the scalar weight %preamble{\usepackage{amsmath}} $s$ to the first volume %preamble{\usepackage{amsmath}} $x_j$.
The first volume is overwritten.
 Volume3D createDirectionalWeights(int dimensions, int[] size, float[] dim, float[] dir, int A, VolumeOperator.FILTER_TYPE t_filt)
          Creates an anisotropic, i.e.
 Volume3D createExponentialDirectionalHighPassFilter(int dimensions, int[] size, float[] dim, float[] dir, int A, float B, float ri, VolumeOperator.FILTER_TYPE t_filt)
          Creates an radially symetric anisotropic quadrature filter according to this definition: %preamble{\\usepackage{amsmath}} \\begin{align*} F(\\mathbf{u}) & = \\left \\{ \\begin{array}{ll} \\left(\\frac{\\displaystyle (\\mathbf{u} \\cdot \\hat{\\mathbf{n}}_k)}{\\displaystyle |\\rho|}\\right)^{2A} \\cdot R(\\rho)& \\text{if}\\qquad\\mathbf{u} \\cdot \\hat{\\mathbf{n}}_k >0 \\\\ 0 & \\text{else} \\end{array} \\right .
 Volume3D createGaussLowPassFilter(int dimensions, int[] size, float[] dim, float alpha)
          Creates an isotropic, i.e.
 Volume3D createHighPassFilter(int dimensions, int[] size, float[] dim, int filt_loop, float lp_upper)
          Creates a directional high pass filter %preamble{\usepackage{amsmath}} $F_{\displaystyle\text{HP}_k}(\mathbf{u}) = 1 - (R_\text{LP}(\rho) \cdot D_k(\mathbf{u}))$.
 Volume3D createLowPassFilter(int dimensions, int[] size, float[] dim, float lp_upper)
          Creates an isotropic low-pass filter, i.e.
 Volume3D createVolume(ij.ImagePlus image, int mirror, int cuty, boolean uneven)
          Creates a Volume3D Object from an ImagePlus.
 Volume3D createVolume(int[] size, float[] dim, int in_dim)
          Creates a new empty volume according to the parameters.
 int divideByVolume(Volume3D vol1, Volume3D vol2)
          Divides the first volume by the second volume element by element
 int divideScalar(Volume3D vol, float realPart, float imagPart)
          Divides the volume by a scalar.
 int fftShift(Volume3D vol)
          Performs the shift required for the FFT, i.e.
static void getFrequencyBoundings(int dimensions, int[] size, float[] dim, float[] f_max, float[] f_delta)
          Computes the maximal frequency for each dimension as %preamble{\usepackage{amsmath}} $f_\text{max} = \frac{\pi}{d_i}$, where $d_i$ is the size of each voxel in dimension $i$ and the size of each frequency bin %preamble{\usepackage{amsmath}} $f_{\delta} = \frac{2 \cdot f_\text{max}}{l_i}$, where $l_i$ is the number of voxels in dimension $i$.
Method currently applies C-style passing of return values.
 int imag(Volume3D vol)
          Maps a complex volume onto its imaginary part.
 void makeComplex(Volume3D vol)
          Makes the volume complex, i.e.
 float max(Volume3D vol)
          Determines the maximum intensity of a pixel in the given volume.
 float mean(Volume3D vol)
          Determines the arithmetic average %preamble{\usepackage{amsmath}} $\mu$ of all %preamble{\usepackage{amsmath}} $N = \prod_i d_i$ voxels %preamble{\usepackage{amsmath}} $x_j$.
 float min(Volume3D vol)
          Determines the minimum intensity of the volume.
 int minOfTwoVolumes(Volume3D vol1, Volume3D vol2)
          Determines the minimal volume element by element.
The output is stored in the first volume.
 int multiply(Volume3D vol1, Volume3D vol2)
          Multiplies two volumes element by element.
 int multiplyScalar(Volume3D vol, float realPart, float imagPart)
          Multiplies the volume by a scalar.
 int real(Volume3D vol)
          Maps volume onto its real part.
 int sigmoid(Volume3D vol, float smoothing, float lowValue, float highValue, float highPassLowerLevel, float highPassUpperLevel)
          Remaps gray values using a sigmoid function.
 Volume3D solveMaximumEigenvalue(Volume3D[][] structureTensor)
          Method to compute the maximal eigenvalue %preamble{\usepackage{amsmath}} $\lambda_1$ per volume element for a volume structure tensor.
 int sqrt(Volume3D vol)
          Replaces every element in the volume with the output of Math.sqrt(), i.e.
 int subtractVolume(Volume3D vol1, Volume3D vol2)
          Subtracts the second volume from the first volume.
 int upperLimit(Volume3D vol, float max)
          Iterates the volume and replaces all entries greater than max with max.
 
Methods inherited from class java.lang.Object
equals, getClass, hashCode, notify, notifyAll, toString, wait, wait, wait
 

Constructor Detail

VolumeOperator

public VolumeOperator()
Method Detail

getFrequencyBoundings

public static void getFrequencyBoundings(int dimensions,
                                         int[] size,
                                         float[] dim,
                                         float[] f_max,
                                         float[] f_delta)
Computes the maximal frequency for each dimension as %preamble{\usepackage{amsmath}} $f_\text{max} = \frac{\pi}{d_i}$, where $d_i$ is the size of each voxel in dimension $i$ and the size of each frequency bin %preamble{\usepackage{amsmath}} $f_{\delta} = \frac{2 \cdot f_\text{max}}{l_i}$, where $l_i$ is the number of voxels in dimension $i$.
Method currently applies C-style passing of return values. An array of float values is passed to the method which is filled throughout the call of the method.

Parameters:
dimensions - the dimension of the volume
size - an array containing $l_i$
dim - an array containing $d_i$
f_max - the return values for %preamble{\usepackage{amsmath}} $f_\text{max}$
f_delta - the return values for %preamble{\usepackage{amsmath}} $f_{\delta}$

solveMaximumEigenvalue

public Volume3D solveMaximumEigenvalue(Volume3D[][] structureTensor)
Method to compute the maximal eigenvalue %preamble{\usepackage{amsmath}} $\lambda_1$ per volume element for a volume structure tensor.

Parameters:
structureTensor - the structure tensor as 2D array of volumes.
Returns:
the values of %preamble{\usepackage{amsmath}} $\lambda_1$ as volume

createDirectionalWeights

public Volume3D createDirectionalWeights(int dimensions,
                                         int[] size,
                                         float[] dim,
                                         float[] dir,
                                         int A,
                                         VolumeOperator.FILTER_TYPE t_filt)
Creates an anisotropic, i.e. %preamble{\usepackage{amsmath}} $D(\mathbf{u})$ is not constant, filter as Volume3D if FILTER_TYPE is QUADRATIC.
%preamble{\usepackage{amsmath}} \begin{align*}
$\rho = |\mathbf{u}|$ and $\hat{\mathbf{n}}_k$ is the direction of the anisotropy. The case selection is only valid in the QUADRATIC case.

Parameters:
dimensions - dimension should be 3
size - number of pixels per dimension
dim - resolution
dir - direction vector $\hat{\mathbf{n}}_k$
A - parameter A
t_filt - filter type
Returns:
the filter

createVolume

public Volume3D createVolume(int[] size,
                             float[] dim,
                             int in_dim)
Creates a new empty volume according to the parameters. Volume is initialized with zero at all elements.

Parameters:
size - the sizes in each direction as array
dim - the dimension of the pixel in mm in each direction
in_dim - the internal dimension (1 = real, 2 = complex)
Returns:
the new volume

createVolume

public Volume3D createVolume(ij.ImagePlus image,
                             int mirror,
                             int cuty,
                             boolean uneven)
Creates a Volume3D Object from an ImagePlus. Parameters match the ImagePlus constructor of Volume3D. If a different type of Volume3D is used (such as CUDAVolume3D) this method is overloaded in the respective inhereted version of of CUDAVolumeOperator. In order to write CUDA compatible code, use this version instead of the Volume3D constructor.

Parameters:
image - the ImagePlus with the image data
mirror - the number of pixels which are mirrored at the boundary
cuty - number of pixels to be cut from the ImagePlus (in order to remove emtpy areas at the top and the bottom of the image
uneven - true if the original number of slices was uneven (may have been altered by cutting y)
Returns:
the volume as Volume3D
See Also:
Volume3D#Volume3D(ImagePlus, double, double, double, int, int, boolean)

createExponentialDirectionalHighPassFilter

public Volume3D createExponentialDirectionalHighPassFilter(int dimensions,
                                                           int[] size,
                                                           float[] dim,
                                                           float[] dir,
                                                           int A,
                                                           float B,
                                                           float ri,
                                                           VolumeOperator.FILTER_TYPE t_filt)
Creates an radially symetric anisotropic quadrature filter according to this definition:
%preamble{\usepackage{amsmath}} \begin{align*}

The filter is suitable to estimate the local orientation in direction $\hat{\mathbf{n}}_k$.

This plot shows $R(\rho)$ with B = 2, A = 1, and $\rho_i$ = 1.5:



Note that the volume is moved to FFT octants with the highest frequencies in the center. It can be directly multiplied to the FFT of a volume.

Parameters:
dimensions - dimension should be 3
size - number of pixels per dimension
dim - resolution
dir - the filter direction $\hat{\mathbf{n}}_k$
A - factor A
B - the relative bandwith of the filter B
ri - $\rho_i$
t_filt - parameter to
Returns:
the filter as volume
See Also:
fftShift(Volume3D)

createHighPassFilter

public Volume3D createHighPassFilter(int dimensions,
                                     int[] size,
                                     float[] dim,
                                     int filt_loop,
                                     float lp_upper)
Creates a directional high pass filter %preamble{\usepackage{amsmath}} $F_{\displaystyle\text{HP}_k}(\mathbf{u}) = 1 - (R_\text{LP}(\rho) \cdot D_k(\mathbf{u}))$. The definition of %preamble{\usepackage{amsmath}} $D_k(\mathbf{u})$ from createDirectionalFilter is used. The low-pass filter is defined as:
%preamble{\usepackage{amsmath}} \begin{align*}

%preamble{\usepackage{amsmath}} $\rho_\text{LP}$ is the upper limit of the low-pass filter.

This plot shows the low-pass filter for %preamble{\usepackage{amsmath}} $\rho_\text{LP}$ = 1.5:


Note that the volume is moved to FFT octants with the highest frequencies in the center. It can be directly multiplied to the FFT of a volume.

Parameters:
dimensions - dimension should be 3
size - number of pixels per dimension
dim - resolution
filt_loop - $k$ for the creation of $\hat{\mathbf{n}}_k$
lp_upper - %preamble{\usepackage{amsmath}} $\rho_\text{LP}$
Returns:
the filter as volume
See Also:
fftShift(Volume3D), createDirectionalWeights(int, int[], float[], float[], int, FILTER_TYPE)

createLowPassFilter

public Volume3D createLowPassFilter(int dimensions,
                                    int[] size,
                                    float[] dim,
                                    float lp_upper)
Creates an isotropic low-pass filter, i.e. %preamble{\usepackage{amsmath}} $D_k(\mathbf{u})$ = 1. The radial shape of the filter is determined by:
%preamble{\usepackage{amsmath}} \begin{align*}

This plot shows the low-pass filter for %preamble{\usepackage{amsmath}} $\rho_\text{LP}$ = 1.5:


Note that the volume is moved to FFT octants with the highest frequencies in the center. It can be directly multiplied to the FFT of a volume.

Parameters:
dimensions - dimension should be 3
size - number of pixels per dimension
dim - resolution
lp_upper - %preamble{\usepackage{amsmath}} $\rho_\text{LP}$
Returns:
the filter as volume
See Also:
fftShift(Volume3D)

createGaussLowPassFilter

public Volume3D createGaussLowPassFilter(int dimensions,
                                         int[] size,
                                         float[] dim,
                                         float alpha)
Creates an isotropic, i.e. %preamble{\usepackage{amsmath}} $D_k(\mathbf{u})$ = 1, 3-D Gaussian filter as Volume. The radial shape is determined as:
%preamble{\usepackage{amsmath}} \begin{align*}


This plot shows the radial shape for %preamble{\usepackage{amsmath}} $\alpha$ = 0.5, 1.0, 1.5:


Note that the volume is moved to FFT octants with the highest frequencies in the center. It can be directly multiplied to the FFT of a volume.

See Also:
fftShift(Volume3D)

mean

public float mean(Volume3D vol)
Determines the arithmetic average %preamble{\usepackage{amsmath}} $\mu$ of all %preamble{\usepackage{amsmath}} $N = \prod_i d_i$ voxels %preamble{\usepackage{amsmath}} $x_j$.
%preamble{\usepackage{amsmath}} \begin{align*}

Parameters:
vol - the volume
Returns:
the mean value %preamble{\usepackage{amsmath}} $\mu$

max

public float max(Volume3D vol)
Determines the maximum intensity of a pixel in the given volume.

Parameters:
vol - the volume
Returns:
the maximal value

min

public float min(Volume3D vol)
Determines the minimum intensity of the volume.

Parameters:
vol - the volume
Returns:
the minimal value

multiplyScalar

public int multiplyScalar(Volume3D vol,
                          float realPart,
                          float imagPart)
Multiplies the volume by a scalar. If the volume is real the imaginary part of the scalar should be set to 0. Any value different from 0 will cause the volume to be converted to complex values.

Parameters:
vol - the volume
realPart - the real part of the scalar
imagPart - the imaginary part of the scalar
Returns:
0, if successful
See Also:
makeComplex(Volume3D)

makeComplex

public void makeComplex(Volume3D vol)
Makes the volume complex, i.e. the internal dimension of the volume is increased to 2, if required and the data is transfered to interleaved complex format.

Parameters:
vol - the volume

divideScalar

public int divideScalar(Volume3D vol,
                        float realPart,
                        float imagPart)
Divides the volume by a scalar. If a real volume is multiplied the imaginary parameter should be set to 0. Otherwise, the volume is made complex.

Parameters:
vol - the volume
realPart - the real part of the scalar
imagPart - the imaginary part of the scalar
Returns:
0, if successful
See Also:
makeComplex(Volume3D)

addScalar

public int addScalar(Volume3D vol,
                     float realPart,
                     float imagPart)
Adds a scalar to the Volume3D. If the volume is real the imaginary part should be set to 0. Otherwise, the volume is made complex.

Parameters:
vol - the volume
realPart - the real part
imagPart - the imaginary part
Returns:
0, if successful.
See Also:
makeComplex(Volume3D)

multiply

public int multiply(Volume3D vol1,
                    Volume3D vol2)
Multiplies two volumes element by element. Volume vol1 is overwritten.

Parameters:
vol1 - the first volume
vol2 - the second volume
Returns:
0, if successful

divideByVolume

public int divideByVolume(Volume3D vol1,
                          Volume3D vol2)
Divides the first volume by the second volume element by element

Parameters:
vol1 - the first volume
vol2 - the second volume
Returns:
0, if successful

addVolume

public int addVolume(Volume3D vol1,
                     Volume3D vol2)
Adds the second volume to the first volume.

Parameters:
vol1 - the first volume
vol2 - the second volume
Returns:
0, if successful

addVolume

public int addVolume(Volume3D vol1,
                     Volume3D vol2,
                     double weight)
Adds the second volume %preamble{\usepackage{amsmath}} $y_j$ multiplied by the scalar weight %preamble{\usepackage{amsmath}} $s$ to the first volume %preamble{\usepackage{amsmath}} $x_j$.
The first volume is overwritten. The second volume remains unchanged:
%preamble{\usepackage{amsmath}} \begin{align*}

Parameters:
vol1 - the first volume %preamble{\usepackage{amsmath}} $x_j$
vol2 - the second volume %preamble{\usepackage{amsmath}} $y_j$
weight - the weighting factor %preamble{\usepackage{amsmath}} $s$
Returns:
0, if successful

subtractVolume

public int subtractVolume(Volume3D vol1,
                          Volume3D vol2)
Subtracts the second volume from the first volume. First volume is overwritten.

Parameters:
vol1 - the first volume
vol2 - the second volume
Returns:
0, if successful

minOfTwoVolumes

public int minOfTwoVolumes(Volume3D vol1,
                           Volume3D vol2)
Determines the minimal volume element by element.
The output is stored in the first volume.

Parameters:
vol1 - the first volume
vol2 - the second volume
Returns:
0, if successful

abs

public int abs(Volume3D vol)
Determines the absolute volume of the input volume.
If the input volume is real Math.abs() is called for each element.
If the input volume is complex the power spectrum is computed.
In any case, the resulting volume is real and has internal dimension 1.
Method works in place and overwrites the old input volume.

Parameters:
vol - the input volume.
Returns:
0, if successful

sqrt

public int sqrt(Volume3D vol)
Replaces every element in the volume with the output of Math.sqrt(), i.e. each element's square root.

Parameters:
vol - the volume to be processed.
Returns:
0, if successful

real

public int real(Volume3D vol)
Maps volume onto its real part.

Parameters:
vol - the volume to be mapped.
Returns:
0, if successful

fftShift

public int fftShift(Volume3D vol)
Performs the shift required for the FFT, i.e. moves the high frequencies to the center and the low frequencies to the periphery. If called again on the same volume, its changes are reversed.

Parameters:
vol - the volume
Returns:
0, if successful

sigmoid

public int sigmoid(Volume3D vol,
                   float smoothing,
                   float lowValue,
                   float highValue,
                   float highPassLowerLevel,
                   float highPassUpperLevel)
Remaps gray values using a sigmoid function. A smoothing factor for the slope of the sigmoid can be set. lowValue and highValue determine the input values into the sigmoid. highPaddLowerLevel and highPassUpperLevel determine the output levels of the sigmoid.
The following figure describes the remapping process:

One configuration in the plot remaps the values 5 to 10 to the range from 0 to 1. If the levels are adjusted, the same range can also be remapped to the range from 0.2 to 1.5.

Parameters:
vol - the volume
smoothing - slope of sigmoid
lowValue - the lower input limit
highValue - the upper input limit
highPassLowerLevel - the lower output value
highPassUpperLevel - the upper output value
Returns:
0, if successful

upperLimit

public int upperLimit(Volume3D vol,
                      float max)
Iterates the volume and replaces all entries greater than max with max.

Parameters:
vol - the volume
max - the maximum
Returns:
0, if successful

imag

public int imag(Volume3D vol)
Maps a complex volume onto its imaginary part.

Parameters:
vol - the volume
Returns:
0, if successful