/* rmatrix.cc */

#include "rmatrix.h"

void dc_rect_matrix::touch( void ) {
  valid = val.Nrows() && val.Ncols();
  T_valid = I_valid = false;
}

dc_rect_matrix::dc_rect_matrix() {
  valid = T_valid = I_valid = false;
}

dc_rect_matrix::dc_rect_matrix( const int m, const int n, const string &label, 
				dc_node *parent = nil) {
  resize( m, n );
  set_both( label, parent );
}

dc_rect_matrix::dc_rect_matrix( const Matrix &M, const string &label, 
				dc_node *parent = nil ) {
  set( M );
  set_both( label, parent );
}

dc_rect_matrix::dc_rect_matrix( const string &label, dc_node *parent = nil ) {
  valid = T_valid = I_valid = false;
  set_both( label, parent );
}

dc_rect_matrix::~dc_rect_matrix() {
  val.CleanUp();
  T.CleanUp();
  I.CleanUp();
}

dc_rect_matrix &dc_rect_matrix::operator=( const Matrix &M ) {
  val = M;
  touch();
  return *this;
}

dc_rect_matrix &dc_rect_matrix::operator=( dc_data &D ) {
  if( D.sub_type() != Rect_Matrix_t ) {
    cerr << "dc_rect_matrix::operator= -- Cannot assign " << D.full_type()
	 << " to " << full_type() << "\n";
    exit( 1 );
  }
  val = ( ( dc_rect_matrix * )&D )->get();
  touch();
  return *this;
}

ostream &dc_rect_matrix::display_c( ostream &stream = cout ) const {
  if( !( is_valid() ) ) {
    return stream << full_type() << "\n";
  } else {
    int n = Nrows();
    int m = Ncols();
    int k = 0;
    Real *d = get_store();
    stream << "( " << n << ", " << m << " ) { ";
    for( int j = 0 ; j < n - 1 ; j++ ) {
      for( int i = 0 ; i < m ; i++ ) {
	stream << d[k++] << ", ";
      }
      stream << "\n  ";
    }
    for( int i = 0 ; i < m - 1 ; i++ ) {
      stream << d[k++] << ", ";
    }
    return stream << d[k] << " }" << units;
  }
}

istream &operator>>( istream &stream, dc_rect_matrix &M ) {
  if( !M.is_valid() ) {
    cerr << "dc_rect_matrix::operator>> -- invalid operation on " 
	 << M.full_type() << "\n";
    exit(1);
  }
  int area = ( ( Matrix )M.get() ).Storage();
  double *n = ( ( Matrix )M.get() ).Store();
  for( int i = 0 ; i < area ; i++ ){
     stream >> n[i];
  }
  
  M.set( n );
  return stream;
}

Real dc_rect_matrix::operator()( const int i, const int j ) const {
  return val( i, j );
}

ReturnMatrix dc_rect_matrix::get( void ) {
  if( !valid ) {
    cerr << "dc_rect_matrix::get -- invalid operation on " << full_type() 
	 << "\n";
    exit(1);
  }

  return val;
}

ReturnMatrix dc_rect_matrix::get_T( void ) {
  if( !valid ) {
    cerr << "dc_rect_matrix::get_T -- invalid operation on " << full_type() 
	 << "\n";
    exit(1);
  }

  if( !T_valid ) {
    T_valid = true;
    T = val.t(); 
  }
  
  return T.ForReturn();
}

ReturnMatrix dc_rect_matrix::get_I( void ) {
  if( !valid ) {
    cerr << "dc_rect_matrix::get_I -- invalid operation on " << full_type() 
	 << "\n";
    exit(1);
  }

  if( !I_valid ) {
    I_valid = true;
    I = val.i(); 
  }
  
  return I.ForReturn();
}

void dc_rect_matrix::get_dim( int &h, int &w ) const {
  if( !valid ) {
    cerr << "dc_rect_matrix::get_dim -- invalid operation on " << full_type() 
	 << "\n";
    exit(1);
  }
  h = val.Nrows();
  w = val.Ncols();
}

void dc_rect_matrix::set( const Matrix &M ) {
  val = M;
  touch();
}

void dc_rect_matrix::set( dc_rect_matrix &D ) {
  val = D.get();
  touch();
}

void dc_rect_matrix::set( double n[] ) {
  val << n;
  touch();
}

void dc_rect_matrix::resize( const int h, const int w ) {
  val.ReDimension( h, w );
  touch();
}

dc_label *dc_rect_matrix::duplicate( dc_node *parent ) const {
  dc_rect_matrix *dupe = new dc_rect_matrix( val );
  dupe->set_units( get_units_c() );

  if( !is_temporary() ) {
    if( dupe->set_both( local_label(), parent ) ) {
      delete( dupe );
      return nil;
    }
  }
  dupe->set_visible( is_visible() );
  return ( dc_label * )dupe;
}

bool operator==( Matrix const &M1, Matrix const &M2 ) {
  if( ( M1.Ncols() != M2.Ncols() ) || ( M1.Nrows() != M2.Nrows() ) ) {
    cerr << "operator== -- cannot compare " << M1.Nrows() << "x" << M1.Ncols()
	 << " matrix with " << M2.Nrows() << "x" << M2.Ncols() << " matrix\n"; 
    exit( 1 );
  }

  double *d1 = M1.Store(); 
  double *d2 = M2.Store();
  for( int a = M1.Storage() - 1 ; a >= 0 ; a -- ) {
    if( d1[a] != d2[a] )
      return false;
  }
  return true;
}

void dc_matrix::get_xtype( xtype &xt ) { 
  xt.T = type();
  xt.subT = sub_type();
  xt.uv = units;
  xt.inf.dim[0] = Nrows();
  xt.inf.dim[1] = Ncols();
}
