cmatrix.c

大地测量专业计算软件

    MTX_ERROR_MSG( "MTX_Malloc returned FALSE." );
    return FALSE;
  }

  for( i = 0; i < M->nrows; i++ )
  {
    Re->data[0][i] = M->cplx[col][i].re;
  }

  return TRUE;
}

BOOL MTX_Imag( const MTX *M, MTX *Im )
{
  unsigned i = 0;
  unsigned j = 0;

  if( MTX_isNull( M ) )
  {
    MTX_ERROR_MSG( "NULL Matrix" );
    return FALSE;
  }

  if( M->isReal )
    return MTX_Calloc( Im, M->nrows, M->ncols, TRUE ); // return a zero matrix

  if( !MTX_Malloc( Im, M->nrows, M->ncols, TRUE ) )
  {
    MTX_ERROR_MSG( "MTX_Malloc returned FALSE." );
    return FALSE;
  }

  for( j = 0; j < M->ncols; j++ )
  {
    for( i = 0; i < M->nrows; i++ )
    {
      Im->data[j][i] = M->cplx[j][i].im;
    }
  }
  return TRUE;
}

BOOL MTX_ImagColumn( const MTX *M, const unsigned col, MTX *Im )
{
  unsigned i = 0;

  if( MTX_isNull( M ) )
  {
    MTX_ERROR_MSG( "NULL Matrix" );
    return FALSE;
  }

  if( col >= M->ncols )
  {
    MTX_ERROR_MSG( "if( col >= M->ncols )" );
    return FALSE;
  }

  if( M->isReal )
    return MTX_Calloc( Im, M->nrows, 1, TRUE ); // return a zero column

  if( !MTX_Malloc( Im, M->nrows, 1, TRUE ) )
  {
    MTX_ERROR_MSG( "MTX_Malloc returned FALSE." );
    return FALSE;
  }

  for( i = 0; i < M->nrows; i++ )
  {
    Im->data[0][i] = M->cplx[col][i].im;
  }
  return TRUE;
}

BOOL MTX_Magnitude( const MTX *M, MTX *Magnitude )
{
  unsigned i = 0;
  unsigned j = 0;
  double re;
  double im;

  if( MTX_isNull( M ) )
  {
    MTX_ERROR_MSG( "NULL Matrix" );
    return FALSE;
  }

  if( M->isReal )
  {
    if( !MTX_Copy( M, Magnitude ) )
    {
      MTX_ERROR_MSG( "MTX_Copy returned FALSE." );
      return FALSE;
    }
    if( !MTX_Abs( Magnitude ) )
    {
      MTX_ERROR_MSG( "MTX_Abs returned FALSE." );
      return FALSE;
    }
    return TRUE;
  }

  if( !MTX_Malloc( Magnitude, M->nrows, M->ncols, TRUE ) )
  {
    MTX_ERROR_MSG( "MTX_Malloc returned FALSE." );
    return FALSE;
  }

  for( j = 0; j < M->ncols; j++ )
  {
    for( i = 0; i < M->nrows; i++ )
    {
      re = M->cplx[j][i].re;
      im = M->cplx[j][i].im;
      Magnitude->data[j][i] = sqrt( re*re + im*im );
    }
  }
  return TRUE;
}

BOOL MTX_Phase( const MTX *M, MTX *Phase )
{
  unsigned i = 0;
  unsigned j = 0;

  if( MTX_isNull( M ) )
  {
    MTX_ERROR_MSG( "NULL Matrix" );
    return FALSE;
  }

  if( !MTX_Calloc( Phase, M->nrows, M->ncols, TRUE ) )
  {
    MTX_ERROR_MSG( "MTX_Calloc returned FALSE." );
    return FALSE;
  }

  if( M->isReal )
  {
    for( j = 0; j < M->ncols; j++ )
    {
      for( i = 0; i < M->nrows; i++ )
      {
        if( M->data[j][i] < 0.0 )
        {
          Phase->data[j][i] = PI;
        }
        else
        {
          Phase->data[j][i] = 0.0;
        }
      }
    }
    return TRUE;
  }

  for( j = 0; j < M->ncols; j++ )
  {
    for( i = 0; i < M->nrows; i++ )
    {
      Phase->data[j][i] = atan2( M->cplx[j][i].im, M->cplx[j][i].re );
    }
  }
  return TRUE;
}


BOOL MTX_Conjugate( MTX *M )
{
  unsigned i=0;
  unsigned j=0;

  if( MTX_isNull( M ) )
  {
    MTX_ERROR_MSG( "NULL Matrix" );
    return FALSE;
  }

  if( M->isReal )
    return TRUE;

  for( j = 0; j < M->ncols; j++ )
  {
    for( i = 0; i < M->nrows; i++ )
    {
      M->cplx[j][i].im = -M->cplx[j][i].im;
    }
  }
  return TRUE;
}

