NetCDF  4.4.0
nc4file.c
Go to the documentation of this file.
1 
12 #include "config.h"
13 #include <errno.h> /* netcdf functions sometimes return system errors */
14 
15 
16 #include "nc.h"
17 #include "nc4internal.h"
18 #include "nc4dispatch.h"
19 
20 /* must be after nc4internal.h */
21 #include <H5DSpublic.h>
22 
23 #ifdef USE_HDF4
24 #include <mfhdf.h>
25 #endif
26 
27 #ifdef USE_DISKLESS
28 #include <hdf5_hl.h>
29 #endif
30 
31 /* This is to track opened HDF5 objects to make sure they are
32  * closed. */
33 #ifdef EXTRA_TESTS
34 extern int num_plists;
35 extern int num_spaces;
36 #endif /* EXTRA_TESTS */
37 
38 #define MIN_DEFLATE_LEVEL 0
39 #define MAX_DEFLATE_LEVEL 9
40 
41 /* These are the special attributes added by the HDF5 dimension scale
42  * API. They will be ignored by netCDF-4. */
43 #define REFERENCE_LIST "REFERENCE_LIST"
44 #define CLASS "CLASS"
45 #define DIMENSION_LIST "DIMENSION_LIST"
46 #define NAME "NAME"
47 
48 /* Define the illegal mode flags */
49 static const int ILLEGAL_OPEN_FLAGS = (NC_MMAP|NC_64BIT_OFFSET);
50 
51 static const int ILLEGAL_CREATE_FLAGS = (NC_NOWRITE|NC_MMAP|NC_INMEMORY|NC_64BIT_OFFSET|NC_CDF5);
52 
59 {
60  hid_t oid; /* HDF5 object ID */
61  char oname[NC_MAX_NAME + 1]; /* Name of object */
62  H5G_stat_t statbuf; /* Information about the object */
63  struct NC4_rec_read_metadata_obj_info *next; /* Pointer to next node in list */
65 
73 {
74  NC4_rec_read_metadata_obj_info_t *grps_head, *grps_tail; /* Pointers to head & tail of list of groups */
75  NC_GRP_INFO_T *grp; /* Pointer to parent group */
77 
78 /* Forward */
79 static int NC4_enddef(int ncid);
80 static int nc4_rec_read_metadata(NC_GRP_INFO_T *grp);
81 static int close_netcdf4_file(NC_HDF5_FILE_INFO_T *h5, int abort);
82 
83 /* These are the default chunk cache sizes for HDF5 files created or
84  * opened with netCDF-4. */
85 size_t nc4_chunk_cache_size = CHUNK_CACHE_SIZE;
86 size_t nc4_chunk_cache_nelems = CHUNK_CACHE_NELEMS;
87 float nc4_chunk_cache_preemption = CHUNK_CACHE_PREEMPTION;
88 
89 /* To turn off HDF5 error messages, I have to catch an early
90  invocation of a netcdf function. */
91 static int virgin = 1;
92 
93 /* For performance, fill this array only the first time, and keep it
94  * in global memory for each further use. */
95 #define NUM_TYPES 12
96 static hid_t h5_native_type_constant_g[NUM_TYPES];
97 static const char nc_type_name_g[NUM_TYPES][NC_MAX_NAME + 1] = {"char", "byte", "short",
98  "int", "float", "double", "ubyte",
99  "ushort", "uint", "int64",
100  "uint64", "string"};
101 static const nc_type nc_type_constant_g[NUM_TYPES] = {NC_CHAR, NC_BYTE, NC_SHORT,
105 static const int nc_type_size_g[NUM_TYPES] = {sizeof(char), sizeof(char), sizeof(short),
106  sizeof(int), sizeof(float), sizeof(double), sizeof(unsigned char),
107  sizeof(unsigned short), sizeof(unsigned int), sizeof(long long),
108  sizeof(unsigned long long), sizeof(char *)};
109 
110 /* Set chunk cache size. Only affects files opened/created *after* it
111  * is called. */
112 int
113 nc_set_chunk_cache(size_t size, size_t nelems, float preemption)
114 {
115  if (preemption < 0 || preemption > 1)
116  return NC_EINVAL;
117  nc4_chunk_cache_size = size;
118  nc4_chunk_cache_nelems = nelems;
119  nc4_chunk_cache_preemption = preemption;
120  return NC_NOERR;
121 }
122 
123 /* Get chunk cache size. Only affects files opened/created *after* it
124  * is called. */
125 int
126 nc_get_chunk_cache(size_t *sizep, size_t *nelemsp, float *preemptionp)
127 {
128  if (sizep)
129  *sizep = nc4_chunk_cache_size;
130 
131  if (nelemsp)
132  *nelemsp = nc4_chunk_cache_nelems;
133 
134  if (preemptionp)
135  *preemptionp = nc4_chunk_cache_preemption;
136  return NC_NOERR;
137 }
138 
139 /* Required for fortran to avoid size_t issues. */
140 int
141 nc_set_chunk_cache_ints(int size, int nelems, int preemption)
142 {
143  if (size <= 0 || nelems <= 0 || preemption < 0 || preemption > 100)
144  return NC_EINVAL;
145  nc4_chunk_cache_size = size;
146  nc4_chunk_cache_nelems = nelems;
147  nc4_chunk_cache_preemption = (float)preemption / 100;
148  return NC_NOERR;
149 }
150 
151 int
152 nc_get_chunk_cache_ints(int *sizep, int *nelemsp, int *preemptionp)
153 {
154  if (sizep)
155  *sizep = (int)nc4_chunk_cache_size;
156  if (nelemsp)
157  *nelemsp = (int)nc4_chunk_cache_nelems;
158  if (preemptionp)
159  *preemptionp = (int)(nc4_chunk_cache_preemption * 100);
160 
161  return NC_NOERR;
162 }
163 
164 /* This will return the length of a netcdf data type in bytes. */
165 int
166 nc4typelen(nc_type type)
167 {
168  switch(type){
169  case NC_BYTE:
170  case NC_CHAR:
171  case NC_UBYTE:
172  return 1;
173  case NC_USHORT:
174  case NC_SHORT:
175  return 2;
176  case NC_FLOAT:
177  case NC_INT:
178  case NC_UINT:
179  return 4;
180  case NC_DOUBLE:
181  case NC_INT64:
182  case NC_UINT64:
183  return 8;
184  }
185  return -1;
186 }
187 
188 /* Given a filename, check to see if it is a HDF5 file. */
189 #define MAGIC_NUMBER_LEN 4
190 #define NC_HDF5_FILE 1
191 #define NC_HDF4_FILE 2
192 static int
193 nc_check_for_hdf(const char *path, int flags, void* parameters, int *hdf_file)
194 {
195  char blob[MAGIC_NUMBER_LEN];
196 #ifdef USE_PARALLEL4
197  int use_parallel = ((flags & NC_MPIIO) == NC_MPIIO);
198  NC_MPI_INFO* mpiinfo = (NC_MPI_INFO*)parameters;
199  MPI_Comm comm = MPI_COMM_WORLD;
200  MPI_Info info = MPI_INFO_NULL;
201 #endif
202  int inmemory = ((flags & NC_INMEMORY) == NC_INMEMORY);
203 #ifdef USE_DISKLESS
204  NC_MEM_INFO* meminfo = (NC_MEM_INFO*)parameters;
205 #endif
206 
207 #ifdef USE_PARALLEL4
208  if(use_parallel) {
209  comm = mpiinfo->comm;
210  info = mpiinfo->info;
211  }
212 #endif
213 
214  assert(hdf_file);
215  LOG((3, "%s: path %s", __func__, path));
216 
217  /* HDF5 function handles possible user block at beginning of file */
218  if(!inmemory && H5Fis_hdf5(path))
219  {
220  *hdf_file = NC_HDF5_FILE;
221  } else {
222 
223 /* Get the 4-byte blob from the beginning of the file. Don't use posix
224  * for parallel, use the MPI functions instead. */
225 #ifdef USE_PARALLEL4
226  if (!inmemory && use_parallel)
227  {
228  MPI_File fh;
229  MPI_Status status;
230  int retval;
231  if ((retval = MPI_File_open(comm, (char *)path, MPI_MODE_RDONLY,
232  info, &fh)) != MPI_SUCCESS)
233  return NC_EPARINIT;
234  if ((retval = MPI_File_read(fh, blob, MAGIC_NUMBER_LEN, MPI_CHAR,
235  &status)) != MPI_SUCCESS)
236  return NC_EPARINIT;
237  if ((retval = MPI_File_close(&fh)) != MPI_SUCCESS)
238  return NC_EPARINIT;
239  }
240  else
241 #endif /* USE_PARALLEL4 */
242  if(!inmemory) {
243  FILE *fp;
244  if (!(fp = fopen(path, "r")) ||
245  fread(blob, MAGIC_NUMBER_LEN, 1, fp) != 1) {
246 
247  if(fp) fclose(fp);
248  return errno;
249  }
250  fclose(fp);
251  } else { /*inmemory*/
252  if(meminfo->size < MAGIC_NUMBER_LEN)
253  return NC_ENOTNC;
254  memcpy(blob,meminfo->memory,MAGIC_NUMBER_LEN);
255  }
256 
257  /* Check for HDF4. */
258  if (memcmp(blob, "\016\003\023\001", 4)==0)
259  *hdf_file = NC_HDF4_FILE;
260  else if (memcmp(blob, "HDF", 3)==0)
261  *hdf_file = NC_HDF5_FILE;
262  else
263  *hdf_file = 0;
264  }
265  return NC_NOERR;
266 }
267 
268 /* Create a HDF5/netcdf-4 file. */
269 
270 static int
271 nc4_create_file(const char *path, int cmode, MPI_Comm comm, MPI_Info info,
272  NC *nc)
273 {
274  hid_t fcpl_id, fapl_id = -1;
275  unsigned flags;
276  FILE *fp;
277  int retval = NC_NOERR;
278  NC_HDF5_FILE_INFO_T* nc4_info = NULL;
279 #ifdef USE_PARALLEL4
280  int comm_duped = 0; /* Whether the MPI Communicator was duplicated */
281  int info_duped = 0; /* Whether the MPI Info object was duplicated */
282 #else /* !USE_PARALLEL4 */
283  int persist = 0; /* Should diskless try to persist its data into file?*/
284 #endif
285 
286  assert(nc);
287 
288  if(cmode & NC_DISKLESS)
289  flags = H5F_ACC_TRUNC;
290  else if(cmode & NC_NOCLOBBER)
291  flags = H5F_ACC_EXCL;
292  else
293  flags = H5F_ACC_TRUNC;
294 
295  LOG((3, "%s: path %s mode 0x%x", __func__, path, cmode));
296  assert(nc && path);
297 
298  /* If this file already exists, and NC_NOCLOBBER is specified,
299  return an error. */
300  if (cmode & NC_DISKLESS) {
301 #ifndef USE_PARALLEL4
302  if(cmode & NC_WRITE)
303  persist = 1;
304 #endif
305  } else if ((cmode & NC_NOCLOBBER) && (fp = fopen(path, "r"))) {
306  fclose(fp);
307  return NC_EEXIST;
308  }
309 
310  /* Add necessary structs to hold netcdf-4 file data. */
311  if ((retval = nc4_nc4f_list_add(nc, path, (NC_WRITE | cmode))))
312  BAIL(retval);
313  nc4_info = NC4_DATA(nc);
314  assert(nc4_info && nc4_info->root_grp);
315 
316 #if 0 /*def USE_PNETCDF*/
317  if (cmode & NC_PNETCDF)
318  return NC_NOERR;
319 #endif
320 
321  /* Need this access plist to control how HDF5 handles open onjects
322  * on file close. (Setting H5F_CLOSE_SEMI will cause H5Fclose to
323  * fail if there are any open objects in the file. */
324  if ((fapl_id = H5Pcreate(H5P_FILE_ACCESS)) < 0)
325  BAIL(NC_EHDFERR);
326 #ifdef EXTRA_TESTS
327  num_plists++;
328 #endif
329 #ifdef EXTRA_TESTS
330  if (H5Pset_fclose_degree(fapl_id, H5F_CLOSE_SEMI))
331  BAIL(NC_EHDFERR);
332 #else
333  if (H5Pset_fclose_degree(fapl_id, H5F_CLOSE_STRONG))
334  BAIL(NC_EHDFERR);
335 #endif /* EXTRA_TESTS */
336 
337 #ifdef USE_PARALLEL4
338  /* If this is a parallel file create, set up the file creation
339  property list. */
340  if ((cmode & NC_MPIIO) || (cmode & NC_MPIPOSIX))
341  {
342  nc4_info->parallel = NC_TRUE;
343  if (cmode & NC_MPIIO) /* MPI/IO */
344  {
345  LOG((4, "creating parallel file with MPI/IO"));
346  if (H5Pset_fapl_mpio(fapl_id, comm, info) < 0)
347  BAIL(NC_EPARINIT);
348  }
349 #ifdef USE_PARALLEL_POSIX
350  else /* MPI/POSIX */
351  {
352  LOG((4, "creating parallel file with MPI/posix"));
353  if (H5Pset_fapl_mpiposix(fapl_id, comm, 0) < 0)
354  BAIL(NC_EPARINIT);
355  }
356 #else /* USE_PARALLEL_POSIX */
357  /* Should not happen! Code in NC4_create/NC4_open should alias the
358  * NC_MPIPOSIX flag to NC_MPIIO, if the MPI-POSIX VFD is not
359  * available in HDF5. -QAK
360  */
361  else /* MPI/POSIX */
362  BAIL(NC_EPARINIT);
363 #endif /* USE_PARALLEL_POSIX */
364 
365  /* Keep copies of the MPI Comm & Info objects */
366  if (MPI_SUCCESS != MPI_Comm_dup(comm, &nc4_info->comm))
367  BAIL(NC_EMPI);
368  comm_duped++;
369  if (MPI_INFO_NULL != info)
370  {
371  if (MPI_SUCCESS != MPI_Info_dup(info, &nc4_info->info))
372  BAIL(NC_EMPI);
373  info_duped++;
374  }
375  else
376  {
377  /* No dup, just copy it. */
378  nc4_info->info = info;
379  }
380  }
381 #else /* only set cache for non-parallel... */
382  if(cmode & NC_DISKLESS) {
383  if (H5Pset_fapl_core(fapl_id, 4096, persist))
384  BAIL(NC_EDISKLESS);
385  }
386  if (H5Pset_cache(fapl_id, 0, nc4_chunk_cache_nelems, nc4_chunk_cache_size,
387  nc4_chunk_cache_preemption) < 0)
388  BAIL(NC_EHDFERR);
389  LOG((4, "%s: set HDF raw chunk cache to size %d nelems %d preemption %f",
390  __func__, nc4_chunk_cache_size, nc4_chunk_cache_nelems, nc4_chunk_cache_preemption));
391 #endif /* USE_PARALLEL4 */
392 
393  if (H5Pset_libver_bounds(fapl_id, H5F_LIBVER_LATEST, H5F_LIBVER_LATEST) < 0)
394  BAIL(NC_EHDFERR);
395 
396  /* Create the property list. */
397  if ((fcpl_id = H5Pcreate(H5P_FILE_CREATE)) < 0)
398  BAIL(NC_EHDFERR);
399 #ifdef EXTRA_TESTS
400  num_plists++;
401 #endif
402 
403  /* RJ: this suppose to be FALSE that is defined in H5 private.h as 0 */
404  if (H5Pset_obj_track_times(fcpl_id,0)<0)
405  BAIL(NC_EHDFERR);
406 
407  /* Set latest_format in access propertly list and
408  * H5P_CRT_ORDER_TRACKED in the creation property list. This turns
409  * on HDF5 creation ordering. */
410  if (H5Pset_link_creation_order(fcpl_id, (H5P_CRT_ORDER_TRACKED |
411  H5P_CRT_ORDER_INDEXED)) < 0)
412  BAIL(NC_EHDFERR);
413  if (H5Pset_attr_creation_order(fcpl_id, (H5P_CRT_ORDER_TRACKED |
414  H5P_CRT_ORDER_INDEXED)) < 0)
415  BAIL(NC_EHDFERR);
416 
417  /* Create the file. */
418  if ((nc4_info->hdfid = H5Fcreate(path, flags, fcpl_id, fapl_id)) < 0)
419  /*Change the return error from NC_EFILEMETADATA to
420  System error EACCES because that is the more likely problem */
421  BAIL(EACCES);
422 
423  /* Open the root group. */
424  if ((nc4_info->root_grp->hdf_grpid = H5Gopen2(nc4_info->hdfid, "/",
425  H5P_DEFAULT)) < 0)
426  BAIL(NC_EFILEMETA);
427 
428  /* Release the property lists. */
429  if (H5Pclose(fapl_id) < 0 ||
430  H5Pclose(fcpl_id) < 0)
431  BAIL(NC_EHDFERR);
432 #ifdef EXTRA_TESTS
433  num_plists--;
434  num_plists--;
435 #endif
436 
437  /* Define mode gets turned on automatically on create. */
438  nc4_info->flags |= NC_INDEF;
439 
440  return NC_NOERR;
441 
442 exit: /*failure exit*/
443 #ifdef USE_PARALLEL4
444  if (comm_duped) MPI_Comm_free(&nc4_info->comm);
445  if (info_duped) MPI_Info_free(&nc4_info->info);
446 #endif
447 #ifdef EXTRA_TESTS
448  num_plists--;
449 #endif
450  if (fapl_id != H5P_DEFAULT) H5Pclose(fapl_id);
451  if(!nc4_info) return retval;
452  close_netcdf4_file(nc4_info,1); /* treat like abort */
453  return retval;
454 }
455 
471 int
472 NC4_create(const char* path, int cmode, size_t initialsz, int basepe,
473  size_t *chunksizehintp, int use_parallel, void *parameters,
474  NC_Dispatch *dispatch, NC* nc_file)
475 {
476  MPI_Comm comm = MPI_COMM_WORLD;
477  MPI_Info info = MPI_INFO_NULL;
478  int res;
479 
480  assert(nc_file && path);
481 
482  LOG((1, "%s: path %s cmode 0x%x comm %d info %d",
483  __func__, path, cmode, comm, info));
484 
485 #ifdef USE_PARALLEL4
486  if (parameters)
487  {
488  comm = ((NC_MPI_INFO *)parameters)->comm;
489  info = ((NC_MPI_INFO *)parameters)->info;
490  }
491 #endif /* USE_PARALLEL4 */
492 
493  /* If this is our first file, turn off HDF5 error messages. */
494  if (virgin)
495  {
496  if (H5Eset_auto(NULL, NULL) < 0)
497  LOG((0, "Couldn't turn off HDF5 error messages!"));
498  LOG((1, "HDF5 error messages have been turned off."));
499  virgin = 0;
500  }
501 
502  /* Check the cmode for validity. */
503  if((cmode & ILLEGAL_CREATE_FLAGS) != 0)
504  return NC_EINVAL;
505 
506  /* Cannot have both */
507  if((cmode & (NC_MPIIO|NC_MPIPOSIX)) == (NC_MPIIO|NC_MPIPOSIX))
508  return NC_EINVAL;
509 
510  /* Currently no parallel diskless io */
511  if((cmode & (NC_MPIIO | NC_MPIPOSIX)) && (cmode & NC_DISKLESS))
512  return NC_EINVAL;
513 
514 #ifndef USE_PARALLEL_POSIX
515 /* If the HDF5 library has been compiled without the MPI-POSIX VFD, alias
516  * the NC_MPIPOSIX flag to NC_MPIIO. -QAK
517  */
518  if(cmode & NC_MPIPOSIX)
519  {
520  cmode &= ~NC_MPIPOSIX;
521  cmode |= NC_MPIIO;
522  }
523 #endif /* USE_PARALLEL_POSIX */
524 
525  cmode |= NC_NETCDF4;
526 
527  /* Apply default create format. */
528  if (nc_get_default_format() == NC_FORMAT_CDF5)
529  cmode |= NC_CDF5;
530  else if (nc_get_default_format() == NC_FORMAT_64BIT_OFFSET)
531  cmode |= NC_64BIT_OFFSET;
532  else if (nc_get_default_format() == NC_FORMAT_NETCDF4_CLASSIC)
533  cmode |= NC_CLASSIC_MODEL;
534 
535  LOG((2, "cmode after applying default format: 0x%x", cmode));
536 
537  nc_file->int_ncid = nc_file->ext_ncid;
538  res = nc4_create_file(path, cmode, comm, info, nc_file);
539 
540 #if 0 /*def USE_PNETCDF*/
541  if (cmode & NC_PNETCDF)
542  {
543  NC_HDF5_FILE_INFO_T* nc4_info;
544  nc4_info = NC4_DATA(nc_file);
545  assert(nc4_info);
546 
547  nc4_info->pnetcdf_file++;
548  res = ncmpi_create(comm, path, cmode, info, &(nc_file->int_ncid));
549  }
550 #endif /* USE_PNETCDF */
551 
552  return res;
553 }
554 
555 /* This function is called by read_dataset when a dimension scale
556  * dataset is encountered. It reads in the dimension data (creating a
557  * new NC_DIM_INFO_T object), and also checks to see if this is a
558  * dimension without a variable - that is, a coordinate dimension
559  * which does not have any coordinate data. */
560 static int
561 read_scale(NC_GRP_INFO_T *grp, hid_t datasetid, const char *obj_name,
562  const H5G_stat_t *statbuf, hsize_t scale_size, hsize_t max_scale_size,
563  NC_DIM_INFO_T **dim)
564 {
565  NC_DIM_INFO_T *new_dim; /* Dimension added to group */
566  char dimscale_name_att[NC_MAX_NAME + 1] = ""; /* Dimscale name, for checking if dim without var */
567  htri_t attr_exists = -1; /* Flag indicating hidden attribute exists */
568  hid_t attid = -1; /* ID of hidden attribute (to store dim ID) */
569  int dimscale_created = 0; /* Remember if a dimension was created (for error recovery) */
570  int initial_grp_ndims = grp->ndims; /* Retain for error recovery */
571  short initial_next_dimid = grp->nc4_info->next_dimid;/* Retain for error recovery */
572  int retval;
573 
574  /* Add a dimension for this scale. */
575  if ((retval = nc4_dim_list_add(&grp->dim, &new_dim)))
576  BAIL(retval);
577  dimscale_created++;
578 
579  /* Does this dataset have a hidden attribute that tells us its
580  * dimid? If so, read it. */
581  if ((attr_exists = H5Aexists(datasetid, NC_DIMID_ATT_NAME)) < 0)
582  BAIL(NC_EHDFERR);
583  if (attr_exists)
584  {
585  if ((attid = H5Aopen_by_name(datasetid, ".", NC_DIMID_ATT_NAME,
586  H5P_DEFAULT, H5P_DEFAULT)) < 0)
587  BAIL(NC_EHDFERR);
588 
589  if (H5Aread(attid, H5T_NATIVE_INT, &new_dim->dimid) < 0)
590  BAIL(NC_EHDFERR);
591 
592  /* Check if scale's dimid should impact the group's next dimid */
593  if (new_dim->dimid >= grp->nc4_info->next_dimid)
594  grp->nc4_info->next_dimid = new_dim->dimid + 1;
595  }
596  else
597  {
598  /* Assign dimid */
599  new_dim->dimid = grp->nc4_info->next_dimid++;
600  }
601 
602  /* Increment number of dimensions. */
603  grp->ndims++;
604 
605  if (!(new_dim->name = strdup(obj_name)))
606  BAIL(NC_ENOMEM);
607  if (SIZEOF_SIZE_T < 8 && scale_size > NC_MAX_UINT)
608  {
609  new_dim->len = NC_MAX_UINT;
610  new_dim->too_long = NC_TRUE;
611  }
612  else
613  new_dim->len = scale_size;
614  new_dim->hdf5_objid.fileno[0] = statbuf->fileno[0];
615  new_dim->hdf5_objid.fileno[1] = statbuf->fileno[1];
616  new_dim->hdf5_objid.objno[0] = statbuf->objno[0];
617  new_dim->hdf5_objid.objno[1] = statbuf->objno[1];
618 
619  /* If the dimscale has an unlimited dimension, then this dimension
620  * is unlimited. */
621  if (max_scale_size == H5S_UNLIMITED)
622  new_dim->unlimited = NC_TRUE;
623 
624  /* If the scale name is set to DIM_WITHOUT_VARIABLE, then this is a
625  * dimension, but not a variable. (If get_scale_name returns an
626  * error, just move on, there's no NAME.) */
627  if (H5DSget_scale_name(datasetid, dimscale_name_att, NC_MAX_NAME) >= 0)
628  {
629  if (!strncmp(dimscale_name_att, DIM_WITHOUT_VARIABLE,
630  strlen(DIM_WITHOUT_VARIABLE)))
631  {
632  if (new_dim->unlimited)
633  {
634  size_t len = 0, *lenp = &len;
635 
636  if ((retval = nc4_find_dim_len(grp, new_dim->dimid, &lenp)))
637  BAIL(retval);
638  new_dim->len = *lenp;
639  }
640 
641  /* Hold open the dataset, since the dimension doesn't have a coordinate variable */
642  new_dim->hdf_dimscaleid = datasetid;
643  H5Iinc_ref(new_dim->hdf_dimscaleid); /* Increment number of objects using ID */
644  }
645  }
646 
647  /* Set the dimension created */
648  *dim = new_dim;
649 
650 exit:
651  /* Close the hidden attribute, if it was opened (error, or no error) */
652  if (attid > 0 && H5Aclose(attid) < 0)
653  BAIL2(NC_EHDFERR);
654 
655  /* On error, undo any dimscale creation */
656  if (retval < 0 && dimscale_created)
657  {
658  /* Delete the dimension */
659  if ((retval = nc4_dim_list_del(&grp->dim, new_dim)))
660  BAIL2(retval);
661 
662  /* Reset the group's information */
663  grp->ndims = initial_grp_ndims;
664  grp->nc4_info->next_dimid = initial_next_dimid;
665  }
666 
667  return retval;
668 }
669 
670 /* This function reads the hacked in coordinates attribute I use for
671  * multi-dimensional coordinates. */
672 static int
673 read_coord_dimids(NC_VAR_INFO_T *var)
674 {
675  hid_t coord_att_typeid = -1, coord_attid = -1, spaceid = -1;
676  hssize_t npoints;
677  int ret = 0;
678 
679  /* There is a hidden attribute telling us the ids of the
680  * dimensions that apply to this multi-dimensional coordinate
681  * variable. Read it. */
682  if ((coord_attid = H5Aopen_name(var->hdf_datasetid, COORDINATES)) < 0) ret++;
683  if (!ret && (coord_att_typeid = H5Aget_type(coord_attid)) < 0) ret++;
684 
685  /* How many dimensions are there? */
686  if (!ret && (spaceid = H5Aget_space(coord_attid)) < 0) ret++;
687 #ifdef EXTRA_TESTS
688  num_spaces++;
689 #endif
690  if (!ret && (npoints = H5Sget_simple_extent_npoints(spaceid)) < 0) ret++;
691 
692  /* Check that the number of points is the same as the number of dimensions
693  * for the variable */
694  if (!ret && npoints != var->ndims) ret++;
695 
696  if (!ret && H5Aread(coord_attid, coord_att_typeid, var->dimids) < 0) ret++;
697  LOG((4, "dimscale %s is multidimensional and has coords", var->name));
698 
699  /* Set my HDF5 IDs free! */
700  if (spaceid >= 0 && H5Sclose(spaceid) < 0) ret++;
701 #ifdef EXTRA_TESTS
702  num_spaces--;
703 #endif
704  if (coord_att_typeid >= 0 && H5Tclose(coord_att_typeid) < 0) ret++;
705  if (coord_attid >= 0 && H5Aclose(coord_attid) < 0) ret++;
706  return ret ? NC_EATTMETA : NC_NOERR;
707 }
708 
709 /* This function is called when reading a file's metadata for each
710  * dimension scale attached to a variable.*/
711 static herr_t
712 dimscale_visitor(hid_t did, unsigned dim, hid_t dsid,
713  void *dimscale_hdf5_objids)
714 {
715  H5G_stat_t statbuf;
716 
717  /* Get more info on the dimscale object.*/
718  if (H5Gget_objinfo(dsid, ".", 1, &statbuf) < 0)
719  return -1;
720 
721  /* Pass this information back to caller. */
722  (*(HDF5_OBJID_T *)dimscale_hdf5_objids).fileno[0] = statbuf.fileno[0];
723  (*(HDF5_OBJID_T *)dimscale_hdf5_objids).fileno[1] = statbuf.fileno[1];
724  (*(HDF5_OBJID_T *)dimscale_hdf5_objids).objno[0] = statbuf.objno[0];
725  (*(HDF5_OBJID_T *)dimscale_hdf5_objids).objno[1] = statbuf.objno[1];
726  return 0;
727 }
728 
729 /* Given an HDF5 type, set a pointer to netcdf type. */
730 static int
731 get_netcdf_type(NC_HDF5_FILE_INFO_T *h5, hid_t native_typeid,
732  nc_type *xtype)
733 {
734  NC_TYPE_INFO_T *type;
735  H5T_class_t class;
736  htri_t is_str, equal = 0;
737 
738  assert(h5 && xtype);
739 
740  if ((class = H5Tget_class(native_typeid)) < 0)
741  return NC_EHDFERR;
742 
743  /* H5Tequal doesn't work with H5T_C_S1 for some reason. But
744  * H5Tget_class will return H5T_STRING if this is a string. */
745  if (class == H5T_STRING)
746  {
747  if ((is_str = H5Tis_variable_str(native_typeid)) < 0)
748  return NC_EHDFERR;
749  if (is_str)
750  *xtype = NC_STRING;
751  else
752  *xtype = NC_CHAR;
753  return NC_NOERR;
754  }
755  else if (class == H5T_INTEGER || class == H5T_FLOAT)
756  {
757  /* For integers and floats, we don't have to worry about
758  * endianness if we compare native types. */
759  if ((equal = H5Tequal(native_typeid, H5T_NATIVE_SCHAR)) < 0)
760  return NC_EHDFERR;
761  if (equal)
762  {
763  *xtype = NC_BYTE;
764  return NC_NOERR;
765  }
766  if ((equal = H5Tequal(native_typeid, H5T_NATIVE_SHORT)) < 0)
767  return NC_EHDFERR;
768  if (equal)
769  {
770  *xtype = NC_SHORT;
771  return NC_NOERR;
772  }
773  if ((equal = H5Tequal(native_typeid, H5T_NATIVE_INT)) < 0)
774  return NC_EHDFERR;
775  if (equal)
776  {
777  *xtype = NC_INT;
778  return NC_NOERR;
779  }
780  if ((equal = H5Tequal(native_typeid, H5T_NATIVE_FLOAT)) < 0)
781  return NC_EHDFERR;
782  if (equal)
783  {
784  *xtype = NC_FLOAT;
785  return NC_NOERR;
786  }
787  if ((equal = H5Tequal(native_typeid, H5T_NATIVE_DOUBLE)) < 0)
788  return NC_EHDFERR;
789  if (equal)
790  {
791  *xtype = NC_DOUBLE;
792  return NC_NOERR;
793  }
794  if ((equal = H5Tequal(native_typeid, H5T_NATIVE_UCHAR)) < 0)
795  return NC_EHDFERR;
796  if (equal)
797  {
798  *xtype = NC_UBYTE;
799  return NC_NOERR;
800  }
801  if ((equal = H5Tequal(native_typeid, H5T_NATIVE_USHORT)) < 0)
802  return NC_EHDFERR;
803  if (equal)
804  {
805  *xtype = NC_USHORT;
806  return NC_NOERR;
807  }
808  if ((equal = H5Tequal(native_typeid, H5T_NATIVE_UINT)) < 0)
809  return NC_EHDFERR;
810  if (equal)
811  {
812  *xtype = NC_UINT;
813  return NC_NOERR;
814  }
815  if ((equal = H5Tequal(native_typeid, H5T_NATIVE_LLONG)) < 0)
816  return NC_EHDFERR;
817  if (equal)
818  {
819  *xtype = NC_INT64;
820  return NC_NOERR;
821  }
822  if ((equal = H5Tequal(native_typeid, H5T_NATIVE_ULLONG)) < 0)
823  return NC_EHDFERR;
824  if (equal)
825  {
826  *xtype = NC_UINT64;
827  return NC_NOERR;
828  }
829  }
830 
831  /* Maybe we already know about this type. */
832  if (!equal)
833  if((type = nc4_rec_find_hdf_type(h5->root_grp, native_typeid)))
834  {
835  *xtype = type->nc_typeid;
836  return NC_NOERR;
837  }
838 
839  *xtype = NC_NAT;
840  return NC_EBADTYPID;
841 }
842 
843 /* Given an HDF5 type, set a pointer to netcdf type_info struct,
844  * either an existing one (for user-defined types) or a newly created
845  * one. */
846 static int
847 get_type_info2(NC_HDF5_FILE_INFO_T *h5, hid_t datasetid,
848  NC_TYPE_INFO_T **type_info)
849 {
850  htri_t is_str, equal = 0;
851  H5T_class_t class;
852  hid_t native_typeid, hdf_typeid;
853  H5T_order_t order;
854  int t;
855 
856  assert(h5 && type_info);
857 
858  /* Because these N5T_NATIVE_* constants are actually function calls
859  * (!) in H5Tpublic.h, I can't initialize this array in the usual
860  * way, because at least some C compilers (like Irix) complain
861  * about calling functions when defining constants. So I have to do
862  * it like this. Note that there's no native types for char or
863  * string. Those are handled later. */
864  if (!h5_native_type_constant_g[1])
865  {
866  h5_native_type_constant_g[1] = H5T_NATIVE_SCHAR;
867  h5_native_type_constant_g[2] = H5T_NATIVE_SHORT;
868  h5_native_type_constant_g[3] = H5T_NATIVE_INT;
869  h5_native_type_constant_g[4] = H5T_NATIVE_FLOAT;
870  h5_native_type_constant_g[5] = H5T_NATIVE_DOUBLE;
871  h5_native_type_constant_g[6] = H5T_NATIVE_UCHAR;
872  h5_native_type_constant_g[7] = H5T_NATIVE_USHORT;
873  h5_native_type_constant_g[8] = H5T_NATIVE_UINT;
874  h5_native_type_constant_g[9] = H5T_NATIVE_LLONG;
875  h5_native_type_constant_g[10] = H5T_NATIVE_ULLONG;
876  }
877 
878  /* Get the HDF5 typeid - we'll need it later. */
879  if ((hdf_typeid = H5Dget_type(datasetid)) < 0)
880  return NC_EHDFERR;
881 
882  /* Get the native typeid. Will be equivalent to hdf_typeid when
883  * creating but not necessarily when reading, a variable. */
884  if ((native_typeid = H5Tget_native_type(hdf_typeid, H5T_DIR_DEFAULT)) < 0)
885  return NC_EHDFERR;
886 
887  /* Is this type an integer, string, compound, or what? */
888  if ((class = H5Tget_class(native_typeid)) < 0)
889  return NC_EHDFERR;
890 
891  /* Is this an atomic type? */
892  if (class == H5T_STRING || class == H5T_INTEGER || class == H5T_FLOAT)
893  {
894  /* Allocate a phony NC_TYPE_INFO_T struct to hold type info. */
895  if (!(*type_info = calloc(1, sizeof(NC_TYPE_INFO_T))))
896  return NC_ENOMEM;
897 
898  /* H5Tequal doesn't work with H5T_C_S1 for some reason. But
899  * H5Tget_class will return H5T_STRING if this is a string. */
900  if (class == H5T_STRING)
901  {
902  if ((is_str = H5Tis_variable_str(native_typeid)) < 0)
903  return NC_EHDFERR;
904  /* Make sure fixed-len strings will work like variable-len strings */
905  if (is_str || H5Tget_size(hdf_typeid) > 1)
906  {
907  /* Set a class for the type */
908  t = NUM_TYPES - 1;
909  (*type_info)->nc_type_class = NC_STRING;
910  }
911  else
912  {
913  /* Set a class for the type */
914  t = 0;
915  (*type_info)->nc_type_class = NC_CHAR;
916  }
917  }
918  else if (class == H5T_INTEGER || class == H5T_FLOAT)
919  {
920  for (t = 1; t < NUM_TYPES - 1; t++)
921  {
922  if ((equal = H5Tequal(native_typeid, h5_native_type_constant_g[t])) < 0)
923  return NC_EHDFERR;
924  if (equal)
925  break;
926  }
927 
928  /* Find out about endianness.
929  * As of HDF 1.8.6, this works with all data types
930  * Not just H5T_INTEGER.
931  *
932  * See https://www.hdfgroup.org/HDF5/doc/RM/RM_H5T.html#Datatype-GetOrder
933  */
934  if((order = H5Tget_order(hdf_typeid)) < 0)
935  return NC_EHDFERR;
936 
937  if(order == H5T_ORDER_LE)
938  (*type_info)->endianness = NC_ENDIAN_LITTLE;
939  else if(order == H5T_ORDER_BE)
940  (*type_info)->endianness = NC_ENDIAN_BIG;
941  else
942  return NC_EBADTYPE;
943 
944  if(class == H5T_INTEGER)
945  (*type_info)->nc_type_class = NC_INT;
946  else
947  (*type_info)->nc_type_class = NC_FLOAT;
948  }
949  (*type_info)->nc_typeid = nc_type_constant_g[t];
950  (*type_info)->size = nc_type_size_g[t];
951  if (!((*type_info)->name = strdup(nc_type_name_g[t])))
952  return NC_ENOMEM;
953  (*type_info)->hdf_typeid = hdf_typeid;
954  (*type_info)->native_hdf_typeid = native_typeid;
955  return NC_NOERR;
956  }
957  else
958  {
959  NC_TYPE_INFO_T *type;
960 
961  /* This is a user-defined type. */
962  if((type = nc4_rec_find_hdf_type(h5->root_grp, native_typeid)))
963  *type_info = type;
964 
965  /* The type entry in the array of user-defined types already has
966  * an open data typeid (and native typeid), so close the ones we
967  * opened above. */
968  if (H5Tclose(native_typeid) < 0)
969  return NC_EHDFERR;
970  if (H5Tclose(hdf_typeid) < 0)
971  return NC_EHDFERR;
972 
973  if (type)
974  return NC_NOERR;
975  }
976 
977  return NC_EBADTYPID;
978 }
979 
980 /* Read an attribute. */
981 static int
982 read_hdf5_att(NC_GRP_INFO_T *grp, hid_t attid, NC_ATT_INFO_T *att)
983 {
984  hid_t spaceid = 0, file_typeid = 0;
985  hsize_t dims[1] = {0}; /* netcdf attributes always 1-D. */
986  int retval = NC_NOERR;
987  size_t type_size;
988  int att_ndims;
989  hssize_t att_npoints;
990  H5T_class_t att_class;
991  int fixed_len_string = 0;
992  size_t fixed_size = 0;
993 
994  assert(att->name);
995  LOG((5, "%s: att->attnum %d att->name %s att->nc_typeid %d att->len %d",
996  __func__, att->attnum, att->name, (int)att->nc_typeid, att->len));
997 
998  /* Get type of attribute in file. */
999  if ((file_typeid = H5Aget_type(attid)) < 0)
1000  return NC_EATTMETA;
1001  if ((att->native_hdf_typeid = H5Tget_native_type(file_typeid, H5T_DIR_DEFAULT)) < 0)
1002  BAIL(NC_EHDFERR);
1003  if ((att_class = H5Tget_class(att->native_hdf_typeid)) < 0)
1004  BAIL(NC_EATTMETA);
1005  if (att_class == H5T_STRING && !H5Tis_variable_str(att->native_hdf_typeid))
1006  {
1007  fixed_len_string++;
1008  if (!(fixed_size = H5Tget_size(att->native_hdf_typeid)))
1009  BAIL(NC_EATTMETA);
1010  }
1011  if ((retval = get_netcdf_type(grp->nc4_info, att->native_hdf_typeid, &(att->nc_typeid))))
1012  BAIL(retval);
1013 
1014 
1015  /* Get len. */
1016  if ((spaceid = H5Aget_space(attid)) < 0)
1017  BAIL(NC_EATTMETA);
1018 #ifdef EXTRA_TESTS
1019  num_spaces++;
1020 #endif
1021  if ((att_ndims = H5Sget_simple_extent_ndims(spaceid)) < 0)
1022  BAIL(NC_EATTMETA);
1023  if ((att_npoints = H5Sget_simple_extent_npoints(spaceid)) < 0)
1024  BAIL(NC_EATTMETA);
1025 
1026  /* If both att_ndims and att_npoints are zero, then this is a
1027  * zero length att. */
1028  if (att_ndims == 0 && att_npoints == 0)
1029  dims[0] = 0;
1030  else if (att->nc_typeid == NC_STRING)
1031  dims[0] = att_npoints;
1032  else if (att->nc_typeid == NC_CHAR)
1033  {
1034  /* NC_CHAR attributes are written as a scalar in HDF5, of type
1035  * H5T_C_S1, of variable length. */
1036  if (att_ndims == 0)
1037  {
1038  if (!(dims[0] = H5Tget_size(file_typeid)))
1039  BAIL(NC_EATTMETA);
1040  }
1041  else
1042  {
1043  /* This is really a string type! */
1044  att->nc_typeid = NC_STRING;
1045  dims[0] = att_npoints;
1046  }
1047  }
1048  else
1049  {
1050  H5S_class_t space_class;
1051 
1052  /* All netcdf attributes are scalar or 1-D only. */
1053  if (att_ndims > 1)
1054  BAIL(NC_EATTMETA);
1055 
1056  /* Check class of HDF5 dataspace */
1057  if ((space_class = H5Sget_simple_extent_type(spaceid)) < 0)
1058  BAIL(NC_EATTMETA);
1059 
1060  /* Check for NULL HDF5 dataspace class (should be weeded out earlier) */
1061  if (H5S_NULL == space_class)
1062  BAIL(NC_EATTMETA);
1063 
1064  /* check for SCALAR HDF5 dataspace class */
1065  if (H5S_SCALAR == space_class)
1066  dims[0] = 1;
1067  else /* Must be "simple" dataspace */
1068  {
1069  /* Read the size of this attribute. */
1070  if (H5Sget_simple_extent_dims(spaceid, dims, NULL) < 0)
1071  BAIL(NC_EATTMETA);
1072  }
1073  }
1074 
1075  /* Tell the user what the length if this attribute is. */
1076  att->len = dims[0];
1077 
1078  /* Allocate some memory if the len is not zero, and read the
1079  attribute. */
1080  if (dims[0])
1081  {
1082  if ((retval = nc4_get_typelen_mem(grp->nc4_info, att->nc_typeid, 0,
1083  &type_size)))
1084  return retval;
1085  if (att_class == H5T_VLEN)
1086  {
1087  if (!(att->vldata = malloc((unsigned int)(att->len * sizeof(hvl_t)))))
1088  BAIL(NC_ENOMEM);
1089  if (H5Aread(attid, att->native_hdf_typeid, att->vldata) < 0)
1090  BAIL(NC_EATTMETA);
1091  }
1092  else if (att->nc_typeid == NC_STRING)
1093  {
1094  if (!(att->stdata = calloc(att->len, sizeof(char *))))
1095  BAIL(NC_ENOMEM);
1096  /* For a fixed length HDF5 string, the read requires
1097  * contiguous memory. Meanwhile, the netCDF API requires that
1098  * nc_free_string be called on string arrays, which would not
1099  * work if one contiguous memory block were used. So here I
1100  * convert the contiguous block of strings into an array of
1101  * malloced strings (each string with its own malloc). Then I
1102  * copy the data and free the contiguous memory. This
1103  * involves copying the data, which is bad, but this only
1104  * occurs for fixed length string attributes, and presumably
1105  * these are small. (And netCDF-4 does not create them - it
1106  * always uses variable length strings. */
1107  if (fixed_len_string)
1108  {
1109  int i;
1110  char *contig_buf, *cur;
1111 
1112  /* Alloc space for the contiguous memory read. */
1113  if (!(contig_buf = malloc(att->len * fixed_size * sizeof(char))))
1114  BAIL(NC_ENOMEM);
1115 
1116  /* Read the fixed-len strings as one big block. */
1117  if (H5Aread(attid, att->native_hdf_typeid, contig_buf) < 0) {
1118  free(contig_buf);
1119  BAIL(NC_EATTMETA);
1120  }
1121 
1122  /* Copy strings, one at a time, into their new home. Alloc
1123  space for each string. The user will later free this
1124  space with nc_free_string. */
1125  cur = contig_buf;
1126  for (i = 0; i < att->len; i++)
1127  {
1128  if (!(att->stdata[i] = malloc(fixed_size))) {
1129  free(contig_buf);
1130  BAIL(NC_ENOMEM);
1131  }
1132  strncpy(att->stdata[i], cur, fixed_size);
1133  cur += fixed_size;
1134  }
1135 
1136  /* Free contiguous memory buffer. */
1137  free(contig_buf);
1138  }
1139  else
1140  {
1141  /* Read variable-length string atts. */
1142  if (H5Aread(attid, att->native_hdf_typeid, att->stdata) < 0)
1143  BAIL(NC_EATTMETA);
1144  }
1145  }
1146  else
1147  {
1148  if (!(att->data = malloc((unsigned int)(att->len * type_size))))
1149  BAIL(NC_ENOMEM);
1150  if (H5Aread(attid, att->native_hdf_typeid, att->data) < 0)
1151  BAIL(NC_EATTMETA);
1152  }
1153  }
1154 
1155  if (H5Tclose(file_typeid) < 0)
1156  BAIL(NC_EHDFERR);
1157  if (H5Sclose(spaceid) < 0)
1158  return NC_EHDFERR;
1159 #ifdef EXTRA_TESTS
1160  num_spaces--;
1161 #endif
1162 
1163  return NC_NOERR;
1164 
1165  exit:
1166  if (H5Tclose(file_typeid) < 0)
1167  BAIL2(NC_EHDFERR);
1168  if (spaceid > 0 && H5Sclose(spaceid) < 0)
1169  BAIL2(NC_EHDFERR);
1170 #ifdef EXTRA_TESTS
1171  num_spaces--;
1172 #endif
1173  return retval;
1174 }
1175 
1176 /* Read information about a user defined type from the HDF5 file, and
1177  * stash it in the group's list of types. */
1178 static int
1179 read_type(NC_GRP_INFO_T *grp, hid_t hdf_typeid, char *type_name)
1180 {
1181  NC_TYPE_INFO_T *type;
1182  H5T_class_t class;
1183  hid_t native_typeid;
1184  size_t type_size;
1185  int retval = NC_NOERR;
1186 
1187  assert(grp && type_name);
1188 
1189  LOG((4, "%s: type_name %s grp->name %s", __func__, type_name, grp->name));
1190 
1191  /* What is the native type for this platform? */
1192  if ((native_typeid = H5Tget_native_type(hdf_typeid, H5T_DIR_DEFAULT)) < 0)
1193  return NC_EHDFERR;
1194 
1195  /* What is the size of this type on this platform. */
1196  if (!(type_size = H5Tget_size(native_typeid)))
1197  return NC_EHDFERR;
1198  LOG((5, "type_size %d", type_size));
1199 
1200  /* Add to the list for this new type, and get a local pointer to it. */
1201  if ((retval = nc4_type_list_add(grp, type_size, type_name, &type)))
1202  return retval;
1203 
1204  /* Remember common info about this type. */
1205  type->committed = NC_TRUE;
1206  type->hdf_typeid = hdf_typeid;
1207  H5Iinc_ref(type->hdf_typeid); /* Increment number of objects using ID */
1208  type->native_hdf_typeid = native_typeid;
1209 
1210  /* What is the class of this type, compound, vlen, etc. */
1211  if ((class = H5Tget_class(hdf_typeid)) < 0)
1212  return NC_EHDFERR;
1213  switch (class)
1214  {
1215  case H5T_STRING:
1216  type->nc_type_class = NC_STRING;
1217  break;
1218 
1219  case H5T_COMPOUND:
1220  {
1221  int nmembers;
1222  unsigned int m;
1223  char* member_name = NULL;
1224 #ifdef JNA
1225  char jna[1001];
1226 #endif
1227 
1228  type->nc_type_class = NC_COMPOUND;
1229 
1230  if ((nmembers = H5Tget_nmembers(hdf_typeid)) < 0)
1231  return NC_EHDFERR;
1232  LOG((5, "compound type has %d members", nmembers));
1233  for (m = 0; m < nmembers; m++)
1234  {
1235  hid_t member_hdf_typeid;
1236  hid_t member_native_typeid;
1237  size_t member_offset;
1238  H5T_class_t mem_class;
1239  nc_type member_xtype;
1240 
1241 
1242  /* Get the typeid and native typeid of this member of the
1243  * compound type. */
1244  if ((member_hdf_typeid = H5Tget_member_type(type->native_hdf_typeid, m)) < 0)
1245  return NC_EHDFERR;
1246 
1247  if ((member_native_typeid = H5Tget_native_type(member_hdf_typeid, H5T_DIR_DEFAULT)) < 0)
1248  return NC_EHDFERR;
1249 
1250  /* Get the name of the member.*/
1251  member_name = H5Tget_member_name(type->native_hdf_typeid, m);
1252  if (!member_name || strlen(member_name) > NC_MAX_NAME) {
1253  retval = NC_EBADNAME;
1254  break;
1255  }
1256 #ifdef JNA
1257  else {
1258  strncpy(jna,member_name,1000);
1259  member_name = jna;
1260  }
1261 #endif
1262 
1263  /* Offset in bytes on *this* platform. */
1264  member_offset = H5Tget_member_offset(type->native_hdf_typeid, m);
1265 
1266  /* Get dimensional data if this member is an array of something. */
1267  if ((mem_class = H5Tget_class(member_hdf_typeid)) < 0)
1268  return NC_EHDFERR;
1269  if (mem_class == H5T_ARRAY)
1270  {
1271  int ndims, dim_size[NC_MAX_VAR_DIMS];
1272  hsize_t dims[NC_MAX_VAR_DIMS];
1273  int d;
1274 
1275  if ((ndims = H5Tget_array_ndims(member_hdf_typeid)) < 0) {
1276  retval = NC_EHDFERR;
1277  break;
1278  }
1279  if (H5Tget_array_dims(member_hdf_typeid, dims, NULL) != ndims) {
1280  retval = NC_EHDFERR;
1281  break;
1282  }
1283  for (d = 0; d < ndims; d++)
1284  dim_size[d] = dims[d];
1285 
1286  /* What is the netCDF typeid of this member? */
1287  if ((retval = get_netcdf_type(grp->nc4_info, H5Tget_super(member_hdf_typeid),
1288  &member_xtype)))
1289  break;
1290 
1291  /* Add this member to our list of fields in this compound type. */
1292  if ((retval = nc4_field_list_add(&type->u.c.field, type->u.c.num_fields++, member_name,
1293  member_offset, H5Tget_super(member_hdf_typeid),
1294  H5Tget_super(member_native_typeid),
1295  member_xtype, ndims, dim_size)))
1296  break;
1297  }
1298  else
1299  {
1300  /* What is the netCDF typeid of this member? */
1301  if ((retval = get_netcdf_type(grp->nc4_info, member_native_typeid,
1302  &member_xtype)))
1303  break;
1304 
1305  /* Add this member to our list of fields in this compound type. */
1306  if ((retval = nc4_field_list_add(&type->u.c.field, type->u.c.num_fields++, member_name,
1307  member_offset, member_hdf_typeid, member_native_typeid,
1308  member_xtype, 0, NULL)))
1309  break;
1310  }
1311 
1312 #ifndef JNA
1313 
1314  /* Free the member name (which HDF5 allocated for us). */
1315  /* On Windows using the microsoft runtime, it is an error
1316  for one library to free memory allocated by a different library.
1317  IF it is available, we should use H5free_memory*/
1318 
1319 #ifdef HDF5_HAS_H5FREE
1320  if(member_name != NULL) H5free_memory(member_name);
1321 #else
1322 #ifndef _MSC_VER
1323  if(member_name != NULL) free(member_name);
1324 #endif
1325 #endif
1326 #endif
1327  member_name = NULL;
1328  }
1329 #ifndef JNA
1330  if(member_name != NULL)
1331  free(member_name);
1332 #endif
1333  if(retval) /* error exit from loop */
1334  return retval;
1335  }
1336  break;
1337 
1338  case H5T_VLEN:
1339  {
1340  htri_t ret;
1341 
1342  /* For conveninence we allow user to pass vlens of strings
1343  * with null terminated strings. This means strings are
1344  * treated slightly differently by the API, although they are
1345  * really just VLENs of characters. */
1346  if ((ret = H5Tis_variable_str(hdf_typeid)) < 0)
1347  return NC_EHDFERR;
1348  if (ret)
1349  type->nc_type_class = NC_STRING;
1350  else
1351  {
1352  hid_t base_hdf_typeid;
1353  nc_type base_nc_type = NC_NAT;
1354 
1355  type->nc_type_class = NC_VLEN;
1356 
1357  /* Find the base type of this vlen (i.e. what is this a
1358  * vlen of?) */
1359  if (!(base_hdf_typeid = H5Tget_super(native_typeid)))
1360  return NC_EHDFERR;
1361 
1362  /* What size is this type? */
1363  if (!(type_size = H5Tget_size(base_hdf_typeid)))
1364  return NC_EHDFERR;
1365 
1366  /* What is the netcdf corresponding type. */
1367  if ((retval = get_netcdf_type(grp->nc4_info, base_hdf_typeid,
1368  &base_nc_type)))
1369  return retval;
1370  LOG((5, "base_hdf_typeid 0x%x type_size %d base_nc_type %d",
1371  base_hdf_typeid, type_size, base_nc_type));
1372 
1373  /* Remember the base types for this vlen */
1374  type->u.v.base_nc_typeid = base_nc_type;
1375  type->u.v.base_hdf_typeid = base_hdf_typeid;
1376  }
1377  }
1378  break;
1379 
1380  case H5T_OPAQUE:
1381  type->nc_type_class = NC_OPAQUE;
1382  break;
1383 
1384  case H5T_ENUM:
1385  {
1386  hid_t base_hdf_typeid;
1387  nc_type base_nc_type = NC_NAT;
1388  void *value;
1389  int i;
1390  char *member_name = NULL;
1391 #ifdef JNA
1392  char jna[1001];
1393 #endif
1394 
1395  type->nc_type_class = NC_ENUM;
1396 
1397  /* Find the base type of this enum (i.e. what is this a
1398  * enum of?) */
1399  if (!(base_hdf_typeid = H5Tget_super(hdf_typeid)))
1400  return NC_EHDFERR;
1401  /* What size is this type? */
1402  if (!(type_size = H5Tget_size(base_hdf_typeid)))
1403  return NC_EHDFERR;
1404  /* What is the netcdf corresponding type. */
1405  if ((retval = get_netcdf_type(grp->nc4_info, base_hdf_typeid,
1406  &base_nc_type)))
1407  return retval;
1408  LOG((5, "base_hdf_typeid 0x%x type_size %d base_nc_type %d",
1409  base_hdf_typeid, type_size, base_nc_type));
1410 
1411  /* Remember the base types for this enum */
1412  type->u.e.base_nc_typeid = base_nc_type;
1413  type->u.e.base_hdf_typeid = base_hdf_typeid;
1414 
1415  /* Find out how many member are in the enum. */
1416  if ((type->u.e.num_members = H5Tget_nmembers(hdf_typeid)) < 0)
1417  return NC_EHDFERR;
1418 
1419  /* Allocate space for one value. */
1420  if (!(value = calloc(1, type_size)))
1421  return NC_ENOMEM;
1422 
1423  /* Read each name and value defined in the enum. */
1424  for (i = 0; i < type->u.e.num_members; i++)
1425  {
1426 
1427  /* Get the name and value from HDF5. */
1428  if (!(member_name = H5Tget_member_name(hdf_typeid, i)))
1429  {
1430  retval = NC_EHDFERR;
1431  break;
1432  }
1433 #ifdef JNA
1434  strncpy(jna,member_name,1000);
1435  member_name = jna;
1436 #endif
1437 
1438  if (strlen(member_name) > NC_MAX_NAME)
1439  {
1440  retval = NC_EBADNAME;
1441  break;
1442  }
1443  if (H5Tget_member_value(hdf_typeid, i, value) < 0)
1444  {
1445  retval = NC_EHDFERR;
1446  break;
1447  }
1448 
1449  /* Insert new field into this type's list of fields. */
1450  if ((retval = nc4_enum_member_add(&type->u.e.enum_member, type->size,
1451  member_name, value)))
1452  {
1453  break;
1454  }
1455 
1456 #ifndef JNA
1457  /* Free the member name (which HDF5 allocated for us). */
1458  if(member_name != NULL) free(member_name);
1459 #endif
1460  member_name = NULL;
1461  }
1462 
1463 #ifndef JNA
1464  if(member_name != NULL)
1465  free(member_name);
1466 #endif
1467  if(value) free(value);
1468  if(retval) /* error exit from loop */
1469  return retval;
1470  }
1471  break;
1472 
1473  default:
1474  LOG((0, "unknown class"));
1475  return NC_EBADCLASS;
1476  }
1477  return retval;
1478 }
1479 
1480 /* This function is called by read_dataset, (which is called by
1481  * nc4_rec_read_metadata) when a netCDF variable is found in the
1482  * file. This function reads in all the metadata about the var,
1483  * including the attributes. */
1484 static int
1485 read_var(NC_GRP_INFO_T *grp, hid_t datasetid, const char *obj_name,
1486  size_t ndims, NC_DIM_INFO_T *dim)
1487 {
1488  NC_VAR_INFO_T *var = NULL;
1489  hid_t access_pid = 0;
1490  int incr_id_rc = 0; /* Whether the dataset ID's ref count has been incremented */
1491  int natts, a, d;
1492 
1493  NC_ATT_INFO_T *att;
1494  hid_t attid = 0;
1495  char att_name[NC_MAX_HDF5_NAME + 1];
1496 
1497 #define CD_NELEMS_ZLIB 1
1498 #define CD_NELEMS_SZIP 4
1499  H5Z_filter_t filter;
1500  int num_filters;
1501  unsigned int cd_values[CD_NELEMS_SZIP];
1502  size_t cd_nelems = CD_NELEMS_SZIP;
1503  hid_t propid = 0;
1504  H5D_fill_value_t fill_status;
1505  H5D_layout_t layout;
1506  hsize_t chunksize[NC_MAX_VAR_DIMS] = {0};
1507  int retval = NC_NOERR;
1508  double rdcc_w0;
1509  int f;
1510 
1511  assert(obj_name && grp);
1512  LOG((4, "%s: obj_name %s", __func__, obj_name));
1513 
1514  /* Add a variable to the end of the group's var list. */
1515  if ((retval = nc4_var_list_add(&grp->var, &var)))
1516  BAIL(retval);
1517 
1518  /* Fill in what we already know. */
1519  var->hdf_datasetid = datasetid;
1520  H5Iinc_ref(var->hdf_datasetid); /* Increment number of objects using ID */
1521  incr_id_rc++; /* Indicate that we've incremented the ref. count (for errors) */
1522  var->varid = grp->nvars++;
1523  var->created = NC_TRUE;
1524  var->ndims = ndims;
1525 
1526  /* We need some room to store information about dimensions for this
1527  * var. */
1528  if (var->ndims)
1529  {
1530  if (!(var->dim = calloc(var->ndims, sizeof(NC_DIM_INFO_T *))))
1531  BAIL(NC_ENOMEM);
1532  if (!(var->dimids = calloc(var->ndims, sizeof(int))))
1533  BAIL(NC_ENOMEM);
1534  }
1535 
1536  /* Get the current chunk cache settings. */
1537  if ((access_pid = H5Dget_access_plist(datasetid)) < 0)
1538  BAIL(NC_EVARMETA);
1539 #ifdef EXTRA_TESTS
1540  num_plists++;
1541 #endif
1542 
1543  /* Learn about current chunk cache settings. */
1544  if ((H5Pget_chunk_cache(access_pid, &(var->chunk_cache_nelems),
1545  &(var->chunk_cache_size), &rdcc_w0)) < 0)
1546  BAIL(NC_EHDFERR);
1547  var->chunk_cache_preemption = rdcc_w0;
1548 
1549  /* Check for a weird case: a non-coordinate variable that has the
1550  * same name as a dimension. It's legal in netcdf, and requires
1551  * that the HDF5 dataset name be changed. */
1552  if (strlen(obj_name) > strlen(NON_COORD_PREPEND) &&
1553  !strncmp(obj_name, NON_COORD_PREPEND, strlen(NON_COORD_PREPEND)))
1554  {
1555  /* Allocate space for the name. */
1556  if (!(var->name = malloc(((strlen(obj_name) - strlen(NON_COORD_PREPEND))+ 1) * sizeof(char))))
1557  BAIL(NC_ENOMEM);
1558 
1559  strcpy(var->name, &obj_name[strlen(NON_COORD_PREPEND)]);
1560 
1561  /* Allocate space for the HDF5 name. */
1562  if (!(var->hdf5_name = malloc((strlen(obj_name) + 1) * sizeof(char))))
1563  BAIL(NC_ENOMEM);
1564 
1565  strcpy(var->hdf5_name, obj_name);
1566  }
1567  else
1568  {
1569  /* Allocate space for the name. */
1570  if (!(var->name = malloc((strlen(obj_name) + 1) * sizeof(char))))
1571  BAIL(NC_ENOMEM);
1572 
1573  strcpy(var->name, obj_name);
1574  }
1575 
1576  /* Find out what filters are applied to this HDF5 dataset,
1577  * fletcher32, deflate, and/or shuffle. All other filters are
1578  * ignored. */
1579  if ((propid = H5Dget_create_plist(datasetid)) < 0)
1580  BAIL(NC_EHDFERR);
1581 #ifdef EXTRA_TESTS
1582  num_plists++;
1583 #endif /* EXTRA_TESTS */
1584 
1585  /* Get the chunking info for non-scalar vars. */
1586  if ((layout = H5Pget_layout(propid)) < -1)
1587  BAIL(NC_EHDFERR);
1588  if (layout == H5D_CHUNKED)
1589  {
1590  if (H5Pget_chunk(propid, NC_MAX_VAR_DIMS, chunksize) < 0)
1591  BAIL(NC_EHDFERR);
1592  if (!(var->chunksizes = malloc(var->ndims * sizeof(size_t))))
1593  BAIL(NC_ENOMEM);
1594  for (d = 0; d < var->ndims; d++)
1595  var->chunksizes[d] = chunksize[d];
1596  }
1597  else if (layout == H5D_CONTIGUOUS || layout == H5D_COMPACT)
1598  var->contiguous = NC_TRUE;
1599 
1600  /* The possible values of filter (which is just an int) can be
1601  * found in H5Zpublic.h. */
1602  if ((num_filters = H5Pget_nfilters(propid)) < 0)
1603  BAIL(NC_EHDFERR);
1604  for (f = 0; f < num_filters; f++)
1605  {
1606  if ((filter = H5Pget_filter2(propid, f, NULL, &cd_nelems,
1607  cd_values, 0, NULL, NULL)) < 0)
1608  BAIL(NC_EHDFERR);
1609  switch (filter)
1610  {
1611  case H5Z_FILTER_SHUFFLE:
1612  var->shuffle = NC_TRUE;
1613  break;
1614 
1615  case H5Z_FILTER_FLETCHER32:
1616  var->fletcher32 = NC_TRUE;
1617  break;
1618 
1619  case H5Z_FILTER_DEFLATE:
1620  var->deflate = NC_TRUE;
1621  if (cd_nelems != CD_NELEMS_ZLIB || cd_values[0] > MAX_DEFLATE_LEVEL)
1622  BAIL(NC_EHDFERR);
1623  var->deflate_level = cd_values[0];
1624  break;
1625 
1626  case H5Z_FILTER_SZIP:
1627  var->szip = NC_TRUE;
1628  if (cd_nelems != CD_NELEMS_SZIP)
1629  BAIL(NC_EHDFERR);
1630  var->options_mask = cd_values[0];
1631  var->pixels_per_block = cd_values[1];
1632  break;
1633 
1634  default:
1635  LOG((1, "Yikes! Unknown filter type found on dataset!"));
1636  break;
1637  }
1638  }
1639 
1640  /* Learn all about the type of this variable. */
1641  if ((retval = get_type_info2(grp->nc4_info, datasetid,
1642  &var->type_info)))
1643  BAIL(retval);
1644 
1645  /* Indicate that the variable has a pointer to the type */
1646  var->type_info->rc++;
1647 
1648  /* Is there a fill value associated with this dataset? */
1649  if (H5Pfill_value_defined(propid, &fill_status) < 0)
1650  BAIL(NC_EHDFERR);
1651 
1652  /* Get the fill value, if there is one defined. */
1653  if (fill_status == H5D_FILL_VALUE_USER_DEFINED)
1654  {
1655  /* Allocate space to hold the fill value. */
1656  if (!var->fill_value)
1657  {
1658  if (var->type_info->nc_type_class == NC_VLEN)
1659  {
1660  if (!(var->fill_value = malloc(sizeof(nc_vlen_t))))
1661  BAIL(NC_ENOMEM);
1662  }
1663  else if (var->type_info->nc_type_class == NC_STRING)
1664  {
1665  if (!(var->fill_value = malloc(sizeof(char *))))
1666  BAIL(NC_ENOMEM);
1667  }
1668  else
1669  {
1670  assert(var->type_info->size);
1671  if (!(var->fill_value = malloc(var->type_info->size)))
1672  BAIL(NC_ENOMEM);
1673  }
1674  }
1675 
1676  /* Get the fill value from the HDF5 property lust. */
1677  if (H5Pget_fill_value(propid, var->type_info->native_hdf_typeid,
1678  var->fill_value) < 0)
1679  BAIL(NC_EHDFERR);
1680  }
1681  else
1682  var->no_fill = NC_TRUE;
1683 
1684  /* If it's a scale, mark it as such. */
1685  if (dim)
1686  {
1687  assert(ndims);
1688  var->dimscale = NC_TRUE;
1689  if (var->ndims > 1)
1690  {
1691  if ((retval = read_coord_dimids(var)))
1692  BAIL(retval);
1693  }
1694  else
1695  {
1696  /* sanity check */
1697  assert(0 == strcmp(var->name, dim->name));
1698 
1699  var->dimids[0] = dim->dimid;
1700  var->dim[0] = dim;
1701  }
1702  dim->coord_var = var;
1703  }
1704  /* If this is not a scale, but has scales, iterate
1705  * through them. (i.e. this is a variable that is not a
1706  * coordinate variable) */
1707  else
1708  {
1709  int num_scales = 0;
1710 
1711  /* Find out how many scales are attached to this
1712  * dataset. H5DSget_num_scales returns an error if there are no
1713  * scales, so convert a negative return value to zero. */
1714  num_scales = H5DSget_num_scales(datasetid, 0);
1715  if (num_scales < 0)
1716  num_scales = 0;
1717 
1718  if (num_scales && ndims)
1719  {
1720  /* Allocate space to remember whether the dimscale has been attached
1721  * for each dimension. */
1722  if (NULL == (var->dimscale_attached = calloc(ndims, sizeof(nc_bool_t))))
1723  BAIL(NC_ENOMEM);
1724 
1725  /* Store id information allowing us to match hdf5
1726  * dimscales to netcdf dimensions. */
1727  if (NULL == (var->dimscale_hdf5_objids = malloc(ndims * sizeof(struct hdf5_objid))))
1728  BAIL(NC_ENOMEM);
1729  for (d = 0; d < var->ndims; d++)
1730  {
1731  if (H5DSiterate_scales(var->hdf_datasetid, d, NULL, dimscale_visitor,
1732  &(var->dimscale_hdf5_objids[d])) < 0)
1733  BAIL(NC_EHDFERR);
1734  var->dimscale_attached[d] = NC_TRUE;
1735  }
1736  }
1737  }
1738 
1739  /* Now read all the attributes of this variable, ignoring the
1740  ones that hold HDF5 dimension scale information. */
1741  if ((natts = H5Aget_num_attrs(var->hdf_datasetid)) < 0)
1742  BAIL(NC_EATTMETA);
1743  for (a = 0; a < natts; a++)
1744  {
1745  /* Close the attribute and try to move on with our
1746  * lives. Like bits through the network port, so
1747  * flows the Days of Our Lives! */
1748  if (attid && H5Aclose(attid) < 0)
1749  BAIL(NC_EHDFERR);
1750 
1751  /* Open the att and get its name. */
1752  if ((attid = H5Aopen_idx(var->hdf_datasetid, (unsigned int)a)) < 0)
1753  BAIL(NC_EATTMETA);
1754  if (H5Aget_name(attid, NC_MAX_HDF5_NAME, att_name) < 0)
1755  BAIL(NC_EATTMETA);
1756  LOG((4, "%s:: a %d att_name %s", __func__, a, att_name));
1757 
1758  /* Should we ignore this attribute? */
1759  if (strcmp(att_name, REFERENCE_LIST) &&
1760  strcmp(att_name, CLASS) &&
1761  strcmp(att_name, DIMENSION_LIST) &&
1762  strcmp(att_name, NAME) &&
1763  strcmp(att_name, COORDINATES) &&
1764  strcmp(att_name, NC_DIMID_ATT_NAME))
1765  {
1766  /* Add to the end of the list of atts for this var. */
1767  if ((retval = nc4_att_list_add(&var->att, &att)))
1768  BAIL(retval);
1769 
1770  /* Fill in the information we know. */
1771  att->attnum = var->natts++;
1772  if (!(att->name = strdup(att_name)))
1773  BAIL(NC_ENOMEM);
1774 
1775  /* Read the rest of the info about the att,
1776  * including its values. */
1777  if ((retval = read_hdf5_att(grp, attid, att)))
1778  {
1779  if (NC_EBADTYPID == retval)
1780  {
1781  if ((retval = nc4_att_list_del(&var->att, att)))
1782  BAIL(retval);
1783  continue;
1784  }
1785  else
1786  BAIL(retval);
1787  }
1788 
1789  att->created = NC_TRUE;
1790  } /* endif not HDF5 att */
1791  } /* next attribute */
1792 
1793  /* Is this a deflated variable with a chunksize greater than the
1794  * current cache size? */
1795  if ((retval = nc4_adjust_var_cache(grp, var)))
1796  BAIL(retval);
1797 
1798 exit:
1799  if (retval)
1800  {
1801  if (incr_id_rc && H5Idec_ref(datasetid) < 0)
1802  BAIL2(NC_EHDFERR);
1803  if (var && nc4_var_list_del(&grp->var, var))
1804  BAIL2(NC_EHDFERR);
1805  }
1806  if (access_pid && H5Pclose(access_pid) < 0)
1807  BAIL2(NC_EHDFERR);
1808 #ifdef EXTRA_TESTS
1809  num_plists--;
1810 #endif
1811  if (propid > 0 && H5Pclose(propid) < 0)
1812  BAIL2(NC_EHDFERR);
1813 #ifdef EXTRA_TESTS
1814  num_plists--;
1815 #endif
1816  if (attid > 0 && H5Aclose(attid) < 0)
1817  BAIL2(NC_EHDFERR);
1818  return retval;
1819 }
1820 
1821 /* This function is called by nc4_rec_read_metadata to read all the
1822  * group level attributes (the NC_GLOBAL atts for this group). */
1823 static int
1824 read_grp_atts(NC_GRP_INFO_T *grp)
1825 {
1826  hid_t attid = 0;
1827  hsize_t num_obj, i;
1828  NC_ATT_INFO_T *att;
1829  NC_TYPE_INFO_T *type;
1830  char obj_name[NC_MAX_HDF5_NAME + 1];
1831  int max_len;
1832  int retval = NC_NOERR;
1833 
1834  num_obj = H5Aget_num_attrs(grp->hdf_grpid);
1835  for (i = 0; i < num_obj; i++)
1836  {
1837  /* Close an attribute from previous loop iteration */
1838  /* (Should be from 'continue' statement, below) */
1839  if (attid && H5Aclose(attid) < 0)
1840  BAIL(NC_EHDFERR);
1841 
1842  if ((attid = H5Aopen_idx(grp->hdf_grpid, (unsigned int)i)) < 0)
1843  BAIL(NC_EATTMETA);
1844  if (H5Aget_name(attid, NC_MAX_NAME + 1, obj_name) < 0)
1845  BAIL(NC_EATTMETA);
1846  LOG((3, "reading attribute of _netCDF group, named %s", obj_name));
1847 
1848  /* This may be an attribute telling us that strict netcdf-3
1849  * rules are in effect. If so, we will make note of the fact,
1850  * but not add this attribute to the metadata. It's not a user
1851  * attribute, but an internal netcdf-4 one. */
1852  if (!strcmp(obj_name, NC3_STRICT_ATT_NAME))
1853  grp->nc4_info->cmode |= NC_CLASSIC_MODEL;
1854  else
1855  {
1856  /* Add an att struct at the end of the list, and then go to it. */
1857  if ((retval = nc4_att_list_add(&grp->att, &att)))
1858  BAIL(retval);
1859 
1860  /* Add the info about this attribute. */
1861  max_len = strlen(obj_name) > NC_MAX_NAME ? NC_MAX_NAME : strlen(obj_name);
1862  if (!(att->name = malloc((max_len + 1) * sizeof(char))))
1863  BAIL(NC_ENOMEM);
1864  strncpy(att->name, obj_name, max_len);
1865  att->name[max_len] = 0;
1866  att->attnum = grp->natts++;
1867  if ((retval = read_hdf5_att(grp, attid, att)))
1868  {
1869  if (NC_EBADTYPID == retval)
1870  {
1871  if ((retval = nc4_att_list_del(&grp->att, att)))
1872  BAIL(retval);
1873  continue;
1874  }
1875  else
1876  BAIL(retval);
1877  }
1878  att->created = NC_TRUE;
1879  if ((retval = nc4_find_type(grp->nc4_info, att->nc_typeid, &type)))
1880  BAIL(retval);
1881  }
1882  }
1883 
1884  exit:
1885  if (attid > 0 && H5Aclose(attid) < 0)
1886  BAIL2(NC_EHDFERR);
1887  return retval;
1888 }
1889 
1890 /* This function is called when nc4_rec_read_metadata encounters an HDF5
1891  * dataset when reading a file. */
1892 static int
1893 read_dataset(NC_GRP_INFO_T *grp, hid_t datasetid, const char *obj_name,
1894  const H5G_stat_t *statbuf)
1895 {
1896  NC_DIM_INFO_T *dim = NULL; /* Dimension created for scales */
1897  hid_t spaceid = 0;
1898  int ndims;
1899  htri_t is_scale;
1900  int retval = NC_NOERR;
1901 
1902  /* Get the dimension information for this dataset. */
1903  if ((spaceid = H5Dget_space(datasetid)) < 0)
1904  BAIL(NC_EHDFERR);
1905 #ifdef EXTRA_TESTS
1906  num_spaces++;
1907 #endif
1908  if ((ndims = H5Sget_simple_extent_ndims(spaceid)) < 0)
1909  BAIL(NC_EHDFERR);
1910 
1911  /* Is this a dimscale? */
1912  if ((is_scale = H5DSis_scale(datasetid)) < 0)
1913  BAIL(NC_EHDFERR);
1914  if (is_scale)
1915  {
1916  hsize_t dims[H5S_MAX_RANK];
1917  hsize_t max_dims[H5S_MAX_RANK];
1918 
1919  /* Query the scale's size & max. size */
1920  if (H5Sget_simple_extent_dims(spaceid, dims, max_dims) < 0)
1921  BAIL(NC_EHDFERR);
1922 
1923  /* Read the scale information. */
1924  if ((retval = read_scale(grp, datasetid, obj_name, statbuf, dims[0],
1925  max_dims[0], &dim)))
1926  BAIL(retval);
1927  }
1928 
1929  /* Add a var to the linked list, and get its metadata,
1930  * unless this is one of those funny dimscales that are a
1931  * dimension in netCDF but not a variable. (Spooky!) */
1932  if (NULL == dim || (dim && !dim->hdf_dimscaleid))
1933  if ((retval = read_var(grp, datasetid, obj_name, ndims, dim)))
1934  BAIL(retval);
1935 
1936 exit:
1937  if (spaceid && H5Sclose(spaceid) <0)
1938  BAIL2(retval);
1939 #ifdef EXTRA_TESTS
1940  num_spaces--;
1941 #endif
1942 
1943  return retval;
1944 }
1945 
1946 static int
1947 nc4_rec_read_metadata_cb_list_add(NC4_rec_read_metadata_obj_info_t **head,
1949  const NC4_rec_read_metadata_obj_info_t *oinfo)
1950 {
1951  NC4_rec_read_metadata_obj_info_t *new_oinfo; /* Pointer to info for object */
1952 
1953  /* Allocate memory for the object's info */
1954  if (!(new_oinfo = calloc(1, sizeof(*new_oinfo))))
1955  return NC_ENOMEM;
1956 
1957  /* Make a copy of the object's info */
1958  memcpy(new_oinfo, oinfo, sizeof(*oinfo));
1959 
1960  if (*tail)
1961  {
1962  assert(*head);
1963  (*tail)->next = new_oinfo;
1964  *tail = new_oinfo;
1965  }
1966  else
1967  {
1968  assert(NULL == *head);
1969  *head = *tail = new_oinfo;
1970  }
1971 
1972  return (NC_NOERR);
1973 }
1974 
1975 static int
1976 nc4_rec_read_metadata_cb(hid_t grpid, const char *name, const H5L_info_t *info,
1977  void *_op_data)
1978 {
1979  NC4_rec_read_metadata_ud_t *udata = (NC4_rec_read_metadata_ud_t *)_op_data; /* Pointer to user data for callback */
1980  NC4_rec_read_metadata_obj_info_t oinfo; /* Pointer to info for object */
1981  int retval = H5_ITER_CONT;
1982 
1983  /* Reset the memory for the object's info */
1984  memset(&oinfo, 0, sizeof(oinfo));
1985 
1986  /* Open this critter. */
1987  if ((oinfo.oid = H5Oopen(grpid, name, H5P_DEFAULT)) < 0)
1988  BAIL(H5_ITER_ERROR);
1989 
1990  /* Get info about the object.*/
1991  if (H5Gget_objinfo(oinfo.oid, ".", 1, &oinfo.statbuf) < 0)
1992  BAIL(H5_ITER_ERROR);
1993 
1994  strncpy(oinfo.oname, name, NC_MAX_NAME);
1995 
1996  /* Add object to list, for later */
1997  switch(oinfo.statbuf.type)
1998  {
1999  case H5G_GROUP:
2000  LOG((3, "found group %s", oinfo.oname));
2001 
2002  /* Defer descending into child group immediately, so that the types
2003  * in the current group can be processed and be ready for use by
2004  * vars in the child group(s).
2005  */
2006  if (nc4_rec_read_metadata_cb_list_add(&udata->grps_head, &udata->grps_tail, &oinfo))
2007  BAIL(H5_ITER_ERROR);
2008  break;
2009 
2010  case H5G_DATASET:
2011  LOG((3, "found dataset %s", oinfo.oname));
2012 
2013  /* Learn all about this dataset, which may be a dimscale
2014  * (i.e. dimension metadata), or real data. */
2015  if ((retval = read_dataset(udata->grp, oinfo.oid, oinfo.oname, &oinfo.statbuf)))
2016  {
2017  /* Allow NC_EBADTYPID to transparently skip over datasets
2018  * which have a datatype that netCDF-4 doesn't undertand
2019  * (currently), but break out of iteration for other
2020  * errors.
2021  */
2022  if(NC_EBADTYPID != retval)
2023  BAIL(H5_ITER_ERROR);
2024  else
2025  retval = H5_ITER_CONT;
2026  }
2027 
2028  /* Close the object */
2029  if (H5Oclose(oinfo.oid) < 0)
2030  BAIL(H5_ITER_ERROR);
2031  break;
2032 
2033  case H5G_TYPE:
2034  LOG((3, "found datatype %s", oinfo.oname));
2035 
2036  /* Process the named datatype */
2037  if (read_type(udata->grp, oinfo.oid, oinfo.oname))
2038  BAIL(H5_ITER_ERROR);
2039 
2040  /* Close the object */
2041  if (H5Oclose(oinfo.oid) < 0)
2042  BAIL(H5_ITER_ERROR);
2043  break;
2044 
2045  default:
2046  LOG((0, "Unknown object class %d in %s!", oinfo.statbuf.type, __func__));
2047  BAIL(H5_ITER_ERROR);
2048  }
2049 
2050 exit:
2051  if (retval)
2052  {
2053  if (oinfo.oid > 0 && H5Oclose(oinfo.oid) < 0)
2054  BAIL2(H5_ITER_ERROR);
2055  }
2056 
2057  return (retval);
2058 }
2059 
2060 /* This is the main function to recursively read all the metadata for the file. */
2061 /* The links in the 'grp' are iterated over and added to the file's metadata
2062  * information. Note that child groups are not immediately processed, but
2063  * are deferred until all the other links in the group are handled (so that
2064  * vars in the child groups are guaranteed to have types that they use in
2065  * a parent group in place).
2066  */
2067 static int
2068 nc4_rec_read_metadata(NC_GRP_INFO_T *grp)
2069 {
2070  NC4_rec_read_metadata_ud_t udata; /* User data for iteration */
2071  NC4_rec_read_metadata_obj_info_t *oinfo; /* Pointer to info for object */
2072  hsize_t idx=0;
2073  hid_t pid = 0;
2074  unsigned crt_order_flags = 0;
2075  H5_index_t iter_index;
2076  int retval = NC_NOERR; /* everything worked! */
2077 
2078  assert(grp && grp->name);
2079  LOG((3, "%s: grp->name %s", __func__, grp->name));
2080 
2081  /* Portably initialize user data for later */
2082  memset(&udata, 0, sizeof(udata));
2083 
2084  /* Open this HDF5 group and retain its grpid. It will remain open
2085  * with HDF5 until this file is nc_closed. */
2086  if (!grp->hdf_grpid)
2087  {
2088  if (grp->parent)
2089  {
2090  if ((grp->hdf_grpid = H5Gopen2(grp->parent->hdf_grpid,
2091  grp->name, H5P_DEFAULT)) < 0)
2092  BAIL(NC_EHDFERR);
2093  }
2094  else
2095  {
2096  if ((grp->hdf_grpid = H5Gopen2(grp->nc4_info->hdfid,
2097  "/", H5P_DEFAULT)) < 0)
2098  BAIL(NC_EHDFERR);
2099  }
2100  }
2101  assert(grp->hdf_grpid > 0);
2102 
2103  /* Get the group creation flags, to check for creation ordering */
2104  pid = H5Gget_create_plist(grp->hdf_grpid);
2105  H5Pget_link_creation_order(pid, &crt_order_flags);
2106  if (H5Pclose(pid) < 0)
2107  BAIL(NC_EHDFERR);
2108 
2109  /* Set the iteration index to use */
2110  if (crt_order_flags & H5P_CRT_ORDER_TRACKED)
2111  iter_index = H5_INDEX_CRT_ORDER;
2112  else
2113  {
2114  NC_HDF5_FILE_INFO_T *h5 = grp->nc4_info;
2115 
2116  /* Without creation ordering, file must be read-only. */
2117  if (!h5->no_write)
2118  BAIL(NC_ECANTWRITE);
2119 
2120  iter_index = H5_INDEX_NAME;
2121  }
2122 
2123  /* Set user data for iteration */
2124  udata.grp = grp;
2125 
2126  /* Iterate over links in this group, building lists for the types,
2127  * datasets and groups encountered
2128  */
2129  if (H5Literate(grp->hdf_grpid, iter_index, H5_ITER_INC, &idx,
2130  nc4_rec_read_metadata_cb, (void *)&udata) < 0)
2131  BAIL(NC_EHDFERR);
2132 
2133  /* Process the child groups found */
2134  /* (Deferred until now, so that the types in the current group get
2135  * processed and are available for vars in the child group(s).)
2136  */
2137  for (oinfo = udata.grps_head; oinfo; oinfo = udata.grps_head)
2138  {
2139  NC_GRP_INFO_T *child_grp;
2140  NC_HDF5_FILE_INFO_T *h5 = grp->nc4_info;
2141 
2142  /* Add group to file's hierarchy */
2143  if ((retval = nc4_grp_list_add(&(grp->children), h5->next_nc_grpid++,
2144  grp, grp->nc4_info->controller, oinfo->oname, &child_grp)))
2145  BAIL(retval);
2146 
2147  /* Recursively read the child group's metadata */
2148  if ((retval = nc4_rec_read_metadata(child_grp)))
2149  BAIL(retval);
2150 
2151  /* Close the object */
2152  if (H5Oclose(oinfo->oid) < 0)
2153  BAIL(NC_EHDFERR);
2154 
2155  /* Advance to next node, free current node */
2156  udata.grps_head = oinfo->next;
2157  free(oinfo);
2158  }
2159 
2160  /* Scan the group for global (i.e. group-level) attributes. */
2161  if ((retval = read_grp_atts(grp)))
2162  BAIL(retval);
2163 
2164 exit:
2165  /* Clean up local information on error, if anything remains */
2166  if (retval)
2167  {
2168  for (oinfo = udata.grps_head; oinfo; oinfo = udata.grps_head)
2169  {
2170  /* Close the object */
2171  if (H5Oclose(oinfo->oid) < 0)
2172  BAIL2(NC_EHDFERR);
2173 
2174  /* Advance to next node, free current node */
2175  udata.grps_head = oinfo->next;
2176  free(oinfo);
2177  }
2178  }
2179 
2180  return retval;
2181 }
2182 
2183 /* Open a netcdf-4 file. Things have already been kicked off in
2184  * ncfunc.c in nc_open, but here the netCDF-4 part of opening a file
2185  * is handled. */
2186 static int
2187 nc4_open_file(const char *path, int mode, void* parameters, NC *nc)
2188 {
2189  hid_t fapl_id = H5P_DEFAULT;
2190  unsigned flags = (mode & NC_WRITE) ?
2191  H5F_ACC_RDWR : H5F_ACC_RDONLY;
2192  int retval;
2193  NC_HDF5_FILE_INFO_T* nc4_info = NULL;
2194  int inmemory = ((mode & NC_INMEMORY) == NC_INMEMORY);
2195 #ifdef USE_DISKLESS
2196  NC_MEM_INFO* meminfo = (NC_MEM_INFO*)parameters;
2197 #endif
2198 #ifdef USE_PARALLEL4
2199  NC_MPI_INFO* mpiinfo = (NC_MPI_INFO*)parameters;
2200  int comm_duped = 0; /* Whether the MPI Communicator was duplicated */
2201  int info_duped = 0; /* Whether the MPI Info object was duplicated */
2202 #endif /* !USE_PARALLEL4 */
2203 
2204  LOG((3, "%s: path %s mode %d", __func__, path, mode));
2205  assert(path && nc);
2206 
2207  /* Add necessary structs to hold netcdf-4 file data. */
2208  if ((retval = nc4_nc4f_list_add(nc, path, mode)))
2209  BAIL(retval);
2210  nc4_info = NC4_DATA(nc);
2211  assert(nc4_info && nc4_info->root_grp);
2212 
2213  /* Need this access plist to control how HDF5 handles open objects
2214  * on file close. (Setting H5F_CLOSE_SEMI will cause H5Fclose to
2215  * fail if there are any open objects in the file. */
2216  if ((fapl_id = H5Pcreate(H5P_FILE_ACCESS)) < 0)
2217  BAIL(NC_EHDFERR);
2218 #ifdef EXTRA_TESTS
2219  num_plists++;
2220 #endif
2221 #ifdef EXTRA_TESTS
2222  if (H5Pset_fclose_degree(fapl_id, H5F_CLOSE_SEMI))
2223  BAIL(NC_EHDFERR);
2224 #else
2225  if (H5Pset_fclose_degree(fapl_id, H5F_CLOSE_STRONG))
2226  BAIL(NC_EHDFERR);
2227 #endif
2228 
2229 #ifdef USE_PARALLEL4
2230  /* If this is a parallel file create, set up the file creation
2231  property list. */
2232  if (mode & NC_MPIIO || mode & NC_MPIPOSIX)
2233  {
2234  nc4_info->parallel = NC_TRUE;
2235  if (mode & NC_MPIIO) /* MPI/IO */
2236  {
2237  LOG((4, "opening parallel file with MPI/IO"));
2238  if (H5Pset_fapl_mpio(fapl_id, mpiinfo->comm, mpiinfo->info) < 0)
2239  BAIL(NC_EPARINIT);
2240  }
2241 #ifdef USE_PARALLEL_POSIX
2242  else /* MPI/POSIX */
2243  {
2244  LOG((4, "opening parallel file with MPI/posix"));
2245  if (H5Pset_fapl_mpiposix(fapl_id, mpiinfo->comm, 0) < 0)
2246  BAIL(NC_EPARINIT);
2247  }
2248 #else /* USE_PARALLEL_POSIX */
2249  /* Should not happen! Code in NC4_create/NC4_open should alias the
2250  * NC_MPIPOSIX flag to NC_MPIIO, if the MPI-POSIX VFD is not
2251  * available in HDF5. -QAK
2252  */
2253  else /* MPI/POSIX */
2254  BAIL(NC_EPARINIT);
2255 #endif /* USE_PARALLEL_POSIX */
2256 
2257  /* Keep copies of the MPI Comm & Info objects */
2258  if (MPI_SUCCESS != MPI_Comm_dup(mpiinfo->comm, &nc4_info->comm))
2259  BAIL(NC_EMPI);
2260  comm_duped++;
2261  if (MPI_INFO_NULL != mpiinfo->info)
2262  {
2263  if (MPI_SUCCESS != MPI_Info_dup(mpiinfo->info, &nc4_info->info))
2264  BAIL(NC_EMPI);
2265  info_duped++;
2266  }
2267  else
2268  {
2269  /* No dup, just copy it. */
2270  nc4_info->info = mpiinfo->info;
2271  }
2272  }
2273 #else /* only set cache for non-parallel. */
2274  if (H5Pset_cache(fapl_id, 0, nc4_chunk_cache_nelems, nc4_chunk_cache_size,
2275  nc4_chunk_cache_preemption) < 0)
2276  BAIL(NC_EHDFERR);
2277  LOG((4, "%s: set HDF raw chunk cache to size %d nelems %d preemption %f",
2278  __func__, nc4_chunk_cache_size, nc4_chunk_cache_nelems, nc4_chunk_cache_preemption));
2279 #endif /* USE_PARALLEL4 */
2280 
2281  /* The NetCDF-3.x prototype contains an mode option NC_SHARE for
2282  multiple processes accessing the dataset concurrently. As there
2283  is no HDF5 equivalent, NC_SHARE is treated as NC_NOWRITE. */
2284  if(inmemory) {
2285  if((nc4_info->hdfid = H5LTopen_file_image(meminfo->memory,meminfo->size,
2286  H5LT_FILE_IMAGE_DONT_COPY|H5LT_FILE_IMAGE_DONT_RELEASE
2287  )) < 0)
2288  BAIL(NC_EHDFERR);
2289  nc4_info->no_write = NC_TRUE;
2290  } else if ((nc4_info->hdfid = H5Fopen(path, flags, fapl_id)) < 0)
2291  BAIL(NC_EHDFERR);
2292 
2293  /* Does the mode specify that this file is read-only? */
2294  if ((mode & NC_WRITE) == 0)
2295  nc4_info->no_write = NC_TRUE;
2296 
2297  /* Now read in all the metadata. Some types and dimscale
2298  * information may be difficult to resolve here, if, for example, a
2299  * dataset of user-defined type is encountered before the
2300  * definition of that type. */
2301  if ((retval = nc4_rec_read_metadata(nc4_info->root_grp)))
2302  BAIL(retval);
2303 
2304  /* Now figure out which netCDF dims are indicated by the dimscale
2305  * information. */
2306  if ((retval = nc4_rec_match_dimscales(nc4_info->root_grp)))
2307  BAIL(retval);
2308 
2309 #ifdef LOGGING
2310  /* This will print out the names, types, lens, etc of the vars and
2311  atts in the file, if the logging level is 2 or greater. */
2312  log_metadata_nc(nc);
2313 #endif
2314 
2315  /* Close the property list. */
2316  if (H5Pclose(fapl_id) < 0)
2317  BAIL(NC_EHDFERR);
2318 #ifdef EXTRA_TESTS
2319  num_plists--;
2320 #endif
2321 
2322  return NC_NOERR;
2323 
2324 exit:
2325 #ifdef USE_PARALLEL4
2326  if (comm_duped) MPI_Comm_free(&nc4_info->comm);
2327  if (info_duped) MPI_Info_free(&nc4_info->info);
2328 #endif
2329 #ifdef EXTRA_TESTS
2330  num_plists--;
2331 #endif
2332  if (fapl_id != H5P_DEFAULT) H5Pclose(fapl_id);
2333  if (!nc4_info) return retval;
2334  close_netcdf4_file(nc4_info,1); /* treat like abort*/
2335  return retval;
2336 }
2337 
2338 #ifdef USE_HDF4
2339 
2352 static int
2353 get_netcdf_type_from_hdf4(NC_HDF5_FILE_INFO_T *h5, int32 hdf4_typeid,
2354  nc_type *xtype, NC_TYPE_INFO_T *type_info)
2355 {
2356  int t = 0;
2357 
2358  /* Added this variable in the course of fixing NCF-332.
2359  * Prior to the fix, all data types were assigned
2360  * NC_ENDIAN_BIG, so I am preserving that here for now.
2361  * Not sure why it wouldn't be NC_ENDIAN_NATIVE, although
2362  * I can hazard a guess or two.
2363  */
2364 
2365  int endianness = NC_ENDIAN_BIG;
2366  assert(h5 && xtype);
2367 
2368  switch(hdf4_typeid)
2369  {
2370  case DFNT_CHAR:
2371  *xtype = NC_CHAR;
2372  t = 0;
2373  break;
2374  case DFNT_UCHAR:
2375  case DFNT_UINT8:
2376  *xtype = NC_UBYTE;
2377  t = 6;
2378  break;
2379  case DFNT_LUINT8:
2380  *xtype = NC_UBYTE;
2381  t = 6;
2382  endianness = NC_ENDIAN_LITTLE;
2383  break;
2384  case DFNT_INT8:
2385  *xtype = NC_BYTE;
2386  t = 1;
2387  break;
2388  case DFNT_LINT8:
2389  *xtype = NC_BYTE;
2390  t = 1;
2391  endianness = NC_ENDIAN_LITTLE;
2392  break;
2393  case DFNT_INT16:
2394  *xtype = NC_SHORT;
2395  t = 2;
2396  break;
2397  case DFNT_LINT16:
2398  *xtype = NC_SHORT;
2399  t = 2;
2400  endianness = NC_ENDIAN_LITTLE;
2401  break;
2402  case DFNT_UINT16:
2403  *xtype = NC_USHORT;
2404  t = 7;
2405  break;
2406  case DFNT_LUINT16:
2407  *xtype = NC_USHORT;
2408  t = 7;
2409  endianness = NC_ENDIAN_LITTLE;
2410  break;
2411  case DFNT_INT32:
2412  *xtype = NC_INT;
2413  t = 3;
2414  break;
2415  case DFNT_LINT32:
2416  *xtype = NC_INT;
2417  t = 3;
2418  endianness = NC_ENDIAN_LITTLE;
2419  break;
2420  case DFNT_UINT32:
2421  *xtype = NC_UINT;
2422  t = 8;
2423  break;
2424  case DFNT_LUINT32:
2425  *xtype = NC_UINT;
2426  t = 8;
2427  endianness = NC_ENDIAN_LITTLE;
2428  break;
2429  case DFNT_FLOAT32:
2430  *xtype = NC_FLOAT;
2431  t = 4;
2432  break;
2433  case DFNT_LFLOAT32:
2434  *xtype = NC_FLOAT;
2435  t = 4;
2436  endianness = NC_ENDIAN_LITTLE;
2437  break;
2438  case DFNT_FLOAT64:
2439  *xtype = NC_DOUBLE;
2440  t = 5;
2441  break;
2442  case DFNT_LFLOAT64:
2443  *xtype = NC_DOUBLE;
2444  t = 5;
2445  endianness = NC_ENDIAN_LITTLE;
2446  break;
2447  default:
2448  *xtype = NC_NAT;
2449  return NC_EBADTYPID;
2450  }
2451 
2452  if (type_info)
2453  {
2454  if (hdf4_typeid == DFNT_FLOAT32)
2455  type_info->nc_type_class = NC_FLOAT;
2456  else if (hdf4_typeid == DFNT_FLOAT64)
2457  type_info->nc_type_class = NC_DOUBLE;
2458  else if (hdf4_typeid == DFNT_CHAR)
2459  type_info->nc_type_class = NC_STRING;
2460  else
2461  type_info->nc_type_class = NC_INT;
2462  type_info->endianness = endianness;
2463  type_info->nc_typeid = *xtype;
2464  type_info->size = nc_type_size_g[t];
2465  if (!(type_info->name = strdup(nc_type_name_g[t])))
2466  return NC_ENOMEM;
2467  }
2468 
2469  return NC_NOERR;
2470 }
2471 
2472 /* Open a HDF4 file. Things have already been kicked off in nc_open,
2473  * but here the netCDF-4 part of opening a file is handled. */
2474 static int
2475 nc4_open_hdf4_file(const char *path, int mode, NC *nc)
2476 {
2477  NC_HDF5_FILE_INFO_T *h5;
2478  NC_GRP_INFO_T *grp;
2479  NC_ATT_INFO_T *att;
2480  int32 num_datasets, num_gatts;
2481  int32 rank;
2482  int v, d, a;
2483  int retval;
2484  NC_HDF5_FILE_INFO_T* nc4_info = NULL;
2485 
2486  LOG((3, "%s: path %s mode %d", __func__, path, mode));
2487  assert(path && nc);
2488 
2489  /* Must be read-only access to hdf4 files. */
2490  if (mode & NC_WRITE)
2491  return NC_EINVAL;
2492 
2493  /* Add necessary structs to hold netcdf-4 file data. */
2494  if ((retval = nc4_nc4f_list_add(nc, path, mode)))
2495  return retval;
2496  nc4_info = NC4_DATA(nc);
2497  assert(nc4_info && nc4_info->root_grp);
2498  h5 = nc4_info;
2499  h5->hdf4 = NC_TRUE;
2500  grp = h5->root_grp;
2501  h5->no_write = NC_TRUE;
2502 
2503  /* Open the file and initialize SD interface. */
2504  if ((h5->sdid = SDstart(path, DFACC_READ)) == FAIL)
2505  return NC_EHDFERR;
2506 
2507  /* Learn how many datasets and global atts we have. */
2508  if (SDfileinfo(h5->sdid, &num_datasets, &num_gatts))
2509  return NC_EHDFERR;
2510 
2511  /* Read the atts. */
2512  for (a = 0; a < num_gatts; a++)
2513  {
2514  int32 att_data_type, att_count;
2515  size_t att_type_size;
2516 
2517  /* Add to the end of the list of atts for this var. */
2518  if ((retval = nc4_att_list_add(&h5->root_grp->att, &att)))
2519  return retval;
2520  att->attnum = grp->natts++;
2521  att->created = NC_TRUE;
2522 
2523  /* Learn about this attribute. */
2524  if (!(att->name = malloc(NC_MAX_HDF4_NAME * sizeof(char))))
2525  return NC_ENOMEM;
2526  if (SDattrinfo(h5->sdid, a, att->name, &att_data_type, &att_count))
2527  return NC_EATTMETA;
2528  if ((retval = get_netcdf_type_from_hdf4(h5, att_data_type,
2529  &att->nc_typeid, NULL)))
2530  return retval;
2531  att->len = att_count;
2532 
2533  /* Allocate memory to hold the data. */
2534  if ((retval = nc4_get_typelen_mem(h5, att->nc_typeid, 0, &att_type_size)))
2535  return retval;
2536  if (!(att->data = malloc(att_type_size * att->len)))
2537  return NC_ENOMEM;
2538 
2539  /* Read the data. */
2540  if (SDreadattr(h5->sdid, a, att->data))
2541  return NC_EHDFERR;
2542  }
2543 
2544  /* Read each dataset. */
2545  for (v = 0; v < num_datasets; v++)
2546  {
2547  NC_VAR_INFO_T *var;
2548  int32 data_type, num_atts;
2549  /* Problem: Number of dims is returned by the call that requires
2550  a pre-allocated array, 'dimsize'.
2551  From SDS_SD website:
2552  http://www.hdfgroup.org/training/HDFtraining/UsersGuide/SDS_SD.fm3.html
2553  The maximum rank is 32, or MAX_VAR_DIMS (as defined in netcdf.h).
2554 
2555  int32 dimsize[MAX_VAR_DIMS];
2556  */
2557  int32 *dimsize = NULL;
2558  size_t var_type_size;
2559  int a;
2560 
2561  /* Add a variable to the end of the group's var list. */
2562  if ((retval = nc4_var_list_add(&grp->var, &var)))
2563  return retval;
2564 
2565  var->varid = grp->nvars++;
2566  var->created = NC_TRUE;
2567  var->written_to = NC_TRUE;
2568 
2569  /* Open this dataset in HDF4 file. */
2570  if ((var->sdsid = SDselect(h5->sdid, v)) == FAIL)
2571  return NC_EVARMETA;
2572 
2573  /* Get shape, name, type, and attribute info about this dataset. */
2574  if (!(var->name = malloc(NC_MAX_HDF4_NAME + 1)))
2575  return NC_ENOMEM;
2576 
2577  /* Invoke SDgetInfo with null dimsize to get rank. */
2578  if (SDgetinfo(var->sdsid, var->name, &rank, NULL, &data_type, &num_atts))
2579  return NC_EVARMETA;
2580 
2581  if(!(dimsize = (int32*)malloc(sizeof(int32)*rank)))
2582  return NC_ENOMEM;
2583 
2584  if (SDgetinfo(var->sdsid, var->name, &rank, dimsize, &data_type, &num_atts)) {
2585  if(dimsize) free(dimsize);
2586  return NC_EVARMETA;
2587  }
2588 
2589  var->ndims = rank;
2590  var->hdf4_data_type = data_type;
2591 
2592  /* Fill special type_info struct for variable type information. */
2593  if (!(var->type_info = calloc(1, sizeof(NC_TYPE_INFO_T)))) {
2594  if(dimsize) free(dimsize);
2595  return NC_ENOMEM;
2596  }
2597 
2598  if ((retval = get_netcdf_type_from_hdf4(h5, data_type, &var->type_info->nc_typeid, var->type_info))) {
2599  if(dimsize) free(dimsize);
2600  return retval;
2601  }
2602 
2603  /* Indicate that the variable has a pointer to the type */
2604  var->type_info->rc++;
2605 
2606  if ((retval = nc4_get_typelen_mem(h5, var->type_info->nc_typeid, 0, &var_type_size))) {
2607  if(dimsize) free(dimsize);
2608  return retval;
2609  }
2610 
2611  var->type_info->size = var_type_size;
2612  LOG((3, "reading HDF4 dataset %s, rank %d netCDF type %d", var->name,
2613  rank, var->type_info->nc_typeid));
2614 
2615  /* Get the fill value. */
2616  if (!(var->fill_value = malloc(var_type_size))) {
2617  if(dimsize) free(dimsize);
2618  return NC_ENOMEM;
2619  }
2620 
2621  if (SDgetfillvalue(var->sdsid, var->fill_value))
2622  {
2623  /* Whoops! No fill value! */
2624  free(var->fill_value);
2625  var->fill_value = NULL;
2626  }
2627 
2628  /* Allocate storage for dimension info in this variable. */
2629  if (var->ndims)
2630  {
2631  if (!(var->dim = malloc(sizeof(NC_DIM_INFO_T *) * var->ndims))) {
2632  if(dimsize) free(dimsize);
2633  return NC_ENOMEM;
2634  }
2635 
2636  if (!(var->dimids = malloc(sizeof(int) * var->ndims))) {
2637  if(dimsize) free(dimsize);
2638  return NC_ENOMEM;
2639  }
2640  }
2641 
2642 
2643  /* Find its dimensions. */
2644  for (d = 0; d < var->ndims; d++)
2645  {
2646  int32 dimid, dim_len, dim_data_type, dim_num_attrs;
2647  char dim_name[NC_MAX_NAME + 1];
2648  NC_DIM_INFO_T *dim;
2649 
2650  if ((dimid = SDgetdimid(var->sdsid, d)) == FAIL) {
2651  if(dimsize) free(dimsize);
2652  return NC_EDIMMETA;
2653  }
2654  if (SDdiminfo(dimid, dim_name, &dim_len, &dim_data_type,
2655  &dim_num_attrs))
2656  {
2657  if(dimsize) free(dimsize);
2658  return NC_EDIMMETA;
2659  }
2660 
2661  /* Do we already have this dimension? HDF4 explicitly uses
2662  * the name to tell. */
2663  for (dim = grp->dim; dim; dim = dim->l.next)
2664  if (!strcmp(dim->name, dim_name))
2665  break;
2666 
2667  /* If we didn't find this dimension, add one. */
2668  if (!dim)
2669  {
2670  LOG((4, "adding dimension %s for HDF4 dataset %s",
2671  dim_name, var->name));
2672  if ((retval = nc4_dim_list_add(&grp->dim, &dim)))
2673  return retval;
2674  grp->ndims++;
2675  dim->dimid = grp->nc4_info->next_dimid++;
2676  if (strlen(dim_name) > NC_MAX_HDF4_NAME)
2677  return NC_EMAXNAME;
2678  if (!(dim->name = strdup(dim_name)))
2679  return NC_ENOMEM;
2680  if (dim_len)
2681  dim->len = dim_len;
2682  else
2683  dim->len = *dimsize;
2684  }
2685 
2686  /* Tell the variable the id of this dimension. */
2687  var->dimids[d] = dim->dimid;
2688  }
2689 
2690  /* Read the atts. */
2691  for (a = 0; a < num_atts; a++)
2692  {
2693  int32 att_data_type, att_count;
2694  size_t att_type_size;
2695 
2696  /* Add to the end of the list of atts for this var. */
2697  if ((retval = nc4_att_list_add(&var->att, &att))) {
2698  if(dimsize) free(dimsize);
2699  return retval;
2700  }
2701  att->attnum = var->natts++;
2702  att->created = NC_TRUE;
2703 
2704  /* Learn about this attribute. */
2705  if (!(att->name = malloc(NC_MAX_HDF4_NAME * sizeof(char)))) {
2706  if(dimsize) free(dimsize);
2707  return NC_ENOMEM;
2708  }
2709  if (SDattrinfo(var->sdsid, a, att->name, &att_data_type, &att_count)) {
2710  if(dimsize) free(dimsize);
2711  return NC_EATTMETA;
2712  }
2713  if ((retval = get_netcdf_type_from_hdf4(h5, att_data_type,
2714  &att->nc_typeid, NULL))) {
2715  if(dimsize) free(dimsize);
2716  return retval;
2717  }
2718 
2719  att->len = att_count;
2720 
2721  /* Allocate memory to hold the data. */
2722  if ((retval = nc4_get_typelen_mem(h5, att->nc_typeid, 0, &att_type_size))) {
2723  if(dimsize) free(dimsize);
2724  return retval;
2725  }
2726  if (!(att->data = malloc(att_type_size * att->len))) {
2727  if(dimsize) free(dimsize);
2728  return NC_ENOMEM;
2729  }
2730 
2731  /* Read the data. */
2732  if (SDreadattr(var->sdsid, a, att->data)) {
2733  if(dimsize) free(dimsize);
2734  return NC_EHDFERR;
2735  }
2736  }
2737  if(dimsize) free(dimsize);
2738 
2739  {
2740  /* HDF4 files can be chunked */
2741  HDF_CHUNK_DEF chunkdefs;
2742  int flag;
2743  if(!SDgetchunkinfo(var->sdsid, &chunkdefs, &flag)) {
2744  if(flag == HDF_NONE)
2745  var->contiguous = NC_TRUE;
2746  else if((flag & HDF_CHUNK) != 0) {
2747  var->contiguous = NC_FALSE;
2748  if (!(var->chunksizes = malloc(var->ndims * sizeof(size_t))))
2749  return NC_ENOMEM;
2750  for (d = 0; d < var->ndims; d++) {
2751  var->chunksizes[d] = chunkdefs.chunk_lengths[d];
2752  }
2753  }
2754  }
2755  }
2756 
2757  } /* next var */
2758 
2759 #ifdef LOGGING
2760  /* This will print out the names, types, lens, etc of the vars and
2761  atts in the file, if the logging level is 2 or greater. */
2762  log_metadata_nc(h5->root_grp->nc4_info->controller);
2763 #endif
2764  return NC_NOERR;
2765  return NC_ENOTBUILT;
2766 }
2767 #endif /* USE_HDF4 */
2768 
2769 int
2770 NC4_open(const char *path, int mode, int basepe, size_t *chunksizehintp,
2771  int use_parallel, void *parameters, NC_Dispatch *dispatch, NC *nc_file)
2772 {
2773  int res;
2774  int hdf_file = 0;
2775 #ifdef USE_PARALLEL4
2776  NC_MPI_INFO mpidfalt = {MPI_COMM_WORLD, MPI_INFO_NULL};
2777 #endif
2778 #if defined USE_PARALLEL4 || defined USE_HDF4
2779  int inmemory = ((mode & NC_INMEMORY) == NC_INMEMORY);
2780 #endif
2781 
2782  assert(nc_file && path);
2783 
2784  LOG((1, "%s: path %s mode %d params %x",
2785  __func__, path, mode, parameters));
2786 
2787 #ifdef USE_PARALLEL4
2788  if (!inmemory && use_parallel && parameters == NULL)
2789  parameters = &mpidfalt;
2790 #endif /* USE_PARALLEL4 */
2791 
2792  /* If this is our first file, turn off HDF5 error messages. */
2793  if (virgin)
2794  {
2795  if (H5Eset_auto(NULL, NULL) < 0)
2796  LOG((0, "Couldn't turn off HDF5 error messages!"));
2797  LOG((1, "HDF5 error messages turned off!"));
2798  virgin = 0;
2799  }
2800 
2801 
2802  /* Check the mode for validity */
2803  if((mode & ILLEGAL_OPEN_FLAGS) != 0)
2804  return NC_EINVAL;
2805 
2806  /* Cannot have both */
2807  if((mode & (NC_MPIIO|NC_MPIPOSIX)) == (NC_MPIIO|NC_MPIPOSIX))
2808  return NC_EINVAL;
2809 
2810 #ifndef USE_PARALLEL_POSIX
2811 /* If the HDF5 library has been compiled without the MPI-POSIX VFD, alias
2812  * the NC_MPIPOSIX flag to NC_MPIIO. -QAK
2813  */
2814  if(mode & NC_MPIPOSIX)
2815  {
2816  mode &= ~NC_MPIPOSIX;
2817  mode |= NC_MPIIO;
2818  }
2819 #endif /* USE_PARALLEL_POSIX */
2820 
2821  /* Figure out if this is a hdf4 or hdf5 file. */
2822  if ((res = nc_check_for_hdf(path, use_parallel, parameters, &hdf_file)))
2823  return res;
2824 
2825  /* Depending on the type of file, open it. */
2826  nc_file->int_ncid = nc_file->ext_ncid;
2827  if (hdf_file == NC_HDF5_FILE)
2828  res = nc4_open_file(path, mode, parameters, nc_file);
2829 #ifdef USE_HDF4
2830  else if (hdf_file == NC_HDF4_FILE && inmemory)
2831  return NC_EDISKLESS;
2832  else if (hdf_file == NC_HDF4_FILE)
2833  res = nc4_open_hdf4_file(path, mode, nc_file);
2834 #endif /* USE_HDF4 */
2835  else
2836  assert(0); /* should never happen */
2837 
2838  return res;
2839 }
2840 
2841 /* Unfortunately HDF only allows specification of fill value only when
2842  a dataset is created. Whereas in netcdf, you first create the
2843  variable and then (optionally) specify the fill value. To
2844  accomplish this in HDF5 I have to delete the dataset, and recreate
2845  it, with the fill value specified. */
2846 /* QAK: This looks completely unused in the code. (?) */
2847 int
2848 NC4_set_fill(int ncid, int fillmode, int *old_modep)
2849 {
2850  NC *nc;
2851  NC_HDF5_FILE_INFO_T* nc4_info;
2852 
2853  LOG((2, "%s: ncid 0x%x fillmode %d", __func__, ncid, fillmode));
2854 
2855  if (!(nc = nc4_find_nc_file(ncid,&nc4_info)))
2856  return NC_EBADID;
2857  assert(nc4_info);
2858 
2859  /* Trying to set fill on a read-only file? You sicken me! */
2860  if (nc4_info->no_write)
2861  return NC_EPERM;
2862 
2863  /* Did you pass me some weird fillmode? */
2864  if (fillmode != NC_FILL && fillmode != NC_NOFILL)
2865  return NC_EINVAL;
2866 
2867  /* If the user wants to know, tell him what the old mode was. */
2868  if (old_modep)
2869  *old_modep = nc4_info->fill_mode;
2870 
2871  nc4_info->fill_mode = fillmode;
2872 
2873 #if 0 /*def USE_PNETCDF*/
2874  /* Take care of files created/opened with parallel-netcdf library. */
2875  if (nc4_info->pnetcdf_file)
2876  return ncmpi_set_fill(nc->int_ncid, fillmode, old_modep);
2877 #endif /* USE_PNETCDF */
2878 
2879 
2880  return NC_NOERR;
2881 }
2882 
2883 /* Put the file back in redef mode. This is done automatically for
2884  * netcdf-4 files, if the user forgets. */
2885 int
2886 NC4_redef(int ncid)
2887 {
2888  NC_HDF5_FILE_INFO_T* nc4_info;
2889 
2890  LOG((1, "%s: ncid 0x%x", __func__, ncid));
2891 
2892  /* Find this file's metadata. */
2893  if (!(nc4_find_nc_file(ncid,&nc4_info)))
2894  return NC_EBADID;
2895  assert(nc4_info);
2896 
2897 #if 0 /*def USE_PNETCDF*/
2898  /* Take care of files created/opened with parallel-netcdf library. */
2899  if (nc4_info->pnetcdf_file)
2900  return ncmpi_redef(nc->int_ncid);
2901 #endif /* USE_PNETCDF */
2902 
2903  /* If we're already in define mode, return an error. */
2904  if (nc4_info->flags & NC_INDEF)
2905  return NC_EINDEFINE;
2906 
2907  /* If the file is read-only, return an error. */
2908  if (nc4_info->no_write)
2909  return NC_EPERM;
2910 
2911  /* Set define mode. */
2912  nc4_info->flags |= NC_INDEF;
2913 
2914  /* For nc_abort, we need to remember if we're in define mode as a
2915  redef. */
2916  nc4_info->redef = NC_TRUE;
2917 
2918  return NC_NOERR;
2919 }
2920 
2921 /* For netcdf-4 files, this just calls nc_enddef, ignoring the extra
2922  * parameters. */
2923 int
2924 NC4__enddef(int ncid, size_t h_minfree, size_t v_align,
2925  size_t v_minfree, size_t r_align)
2926 {
2927  if (nc4_find_nc_file(ncid,NULL) == NULL)
2928  return NC_EBADID;
2929 
2930  return NC4_enddef(ncid);
2931 }
2932 
2933 /* Take the file out of define mode. This is called automatically for
2934  * netcdf-4 files, if the user forgets. */
2935 static int NC4_enddef(int ncid)
2936 {
2937  NC *nc;
2938  NC_HDF5_FILE_INFO_T* nc4_info;
2939 
2940  LOG((1, "%s: ncid 0x%x", __func__, ncid));
2941 
2942  if (!(nc = nc4_find_nc_file(ncid,&nc4_info)))
2943  return NC_EBADID;
2944  assert(nc4_info);
2945 
2946 #if 0 /*def USE_PNETCDF*/
2947  if (nc4_info->pnetcdf_file)
2948  {
2949  int res;
2950  res = ncmpi_enddef(nc->int_ncid);
2951  if (!res)
2952  {
2953  if (nc4_info->pnetcdf_access_mode == NC_INDEPENDENT)
2954  res = ncmpi_begin_indep_data(nc->int_ncid);
2955  }
2956  return res;
2957  }
2958 #endif /* USE_PNETCDF */
2959 
2960  return nc4_enddef_netcdf4_file(nc4_info);
2961 }
2962 
2963 /* This function will write all changed metadata, and (someday) reread
2964  * all metadata from the file. */
2965 static int
2966 sync_netcdf4_file(NC_HDF5_FILE_INFO_T *h5)
2967 {
2968  int retval;
2969 
2970  assert(h5);
2971  LOG((3, "%s", __func__));
2972 
2973  /* If we're in define mode, that's an error, for strict nc3 rules,
2974  * otherwise, end define mode. */
2975  if (h5->flags & NC_INDEF)
2976  {
2977  if (h5->cmode & NC_CLASSIC_MODEL)
2978  return NC_EINDEFINE;
2979 
2980  /* Turn define mode off. */
2981  h5->flags ^= NC_INDEF;
2982 
2983  /* Redef mode needs to be tracked separately for nc_abort. */
2984  h5->redef = NC_FALSE;
2985  }
2986 
2987 #ifdef LOGGING
2988  /* This will print out the names, types, lens, etc of the vars and
2989  atts in the file, if the logging level is 2 or greater. */
2990  log_metadata_nc(h5->root_grp->nc4_info->controller);
2991 #endif
2992 
2993  /* Write any metadata that has changed. */
2994  if (!(h5->cmode & NC_NOWRITE))
2995  {
2996  nc_bool_t bad_coord_order = NC_FALSE; /* if detected, propagate to all groups to consistently store dimids */
2997 
2998  if ((retval = nc4_rec_write_groups_types(h5->root_grp)))
2999  return retval;
3000  if ((retval = nc4_rec_detect_need_to_preserve_dimids(h5->root_grp, &bad_coord_order)))
3001  return retval;
3002  if ((retval = nc4_rec_write_metadata(h5->root_grp, bad_coord_order)))
3003  return retval;
3004  }
3005 
3006  if (H5Fflush(h5->hdfid, H5F_SCOPE_GLOBAL) < 0)
3007  return NC_EHDFERR;
3008 
3009  return retval;
3010 }
3011 
3012 /* Flushes all buffers associated with the file, after writing all
3013  changed metadata. This may only be called in data mode. */
3014 int
3015 NC4_sync(int ncid)
3016 {
3017  NC *nc;
3018  int retval;
3019  NC_HDF5_FILE_INFO_T* nc4_info;
3020 
3021  LOG((2, "%s: ncid 0x%x", __func__, ncid));
3022 
3023  if (!(nc = nc4_find_nc_file(ncid,&nc4_info)))
3024  return NC_EBADID;
3025  assert(nc4_info);
3026 
3027 #if 0 /*def USE_PNETCDF*/
3028  /* Take care of files created/opened with parallel-netcdf library. */
3029  if (nc4_info->pnetcdf_file)
3030  return ncmpi_sync(nc->int_ncid);
3031 #endif /* USE_PNETCDF */
3032 
3033  /* If we're in define mode, we can't sync. */
3034  if (nc4_info && nc4_info->flags & NC_INDEF)
3035  {
3036  if (nc4_info->cmode & NC_CLASSIC_MODEL)
3037  return NC_EINDEFINE;
3038  if ((retval = NC4_enddef(ncid)))
3039  return retval;
3040  }
3041 
3042  return sync_netcdf4_file(nc4_info);
3043 }
3044 
3045 /* This function will free all allocated metadata memory, and close
3046  the HDF5 file. The group that is passed in must be the root group
3047  of the file. */
3048 static int
3049 close_netcdf4_file(NC_HDF5_FILE_INFO_T *h5, int abort)
3050 {
3051  int retval = NC_NOERR;
3052 
3053  assert(h5 && h5->root_grp);
3054  LOG((3, "%s: h5->path %s abort %d", __func__, h5->controller->path, abort));
3055 
3056  /* According to the docs, always end define mode on close. */
3057  if (h5->flags & NC_INDEF)
3058  h5->flags ^= NC_INDEF;
3059 
3060  /* Sync the file, unless we're aborting, or this is a read-only
3061  * file. */
3062  if (!h5->no_write && !abort)
3063  if ((retval = sync_netcdf4_file(h5)))
3064  goto exit;
3065 
3066  /* Delete all the list contents for vars, dims, and atts, in each
3067  * group. */
3068  if ((retval = nc4_rec_grp_del(&h5->root_grp, h5->root_grp)))
3069  goto exit;
3070 
3071  /* Close hdf file. */
3072 #ifdef USE_HDF4
3073  if (h5->hdf4)
3074  {
3075  if (SDend(h5->sdid))
3076  BAIL_QUIET(NC_EHDFERR);
3077  }
3078  else
3079 #endif /* USE_HDF4 */
3080  {
3081 #ifdef USE_PARALLEL4
3082  /* Free the MPI Comm & Info objects, if we opened the file in parallel */
3083  if(h5->parallel)
3084  {
3085  MPI_Comm_free(&h5->comm);
3086  if(MPI_INFO_NULL != h5->info)
3087  MPI_Info_free(&h5->info);
3088  }
3089 #endif
3090  if (H5Fclose(h5->hdfid) < 0)
3091  {
3092  int nobjs;
3093 
3094  nobjs = H5Fget_obj_count(h5->hdfid, H5F_OBJ_ALL);
3095  /* Apparently we can get an error even when nobjs == 0 */
3096  if(nobjs < 0) {
3097  BAIL_QUIET(NC_EHDFERR);
3098  } else if(nobjs > 0) {
3099 #ifdef LOGGING
3100  /* If the close doesn't work, probably there are still some HDF5
3101  * objects open, which means there's a bug in the library. So
3102  * print out some info on to help the poor programmer figure it
3103  * out. */
3104  LOG((0, "There are %d HDF5 objects open!", nobjs));
3105 #endif
3106  BAIL_QUIET(NC_EHDFERR);
3107  }
3108  }
3109  }
3110 
3111 exit:
3112  /* Free the nc4_info struct; above code should have reclaimed
3113  everything else */
3114  if(h5 != NULL)
3115  free(h5);
3116 
3117  return retval;
3118 }
3119 
3120 /* From the netcdf-3 docs: The function nc_abort just closes the
3121  netCDF dataset, if not in define mode. If the dataset is being
3122  created and is still in define mode, the dataset is deleted. If
3123  define mode was entered by a call to nc_redef, the netCDF dataset
3124  is restored to its state before definition mode was entered and the
3125  dataset is closed. */
3126 int
3127 NC4_abort(int ncid)
3128 {
3129  NC *nc;
3130  int delete_file = 0;
3131  char path[NC_MAX_NAME + 1];
3132  int retval = NC_NOERR;
3133  NC_HDF5_FILE_INFO_T* nc4_info;
3134 
3135  LOG((2, "%s: ncid 0x%x", __func__, ncid));
3136 
3137  /* Find metadata for this file. */
3138  if (!(nc = nc4_find_nc_file(ncid,&nc4_info)))
3139  return NC_EBADID;
3140 
3141  assert(nc4_info);
3142 
3143 #if 0 /*def USE_PNETCDF*/
3144  /* Take care of files created/opened with parallel-netcdf library. */
3145  if (nc4_info->pnetcdf_file)
3146  return ncmpi_abort(nc->int_ncid);
3147 #endif /* USE_PNETCDF */
3148 
3149  /* If we're in define mode, but not redefing the file, delete it. */
3150  if (nc4_info->flags & NC_INDEF && !nc4_info->redef)
3151  {
3152  delete_file++;
3153  strncpy(path, nc->path,NC_MAX_NAME);
3154  }
3155 
3156  /* Free any resources the netcdf-4 library has for this file's
3157  * metadata. */
3158  if ((retval = close_netcdf4_file(nc4_info, 1)))
3159  return retval;
3160 
3161  /* Delete the file, if we should. */
3162  if (delete_file)
3163  if (remove(path) < 0)
3164  return NC_ECANTREMOVE;
3165 
3166  return retval;
3167 }
3168 
3169 /* Close the netcdf file, writing any changes first. */
3170 int
3171 NC4_close(int ncid)
3172 {
3173  NC_GRP_INFO_T *grp;
3174  NC *nc;
3175  NC_HDF5_FILE_INFO_T *h5;
3176  int retval;
3177 
3178  LOG((1, "%s: ncid 0x%x", __func__, ncid));
3179 
3180  /* Find our metadata for this file. */
3181  if ((retval = nc4_find_nc_grp_h5(ncid, &nc, &grp, &h5)))
3182  return retval;
3183 
3184  assert(nc && h5 && grp);
3185 
3186  /* This must be the root group. */
3187  if (grp->parent)
3188  return NC_EBADGRPID;
3189 
3190 #if 0 /*def USE_PNETCDF*/
3191  /* Take care of files created/opened with parallel-netcdf library. */
3192  if (h5->pnetcdf_file)
3193  return ncmpi_close(nc->int_ncid);
3194 #endif /* USE_PNETCDF */
3195 
3196  /* Call the nc4 close. */
3197  if ((retval = close_netcdf4_file(grp->nc4_info, 0)))
3198  return retval;
3199 
3200  return NC_NOERR;
3201 }
3202 
3203 /* It's possible for any of these pointers to be NULL, in which case
3204  don't try to figure out that value. */
3205 int
3206 NC4_inq(int ncid, int *ndimsp, int *nvarsp, int *nattsp, int *unlimdimidp)
3207 {
3208  NC *nc;
3209  NC_HDF5_FILE_INFO_T *h5;
3210  NC_GRP_INFO_T *grp;
3211  NC_DIM_INFO_T *dim;
3212  NC_ATT_INFO_T *att;
3213  NC_VAR_INFO_T *var;
3214  int retval;
3215 
3216  LOG((2, "%s: ncid 0x%x", __func__, ncid));
3217 
3218  /* Find file metadata. */
3219  if ((retval = nc4_find_nc_grp_h5(ncid, &nc, &grp, &h5)))
3220  return retval;
3221 
3222  assert(h5 && grp && nc);
3223 
3224 #if 0 /*def USE_PNETCDF*/
3225  /* Take care of files created/opened with parallel-netcdf library. */
3226  if (h5->pnetcdf_file)
3227  return ncmpi_inq(nc->int_ncid, ndimsp, nvarsp, nattsp, unlimdimidp);
3228 #endif /* USE_PNETCDF */
3229 
3230  /* Count the number of dims, vars, and global atts. */
3231  if (ndimsp)
3232  {
3233  *ndimsp = 0;
3234  for (dim = grp->dim; dim; dim = dim->l.next)
3235  (*ndimsp)++;
3236  }
3237  if (nvarsp)
3238  {
3239  *nvarsp = 0;
3240  for (var = grp->var; var; var= var->l.next)
3241  (*nvarsp)++;
3242  }
3243  if (nattsp)
3244  {
3245  *nattsp = 0;
3246  for (att = grp->att; att; att = att->l.next)
3247  (*nattsp)++;
3248  }
3249 
3250  if (unlimdimidp)
3251  {
3252  /* Default, no unlimited dimension */
3253  *unlimdimidp = -1;
3254 
3255  /* If there's more than one unlimited dim, which was not possible
3256  with netcdf-3, then only the last unlimited one will be reported
3257  back in xtendimp. */
3258  /* Note that this code is inconsistent with nc_inq_unlimid() */
3259  for (dim = grp->dim; dim; dim = dim->l.next)
3260  if (dim->unlimited)
3261  {
3262  *unlimdimidp = dim->dimid;
3263  break;
3264  }
3265  }
3266 
3267  return NC_NOERR;
3268 }
3269 
3270 #if 0
3271 int
3272 NC4_set_content(int ncid, size_t size, void* memory)
3273 {
3274  int retval = NC_NOERR;
3275  herr_t herr;
3276  NC *nc;
3277  NC_HDF5_FILE_INFO_T *h5;
3278  NC_GRP_INFO_T *grp;
3279 
3280  LOG((4,"%s: ncid 0x%x size %ld memory 0x%x", __func__, ncid, size, memory));
3281 
3282  /* Find file metadata. */
3283  if ((retval = nc4_find_nc_grp_h5(ncid, &nc, &grp, &h5)))
3284  return retval;
3285  assert(h5 && grp && nc);
3286 
3287 #ifdef USE_DISKLESS
3288  herr = H5Pset_file_image(h5->hdfid,memory,size);
3289  if(herr)
3290  BAIL(NC_EHDFERR);
3291 #else
3292  retval = NC_EDISKLESS;
3293 #endif
3294 
3295 done:
3296  return retval;
3297 }
3298 #endif
3299 
3300 /* This function will do the enddef stuff for a netcdf-4 file. */
3301 int
3302 nc4_enddef_netcdf4_file(NC_HDF5_FILE_INFO_T *h5)
3303 {
3304  assert(h5);
3305  LOG((3, "%s", __func__));
3306 
3307  /* If we're not in define mode, return an error. */
3308  if (!(h5->flags & NC_INDEF))
3309  return NC_ENOTINDEFINE;
3310 
3311  /* Turn define mode off. */
3312  h5->flags ^= NC_INDEF;
3313 
3314  /* Redef mode needs to be tracked separately for nc_abort. */
3315  h5->redef = NC_FALSE;
3316 
3317  return sync_netcdf4_file(h5);
3318 }
3319 
3320 #ifdef EXTRA_TESTS
3321 int
3322 nc_exit()
3323 {
3324  if (num_plists || num_spaces)
3325  return NC_EHDFERR;
3326 
3327  return NC_NOERR;
3328 }
3329 #endif /* EXTRA_TESTS */
3330 
3331 #ifdef USE_PARALLEL4
3332 int
3333 nc_use_parallel_enabled()
3334 {
3335  return 0;
3336 }
3337 #endif /* USE_PARALLEL4 */
#define NC_PNETCDF
Use parallel-netcdf library; alias for NC_MPIIO.
Definition: netcdf.h:165
#define NC_ENOMEM
Memory allocation (malloc) failure.
Definition: netcdf.h:395
#define NC_CHAR
ISO/ASCII character.
Definition: netcdf.h:39
#define NC_ECANTWRITE
Can&#39;t write.
Definition: netcdf.h:428
int NC4_create(const char *path, int cmode, size_t initialsz, int basepe, size_t *chunksizehintp, int use_parallel, void *parameters, NC_Dispatch *dispatch, NC *nc_file)
Create a netCDF-4/HDF5 file.
Definition: nc4file.c:472
#define NC_UBYTE
unsigned 1 byte int
Definition: netcdf.h:45
#define NC_CLASSIC_MODEL
Enforce classic model on netCDF-4.
Definition: netcdf.h:141
#define NC_MAX_VAR_DIMS
max per variable dimensions
Definition: netcdf.h:269
#define NC_UINT
unsigned 4-byte int
Definition: netcdf.h:47
#define NC_NOCLOBBER
Don&#39;t destroy existing file.
Definition: netcdf.h:133
#define NC_INMEMORY
Read from memory.
Definition: netcdf.h:163
#define NC_EHDFERR
Error at HDF5 layer.
Definition: netcdf.h:426
#define NC_OPAQUE
opaque types
Definition: netcdf.h:57
#define NC_MPIIO
Turn on MPI I/O.
Definition: netcdf.h:158
#define NC_INT64
signed 8-byte int
Definition: netcdf.h:48
#define NC_STRING
string
Definition: netcdf.h:50
#define NC_ENOTINDEFINE
Operation not allowed in data mode.
Definition: netcdf.h:331
#define NC_DOUBLE
double precision floating point number
Definition: netcdf.h:44
#define NC_EBADCLASS
Bad class.
Definition: netcdf.h:447
int nc_type
The nc_type type is just an int.
Definition: netcdf.h:28
#define NC_64BIT_OFFSET
Use large (64-bit) file offsets.
Definition: netcdf.h:142
#define NC_NOWRITE
Set read-only access for nc_open().
Definition: netcdf.h:130
#define NC_BYTE
signed 1 byte integer
Definition: netcdf.h:38
#define NC_EINDEFINE
Operation not allowed in define mode.
Definition: netcdf.h:340
#define NC_ENOTNC
Not a netcdf file.
Definition: netcdf.h:371
#define NC_FORMAT_CDF5
Format specifier for nc_set_default_format() and returned by nc_inq_format.
Definition: netcdf.h:187
#define NC_ENOTBUILT
Attempt to use feature that was not turned on when netCDF was built.
Definition: netcdf.h:455
#define NC_ENDIAN_LITTLE
In HDF5 files you can set the endianness of variables with nc_def_var_endian().
Definition: netcdf.h:279
#define NC_EATTMETA
Problem with attribute metadata.
Definition: netcdf.h:432
#define NC_VLEN
vlen (variable-length) types
Definition: netcdf.h:56
#define NC_EMPI
MPI operation failed.
Definition: netcdf.h:458
#define NC_EDISKLESS
Error in using diskless access.
Definition: netcdf.h:456
#define NC_EFILEMETA
Problem with file metadata.
Definition: netcdf.h:430
This is the type of arrays of vlens.
Definition: netcdf.h:666
#define NC_EBADTYPE
Not a netcdf data type.
Definition: netcdf.h:357
#define NC_EBADNAME
Attribute or variable name contains illegal characters.
Definition: netcdf.h:387
#define NC_EDIMMETA
Problem with dimension metadata.
Definition: netcdf.h:431
#define NC_EINVAL
Invalid Argument.
Definition: netcdf.h:325
#define NC_INT
signed 4 byte integer
Definition: netcdf.h:41
#define NC_EBADGRPID
Bad group ID.
Definition: netcdf.h:443
#define NC_NOFILL
Argument to nc_set_fill() to turn off filling of data.
Definition: netcdf.h:120
#define NC_ENDIAN_BIG
In HDF5 files you can set the endianness of variables with nc_def_var_endian().
Definition: netcdf.h:280
#define NC_MAX_NAME
Maximum for classic library.
Definition: netcdf.h:268
#define NC_ECANTREMOVE
Can&#39;t remove file.
Definition: netcdf.h:418
#define NC_NAT
Not A Type.
Definition: netcdf.h:37
#define NC_EBADTYPID
Bad type ID.
Definition: netcdf.h:444
User data struct for call to H5Literate() in nc4_rec_read_metadata()
Definition: nc4file.c:72
#define NC_USHORT
unsigned 2-byte int
Definition: netcdf.h:46
#define NC_EPARINIT
Error initializing for parallel access.
Definition: netcdf.h:442
#define NC_NETCDF4
Use netCDF-4/HDF5 format.
Definition: netcdf.h:154
#define NC_EEXIST
netcdf file exists && NC_NOCLOBBER
Definition: netcdf.h:324
#define NC_FORMAT_NETCDF4_CLASSIC
Format specifier for nc_set_default_format() and returned by nc_inq_format.
Definition: netcdf.h:183
#define NC_EBADID
Not a netcdf id.
Definition: netcdf.h:322
#define NC_EVARMETA
Problem with variable metadata.
Definition: netcdf.h:433
struct NC4_rec_read_metadata_ud NC4_rec_read_metadata_ud_t
User data struct for call to H5Literate() in nc4_rec_read_metadata()
#define NC_MAX_UINT
Max or min values for a type.
Definition: netcdf.h:104
#define NC_SHORT
signed 2 byte integer
Definition: netcdf.h:40
#define NC_CDF5
Alias NC_CDF5 to NC_64BIT_DATA.
Definition: netcdf.h:139
#define NC_WRITE
Set read-write access for nc_open().
Definition: netcdf.h:131
#define NC_EMAXNAME
NC_MAX_NAME exceeded.
Definition: netcdf.h:373
#define NC_EPERM
Write to read only.
Definition: netcdf.h:326
#define NC_MAX_HDF4_NAME
This is the max size of an SD dataset name in HDF4 (from HDF4 documentation).
Definition: netcdf.h:273
#define NC_NOERR
No Error.
Definition: netcdf.h:315
#define NC_ENUM
enum types
Definition: netcdf.h:58
#define NC_DISKLESS
Use diskless file.
Definition: netcdf.h:135
Struct to track information about objects in a group, for nc4_rec_read_metadata() ...
Definition: nc4file.c:58
struct NC4_rec_read_metadata_obj_info NC4_rec_read_metadata_obj_info_t
Struct to track information about objects in a group, for nc4_rec_read_metadata() ...
#define NC_COMPOUND
compound types
Definition: netcdf.h:59
#define NC_FORMAT_64BIT_OFFSET
Format specifier for nc_set_default_format() and returned by nc_inq_format.
Definition: netcdf.h:180
#define NC_MMAP
Use diskless file with mmap.
Definition: netcdf.h:136
#define NC_FILL
Argument to nc_set_fill() to clear NC_NOFILL.
Definition: netcdf.h:119
#define NC_FLOAT
single precision floating point number
Definition: netcdf.h:43
#define NC_UINT64
unsigned 8-byte int
Definition: netcdf.h:49
#define NC_MPIPOSIX
Turn on MPI POSIX I/O.
Definition: netcdf.h:161

Return to the Main Unidata NetCDF page.
Generated on Sun Mar 27 2016 11:43:30 for NetCDF. NetCDF is a Unidata library.