/////////////////////////////////////////////////////////////////////////////////////////////////
// 
// FilterStd.cc: 2-dimensional grid simulation program
// to support HighPerformance Computing benchmark study.
// Copyright 2006 National Center for Ecological Analysis 
// and Synthesis (NCEAS) 
// Author: Rick Reeves, NCEAS Scientific Programmer
//
// Starting Point: this program defines a 2-dimensional grid,
// initializes the grid to a set of random numbers, and then
// replaces each random value with the integer average of the 
// cell and its surrounding FILTER_DIM * FILTER_DIM  neighbors. 
// The program employs matrix 'wrapping' logic to compute the
// average for cells along the outer edges of the matrix. 
// the program also uses a Clock object to generate random number 
// series based on system time and also to time the program execution.
//
// This program's results best used in conjunction with the 
// companion program SimHpc.cc, which solves the same problem
// using multiple processors and shared memory via the MPI
// applications programming interface and UNIX Shared Memory
// programming, respectively. 
//
/////////////////////////////////////////////////////////////////////////////////////////////////
//
#include <stdio.h>
#include <stdlib.h>
#include "Clock.h" 

#define GRID_DIM   1000
#define FILTER_DIM    9 
#define FILTER_SIZE (FILTER_DIM * FILTER_DIM) 

#define MAX_FILTER_ITERATIONS 100

#define MAX_PRINT_GRID_ROW 25
#define MAX_PRINT_GRID_COL 25

#define TRUE  1
#define FALSE 0

#undef  TIME_INPUT_MATRIX_GEN
#define TIME_INPUT_MATRIX_GEN 1

#undef PRINT_MATRIX 
#define PRINT_MATRIX 1

void PrintGrid(int (*matptr)[GRID_DIM]);

int iGrid[GRID_DIM][GRID_DIM];
int iGrid2[GRID_DIM][GRID_DIM];

