Matrix< Real > Class Template Reference

A class for storing matrices. More...

#include <kaldi-matrix.h>

Inheritance diagram for Matrix< Real >:
Collaboration diagram for Matrix< Real >:

Public Member Functions

 Matrix ()
 Empty constructor. More...
 
 Matrix (const MatrixIndexT r, const MatrixIndexT c, MatrixResizeType resize_type=kSetZero, MatrixStrideType stride_type=kDefaultStride)
 Basic constructor. More...
 
template<typename OtherReal >
 Matrix (const CuMatrixBase< OtherReal > &cu, MatrixTransposeType trans=kNoTrans)
 Copy constructor from CUDA matrix This is defined in ../cudamatrix/cu-matrix.h. More...
 
void Swap (Matrix< Real > *other)
 Swaps the contents of *this and *other. Shallow swap. More...
 
void Swap (CuMatrix< Real > *mat)
 Defined in ../cudamatrix/cu-matrix.cc. More...
 
 Matrix (const MatrixBase< Real > &M, MatrixTransposeType trans=kNoTrans)
 Constructor from any MatrixBase. More...
 
 Matrix (const Matrix< Real > &M)
 Same as above, but need to avoid default copy constructor. More...
 
template<typename OtherReal >
 Matrix (const MatrixBase< OtherReal > &M, MatrixTransposeType trans=kNoTrans)
 Copy constructor: as above, but from another type. More...
 
template<typename OtherReal >
 Matrix (const SpMatrix< OtherReal > &M)
 Copy constructor taking SpMatrix... More...
 
 Matrix (const CompressedMatrix &C)
 Constructor from CompressedMatrix. More...
 
template<typename OtherReal >
 Matrix (const TpMatrix< OtherReal > &M, MatrixTransposeType trans=kNoTrans)
 Copy constructor taking TpMatrix... More...
 
void Read (std::istream &in, bool binary, bool add=false)
 read from stream. More...
 
void RemoveRow (MatrixIndexT i)
 Remove a specified row. More...
 
void Transpose ()
 Transpose the matrix. More...
 
 ~Matrix ()
 Distructor to free matrices. More...
 
void Resize (const MatrixIndexT r, const MatrixIndexT c, MatrixResizeType resize_type=kSetZero, MatrixStrideType stride_type=kDefaultStride)
 Sets matrix to a specified size (zero is OK as long as both r and c are zero). More...
 
Matrix< Real > & operator= (const MatrixBase< Real > &other)
 Assignment operator that takes MatrixBase. More...
 
Matrix< Real > & operator= (const Matrix< Real > &other)
 Assignment operator. Needed for inclusion in std::vector. More...
 
- Public Member Functions inherited from MatrixBase< Real >
MatrixIndexT NumRows () const
 Returns number of rows (or zero for empty matrix). More...
 
MatrixIndexT NumCols () const
 Returns number of columns (or zero for empty matrix). More...
 
MatrixIndexT Stride () const
 Stride (distance in memory between each row). Will be >= NumCols. More...
 
size_t SizeInBytes () const
 Returns size in bytes of the data held by the matrix. More...
 
const Real * Data () const
 Gives pointer to raw data (const). More...
 
Real * Data ()
 Gives pointer to raw data (non-const). More...
 
Real * RowData (MatrixIndexT i)
 Returns pointer to data for one row (non-const) More...
 
const Real * RowData (MatrixIndexT i) const
 Returns pointer to data for one row (const) More...
 
Real & operator() (MatrixIndexT r, MatrixIndexT c)
 Indexing operator, non-const (only checks sizes if compiled with -DKALDI_PARANOID) More...
 
Real & Index (MatrixIndexT r, MatrixIndexT c)
 Indexing operator, provided for ease of debugging (gdb doesn't work with parenthesis operator). More...
 
const Real operator() (MatrixIndexT r, MatrixIndexT c) const
 Indexing operator, const (only checks sizes if compiled with -DKALDI_PARANOID) More...
 
void SetZero ()
 Sets matrix to zero. More...
 
void Set (Real)
 Sets all elements to a specific value. More...
 
void SetUnit ()
 Sets to zero, except ones along diagonal [for non-square matrices too]. More...
 
void SetRandn ()
 Sets to random values of a normal distribution. More...
 
void SetRandUniform ()
 Sets to numbers uniformly distributed on (0, 1) More...
 
template<typename OtherReal >
void CopyFromMat (const MatrixBase< OtherReal > &M, MatrixTransposeType trans=kNoTrans)
 Copy given matrix. (no resize is done). More...
 
void CopyFromMat (const CompressedMatrix &M)
 Copy from compressed matrix. More...
 
template<typename OtherReal >
void CopyFromSp (const SpMatrix< OtherReal > &M)
 Copy given spmatrix. (no resize is done). More...
 
template<typename OtherReal >
void CopyFromTp (const TpMatrix< OtherReal > &M, MatrixTransposeType trans=kNoTrans)
 Copy given tpmatrix. (no resize is done). More...
 
template<typename OtherReal >
void CopyFromMat (const CuMatrixBase< OtherReal > &M, MatrixTransposeType trans=kNoTrans)
 Copy from CUDA matrix. Implemented in ../cudamatrix/cu-matrix.h. More...
 
void CopyRowsFromVec (const VectorBase< Real > &v)
 This function has two modes of operation. More...
 
void CopyRowsFromVec (const CuVectorBase< Real > &v)
 This version of CopyRowsFromVec is implemented in ../cudamatrix/cu-vector.cc. More...
 
template<typename OtherReal >
void CopyRowsFromVec (const VectorBase< OtherReal > &v)
 
void CopyColsFromVec (const VectorBase< Real > &v)
 Copies vector into matrix, column-by-column. More...
 
void CopyColFromVec (const VectorBase< Real > &v, const MatrixIndexT col)
 Copy vector into specific column of matrix. More...
 
void CopyRowFromVec (const VectorBase< Real > &v, const MatrixIndexT row)
 Copy vector into specific row of matrix. More...
 
void CopyDiagFromVec (const VectorBase< Real > &v)
 Copy vector into diagonal of matrix. More...
 
const SubVector< Real > Row (MatrixIndexT i) const
 Return specific row of matrix [const]. More...
 
SubVector< Real > Row (MatrixIndexT i)
 Return specific row of matrix. More...
 
SubMatrix< Real > Range (const MatrixIndexT row_offset, const MatrixIndexT num_rows, const MatrixIndexT col_offset, const MatrixIndexT num_cols) const
 Return a sub-part of matrix. More...
 
SubMatrix< Real > RowRange (const MatrixIndexT row_offset, const MatrixIndexT num_rows) const
 
SubMatrix< Real > ColRange (const MatrixIndexT col_offset, const MatrixIndexT num_cols) const
 
Real Sum () const
 Returns sum of all elements in matrix. More...
 
Real Trace (bool check_square=true) const
 Returns trace of matrix. More...
 
Real Max () const
 Returns maximum element of matrix. More...
 
Real Min () const
 Returns minimum element of matrix. More...
 
void MulElements (const MatrixBase< Real > &A)
 Element by element multiplication with a given matrix. More...
 
void DivElements (const MatrixBase< Real > &A)
 Divide each element by the corresponding element of a given matrix. More...
 
void Scale (Real alpha)
 Multiply each element with a scalar value. More...
 
void Max (const MatrixBase< Real > &A)
 Set, element-by-element, *this = max(*this, A) More...
 
void Min (const MatrixBase< Real > &A)
 Set, element-by-element, *this = min(*this, A) More...
 
void MulColsVec (const VectorBase< Real > &scale)
 Equivalent to (*this) = (*this) * diag(scale). More...
 
void MulRowsVec (const VectorBase< Real > &scale)
 Equivalent to (*this) = diag(scale) * (*this). More...
 
void MulRowsGroupMat (const MatrixBase< Real > &src)
 Divide each row into src.NumCols() equal groups, and then scale i'th row's j'th group of elements by src(i, j). More...
 
Real LogDet (Real *det_sign=NULL) const
 Returns logdet of matrix. More...
 
void Invert (Real *log_det=NULL, Real *det_sign=NULL, bool inverse_needed=true)
 matrix inverse. More...
 
void InvertDouble (Real *LogDet=NULL, Real *det_sign=NULL, bool inverse_needed=true)
 matrix inverse [double]. More...
 
void InvertElements ()
 Inverts all the elements of the matrix. More...
 
void Transpose ()
 Transpose the matrix. More...
 
void CopyCols (const MatrixBase< Real > &src, const MatrixIndexT *indices)
 Copies column r from column indices[r] of src. More...
 
void CopyRows (const MatrixBase< Real > &src, const MatrixIndexT *indices)
 Copies row r from row indices[r] of src (does nothing As a special case, if indexes[i] == -1, sets row i to zero. More...
 
void AddCols (const MatrixBase< Real > &src, const MatrixIndexT *indices)
 Add column indices[r] of src to column r. More...
 
void CopyRows (const Real *const *src)
 Copies row r of this matrix from an array of floats at the location given by src[r]. More...
 
void CopyToRows (Real *const *dst) const
 Copies row r of this matrix to the array of floats at the location given by dst[r]. More...
 
void AddRows (Real alpha, const MatrixBase< Real > &src, const MatrixIndexT *indexes)
 Does for each row r, this.Row(r) += alpha * src.row(indexes[r]). More...
 
void AddRows (Real alpha, const Real *const *src)
 Does for each row r, this.Row(r) += alpha * src[r], treating src[r] as the beginning of a region of memory representing a vector of floats, of the same length as this.NumCols(). More...
 
void AddToRows (Real alpha, Real *const *dst) const
 For each row r of this matrix, adds it (times alpha) to the array of floats at the location given by dst[r]. More...
 
void AddToRows (Real alpha, const MatrixIndexT *indexes, MatrixBase< Real > *dst) const
 For each row i of *this, adds this->Row(i) to dst->Row(indexes(i)) if indexes(i) >= 0, else do nothing. More...
 
void ApplyPow (Real power)
 
void ApplyPowAbs (Real power, bool include_sign=false)
 
void ApplyHeaviside ()
 
void ApplyFloor (Real floor_val)
 
void ApplyCeiling (Real ceiling_val)
 
void ApplyExp ()
 
void ApplyExpSpecial ()
 
void ApplyExpLimited (Real lower_limit, Real upper_limit)
 
void ApplyLog ()
 
void Eig (MatrixBase< Real > *P, VectorBase< Real > *eigs_real, VectorBase< Real > *eigs_imag) const
 Eigenvalue Decomposition of a square NxN matrix into the form (*this) = P D P^{-1}. More...
 
bool Power (Real pow)
 The Power method attempts to take the matrix to a power using a method that works in general for fractional and negative powers. More...
 
void DestructiveSvd (VectorBase< Real > *s, MatrixBase< Real > *U, MatrixBase< Real > *Vt)
 Singular value decomposition Major limitations: For nonsquare matrices, we assume m>=n (NumRows >= NumCols), and we return the "skinny" Svd, i.e. More...
 
void Svd (VectorBase< Real > *s, MatrixBase< Real > *U, MatrixBase< Real > *Vt) const
 Compute SVD (*this) = U diag(s) Vt. More...
 
void Svd (VectorBase< Real > *s) const
 Compute SVD but only retain the singular values. More...
 
Real MinSingularValue () const
 Returns smallest singular value. More...
 
void TestUninitialized () const
 
Real Cond () const
 Returns condition number by computing Svd. More...
 
bool IsSymmetric (Real cutoff=1.0e-05) const
 Returns true if matrix is Symmetric. More...
 
bool IsDiagonal (Real cutoff=1.0e-05) const
 Returns true if matrix is Diagonal. More...
 
bool IsUnit (Real cutoff=1.0e-05) const
 Returns true if the matrix is all zeros, except for ones on diagonal. More...
 
bool IsZero (Real cutoff=1.0e-05) const
 Returns true if matrix is all zeros. More...
 
Real FrobeniusNorm () const
 Frobenius norm, which is the sqrt of sum of square elements. More...
 
bool ApproxEqual (const MatrixBase< Real > &other, float tol=0.01) const
 Returns true if ((*this)-other).FrobeniusNorm() <= tol * (*this).FrobeniusNorm(). More...
 
bool Equal (const MatrixBase< Real > &other) const
 Tests for exact equality. It's usually preferable to use ApproxEqual. More...
 
Real LargestAbsElem () const
 largest absolute value. More...
 
Real LogSumExp (Real prune=-1.0) const
 Returns log(sum(exp())) without exp overflow If prune > 0.0, it uses a pruning beam, discarding terms less than (max - prune). More...
 
Real ApplySoftMax ()
 Apply soft-max to the collection of all elements of the matrix and return normalizer (log sum of exponentials). More...
 
void Sigmoid (const MatrixBase< Real > &src)
 Set each element to the sigmoid of the corresponding element of "src". More...
 
void Heaviside (const MatrixBase< Real > &src)
 Sets each element to the Heaviside step function (x > 0 ? 1 : 0) of the corresponding element in "src". More...
 
void Exp (const MatrixBase< Real > &src)
 
void Pow (const MatrixBase< Real > &src, Real power)
 
void Log (const MatrixBase< Real > &src)
 
void PowAbs (const MatrixBase< Real > &src, Real power, bool include_sign=false)
 Apply power to the absolute value of each element. More...
 
void Floor (const MatrixBase< Real > &src, Real floor_val)
 
void Ceiling (const MatrixBase< Real > &src, Real ceiling_val)
 
void ExpSpecial (const MatrixBase< Real > &src)
 For each element x of the matrix, set it to (x < 0 ? exp(x) : x + 1). More...
 
void ExpLimited (const MatrixBase< Real > &src, Real lower_limit, Real upper_limit)
 This is equivalent to running: Floor(src, lower_limit); Ceiling(src, upper_limit); Exp(src) More...
 
void SoftHinge (const MatrixBase< Real > &src)
 Set each element to y = log(1 + exp(x)) More...
 
void GroupPnorm (const MatrixBase< Real > &src, Real power)
 Apply the function y(i) = (sum_{j = i*G}^{(i+1)*G-1} x_j^(power))^(1 / p). More...
 
void GroupPnormDeriv (const MatrixBase< Real > &input, const MatrixBase< Real > &output, Real power)
 Calculate derivatives for the GroupPnorm function above... More...
 
void GroupMax (const MatrixBase< Real > &src)
 Apply the function y(i) = (max_{j = i*G}^{(i+1)*G-1} x_j Requires src.NumRows() == this->NumRows() and src.NumCols() % this->NumCols() == 0. More...
 
void GroupMaxDeriv (const MatrixBase< Real > &input, const MatrixBase< Real > &output)
 Calculate derivatives for the GroupMax function above, where "input" is the input to the GroupMax function above (i.e. More...
 
void Tanh (const MatrixBase< Real > &src)
 Set each element to the tanh of the corresponding element of "src". More...
 
void DiffSigmoid (const MatrixBase< Real > &value, const MatrixBase< Real > &diff)
 
void DiffTanh (const MatrixBase< Real > &value, const MatrixBase< Real > &diff)
 
void SymPosSemiDefEig (VectorBase< Real > *s, MatrixBase< Real > *P, Real check_thresh=0.001)
 Uses Svd to compute the eigenvalue decomposition of a symmetric positive semi-definite matrix: (*this) = rP * diag(rS) * rP^T, with rP an orthogonal matrix so rP^{-1} = rP^T. More...
 
void Add (const Real alpha)
 Add a scalar to each element. More...
 
void AddToDiag (const Real alpha)
 Add a scalar to each diagonal element. More...
 
template<typename OtherReal >
void AddVecVec (const Real alpha, const VectorBase< OtherReal > &a, const VectorBase< OtherReal > &b)
 *this += alpha * a * b^T More...
 
template<typename OtherReal >
void AddVecToRows (const Real alpha, const VectorBase< OtherReal > &v)
 [each row of *this] += alpha * v More...
 
template<typename OtherReal >
void AddVecToCols (const Real alpha, const VectorBase< OtherReal > &v)
 [each col of *this] += alpha * v More...
 
void AddMat (const Real alpha, const MatrixBase< Real > &M, MatrixTransposeType transA=kNoTrans)
 *this += alpha * M [or M^T] More...
 
void AddSmat (Real alpha, const SparseMatrix< Real > &A, MatrixTransposeType trans=kNoTrans)
 *this += alpha * A [or A^T]. More...
 
void AddSmatMat (Real alpha, const SparseMatrix< Real > &A, MatrixTransposeType transA, const MatrixBase< Real > &B, Real beta)
 (*this) = alpha * op(A) * B + beta * (*this), where A is sparse. More...
 
void AddMatSmat (Real alpha, const MatrixBase< Real > &A, const SparseMatrix< Real > &B, MatrixTransposeType transB, Real beta)
 (*this) = alpha * A * op(B) + beta * (*this), where B is sparse and op(B) is either B or trans(B) depending on the 'transB' argument. More...
 
void SymAddMat2 (const Real alpha, const MatrixBase< Real > &M, MatrixTransposeType transA, Real beta)
 *this = beta * *this + alpha * M M^T, for symmetric matrices. More...
 
void AddDiagVecMat (const Real alpha, const VectorBase< Real > &v, const MatrixBase< Real > &M, MatrixTransposeType transM, Real beta=1.0)
 *this = beta * *this + alpha * diag(v) * M [or M^T]. More...
 
void AddMatDiagVec (const Real alpha, const MatrixBase< Real > &M, MatrixTransposeType transM, VectorBase< Real > &v, Real beta=1.0)
 *this = beta * *this + alpha * M [or M^T] * diag(v) The same as adding M but scaling each column M_j by v(j). More...
 
void AddMatMatElements (const Real alpha, const MatrixBase< Real > &A, const MatrixBase< Real > &B, const Real beta)
 *this = beta * *this + alpha * A .* B (.* element by element multiplication) More...
 
template<typename OtherReal >
void AddSp (const Real alpha, const SpMatrix< OtherReal > &S)
 *this += alpha * S More...
 
void AddMatMat (const Real alpha, const MatrixBase< Real > &A, MatrixTransposeType transA, const MatrixBase< Real > &B, MatrixTransposeType transB, const Real beta)
 
void SetMatMatDivMat (const MatrixBase< Real > &A, const MatrixBase< Real > &B, const MatrixBase< Real > &C)
 *this = a * b / c (by element; when c = 0, *this = a) More...
 
void AddMatSmat (const Real alpha, const MatrixBase< Real > &A, MatrixTransposeType transA, const MatrixBase< Real > &B, MatrixTransposeType transB, const Real beta)
 A version of AddMatMat specialized for when the second argument contains a lot of zeroes. More...
 
void AddSmatMat (const Real alpha, const MatrixBase< Real > &A, MatrixTransposeType transA, const MatrixBase< Real > &B, MatrixTransposeType transB, const Real beta)
 A version of AddMatMat specialized for when the first argument contains a lot of zeroes. More...
 
void AddMatMatMat (const Real alpha, const MatrixBase< Real > &A, MatrixTransposeType transA, const MatrixBase< Real > &B, MatrixTransposeType transB, const MatrixBase< Real > &C, MatrixTransposeType transC, const Real beta)
 this <– beta*this + alpha*A*B*C. More...
 
void AddSpMat (const Real alpha, const SpMatrix< Real > &A, const MatrixBase< Real > &B, MatrixTransposeType transB, const Real beta)
 this <– beta*this + alpha*SpA*B. More...
 
void AddTpMat (const Real alpha, const TpMatrix< Real > &A, MatrixTransposeType transA, const MatrixBase< Real > &B, MatrixTransposeType transB, const Real beta)
 this <– beta*this + alpha*A*B. More...
 
void AddMatSp (const Real alpha, const MatrixBase< Real > &A, MatrixTransposeType transA, const SpMatrix< Real > &B, const Real beta)
 this <– beta*this + alpha*A*B. More...
 
void AddSpMatSp (const Real alpha, const SpMatrix< Real > &A, const MatrixBase< Real > &B, MatrixTransposeType transB, const SpMatrix< Real > &C, const Real beta)
 this <– beta*this + alpha*A*B*C. More...
 
void AddMatTp (const Real alpha, const MatrixBase< Real > &A, MatrixTransposeType transA, const TpMatrix< Real > &B, MatrixTransposeType transB, const Real beta)
 this <– beta*this + alpha*A*B. More...
 
void AddTpTp (const Real alpha, const TpMatrix< Real > &A, MatrixTransposeType transA, const TpMatrix< Real > &B, MatrixTransposeType transB, const Real beta)
 this <– beta*this + alpha*A*B. More...
 
void AddSpSp (const Real alpha, const SpMatrix< Real > &A, const SpMatrix< Real > &B, const Real beta)
 this <– beta*this + alpha*A*B. More...
 
void CopyLowerToUpper ()
 Copy lower triangle to upper triangle (symmetrize) More...
 
void CopyUpperToLower ()
 Copy upper triangle to lower triangle (symmetrize) More...
 
void OrthogonalizeRows ()
 This function orthogonalizes the rows of a matrix using the Gram-Schmidt process. More...
 
void Read (std::istream &in, bool binary, bool add=false)
 stream read. More...
 
void Write (std::ostream &out, bool binary) const
 write to stream. More...
 
void LapackGesvd (VectorBase< Real > *s, MatrixBase< Real > *U, MatrixBase< Real > *Vt)
 
template<>
void AddVecVec (const float alpha, const VectorBase< float > &ra, const VectorBase< float > &rb)
 
template<>
void AddVecVec (const double alpha, const VectorBase< double > &ra, const VectorBase< double > &rb)
 
template<>
void AddVecVec (const float alpha, const VectorBase< float > &a, const VectorBase< float > &rb)
 
template<>
void AddVecVec (const double alpha, const VectorBase< double > &a, const VectorBase< double > &rb)
 
template<>
void CopyFromSp (const SpMatrix< float > &M)
 
template<>
void CopyFromSp (const SpMatrix< double > &M)
 

Private Member Functions

void Destroy ()
 Deallocates memory and sets to empty matrix (dimension 0, 0). More...
 
void Init (const MatrixIndexT r, const MatrixIndexT c, const MatrixStrideType stride_type)
 Init assumes the current class contents are invalid (i.e. More...
 

Additional Inherited Members

- Protected Member Functions inherited from MatrixBase< Real >
 MatrixBase (Real *data, MatrixIndexT cols, MatrixIndexT rows, MatrixIndexT stride)
 Initializer, callable only from child. More...
 
 MatrixBase ()
 Initializer, callable only from child. More...
 
 ~MatrixBase ()
 
Real * Data_workaround () const
 A workaround that allows SubMatrix to get a pointer to non-const data for const Matrix. More...
 
- Protected Attributes inherited from MatrixBase< Real >
Real * data_
 data memory area More...
 
MatrixIndexT num_cols_
 these attributes store the real matrix size as it is stored in memory including memalignment More...
 
MatrixIndexT num_rows_
 < Number of columns More...
 
MatrixIndexT stride_
 < Number of rows More...
 

Detailed Description

template<typename Real>
class kaldi::Matrix< Real >

A class for storing matrices.

Definition at line 823 of file kaldi-matrix.h.

Constructor & Destructor Documentation

◆ Matrix() [1/9]

Matrix ( )

Empty constructor.

Definition at line 29 of file kaldi-matrix-inl.h.

29 : MatrixBase<Real>(NULL, 0, 0, 0) { }

◆ Matrix() [2/9]

Matrix ( const MatrixIndexT  r,
const MatrixIndexT  c,
MatrixResizeType  resize_type = kSetZero,
MatrixStrideType  stride_type = kDefaultStride 
)
inline

Basic constructor.

Definition at line 830 of file kaldi-matrix.h.

832  :
833  MatrixBase<Real>() { Resize(r, c, resize_type, stride_type); }
void Resize(const MatrixIndexT r, const MatrixIndexT c, MatrixResizeType resize_type=kSetZero, MatrixStrideType stride_type=kDefaultStride)
Sets matrix to a specified size (zero is OK as long as both r and c are zero).

◆ Matrix() [3/9]

Matrix ( const CuMatrixBase< OtherReal > &  cu,
MatrixTransposeType  trans = kNoTrans 
)
explicit

Copy constructor from CUDA matrix This is defined in ../cudamatrix/cu-matrix.h.

Definition at line 966 of file cu-matrix.h.

967  {
968  if (trans == kNoTrans) Init(M.NumRows(), M.NumCols(), kDefaultStride);
969  else Init(M.NumCols(), M.NumRows(), kDefaultStride);
970  M.CopyToMat(this, trans);
971 }
void Init(const MatrixIndexT r, const MatrixIndexT c, const MatrixStrideType stride_type)
Init assumes the current class contents are invalid (i.e.

◆ Matrix() [4/9]

Matrix ( const MatrixBase< Real > &  M,
MatrixTransposeType  trans = kNoTrans 
)
explicit

Constructor from any MatrixBase.

Can also copy with transpose. Allocates new memory.

Definition at line 740 of file kaldi-matrix.cc.

742  : MatrixBase<Real>() {
743  if (trans == kNoTrans) {
744  Resize(M.num_rows_, M.num_cols_);
745  this->CopyFromMat(M);
746  } else {
747  Resize(M.num_cols_, M.num_rows_);
748  this->CopyFromMat(M, kTrans);
749  }
750 }
void CopyFromMat(const MatrixBase< OtherReal > &M, MatrixTransposeType trans=kNoTrans)
Copy given matrix. (no resize is done).
void Resize(const MatrixIndexT r, const MatrixIndexT c, MatrixResizeType resize_type=kSetZero, MatrixStrideType stride_type=kDefaultStride)
Sets matrix to a specified size (zero is OK as long as both r and c are zero).

◆ Matrix() [5/9]

Matrix ( const Matrix< Real > &  M)

Same as above, but need to avoid default copy constructor.

Definition at line 754 of file kaldi-matrix.cc.

754  :
755  MatrixBase<Real>() {
756  Resize(M.num_rows_, M.num_cols_);
757  this->CopyFromMat(M);
758 }
void CopyFromMat(const MatrixBase< OtherReal > &M, MatrixTransposeType trans=kNoTrans)
Copy given matrix. (no resize is done).
void Resize(const MatrixIndexT r, const MatrixIndexT c, MatrixResizeType resize_type=kSetZero, MatrixStrideType stride_type=kDefaultStride)
Sets matrix to a specified size (zero is OK as long as both r and c are zero).

◆ Matrix() [6/9]

Matrix ( const MatrixBase< OtherReal > &  M,
MatrixTransposeType  trans = kNoTrans 
)
explicit

Copy constructor: as above, but from another type.

Copy constructor from another type.

Definition at line 763 of file kaldi-matrix.cc.

764  : MatrixBase<Real>() {
765  if (trans == kNoTrans) {
766  Resize(M.NumRows(), M.NumCols());
767  this->CopyFromMat(M);
768  } else {
769  Resize(M.NumCols(), M.NumRows());
770  this->CopyFromMat(M, kTrans);
771  }
772 }
void CopyFromMat(const MatrixBase< OtherReal > &M, MatrixTransposeType trans=kNoTrans)
Copy given matrix. (no resize is done).
void Resize(const MatrixIndexT r, const MatrixIndexT c, MatrixResizeType resize_type=kSetZero, MatrixStrideType stride_type=kDefaultStride)
Sets matrix to a specified size (zero is OK as long as both r and c are zero).

◆ Matrix() [7/9]

Matrix ( const SpMatrix< OtherReal > &  M)
inlineexplicit

Copy constructor taking SpMatrix...

It is symmetric, so no option for transpose, and NumRows == Cols

Definition at line 864 of file kaldi-matrix.h.

864  : MatrixBase<Real>() {
865  Resize(M.NumRows(), M.NumRows(), kUndefined);
866  this->CopyFromSp(M);
867  }
void CopyFromSp(const SpMatrix< OtherReal > &M)
Copy given spmatrix. (no resize is done).
void Resize(const MatrixIndexT r, const MatrixIndexT c, MatrixResizeType resize_type=kSetZero, MatrixStrideType stride_type=kDefaultStride)
Sets matrix to a specified size (zero is OK as long as both r and c are zero).

◆ Matrix() [8/9]

Matrix ( const CompressedMatrix< Real > &  C)
explicit

Constructor from CompressedMatrix.

Definition at line 2062 of file kaldi-matrix.cc.

2062  : MatrixBase<Real>() {
2063  Resize(M.NumRows(), M.NumCols(), kUndefined);
2064  M.CopyToMat(this);
2065 }
void Resize(const MatrixIndexT r, const MatrixIndexT c, MatrixResizeType resize_type=kSetZero, MatrixStrideType stride_type=kDefaultStride)
Sets matrix to a specified size (zero is OK as long as both r and c are zero).

◆ Matrix() [9/9]

Matrix ( const TpMatrix< OtherReal > &  M,
MatrixTransposeType  trans = kNoTrans 
)
inlineexplicit

Copy constructor taking TpMatrix...

Definition at line 874 of file kaldi-matrix.h.

875  : MatrixBase<Real>() {
876  if (trans == kNoTrans) {
877  Resize(M.NumRows(), M.NumCols(), kUndefined);
878  this->CopyFromTp(M);
879  } else {
880  Resize(M.NumCols(), M.NumRows(), kUndefined);
881  this->CopyFromTp(M, kTrans);
882  }
883  }
void CopyFromTp(const TpMatrix< OtherReal > &M, MatrixTransposeType trans=kNoTrans)
Copy given tpmatrix. (no resize is done).
void Resize(const MatrixIndexT r, const MatrixIndexT c, MatrixResizeType resize_type=kSetZero, MatrixStrideType stride_type=kDefaultStride)
Sets matrix to a specified size (zero is OK as long as both r and c are zero).

◆ ~Matrix()

~Matrix ( )
inline

Distructor to free matrices.

Definition at line 897 of file kaldi-matrix.h.

897 { Destroy(); }
void Destroy()
Deallocates memory and sets to empty matrix (dimension 0, 0).

Member Function Documentation

◆ Destroy()

void Destroy ( )
private

Deallocates memory and sets to empty matrix (dimension 0, 0).

Definition at line 1128 of file kaldi-matrix.cc.

1128  {
1129  // we need to free the data block if it was defined
1130  if (NULL != MatrixBase<Real>::data_)
1132  MatrixBase<Real>::data_ = NULL;
1135 }
MatrixIndexT stride_
< Number of rows
Definition: kaldi-matrix.h:816
Real * data_
data memory area
Definition: kaldi-matrix.h:808
MatrixIndexT num_rows_
< Number of columns
Definition: kaldi-matrix.h:813
#define KALDI_MEMALIGN_FREE(x)
Definition: kaldi-utils.h:60
MatrixIndexT num_cols_
these attributes store the real matrix size as it is stored in memory including memalignment ...
Definition: kaldi-matrix.h:812

◆ Init()

void Init ( const MatrixIndexT  r,
const MatrixIndexT  c,
const MatrixStrideType  stride_type 
)
inlineprivate

Init assumes the current class contents are invalid (i.e.

junk or have already been freed), and it sets the matrix to newly allocated memory with the specified number of rows and columns. r == c == 0 is acceptable. The data memory contents will be undefined.

Definition at line 783 of file kaldi-matrix.cc.

785  {
786  if (rows * cols == 0) {
787  KALDI_ASSERT(rows == 0 && cols == 0);
788  this->num_rows_ = 0;
789  this->num_cols_ = 0;
790  this->stride_ = 0;
791  this->data_ = NULL;
792  return;
793  }
794  KALDI_ASSERT(rows > 0 && cols > 0);
795  MatrixIndexT skip, stride;
796  size_t size;
797  void *data; // aligned memory block
798  void *temp; // memory block to be really freed
799 
800  // compute the size of skip and real cols
801  skip = ((16 / sizeof(Real)) - cols % (16 / sizeof(Real)))
802  % (16 / sizeof(Real));
803  stride = cols + skip;
804  size = static_cast<size_t>(rows) * static_cast<size_t>(stride)
805  * sizeof(Real);
806 
807  // allocate the memory and set the right dimensions and parameters
808  if (NULL != (data = KALDI_MEMALIGN(16, size, &temp))) {
809  MatrixBase<Real>::data_ = static_cast<Real *> (data);
812  MatrixBase<Real>::stride_ = (stride_type == kDefaultStride ? stride : cols);
813  } else {
814  throw std::bad_alloc();
815  }
816 }
MatrixIndexT stride_
< Number of rows
Definition: kaldi-matrix.h:816
Real * data_
data memory area
Definition: kaldi-matrix.h:808
int32 MatrixIndexT
Definition: matrix-common.h:98
MatrixIndexT num_rows_
< Number of columns
Definition: kaldi-matrix.h:813
#define KALDI_MEMALIGN(align, size, pp_orig)
Definition: kaldi-utils.h:58
MatrixIndexT num_cols_
these attributes store the real matrix size as it is stored in memory including memalignment ...
Definition: kaldi-matrix.h:812
#define KALDI_ASSERT(cond)
Definition: kaldi-error.h:185

◆ operator=() [1/2]

Matrix<Real>& operator= ( const MatrixBase< Real > &  other)
inline

Assignment operator that takes MatrixBase.

Definition at line 917 of file kaldi-matrix.h.

917  {
918  if (MatrixBase<Real>::NumRows() != other.NumRows() ||
919  MatrixBase<Real>::NumCols() != other.NumCols())
920  Resize(other.NumRows(), other.NumCols(), kUndefined);
922  return *this;
923  }
MatrixIndexT NumCols() const
Returns number of columns (or zero for empty matrix).
Definition: kaldi-matrix.h:67
void CopyFromMat(const MatrixBase< OtherReal > &M, MatrixTransposeType trans=kNoTrans)
Copy given matrix. (no resize is done).
MatrixIndexT NumRows() const
Returns number of rows (or zero for empty matrix).
Definition: kaldi-matrix.h:64
void Resize(const MatrixIndexT r, const MatrixIndexT c, MatrixResizeType resize_type=kSetZero, MatrixStrideType stride_type=kDefaultStride)
Sets matrix to a specified size (zero is OK as long as both r and c are zero).

◆ operator=() [2/2]

Matrix<Real>& operator= ( const Matrix< Real > &  other)
inline

Assignment operator. Needed for inclusion in std::vector.

Definition at line 926 of file kaldi-matrix.h.

926  {
927  if (MatrixBase<Real>::NumRows() != other.NumRows() ||
928  MatrixBase<Real>::NumCols() != other.NumCols())
929  Resize(other.NumRows(), other.NumCols(), kUndefined);
931  return *this;
932  }
MatrixIndexT NumCols() const
Returns number of columns (or zero for empty matrix).
Definition: kaldi-matrix.h:67
void CopyFromMat(const MatrixBase< OtherReal > &M, MatrixTransposeType trans=kNoTrans)
Copy given matrix. (no resize is done).
MatrixIndexT NumRows() const
Returns number of rows (or zero for empty matrix).
Definition: kaldi-matrix.h:64
void Resize(const MatrixIndexT r, const MatrixIndexT c, MatrixResizeType resize_type=kSetZero, MatrixStrideType stride_type=kDefaultStride)
Sets matrix to a specified size (zero is OK as long as both r and c are zero).

◆ Read()

void Read ( std::istream &  in,
bool  binary,
bool  add = false 
)

read from stream.

Definition at line 1450 of file kaldi-matrix.cc.

Referenced by OnlineUdpInput::Compute(), main(), kaldi::operator>>(), AffineXformStats::Read(), CompressedAffineXformStats::Read(), LogisticRegression::Read(), AccumFullGmm::Read(), LdaEstimate::Read(), Sgmm2FmllrGlobalParams::Read(), AccumDiagGmm::Read(), MleAmSgmm2Accs::Read(), FullGmm::Read(), Plda::Read(), CompressedMatrix::Read(), DiagGmm::Read(), Fmpe::Read(), OnlineGmmAdaptationState::Read(), OnlineCmvnState::Read(), IvectorExtractorStats::Read(), MatrixBase< float >::Read(), CuMatrix< float >::Read(), Matrix< BaseFloat >::Read(), KlHmm::ReadData(), kaldi::ReadKaldiObject(), kaldi::UnitTestCompressedMatrix(), kaldi::UnitTestComputeGPE(), UnitTestEndless1(), UnitTestEndless2(), kaldi::UnitTestIo(), kaldi::UnitTestIoCross(), UnitTestMono22K(), and UnitTestStereo8K().

1450  {
1451  if (add) {
1452  Matrix<Real> tmp;
1453  tmp.Read(is, binary, false); // read without adding.
1454  if (this->num_rows_ == 0) this->Resize(tmp.num_rows_, tmp.num_cols_);
1455  else {
1456  if (this->num_rows_ != tmp.num_rows_ || this->num_cols_ != tmp.num_cols_) {
1457  if (tmp.num_rows_ == 0) return; // do nothing in this case.
1458  else KALDI_ERR << "Matrix::Read, size mismatch "
1459  << this->num_rows_ << ", " << this->num_cols_
1460  << " vs. " << tmp.num_rows_ << ", " << tmp.num_cols_;
1461  }
1462  }
1463  this->AddMat(1.0, tmp);
1464  return;
1465  }
1466 
1467  // now assume add == false.
1468  MatrixIndexT pos_at_start = is.tellg();
1469  std::ostringstream specific_error;
1470 
1471  if (binary) { // Read in binary mode.
1472  int peekval = Peek(is, binary);
1473  if (peekval == 'C') {
1474  // This code enables us to read CompressedMatrix as a regular matrix.
1475  CompressedMatrix compressed_mat;
1476  compressed_mat.Read(is, binary); // at this point, add == false.
1477  this->Resize(compressed_mat.NumRows(), compressed_mat.NumCols());
1478  compressed_mat.CopyToMat(this);
1479  return;
1480  }
1481  const char *my_token = (sizeof(Real) == 4 ? "FM" : "DM");
1482  char other_token_start = (sizeof(Real) == 4 ? 'D' : 'F');
1483  if (peekval == other_token_start) { // need to instantiate the other type to read it.
1484  typedef typename OtherReal<Real>::Real OtherType; // if Real == float, OtherType == double, and vice versa.
1485  Matrix<OtherType> other(this->num_rows_, this->num_cols_);
1486  other.Read(is, binary, false); // add is false at this point anyway.
1487  this->Resize(other.NumRows(), other.NumCols());
1488  this->CopyFromMat(other);
1489  return;
1490  }
1491  std::string token;
1492  ReadToken(is, binary, &token);
1493  if (token != my_token) {
1494  if (token.length() > 20) token = token.substr(0, 17) + "...";
1495  specific_error << ": Expected token " << my_token << ", got " << token;
1496  goto bad;
1497  }
1498  int32 rows, cols;
1499  ReadBasicType(is, binary, &rows); // throws on error.
1500  ReadBasicType(is, binary, &cols); // throws on error.
1501  if ((MatrixIndexT)rows != this->num_rows_ || (MatrixIndexT)cols != this->num_cols_) {
1502  this->Resize(rows, cols);
1503  }
1504  if (this->Stride() == this->NumCols() && rows*cols!=0) {
1505  is.read(reinterpret_cast<char*>(this->Data()),
1506  sizeof(Real)*rows*cols);
1507  if (is.fail()) goto bad;
1508  } else {
1509  for (MatrixIndexT i = 0; i < (MatrixIndexT)rows; i++) {
1510  is.read(reinterpret_cast<char*>(this->RowData(i)), sizeof(Real)*cols);
1511  if (is.fail()) goto bad;
1512  }
1513  }
1514  if (is.eof()) return;
1515  if (is.fail()) goto bad;
1516  return;
1517  } else { // Text mode.
1518  std::string str;
1519  is >> str; // get a token
1520  if (is.fail()) { specific_error << ": Expected \"[\", got EOF"; goto bad; }
1521  // if ((str.compare("DM") == 0) || (str.compare("FM") == 0)) { // Back compatibility.
1522  // is >> str; // get #rows
1523  // is >> str; // get #cols
1524  // is >> str; // get "["
1525  // }
1526  if (str == "[]") { Resize(0, 0); return; } // Be tolerant of variants.
1527  else if (str != "[") {
1528  if (str.length() > 20) str = str.substr(0, 17) + "...";
1529  specific_error << ": Expected \"[\", got \"" << str << '"';
1530  goto bad;
1531  }
1532  // At this point, we have read "[".
1533  std::vector<std::vector<Real>* > data;
1534  std::vector<Real> *cur_row = new std::vector<Real>;
1535  while (1) {
1536  int i = is.peek();
1537  if (i == -1) { specific_error << "Got EOF while reading matrix data"; goto cleanup; }
1538  else if (static_cast<char>(i) == ']') { // Finished reading matrix.
1539  is.get(); // eat the "]".
1540  i = is.peek();
1541  if (static_cast<char>(i) == '\r') {
1542  is.get();
1543  is.get(); // get \r\n (must eat what we wrote)
1544  } else if (static_cast<char>(i) == '\n') { is.get(); } // get \n (must eat what we wrote)
1545  if (is.fail()) {
1546  KALDI_WARN << "After end of matrix data, read error.";
1547  // we got the data we needed, so just warn for this error.
1548  }
1549  // Now process the data.
1550  if (!cur_row->empty()) data.push_back(cur_row);
1551  else delete(cur_row);
1552  cur_row = NULL;
1553  if (data.empty()) { this->Resize(0, 0); return; }
1554  else {
1555  int32 num_rows = data.size(), num_cols = data[0]->size();
1556  this->Resize(num_rows, num_cols);
1557  for (int32 i = 0; i < num_rows; i++) {
1558  if (static_cast<int32>(data[i]->size()) != num_cols) {
1559  specific_error << "Matrix has inconsistent #cols: " << num_cols
1560  << " vs." << data[i]->size() << " (processing row"
1561  << i << ")";
1562  goto cleanup;
1563  }
1564  for (int32 j = 0; j < num_cols; j++)
1565  (*this)(i, j) = (*(data[i]))[j];
1566  delete data[i];
1567  data[i] = NULL;
1568  }
1569  }
1570  return;
1571  } else if (static_cast<char>(i) == '\n' || static_cast<char>(i) == ';') {
1572  // End of matrix row.
1573  is.get();
1574  if (cur_row->size() != 0) {
1575  data.push_back(cur_row);
1576  cur_row = new std::vector<Real>;
1577  cur_row->reserve(data.back()->size());
1578  }
1579  } else if ( (i >= '0' && i <= '9') || i == '-' ) { // A number...
1580  Real r;
1581  is >> r;
1582  if (is.fail()) {
1583  specific_error << "Stream failure/EOF while reading matrix data.";
1584  goto cleanup;
1585  }
1586  cur_row->push_back(r);
1587  } else if (isspace(i)) {
1588  is.get(); // eat the space and do nothing.
1589  } else { // NaN or inf or error.
1590  std::string str;
1591  is >> str;
1592  if (!KALDI_STRCASECMP(str.c_str(), "inf") ||
1593  !KALDI_STRCASECMP(str.c_str(), "infinity")) {
1594  cur_row->push_back(std::numeric_limits<Real>::infinity());
1595  KALDI_WARN << "Reading infinite value into matrix.";
1596  } else if (!KALDI_STRCASECMP(str.c_str(), "nan")) {
1597  cur_row->push_back(std::numeric_limits<Real>::quiet_NaN());
1598  KALDI_WARN << "Reading NaN value into matrix.";
1599  } else {
1600  if (str.length() > 20) str = str.substr(0, 17) + "...";
1601  specific_error << "Expecting numeric matrix data, got " << str;
1602  goto cleanup;
1603  }
1604  }
1605  }
1606  // Note, we never leave the while () loop before this
1607  // line (we return from it.)
1608  cleanup: // We only reach here in case of error in the while loop above.
1609  if(cur_row != NULL)
1610  delete cur_row;
1611  for (size_t i = 0; i < data.size(); i++)
1612  if(data[i] != NULL)
1613  delete data[i];
1614  // and then go on to "bad" below, where we print error.
1615  }
1616 bad:
1617  KALDI_ERR << "Failed to read matrix from stream. " << specific_error.str()
1618  << " File position at start is "
1619  << pos_at_start << ", currently " << is.tellg();
1620 }
MatrixIndexT NumCols() const
Returns number of columns (or zero for empty matrix).
Definition: kaldi-matrix.h:67
const Real * Data() const
Gives pointer to raw data (const).
Definition: kaldi-matrix.h:79
void ReadBasicType(std::istream &is, bool binary, T *t)
ReadBasicType is the name of the read function for bool, integer types, and floating-point types...
Definition: io-funcs-inl.h:55
Real * RowData(MatrixIndexT i)
Returns pointer to data for one row (non-const)
Definition: kaldi-matrix.h:87
void AddMat(const Real alpha, const MatrixBase< Real > &M, MatrixTransposeType transA=kNoTrans)
*this += alpha * M [or M^T]
kaldi::int32 int32
void ReadToken(std::istream &is, bool binary, std::string *str)
ReadToken gets the next token and puts it in str (exception on failure).
Definition: io-funcs.cc:154
void CopyFromMat(const MatrixBase< OtherReal > &M, MatrixTransposeType trans=kNoTrans)
Copy given matrix. (no resize is done).
int Peek(std::istream &is, bool binary)
Peek consumes whitespace (if binary == false) and then returns the peek() value of the stream...
Definition: io-funcs.cc:145
MatrixIndexT Stride() const
Stride (distance in memory between each row). Will be >= NumCols.
Definition: kaldi-matrix.h:70
int32 MatrixIndexT
Definition: matrix-common.h:98
MatrixIndexT num_rows_
< Number of columns
Definition: kaldi-matrix.h:813
#define KALDI_ERR
Definition: kaldi-error.h:147
MatrixIndexT num_cols_
these attributes store the real matrix size as it is stored in memory including memalignment ...
Definition: kaldi-matrix.h:812
#define KALDI_WARN
Definition: kaldi-error.h:150
friend class Matrix< Real >
Definition: kaldi-matrix.h:52
#define KALDI_STRCASECMP
Definition: kaldi-utils.h:147
void Resize(const MatrixIndexT r, const MatrixIndexT c, MatrixResizeType resize_type=kSetZero, MatrixStrideType stride_type=kDefaultStride)
Sets matrix to a specified size (zero is OK as long as both r and c are zero).

◆ RemoveRow()

void RemoveRow ( MatrixIndexT  i)

Remove a specified row.

Definition at line 1118 of file kaldi-matrix.cc.

Referenced by FullGmm::Merge(), DiagGmm::Merge(), DiagGmm::RemoveComponent(), FullGmm::RemoveComponent(), and kaldi::UnitTestRemoveRow().

1118  {
1119  KALDI_ASSERT(static_cast<UnsignedMatrixIndexT>(i) <
1120  static_cast<UnsignedMatrixIndexT>(MatrixBase<Real>::num_rows_)
1121  && "Access out of matrix");
1122  for (MatrixIndexT j = i + 1; j < MatrixBase<Real>::num_rows_; j++)
1123  MatrixBase<Real>::Row(j-1).CopyFromVec( MatrixBase<Real>::Row(j));
1125 }
int32 MatrixIndexT
Definition: matrix-common.h:98
const SubVector< Real > Row(MatrixIndexT i) const
Return specific row of matrix [const].
Definition: kaldi-matrix.h:188
MatrixIndexT num_rows_
< Number of columns
Definition: kaldi-matrix.h:813
#define KALDI_ASSERT(cond)
Definition: kaldi-error.h:185

◆ Resize()

void Resize ( const MatrixIndexT  r,
const MatrixIndexT  c,
MatrixResizeType  resize_type = kSetZero,
MatrixStrideType  stride_type = kDefaultStride 
)

Sets matrix to a specified size (zero is OK as long as both r and c are zero).

The value of the new data depends on resize_type: -if kSetZero, the new data will be zero -if kUndefined, the new data will be undefined -if kCopyData, the new data will be the same as the old data in any shared positions, and zero elsewhere.

You can set stride_type to kStrideEqualNumCols to force the stride to equal the number of columns; by default it is set so that the stride in bytes is a multiple of 16.

This function takes time proportional to the number of data elements.

Definition at line 819 of file kaldi-matrix.cc.

Referenced by DecodableMatrixMappedOffset::AcceptLoglikes(), OnlinePitchFeatureImpl::AcceptWaveform(), LdaEstimate::AddMeanOffset(), kaldi::AddToConfusionMatrix(), DecodableNnetLoopedOnlineBase::AdvanceChunk(), DecodableNnetSimpleLooped::AdvanceChunk(), kaldi::nnet2::AppendDiscriminativeExamples(), kaldi::AppendFeats(), OnlineDeltaInput::AppendFrames(), kaldi::AppendPostToFeats(), kaldi::AppendVectorToFeats(), kaldi::ApplyPca(), Sgmm2Project::ApplyProjection(), Plda::ApplyTransform(), BatchedXvectorComputer::BatchedXvectorComputer(), AmSgmm2::ComponentLogLikes(), kaldi::ComposeTransforms(), OnlineMatrixInput::Compute(), NnetComputerFromEg::Compute(), OnlineCmnInput::Compute(), DecodableAmNnetParallel::Compute(), OfflineFeatureTpl< F >::Compute(), OnlineUdpInput::Compute(), OnlineLdaInput::Compute(), OnlineDeltaInput::Compute(), OnlineFeInput< E >::Compute(), kaldi::ComputeAmGmmFeatureDeriv(), kaldi::ComputeAndProcessKaldiPitch(), kaldi::ComputeDeltas(), IvectorExtractor::ComputeDerivedVars(), kaldi::ComputeFeatureNormalizingTransform(), Fmpe::ComputeFeatures(), AmSgmm2::ComputeFmllrPreXform(), DecodableNnet2Online::ComputeForFrame(), OnlineCmnInput::ComputeInternal(), kaldi::ComputeKaldiPitch(), kaldi::ComputeKaldiPitchFirstPass(), Sgmm2Project::ComputeLdaTransform(), OnlineLdaInput::ComputeNextRemainder(), kaldi::ComputePca(), AmSgmm2::ComputePerSpkDerivedVars(), Sgmm2Project::ComputeProjection(), kaldi::ComputeShiftedDeltas(), Fmpe::ComputeStddevs(), BasisFmllrEstimate::ComputeTransform(), kaldi::ConcatFeats(), FmllrRawAccs::ConvertToPerRowStats(), CompressedAffineXformStats::CopyFromAffineXformStats(), SvdApplier::DecomposeComponent(), OnlineDeltaInput::DeltaComputation(), DiscriminativeExampleSplitter::DoExcise(), LdaEstimate::Estimate(), FeatureTransformEstimateMulti::Estimate(), FeatureTransformEstimate::EstimateInternal(), kaldi::EstimateIvectorsOnline(), FeatureTransformEstimateMulti::EstimateTransformPart(), kaldi::EstPca(), kaldi::ExtractObjectRange(), kaldi::FilterCompressedMatrixRows(), kaldi::FilterMatrixRows(), FmllrRawAccs::FmllrRawAccs(), Fmpe::Fmpe(), kaldi::nnet2::FormatNnetInput(), NnetBatchComputer::FormatOutputs(), OnlineFeaturePipeline::GetAsMatrix(), OnlineCacheInput::GetCachedData(), FullGmm::GetCovarsAndMeans(), kaldi::GetFeatDeriv(), OnlineCmvn::GetFrame(), LogisticRegression::GetLogPosteriors(), GeneralMatrix::GetMatrix(), FullGmm::GetMeans(), DiagGmm::GetMeans(), IvectorExtractorStats::GetOrthogonalIvectorTransform(), kaldi::GetOutput(), PldaEstimator::GetOutput(), OnlineCmvn::GetState(), DiagGmm::GetVars(), RegtreeFmllrDiagGmm::GetXformMatrix(), AmSgmm2::IncreasePhoneSpaceDim(), kaldi::IncreaseTransformDimension(), AffineXformStats::Init(), LdaEstimate::Init(), OnlineFeaturePipeline::Init(), kaldi::InitCmvnStats(), OnlinePreconditionerSimple::InitDefault(), OnlineNaturalGradientSimple::InitDefault(), kaldi::InitFmllr(), AmSgmm2::InitializeMw(), kaldi::InitIdftBases(), IvectorExtractor::IvectorExtractor(), IvectorExtractorStats::IvectorExtractorStats(), MatrixBase< float >::LapackGesvd(), kaldi::nnet2::LatticeToDiscriminativeExample(), AffineComponent::LimitRank(), DiagGmm::LogLikelihoods(), main(), FullGmm::Merge(), DiagGmm::Merge(), FullGmm::MergePreselect(), kaldi::nnet3::MergeTaskOutput(), MfccComputer::MfccComputer(), LogisticRegression::MixUp(), OnlineLdaInput::OnlineLdaInput(), LimitRankClass::operator()(), kaldi::PosteriorToMatrix(), kaldi::PosteriorToPdfMatrix(), kaldi::ProcessPitch(), AffineXformStats::Read(), kaldi::ReadData(), NnetExampleBackgroundReader::ReadExamples(), kaldi::ReadHtk(), kaldi::nnet3::ReduceRankOfComponents(), DiagGmmNormal::Resize(), FullGmm::Resize(), FullGmmNormal::Resize(), DiagGmm::Resize(), AccumFullGmm::Resize(), AccumDiagGmm::Resize(), Sgmm2PerFrameDerivedVars::Resize(), MleAmSgmm2Accs::ResizeAccumulators(), kaldi::ReverseFrames(), SingleUtteranceNnet2DecoderThreaded::RunNnetEvaluationInternal(), KlHmm::SetStats(), LogisticRegression::SetWeights(), kaldi::SpliceFrames(), OnlineLdaInput::SpliceFrames(), FullGmm::Split(), DiagGmm::Split(), CuMatrix< float >::Swap(), TestSgmm2Fmllr(), LogisticRegression::Train(), kaldi::TransformIvectors(), OnlineLdaInput::TransformToOutput(), kaldi::TypeThreeUsage(), kaldi::UnitTestAddMat2(), kaldi::UnitTestDelay(), UnitTestEstimateSgmm2(), kaldi::UnitTestIo(), UnitTestMleAmDiagGmm(), kaldi::UnitTestPieces(), kaldi::UnitTestResize(), kaldi::UnitTestResizeCopyDataDifferentStrideType(), kaldi::UnitTestTransposeScatter(), and kaldi::nnet2::UpdateHash().

822  {
823  // the next block uses recursion to handle what we have to do if
824  // resize_type == kCopyData.
825  if (resize_type == kCopyData) {
826  if (this->data_ == NULL || rows == 0) resize_type = kSetZero; // nothing to copy.
827  else if (rows == this->num_rows_ && cols == this->num_cols_ &&
828  (stride_type == kDefaultStride || this->stride_ == this->num_cols_)) { return; } // nothing to do.
829  else {
830  // set tmp to a matrix of the desired size; if new matrix
831  // is bigger in some dimension, zero it.
832  MatrixResizeType new_resize_type =
833  (rows > this->num_rows_ || cols > this->num_cols_) ? kSetZero : kUndefined;
834  Matrix<Real> tmp(rows, cols, new_resize_type, stride_type);
835  MatrixIndexT rows_min = std::min(rows, this->num_rows_),
836  cols_min = std::min(cols, this->num_cols_);
837  tmp.Range(0, rows_min, 0, cols_min).
838  CopyFromMat(this->Range(0, rows_min, 0, cols_min));
839  tmp.Swap(this);
840  // and now let tmp go out of scope, deleting what was in *this.
841  return;
842  }
843  }
844  // At this point, resize_type == kSetZero or kUndefined.
845 
846  if (MatrixBase<Real>::data_ != NULL) {
847  if (rows == MatrixBase<Real>::num_rows_
848  && cols == MatrixBase<Real>::num_cols_) {
849  if (resize_type == kSetZero)
850  this->SetZero();
851  return;
852  }
853  else
854  Destroy();
855  }
856  Init(rows, cols, stride_type);
857  if (resize_type == kSetZero) MatrixBase<Real>::SetZero();
858 }
MatrixResizeType
Definition: matrix-common.h:37
MatrixIndexT stride_
< Number of rows
Definition: kaldi-matrix.h:816
Real * data_
data memory area
Definition: kaldi-matrix.h:808
void Init(const MatrixIndexT r, const MatrixIndexT c, const MatrixStrideType stride_type)
Init assumes the current class contents are invalid (i.e.
void CopyFromMat(const MatrixBase< OtherReal > &M, MatrixTransposeType trans=kNoTrans)
Copy given matrix. (no resize is done).
int32 MatrixIndexT
Definition: matrix-common.h:98
MatrixIndexT num_rows_
< Number of columns
Definition: kaldi-matrix.h:813
MatrixIndexT num_cols_
these attributes store the real matrix size as it is stored in memory including memalignment ...
Definition: kaldi-matrix.h:812
void SetZero()
Sets matrix to zero.
void Destroy()
Deallocates memory and sets to empty matrix (dimension 0, 0).
friend class Matrix< Real >
Definition: kaldi-matrix.h:52
SubMatrix< Real > Range(const MatrixIndexT row_offset, const MatrixIndexT num_rows, const MatrixIndexT col_offset, const MatrixIndexT num_cols) const
Return a sub-part of matrix.
Definition: kaldi-matrix.h:202

◆ Swap() [1/2]

void Swap ( Matrix< Real > *  other)

Swaps the contents of *this and *other. Shallow swap.

Definition at line 2255 of file kaldi-matrix.cc.

Referenced by DecodableNnetLoopedOnlineBase::AdvanceChunk(), DecodableNnetSimpleLooped::AdvanceChunk(), OnlineCmnInput::ComputeInternal(), DecodableAmNnet::DecodableAmNnet(), NnetExampleBackgroundReader::GetNextMinibatch(), Matrix< BaseFloat >::Resize(), SingleUtteranceNnet2DecoderThreaded::RunNnetEvaluationInternal(), CuMatrix< float >::Swap(), GeneralMatrix::SwapFullMatrix(), and kaldi::UnitTestIvectorExtractor().

2255  {
2256  std::swap(this->data_, other->data_);
2257  std::swap(this->num_cols_, other->num_cols_);
2258  std::swap(this->num_rows_, other->num_rows_);
2259  std::swap(this->stride_, other->stride_);
2260 }
MatrixIndexT stride_
< Number of rows
Definition: kaldi-matrix.h:816
Real * data_
data memory area
Definition: kaldi-matrix.h:808
void swap(basic_filebuf< CharT, Traits > &x, basic_filebuf< CharT, Traits > &y)
MatrixIndexT num_rows_
< Number of columns
Definition: kaldi-matrix.h:813
MatrixIndexT num_cols_
these attributes store the real matrix size as it is stored in memory including memalignment ...
Definition: kaldi-matrix.h:812

◆ Swap() [2/2]

void Swap ( CuMatrix< Real > *  mat)

Defined in ../cudamatrix/cu-matrix.cc.

Definition at line 3160 of file cu-matrix.cc.

3160 { mat->Swap(this); }

◆ Transpose()

void Transpose ( )

Transpose the matrix.

Works for non-square matrices as well as square ones.

Definition at line 2091 of file kaldi-matrix.cc.

Referenced by kaldi::ComputePca(), BasisFmllrEstimate::ComputeTransform(), kaldi::CuVectorUnitTestCopyDiagFromMat(), BasisFmllrEstimate::EstimateFmllrBasis(), kaldi::UnitTestAddDiagVecMat(), kaldi::UnitTestAddMat2Sp(), kaldi::UnitTestAddMatSmat(), kaldi::UnitTestCuMatrixTranspose(), kaldi::UnitTestSymAddMat2(), kaldi::UnitTestTranspose(), kaldi::UnitTestVecmul(), and MlltAccs::Update().

2091  {
2092  if (this->num_rows_ != this->num_cols_) {
2093  Matrix<Real> tmp(*this, kTrans);
2094  Resize(this->num_cols_, this->num_rows_);
2095  this->CopyFromMat(tmp);
2096  } else {
2097  (static_cast<MatrixBase<Real>&>(*this)).Transpose();
2098  }
2099 }
void Transpose()
Transpose the matrix.
void CopyFromMat(const MatrixBase< OtherReal > &M, MatrixTransposeType trans=kNoTrans)
Copy given matrix. (no resize is done).
MatrixIndexT num_rows_
< Number of columns
Definition: kaldi-matrix.h:813
MatrixIndexT num_cols_
these attributes store the real matrix size as it is stored in memory including memalignment ...
Definition: kaldi-matrix.h:812
friend class Matrix< Real >
Definition: kaldi-matrix.h:52
void Resize(const MatrixIndexT r, const MatrixIndexT c, MatrixResizeType resize_type=kSetZero, MatrixStrideType stride_type=kDefaultStride)
Sets matrix to a specified size (zero is OK as long as both r and c are zero).

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