path: root/src/DumpVar.c

/*@@
   @file      DumpVar.c
   @date      Fri May 21 1999
   @author    Thomas Radke / Gerd Lanfermann
   @desc 
   Do the actual writing of a variable.
   @enddesc 
   @history
   @hauthor Thomas Radke @hdate May 21 1999
   @hdesc Just copied from thorn FlexIO.
   @hendhistory
 @@*/


#include <stdio.h>
#include <stdlib.h>

#include "cctk.h"
#include "cctk_Parameters.h"
#include "CactusBase/IOUtil/src/ioGH.h"
#include "CactusPUGH/PUGHSlab/src/PUGHSlab.h"
#include "StreamedHDF5GH.h"


/*#define StreamedHDF5_DEBUG    1 */

/* local function prototypes */
static void StreamedHDF5_Dump (cGH *GH, int index, int timelevel, void *outme,
                               int hdim, int *hsizes, int hdf5_type, hid_t fid);
static void StreamedHDF5_AddCommonAttributes (cGH *GH, int index, int timelevel,
                                              hid_t dataset);
static int GetHDF5type (StreamedHDF5GH *myGH, int vtype, hid_t *hdf5_type);

/*@@
  @routine    StreamedHDF5_DumpVar
  @date       16 Apr 1999
  @author     Thomas Radke
  @desc
              Generic dump routine, just calls the appropriate dump routine
	      Extracts the data with Hyperslab.
  @enddesc
  @calledby   StreamedHDF5_DumpGH IOHDF5_Write3D
  @history
  @endhistory
  @var        GH
  @vdesc      Pointer to CCTK grid hierarchy
  @vtype      cGH
  @vio        in
  @endvar
  @var        index
  @vdesc      global index of the variable to be dumped
  @vtype      int
  @vio        in
  @endvar
  @var        timelevel
  @vdesc      the timelevel to store
  @vtype      int
  @vio        in
  @endvar
  @var        fid
  @vdesc      the HDF5 file handle
  @vtype      hid_t
  @vio        in
  @endvar
@@*/
int StreamedHDF5_DumpVar (cGH *GH, int index, int timelevel, hid_t fid)
{
  DECLARE_CCTK_PARAMETERS
  void *data;
  int *hsizes;
  hid_t hdf5_type;
  StreamedHDF5GH *myGH;
  StreamGeo_t geo;
  int vtype;
  
  int sdir[3] = {0,0,0};


  /* Get the handle for StreamedHDF5 extensions */
  myGH = (StreamedHDF5GH *) GH->extensions [CCTK_GHExtensionHandle ("StreamedHDF5")];

  vtype = CCTK_VarTypeI (index);
  if (GetHDF5type (myGH, vtype, &hdf5_type) < 0)
  {
    CCTK_VWarn (1, __LINE__, __FILE__, CCTK_THORNSTRING,
                "Unsupported variable type %d", vtype);
    return (-1);
  }
  /* shortcut to the slab geometry */
  geo = myGH->geo_output[index];

  /* get the dimension of the variable */

  hsizes = (int*)malloc(geo.sdim*sizeof(int));

  data = NULL;

  /* TEMPORARY FIX DUE TO INCONSISTENT DIR SPECIFICATION */
  /* direction vector of 1d line in 3d volume */
  if (geo.sdim==1) {
    sdir[0] = (geo.direction[0]==0) ? 1:0;
    sdir[1] = (geo.direction[0]==1) ? 1:0;
    sdir[2] = (geo.direction[0]==2) ? 1:0;
  } 
  /* norm vector of 2d surface in 3d volume */
  else if (geo.sdim==2)
  {
    sdir[0] = (geo.direction[1]&&geo.direction[2]) ? 1:0;
    sdir[1] = (geo.direction[0]&&geo.direction[2]) ? 1:0;
    sdir[2] = (geo.direction[0]&&geo.direction[1]) ? 1:0;
  }
  /* spanning directions for 3d */
  else if (geo.sdim==3)
  {
    sdir[0] = geo.direction[0];
    sdir[1] = geo.direction[1];
    sdir[2] = geo.direction[2];
  }

  if (Hyperslab_GetHyperslab (GH, 0, index, timelevel, geo.sdim, geo.origin,
                              sdir, geo.length, geo.downs,
                              &data, hsizes) < 0)
  {
    char *fullname = CCTK_FullName (index);


    CCTK_VWarn (1, __LINE__, __FILE__, CCTK_THORNSTRING,
                "Failed to extract hyperslab for variable '%s'", fullname);
    free (fullname);
    return (-1);
  }

  if (fid >= 0)
  {
    StreamedHDF5_Dump (GH, index, timelevel, data, geo.sdim, hsizes, hdf5_type, fid);
  }

  if (data)
    free (data);

  return (0);
}


