/*@@
   @file      Write.c
   @date      Thu 18 April 2002
   @author    Thomas Radke
   @desc
              Output two-dimensional slices in Jpeg image format.
   @enddesc
   @version   $Id$
 @@*/


#include <stdlib.h>
#include <string.h>

#include "cctk.h"
#include "cctk_Parameters.h"
#include "CactusBase/IOUtil/src/ioutil_AdvertisedFiles.h"
#include "ioJpegGH.h"

/* the rcs ID and its dummy function to use it */
static const char *rcsid = "$Header$";
CCTK_FILEVERSION(CactusIO_IOJpeg_Write_c)


/********************************************************************
 ********************    Internal Routines   ************************
 ********************************************************************/
static void WriteData (const cGH *GH, int vindex, const char *alias, int dim,
                       int dir, CCTK_REAL min, CCTK_REAL max,
                       const CCTK_INT hsize[2], const void *hdata);


/*@@
   @routine   IOJpeg_Write
   @date      Thu 18 April 2002
   @author    Thomas Radke
   @desc
              Writes the slices of a variable into separate Jpeg files.
   @enddesc
   @calls     Hyperslab_DefineGlobalMappingByIndex
              Hyperslab_FreeMapping
              Hyperslab_GetList
              OpenFile
              WriteData

   @var       GH
   @vdesc     pointer to CCTK GH
   @vtype     const cGH *
   @vio       in
   @endvar
   @var       vindex
   @vdesc     index of variable to output
   @vtype     CCTK_INT
   @vio       in
   @endvar
   @var       alias
   @vdesc     alias name of variable to output
   @vtype     const char *
   @vio       in
   @endvar

   @returntype int
   @returndesc
                0 for success, or<BR>
               -1 if variable has no storage assigned<BR>
               -2 if output file couldn't be opened<BR>
               -3 if hyperslab for coordinates and/or variable couldn't be
                  extracted
   @endreturndesc
@@*/
int IOJpeg_Write (const cGH *GH, CCTK_INT vindex, const char *alias)
{
  const ioJpegGH *myGH;
  int i, mapping, total_hsize;
  int dir, dir_i, dir_j, maxdir, myproc, groupindex;
  cGroup gdata;
  char *fullname;
  int extent_int[3];
  CCTK_INT origin[3], extent[2], direction[6], hsize[2];
  void *hdata;
  CCTK_REAL min, max;
  const CCTK_INT htype = CCTK_VARIABLE_REAL;
  /*** FIXME: can CCTK_Reduce() have a 'const cGH *' parameter ?? ***/
  union
  {
    const cGH *const_ptr;
    cGH *non_const_ptr;
  } GH_fake_const;
  DECLARE_CCTK_PARAMETERS


  /* get the variable name and group information */
  fullname = CCTK_FullName (vindex);
  groupindex = CCTK_GroupIndexFromVarI (vindex);
  CCTK_GroupData (groupindex, &gdata);

  /* check if variable has storage assigned */
  if (! CCTK_QueryGroupStorageI (GH, groupindex))
  {
    CCTK_VWarn (2, __LINE__, __FILE__, CCTK_THORNSTRING,
                "No IOJpeg output for '%s' (no storage)", fullname);
    free (fullname);
    return (-1);
  }

  /* get the handle for IOJpeg extensions */
  myproc = CCTK_MyProc (GH);
  myGH = (const ioJpegGH *) CCTK_GHExtension (GH, "IOJpeg");

  /* set the minimum/maximum for the colormap */
  if (CCTK_Equals (colormap, "custom"))
  {
    min = colormap_min;
    max = colormap_max;
  }
  else
  {
    GH_fake_const.const_ptr = GH;
    i = CCTK_ReductionHandle ("minimum");
    CCTK_Reduce (GH_fake_const.non_const_ptr, 0, i, 1, CCTK_VARIABLE_REAL,
                 &min, 1, vindex);
    i = CCTK_ReductionHandle ("maximum");
    CCTK_Reduce (GH_fake_const.non_const_ptr, 0, i, 1, CCTK_VARIABLE_REAL,
                 &max, 1, vindex);
  }

  /* get the number of slices to output */
  /* in general: maxdir = gdata.dim * (gdata.dim - 1) / 2; */
  maxdir = gdata.dim == 2 ? 1 : 3;

  /* get the extents of the variable */
  CCTK_GroupgshVI (GH, 3, extent_int, vindex);

  /* now do the actual I/O looping over all directions */
  for (dir = 0; dir < maxdir; dir++)
  {
    /* get the directions to span the hyperslab */
    if (dir == 0)
    {
      dir_i = 0; dir_j = 1;   /* xy */
    }
    else if (dir == 1)
    {
      dir_i = 0; dir_j = 2;   /* xz */
    }
    else
    {
      dir_i = 1; dir_j = 2;   /* yz */
    }

    /* set the extent vector */
    extent[0] = extent_int[dir_i];
    extent[1] = extent_int[dir_j];

    /* set the origin using the slice center from IOUtil */
    memset (origin, 0, sizeof (origin));
    if (gdata.grouptype == CCTK_GF)
    {
      origin[maxdir-dir-1] = myGH->sp2xyz[gdata.dim - 1][dir];
    }

    /* set the direction vector */
    memset (direction, 0, sizeof (direction));
    direction[dir_i] = direction[gdata.dim + dir_j] = 1;

    mapping = Hyperslab_DefineGlobalMappingByIndex (GH, vindex, 2,
                                                    direction, origin, extent,
                                                    NULL, -1, NULL, hsize);
    if (mapping < 0)
    {
      CCTK_VWarn (1, __LINE__, __FILE__, CCTK_THORNSTRING,
                  "Failed to define hyperslab mapping for variable '%s'",
                  fullname);
      continue;
    }
    total_hsize = hsize[0] * hsize[1] * CCTK_VarTypeSize (gdata.vartype);
    if (total_hsize <= 0)
    {
      CCTK_VWarn (1, __LINE__, __FILE__, CCTK_THORNSTRING,
                  "Selected hyperslab has zero size for variable '%s' "
                  "direction %d", fullname, dir);
      Hyperslab_FreeMapping (mapping);
      continue;
    }

    /* allocate hyperslab buffer */
    hdata = myproc == 0 ? malloc (total_hsize) : NULL;

    /* get the hyperslab */
    i = Hyperslab_Get (GH, mapping, 0, vindex, 0, htype, hdata);
    if (i)
    {
      CCTK_VWarn (1, __LINE__, __FILE__, CCTK_THORNSTRING,
                  "Failed to extract hyperslab for variable '%s'", fullname);
    }

    /* release the mapping structure */
    Hyperslab_FreeMapping (mapping);

    /* and dump the data to file */
    if (myproc == 0)
    {
      if (i == 0)
      {
        WriteData (GH, vindex, alias, gdata.dim, dir, min, max, hsize, hdata);
      }

      /* clean up */
      free (hdata);
    }
  } /* end of looping through xyz directions */

  free (fullname);

  return (0);
}


