public class DenseDoubleMatrix3D extends DoubleMatrix3D
Implementation:
Internally holds one single contiguous one-dimensional array, addressed in (in decreasing order of significance): slice major, row major, column major. Note that this implementation is not synchronized.
Time complexity:
O(1) (i.e. constant time) for the basic operations get, getQuick, set, setQuick and size,
Applications demanding utmost speed can exploit knowledge about the internal addressing. Setting/getting values in a loop slice-by-slice, row-by-row, column-by-column is quicker than, for example, column-by-column, row-by-row, slice-by-slice. Thus
for (int slice = 0; slice < slices; slice++) {
for (int row = 0; row < rows; row++) {
for (int column = 0; column < columns; column++) {
matrix.setQuick(slice, row, column, someValue);
}
}
}
is quicker than
for (int column = 0; column < columns; column++) {
for (int row = 0; row < rows; row++) {
for (int slice = 0; slice < slices; slice++) {
matrix.setQuick(slice, row, column, someValue);
}
}
}
| Constructor and Description |
|---|
DenseDoubleMatrix3D(double[][][] values)
Constructs a matrix with a copy of the given values.
|
DenseDoubleMatrix3D(int slices,
int rows,
int columns)
Constructs a matrix with a given number of slices, rows and columns.
|
DenseDoubleMatrix3D(int slices,
int rows,
int columns,
double[] elements,
int sliceZero,
int rowZero,
int columnZero,
int sliceStride,
int rowStride,
int columnStride,
boolean isView)
Constructs a matrix with the given parameters.
|
| Modifier and Type | Method and Description |
|---|---|
double |
aggregate(DoubleDoubleFunction aggr,
DoubleFunction f)
Applies a function to each cell and aggregates the results.
|
double |
aggregate(DoubleDoubleFunction aggr,
DoubleFunction f,
DoubleProcedure cond)
Applies a function to each cell that satisfies a condition and aggregates
the results.
|
double |
aggregate(DoubleDoubleFunction aggr,
DoubleFunction f,
IntArrayList sliceList,
IntArrayList rowList,
IntArrayList columnList)
Applies a function to all cells with a given indexes and aggregates the
results.
|
double |
aggregate(DoubleMatrix3D other,
DoubleDoubleFunction aggr,
DoubleDoubleFunction f)
Applies a function to each corresponding cell of two matrices and
aggregates the results.
|
DoubleMatrix3D |
assign(double value)
Sets all cells to the state specified by value.
|
DoubleMatrix3D |
assign(double[] values)
Sets all cells to the state specified by values.
|
DoubleMatrix3D |
assign(double[][][] values)
Sets all cells to the state specified by values.
|
DoubleMatrix3D |
assign(DoubleFunction function)
Assigns the result of a function to each cell;
x[slice,row,col] = function(x[slice,row,col]).
|
DoubleMatrix3D |
assign(DoubleMatrix3D source)
Replaces all cell values of the receiver with the values of another
matrix.
|
DoubleMatrix3D |
assign(DoubleMatrix3D y,
DoubleDoubleFunction function)
Assigns the result of a function to each cell;
x[row,col] = function(x[row,col],y[row,col]).
|
DoubleMatrix3D |
assign(DoubleMatrix3D y,
DoubleDoubleFunction function,
IntArrayList sliceList,
IntArrayList rowList,
IntArrayList columnList)
Assigns the result of a function to all cells with a given indexes
|
DoubleMatrix3D |
assign(DoubleProcedure cond,
double value)
Assigns a value to all cells that satisfy a condition.
|
DoubleMatrix3D |
assign(DoubleProcedure cond,
DoubleFunction f)
Assigns the result of a function to all cells that satisfy a condition.
|
int |
cardinality()
Returns the number of cells having non-zero values; ignores tolerance.
|
void |
dct2Slices(boolean scale)
Computes the 2D discrete cosine transform (DCT-II) of each slice of this
matrix.
|
void |
dct3(boolean scale)
Computes the 3D discrete cosine transform (DCT-II) of this matrix.
|
void |
dht2Slices()
Computes the 2D discrete Hartley transform (DHT) of each slice of this
matrix.
|
void |
dht3()
Computes the 3D discrete Hartley transform (DHT) of this matrix.
|
void |
dst2Slices(boolean scale)
Computes the 2D discrete sine transform (DST-II) of each slice of this
matrix.
|
void |
dst3(boolean scale)
Computes the 3D discrete sine transform (DST-II) of this matrix.
|
double[] |
elements()
Returns the elements of this matrix.
|
void |
fft3()
Computes the 3D discrete Fourier transform (DFT) of this matrix.
|
DenseDComplexMatrix3D |
getFft2Slices()
Returns new complex matrix which is the 2D discrete Fourier transform
(DFT) of each slice of this matrix.
|
DenseDComplexMatrix3D |
getFft3()
Returns new complex matrix which is the 3D discrete Fourier transform
(DFT) of this matrix.
|
DenseDComplexMatrix3D |
getIfft2Slices(boolean scale)
Returns new complex matrix which is the 2D inverse of the discrete
Fourier transform (IDFT) of each slice of this matrix.
|
DenseDComplexMatrix3D |
getIfft3(boolean scale)
Returns new complex matrix which is the 3D inverse of the discrete
Fourier transform (IDFT) of this matrix.
|
double[] |
getMaxLocation()
Return maximum value of this matrix together with its location
|
double[] |
getMinLocation()
Returns minimum value of this matrix together with its location
|
void |
getNegativeValues(IntArrayList sliceList,
IntArrayList rowList,
IntArrayList columnList,
DoubleArrayList valueList)
Fills the coordinates and values of cells having negative values into the
specified lists.
|
void |
getNonZeros(IntArrayList sliceList,
IntArrayList rowList,
IntArrayList columnList,
DoubleArrayList valueList)
Fills the coordinates and values of cells having non-zero values into the
specified lists.
|
void |
getPositiveValues(IntArrayList sliceList,
IntArrayList rowList,
IntArrayList columnList,
DoubleArrayList valueList)
Fills the coordinates and values of cells having positive values into the
specified lists.
|
double |
getQuick(int slice,
int row,
int column)
Returns the matrix cell value at coordinate [slice,row,column].
|
void |
idct2Slices(boolean scale)
Computes the 2D inverse of the discrete cosine transform (DCT-III) of
each slice of this matrix.
|
void |
idct3(boolean scale)
Computes the 3D inverse of the discrete cosine transform (DCT-III) of
this matrix.
|
void |
idht2Slices(boolean scale)
Computes the 2D inverse of the discrete Hartley transform (IDHT) of each
slice of this matrix.
|
void |
idht3(boolean scale)
Computes the 3D inverse of the discrete Hartley transform (IDHT) of this
matrix.
|
void |
idst2Slices(boolean scale)
Computes the 2D inverse of the discrete sine transform (DST-III) of each
slice of this matrix.
|
void |
idst3(boolean scale)
Computes the 3D inverse of the discrete sine transform (DST-III) of this
matrix.
|
void |
ifft3(boolean scale)
Computes the 3D inverse of the discrete Fourier transform (IDFT) of this
matrix.
|
long |
index(int slice,
int row,
int column)
Returns the position of the given coordinate within the (virtual or
non-virtual) internal 1-dimensional array.
|
DoubleMatrix3D |
like(int slices,
int rows,
int columns)
Construct and returns a new empty matrix of the same dynamic type
as the receiver, having the specified number of slices, rows and columns.
|
DoubleMatrix2D |
like2D(int rows,
int columns)
Construct and returns a new 2-d matrix of the corresponding dynamic
type, sharing the same cells.
|
void |
setQuick(int slice,
int row,
int column,
double value)
Sets the matrix cell at coordinate [slice,row,column] to the
specified value.
|
double[][][] |
toArray()
Constructs and returns a 2-dimensional array containing the cell values.
|
DoubleMatrix1D |
vectorize()
Returns a vector obtained by stacking the columns of each slice of the
matrix on top of one another.
|
void |
zAssign27Neighbors(DoubleMatrix3D B,
Double27Function function)
27 neighbor stencil transformation.
|
double |
zSum()
Returns the sum of all cells; Sum( x[i,j,k] ).
|
copy, equals, equals, get, like, normalize, set, toString, viewColumn, viewColumnFlip, viewDice, viewPart, viewRow, viewRowFlip, viewSelection, viewSelection, viewSlice, viewSliceFlip, viewSorted, viewStridescheckShape, checkShape, columns, columnStride, rows, rowStride, size, slices, sliceStride, toStringShortensureCapacity, isView, trimToSizeclonepublic DenseDoubleMatrix3D(double[][][] values)
The values are copied. So subsequent changes in values are not reflected in the matrix, and vice-versa.
values - The values to be filled into the new matrix.IllegalArgumentException - if
for any 1 <= slice < values.length: values[slice].length != values[slice-1].length
.IllegalArgumentException - if
for any 1 <= row < values[0].length: values[slice][row].length != values[slice][row-1].length
.public DenseDoubleMatrix3D(int slices,
int rows,
int columns)
slices - the number of slices the matrix shall have.rows - the number of rows the matrix shall have.columns - the number of columns the matrix shall have.IllegalArgumentException - if (double)slices*columns*rows > Integer.MAX_VALUE.IllegalArgumentException - if slices<0 || rows<0 || columns<0.public DenseDoubleMatrix3D(int slices,
int rows,
int columns,
double[] elements,
int sliceZero,
int rowZero,
int columnZero,
int sliceStride,
int rowStride,
int columnStride,
boolean isView)
slices - the number of slices the matrix shall have.rows - the number of rows the matrix shall have.columns - the number of columns the matrix shall have.elements - the cells.sliceZero - the position of the first element.rowZero - the position of the first element.columnZero - the position of the first element.sliceStride - the number of elements between two slices, i.e.
index(k+1,i,j)-index(k,i,j).rowStride - the number of elements between two rows, i.e.
index(k,i+1,j)-index(k,i,j).columnStride - the number of elements between two columns, i.e.
index(k,i,j+1)-index(k,i,j).isView - if true then a matrix view is constructedIllegalArgumentException - if (double)slices*columns*rows > Integer.MAX_VALUE.IllegalArgumentException - if slices<0 || rows<0 || columns<0.public double aggregate(DoubleDoubleFunction aggr, DoubleFunction f)
DoubleMatrix3DExample:
cern.jet.math.Functions F = cern.jet.math.Functions.functions;
2 x 2 x 2 matrix
0 1
2 3
4 5
6 7
// Sum( x[slice,row,col]*x[slice,row,col] )
matrix.aggregate(F.plus,F.square);
--> 140
For further examples, see the package doc.aggregate in class DoubleMatrix3Daggr - an aggregation function taking as first argument the current
aggregation and as second argument the transformed current
cell value.f - a function transforming the current cell value.DoubleFunctionspublic double aggregate(DoubleDoubleFunction aggr, DoubleFunction f, DoubleProcedure cond)
DoubleMatrix3Daggregate in class DoubleMatrix3Daggr - an aggregation function taking as first argument the current
aggregation and as second argument the transformed current
cell value.f - a function transforming the current cell value.cond - a condition.DoubleFunctionspublic double aggregate(DoubleDoubleFunction aggr, DoubleFunction f, IntArrayList sliceList, IntArrayList rowList, IntArrayList columnList)
DoubleMatrix3Daggregate in class DoubleMatrix3Daggr - an aggregation function taking as first argument the current
aggregation and as second argument the transformed current
cell value.f - a function transforming the current cell value.sliceList - slice indexes.rowList - row indexes.columnList - column indexes.DoubleFunctionspublic double aggregate(DoubleMatrix3D other, DoubleDoubleFunction aggr, DoubleDoubleFunction f)
DoubleMatrix3DExample:
cern.jet.math.Functions F = cern.jet.math.Functions.functions;
x = 2 x 2 x 2 matrix
0 1
2 3
4 5
6 7
y = 2 x 2 x 2 matrix
0 1
2 3
4 5
6 7
// Sum( x[slice,row,col] * y[slice,row,col] )
x.aggregate(y, F.plus, F.mult);
--> 140
// Sum( (x[slice,row,col] + y[slice,row,col])ˆ2 )
x.aggregate(y, F.plus, F.chain(F.square,F.plus));
--> 560
For further examples, see the package doc.aggregate in class DoubleMatrix3Daggr - an aggregation function taking as first argument the current
aggregation and as second argument the transformed current
cell values.f - a function transforming the current cell values.DoubleFunctionspublic DoubleMatrix3D assign(DoubleFunction function)
DoubleMatrix3DExample:
matrix = 1 x 2 x 2 matrix
0.5 1.5
2.5 3.5
// change each cell to its sine
matrix.assign(cern.jet.math.Functions.sin);
-->
1 x 2 x 2 matrix
0.479426 0.997495
0.598472 -0.350783
For further examples, see the package doc.assign in class DoubleMatrix3Dfunction - a function object taking as argument the current cell's value.DoubleFunctionspublic DoubleMatrix3D assign(DoubleProcedure cond, DoubleFunction f)
DoubleMatrix3Dassign in class DoubleMatrix3Dcond - a condition.f - a function object.DoubleFunctionspublic DoubleMatrix3D assign(DoubleProcedure cond, double value)
DoubleMatrix3Dassign in class DoubleMatrix3Dcond - a condition.value - a value.public DoubleMatrix3D assign(double value)
DoubleMatrix3Dassign in class DoubleMatrix3Dvalue - the value to be filled into the cells.public DoubleMatrix3D assign(double[] values)
DoubleMatrix3DThe values are copied. So subsequent changes in values are not reflected in the matrix, and vice-versa.
assign in class DoubleMatrix3Dvalues - the values to be filled into the cells.public DoubleMatrix3D assign(double[][][] values)
DoubleMatrix3DThe values are copied. So subsequent changes in values are not reflected in the matrix, and vice-versa.
assign in class DoubleMatrix3Dvalues - the values to be filled into the cells.public DoubleMatrix3D assign(DoubleMatrix3D source)
DoubleMatrix3Dassign in class DoubleMatrix3Dsource - the source matrix to copy from (may be identical to the
receiver).public DoubleMatrix3D assign(DoubleMatrix3D y, DoubleDoubleFunction function)
DoubleMatrix3DExample:
// assign x[row,col] = x[row,col]<sup>y[row,col]</sup>
m1 = 1 x 2 x 2 matrix
0 1
2 3
m2 = 1 x 2 x 2 matrix
0 2
4 6
m1.assign(m2, cern.jet.math.Functions.pow);
-->
m1 == 1 x 2 x 2 matrix
1 1
16 729
For further examples, see the package doc.assign in class DoubleMatrix3Dy - the secondary matrix to operate on.function - a function object taking as first argument the current cell's
value of this, and as second argument the current
cell's value of y,DoubleFunctionspublic DoubleMatrix3D assign(DoubleMatrix3D y, DoubleDoubleFunction function, IntArrayList sliceList, IntArrayList rowList, IntArrayList columnList)
DoubleMatrix3Dassign in class DoubleMatrix3Dy - the secondary matrix to operate on.function - a function object taking as first argument the current cell's
value of this, and as second argument the current
cell's value of y, *sliceList - slice indexes.rowList - row indexes.columnList - column indexes.DoubleFunctionspublic int cardinality()
DoubleMatrix3Dcardinality in class DoubleMatrix3Dpublic void dct2Slices(boolean scale)
scale - if true then scaling is performedpublic void dct3(boolean scale)
scale - if true then scaling is performedpublic void dht2Slices()
public void dht3()
public void dst2Slices(boolean scale)
scale - if true then scaling is performedpublic void dst3(boolean scale)
scale - if true then scaling is performedpublic double[] elements()
DoubleMatrix3Delements in class DoubleMatrix3Dpublic void fft3()
this[k1][k2][2*k3] = Re[k1][k2][k3]
= Re[(n1-k1)%n1][(n2-k2)%n2][n3-k3],
this[k1][k2][2*k3+1] = Im[k1][k2][k3]
= -Im[(n1-k1)%n1][(n2-k2)%n2][n3-k3],
0<=k1<n1, 0<=k2<n2, 0<k3<n3/2,
this[k1][k2][0] = Re[k1][k2][0]
= Re[(n1-k1)%n1][n2-k2][0],
this[k1][k2][1] = Im[k1][k2][0]
= -Im[(n1-k1)%n1][n2-k2][0],
this[k1][n2-k2][1] = Re[(n1-k1)%n1][k2][n3/2]
= Re[k1][n2-k2][n3/2],
this[k1][n2-k2][0] = -Im[(n1-k1)%n1][k2][n3/2]
= Im[k1][n2-k2][n3/2],
0<=k1<n1, 0<k2<n2/2,
this[k1][0][0] = Re[k1][0][0]
= Re[n1-k1][0][0],
this[k1][0][1] = Im[k1][0][0]
= -Im[n1-k1][0][0],
this[k1][n2/2][0] = Re[k1][n2/2][0]
= Re[n1-k1][n2/2][0],
this[k1][n2/2][1] = Im[k1][n2/2][0]
= -Im[n1-k1][n2/2][0],
this[n1-k1][0][1] = Re[k1][0][n3/2]
= Re[n1-k1][0][n3/2],
this[n1-k1][0][0] = -Im[k1][0][n3/2]
= Im[n1-k1][0][n3/2],
this[n1-k1][n2/2][1] = Re[k1][n2/2][n3/2]
= Re[n1-k1][n2/2][n3/2],
this[n1-k1][n2/2][0] = -Im[k1][n2/2][n3/2]
= Im[n1-k1][n2/2][n3/2],
0<k1<n1/2,
this[0][0][0] = Re[0][0][0],
this[0][0][1] = Re[0][0][n3/2],
this[0][n2/2][0] = Re[0][n2/2][0],
this[0][n2/2][1] = Re[0][n2/2][n3/2],
this[n1/2][0][0] = Re[n1/2][0][0],
this[n1/2][0][1] = Re[n1/2][0][n3/2],
this[n1/2][n2/2][0] = Re[n1/2][n2/2][0],
this[n1/2][n2/2][1] = Re[n1/2][n2/2][n3/2]
This method computes only half of the elements of the real transform. The
other half satisfies the symmetry condition. If you want the full real
forward transform, use getFft3. To get back the original
data, use ifft3.IllegalArgumentException - if the slice size or the row size or the column size of this
matrix is not a power of 2 number.public DenseDComplexMatrix3D getFft2Slices()
public DenseDComplexMatrix3D getFft3()
public DenseDComplexMatrix3D getIfft2Slices(boolean scale)
scale - if true then scaling is performedpublic DenseDComplexMatrix3D getIfft3(boolean scale)
scale - if true then scaling is performedpublic double[] getMaxLocation()
DoubleMatrix3DgetMaxLocation in class DoubleMatrix3Dpublic double[] getMinLocation()
DoubleMatrix3DgetMinLocation in class DoubleMatrix3Dpublic void getNegativeValues(IntArrayList sliceList, IntArrayList rowList, IntArrayList columnList, DoubleArrayList valueList)
DoubleMatrix3DgetNegativeValues in class DoubleMatrix3DsliceList - the list to be filled with slice indexes, can have any size.rowList - the list to be filled with row indexes, can have any size.columnList - the list to be filled with column indexes, can have any size.valueList - the list to be filled with values, can have any size.public void getNonZeros(IntArrayList sliceList, IntArrayList rowList, IntArrayList columnList, DoubleArrayList valueList)
DoubleMatrix3D
In general, fill order is unspecified. This implementation fill
like:
for (slice = 0..slices-1) for (row = 0..rows-1) for (column = 0..colums-1) do ...
. However, subclasses are free to us any other order, even an order that
may change over time as cell values are changed. (Of course, result lists
indexes are guaranteed to correspond to the same cell). For an example,
see
DoubleMatrix3D.getNonZeros(IntArrayList,IntArrayList,IntArrayList,DoubleArrayList).
getNonZeros in class DoubleMatrix3DsliceList - the list to be filled with slice indexes, can have any size.rowList - the list to be filled with row indexes, can have any size.columnList - the list to be filled with column indexes, can have any size.valueList - the list to be filled with values, can have any size.public void getPositiveValues(IntArrayList sliceList, IntArrayList rowList, IntArrayList columnList, DoubleArrayList valueList)
DoubleMatrix3DgetPositiveValues in class DoubleMatrix3DsliceList - the list to be filled with slice indexes, can have any size.rowList - the list to be filled with row indexes, can have any size.columnList - the list to be filled with column indexes, can have any size.valueList - the list to be filled with values, can have any size.public double getQuick(int slice,
int row,
int column)
DoubleMatrix3DProvided with invalid parameters this method may return invalid objects without throwing any exception. You should only use this method when you are absolutely sure that the coordinate is within bounds. Precondition (unchecked): slice<0 || slice>=slices() || row<0 || row>=rows() || column<0 || column>=column().
getQuick in class DoubleMatrix3Dslice - the index of the slice-coordinate.row - the index of the row-coordinate.column - the index of the column-coordinate.public void idct2Slices(boolean scale)
scale - if true then scaling is performedpublic void idct3(boolean scale)
scale - if true then scaling is performedpublic void idht2Slices(boolean scale)
scale - if true then scaling is performedpublic void idht3(boolean scale)
scale - if true then scaling is performedpublic void idst2Slices(boolean scale)
scale - if true then scaling is performedpublic void idst3(boolean scale)
scale - if true then scaling is performedpublic void ifft3(boolean scale)
this[k1][k2][2*k3] = Re[k1][k2][k3]
= Re[(n1-k1)%n1][(n2-k2)%n2][n3-k3],
this[k1][k2][2*k3+1] = Im[k1][k2][k3]
= -Im[(n1-k1)%n1][(n2-k2)%n2][n3-k3],
0<=k1<n1, 0<=k2<n2, 0<k3<n3/2,
this[k1][k2][0] = Re[k1][k2][0]
= Re[(n1-k1)%n1][n2-k2][0],
this[k1][k2][1] = Im[k1][k2][0]
= -Im[(n1-k1)%n1][n2-k2][0],
this[k1][n2-k2][1] = Re[(n1-k1)%n1][k2][n3/2]
= Re[k1][n2-k2][n3/2],
this[k1][n2-k2][0] = -Im[(n1-k1)%n1][k2][n3/2]
= Im[k1][n2-k2][n3/2],
0<=k1<n1, 0<k2<n2/2,
this[k1][0][0] = Re[k1][0][0]
= Re[n1-k1][0][0],
this[k1][0][1] = Im[k1][0][0]
= -Im[n1-k1][0][0],
this[k1][n2/2][0] = Re[k1][n2/2][0]
= Re[n1-k1][n2/2][0],
this[k1][n2/2][1] = Im[k1][n2/2][0]
= -Im[n1-k1][n2/2][0],
this[n1-k1][0][1] = Re[k1][0][n3/2]
= Re[n1-k1][0][n3/2],
this[n1-k1][0][0] = -Im[k1][0][n3/2]
= Im[n1-k1][0][n3/2],
this[n1-k1][n2/2][1] = Re[k1][n2/2][n3/2]
= Re[n1-k1][n2/2][n3/2],
this[n1-k1][n2/2][0] = -Im[k1][n2/2][n3/2]
= Im[n1-k1][n2/2][n3/2],
0<k1<n1/2,
this[0][0][0] = Re[0][0][0],
this[0][0][1] = Re[0][0][n3/2],
this[0][n2/2][0] = Re[0][n2/2][0],
this[0][n2/2][1] = Re[0][n2/2][n3/2],
this[n1/2][0][0] = Re[n1/2][0][0],
this[n1/2][0][1] = Re[n1/2][0][n3/2],
this[n1/2][n2/2][0] = Re[n1/2][n2/2][0],
this[n1/2][n2/2][1] = Re[n1/2][n2/2][n3/2]
This method computes only half of the elements of the real transform. The
other half satisfies the symmetry condition. If you want the full real
inverse transform, use getIfft3.scale - if true then scaling is performedIllegalArgumentException - if the slice size or the row size or the column size of this
matrix is not a power of 2 number.public long index(int slice,
int row,
int column)
AbstractMatrix3Dindex in class AbstractMatrix3Dslice - the index of the slice-coordinate.row - the index of the row-coordinate.column - the index of the third-coordinate.public DoubleMatrix3D like(int slices, int rows, int columns)
DoubleMatrix3Dlike in class DoubleMatrix3Dslices - the number of slices the matrix shall have.rows - the number of rows the matrix shall have.columns - the number of columns the matrix shall have.public DoubleMatrix2D like2D(int rows, int columns)
DoubleMatrix3Dlike2D in class DoubleMatrix3Drows - the number of rows the matrix shall have.columns - the number of columns the matrix shall have.public void setQuick(int slice,
int row,
int column,
double value)
DoubleMatrix3DProvided with invalid parameters this method may access illegal indexes without throwing any exception. You should only use this method when you are absolutely sure that the coordinate is within bounds. Precondition (unchecked): slice<0 || slice>=slices() || row<0 || row>=rows() || column<0 || column>=column().
setQuick in class DoubleMatrix3Dslice - the index of the slice-coordinate.row - the index of the row-coordinate.column - the index of the column-coordinate.value - the value to be filled into the specified cell.public double[][][] toArray()
DoubleMatrix3DThe values are copied. So subsequent changes in values are not reflected in the matrix, and vice-versa.
toArray in class DoubleMatrix3Dpublic DoubleMatrix1D vectorize()
DoubleMatrix3Dvectorize in class DoubleMatrix3Dpublic void zAssign27Neighbors(DoubleMatrix3D B, Double27Function function)
DoubleMatrix3D
B[k,i,j] = function.apply(
A[k-1,i-1,j-1], A[k-1,i-1,j], A[k-1,i-1,j+1],
A[k-1,i, j-1], A[k-1,i, j], A[k-1,i, j+1],
A[k-1,i+1,j-1], A[k-1,i+1,j], A[k-1,i+1,j+1],
A[k ,i-1,j-1], A[k ,i-1,j], A[k ,i-1,j+1],
A[k ,i, j-1], A[k ,i, j], A[k ,i, j+1],
A[k ,i+1,j-1], A[k ,i+1,j], A[k ,i+1,j+1],
A[k+1,i-1,j-1], A[k+1,i-1,j], A[k+1,i-1,j+1],
A[k+1,i, j-1], A[k+1,i, j], A[k+1,i, j+1],
A[k+1,i+1,j-1], A[k+1,i+1,j], A[k+1,i+1,j+1]
)
x x x - - x x x - - - -
x o x - - x o x - - - -
x x x - - x x x ... - x x x
- - - - - - - - - x o x
- - - - - - - - - x x x
Make sure that cells of this and B do not overlap. In
case of overlapping views, behaviour is unspecified.
Example:
final double alpha = 0.25; final double beta = 0.75;
cern.colt.function.Double27Function f = new
cern.colt.function.Double27Function() { public final
double apply( double a000, double
a001, double a002, double a010,
double a011, double a012, double
a020, double a021, double a022,
double a100, double a101, double
a102, double a110, double a111,
double a112, double a120, double
a121, double a122,
double a200, double a201, double
a202, double a210, double a211,
double a212, double a220, double
a221, double a222) {
return beta*a111 +
alpha*(a000 + ... + a222); } };
A.zAssign27Neighbors(B,f);
zAssign27Neighbors in class DoubleMatrix3DB - the matrix to hold the results.function - the function to be applied to the 27 cells.public double zSum()
DoubleMatrix3DzSum in class DoubleMatrix3DJump to the Parallel Colt Homepage