/*
 * EasyWave - A realtime tsunami simulation program with GPU support.
 * Copyright (C) 2014  Andrey Babeyko, Johannes Spazier
 * GFZ German Research Centre for Geosciences (http://www.gfz-potsdam.de)
 *
 * Parts of this program (especially the GPU extension) were developed
 * within the context of the following publicly funded project:
 * - TRIDEC, EU 7th Framework Programme, Grant Agreement 258723
 *   (http://www.tridec-online.eu)
 *
 * This program is free software: you can redistribute it and/or modify
 * it under the terms of the GNU Affero General Public License as
 * published by the Free Software Foundation, either version 3 of the
 * License, or (at your option) any later version.
 *
 * This program is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU Affero General Public License for more details.
 *
 * You should have received a copy of the GNU Affero General Public License
 * along with this program.  If not, see <http://www.gnu.org/licenses/>.
 */

#include <stdio.h>
#include <string.h>
#include <utilits.h>
#include "easywave.h"


static char *IndexFile;
static int Nrec2DOutput;


int ewStart2DOutput()
{
  FILE *fp;
  char buf[64];


  // start index file
  sprintf( buf, "%s.2D.idx", Par.modelName );
  IndexFile = strdup(buf);

  fp = fopen( IndexFile, "wt" );

  fprintf( fp, "%g %g %d %g %g %d\n", LonMin, LonMax, NLon, LatMin, LatMax, NLat );

  fclose( fp );

  Nrec2DOutput = 0;

  return 0;
}



int ewOut2D()
{
  FILE *fp;
  short nOutI,nOutJ;
  int i,j,m;
  float ftmp;
  double dtmp,lonOutMin,lonOutMax,latOutMin,latOutMax;
  char record[128];

  CNode& Node = *gNode;

  Nrec2DOutput++;

  nOutI = Imax-Imin+1;
  lonOutMin = getLon(Imin); lonOutMax = getLon(Imax);
  nOutJ = Jmax-Jmin+1;
  latOutMin = getLat(Jmin); latOutMax = getLat(Jmax);

  // write ssh
  sprintf( record, "%s.2D.%5.5d.ssh", Par.modelName, Par.time );
  fp = fopen( record, "wb" );
  fwrite( "DSBB", 4, 1, fp );
  fwrite( &nOutI, sizeof(short), 1, fp );
  fwrite( &nOutJ, sizeof(short), 1, fp );
  fwrite( &lonOutMin, sizeof(double), 1, fp );
  fwrite( &lonOutMax, sizeof(double), 1, fp );
  fwrite( &latOutMin, sizeof(double), 1, fp );
  fwrite( &latOutMax, sizeof(double), 1, fp );
  dtmp = -1.; fwrite( &dtmp, sizeof(double), 1, fp );
  dtmp = +1.; fwrite( &dtmp, sizeof(double), 1, fp );
  for( j=Jmin; j<=Jmax; j++ ) {
    for( i=Imin; i<=Imax; i++ ) {
      m = idx(j,i);
      if( fabs(Node(m, iH)) < Par.sshTransparencyThreshold )
        ftmp = (float)9999;
      else
        ftmp = (float)Node(m, iH);
      fwrite( &ftmp, sizeof(float), 1, fp );
    }
  }
  fclose( fp );

  // updating contents file
  fp = fopen( IndexFile, "at" );
  fprintf( fp, "%3.3d %s %d %d %d %d\n", Nrec2DOutput, utlTimeSplitString(Par.time), Imin, Imax, Jmin, Jmax );
  fclose( fp );

  return 0;
}


int ewDump2D()
{
  FILE *fp;
  short nOutI,nOutJ;
  int i,j,m;
  float ftmp;
  double dtmp,lonOutMin,lonOutMax,latOutMin,latOutMax;
  char record[128];

  CNode& Node = *gNode;

  nOutI = Imax-Imin+1;
  lonOutMin = getLon(Imin); lonOutMax = getLon(Imax);
  nOutJ = Jmax-Jmin+1;
  latOutMin = getLat(Jmin); latOutMax = getLat(Jmax);

  // write ssh max
  sprintf( record, "%s.2D.sshmax", Par.modelName );
  fp = fopen( record, "wb" );
  fwrite( "DSBB", 4, 1, fp );
  fwrite( &nOutI, sizeof(short), 1, fp );
  fwrite( &nOutJ, sizeof(short), 1, fp );
  fwrite( &lonOutMin, sizeof(double), 1, fp );
  fwrite( &lonOutMax, sizeof(double), 1, fp );
  fwrite( &latOutMin, sizeof(double), 1, fp );
  fwrite( &latOutMax, sizeof(double), 1, fp );
  dtmp = 0.; fwrite( &dtmp, sizeof(double), 1, fp );
  dtmp = 1.; fwrite( &dtmp, sizeof(double), 1, fp );
  for( j=Jmin; j<=Jmax; j++ ) {
    for( i=Imin; i<=Imax; i++ ) {
      ftmp = (float)Node(idx(j,i), iHmax);
      fwrite( &ftmp, sizeof(float), 1, fp );
    }
  }
  fclose( fp );

  // write arrival times
  sprintf( record, "%s.2D.time", Par.modelName );
  fp = fopen( record, "wb" );
  fwrite( "DSBB", 4, 1, fp );
  fwrite( &nOutI, sizeof(short), 1, fp );
  fwrite( &nOutJ, sizeof(short), 1, fp );
  fwrite( &lonOutMin, sizeof(double), 1, fp );
  fwrite( &lonOutMax, sizeof(double), 1, fp );
  fwrite( &latOutMin, sizeof(double), 1, fp );
  fwrite( &latOutMax, sizeof(double), 1, fp );
  dtmp = 0.; fwrite( &dtmp, sizeof(double), 1, fp );
  dtmp = 1.; fwrite( &dtmp, sizeof(double), 1, fp );
  for( j=Jmin; j<=Jmax; j++ ) {
    for( i=Imin; i<=Imax; i++ ) {
      ftmp = (float)Node(idx(j,i), iTime) / 60;
      fwrite( &ftmp, sizeof(float), 1, fp );
    }
  }
  fclose( fp );

  return 0;
}