src/hdf5/parallel/H5Part.cc

Go to the documentation of this file.

#include <stdio.h>
#include <stdlib.h>
#include <hdf5.h>
#include "H5Part.hh"

/********* Private Function Declarations *************/
herr_t H5PartIOcounter(hid_t group_id,
                   const char *member_name,
                   void *operator_data);


/********** Definitions of Functions/Subroutines ******/

/*========== File Opening/Closing ===============*/
/*
  H5PartOpenFile:  Opens file with specified filename. 
    If you open with flag H5PART_WRITE, it will truncate any
    file with the specified filename and start writing to it.
    If you open with H5PART_READ, then it will open the file
    readonly.  I could probably do an APPEND mode as well, but
    just haven't done so yet.  APPEND would need to check the
    file to determine its current state and modify the 
    H5PartFile datastructure accordingly.

    H5PartFile should be treated as an essentially opaque
    datastructure.  It acts as the file handle, but internally
    it maintains several key state variables associated with 
    the file.
*/

H5PartFile *H5PartOpenFileParallel(const char *filename,unsigned flags
#ifdef PARALLEL_IO
                ,                  MPI_Comm comm
#endif
){
  H5PartFile *f=(H5PartFile*)malloc(sizeof(H5PartFile));
#ifdef PARALLEL_IO
  MPI_Info info = MPI_INFO_NULL; /* for the SP2... perhaps different for linux */
#endif
  f->xfer_prop = f->create_prop = f->access_prop=H5P_DEFAULT;
#ifdef PARALLEL_IO
  MPI_Comm_size(comm,&f->nprocs);
  MPI_Comm_rank(comm,&f->myproc);
  f->pnparticles=(long long*)malloc(sizeof(long long)*f->nprocs);
  /* access_prop: properties like chunking or parallel I/O */
  f->access_prop = H5Pcreate(H5P_FILE_ACCESS);
  if(H5Pset_fapl_mpio(f->access_prop,comm,info)<0){
        fprintf(stderr,"Total failure trying to setup mpi-io for access\n");
        exit(0);
  }
  /* create_prop: tunable parameters like blocksize and btree sizes */
  /* f->create_prop = H5Pcreate(H5P_FILE_CREATE); */
  f->create_prop = H5P_DEFAULT;
  /* currently create_prop is empty */
  /* xfer_prop:  also used for parallel I/O, during actual writes
          rather than the access_prop which is for file creation. */
  f->xfer_prop = H5Pcreate(H5P_DATASET_XFER);
  H5Pset_dxpl_mpio(f->xfer_prop,H5FD_MPIO_COLLECTIVE);
  f->comm = comm;
#endif
  if(flags==H5PART_READ){
    /*  f->file=H5Fopen(filename,H5F_ACC_RDONLY,f->access_prop); */
    /* just do this serially */
    f->file = H5Fopen(filename,H5F_ACC_RDONLY,H5P_DEFAULT);
  }
  else if(flags==H5PART_WRITE){
    f->file=H5Fcreate(filename,H5F_ACC_TRUNC,f->create_prop,f->access_prop);
  }
  else {
    f->file=-1;
    fprintf(stderr,"H5Open: Invalid file access type %d\n",flags);
  }
  if(f->file<=0) {
    free(f);
    fprintf(stderr,"H5Open: File open failed for file [%s]\n",
            filename);
        exit(0); 
        return 0;
  }
  f->mode=flags;
  f->maxstep=0;
  f->timegroup=0;
  f->timestep=0;
  f->shape=0;
  f->diskshape=H5S_ALL;
  f->memshape=H5S_ALL;
  return f;
}

