mmat_10.cc

这是一个从音频信号里提取特征参量的程序

(const MMatrix<ISIP_TEMPLATE_TARGET>&, MVector<ISIP_TEMPLATE_TARGET>&, long);

// method: getColumn
//
// arguments:
//  const MMatrix<TScalar, TIntegral>& this: (output) class operand
//  MVector<TScalar, TIntegral>& vector: (output) a vector at the given column
//  long col_index: (input) column index
//
// return: a boolean value indicating status
//
template<class TScalar, class TIntegral>
boolean MMatrixMethods::getColumn(const MMatrix<TScalar, TIntegral>& this_a, 
				  MVector<TScalar, TIntegral>& vector_a,
				  long col_index_a) {

  // check the argument
  //
  if ((col_index_a < 0) || (col_index_a > this_a.ncols_d)) {
    this_a.debug(L"this_a");
    vector_a.debug(L"vector_a");    
    return Error::handle(name(), L"getColumn", Error::ARG, __FILE__, __LINE__);
  }

  // create space for the output vector
  //
  vector_a.setLength(this_a.nrows_d, false);

  // since most of our matrices are stored in a row-major format, the
  // easiest way to do this is to simply loop over the elements and
  // copy them one by one
  //
  for (long row = 0; row < this_a.nrows_d; row++) {
    vector_a(row) = this_a.getValue(row, col_index_a);
  }
  
  // exit gracefully
  //
  return true;
}

// explicit instantiations
//
template boolean
MMatrixMethods::getColumn<ISIP_TEMPLATE_TARGET>
(const MMatrix<ISIP_TEMPLATE_TARGET>&, MVector<ISIP_TEMPLATE_TARGET>&, long);

