#include "Grid.h"

Grid::Grid()
{
  latitude = NULL;
}

Grid::Grid(const string &inGridName)
{
  latitude = NULL;
  setGrid(inGridName);
}

Grid::~Grid(){
  if (latitude){
    delete [] latitude;
    latitude = NULL;
  }
}

/********************************************************************//**
 * Setup the MIX grid
 *
 * \param[in] gridName:  String describing the grid.
 *
 * FIXME: If "gridName" is not understood, we should read that grid
 * from a file.
 * 
 * Example: 
 *    >>> setGrid("2x2_44")
 *    latitude = [ 0.0, 0.1, 0.25, 0.5, 1.0, 2, 4, 6, 8, ..., 42, 44 ]
 * ************************************************************************/
void Grid::setGrid(const string &inGridName)
{
  gridName = inGridName;
  int lowLatBoundary = 0;
  
  if (gridName == "1x1_44"){
    // 1 degree lat/lon spacing, down to colatitude of 44 degrees.
    gridSpacing = 1;
    lowLatBoundary = 44;
    nLat = 48;
    
    forceSymmetricPotential = false;
  }
  else if (gridName == "1x1_60"){
    // 1 degree lat/lon spacing, down to colatitude of 60 degrees.
    gridSpacing = 1;
    lowLatBoundary = 60;
    nLat = 64;

    forceSymmetricPotential = true;
  }
  else if (gridName == string("2x2_44")){
    // 2 degree lat/lon spacing, down to colatitude of 44 degrees
    gridSpacing = 2;
    lowLatBoundary = 44;
    nLat = 27;
 
    forceSymmetricPotential = false;
  }
  else if (gridName == "2x2_60"){
    // 2 degree lat/lon spacing, down to a colatitude of 60 degrees.
    gridSpacing = 2;
    lowLatBoundary = 60;
    nLat = 35;

    forceSymmetricPotential = true;    
  }
  else if (gridName == "2x2_30"){
    // 2 degree lat/lon spacing, down to a colatitude of 30 degrees.
    gridSpacing = 2;
    lowLatBoundary = 30;
    nLat = 20;

    forceSymmetricPotential = false;    
  }
  else if (gridName == "1x1_30"){
    // 1 degree lat/lon spacing, down to a colatitude of 30 degrees.
    gridSpacing = 1;
    lowLatBoundary = 30;
    nLat = 34;

    forceSymmetricPotential = false;    
  }
  else {
    //FIXME: If "gridName" is not understood, we should read that grid from a file.
    cerr << "ERROR:  Did not understand grid type: \"" << gridName << "\"." << endl;
    exit(1);
  }

  // Now that we have some basic grid information, define the latitude array:
  if (latitude != NULL){
    delete [] latitude;
    latitude = NULL;
  }
  latitude = new double[nLat+1];
  
  // Dense grid near the pole
  latitude[0] = 0.0;
  latitude[1] = 0.1;
  latitude[2] = 0.25;
  latitude[3] = 0.5;
  latitude[4] = 1.0;
  // Set evenly spaced grid beyond 1 degree colatitude:
  fillLatitudes(5, 2, lowLatBoundary, gridSpacing);
}

/********************************************************************//**
 * Fill in the latitude array with constant spacing up to stopLat.
 *
 * \param[in] startIndex:  First index that should be filled in  
 * \param[in] startLat:  Latitude at first index
 * \param[in] stopLat:  fill in latitude array until reaching stopLat
 * \param[in] spacing:  constant spacing between startLat & stopLat
 *
 * Example: 
 *    >>> fillLatitudes(5,2,12,2)
 *    latitude = [ 0.0, 0.1, 0.25, 0.5, 1.0, 2, 4, 6, 8, 10, 12 ]
 * ************************************************************************/
void Grid::fillLatitudes(const int &startIndex, const double &startLat, const double &stopLat, const double &spacing)
{
  int idx = startIndex;
  do{
    latitude[idx] = startLat + double(idx-startIndex) * spacing;
    idx++;
  }
  while ( latitude[idx-1] < stopLat-1e-6 );
}