H5PartFile *H5PartOpenFile(const char *filename,unsigned flags){
  H5PartFile *f=(H5PartFile*)malloc(sizeof(H5PartFile));
  f->xfer_prop = f->create_prop = f->access_prop=H5P_DEFAULT;
#ifdef PARALLEL_IO
  f->pnparticles=0;
  f->comm = MPI_COMM_WORLD;
  f->nprocs=1;
  f->myproc=0;
#endif
  if(flags==H5PART_READ){
    /*  f->file=H5Fopen(filename,H5F_ACC_RDONLY,f->access_prop); */
    /* just do this serially */
    f->file = H5Fopen(filename,H5F_ACC_RDONLY,H5P_DEFAULT);
  }
  else if(flags==H5PART_WRITE){
    f->file=H5Fcreate(filename,H5F_ACC_TRUNC,H5P_DEFAULT,H5P_DEFAULT);
  }
  else {
    f->file=-1;
    fprintf(stderr,"H5Open: Invalid file access type %d\n",flags);
  }
  if(f->file<0) {
    free(f);
    fprintf(stderr,"H5Open: File open failed for file [%s]\n",
            filename);
    return 0;
  }
  f->mode=flags;
  f->maxstep=0;
  f->timegroup=0;
  f->timestep=0;
  f->shape=0;
  f->diskshape=H5S_ALL;
  f->memshape=H5S_ALL;
  return f;
}
/*
  H5PartCloseFile:  Does what it says it does.
*/
void H5PartCloseFile(H5PartFile *f){
  if(f->shape>0) H5Sclose(f->shape); f->shape=0;
  if(f->timegroup>0) H5Gclose(f->timegroup); f->timegroup=0;
  if(f->diskshape!=H5S_ALL) H5Sclose(f->diskshape);
#ifdef PARALLEL_IO
  if(f->pnparticles) free(f->pnparticles);
#endif
  if(f->xfer_prop!=H5P_DEFAULT) {
    H5Pclose(f->xfer_prop); f->xfer_prop=H5P_DEFAULT;
  }
  if(f->access_prop!=H5P_DEFAULT) {
    H5Pclose(f->access_prop); f->access_prop=H5P_DEFAULT;
  }  
  if(f->create_prop!=H5P_DEFAULT) {
    H5Pclose(f->create_prop); f->create_prop=H5P_DEFAULT;
  }
  if(f->file) H5Fclose(f->file); f->file=0;
  free(f);
}


/*============== File Writing Functions ==================== */
/*
  H5PartSetNumParticles:  This is a private function used by the 
  other functions that write arrays.  Its sole purpose is to 
  prevent needless creation of new HDF5 DataSpace handles if
  the number of particles is invariant throughout the sim.
  That's its only reason for existence.
*/
void H5PartSetNumParticles(H5PartFile *f,long long nparticles){
  hssize_t start[1];
  hsize_t stride[1],count[1],total;
  register int i,r;
  /*  printf("Set Number of Particles to %lld\n",nparticles); */
  if(f->nparticles==nparticles) return; // no change
  if(f->diskshape!=H5S_ALL) H5Sclose(f->diskshape);
  if(f->memshape!=H5S_ALL) H5Sclose(f->memshape);
  f->memshape=f->diskshape = H5S_ALL;
  if(f->shape) H5Sclose(f->shape);
  f->nparticles=(hsize_t)nparticles;
#ifndef PARALLEL_IO
  f->shape=H5Screate_simple(1, /* 1 dimensional shape */
                            &(f->nparticles), /* with nparticles elements */
                            NULL); /* ignore this :-) */
#else /* PARALLEL_IO */
  /* The Gameplan here is to declare the overall size of the on-disk
     data structure the same way we do for the serial case.  But
     then we must have additional "DataSpace" structures to define
     our in-memory layout of our domain-decomposed portion of the particle
     list as well as a "selection" of a subset of the on-disk 
     data layout that will be written in parallel to mutually exclusive
     regions by all of the processors during a parallel I/O operation.
     These are f->shape, f->memshape and f->diskshape respectively. */
  /* acquire the number of particles to be written from each MPI process */
  MPI_Allgather(&nparticles,1,MPI_LONG_LONG,f->pnparticles,1,MPI_LONG_LONG,f->comm);

  //  if(f->myproc==0){
  //    puts("AllGather:  Particle offsets:\n");
  //    for(i=0;i<f->nprocs;i++) 
  //            printf("\tnp=%u\n",f->pnparticles[i]); /* compute total nparticles */
  //}
  /* should I create a selection here? */
  stride[0]=1;
  start[0]=0;
  for(i=0;i<f->myproc;i++) start[0]+=f->pnparticles[i]; /* compute start offsets */
  total=0;
  for(i=0;i<f->nprocs;i++) total+=f->pnparticles[i]; /* compute total nparticles */
  f->shape = H5Screate_simple(1,&total,NULL); /* declare overall datasize */
  f->diskshape = H5Screate_simple(1,&total,NULL); /* declare overall data size  but then will select a subset */
  f->memshape = H5Screate_simple(1,&(f->nparticles),NULL); /* declare local memory datasize */

  if(f->shape<0 || f->memshape<0 || f->diskshape<0){
        fprintf(stderr,"Abort:  shape construction failed\n");
        exit(0);
}
  count[0]=nparticles; /* based on local nparticles (for the selection */
  /* and then set the subset of the data you will write to */
  r=H5Sselect_hyperslab(f->diskshape,H5S_SELECT_SET,start,stride,count,NULL);
  if(r<0){
        fprintf(stderr,"Abort: Selection Failed!\n");
        exit(0);
  }
#endif
}

