// tsncom05.cc
// Nick Macks (macks@fnal.gov), Summer 2002

// CovMatrices: CovMatrixX subtraction

#include "CovMatrices/CovMatrix5.h"
#include "CovMatrices/CovMatrix6.h"

using namespace std;

using CovMatrices::CovMatrix5;
using CovMatrices::CovMatrix6;

// Function declarations
void Identification();
void GetDimension(int&);
void Describe(const int&);
void BuildMatrix(CovMatrix5&, CovMatrix6&, const int&);
void BuildMatrix2(CovMatrix5&, CovMatrix6&, const int&);
void Show(const CovMatrix5&, const CovMatrix6&, const int&);
void GetCalcs(int&);
void GetRuns(int&);
void Calculate(CovMatrix5&, CovMatrix5&, CovMatrix6&, CovMatrix6&,
               const int&, const int&, const int&);


// Package ID
void Identification()
{
  cout << "Using CovMatrices\n";
}

// Ask for the desired dimension
void GetDimension(int& dimension)
{
  cout << "Enter the desired dimension(5-6): ";
  cin >> dimension;

  while ( (dimension < 5) || (dimension > 6) )
  {
    cout << "Enter 5 or 6: ";
    cin >> dimension;
  }
}

// Describe how the matrices were defined
void Describe(const int& dimension)
{
  cout << "Matrices were defined 'CovMatrix" << dimension << "'\n";
}

// Build the matrix (large values)
void BuildMatrix(CovMatrix5& matrix5, CovMatrix6& matrix6, const int& dimension)
{
  if (dimension == 5)
  {
    for (int i=0; i<dimension; ++i)
    {
      for (int j=0; j<dimension; ++j)
      {
        if (i == j)
          matrix5(i,j) = 1000000 + i*100000;
        else if (i > j)
          matrix5(i,j) = 2000 + i*2000 + j*(dimension*2000);
      }
    }
  }
  else if (dimension == 6)
  {
    for (int i=0; i<dimension; ++i)
    {
      for (int j=0; j<dimension; ++j)
      {
        if (i == j)
          matrix6(i,j) = 1000 + i*100;
        else if (i > j)
          matrix6(i,j) = 2 + i*2 + j*(dimension*2);
      }
    }
  }
}

// Build the matrix (small values)
void BuildMatrix2(CovMatrix5& matrix5, CovMatrix6& matrix6,
                  const int& dimension)
{
  if (dimension == 5)
  {
    for (int i=0; i<dimension; ++i)
    {
      for (int j=0; j<dimension; ++j)
      {
        if (i == j)
          matrix5(i,j) = 1000 + i*100;
        else if (i > j)
          matrix5(i,j) = 2 + i*2 + j*(dimension*2);
       }
    }
  }
  else if (dimension == 6)
  {
    for (int i=0; i<dimension; ++i)
    {
      for (int j=0; j<dimension; ++j)
      {
        if (i == j)
          matrix6(i,j) = 1000 + i*100;
        else if (i > j)
          matrix6(i,j) = 2 + i*2 + j*(dimension*2);
      }
    }
  }
}

// Show the matrix
void Show(CovMatrix5& matrix5, CovMatrix6& matrix6, int dimension)
{
  if (dimension == 5)
    cout << "Using matrix =\n" << matrix5 << "\n";
  else if (dimension == 6)
    cout << "Using matrix =\n" << matrix6 << "\n";
}

// Get the number of calculations
void GetCalcs(int& N)
{
  cout << "How many subtractions? ";
  cin >> N;
}

// Get the number of runs
void GetRuns(int& run)
{
  cout << "How many runs? ";
  cin >> run;
  cout << "\n";
}

// Subtract the matrices
void Calculate(CovMatrix5& matrix51, CovMatrix5& matrix52,
               CovMatrix6& matrix61, CovMatrix6& matrix62,
               const int& dimension, const int& N, const int& runs)
{
  // Define a counter
  int counter = 0;

  // Define two matrices to hold the result
  CovMatrix5 difference5;
  CovMatrix6 difference6;

  if (dimension == 5)
  {
    // Loop "run" times
    for (int r=0; r<runs; ++r)
    {
      // Do the subtractions
      for (int i=0; i<N; ++i)
      {
        difference5 = matrix51 - matrix52;
	matrix51 = difference5;
        ++counter;
      }

      // Write the output
      cout << counter << " subtractions took place\n";
      cout << "Difference =\n" << difference5 << "\n";

      // Reset the counter
      counter = 0;
    }
  }
  else if(dimension == 6)
  {
    // Loop "run" times
    for (int r=0; r<runs; ++r)
    {
      // Do the subtractions
      for (int i=0; i<N; ++i)
      {
        difference6 = matrix61 - matrix62;
	matrix61 = difference6;
        ++counter;
      }

      // Write the output
      cout << counter << " subtractions took place\n";  
      cout << "Difference =\n" << difference6 << "\n";

      // Reset the counter
      counter = 0;
    }
  }
}


int main()
{
  // What package are we using?
  Identification();

  // Get the dim(-ension)
  int dim = 0; GetDimension(dim);

  // Define two pairs of CovMatrices
  CovMatrix5 m51, m52;
  CovMatrix6 m61, m62;

  // Describe how the matrices were defined
  Describe (dim);

  // Build the matrices
  BuildMatrix(m51, m61, dim);
  BuildMatrix2(m52, m62, dim);

  // Show the actual matrices
  Show(m51, m61, dim);
  Show(m52, m62, dim);

  // Calculation parameters
  int tim = 0; GetCalcs(tim); // times to perform the calculation
  int run = 0; GetRuns(run);  // number of runs

  // Do the matrix subtractions
  Calculate(m51, m52, m61, m62, dim, tim, run);

  return 0;
}