//$$ submat.cxx                         submatrices

// Copyright (C) 1991: R B Davies and DSIR

#include "include.hxx"

#include "newmat.hxx"
#include "newmatrc.hxx"


//#define REPORT { static ExeCounter ExeCount(__LINE__,6); ExeCount++; }

#define REPORT {}


/****************************** submatrices *********************************/

GetSubMatrix BaseMatrix::SubMatrix(int first_row, int last_row, int first_col,
   int last_col)
{
   REPORT
   int a = first_row - 1; int b = last_row - first_row + 1;
   int c = first_col - 1; int d = last_col - first_col + 1;
   if (a<0 || b<=0 || c<0 || d<=0) MatrixError("SubMatrix dimension error");
   return GetSubMatrix(this, a, b, c, d, Type().sub());
}

GetSubMatrix BaseMatrix::SymSubMatrix(int first_row, int last_row)
{
   REPORT
   int a = first_row - 1; int b = last_row - first_row + 1;
   if (a<0 || b<=0) MatrixError("SubMatrix dimension error");
   return GetSubMatrix( this, a, b, a, b, Type().ssub());
}

GetSubMatrix BaseMatrix::Row(int first_row)
{
   REPORT
   int a = first_row - 1;
   if (a<0) MatrixError("SubMatrix dimension error");
   return GetSubMatrix(this, a, 1, 0, NcolsV(), MatrixType::RowV);
}

GetSubMatrix BaseMatrix::Rows(int first_row, int last_row)
{
   REPORT
   int a = first_row - 1; int b = last_row - first_row + 1;
   if (a<0 || b<=0) MatrixError("SubMatrix dimension error");
   return GetSubMatrix(this, a, b, 0, NcolsV(), Type().sub());
}

GetSubMatrix BaseMatrix::Column(int first_col)
{
   REPORT
   int c = first_col - 1;
   if (c<0) MatrixError("SubMatrix dimension error");
   return GetSubMatrix(this, 0, NrowsV(), c, 1, MatrixType::ColV);
}

GetSubMatrix BaseMatrix::Columns(int first_col, int last_col)
{
   REPORT
   int c = first_col - 1; int d = last_col - first_col + 1;
   if (c<0 || d<=0) MatrixError("SubMatrix dimension error");
   return GetSubMatrix(this, 0, NrowsV(), c, d, Type().sub());
}

void GetSubMatrix::operator<<(BaseMatrix& bmx)
{
   REPORT
   GeneralMatrix* gmx = bmx.Evaluate();
   GeneralMatrix* gm = bm->Evaluate();
   if ((BaseMatrix*)gm!=bm) MatrixError("Illegal argument for SubMatrix");
   if (row_number != gmx->Nrows() || col_number != gmx->Ncols())
      MatrixError("Dimensions don't match");
   if (row_skip+row_number > gm->Nrows() || col_skip+col_number > gm->Ncols())
      MatrixError("SubMatrix dimension error");
   int i = row_number;
   MatrixRow mrx(gmx, LoadOnEntry); 
   MatrixRow mr(gm, LoadOnEntry+StoreOnExit+DirectPart, row_skip);
                                  // do need LoadOnEntry
   MatrixRowCol sub;
   while (i--)
   {
      mr.SubRowCol(sub, col_skip, col_number);   // put values in sub
      sub.Copy(mrx); mr.Next(); mrx.Next();
   }
   gmx->tDelete();
}   

void GetSubMatrix::operator<<(const real* r)
{
   REPORT
   GeneralMatrix* gm = bm->Evaluate();
   if ((BaseMatrix*)gm!=bm) MatrixError("Illegal argument for SubMatrix");
   if (row_skip+row_number > gm->Nrows() || col_skip+col_number > gm->Ncols())
      MatrixError("SubMatrix dimension error");
   int i = row_number;
   MatrixRow mr(gm, LoadOnEntry+StoreOnExit+DirectPart, row_skip);
                                  // do need LoadOnEntry
   MatrixRowCol sub;
   while (i--)
   {
      mr.SubRowCol(sub, col_skip, col_number);   // put values in sub
      sub.Copy(r); mr.Next();
   }
}   

void GetSubMatrix::operator<<(real r)
{
   REPORT
   GeneralMatrix* gm = bm->Evaluate();
   if ((BaseMatrix*)gm!=bm) MatrixError("Illegal argument for SubMatrix");
   if (row_skip+row_number > gm->Nrows() || col_skip+col_number > gm->Ncols())
      MatrixError("SubMatrix dimension error");
   int i = row_number;
   MatrixRow mr(gm, LoadOnEntry+StoreOnExit+DirectPart, row_skip);
                                  // do need LoadOnEntry
   MatrixRowCol sub;
   while (i--)
   {
      mr.SubRowCol(sub, col_skip, col_number);   // put values in sub
      sub.Copy(r); mr.Next();
   }
}   

void GetSubMatrix::Inject(const GeneralMatrix& gmx)
{
   REPORT
   GeneralMatrix* gm = bm->Evaluate();
   if ((BaseMatrix*)gm!=bm) MatrixError("Illegal argument for SubMatrix");
   if (row_number != gmx.Nrows() || col_number != gmx.Ncols())
      MatrixError("Dimensions don't match");
   if (row_skip+row_number > gm->Nrows() || col_skip+col_number > gm->Ncols())
      MatrixError("SubMatrix dimension error");
   int i = row_number;
   MatrixRow mrx((GeneralMatrix*)(&gmx), LoadOnEntry);
   MatrixRow mr(gm, LoadOnEntry+StoreOnExit+DirectPart, row_skip);
                                  // do need LoadOnEntry
   MatrixRowCol sub;
   while (i--)
   {
      mr.SubRowCol(sub, col_skip, col_number);   // put values in sub
      sub.Inject(mrx); mr.Next(); mrx.Next();
   }
}

