#ifdef  HAVE_CONFIG_H
#include "../src/config.h"  /* The size of off_t depends on _FILE_OFFSET_BITS (on SUN) */
#endif

#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <math.h>
#include <limits.h>

#include "cgribex.h"
#include "namelist.h"

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

#define  UNDEFID  -1
#define  JP24SET    0xFFFFFF  // 2**24     (---> 16777215)

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

int gribTimer = 0;
int timer_total, timer_read, timer_write, timer_decode, timer_encode;
int timer_szip, timer_sunzip;

static char *Progname;

static
void version(void)
{
  fprintf(stderr, "%s version 1.7.0\n", Progname);
  fprintf(stderr, "  CGRIBEX library version : %s\n", cgribexLibraryVersion());
  fprintf(stderr, "\n");
/*
  1.0.0   6 Dec 2000 : initial version
  1.0.1  10 Jan 2003 : update to libgrib.a
  1.0.2  15 Jan 2003 : create output buffer for copy
  1.0.3   2 Feb 2003 : change macro MAX/MIN to if statement
  1.0.4  27 Jul 2003 : handle missing values
  1.1.0  22 Mar 2005 : change level handling
  1.1.1  19 Apr 2005 : change date/time format to ISO
  1.2.0  12 Feb 2006 : SZ compression with -szip and -sunzip
  1.3.0  12 Feb 2006 : add parameter select
  1.3.1  20 Mar 2006 : add select parameter - bits
  1.3.2  25 Aug 2006 : use ISEC1 only if ljdecode is TRUE (bug fix)
  1.3.3  11 Jan 2007 : add select parameter - fcnum
  1.3.4  25 Jun 2007 : replace decodeDate/decodeTime
  1.3.5  27 Dec 2007 : add option read
  1.3.6  03 Mar 2008 : add select parameter - ltype
  1.3.7  07 Jul 2010 : check1 - print 24 bit data with compression ratio 1 
  1.4.0   3 Nov 2011 : added option -p to print paramstr 
  1.5.0  27 Feb 2012 : added option -T  
  1.6.0  29 Jan 2013 : added option -aec/-unaec  
  1.7.0   8 Jul 2021 : fix timer bug  
*/
}

static
void usage(void)
{
  fprintf(stderr, "usage : %s [-options]  ifile  [ofile]\n", Progname);
  fprintf(stderr, "with:\n");

  fprintf(stderr, "   -bits <num>    pack with num bits\n");
  fprintf(stderr, "   -info          print GRIB contents\n");
  fprintf(stderr, "   -ginfo         print GRIB information\n");
  fprintf(stderr, "   -linfo         print long information (gribex section)\n");
  fprintf(stderr, "   -gdsinfo       print grid description section\n");
  fprintf(stderr, "   -pdsinfo       print product definition section\n");
  fprintf(stderr, "   -bmsinfo       print bit-map section\n");
  fprintf(stderr, "   -bdsinfo       print binary data section\n");

  fprintf(stderr, "   -copy          decode and encode GRIB records\n");
  fprintf(stderr, "   -read <num>    read file <num> times (performance test)\n");
  fprintf(stderr, "   -regular       regular grid (HOPER='R')\n");

  fprintf(stderr, "   -select        select records via namelist:\n");
  fprintf(stderr, "                    bits   INTEGER  List of bits numbers\n");
  fprintf(stderr, "                    code   INTEGER  List of code numbers\n");
  fprintf(stderr, "                    fcnum  INTEGER  List of forecast numbers\n");
  fprintf(stderr, "                    day    INTEGER  List of days\n");
  fprintf(stderr, "                    hour   INTEGER  List of hours\n");
  fprintf(stderr, "                    level  INTEGER  List of levels\n");
  fprintf(stderr, "                    ltype  INTEGER  List of level types\n");
  fprintf(stderr, "                    minute INTEGER  List of minutes\n");
  fprintf(stderr, "                    month  INTEGER  List of months\n");
  fprintf(stderr, "                    rec    INTEGER  Start and end record number\n");
  fprintf(stderr, "                    record INTEGER  List of records\n");
  fprintf(stderr, "                    rdate  INTEGER  List of reference dates [YYYYMMDD]\n");
  fprintf(stderr, "                    rtime  INTEGER  List of reference times [hhmmss]\n");
  fprintf(stderr, "                    table  INTEGER  List of table numbers\n");
  fprintf(stderr, "                    vdate  INTEGER  List of verification dates [YYYYMMDD]\n");
  fprintf(stderr, "                    vtime  INTEGER  List of verification times [hhmmss]\n");
  fprintf(stderr, "                    year   INTEGER  List of years\n");

  fprintf(stderr, "   -szip          szip compression\n");
  fprintf(stderr, "   -sunzip        szip decompression\n");
  fprintf(stderr, "   -aec           aec compression\n");
  fprintf(stderr, "   -unaec         aec decompression\n");

  fprintf(stderr, "   -float         single precision\n");

  fprintf(stderr, "   -debug         debugging mode\n");
  fprintf(stderr, "   -check1        print 24 bit data with compression ratio 1\n");
  fprintf(stderr, "   -Version       version number\n");
}

