#include <stdio.h>
#include <string.h>
#include <stdlib.h>
#include <stdbool.h>
#include <unistd.h>
#include <math.h>

int  timer_new(const char *text);
void timer_report(void);
void timer_start(int it);
void timer_stop(int it);

//#undef _GET_IBM_COUNTER

#ifdef _GET_IBM_COUNTER
#include <libhpc.h>
#endif

#include "cgribex.h"

#define  NINT(x)  ((int)(lround((x))))

#define  NLEV    1 // Number of levels    
#define  WITH_IO 1

static
void cgribexDefaultSec0(int *isec0)
{
  ISEC0_GRIB_Len     = 0;
  ISEC0_GRIB_Version = 0;
}

static
void cgribexDefaultSec1(int *isec1)
{
  ISEC1_CenterID    = 0;
  ISEC1_SubCenterID = 0;
  ISEC1_LocalFLag   = 0;
}

static
void cgribexDefaultSec4(int *isec4)
{
  for ( long i = 2; i <= 10; ++i ) isec4[i] = 0;
}

static
void cgribexDefGridLonLat(int *isec1, int *isec2, int *isec4, int nlon, int nlat, double di)
{
  (void)isec4;
  memset(isec2, 0, 16*sizeof(int));

  double xfirst = 0;
  double xlast  = 360. - di/2;
  double yfirst = -90. + di/2;
  double ylast  =  90. - di/2;

  ISEC1_Sec2Or3Flag = 128;
  ISEC1_GridDefinition = 255;

  ISEC2_GridType = GRIB1_GTYPE_LATLON;
  ISEC2_NumLon   = nlon;
  ISEC2_NumLat   = nlat;
  ISEC2_FirstLat = NINT(yfirst*1000);
  ISEC2_LastLat  = NINT(ylast*1000);
  ISEC2_LatIncr  = NINT(di*1000);
  ISEC2_FirstLon = NINT(xfirst*1000);
  ISEC2_LastLon  = NINT(xlast*1000);
  ISEC2_LonIncr  = NINT(di*1000);

  ISEC2_Reduced  = 0;
  ISEC2_ResFlag  = 128;
  ISEC2_ScanFlag = 128;
  ISEC2_ScanFlag += 64;
}

static
void cgribexDefLevelSurface(int *isec1, int *isec2, double *fsec2)
{
  (void)fsec2;
  ISEC2_NumVCP    = 0;
  ISEC1_LevelType = GRIB1_LTYPE_SURFACE;
  ISEC1_Level1    = 0;
  ISEC1_Level2    = 0;
}

static
void cgribexDefTime(int *isec1, int time)
{
  int year = 2000, month = 1, day = 1;

  int hour = time/10000;
  int minute = (time - hour*10000) / 100;

  int century =  year / 100;

  ISEC1_Year = year - century*100;
  
  if ( year < 0 )
    {
      century = -century;
      ISEC1_Year = -ISEC1_Year;
    }

  if ( ISEC1_Year == 0 )
    {
      century -= 1;
      ISEC1_Year = 100;
    }

  century += 1;
  if ( year < 0 ) century = -century;

  ISEC1_Month  = month;
  ISEC1_Day    = day;
  ISEC1_Hour   = hour;
  ISEC1_Minute = minute;

  ISEC1_Century = century;

  ISEC1_TimeUnit = ISEC1_TABLE4_MINUTE;
  ISEC1_TimeRange   = 10;

  ISEC1_TimePeriod1 = 0;
  ISEC1_TimePeriod2 = 0;
  ISEC1_AvgNum         = 0;
  ISEC1_AvgMiss        = 0;
  ISEC1_DecScaleFactor = 0;
}

void usage(void)
{
  fprintf(stderr,"\ngrib_encode [-b <bits per sample> -d <grid direction increment> -n <num recs> -z -t] <filename>\n\n");
  exit(-1);
}