/********************************************************************
 ********************    Internal Routines   ************************
 ********************************************************************/
/*@@
   @routine    WriteData
   @date       Thu 18 April 2002
   @author     Thomas Radke
   @desc
	       Writes the given hyperslab into a Jpeg output file.
   @enddesc
 @@*/
static void WriteData (const cGH *GH, int vindex, const char *alias, int dim,
                       int dir, CCTK_REAL min, CCTK_REAL max,
                       const CCTK_INT hsize[2], const void *hdata)
{
  FILE *file;
  unsigned char *dataout;
  const ioJpegGH *myGH;
  ioAdvertisedFileDesc advertised_file;
  char *filename, *fullname;
  char slicename[30];
  const char *extensions[] = {"xy", "xz", "yz"};
  DECLARE_CCTK_PARAMETERS


  myGH = (const ioJpegGH *) CCTK_GHExtension (GH, "IOJpeg");

  /* allocate the RGB image buffer */
  dataout = (unsigned char *) malloc (3 * hsize[0] * hsize[1]);

  AutoColorDataSlice (hsize[0], hsize[1], hdata, dataout, min, max,
                      colormap_bias, colormap_factor);

  /* open the file */
  if (dim == 2)
  {
    strcpy (slicename, "2D");
  }
  else
  {
    /* give the slice origin as range [1 .. n] */
    sprintf (slicename, "%s_[%d]", extensions[dir], myGH->sp2xyz[dim-1][dir]);
  }

  filename = (char *) malloc (strlen (myGH->out_dir) + strlen (alias) +
                              sizeof (slicename) + 20);

  if (CCTK_Equals (mode, "remove"))
  {
    sprintf (filename, "%s%s_%s.jpeg", myGH->out_dir, alias, slicename);
  }
  else
  {
    sprintf (filename, "%s%s_%s.it_%d.jpeg", myGH->out_dir, alias, slicename,
             GH->cctk_iteration);
  }

  file = fopen (filename, "w");
  if (file)
  {
    /* write the data */
    WriteJPEGToFileRGB (hsize[0], hsize[1], dataout, colormap_quality, file);

    /* advertise the file for downloading */
    if (CCTK_Equals (mode, "remove") && myGH->out_last[vindex] < 0)
    {
      fullname = CCTK_FullName (vindex);
      advertised_file.slice = slicename;
      advertised_file.thorn = CCTK_THORNSTRING;
      advertised_file.varname = fullname;
      advertised_file.description = "Jpegs of slices";
      advertised_file.mimetype = "image/jpeg";

      IOUtil_AdvertiseFile (GH, filename, &advertised_file);

      free (fullname);
    }

    /* close the file */
    fclose (file);
  }
  else
  {
    CCTK_VWarn (1, __LINE__, __FILE__, CCTK_THORNSTRING,
                "Cannot open IOJpeg output file '%s'", filename);
  }

  /* clean up */
  free (dataout);
  free (filename);
}