/*
  H5PartWriteArray:  This writes an individual array of data
    for a given timestep.  It probably would be convenient if this
    were public, but its private because it there are some order/
    state dependencies that get kind of nasty to manage (and therefore
    are not managed).  It is otherwise assuming that the timegroup
    has already been created and other such details.
    It is called by the mongo H5PartWriteStep() routine.
 */
int H5PartWriteDataFloat64(H5PartFile *f,char *name,double *array){
  register int r;
  hid_t dataset;
  //fprintf(stderr,"Create a dataset[%s]  mounted on the timegroup %u\n",
  //     name,f->timestep);
  if(f->timegroup<=0) {fprintf(stderr,"Eeeerrrrroooorrrr!!!!\n");}
  dataset = H5Dcreate(f->timegroup,name,H5T_NATIVE_DOUBLE,f->shape,H5P_DEFAULT);
  if(dataset<0)
    fprintf(stderr,"Dataset open failed for name=[%s] step=%u!\n",
            name,f->timestep);
  r=H5Dwrite(dataset,H5T_NATIVE_DOUBLE,f->memshape,f->diskshape,H5P_DEFAULT,array);
  H5Dclose(dataset);
  return r;
}
int H5PartWriteDataInt64(H5PartFile *f,char *name,long long *array){
  register int r;
  hid_t dataset;
  //fprintf(stderr,"Create a dataset[%s]  mounted on the timegroup %u\n",
  //     name,f->timestep);
  if(f->timegroup<=0) {fprintf(stderr,"Eeeerrrrroooorrrr!!!!\n");}
  dataset = H5Dcreate(f->timegroup,name,H5T_NATIVE_INT64,f->shape,H5P_DEFAULT);
  if(dataset<0){
    fprintf(stderr,"Dataset open failed for name=[%s] step=%u!\n",
            name,f->timestep);
    exit(0);
  }
  r=H5Dwrite(dataset,H5T_NATIVE_INT64,f->memshape,f->diskshape,H5P_DEFAULT,array);
  if(r<0) {
        fprintf(stderr,"Attempt to write dataset failed!\n");
        exit(0);
  }
  H5Dclose(dataset);
  return r;
}

int H5PartWriteStepAttrib(H5PartFile *f,char *name,
                          hid_t type,void *value,int nelem){
  register int r;
  hid_t attrib;
  hid_t space;
  hsize_t len;
  //fprintf(stderr,"Create a attribute[%s]  mounted on the timegroup %u\n",
  //     name,f->timestep);
  if(f->timegroup<=0) {fprintf(stderr,"Eeeerrrrroooorrrr!!!!\n");}
  len = nelem;
  space = H5Screate_simple(1,&len,NULL);
  attrib = H5Acreate(f->timegroup,name,type,space,H5P_DEFAULT);
  r=H5Awrite(attrib,type,value);
  H5Aclose(attrib);
  H5Sclose(space);
  return r;
}

int H5PartWriteAttrib(H5PartFile *f,char *name,
    hid_t type,void *value,int nelem){
    return H5PartWriteStepAttrib(f,name,type,value,nelem);
}