// method: setValue
//
// arguments:
//  MMatrix<TScalar, TIntegral>& this: (output) class operand
//  long row_index: (input) row_index
//  long col_index: (input) col_index
//  TIntegral value: (input) input value
//
// return: a boolean value indicating status
//
// this method sets the specified position of the matrix to the input value
//
template<class TScalar, class TIntegral>
boolean MMatrixMethods::setValue(MMatrix<TScalar, TIntegral>& this_a, 
				 long row_index_a, long col_index_a, 
				 TIntegral value_a) {

  // check the arguments
  //
  if ((row_index_a >= this_a.getNumRows()) ||
      (col_index_a >= this_a.getNumColumns())) {
    this_a.debug(L"this_a");    
    return Error::handle(name(), L"setValue()", Error::ARG, 
                         __FILE__, __LINE__);
  }

  // type: DIAGONAL
  //  this general idea here is that it is only allowed to set non-zero value
  //  for diagonal element, and zero for any element.
  //
  if (this_a.type_d == Integral::DIAGONAL) {

    // setting the non-zero value for off-diagonal element will change
    // the type of the matrix, so error
    //
    if (row_index_a != col_index_a ) {
      if (value_a != 0) {
	this_a.debug(L"this_a");	
        return Error::handle(name(), L"setValue()", 
			     MMatrix<TScalar, TIntegral>::ERR_OPTYPE, 
                             __FILE__, __LINE__);
      }
    }

    // set the value
    //
    else {
      this_a.m_d(row_index_a) = value_a;
    }
  }
  
  // type: FULL
  //
  else if (this_a.type_d == Integral::FULL) {
    (this_a.m_d(this_a.index(row_index_a, col_index_a))) = value_a;
  }

  // type: SYMMETRIC
  //  as we only keep elements of symmetric matrix in lower triangular, for
  //  upper triangular elements, we have a(i, j) = a(j, i) if i < j. so the
  //  general idea here is that we first check whether this position is lower
  //  or upper triangular, and then set the corresponding value in lower
  //  triangle.
  //
  else if (this_a.type_d == Integral::SYMMETRIC) {

    // if the position is in lower triangular part of the matrix, set the value
    //
    if (row_index_a > col_index_a) {
      (this_a.m_d(this_a.index(row_index_a, col_index_a))) = value_a;      
    }

    // else set the corresponding value in lower triangle
    //
    else {
      (this_a.m_d(this_a.index(col_index_a, row_index_a))) = value_a;
    }
  }

  // type: LOWER_TRIANGULAR
  //  this general idea here is that it is only allowed to set non-zero value
  //  for lower_triangular element, and zero for any element.
  //
  else if (this_a.type_d == Integral::LOWER_TRIANGULAR) {

    // set the value if the position is in lower triangular portion of matrix
    //
    if (row_index_a >= col_index_a) {
      (this_a.m_d(this_a.index(row_index_a, col_index_a))) = value_a;      
    }

    // error if a non-zero value is to be set in upper triangular
    // portion of the matrix, as it will change the type of the matrix
    //
    else if (value_a != 0) {
      this_a.debug(L"this_a");      
      return Error::handle(name(), L"setValue()", 
			   MMatrix<TScalar, TIntegral>::ERR_OPTYPE, 
			   __FILE__, __LINE__);
    }
  }

  // type: UPPER_TRIANGULAR
  //  this general idea here is that it is only allowed to set non-zero value
  //  for upper_triangular element, and zero for any element.
  //  
  else if (this_a.type_d == Integral::UPPER_TRIANGULAR) {

    // set the value if the position is in upper triangular portion of matrix
    //    
    if (row_index_a <= col_index_a) {
      (this_a.m_d(this_a.index(row_index_a, col_index_a))) = value_a;      
      return true;
    }

    // error if a non-zero value is to be set in lower triangular
    // portion of the matrix, as it will change the type of the matrix
    //    
    else if (value_a != 0) {
      this_a.debug(L"this_a");      
      return Error::handle(name(), L"setValue", 
			   MMatrix<TScalar, TIntegral>::ERR_OPTYPE, 
			   __FILE__, __LINE__);
    }
  }
  
  // type: SPARSE
  //  first we need to check whether this element exists or not. the general
  //  idea is we replace an existing value or insert a new value.
  //
  else if (this_a.type_d == Integral::SPARSE) {

    // find the position in the matrix where we can insert the value:
    //  the matrix elements are ordered by row (all elements for the
    //  first row followed by all elements for the second row ...). the
    //  value returned by this code, index, is the place at which the
    //  value should replace an existing value or insert a new value.
    // 
    //
    long index = 0;
    long num_elements = this_a.m_d.length();

    while ((index < this_a.col_index_d.length()) &&
           ((row_index_a > this_a.row_index_d(index)) ||
	    ((row_index_a == this_a.row_index_d(index)) &&
	     (col_index_a > this_a.col_index_d(index))))) {
      index++;
    }
    
    // if the element doesn't exist, or the matrix is empty, add it
    //
    if ((index >= num_elements) ||
	(row_index_a != this_a.row_index_d(index)) ||
	(col_index_a != this_a.col_index_d(index))) {
      
      // no need to add zero into sparse matrix
      //
      if (value_a == 0) {

	// exit gracefully
	//
	return true;
      }

      // increase the length by one
      //
      long new_len = num_elements + 1;
      this_a.m_d.setLength(new_len);
      this_a.row_index_d.setLength(new_len);
      this_a.col_index_d.setLength(new_len);

      // shift the data up by one
      //
      for (long i = num_elements; i > index; i--) {
	this_a.m_d(i) = this_a.m_d(i - 1);
	this_a.row_index_d(i) = this_a.row_index_d(i - 1);
	this_a.col_index_d(i) = this_a.col_index_d(i - 1);
      }
    
      // assign the new value
      //
      this_a.m_d(index) = value_a;
      this_a.row_index_d(index) = row_index_a;
      this_a.col_index_d(index) = col_index_a;
    }

    // if this element already exists, modify it
    //
    else {
      
      // if the input value is set to zero, remove this value
      //
      if (value_a == 0) {

	// shift the data down by one
	//
        for (long i = index; i < num_elements - 1; i++) {
          this_a.m_d(i) = this_a.m_d(i + 1);
          this_a.row_index_d(i) = this_a.row_index_d(i + 1);
          this_a.col_index_d(i) = this_a.col_index_d(i + 1);
        }

	// decrease the length by one
	//
        long new_len = num_elements - 1;
        this_a.m_d.setLength(new_len);
        this_a.row_index_d.setLength(new_len);
        this_a.col_index_d.setLength(new_len);
      }
      
      // the value is not zero, and the element can be replaced
      //
      else {
        this_a.m_d(index) = value_a;
      }
    }
  }

  // exit gracefully
  //
  return true;
}

// explicit instantiations
//
template boolean
MMatrixMethods::setValue<ISIP_TEMPLATE_TARGET>
(MMatrix<ISIP_TEMPLATE_TARGET>&, long, long, ISIP_TEMPLATE_T1);