/***************************** local functions ****************************/

static int GetHDF5type (StreamedHDF5GH *myGH, int vtype, hid_t *hdf5_type)
{
  int retval = 0;


  switch (vtype)
  {
    case CCTK_VARIABLE_CHAR:
      *hdf5_type = IOHDF5_CHAR;
      break;
    case CCTK_VARIABLE_INT:
      *hdf5_type = IOHDF5_INT;
      break;
    case CCTK_VARIABLE_REAL: 
      *hdf5_type = IOHDF5_REAL;
      break;
    case CCTK_VARIABLE_COMPLEX:
      *hdf5_type = myGH->IOHDF5_COMPLEX;
      break;
    default:
      retval = -1;
      break;
  }

  return (retval);
}


/*@@
  @routine    StreamedHDF5_AddCommonAttributes
  @date       May 21 1999
  @author     Thomas Radke
  @desc 
  Add "Common" attributes, these are:
  <ul>
    <li> the variable's groupname
    <li> the grouptype
    <li> number of timelevels
    <li> the current date
         Note that the datestamp should be turned of if you are byte 
         comparing two output files.
    <li> simulation time
    <li> origin
    <li> bounding box
    <li> gridspacings (both downsampled and evolution)
    <li> global grid size
  </ul>
  @enddesc 
  @calledby   StreamedHDF5_DumpGS, IOHDF5_collectiveDump, IOHDF5_procDump
  @history 
  @endhistory 
  @var        GH
  @vdesc      Pointer to CCTK grid hierarchy
  @vtype      cGH
  @vio        in
  @endvar
  @var        index
  @vdesc      global index of the variable to be dumped
  @vtype      int
  @vio        in
  @endvar
  @var        timelevel
  @vdesc      the timelevel to store
  @vtype      int
  @vio        in
  @endvar
  @var        dataset
  @vdesc      the dataset handle where the attributes should be attached to
  @vtype      hid_t
  @vio        in
  @endvar
@@*/
static void StreamedHDF5_AddCommonAttributes (cGH *GH, int index, int timelevel,
                                              hid_t dataset)
{
  DECLARE_CCTK_PARAMETERS
#if 0
  int dim;
  CCTK_INT intAttr;
  char *groupname;
#endif
  CCTK_REAL realAttr [6];
  StreamedHDF5GH *myGH;
  CCTK_REAL dummy;

  /* Get the handles for IOUtil and StreamedHDF5 extensions */
  myGH = (StreamedHDF5GH *) GH->extensions [CCTK_GHExtensionHandle ("StreamedHDF5")];

#if 0
  /* get the dimension of the variable */
  dim = CCTK_GroupDimI (CCTK_GroupIndexFromVarI (index));

  groupname = CCTK_GroupNameFromVarI (index);
  WRITE_ATTRIBUTE ("groupname", groupname,
                   dataset, myGH->scalarDataspace, 0, myGH->IOHDF5_STRING);
  free (groupname);
  intAttr = CCTK_GroupTypeFromVarI (index);
  WRITE_ATTRIBUTE ("grouptype", &intAttr, dataset, myGH->scalarDataspace, 0,
                   IOHDF5_INT);
  intAttr = CCTK_NumTimeLevelsFromVarI (index);
  WRITE_ATTRIBUTE ("ntimelevels", &intAttr, dataset, myGH->scalarDataspace, 0,
                   IOHDF5_INT);
#endif

  WRITE_ATTRIBUTE ("time", &GH->cctk_time, dataset, myGH->scalarDataspace, 0,
                   IOHDF5_REAL);

  /* NOTE: the attributes "origin", "min_ext", "max_ext", and "delta"
           are always stored as 3-dimensional points */
  CCTK_CoordRange (GH, &realAttr [0], &dummy, -1, "x", "cart3d");
  CCTK_CoordRange (GH, &realAttr [1], &dummy, -1, "y", "cart3d");
  CCTK_CoordRange (GH, &realAttr [2], &dummy, -1, "z", "cart3d");

  WRITE_ATTRIBUTE ("origin", realAttr, dataset, myGH->arrayDataspace, 3,
                   IOHDF5_REAL);

  CCTK_CoordRange (GH, &realAttr [0], &realAttr [3], -1, "x","cart3d");
  CCTK_CoordRange (GH, &realAttr [1], &realAttr [4], -1, "y","cart3d");
  CCTK_CoordRange (GH, &realAttr [2], &realAttr [5], -1, "z","cart3d");
  WRITE_ATTRIBUTE ("min_ext", realAttr, dataset, myGH->arrayDataspace, 3,
                   IOHDF5_REAL);
  WRITE_ATTRIBUTE ("max_ext", realAttr + 3, dataset, myGH->arrayDataspace, 3,
                   IOHDF5_REAL);
 
  WRITE_ATTRIBUTE ("delta", GH->cctk_delta_space, dataset, myGH->arrayDataspace,
                   3, IOHDF5_REAL);
}