int H5PartWriteFileAttrib(H5PartFile *f,char *name,
                          hid_t type,void *value,int nelem){
  register int r;
  hid_t attrib;
  hid_t space;
  hsize_t len;
  //fprintf(stderr,"Create a attribute[%s]  mounted on the timegroup %u\n",
  //     name,f->timestep);
  if(f->file<=0) {fprintf(stderr,"Eeeerrrrroooorrrr!!!!\n");}
  len = nelem;
  space = H5Screate_simple(1,&len,NULL);
  attrib = H5Acreate(f->file,name,type,space,H5P_DEFAULT);
  r=H5Awrite(attrib,type,value);
  H5Aclose(attrib);
  H5Sclose(space);
  return r;
}
herr_t H5PartAttribcounter(hid_t group_id,
                   const char *member_name,
                   void *operator_data){
  int *count = (int*)operator_data;
  (*count)++;
  return 0;
}
int H5PartGetNumStepAttribs(H5PartFile *f){ /* for current filestep */
  return H5Aget_num_attrs(f->timegroup);
}
int H5PartGetNumFileAttribs(H5PartFile *f){
  /* must walk the attributes */
  /* iterate to get numsteps */
  int nattribs;
  hid_t rootgroup=H5Gopen(f->file,"/");;
  nattribs=H5Aget_num_attrs(rootgroup); /* open / for the file?
                                       or is there a problem with 
                                       file attributes? */
  H5Gclose(rootgroup);
  return nattribs;
}

hid_t H5PartNormType(hid_t type){
  H5T_class_t tclass = H5Tget_class(type);
  int size = H5Tget_size(type);
  switch(tclass){
  case H5T_INTEGER:
    if(size==8) return H5T_NATIVE_INT64;
    else if(size==1) return H5T_NATIVE_CHAR;
    else fprintf(stderr,"Unknown type %u.  Cannot infer normalized type\n",type);
    break;
  case H5T_FLOAT:
    return H5T_NATIVE_DOUBLE;
    break;
  }
  return -1;
}

void H5PartGetStepAttribInfo(H5PartFile *f,int idx,
                         char *name,size_t maxname,hid_t *type,int *nelem){
  hid_t attrib=H5Aopen_idx(f->timegroup,idx);
  hid_t mytype=H5Aget_type(attrib);
  hid_t space = H5Aget_space(attrib);
  if(nelem)
    *nelem=H5Sget_simple_extent_npoints(space);
  if(name)
    H5Aget_name(attrib,maxname,name);
  if(type)
    *type=H5PartNormType(mytype); /* normalize the type to be machine-native
                                 this means one of 
                                H5T_NATIVE_DOUBLE
                                H5T_NATIVE_INT64
                                H5T_NATIVE_CHAR */
  H5Sclose(space);
  H5Tclose(mytype);
  H5Aclose(attrib);
}
void H5PartGetFileAttribInfo(H5PartFile *f,int idx,
                         char *name,size_t maxname,hid_t *type,int *nelem){
  hid_t attrib=H5Aopen_idx(f->file,idx);
  hid_t mytype=H5Aget_type(attrib);
  hid_t space = H5Aget_space(attrib);
  if(nelem)
    *nelem=H5Sget_simple_extent_npoints(space);
  if(name)
    H5Aget_name(attrib,maxname,name);
  if(type)
    *type=H5PartNormType(mytype); /* normalize the type to be machine-native
                                 this means one of 
                                H5T_NATIVE_DOUBLE
                                H5T_NATIVE_INT64
                                H5T_NATIVE_CHAR */
  H5Sclose(space);
  H5Tclose(mytype);
  H5Aclose(attrib);
}
void H5ReadStepAttrib(H5PartFile *f,char *name,void *data){
  /* use the open attribute by name mode of operation */
  hid_t attrib=H5Aopen_name(f->timegroup,name);
  hid_t mytype=H5Aget_type(attrib);
  hid_t space = H5Aget_space(attrib);
  hsize_t nelem=H5Sget_simple_extent_npoints(space);
  hid_t type=H5PartNormType(mytype); /* normalize the type to be machine-native
                                 this means one of 
                                H5T_NATIVE_DOUBLE
                                H5T_NATIVE_INT64
                                H5T_NATIVE_CHAR */
  H5Aread(attrib,type,data);
  H5Sclose(space);
  H5Tclose(mytype);
  H5Aclose(attrib);
}