BOOL MTX_RemoveColumn( MTX *M, const unsigned col )
{
  unsigned j = 0;
  unsigned k = 0;
  double **dptr = NULL;
  stComplex **cptr = NULL;

  if( MTX_isNull( M ) )
  {
    MTX_ERROR_MSG( "NULL Matrix" );
    return FALSE;
  }

  if( col >= M->ncols )
  {
    MTX_ERROR_MSG( "if( col >= M->ncols )" );
    return FALSE;
  }

  // special case
  if( M->ncols == 1 )
  {
    return MTX_Free( M );
  }

  // allocate a new array of column vectors
  if( M->isReal )
  {
    dptr = (double**)malloc( (M->ncols-1)*sizeof(double*) );
    if( !dptr )
    {
      MTX_ERROR_MSG( "malloc returned NULL." );
      return FALSE;
    }
  }
  else
  {
    cptr = (stComplex**)malloc( (M->ncols-1)*sizeof(stComplex*) );
    if( !cptr )
    {
      MTX_ERROR_MSG( "malloc returned NULL." );
      return FALSE;
    }
  }

  // copy the previous array of pointers
  // except the one to remove
  k = 0;
  for( j = 0; j < M->ncols; j++ )
  {
    if( j != col )
    {
      if( M->isReal )
        dptr[k] = M->data[j];
      else
        cptr[k] = M->cplx[j];
      k++;
    }
    else
    {
      if( M->isReal )
        free( M->data[j] );
      else
        free( M->cplx[j] );
    }
  }

  // free the old column array, and copy the new
  if( M->isReal )
  {
    free( M->data );
    M->data = dptr;
  }
  else
  {
    free( M->cplx );
    M->cplx = cptr;
  }
  M->ncols--;

  return TRUE;
}

BOOL MTX_RemoveColumnsAfterIndex( MTX *dst, const unsigned col )
{
  unsigned ncols;
  unsigned j = 0;
  double **dptr = NULL;
  stComplex **cptr = NULL;

  if( MTX_isNull( dst ) )
  {
    MTX_ERROR_MSG( "NULL Matrix" );
    return FALSE;
  }

  if( col >= dst->ncols )
  {
    MTX_ERROR_MSG( "if( col >= dst->ncols )" );
    return FALSE;
  }

  // special case
  if( dst->ncols == 1 )
  {
    return MTX_Free( dst );
  }

  ncols = col+1;

  // allocate a new array of column vectors
  if( dst->isReal )
  {
    dptr = (double**)malloc( ncols*sizeof(double*) );
    if( !dptr )
    {
      MTX_ERROR_MSG( "malloc returned NULL." );
      return FALSE;
    }
  }
  else
  {
    cptr = (stComplex**)malloc( ncols*sizeof(stComplex*) );
    if( !cptr )
    {
      MTX_ERROR_MSG( "malloc returned NULL." );
      return FALSE;
    }
  }

  for( j = 0; j < dst->ncols; j++ )
  {
    if( j < ncols )
    {
      if( dst->isReal )
        dptr[j] = dst->data[j];
      else
        cptr[j] = dst->cplx[j];
    }
    else
    {
      if( dst->isReal )
        free( dst->data[j] );
      else
        free( dst->cplx[j] );
    }
  }

  // free the old column array, and copy the new
  if( dst->isReal )
  {
    free( dst->data );
    dst->data = dptr;
  }
  else
  {
    free( dst->cplx );
    dst->cplx = cptr;
  }
  dst->ncols = ncols;

  return TRUE;
}