void cdiDecodeDate(int date, int *year, int *month, int *day);
void cdiDecodeTime(int time, int *hour, int *minute, int *second);

static
void printHead(int rec, int date, int time, int table, double level, int datasize, int nmiss)
{
  if ( ! rec )
  {
    if ( table == -1 )
      fprintf(stdout, 
	      "   Rec :       Date     Time   Level Gridsize    Miss :     Minimum        Mean     Maximum : Code number\n");
    else
      fprintf(stdout, 
              "   Rec :       Date     Time   Level Gridsize    Miss :     Minimum        Mean     Maximum : Parameter ID\n");
    //         ----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+
  }

  int year, month, day, hour, minute, second;
  cdiDecodeDate(date, &year, &month, &day);
  cdiDecodeTime(time, &hour, &minute, &second);

  fprintf(stdout, "%6d : %4.4d-%2.2d-%2.2d %2.2d:%2.2d:%2.2d %7g ", rec+1, year, month, day, hour, minute, second, level);
  fprintf(stdout, "%8d ", datasize);
  fprintf(stdout, "%7d :", nmiss);
}

static
void printTail(int table, int code)
{
  char paramstr[32];

  if ( table == 255 || table == 0 )
    sprintf(paramstr, "%d", code);
  else
    sprintf(paramstr, "%d.%d", code, table);

 if ( table == -1 )
   fprintf(stdout, " : %4d\n", code);
 else
   fprintf(stdout, " : %-11s\n", paramstr);
}

static
void printInfoDP(int date, int time, int table, int code, double level, int datasize, int nmiss, double *data)
{
  static int rec = 0;
  int i, ivals = 0, imiss = 0;
  double arrmean, arrmin, arrmax;

  printHead(rec++, date, time, table, level, datasize, nmiss);

    {
      if ( nmiss > 0 )
	{
	  arrmean = 0;
	  arrmin  =  1.e300;
	  arrmax  = -1.e300;
	  for ( i = 0; i < datasize; i++ )
	    {
	      if ( fabs(data[i]-GRIB_MISSVAL) > 0 )
		{
		  if ( data[i] < arrmin ) arrmin = data[i];
		  if ( data[i] > arrmax ) arrmax = data[i];
		  arrmean += data[i];
		  ivals++;
		}
	    }
	  imiss = datasize - ivals;
	  datasize = ivals;
	}
      else
	{
	  arrmean = 0;
	  arrmin  =  1.e300;
	  arrmax  = -1.e300;
	  for ( i = 0; i < datasize; i++ )
	    {
	      if ( data[i] < arrmin ) arrmin = data[i];
	      if ( data[i] > arrmax ) arrmax = data[i];
	      arrmean += data[i];
	    }
	}

      if ( datasize > 0 ) arrmean /= datasize;

      fprintf(stdout, "%#12.5g%#12.5g%#12.5g", arrmin, arrmean, arrmax);
    }
    
  printTail(table, code);

  if ( imiss != nmiss && nmiss > 0 )
    fprintf(stdout, "Found %d of %d missing values!\n", imiss, nmiss);
}

static
void printInfoSP(int date, int time, int table, int code, double level, int datasize, int nmiss, float *data)
{
  static int rec = 0;
  int i, ivals = 0, imiss = 0;
  double arrmean, arrmin, arrmax;

  printHead(rec++, date, time, table, level, datasize, nmiss);

    {
      if ( nmiss > 0 )
	{
	  arrmean = 0;
	  arrmin  =  1.e300;
	  arrmax  = -1.e300;
	  for ( i = 0; i < datasize; i++ )
	    {
	      if ( fabs(data[i]-GRIB_MISSVAL) > 0 )
		{
		  if ( data[i] < arrmin ) arrmin = data[i];
		  if ( data[i] > arrmax ) arrmax = data[i];
		  arrmean += data[i];
		  ivals++;
		}
	    }
	  imiss = datasize - ivals;
	  datasize = ivals;
	}
      else
	{
	  arrmean = 0;
	  arrmin  =  1.e300;
	  arrmax  = -1.e300;
	  for ( i = 0; i < datasize; i++ )
	    {
	      if ( data[i] < arrmin ) arrmin = data[i];
	      if ( data[i] > arrmax ) arrmax = data[i];
	      arrmean += data[i];
	    }
	}

      if ( datasize > 0 ) arrmean /= datasize;

      fprintf(stdout, "%#12.5g%#12.5g%#12.5g", arrmin, arrmean, arrmax);
    }
    
  printTail(table, code);

  if ( imiss != nmiss && nmiss > 0 )
    fprintf(stdout, "Found %d of %d missing values!\n", imiss, nmiss);
}