void H5PartReadAttrib(H5PartFile *f,char *name,void *data){
    /* use the open attribute by name mode of operation */
    H5ReadStepAttrib(f,name,data);
}

void H5ReadFileAttrib(H5PartFile *f,char *name,void *data){
  hid_t attrib=H5Aopen_name(f->file,name);
  hid_t mytype=H5Aget_type(attrib);
  hid_t space = H5Aget_space(attrib);
  hsize_t nelem=H5Sget_simple_extent_npoints(space);
  hid_t type=H5PartNormType(mytype); /* normalize the type to be machine-native
                                 this means one of 
                                H5T_NATIVE_DOUBLE
                                H5T_NATIVE_INT64
                                H5T_NATIVE_CHAR */
  H5Aread(attrib,type,data);
  H5Sclose(space);
  H5Tclose(mytype);
  H5Aclose(attrib);
}


/*================== File Reading Routines =================*/
/*
  H5PartSetStep:  This routine is *only* useful when you are reading
  data.  Using it while you are writing will generate nonsense results!
  (This file format is only half-baked... robustness is not std equipment!)
  So you use this to random-access the file for a particular timestep.
  Failure to explicitly set the timestep on each read will leave you
  stuck on the same timestep for *all* of your reads.  That is to say
  the writes auto-advance the file pointer, but the reads do not
  (they require explicit advancing by selecting a particular timestep).
*/
void H5PartSetStep(H5PartFile *f,int step){
  char name[128];
  /*printf("Proc[%u] SetStep to %u for file %u\n",f->myproc,step,f->file);  */
  f->timestep=step;
  if(f->timegroup>0) {
        H5Gclose(f->timegroup); 
  }
  f->timegroup=-1;
  sprintf(name,"Particles#%u",f->timestep);
  /*fprintf(stderr,"P[%u] Open timegroup [%s] for file %d\n",
        f->myproc,name,f->file);*/
  if(f->mode==H5PART_READ) { /* kludge to prevent error messages */
/*      fprintf(stderr,"P[%u] open group\n"); */
    f->timegroup=H5Gopen(f->file,name);
  }
  if(f->timegroup<=0){
    /* timegroup doesn't exist, so we need to create one */
        /* fprintf(stderr,"P[%u] create group\n"); */
    f->timegroup = H5Gcreate(f->file,name,0);
    if(f->timegroup<0) 
        fprintf(stderr,"p[%u] timegroup creation failed!\n",
                f->myproc);
  }
}

/*
  H5PartIOcounter:  This is an entirely internal callback function 
   which is used in conjunction with HDF5 iterators.  The HDF5 Group
   iterator will call this repeatedly in order to count how many
   timesteps of data have been stored in a particular file.
   This is used by H5PartGetNumSteps().
*/
herr_t H5PartIOcounter(hid_t group_id,
                   const char *member_name,
                   void *operator_data){
  int *count = (int*)operator_data;
  H5G_stat_t objinfo;
  /* only count the particle groups... ignore all others */
  if(!strncmp(member_name,"Particles",9)) (*count)++;
  return 0;
}
herr_t H5PartDScounter(hid_t group_id,
                   const char *member_name,
                   void *operator_data){
  int *count = (int*)operator_data;
  H5G_stat_t objinfo;
  /* only count the datasets... ignore all others */
  if(H5Gget_objinfo(group_id, 
                 member_name,
                 1 /* follow links */, 
                    &objinfo)<0) {
    return 0; /* error (probably bad symlink) */
  }
  if(objinfo.type==H5G_DATASET)    (*count)++;
  return 0;
}
struct H5IO_getname_t {
  int index,count;
  char *name;
};
herr_t H5IOgetname(hid_t group_id,
                          const char *member_name, 
                          void *operator_data){
  H5IO_getname_t *getn = (H5IO_getname_t*)operator_data;
  // check type first (only respond if it is a dataset)
  H5G_stat_t objinfo;
  // request info about the type of objects in root group
  if(H5Gget_objinfo(group_id, 
                    member_name,
                    1 /* follow links */, 
                    &objinfo)<0) return 0; // error (probably bad symlink)
  // only count objects that are datasets (not subgroups)
  if(objinfo.type!=H5G_DATASET) 
    return 0; // do not increment count if it isn't a dataset.
  if(getn->index==getn->count){
    strcpy(getn->name,member_name);
    return 1; // success
  }
  getn->count++;
  return 0;
}
/* 
   H5PartGetNumSteps:  This reads the number of datasteps that are 
     currently stored in the datafile.  (simple return of an int).
     It works for both reading and writing of files, but is probably
     only typically used when you are reading.
*/
int H5PartGetNumSteps(H5PartFile *f){
  /* iterate to get numsteps */
  int count=0;
  int idx=0;
  while(H5Giterate(f->file, /* hid_t loc_id, */
                   "/", /*const char *name, */
                   &idx, /* int *idx, */
                   H5PartIOcounter,
                   &count)<0){}  
  return count;
}
int H5PartGetNumDatasets(H5PartFile *f){
  int count=0;
  int idx=0;
  char name[128];
  sprintf(name,"Particles#%u",f->timestep);
  while(H5Giterate(f->file, /* hid_t loc_id */
                   name,
                   &idx,
                   H5PartDScounter,
                   &count)<0){}
  return count;   
}
int H5PartGetDatasetName(H5PartFile *f,
                         int _index, /* index of the dataset (0 to N-1) */
                         char *name, /* buffer to store name of dataset */
                         int maxlen){ /* max size of the "name" buffer */
  H5IO_getname_t getn;
  int idx=_index;
  char gname[128];
  sprintf(gname,"Particles#%u",f->timestep);
  getn.index=_index; getn.name=name; getn.count=_index;
  while(H5Giterate(f->file, /* hid_t loc_id, */
                   gname, /*const char *name, */
                   &idx, /* int *idx, */
                   H5IOgetname,
                   &getn)<0){}
  return 1;
}

/*
  H5PartGetNumParticles:  This gets the number of particles 
    that are stored in the current timestep.  It will arbitrarily
    select a timestep if you haven't already set the timestep
    with H5PartSetStep().
 */
long long H5PartGetNumParticles(H5PartFile *f){
  hid_t space,dataset;
  hsize_t nparticles;
  if(f->timegroup<=0) H5PartSetStep(f,f->timestep); /* choose a step */
  /* lets open up a dataset in this group... we know "x" is there */
  dataset=H5Dopen(f->timegroup,"x");/* choice of X is a kludge of sorts */
  space = H5Dget_space(dataset);
  nparticles= H5Sget_simple_extent_npoints(space);
  H5Sclose(space); /* release data shape handle */
  H5Dclose(dataset); /* release the dataset handle */
  return (long long)nparticles;
}


/*
  H5PartReadArray:  This reads in an individual array from a
    particlar timestep.  Unlike its "Write" mode sibling, this routine
    is completely public and will not have any wacky state dependencies.
    If you haven't selected a particular timestep, it will pick
    the current one.
*/
int H5PartReadDataFloat64(H5PartFile *f,char *name,double *array){
  hid_t space,dataset,datatype;
  if(!f->timegroup) H5PartSetStep(f,f->timestep); /* choose a step */
  dataset=H5Dopen(f->timegroup,name);
  space = H5Dget_space(dataset);
  /*  datatype=H5Dget_type(dataset); */
  H5Dread(dataset, /* handle for the dataset */
          H5T_NATIVE_DOUBLE, /* the datatype we use in memory 
                              you can change it to FLOAT if you want */
          H5S_ALL, /* shape/size of data in memory (currently get all) */
          H5S_ALL, /* shape/size of data on disk  (currently get all) */
          H5P_DEFAULT,/* ignore... its for parallel reads */
          array); /* the data array we are reading into */
  H5Sclose(space); /* release data shape handle */
  H5Dclose(dataset); /* release the dataset handle */
  return 1;
}

int H5PartReadDataInt64(H5PartFile *f,char *name,long long *array){
  hid_t space,dataset,datatype;
  if(!f->timegroup) H5PartSetStep(f,f->timestep); /* choose a step */
  dataset=H5Dopen(f->timegroup,name);
  space = H5Dget_space(dataset);
  /*  datatype=H5Dget_type(dataset); */
  H5Dread(dataset, /* handle for the dataset */
          H5T_NATIVE_INT64, /* the datatype we use in memory 
                              you can change it to FLOAT if you want */
          H5S_ALL, /* shape/size of data in memory (currently get all) */
          H5S_ALL, /* shape/size of data on disk  (currently get all) */
          H5P_DEFAULT,/* ignore... its for parallel reads */
          array); /* the data array we are reading into */
  H5Sclose(space); /* release data shape handle */
  H5Dclose(dataset); /* release the dataset handle */
  return 1; /* totally bogus */
}

/*
  H5PartReadParticleStep:  This is the mongo read function that
    pulls in all of the data for a given timestep in one shot.
    It also takes the timestep as an argument and will call
    H5PartSetStep() internally so that you don't have to 
    make that call separately.
    See also: H5PartReadArray() if you want to just
    read in one of the many arrays.
*/
int H5PartReadParticleStep(H5PartFile *f,
                           int step,
                           double *x,double *y,double *z,
                           long long *id){
  H5PartSetStep(f,step);
  /* or smuggle it into the array names */
  H5PartReadDataFloat64(f,"x",x);
  H5PartReadDataFloat64(f,"y",y);
  H5PartReadDataFloat64(f,"z",z);
  H5PartReadDataInt64(f,"id",id);
  
  return 1;
}


#ifdef REGRESSIONTEST

/*
  A simple regression test that shows how you use this API
  to write and read multi-timestep files of particle data.
*/
#ifdef PARALLEL_IO

int main(int argc,char *argv[]){
  int sz=5;
  double *x,*y,*z;
  long long *id;
  char name[64];
  H5PartFile *file;
  int i,t,nt,nds;
  int nprocs,myproc;
  hid_t gid;
  MPI_Comm comm=MPI_COMM_WORLD;
  
  MPI_Init(&argc,&argv);
  MPI_Comm_size(comm,&nprocs);
  MPI_Comm_rank(comm,&myproc);

  x=(double*)malloc(sz*nprocs*sizeof(double));
  y=(double*)malloc(sz*nprocs*sizeof(double));
  z=(double*)malloc(sz*nprocs*sizeof(double));
  id=(long long*)malloc(sz*nprocs*sizeof(long long));
  /* parallel file creation */
  file=H5PartOpenFileParallel("parttest.h5",H5PART_WRITE,comm);
  if(!file) {
    perror("File open failed:  exiting!");
    exit(0);
  }
  
  for(t=0;t<5;t++){
    MPI_Barrier(comm);
    for(i=0;i<sz;i++) {
      x[i]=(double)(i+t)+10.0*(double)myproc;
      y[i]=0.1 + (double)(i+t);
      z[i]=0.2 + (double)(i+t*10);
      id[i]=i+sz*myproc;
    }
    printf("Proc[%u] Writing timestep %u file=%u\n",myproc,t,file->file);
    H5PartSetStep(file,t); /* must set the current timestep in file */
    H5PartSetNumParticles(file,sz); /* then set number of particles to store */
    /* now write different tuples of data into this timestep of the file */
    H5PartWriteDataFloat64(file,"x",x); 
    H5PartWriteDataFloat64(file,"y",y);
    H5PartWriteDataFloat64(file,"z",z);
    H5PartWriteDataInt64(file,"id",id);
  }
  printf("AllDone p[%u]\n",myproc);
  H5PartCloseFile(file);
    MPI_Barrier(comm);
  printf("p[%u:%u] : OK, close file and reopen for reading \n",myproc,nprocs);
  if(myproc==0){ /* now only proc 0 reads the file serially */
    file= H5PartOpenFileSerial("parttest.h5",H5PART_READ);
    nt=H5PartGetNumSteps(file); /* get number of steps in file */
    H5PartSetStep(file,0);
    nds=H5PartGetNumDatasets(file); /* get number of datasets in timestep 0 */
    
    puts("\n\n===============================");
    for(i=0;i<nds;i++){ /* and print out those names */
      H5PartGetDatasetName(file,i,name,64);
      printf("\tDataset[%u] name=[%s]\n",
             i,name);
    }
    puts("===============================\n\n");
    printf("Number of datasteps in file is %u\n",nt);
    for(t=0;t<nt;t++){
      int nparticles;
      H5PartSetStep(file,t); /* select a timestep */
      nparticles=(int)H5PartGetNumParticles(file);
      printf("Step[%u] nparticles this step=%u\n",
             t,nparticles); /* get num particles this step */
      H5PartReadParticleStep(file,t,/* do a mongo read of all data this step */
                             x,y,z,id);
      printf("\tid\t\tx\t\ty\t\tz\n");
      puts("\t----------------------------------------------------");
      for(i=0;i<nparticles;i++) {
        printf("\t%llu\t%f\t%f\t%f\n\n",id[i],x[i],y[i],z[i]);
      }
    }
    H5PartCloseFile(file);
  }
  if(x) free(x); 
  if(y) free(y);
  if(z) free(z);
  if(id) free(id);
  MPI_Barrier(comm);
  fprintf(stderr,"proc[%u]:  done\n",myproc);
  return MPI_Finalize();
}

