Isis Developer Reference
Isis::SparseBlockMatrix Class Reference

SparseBlockMatrix. More...

#include <SparseBlockMatrix.h>

Inheritance diagram for Isis::SparseBlockMatrix:
Inheritance graph
Collaboration diagram for Isis::SparseBlockMatrix:
Collaboration graph

Public Member Functions

 SparseBlockMatrix ()
 
 ~SparseBlockMatrix ()
 Destructor.
 
 SparseBlockMatrix (const SparseBlockMatrix &src)
 Copy constructor.
 
SparseBlockMatrixoperator= (const SparseBlockMatrix &src)
 "Equals" operator.
 
void wipe ()
 Deletes all pointer elements and removes them from the list.
 
void copy (const SparseBlockMatrix &src)
 Copy method.
 
bool setNumberOfColumns (int n)
 Initializes number of columns (SparseBlockColumnMatrix).
 
void zeroBlocks ()
 Sets all elements of all matrix blocks to zero.
 
bool insertMatrixBlock (int nColumnBlock, int nRowBlock, int nRows, int nCols)
 Inserts a "newed" boost LinearAlgebra::Matrix pointer of size (nRows, nCols) into the matrix at nColumnBlock, nRowBlock.
 
LinearAlgebra::MatrixgetBlock (int column, int row)
 Returns pointer to boost matrix at position (column, row).
 
int numberOfBlocks ()
 Returns total number of blocks in matrix.
 
int numberOfDiagonalBlocks ()
 Returns number of diagonal matrix blocks (equivalent to size - there is one per column).
 
int numberOfOffDiagonalBlocks ()
 Returns number of off-diagonal matrix blocks.
 
int numberOfElements ()
 Returns number of matrix elements in matrix.
 
void print (std::ostream &outstream)
 Prints matrix blocks to std output stream out for debugging.
 
void printClean (std::ostream &outstream)
 Prints matrix blocks to std output stream out for debugging.
 
bool write (std::ofstream &fp_out, bool binary=true)
 
int getLeadingColumnsForBlock (int nblockColumn)
 Sums and returns the number of columns in each matrix block prior to nblockColumn.
 
int getLeadingRowsForBlock (int nblockRow)
 Sums and returns the number of rows in each matrix block prior to nblockRow.
 

Detailed Description

SparseBlockMatrix.

The CHOLMOD (Cholesky decomposition) package uses the compressed column storage (CCS) matrix format which is efficient in memory storage but inefficient to fill in an arbitrary manner because the insertion of every non-zero entry requires that all succeeding entries be shifted.

To build the reduced normal equations matrix, an interim sparse matrix structure is required that can be efficiently populated in a random fashion and can be traversed by column in row order to subsequently construct the CCS matrix required by CHOLMOD. We use a type of Block Compressed Column Storage (BCCS) which consists of an array of map containers (QList of SparseBlockColumnMatrices), each representing a column of square matrix blocks in the reduced normal equations. The key into each column map is the block’s row index.The value at each key is a square dense matrix (Boost matrix) with a dimension equivalent to the number of exterior orientation parameters used for the image. Zero blocks are not stored. The BCCS matrix is created only in the first iteration of the bundle adjustment; in subsequent iterations it need only be repopulated. As the normal equations matrix is symmetric, only the upper triangular portion is stored in memory.

Author
2011-07-29 Ken Edmundson

Constructor & Destructor Documentation

◆ SparseBlockMatrix() [1/2]

Isis::SparseBlockMatrix::SparseBlockMatrix ( )
inline

◆ ~SparseBlockMatrix()

Isis::SparseBlockMatrix::~SparseBlockMatrix ( )

Destructor.

References wipe().

◆ SparseBlockMatrix() [2/2]

Isis::SparseBlockMatrix::SparseBlockMatrix ( const SparseBlockMatrix & src)

Copy constructor.

Calls copy method immediately below.

Parameters
srcSparseBlockMatrix to copy

References copy().

Member Function Documentation

◆ copy()

void Isis::SparseBlockMatrix::copy ( const SparseBlockMatrix & src)

Copy method.

Parameters
srcSparseBlockMatrix to copy

References wipe().

Referenced by operator=(), and SparseBlockMatrix().

◆ getBlock()

