#include "stdlib.h"
#include "stdio.h"
#include "netcdf.h"
#include "errno.h"
#include "string.h"
#include <sys/types.h>
#include <sys/stat.h>
#include <fcntl.h>


/* gcc main.c -I/usr/include/netcdf-3 -L/usr/lib/netcdf-3 -lnetcdf */

main(int argc, char *argv[])
{
 int ncfid;
 int rc;
 int tile_number, sNx, sNy, nSx, nSy, nPx, nPy, Nx, Ny, Nr, iLo, jLo;
 int vId, dId;
 int dimT, dimX, dimY, dimZ, dError;
 float *arrIn, *arrOut;
 char  *fNam, *vNam,fNam1[1024],fNam2[1024];
 int i, j, k, t, iOffIn, iOffOut, iOffFile;
 long iSeek;
 int ti, tj;
 int fid;
 
 int dbLev;

 dbLev=0;

 /* Open a tile file */
  fNam=argv[1];
  rc=nc_open( fNam, 0, &ncfid );
  vNam = argv[2];

 /* Get global attributes */
 /* o tile number       */
 rc=nc_get_att_int( ncfid, NC_GLOBAL,"tile_number", &tile_number );
 if( dbLev > 0 ) printf( "tile_number == %d\n", tile_number );
 /* o sNx, sNy, nSx, nSy, nPx, nPy, Nx, Ny, Nr */
 rc=nc_get_att_int(ncfid,NC_GLOBAL,"sNx",&sNx);
 if ( dbLev > 2 ) printf("sNx == %d\n",sNx);
 rc=nc_get_att_int(ncfid,NC_GLOBAL,"sNy",&sNy);
 if ( dbLev > 2 ) printf("sNy == %d\n",sNy);
 rc=nc_get_att_int(ncfid,NC_GLOBAL,"nSx",&nSx);
 if ( dbLev > 2 ) printf("nSx == %d\n",nSx);
 rc=nc_get_att_int(ncfid,NC_GLOBAL,"nSy",&nSy);
 if ( dbLev > 2 ) printf("nSy == %d\n",nSy);
 rc=nc_get_att_int(ncfid,NC_GLOBAL,"nPx",&nPx);
 if ( dbLev > 2 ) printf("nPx == %d\n",nPx);
 rc=nc_get_att_int(ncfid,NC_GLOBAL,"nPy",&nPy);
 if ( dbLev > 2 ) printf("nPy == %d\n",nPy);
 rc=nc_get_att_int(ncfid,NC_GLOBAL,"Nx",&Nx);
 if ( dbLev > 2 ) printf("Nx == %d\n",Nx);
 rc=nc_get_att_int(ncfid,NC_GLOBAL,"Ny",&Ny);
 if ( dbLev > 2 ) printf("Ny == %d\n",Ny);
 rc=nc_get_att_int(ncfid,NC_GLOBAL,"Nr",&Nr);
 if ( dbLev > 2 ) printf("Nr == %d\n",Nr);

/* Check dimensions - for now these should match */
 /* sNx, sNy, Nr and have T = 1                   */
 rc=nc_inq_dimid( ncfid, "X", &dId );
 rc=nc_inq_dimlen( ncfid, dId, &dimX );
 if ( dbLev > 2 ) printf( "dimX == %d\n", dimX );
 rc=nc_inq_dimid( ncfid, "Y", &dId );
 rc=nc_inq_dimlen( ncfid, dId, &dimY );
 if( dbLev > 2 ) printf( "dimY == %d\n", dimY );
 rc=nc_inq_dimid( ncfid, "Z", &dId );
 rc=nc_inq_dimlen( ncfid, dId, &dimZ );
 if( dbLev > 2 ) printf( "dimZ == %d\n", dimZ );
 rc=nc_inq_dimid( ncfid, "T", &dId );
 rc=nc_inq_dimlen( ncfid, dId, &dimT );
 if( dbLev > 2 ) printf( "dimT == %d\n", dimT );

 dError = 0;
 if( dimX != sNx ) {dError=1;};
 if( dimY != sNy ) {dError=1;};
 if( dimZ != Nr  ) {dError=1;};
 if( dimT !=  1  ) {dError=1;};
 if( dError != 0 ) {
  printf(
  "Error: Inconsistent dimension variable sizes and global variable attribute settings\n");
 }


 /* Get a piece of the data */
 arrIn  = calloc( sNx*sNy*Nr, sizeof(*arrIn) );
 arrOut = calloc( sNx*sNy*Nr, sizeof(*arrOut) );
 rc = nc_inq_varid( ncfid, vNam, &vId );
 rc = nc_get_var_float( ncfid, vId, arrIn );
 /* Reorder from T Z Y X to X Y Z */
 iOffIn=0;
 for( t=1; t<2; t++ ){
  for( k=1; k<=Nr; k++ ){
   for( j=1; j<=sNy; j++ ){
    for( i=1; i<=sNx; i++ ){
     iOffOut=(k-1)*sNx*sNy+(j-1)*sNx+i-1;
     arrOut[iOffOut]=arrIn[iOffIn];
     iOffIn++;
    }
   }
  }
 }


 /* Put my piece in the global binary file */
 /* tRank=(tj-1)*nSx+ti;                   */
 /* ti=mod(tRank,nSx)+1;                   */
 /* tj=tRank/nSx+1;                        */
 ti=(tile_number-1)%(nSx*nPx)+1;
 iLo=(ti-1)*sNx+1;
 if( dbLev > 2 ) printf( "ti == %d\n", ti );
 if( dbLev > 2 ) printf( "iLo == %d\n", iLo );
 tj=(tile_number-1)/(nSx*nPx)+1;
 jLo=(tj-1)*sNy+1;
 if( dbLev > 2 ) printf( "tj == %d\n", tj );
 if( dbLev > 2 ) printf( "jLo == %d\n", jLo );
 /* Write data to global file one row at a time */
 sscanf( fNam, "%[^.].%[^.].", fNam1, fNam2 );
 sprintf( fNam1, "%s.%s.data", vNam, fNam2 );
 printf( "Writing to %s\n", fNam1 );
 fid=open( fNam1, O_CREAT | O_RDWR , S_IRWXU );
 if( fid <= 0 ){
  perror( "open error" );
  exit( 0 );
 }
 for( t=1; t<2; t++ ){
  for( k=1; k<=Nr; k++ ){
   for( j=1; j<=sNy; j++ ){
    iOffOut=(k-1)*sNx*sNy+(j-1)*sNx;
    iOffFile=(k-1)*Nx*Ny+(jLo-1)*Nx+(j-1)*Nx+iLo-1;
    iSeek=iOffFile*sizeof(*arrOut);
    /* printf( "Seeking to %d\n", iSeek ); */
    rc=lseek( fid, (off_t)iSeek, SEEK_SET );
    if ( rc != iSeek ){
     perror( "lseek error" );
     printf( "rc = %d\n", rc );
     exit( 0 );
    }
    rc=write( fid, arrOut+iOffOut, sizeof(*arrOut)*sNx );
   }
  }
 }
 close( fid );

}