int main(int argc,char *argv[])
{
   int iCtr,jCtr,kCtr,iCtr2,jCtr2,iPtr,jPtr,iSum,iHalfFilterDim;
   unsigned long lElapsedTimeSec,lElapsedTimeMin; 
   float fAvg;
   Clock myClock;

   char sMsgBuf[180];
//
// get the aproximate execution start time from the O/S
//
#ifdef TIME_INPUT_MATRIX_GEN
   myClock.GetStartingTime();
#endif
//
// Initialize the matrix to random number series seeded and generated 
// by the Clock objects. Note: the Clock object automatically maintains
// a vector cache of random numbes, replenishing it with a new series 
// when necessary.
//
   printf("populating the  input matrix...\n");
   for (iCtr = 0; iCtr < GRID_DIM; iCtr++)
   { 
      if (!(iCtr % 100))
         printf("generating input matrix row %d...\n",iCtr);
      for (jCtr = 0; jCtr < GRID_DIM; jCtr++)
      {
         iGrid[iCtr][jCtr] = int (GRID_DIM * (myClock.GetNextRandValFromCache()/(RAND_MAX + 1.0))); 
      }
   }
//
// for debug purposes, print the randomly generated 
// input matrix.
// Note that the results of the filtering process
// will be written to a new matrix, so that the input
// matrix is not corrupted by the filtering..
//
#ifdef PRINT_MATRIX
printf("printing INPUT matrix...\n");
   PrintGrid(iGrid);          
#endif
#ifdef TIME_INPUT_MATRIX_GEN
//
// Elapsed time calculation: uncomment this to see how long random # matrix takes to generate
// 
   myClock.GetCurrentTime();
   lElapsedTimeSec = myClock.GetElapsedTimeInSec();
   lElapsedTimeSec = lElapsedTimeSec;
   lElapsedTimeMin = lElapsedTimeSec / 60;
   if (lElapsedTimeMin > 0)
   lElapsedTimeSec = lElapsedTimeSec % lElapsedTimeMin;

   sprintf(sMsgBuf,"(Std Version) elapsed time to generate input matrix: %ld minutes %ld seconds.\n",
           lElapsedTimeMin,lElapsedTimeSec);
#endif
//
// Now, do some computing across the grid. For the demonstration, run a simple
// low-pass filter over the grid, replacing each cell with the average of 
// the FILTER_DIM * FILTER_DIM elements surrounding each cell. We selected
// this filter because it is straightforward to understand, easily coded, 
// and for large grid or filter kernal sizes, can consume lots of CPU-seconds.
//
// Scientists can substitute their own grid processing function here.
// if you timed the random number generation, uncomment the following line
// get a second execution start time to clock the second half.
// Reset the system clock to a new starting time for the grid calculations
//
   myClock.GetStartingTime();
//
// Note: to increase the consumption of execution time, repeat the filter operation
// multiple times (MAX_FILTER_ITERATIONS, to be exact).
//
// Filter is coded inline (rather than wrapped inside a function) to minimize overhead.
//
   iHalfFilterDim = ((FILTER_DIM - 1) / 2);
   for (kCtr = 0; kCtr < MAX_FILTER_ITERATIONS; kCtr++)
   {
      if (!(kCtr % 100))
        printf("generating FILTERED matrix iteration %d..\n",kCtr);
      for (iCtr = 0; iCtr < GRID_DIM; iCtr++)
      {
         for (jCtr = 0; jCtr < GRID_DIM; jCtr++)
         {
//
// the convolution is done for each cell: each cell replaced by
// the average of its FILTER_DIM x FILTER_DIM neighborhood.
//
            iSum = 0;
            fAvg = 0.0;    
//
// iPtr is the ULC of the local neighborhood
// 
            iPtr = iCtr - iHalfFilterDim;
            if (iPtr < 0)
               iPtr = GRID_DIM - iHalfFilterDim;
            for (iCtr2 = 1; iCtr2 <= FILTER_DIM; iCtr2++)
            {
               jPtr = jCtr - iHalfFilterDim;
               if (jPtr < 0)
                  jPtr = GRID_DIM - iHalfFilterDim;
               for (jCtr2 = 1; jCtr2 <= FILTER_DIM; jCtr2++)
               {
                  iSum += iGrid[iPtr][jPtr];
                  jPtr++;
                  if (jPtr == GRID_DIM)
                     jPtr = 0;
               }
               iPtr++;
               if (iPtr == GRID_DIM)
                 iPtr = 0;
            }
            iGrid2[iCtr][jCtr] = int(iSum / (float)FILTER_SIZE);
         }
      }
   }
#ifdef PRINT_MATRIX
//
// print the output matrix
//
printf("printing OUTPUT matrix...\n");
   PrintGrid(iGrid2);          
#endif
//
// elapsed time calculation
// 
   myClock.GetCurrentTime();
   lElapsedTimeSec = myClock.GetElapsedTimeInSec();
   lElapsedTimeMin = lElapsedTimeSec / 60;
   if (lElapsedTimeMin > 0)
   lElapsedTimeSec = lElapsedTimeSec % lElapsedTimeMin;

   printf("(Std Version) Program is done.\n");
//
#ifdef TIME_INPUT_MATRIX_GEN
    printf("%s\n",sMsgBuf);
#endif   
    printf("Elapsed time to generate ouptut matrix is %ld minutes %ld seconds.\n",
           lElapsedTimeMin,lElapsedTimeSec);

    return(0);
}
/////////////////////////////////////////////////////////////////////////////////////////////////
// 
// PrintGrid() function prints a 2-dimensional grid to the terminal screen.
// Copyright 2006 National Center for Ecological Analysis 
// and Synthesis (NCEAS) 
// Author: Rick Reeves, NCEAS Scientific Programmer
//
// Input argument: a pointer (scaled to a one-dimenisional array) to the 
// start of the desired grid.
// Note: sometimes the grid will be larger than can be displayed on the 
// screen, so we print out only the upper left-hand MAX_PRINT_GRID_ROW(COL)
// portion of the matrix, or the entire matrix, whichever is smaller.
//
/////////////////////////////////////////////////////////////////////////////////////////////////
//
void PrintGrid(int (*matptr)[GRID_DIM])
{
   int iCtr, jCtr,iPrintRowLim,iPrintColLim;
//
// set the upper bound on grid output.
//
   if (GRID_DIM > MAX_PRINT_GRID_ROW)
      iPrintRowLim = MAX_PRINT_GRID_ROW;
   else
      iPrintRowLim = GRID_DIM;
   if (GRID_DIM > MAX_PRINT_GRID_COL)
      iPrintColLim = MAX_PRINT_GRID_COL;
   else
      iPrintColLim = GRID_DIM;

   fprintf(stdout,"PrintGrid printing the matrix:\n");
   for (iCtr = 0; iCtr < iPrintRowLim; iCtr++)
   {
      for (jCtr = 0; jCtr < iPrintColLim; jCtr++)
      {
        fprintf(stdout,"%3d ",matptr[iCtr][jCtr]);
      }
      fprintf(stdout,"\n");
   }
   fprintf(stdout,"PrintGrid finished printing the matrix\n");
}