LinearAlgebra::Matrix * Isis::SparseBlockMatrix::getBlock ( int column,
int row )

Returns pointer to boost matrix at position (column, row).

Parameters
columnblock column number
rowblock row number
Returns
LinearAlgebra::Matrix Pointer to Boost matrix at position (column, row)

◆ getLeadingColumnsForBlock()

int Isis::SparseBlockMatrix::getLeadingColumnsForBlock ( int nblockColumn)

Sums and returns the number of columns in each matrix block prior to nblockColumn.

Parameters
nblockColumn
Returns
int Number of leading column elements for block at nblockColumn

◆ getLeadingRowsForBlock()

int Isis::SparseBlockMatrix::getLeadingRowsForBlock ( int nblockRow)

Sums and returns the number of rows in each matrix block prior to nblockRow.

Parameters
nblockRow
Returns
int Number of leading row elements for block at nblockRow

◆ insertMatrixBlock()

bool Isis::SparseBlockMatrix::insertMatrixBlock ( int nColumnBlock,
int nRowBlock,
int nRows,
int nCols )

Inserts a "newed" boost LinearAlgebra::Matrix pointer of size (nRows, nCols) into the matrix at nColumnBlock, nRowBlock.

The inserted matrix elements are initialized to zero. If an entry exists at nColumnBlock, RowBlock, no insertion is made.

Parameters
nColumnBlockblock column number of inserted matrix (QList index)
nRowBlockblock row number of inserted matrix (key into map)
nRowsnumber of rows in matrix to be inserted
nColsnumber of columns in matrix to be inserted
Returns
bool Returns result of SparseBlockColumnMatrix::insertMatrixBlock()

◆ numberOfBlocks()

int Isis::SparseBlockMatrix::numberOfBlocks ( )

Returns total number of blocks in matrix.

Returns
int Total number of blocks in matrix

Referenced by numberOfOffDiagonalBlocks().

◆ numberOfDiagonalBlocks()

int Isis::SparseBlockMatrix::numberOfDiagonalBlocks ( )

Returns number of diagonal matrix blocks (equivalent to size - there is one per column).

Returns
int Number of diagnonal blocks in matrix

Referenced by numberOfOffDiagonalBlocks().

◆ numberOfElements()

int Isis::SparseBlockMatrix::numberOfElements ( )

Returns number of matrix elements in matrix.

Returns
int Total number of matrix elements

◆ numberOfOffDiagonalBlocks()

int Isis::SparseBlockMatrix::numberOfOffDiagonalBlocks ( )

Returns number of off-diagonal matrix blocks.

Returns
int Number of off-diagonal blocks in matrix

References numberOfBlocks(), and numberOfDiagonalBlocks().

◆ operator=()

SparseBlockMatrix & Isis::SparseBlockMatrix::operator= ( const SparseBlockMatrix & src)

"Equals" operator.

Parameters
srcSparseBlockMatrix to check against

References copy().

◆ print()

void Isis::SparseBlockMatrix::print ( std::ostream & outstream)

Prints matrix blocks to std output stream out for debugging.

Parameters
outstreamoutput stream

References Isis::SparseBlockColumnMatrix::print().

◆ printClean()

void Isis::SparseBlockMatrix::printClean ( std::ostream & outstream)

Prints matrix blocks to std output stream out for debugging.

Parameters
outstreamoutput stream

References Isis::SparseBlockColumnMatrix::printClean().

◆ setNumberOfColumns()

bool Isis::SparseBlockMatrix::setNumberOfColumns ( int n)

Initializes number of columns (SparseBlockColumnMatrix).

Parameters
nnumber of columns to insert
Returns
bool Always returns true TODO: why bother returning bool? should there be a false condition?

◆ wipe()

void Isis::SparseBlockMatrix::wipe ( )

Deletes all pointer elements and removes them from the list.

Effectively, a destructor, and in fact, called by the ~SparseBlockMatrix above.

Referenced by copy(), and ~SparseBlockMatrix().

◆ write()

bool Isis::SparseBlockMatrix::write ( std::ofstream & fp_out,
bool binary = true )

◆ zeroBlocks()

void Isis::SparseBlockMatrix::zeroBlocks ( )

Sets all elements of all matrix blocks to zero.

Referenced by Isis::BundleAdjust::solveCholesky().


The documentation for this class was generated from the following files: