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