/*@@ @file hdf5_bitant_to_fullmode.c @date Thu 10 Jan 2002 @author Thomas Radke @desc This utility program copies a Cactus HDF5 datafile turning bitant mode datasets into full mode. @enddesc @version $Id$ @@*/ #include "cctk.h" #include #include #include #include #include "BetaThorns/IOHDF5Util/src/ioHDF5UtilGH.h" /* the rcs ID and its dummy function to use it */ static const char *rcsid = "$Header$"; CCTK_FILEVERSION(AlphaThorns_IOHDF5_hdf5_bitant_to_fullmode_c) /*****************************************************************************/ /* macro definitions */ /*****************************************************************************/ /* uncomment the following to get some debugging output */ /* #define IOHDF5_DEBUG 1 */ /* default stencil width if no '-stencil' option was given */ #define STENCIL 2 /* macro to do an HDF5 call, check its return code, and print a warning in case of an error */ #define CHECK_ERROR(hdf5_call) \ do \ { \ int _error_code = hdf5_call; \ \ \ if (_error_code < 0) \ { \ fprintf (stderr, "WARNING: line %d: HDF5 call '%s' returned " \ "error code %d\n", \ __LINE__, #hdf5_call, _error_code); \ nerrors++; \ } \ } while (0) /*****************************************************************************/ /* global variables */ /*****************************************************************************/ /* NOTE: although it isn't good programming practice we make these variables global for convenience since they are accessed from recursively or indirectly called routines which only get passed a single user-supplied argument */ static char *pathname = NULL; /* pathname of the current object */ static unsigned int nerrors = 0; /* global error counter */ static int stencil = STENCIL; /* stencil width for reflection */ /*****************************************************************************/ /* local function prototypes */ /*****************************************************************************/ static herr_t CopyObject (hid_t copy_from, const char *objectname, void *arg); static herr_t CopyAttribute (hid_t src, const char *attr_name, void *arg); static void FixBoundingBox (hid_t object); /*@@ @routine main @date Thu 10 Jan 2002 @author Thomas Radke @desc Main routine of the HDF5 bitant-to-full converter @enddesc @calls CopyObject @var argc @vdesc number of command line arguments @vtype int @vio in @endvar @var argv @vdesc command line arguments @vtype char *[] @vio in @endvar @returntype int @returndesc 0 for success, negative return values indicate an error @endreturndesc @@*/ int main (int argc, char *argv[]) { int i; hid_t infile, outfile; /* check for the '-stencil' option */ if (argc > 1 && strncmp (argv[1], "-stencil=", 9) == 0) { stencil = atoi (argv[1] + 9); for (i = 2; i < argc; i++) { argv[i-1] = argv[i]; } argc--; } /* give some help if called with incorrect number of parameters */ if (argc != 3) { fprintf (stderr, "Usage: %s [-stencil=] " "\n", argv[0]); fprintf (stderr, " eg, %s alp.h5 alp_fullmode.h5\n\n", argv[0]); return (0); } H5E_BEGIN_TRY { /* open the bitant mode file */ infile = H5Fopen (argv[1], H5F_ACC_RDONLY, H5P_DEFAULT); if (infile < 0) { fprintf (stderr, "ERROR: Cannot open HDF5 input file '%s' !\n\n",argv[1]); return (-1); } /* create output file */ outfile = H5Fcreate (argv[2], H5F_ACC_TRUNC, H5P_DEFAULT, H5P_DEFAULT); if (outfile < 0) { fprintf (stderr, "ERROR: Cannot create HDF5 output file '%s' !\n\n", argv[2]); return (-1); } } H5E_END_TRY printf ("\n ------------------------------------------\n" " Cactus 4 HDF5 bitant-to-fullmode Converter\n" " ------------------------------------------\n" " Using stencil width of %d for reflection...\n", stencil); /* do the copying by iterating over all objects */ pathname = ""; CHECK_ERROR (H5Giterate (infile, "/", NULL, CopyObject, &outfile)); /* finally, close all open files */ CHECK_ERROR (H5Fclose (infile)); CHECK_ERROR (H5Fclose (outfile)); /* report status */ if (nerrors == 0) { printf ("\n\n *** All Cactus data successfully converted. ***\n\n"); } else { fprintf (stderr, "\n\n *** WARNING: %d errors occured during " "conversion. ***\n\n", nerrors); } return (0); } /*****************************************************************************/ /* local routines */ /*****************************************************************************/ /*@@ @routine CopyObject @date Thu 10 Jan 2002 @author Thomas Radke @desc Iterator recursively called by H5Giterate() for every object in the input file It copies the current object or, if the object is a 3D dataset, reflects the data in the z-plane before copying. @enddesc @calls CopyAttribute @var from @vdesc identifier for the group the current object belongs to @vtype hid_t @vio in @endvar @var objectname @vdesc name of the current object @vtype const char * @vio in @endvar @var _to @vdesc user-supplied argument indicating the output object identifier @vtype hid_t @vio in @endvar @returntype int @returndesc 0 - continue the iteration for following group objects 1 - short-curcuit, no further iteration of this group @endreturndesc @@*/ static herr_t CopyObject (hid_t from, const char *objectname, void *_to) { hid_t to, datatype, dataspace_from, dataspace_to; H5G_stat_t objectinfo; char *current_pathname; int i, offset, xyplane_size, ndims, is_simple; size_t objectsize, datatypesize; hsize_t dims[4]; char *data; /* build the full pathname for the current to object to process */ current_pathname = pathname; pathname = (char *) malloc (strlen (current_pathname) + strlen (objectname) + 2); sprintf (pathname, "%s/%s", current_pathname, objectname); /* get the output object identifier */ to = *(hid_t *) _to; /* check the type of the current object */ CHECK_ERROR (H5Gget_objinfo (from, objectname, 0, &objectinfo)); if (objectinfo.type == H5G_GROUP) { printf (" copying group '%s'\n", pathname); CHECK_ERROR (from = H5Gopen (from, objectname)); CHECK_ERROR (to = H5Gcreate (to, objectname, 0)); /* iterate over all objects in the (first) input file */ CHECK_ERROR (H5Giterate (from, ".", NULL, CopyObject, &to)); CHECK_ERROR (H5Aiterate (from, NULL, CopyAttribute, &to)); FixBoundingBox (to); CHECK_ERROR (H5Gclose (to)); CHECK_ERROR (H5Gclose (from)); } else if (objectinfo.type == H5G_DATASET) { CHECK_ERROR (from = H5Dopen (from, objectname)); CHECK_ERROR (datatype = H5Dget_type (from)); CHECK_ERROR (dataspace_from = H5Dget_space (from)); CHECK_ERROR (is_simple = H5Sis_simple (dataspace_from)); CHECK_ERROR (ndims = H5Sget_simple_extent_ndims (dataspace_from)); if (is_simple && ndims == 3) { printf (" reflecting dataset '%s'\n", pathname); CHECK_ERROR (H5Sget_simple_extent_dims (dataspace_from, dims, NULL)); /* save the bitant zmax and compute the fullgrid zmax */ dims[3] = dims[0] - stencil; dims[0] = 2 * dims[3]; CHECK_ERROR (dataspace_to = H5Screate_simple (3, dims, NULL)); } else { printf (" copying dataset '%s'\n", pathname); CHECK_ERROR (dataspace_to = H5Scopy (dataspace_from)); } CHECK_ERROR (to = H5Dcreate (to, objectname, datatype, dataspace_to, H5P_DEFAULT)); datatypesize = objectsize = H5Tget_size (datatype); xyplane_size = 0; if (is_simple) { objectsize *= H5Sget_simple_extent_npoints (dataspace_from); } if (objectsize > 0) { if (is_simple && ndims == 3) { xyplane_size = dims[1] * dims[2] * datatypesize; offset = (dims[3] - stencil) * xyplane_size; objectsize += offset; } else { offset = 0; } data = (char *) malloc (objectsize); CHECK_ERROR (H5Dread (from, datatype, H5S_ALL, H5S_ALL, H5P_DEFAULT, data + offset)); if (is_simple && ndims == 3) { /* copy the +z xy planes into the -z xy planes */ for (i = stencil + 1; i <= (int) dims[3]; i++) { memcpy (data + offset - (i-stencil)*xyplane_size, data + offset + (i+stencil-1)*xyplane_size, xyplane_size); } } #if 0 /* flip the sign of the reflected elements */ if (flipsign) { if (H5Tequal (datatype, H5T_NATIVE_DOUBLE) > 0) { for (i = offset/sizeof (double) - 1; i >= 0; i--) { ((double *) data)[i] = -((double *) data)[i]; } } else { fprintf (stderr, "Warning: can only flip sign for datasets of " "datatype 'native double'\n"); } } #endif CHECK_ERROR (H5Dwrite (to, datatype, H5S_ALL, H5S_ALL, H5P_DEFAULT,data)); free (data); } CHECK_ERROR (H5Aiterate (from, NULL, CopyAttribute, &to)); if (is_simple && ndims == 3) { FixBoundingBox (to); } CHECK_ERROR (H5Dclose (to)); CHECK_ERROR (H5Dclose (from)); CHECK_ERROR (H5Sclose (dataspace_to)); CHECK_ERROR (H5Sclose (dataspace_from)); CHECK_ERROR (H5Tclose (datatype)); } else { fprintf (stderr, "WARNING: Found object '%s' which is not neither a " "group nor a dataset ! Object will not be copied.\n", pathname); nerrors++; } /* reset the pathname */ free (pathname); pathname = current_pathname; return (0); } /*@@ @routine CopyAttribute @date Thu 10 Jan 2002 @author Thomas Radke @desc Iterator recursively called by H5Aiterate() for every attribute of an object (dataset or group) @enddesc @var from @vdesc identifier for the group or dataset to read the attribute from @vtype hid_t @vio in @endvar @var attrname @vdesc name of the current attribute @vtype const char * @vio in @endvar @var _to @vdesc user-supplied argument indicating the group or dataset to copy the attribute to @vtype hid_t @vio in @endvar @returntype int @returndesc 0 - continue the iteration for following attributes @endreturndesc @@*/ static herr_t CopyAttribute (hid_t from, const char *attrname, void *_to) { hid_t attr, datatype, dataspace, to; size_t attrsize; void *value; /* get the target group/dataset identifier */ to = *(hid_t *) _to; /* open the attribute given by its name, get type, dataspace, and value and just copy it */ CHECK_ERROR (attr = H5Aopen_name (from, attrname)); CHECK_ERROR (datatype = H5Aget_type (attr)); CHECK_ERROR (dataspace = H5Aget_space (attr)); attrsize = H5Tget_size (datatype); if (H5Sis_simple (dataspace) > 0) { attrsize *= H5Sget_simple_extent_npoints (dataspace); } if (attrsize > 0) { value = malloc (attrsize); CHECK_ERROR (H5Aread (attr, datatype, value)); CHECK_ERROR (H5Aclose (attr)); CHECK_ERROR (attr = H5Acreate (to, attrname, datatype, dataspace, H5P_DEFAULT)); CHECK_ERROR (H5Awrite (attr, datatype, value)); free (value); } CHECK_ERROR (H5Aclose (attr)); CHECK_ERROR (H5Sclose (dataspace)); CHECK_ERROR (H5Tclose (datatype)); return (0); } /*@@ @routine FixBoundingBox @date Mon 18 Feb 2002 @author Thomas Radke @desc Fix the "global_size", "min_ext", and "origin" attributes if found in the object (dataset or group) given. @enddesc @var object @vdesc identifier for the group or dataset to fix the attributes @vtype hid_t @vio inout @endvar @@*/ static void FixBoundingBox (hid_t object) { hid_t attr, datatype, dataspace; double delta[3], min_ext[3]; int global_size[3], do_it; /* try to read the "delta" and "global_size" attributes */ H5E_BEGIN_TRY { attr = H5Aopen_name (object, "delta"); do_it = attr >= 0; if (do_it) { CHECK_ERROR (datatype = H5Aget_type (attr)); CHECK_ERROR (dataspace = H5Aget_space (attr)); if (H5Tget_class (datatype) == H5T_FLOAT && H5Sis_simple (dataspace) > 0 && H5Sget_simple_extent_npoints (dataspace) == 3) { CHECK_ERROR (H5Aread (attr, H5T_NATIVE_DOUBLE, delta)); } else { do_it = 0; } CHECK_ERROR (H5Aclose (attr)); CHECK_ERROR (H5Sclose (dataspace)); CHECK_ERROR (H5Tclose (datatype)); } if (do_it) { attr = H5Aopen_name (object, "global_size"); do_it = attr >= 0; } if (do_it) { CHECK_ERROR (datatype = H5Aget_type (attr)); CHECK_ERROR (dataspace = H5Aget_space (attr)); if (H5Tget_class (datatype) == H5T_INTEGER && H5Sis_simple (dataspace) > 0 && H5Sget_simple_extent_npoints (dataspace) == 3) { CHECK_ERROR (H5Aread (attr, H5T_NATIVE_INT, global_size)); global_size[2] += global_size[2] - 1; CHECK_ERROR (H5Awrite (attr, H5T_NATIVE_INT, global_size)); } else { do_it = 0; } CHECK_ERROR (H5Aclose (attr)); CHECK_ERROR (H5Sclose (dataspace)); CHECK_ERROR (H5Tclose (datatype)); } } H5E_END_TRY; /* now fix the "min_ext" and "origin" attributes */ if (do_it) { CHECK_ERROR (attr = H5Aopen_name (object, "min_ext")); if (attr >= 0) { CHECK_ERROR (datatype = H5Aget_type (attr)); CHECK_ERROR (dataspace = H5Aget_space (attr)); if (H5Tget_class (datatype) == H5T_FLOAT && H5Sis_simple (dataspace) > 0 && H5Sget_simple_extent_npoints (dataspace) == 3) { CHECK_ERROR (H5Aread (attr, H5T_NATIVE_DOUBLE, min_ext)); min_ext[2] -= delta[2] * ((global_size[2] + 1) / 2); CHECK_ERROR (H5Awrite (attr, H5T_NATIVE_DOUBLE, min_ext)); } CHECK_ERROR (H5Aclose (attr)); CHECK_ERROR (H5Sclose (dataspace)); CHECK_ERROR (H5Tclose (datatype)); } CHECK_ERROR (attr = H5Aopen_name (object, "origin")); if (attr >= 0) { CHECK_ERROR (datatype = H5Aget_type (attr)); CHECK_ERROR (dataspace = H5Aget_space (attr)); if (H5Tget_class (datatype) == H5T_FLOAT && H5Sis_simple (dataspace) > 0 && H5Sget_simple_extent_npoints (dataspace) == 3) { CHECK_ERROR (H5Aread (attr, H5T_NATIVE_DOUBLE, min_ext)); min_ext[2] -= delta[2] * ((global_size[2] + 1) / 2); CHECK_ERROR (H5Awrite (attr, H5T_NATIVE_DOUBLE, min_ext)); } CHECK_ERROR (H5Aclose (attr)); CHECK_ERROR (H5Sclose (dataspace)); CHECK_ERROR (H5Tclose (datatype)); } } }