static
bool notinlist_int(int nval, int *vallist, int val)
{
  bool found = (nval == 0);

  for ( int i = 0; i < nval; i++ )
    if ( val == vallist[i] ) { found = true; break; }

  return !found;
}

static
bool notinlist_flt(int nval, double *vallist, double val)
{
  bool found = (nval == 0);

  for ( int i = 0; i < nval; i++ )
    if ( ! (fabs(val - vallist[i]) > 0.0) ) { found = true; break; }

  return !found;
}

#define GRIB_BUFSIZE_EXTRA  4096 

int main(int argc, char *argv[])
{
  namelist_t *nml = NULL;
  NML_DEF_INT(bits,     32, 0);
  NML_DEF_INT(code,   1024, 0);
  NML_DEF_INT(day,      32, 0);
  NML_DEF_INT(fcnum,  1024, 0);
  NML_DEF_INT(hour,     25, 0);
  NML_DEF_INT(level,  1024, 0);
  NML_DEF_INT(ltype,   256, 0);
  NML_DEF_INT(minute,   61, 0);
  NML_DEF_INT(month,    16, 0);
  NML_DEF_FLT(param,  1024, 0);
  NML_DEF_INT(rec,       2, 0);
  /* no idea why it's here */
  (void)rec;
  NML_DEF_INT(record, 1024, 0);
  NML_DEF_INT(rdate,  1024, 0);
  NML_DEF_INT(rtime,    64, 0);
  NML_DEF_INT(table,  1024, 0);
  NML_DEF_INT(vdate,  1024, 0);
  NML_DEF_INT(vtime,    64, 0);
  NML_DEF_INT(year,   1024, 0);

  int repair1 = 0;
  int check1 = 0;
  int second;
  bool select = false;
  char *ifile = NULL;
  char *ofile = NULL;
  int nrec = 0;
  int ifileID = 0, ofileID = 0;
  int debug = 0;
  int sp_mode = 0;
  bool copy = 0, move = 0;
  int info = 0, ginfo = 0;
  int pdsinfo = 0, gdsinfo = 0, bmsinfo = 0, bdsinfo = 0;
  int paraminfo = 0;
  int linfo = 0;
  int dinfo = 0;
  int nbits = 0;
  int nread = 0;
  long recpos = 0;
  unsigned char *igribbuffer = NULL;
  unsigned char *itmpbuffer = NULL;
  unsigned char *ogribbuffer = NULL;
  unsigned char *buffer = NULL;
  size_t igribbuffersize = 0;
  size_t itmpbuffersize = 0;
  size_t ogribbuffersize = 0;
  size_t buffersize = 0;
  int iret = 0, iword = 0;
  size_t datasize = 0;
  int isec0[2], isec1[4096], isec2[4096], isec3[2], isec4[512];
  double fsec2[512], fsec3[2], *data = NULL;
  float fsec2f[512], fsec3f[2], *dataf = NULL;
  int date = 0, time = 0;
  int level1, level2;
  double dlevel = 0;
  size_t allocdatasize = 0;
  int regular = 0;
  int gribversion = -1;
  int izip;
  bool szip = false, sunzip = false;
  int i;
  size_t nbytes;
  size_t unzipsize;

  memset(isec0, 0, sizeof(isec0));
  memset(isec1, 0, sizeof(isec1));
  memset(isec2, 0, sizeof(isec2));
  memset(isec3, 0, sizeof(isec3));
  memset(isec4, 0, sizeof(isec4));
  memset(fsec2, 0, sizeof(fsec2));
  memset(fsec3, 0, sizeof(fsec3));
  memset(fsec2f, 0, sizeof(fsec2f));
  memset(fsec3f, 0, sizeof(fsec3f));

  Progname = strrchr(argv[0], '/');
  if ( Progname == NULL ) Progname = argv[0];
  else                    Progname++;

  ISEC3_MissVal = 0;
  FSEC3_MissVal = GRIB_MISSVAL;

  while ( --argc && (*++argv)[0] == '-' )
    {
      int c = *++argv[0];
      char *cstring = argv[0];
      size_t len = strlen(cstring);
      switch (c)
      {
      case 'a':
        if ( !strncmp(cstring, "aec", len) ) szip = true;
        break;
      case 'b':
        if      ( !strncmp(cstring, "bits", len) ) { nbits = atoi((++argv)[0]); --argc; }
        else if ( !strncmp(cstring, "bdsinfo", len) ) bdsinfo = 1;
        else if ( !strncmp(cstring, "bmsinfo", len) ) bmsinfo = 1;
        break;
      case 'c':
        if      ( !strncmp(cstring, "copy", len) ) copy = true;
        else if ( !strncmp(cstring, "check1", len) ) check1 = 1;
        break;   
      case 'd':
        if      ( !strncmp(cstring, "debug", len) ) debug = 1;
        else if ( !strncmp(cstring, "dinfo", len) ) dinfo = 1;
        break;
      case 'f':
        if ( !strncmp(cstring, "float", len) ) sp_mode = 1;
        break;
      case 'g':
        if      ( !strncmp(cstring, "ginfo", len) ) ginfo = 1;
        else if ( !strncmp(cstring, "gdsinfo", len) ) gdsinfo = 1;
        break;
      case 'h':
        if ( !strncmp(cstring, "help", len) ) { usage (); return 0; }
        break;
      case 'i':
        if      ( !strncmp(cstring, "inputfile", len) ) { ifile = (++argv)[0]; --argc; }
        else if ( !strncmp(cstring, "info", 4) ) { info = 1; }
        break;
      case 'l':
        if ( !strncmp(cstring, "linfo", len) ) linfo = 1;
        break;
      case 'o':
        if ( !strncmp(cstring, "outputfile", len) ) { ofile = (++argv)[0];  --argc; }
        break;
      case 'p':
        if      ( len == 1 && *cstring == 'p' )  paraminfo = 1;
        else if ( !strncmp(cstring, "pdsinfo", len) ) pdsinfo = 1;
        break;
      case 'r':
        if      ( !strncmp(cstring, "read", len) )    { nread = atoi((++argv)[0]); --argc; }
        else if ( !strncmp(cstring, "regular", len) ) { regular = 1; }
        else if ( !strncmp(cstring, "repair1", len) ) { repair1 = 1; }
        break;
      case 's':
        if      ( !strncmp(cstring, "szip", len) ) szip = true;
        else if ( !strncmp(cstring, "sunzip", len) ) sunzip = true;
        else if ( !strncmp(cstring, "select", len) ) select = true;
        break;
      case 'T':
        gribTimer = 1;
        break;
      case 'V':
        if ( !strncmp(cstring, "Version", len) ) { version(); return 0; }
        break;
      case 'u':
        if ( !strncmp(cstring, "unaec", len) ) sunzip = true;
        break;
      default:
        usage();
        fprintf(stderr, "illegal option %s\n", cstring);
        return -1;
      }
    }

  if ( argc && !ifile )
    {
      ifile = argv[0];
      argc--;
    }
  if ( argc && !ofile )
    {
      ofile = (++argv)[0];
      argc--;
    }

  if ( !ifile || ( !ofile && copy ) )
    {
      usage();
      return -1;
    }

  if ( debug )
    {
      fprintf(stderr, "\n");
      if ( ifile )
	fprintf(stderr, "  < %s  \n", ifile);
      if ( ofile )
	fprintf(stderr, "  > %s\n\n", ofile);

       gribSetDebug(debug);
    }

  if ( ofile && !copy ) move = true;

  if ( !move && !copy && !check1 && !bdsinfo && !bmsinfo && !ginfo && !gdsinfo &&
       !linfo && !pdsinfo && !info && !dinfo && !nread ) info = 1;

  if ( gribTimer )
    {
      timer_total    = timer_new("total");
      timer_read     = timer_new("read");
      timer_write    = timer_new("write");
      timer_decode   = timer_new("decode");
      timer_encode   = timer_new("encode");
      timer_sunzip   = timer_new("sunzip");
      timer_szip     = timer_new("szip");
    }

  if ( gribTimer ) timer_start(timer_total);

 LABEL_READ:

  ifileID = gribOpen(ifile, "r");
  if ( ifileID == UNDEFID )
    {
      fprintf(stderr, "Open failed on %s\n", ifile);
      return 1;
    }

  if ( ofile )
    {
      ofileID = gribOpen(ofile, "w");
      if ( ofileID == UNDEFID )
	{
	  fprintf(stderr, "Open failed on %s\n", ofile);
	  return 1;
	}
    }

  /* read namelist */
  if ( select )
    {      
      nml = namelistNew("SELECT");

      NML_ADD_INT(nml, bits);
      NML_ADD_INT(nml, code);
      NML_ADD_INT(nml, fcnum);
      NML_ADD_INT(nml, day);
      NML_ADD_INT(nml, hour);
      NML_ADD_INT(nml, level);
      NML_ADD_INT(nml, ltype);
      NML_ADD_INT(nml, minute);
      NML_ADD_INT(nml, month);
      NML_ADD_FLT(nml, param);
      NML_ADD_INT(nml, rec);
      NML_ADD_INT(nml, record);
      NML_ADD_INT(nml, rdate);
      NML_ADD_INT(nml, rtime);
      NML_ADD_INT(nml, table);
      NML_ADD_INT(nml, vdate);
      NML_ADD_INT(nml, vtime);
      NML_ADD_INT(nml, year);

      namelistRead(stdin, nml);
      /* namelistPrint(nml); */

      NML_NUM(nml, bits);
      NML_NUM(nml, code);
      NML_NUM(nml, fcnum);
      NML_NUM(nml, day);
      NML_NUM(nml, hour);
      NML_NUM(nml, level);
      NML_NUM(nml, ltype);
      NML_NUM(nml, minute);
      NML_NUM(nml, month);
      NML_NUM(nml, param);
      NML_NUM(nml, rec);
      NML_NUM(nml, record);
      NML_NUM(nml, rdate);
      NML_NUM(nml, rtime);
      NML_NUM(nml, table);
      NML_NUM(nml, vdate);
      NML_NUM(nml, vtime);
      NML_NUM(nml, year);

      /*
      if ( nselparam )
	{
	  int itable, icode;
	  for ( i = 0; i < nselparam; i++ )
	    {
	      itable = (int) selparam[i];
              icode  = (int) (1000*(selparam[i] - itable));
	    }
	}
      */
    }

  bool ldecode  = ( info || linfo || dinfo || copy );
  bool ljdecode = ( info || linfo || dinfo || copy || select || repair1 );
  bool lencode  = copy;
  bool lzip     = szip;
  bool lunzip   = sunzip || lencode || info;
  /*
  fprintf(stderr, "lzip     %d\n", lzip);
  fprintf(stderr, "lunzip   %d\n", lunzip);
  fprintf(stderr, "lencode  %d\n", lencode);
  fprintf(stderr, "ldecode  %d\n", ldecode);
  fprintf(stderr, "ljdecode %d\n", ljdecode);
  */
  size_t irecsize = 0;
  while ( 1 )
    {
      off_t offset   = -(recpos + (off_t)irecsize);
      irecsize = gribGetSize(ifileID);
      recpos   = gribGetPos(ifileID);
      offset  += recpos;

      if ( irecsize == 0 ) break;

      if ( irecsize > UINT_MAX )
	{
	  fprintf(stderr, "GRIB record size limit exceeded (%zu)\n", irecsize);
	  return 1;
	}

      if ( igribbuffersize < (irecsize+GRIB_BUFSIZE_EXTRA) )
	{
	  igribbuffersize = irecsize+GRIB_BUFSIZE_EXTRA;
	  igribbuffer = (unsigned char *) realloc(igribbuffer, igribbuffersize);
	}

      size_t readsize = irecsize;
      if ( gribTimer ) timer_start(timer_read);
      int ierr = gribRead(ifileID, igribbuffer, &readsize);
      if ( gribTimer ) timer_stop(timer_read);
      if ( ierr != 0 ) break;
      nrec++;

      gribversion = gribVersion(igribbuffer, readsize);

      if ( lunzip && (izip = gribGetZip(irecsize, igribbuffer, &unzipsize)) > 0 )
	{
	  // fprintf(stderr, "compressed data: %d %d %ld\n", izip, irecsize, unzipsize);
	  if ( izip == 128 || izip == 130 ) /* szip or aec */
	    {
	      if ( unzipsize < irecsize )
		{
		  fprintf(stderr, "Decompressed size smaller than compressed size (in %zu; out %zu)!\n",
			  irecsize, unzipsize);
		  continue;
		}

	      if ( itmpbuffersize < irecsize )
		{
		  itmpbuffersize = irecsize;
		  itmpbuffer = (unsigned char *) realloc(itmpbuffer, itmpbuffersize);
		}

	      memcpy(itmpbuffer, igribbuffer, itmpbuffersize);

	      unzipsize += 100; // need 0 to 1 bytes for rounding of bds
	      if ( igribbuffersize < unzipsize )
		{
		  igribbuffersize = unzipsize;
		  igribbuffer = (unsigned char *) realloc(igribbuffer, igribbuffersize);
		}

	      if ( gribTimer ) timer_start(timer_sunzip);
	      irecsize = (size_t)gribUnzip(igribbuffer, (long)unzipsize, itmpbuffer, (long)irecsize);
	      if ( gribTimer ) timer_stop(timer_sunzip);
	    }
	  else
	    {
	      fprintf(stderr, "Decompression for %d not implemented!\n", izip);
	    }
	}

      if ( ljdecode )
	{
	  if ( gribTimer ) timer_start(timer_decode);
	  if ( sp_mode )
	    {
	      gribExSP(isec0, isec1, isec2, fsec2f, isec3, fsec3f, isec4, dataf,
		       (int)datasize, (int *) igribbuffer, (int)igribbuffersize, &iword, "J", &iret);
	    }
	  else
	    {
	      gribExDP(isec0, isec1, isec2, fsec2, isec3, fsec3, isec4, data,
		       (int)datasize, (int *) igribbuffer, (int)igribbuffersize, &iword, "J", &iret);
	    }
	  if ( gribTimer ) timer_stop(timer_decode);
	}

      if ( ldecode )
	{
	  if ( regular && ISEC2_Reduced && iret != -801 )
	    {
	      int nlon = 0;
	      // datasize = ISEC2_NumLat*ISEC2_NumLat*2;
	      for ( int ilat = 0; ilat < ISEC2_NumLat; ++ilat )
		if ( ISEC2_ReducedPoints(ilat) > nlon ) nlon = ISEC2_ReducedPoints(ilat);
	      datasize = (size_t)ISEC2_NumLat*(size_t)nlon;
	    }
	  else
	    datasize = (size_t)ISEC4_NumValues;

	  if ( datasize > UINT_MAX )
	    {
	      fprintf(stderr, "GRIB data size limit exceeded (%zu)\n", datasize);
	      return 1;
	    }

	  if ( allocdatasize < datasize )
	    {
	      allocdatasize = datasize;
	      data  = (double *) realloc(data, datasize*sizeof(double));
	      dataf = (float *) realloc(dataf, datasize*sizeof(float));
	    }

	  if ( gribTimer ) timer_start(timer_decode);

	  const char *oper = regular ? "R" : "D";

	  if ( sp_mode )
	    {
	      gribExSP(isec0, isec1, isec2, fsec2f, isec3, fsec3f, isec4, dataf,
		       (int)datasize, (int *) igribbuffer, (int)igribbuffersize, &iword, oper, &iret);
	    }
	  else
	    {
	      gribExDP(isec0, isec1, isec2, fsec2, isec3, fsec3, isec4, data,
		       (int)datasize, (int *) igribbuffer, (int)igribbuffersize, &iword, oper, &iret);
	    }
	  if ( gribTimer ) timer_stop(timer_decode);
	}

      if ( ljdecode )
	{
	  table  = ISEC1_CodeTable;
	  code   = ISEC1_Parameter;
	  param  = table + ((double)code)/1000;
	  level1 = ISEC1_Level1;
	  level2 = ISEC1_Level2;
	  dlevel = level1;
	  if ( ISEC1_LevelType != GRIB1_LTYPE_HYBRID       &&
	       ISEC1_LevelType != GRIB1_LTYPE_SIGMA        &&
	       ISEC1_LevelType != GRIB1_LTYPE_99           &&
	       ISEC1_LevelType != GRIB1_LTYPE_ISOBARIC     &&
	       ISEC1_LevelType != GRIB1_LTYPE_ISOBARIC_PA  &&
	       ISEC1_LevelType != GRIB1_LTYPE_HEIGHT       &&
	       ISEC1_LevelType != GRIB1_LTYPE_LANDDEPTH    &&
	       ISEC1_LevelType != GRIB1_LTYPE_SEADEPTH )
	    {
	      dlevel = (level1+level2)/2.;
	    }
      
	  if ( ISEC1_LevelType == GRIB1_LTYPE_ISOBARIC ) dlevel *= 100;
	  
	  gribDateTime(isec1, &date, &time);
	  
	  vdate = date;
	  vtime = time;

	  rdate = gribRefDate(isec1);
	  rtime = gribRefTime(isec1);
	}

      if ( select )
	{
	  record++;

	  if ( nselrec >= 1 && record < selrec[0] ) continue;
	  if ( nselrec == 2 && record > selrec[1] ) break;

	  cdiDecodeDate(date, &year, &month, &day);
	  cdiDecodeTime(time, &hour, &minute, &second);

	  ltype = ISEC1_LevelType;
	  level = (int) dlevel;
	  bits  = ISEC4_NumBits;

	  fcnum = 0;
	  
	  if ( (ISEC1_CenterID    == 98 && ISEC1_LocalFLag ==  1) ||
	       (ISEC1_SubCenterID == 98 && ISEC1_LocalFLag ==  1) ||
	       (ISEC1_CenterID    ==  7 && ISEC1_SubCenterID == 98) )
	    if ( isec1[36] == 1 )
	      fcnum = isec1[41];

	  if ( notinlist_int(NML_PAR(bits))   ||
	       notinlist_int(NML_PAR(code))   ||
	       notinlist_int(NML_PAR(fcnum))  ||
	       notinlist_int(NML_PAR(day))    ||
	       notinlist_int(NML_PAR(hour))   ||
	       notinlist_int(NML_PAR(level))  ||
	       notinlist_int(NML_PAR(ltype))  ||
	       notinlist_int(NML_PAR(minute)) ||
	       notinlist_int(NML_PAR(month))  ||
	       notinlist_flt(NML_PAR(param))  ||
	       notinlist_int(NML_PAR(record)) ||
	       notinlist_int(NML_PAR(rdate))  ||
	       notinlist_int(NML_PAR(rtime))  ||
	       notinlist_int(NML_PAR(table))  ||
	       notinlist_int(NML_PAR(vdate))  ||
	       notinlist_int(NML_PAR(vtime))  ||
	       notinlist_int(NML_PAR(year)) )
	    continue;
	}

      if ( info )
	{
          int nmiss = (ISEC1_Sec2Or3Flag & 64) ? ((int)datasize - ISEC4_NumNonMissValues) : 0;

	  // printf("nmiss %d %d %d\n", datasize, ISEC4_NumNonMissValues, nmiss);
	  if ( paraminfo == 0 ) table = -1;
	  if ( sp_mode )
	    printInfoSP(date, time, table, code, dlevel, (int)datasize, nmiss, dataf);
	  else
	    printInfoDP(date, time, table, code, dlevel, (int)datasize, nmiss, data);
	}
      else if ( linfo )
	{
	  gribPrintSec0(isec0);
	  gribPrintSec1(isec0, isec1);

	  if ( sp_mode )
	    {
	      if ( (isec1[4] == 0) || (isec1[4] == 64) )
		printf(" No section 2 in GRIB message.\n");
	      else
		gribPrintSec2SP(isec0, isec2, fsec2f);
	      if ( (isec1[4] == 0) || (isec1[4] == 128) )
		printf(" No section 3 in GRIB message.\n");
	      else
		gribPrintSec3SP(isec0, isec3, fsec3f);
	      
	      gribPrintSec4SP(isec0, isec4, dataf);
	    }
	  else
	    {
	      if ( (isec1[4] == 0) || (isec1[4] == 64) )
		printf(" No section 2 in GRIB message.\n");
	      else
		gribPrintSec2DP(isec0, isec2, fsec2);
	      if ( (isec1[4] == 0) || (isec1[4] == 128) )
		printf(" No section 3 in GRIB message.\n");
	      else
		gribPrintSec3DP(isec0, isec3, fsec3);
	      
	      gribPrintSec4DP(isec0, isec4, data);
	    }
	}
      else if ( dinfo )
	{
	  printf("Code %d   Level %d\n", ISEC1_Parameter, ISEC1_Level1);
	  for ( i = 0; i < isec4[0]; i++ )
	    {
	      printf(" %6d  %g\n", i+1, data[i]);
	    }
	}
      else if ( ginfo   )  gribPrintALL(nrec, (long)offset, recpos, (long)irecsize, igribbuffer);
      else if ( pdsinfo )  gribPrintPDS(nrec, recpos, (long)irecsize, igribbuffer);
      else if ( gdsinfo )  gribPrintGDS(nrec, recpos, (long)irecsize, igribbuffer);
      else if ( bmsinfo )  gribPrintBMS(nrec, recpos, (long)irecsize, igribbuffer);
      else if ( bdsinfo )  gribPrintBDS(nrec, recpos, (long)irecsize, igribbuffer);
      else if ( check1  )  gribCheck1(nrec, recpos, (long)irecsize, igribbuffer);
      else if ( repair1 )  gribRepair1(nrec, (long)irecsize, igribbuffer);

      if ( ofile )
	{
	  if ( lencode )
	    {
	      if ( nbits > 0 ) isec4[1] = nbits;
	      size_t orecsize = (size_t)isec4[0]*(size_t)isec4[1]/8 + (size_t)isec4[0]/8 + 3000; /* datasize + mask + header */
	      if ( ogribbuffersize < orecsize )
		{
		  ogribbuffersize = orecsize;
		  ogribbuffer = (unsigned char *) realloc(ogribbuffer, ogribbuffersize);
		}
	      /*
	      printf("datasize = %d\n", isec4[0]);
	      printf("nbits = %d\n", isec4[1]);
	      printf("igribbuffersize = %d\n", igribbuffersize);
	      printf("ogribbuffersize = %d\n", ogribbuffersize);
	      */
	      if ( ISEC1_CodeTable == 0 ) ISEC1_CodeTable = 128;   // bug fix
	      if ( ISEC2_RepMode == 128 &&
		   ISEC2_GridType == 50 ) ISEC2_RepMode = 1;       // bug fix
	      /*
	      fsec3[1] = 999999.9;
              isec1[4] = 192;
	      */

	      if ( gribTimer ) timer_start(timer_encode);
	      if ( sp_mode )
		{
		  /*
		  for ( i = 0; i < (int)(sizeof(fsec2)/sizeof(double)); ++i ) fsec2f[i] = (float) fsec2[i];
		  for ( i = 0; i < (int)(sizeof(fsec3)/sizeof(double)); ++i ) fsec3f[i] = (float) fsec3[i];
		  for ( i = 0; i < (int)datasize; ++i ) dataf[i] = (float) data[i];
		  */
		  gribExSP(isec0, isec1, isec2, fsec2f, isec3, fsec3f, isec4, dataf,
			   (int)datasize, (int *) ogribbuffer, (int)ogribbuffersize, &iword, "C", &iret);
		}
	      else
		{
		  gribExDP(isec0, isec1, isec2, fsec2, isec3, fsec3, isec4, data,
			   (int)datasize, (int *) ogribbuffer, (int)ogribbuffersize, &iword, "C", &iret);
		}
	      if ( gribTimer ) timer_stop(timer_encode);

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

	      if ( lzip && gribversion == 1 )
		{
		  if ( buffersize < ogribbuffersize )
		    {
		      buffersize = ogribbuffersize;
		      buffer = (unsigned char *) realloc(buffer, buffersize);
		    }
		  // memcpy(buffer, ogribbuffer, ogribbuffersize);
		  if ( gribTimer ) timer_start(timer_szip);
		  nbytes = (size_t)gribZip(ogribbuffer, (long)ogribbuffersize, buffer, (long)buffersize);
		  if ( gribTimer ) timer_stop(timer_szip);
		}

	      while ( nbytes & 7 ) ogribbuffer[nbytes++] = 0;

	      if ( gribTimer ) timer_start(timer_write);
	      (void) gribWrite(ofileID, ogribbuffer, nbytes);
	      if ( gribTimer ) timer_stop(timer_write);
	    }
	  else
	    {
	      nbytes = irecsize;

	      if ( lzip && gribversion == 1 )
		{
		  if ( buffersize < igribbuffersize )
		    {
		      buffersize = igribbuffersize;
		      buffer = (unsigned char *) realloc(buffer, buffersize);
		    }
		  // memcpy(buffer, igribbuffer, igribbuffersize);
		  if ( gribTimer ) timer_start(timer_szip);
		  nbytes = (size_t)gribZip(igribbuffer, (long)igribbuffersize, buffer, (long)buffersize);
		  if ( gribTimer ) timer_stop(timer_szip);
		}

	      if ( nbytes > JP24SET )
		while ( nbytes%120 ) igribbuffer[nbytes++] = 0;
	      else
		while ( nbytes & 7 ) igribbuffer[nbytes++] = 0;
	  
	      if ( gribTimer ) timer_start(timer_write);
	      (void) gribWrite(ofileID, igribbuffer, nbytes);
	      if ( gribTimer ) timer_stop(timer_write);
	    }
	}
    }

  if ( ofile ) gribClose(ofileID);

  gribClose(ifileID);

  if ( --nread > 0 ) goto LABEL_READ;

  if ( igribbuffer ) free(igribbuffer);
  if ( itmpbuffer  ) free(itmpbuffer);
  if ( ogribbuffer ) free(ogribbuffer);
  if ( buffer ) free(buffer);
  if ( data ) free(data);
  if ( dataf ) free(dataf);

  if ( select ) namelistDelete(nml);
 
  if ( gribTimer ) timer_stop(timer_total);

  if ( gribTimer ) timer_report();

  return 0;
}