#else
  
  int main(int argc,char *argv){
    const int sz=5;
    double x[sz],y[sz],z[sz];
  long long id[sz];
  char name[64];
  H5PartFile *file;
  int i,t,nt,nds,myproc;
  int nfattribs,nsattribs;
  

  if (argc==1) {

    file=H5PartOpenFile("parttest.h5",H5PART_WRITE);
    if(!file) {
      perror("File open failed:  exiting!");
      exit(0);
    }
    for(t=0;t<5;t++){
      long long step=t;
      printf("Writing timestep %u\n",t);
      for(i=0;i<sz;i++) {
        x[i]=(double)(i+t);
        y[i]=0.1 + (double)(i+t);
        z[i]=0.2 + (double)(i+t);
        id[i]=i;
      }
      H5PartSetStep(file,t); /* must set the current timestep in file */
      H5PartSetNumParticles(file,sz); /* then set number of particles to store */
      /* now write different tuples of data into this timestep of the file */
      H5PartWriteDataFloat64(file,"x",x); 
      H5PartWriteDataFloat64(file,"y",y);
      H5PartWriteDataFloat64(file,"z",z);
      H5PartWriteDataInt64(file,"id",id);
      H5PartWriteStepAttrib(file,"Step",H5T_NATIVE_INT64,&step,1);
    }
    H5PartCloseFile(file);
    printf("OK, close file and reopen for reading\n");
    file= H5PartOpenFile("parttest.h5",H5PART_READ);
    nt=H5PartGetNumSteps(file); /* get number of steps in file */
    H5PartSetStep(file,0);
    nds=H5PartGetNumDatasets(file); /* get number of datasets in timestep 0 */
    
    puts("\n\n===============================");
    for(i=0;i<nds;i++){ /* and print out those names */
      H5PartGetDatasetName(file,i,name,64);
      printf("\tDataset[%u] name=[%s]\n",
             i,name);
    }
    puts("===============================\n\n");
    
    nfattribs=H5PartGetNumFileAttribs(file);
    printf("Number of datasteps in file is %u num file attribs=%d\n",
           nt,nfattribs);
    for(t=0;t<nt;t++){
      int nparticles;
      H5PartSetStep(file,t); /* select a timestep */
      nparticles=(int)H5PartGetNumParticles(file);
      nsattribs=H5PartGetNumStepAttribs(file);
      printf("Step[%u] nparticles this step=%u stepattribs=%u\n",
             t,nparticles,nsattribs); /* get num particles this step */
      if(nsattribs>0){
        char attrname[32];
        H5PartGetStepAttribInfo(file,0,attrname,32,0,0);
        printf("First Attrib name is [%s]\n",attrname);
      }
      H5PartReadParticleStep(file,t,/* do a mongo read of all data this step */
                             x,y,z,id);
      printf("\tid\t\tx\t\ty\t\tz\n");
      puts("\t----------------------------------------------------");
      for(i=0;i<sz;i++) {
        printf("\t%llu\t%f\t%f\t%f\n\n",id[i],x[i],y[i],z[i]);
      }
    }
    H5PartCloseFile(file);
  }
  }
#endif

#endif

---