BOOL MTX_InsertColumn( MTX *dst, const MTX *src, const unsigned dst_col, const unsigned src_col )
{
  unsigned i = 0;
  unsigned j = 0;
  unsigned k = 0;
  unsigned m = 0;
  double **dptr = NULL;
  stComplex **cptr = NULL;

  if( MTX_isNull( dst ) )
  {
    MTX_ERROR_MSG( "NULL Matrix" );
    return FALSE;
  }

  if( MTX_isNull( src ) )
  {
    MTX_ERROR_MSG( "NULL Matrix" );
    return FALSE;
  }

  if( dst->nrows != src->nrows )
  {
    MTX_ERROR_MSG( "if( dst->nrows != src->nrows )" );
    return FALSE;
  }

  // note the missing '=' here, a column can be inserted at the end (i.e. AddColumn )
  if( dst_col > dst->ncols )
  {
    MTX_ERROR_MSG( "if( dst_col > dst->ncols )" );
    return FALSE;
  }

  if( src_col >= src->ncols )
  {
    MTX_ERROR_MSG( "if( src_col >= src->ncols )" );
    return FALSE;
  }

  if( !src->isReal && dst->isReal )
  {
    // convert the destination matrix to complex
    if( !MTX_ConvertRealToComplex( dst ) )
    {
      MTX_ERROR_MSG( "MTX_ConvertRealToComplex returned FALSE." );
      return FALSE;
    }
  }

  // allocate a new array of column vectors
  if( dst->isReal )
  {
    dptr = (double**)malloc( (dst->ncols+1)*sizeof(double*) );
    if( !dptr )
    {
      MTX_ERROR_MSG( "malloc returned NULL." );
      return FALSE;
    }
  }
  else
  {
    cptr = (stComplex**)malloc( (dst->ncols+1)*sizeof(stComplex*) );
    if( !cptr )
    {
      MTX_ERROR_MSG( "malloc returned NULL." );
      return FALSE;
    }
  }

  // copy the previous array of pointers
  // and add the new column
  k = 0;
  for( j = 0; j <= dst->ncols; j++ )
  {
    if( j == dst_col )
    {
      // allocate a new column vector
      if( dst->isReal )
        dptr[k] = (double*)malloc( dst->nrows*sizeof(double) );
      else
        cptr[k] = (stComplex*)malloc( dst->nrows*sizeof(stComplex) );
      if( dst->isReal )
      {
        if( !dptr[k] )
        {
          // this is most likely to occur if allocating more memory than available
          for( m = 0; m < k; m++ )
          {
            free( dptr[m] );
          }

          MTX_ERROR_MSG( "malloc returned NULL." );
          free( dptr );
          return FALSE;
        }
      }
      else
      {
        if( !cptr[k] )
        {
          // this is most likely to occur if allocating more memory than available
          for( m = 0; m < k; m++ )
          {
            free( cptr[m] );
          }
          MTX_ERROR_MSG( "malloc returned NULL." );
          free( cptr );
          return FALSE;
        }
      }

      // copy the src column vector
      for( i = 0; i < dst->nrows; i++ )
      {
        if( dst->isReal )
        {
          dptr[k][i] = src->data[src_col][i];
        }
        else
        {
          if( src->isReal )
          {
            cptr[k][i].re = src->data[src_col][i];
            cptr[k][i].im = 0;
          }
          else
          {
            cptr[k][i] = src->cplx[src_col][i];
          }
        }
      }
      // copy the data that was at the insertion index
      // unless this is the after the last column
      if( j != dst->ncols )
      {
        k++;
        if( dst->isReal )
          dptr[k] = dst->data[j];
        else
          cptr[k] = dst->cplx[j];
      }
    }
    else
    {
      if( j != dst->ncols )
      {
        if( dst->isReal )
          dptr[k] = dst->data[j];
        else
          cptr[k] = dst->cplx[j];
      }
    }
    k++;
  }

  // free the old column array, and copy the new
  if( dst->isReal )
  {
    free( dst->data );
    dst->data = dptr;
  }
  else
  {
    free( dst->cplx );
    dst->cplx = cptr;
  }
  dst->ncols++;

  return TRUE;
}

BOOL MTX_AddColumn( MTX *dst, const MTX *src, const unsigned src_col )
{
  return MTX_InsertColumn( dst, src, dst->ncols, src_col );
}

BOOL MTX_Concatonate( MTX *dst, const MTX *src )
{
  unsigned i = 0;
  unsigned j = 0;
  unsigned ncols;
  unsigned m = 0;
  double **dptr = NULL;
  stComplex **cptr = NULL;

  if( dst == NULL )
  {
    MTX_ERROR_MSG( "dst is a NULL Matrix" );
    return FALSE;
  }

  if( MTX_isNull( src ) )
  {
    MTX_ERROR_MSG( "NULL Matrix" );
    return FALSE;
  }

  if( dst->nrows == 0 && dst->ncols == 0 )
  {
    return MTX_Copy( src, dst );
  }
  else if( dst->nrows != src->nrows )
  {
    MTX_ERROR_MSG( "if( dst->nrows != src->nrows )" );
    return FALSE;
  }

  ncols = dst->ncols+src->ncols;

  if( dst->isReal && !src->isReal )
  {
    // Convert dst to complex
    if( !MTX_ConvertRealToComplex( dst ) )
    {
      MTX_ERROR_MSG( "MTX_ConvertRealToComplex returned FALSE." );
      return FALSE;
    }
  }

  // allocate a new array of column vectors
  if( dst->isReal )
  {
    dptr = (double**)malloc( ncols*sizeof(double*) );
    if( !dptr )
    {