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

Return to the Main Unidata NetCDF page.
Generated on Sat Dec 30 2017 10:59:33 for NetCDF. NetCDF is a Unidata library.