// method: setRow
//
// arguments:
//  MMatrix<TScalar, TIntegral>& this: (output) class operand
//  long row_index: (input) row index for the row which need to be set
//  const MVector<TScalar, TIntegral>& vector: (input) a row vector
//
// return: a boolean value indicating status
//
// this method sets the values of the row specified by the row index
//
template<class TScalar, class TIntegral>
boolean MMatrixMethods::setRow(MMatrix<TScalar, TIntegral>& this_a, 
			       long row_index_a,
			       const MVector<TScalar, TIntegral>& vector_a) { 

  // check the argument
  //
  if ((row_index_a < 0) || (row_index_a >= this_a.nrows_d) ||
      (this_a.ncols_d != vector_a.length())) {
    this_a.debug(L"this_a");
    vector_a.debug(L"vector_a");    
    return Error::handle(name(), L"setRow", Error::ARG, __FILE__, __LINE__);
  }
  
  // get the number of columns of the matrix
  //
  long ncols = this_a.getNumColumns();

  // type: FULL
  //  the general idea is that we loop through the elements of the specified
  //  row and assign the values.
  //
  if (this_a.type_d == Integral::FULL) {

    // get the starting index
    //
    long index = this_a.index(row_index_a, 0);

    // copy the elements
    //
    for (long col = 0; col < ncols; col++) {  
      this_a.m_d(index++) = vector_a(col);
    }
  }

  // type: DIAGONAL
  //  user shouldn't use setRow for diagonal matrix as it will probably
  //  break its diagonal property.
  //
  else if (this_a.type_d == Integral::DIAGONAL) {
    this_a.debug(L"this_a");    
    return Error::handle(name(), L"setRow", Error::ARG, __FILE__, __LINE__);
  }

  // type: SYMMETRIC or LOWER_TRIANGULAR
  //  note that we only need to set the lower triangular part of the row  
  //
  else if ((this_a.type_d == Integral::SYMMETRIC) ||
	   (this_a.type_d == Integral::LOWER_TRIANGULAR)) {

    // get the starting index
    //
    long index = this_a.index(row_index_a, 0);

    // copy the elements
    //
    for (long col = 0; col <= row_index_a; col++) {  
      this_a.m_d(index++) = vector_a(col);
    }
  }

  // type: UPPER_TRIANGULAR
  //  note that we only need to set the upper triangular part of the row  
  //  since this format is stored in a column-oriented order, we use
  //  the indexing function.
  //
  else if (this_a.type_d == Integral::UPPER_TRIANGULAR) {

    // copy each element
    //    
    for (long col = row_index_a; col < ncols; col++) {
      this_a.m_d(this_a.index(row_index_a, col)) = vector_a(col);
    }
  }

  // type: SPARSE
  //  we try to do this efficiently by copying moving the elements in blocks
  //
  else {

    // get the length of non-zero elements in sparse matrix
    //
    long num_elements = this_a.m_d.length();

    // loop through each element to search the start index to insert
    // the input vector. 
    // 
    long start_index = Integral::NO_POS;

    for (long index = 0; index < num_elements; index++) {
      if (this_a.row_index_d(index) >= row_index_a) {
	start_index = index;
	break;
      }
    }

    if (start_index == Integral::NO_POS) {

      // will insert the input vector at the end of storing vectors of matrix
      //
      start_index = num_elements;
    }

    // copy the existing data
    //
    MVector<TScalar, TIntegral> tmp_m(this_a.m_d);

    VectorLong tmp_ri;
    tmp_ri.assign(this_a.row_index_d);

    VectorLong tmp_ci;
    tmp_ci.assign(this_a.col_index_d);    

    // get the maximum length of the number of non-zero elements in the
    // new sparse matrix
    //
    long max_num_elem = num_elements + vector_a.numNotEqual(0);

    // increase the capacity in preparation for the data copies
    //
    this_a.m_d.setLength(max_num_elem);
    this_a.row_index_d.setLength(max_num_elem);
    this_a.col_index_d.setLength(max_num_elem);

    // loop over all elements in the new vector to copy the non-zero elements
    // to current matrix
    //
    long index = start_index;

    for (long col = 0; col < ncols; col++) {

      // increase the length and copy
      //
      if ((TIntegral)(vector_a(col)) != 0) {
	this_a.m_d(index) = vector_a(col);
	this_a.row_index_d(index) = row_index_a;
	this_a.col_index_d(index++) = col;
      }
    }

    // from start_index, search the starting index of the next row
    // 
    long next_index = Integral::NO_POS;
    long next_row = row_index_a + 1;

    for (long i = start_index; i < num_elements; i++) {
      if (tmp_ri(i) >= next_row) {
	next_index = i;
	break;
      }
    }

    // is there any elements remaining to be copied?
    //
    long num_remain;

    if (next_index != Integral::NO_POS) {

      // compute the number of the rest elements
      //
      num_remain = num_elements - next_index;
    
      // copy those elements from [next_index to num_elements]
      //
      this_a.m_d.move(tmp_m, num_remain, next_index, index);
      this_a.row_index_d.move(tmp_ri, num_remain, next_index, index);
      this_a.col_index_d.move(tmp_ci, num_remain, next_index, index);
    }

    // no elements remaining to be copied
    //
    else {
      num_remain = 0;