/*@@
  @routine StreamedHDF5_procDump
  @author  Thomas Radke
  @date    May 21 1999
  @desc    All IO processors dump the data of their group into the file,
           either as one chunk per processor or unchunked
           (depending on parameter "unchunked").
  @enddesc
  @calledby   StreamedHDF5_DumpGA
  @history 
  @endhistory 
  @var        GH
  @vdesc      Pointer to CCTK grid hierarchy
  @vtype      cGH
  @vio        in
  @endvar
  @var        index
  @vdesc      global index of the variable to be dumped
  @vtype      int
  @vio        in
  @endvar
  @var        timelevel
  @vdesc      the timelevel to store
  @vtype      int
  @vio        in
  @endvar
  @var        outme
  @vdesc      pointer to the chunk to dump
  @vtype      void *
  @vio        in
  @endvar
  @var        geom
  @vdesc      bounds and sizes of the output
  @vtype      CCTK_INT [3*dim]
  @vio        in
  @vcomment   geom [0*dim..1*dim-1]: lower bounds
              geom [1*dim..2*dim-1]: (local) data sizes
              geom [2*dim..3*dim-1]: global space
  @endvar 
  @var        proc
  @vdesc      the processor which's chunk is to be dumped
  @vtype      int
  @vio        in
  @endvar
  @var        hdf5io_type
  @vdesc      the HDF5 datatype identifier
  @vtype      int
  @vio        in
  @endvar
  @var        fid
  @vdesc      the HDF5 file handle
  @vtype      hid_t
  @vio        in
  @endvar
@@*/
static void StreamedHDF5_Dump (cGH *GH, int index, int timelevel, void *outme,
                               int hdim, int *hsizes, int hdf5_type, hid_t fid)
{
  DECLARE_CCTK_PARAMETERS
  int i;
  hid_t dataset, memspace;
  char *name, *datasetname;
  hsize_t *dims;


  dims = (hsize_t *)  malloc (hdim * sizeof (hsize_t));

  /* HDF5 needs it in reverse order */
  for (i = 0; i < hdim; i++)
  {
    dims [i] = hsizes [hdim - 1 - i];
  }

  /* build the unique dataset name */
  name = CCTK_FullName (index);
  datasetname = (char *) malloc (strlen (name) + 30);  /* 30 chars should be
                                                          enough for 2 ints */
  sprintf (datasetname, "%s@%d@%d", name, GH->cctk_iteration, timelevel);
  free (name);

  /* create the memspace according to chunk dims */
  CACTUS_IOHDF5_ERROR (memspace = H5Screate_simple (hdim, dims, NULL));

  /* create the dataset */
  CACTUS_IOHDF5_ERROR (dataset = H5Dcreate (fid, datasetname, hdf5_type,
                                            memspace, H5P_DEFAULT));
  /* write the data */
  CACTUS_IOHDF5_ERROR (H5Dwrite (dataset, hdf5_type, H5S_ALL, H5S_ALL,
                                 H5P_DEFAULT, outme));

  StreamedHDF5_AddCommonAttributes (GH, index, timelevel, dataset);

  /* close the dataset and the memspace */
  CACTUS_IOHDF5_ERROR (H5Dclose (dataset));
  CACTUS_IOHDF5_ERROR (H5Sclose (memspace));

  free (datasetname);
  free (dims);
}