int main(int argc, char *argv[])
{
  int isec0[2], isec1[4096], isec2[4096], isec3[2], isec4[512];
  double fsec2[512], fsec3[2];
  int nbits = 16;
  int c, iret, iword;
  int fileID;
  int nts = 100;
  double di = 0.1; // grid direction increment
  bool lzip = false;
  bool gribTimer = false;

  while ((c = getopt(argc, argv, "b:d:n:azt")) != -1) 
    {
      switch (c) 
	{
	case 'a':
	case 'z':
          lzip = true;
	  break;
	case 'b':
	  nbits = atoi(optarg);
	  break;
	case 'd':
	  di = atof(optarg);
	  break;
	case 'n':
	  nts = atoi(optarg);
	  break;
        case 't':
          gribTimer = true;
          break;
	case '?':
	default:
	  usage();
	}
    }
  argc -= optind;
  argv += optind;

  if (nbits < 1) usage(); 
  if ((argc == 0) | (argc > 1)) usage();

  char *filename = (char *) malloc(sizeof(char)*(strlen(argv[0])+1));
  strcpy(filename, argv[0]);

#ifdef _GET_IBM_COUNTER
  hpmInit(0,"main");
#endif

  int nlat = (int) (180./di);
  int nlon = nlat*2;
  size_t gridsize = (size_t)nlon*(size_t)nlat;
  
  double *data = (double *) malloc(gridsize*sizeof(double));
  for ( size_t i = 0; i < gridsize; ++i )
    {
      // data[i] = 273. + 50*(rand()/(RAND_MAX+1.0));
      int ilat = (int)(i%gridsize)/ nlon;
      int ilon = (int)i%nlon;
      data[i] = cos(2.0 * M_PI * (0          + ilon*di)/360)
	      * sin(2.0 * M_PI * (-90 + di/2 + ilat*di)/180);
    }

  memset(isec1, 0, 32*sizeof(int));
  fsec2[0] = 0; fsec2[1] = 0;

  cgribexDefaultSec0(isec0);
  cgribexDefaultSec1(isec1);
  cgribexDefaultSec4(isec4);

  ISEC1_CodeTable = 201;
  ISEC1_Parameter = 1;

  cgribexDefGridLonLat(isec1, isec2, isec4, nlon, nlat, di);
  cgribexDefLevelSurface(isec1, isec2, fsec2);

  if ( nbits > 0 ) isec4[1] = nbits;

  ISEC4_NumValues = (int)gridsize;
  ISEC4_NumBits   = nbits;

  size_t gribbuffersize = (size_t)isec4[0]*(size_t)isec4[1]/8 + (size_t)isec4[0]/8 + 3000; // datasize + mask + header
  long gribsize = (long)(gribbuffersize / sizeof(int));
  unsigned char *gribbuffer = (unsigned char *) malloc(gribbuffersize);
  size_t buffersize = 0;
  unsigned char *buffer = NULL;

  int timer_total, timer_write, timer_encode, timer_szip;
  if ( gribTimer )
    {
      timer_total      = timer_new("total");
      timer_write      = timer_new("write");
      timer_encode     = timer_new("encode");
      timer_szip       = timer_new("szip");
    }

  if ( gribTimer ) timer_start(timer_total);

  if ( WITH_IO ) fileID = gribOpen(filename, "w");

  for ( int i = 0; i < nts; ++i )
    {
      cgribexDefTime(isec1, 120000+i*100);

      if ( gribTimer ) timer_start(timer_encode);
      gribExDP(isec0, isec1, isec2, fsec2, isec3, fsec3, isec4, (double*) data,
	       (int)gridsize, (int *) gribbuffer, (int)gribsize, &iword, "C", &iret);
      if ( gribTimer ) timer_stop(timer_encode);

      if ( iret )
	{
	  fprintf(stderr, "Problem during GRIB encode (errno = %d)!\n", iret);
	  break;
	}

      size_t nbytes = (size_t)iword*sizeof(int);

      if ( lzip )
        {
          if ( buffersize < gribbuffersize )
            {
              buffersize = gribbuffersize;
              buffer = (unsigned char *) realloc(buffer, buffersize);
            }

          if ( gribTimer ) timer_start(timer_szip);
          nbytes = (size_t)gribZip(gribbuffer, (long)gribbuffersize, buffer, (long)buffersize);
          if ( gribTimer ) timer_stop(timer_szip);
        }
 
      if ( gribTimer ) timer_start(timer_write);
      if ( WITH_IO ) (void) gribWrite(fileID, gribbuffer, nbytes);
      if ( gribTimer ) timer_stop(timer_write);
    }

  if ( WITH_IO ) gribClose(fileID);

  free(gribbuffer);
  free(data);
  free(filename);

#ifdef _GET_IBM_COUNTER
  hpmTerminate(0);
#endif

  if ( gribTimer ) timer_stop(timer_total);
  if ( gribTimer ) timer_report();

  return 0;
}
