NetCDF  4.6.1
 All Data Structures Files Functions Variables Typedefs Macros Modules Pages
nc4hdf.c
Go to the documentation of this file.
1 
18 #include "config.h"
19 #include "nc4internal.h"
20 #include "nc4dispatch.h"
21 #include <H5DSpublic.h>
22 #include <math.h>
23 
24 #ifdef HAVE_INTTYPES_H
25 #define __STDC_FORMAT_MACROS
26 #include <inttypes.h>
27 #endif
28 
29 #ifdef USE_PARALLEL
30 #include "netcdf_par.h"
31 #endif
32 
33 #define NC3_STRICT_ATT_NAME "_nc3_strict"
35 #define NC_HDF5_MAX_NAME 1024
37 #define MAXNAME 1024
40 static unsigned int OTYPES[5] = {H5F_OBJ_FILE, H5F_OBJ_DATASET, H5F_OBJ_GROUP,
41  H5F_OBJ_DATATYPE, H5F_OBJ_ATTR};
42 
50 static int
51 flag_atts_dirty(NC_ATT_INFO_T **attlist) {
52 
53  NC_ATT_INFO_T *att = NULL;
54 
55  if(attlist == NULL) {
56  return NC_NOERR;
57  }
58 
59  for(att = *attlist; att; att = att->l.next) {
60  att->dirty = NC_TRUE;
61  }
62 
63  return NC_NOERR;
64 
65 }
66 
83 int
84 rec_reattach_scales(NC_GRP_INFO_T *grp, int dimid, hid_t dimscaleid)
85 {
86  NC_VAR_INFO_T *var;
87  NC_GRP_INFO_T *child_grp;
88  int d, i;
89  int retval;
90 
91  assert(grp && grp->name && dimid >= 0 && dimscaleid >= 0);
92  LOG((3, "%s: grp->name %s", __func__, grp->name));
93 
94  /* If there are any child groups, attach dimscale there, if needed. */
95  for (child_grp = grp->children; child_grp; child_grp = child_grp->l.next)
96  if ((retval = rec_reattach_scales(child_grp, dimid, dimscaleid)))
97  return retval;
98 
99  /* Find any vars that use this dimension id. */
100  for (i=0; i < grp->vars.nelems; i++)
101  {
102  var = grp->vars.value[i];
103  if (!var) continue;
104  for (d = 0; d < var->ndims; d++)
105  if (var->dimids[d] == dimid && !var->dimscale)
106  {
107  LOG((2, "%s: attaching scale for dimid %d to var %s",
108  __func__, var->dimids[d], var->name));
109  if (var->created)
110  {
111  if (H5DSattach_scale(var->hdf_datasetid, dimscaleid, d) < 0)
112  return NC_EHDFERR;
113  var->dimscale_attached[d] = NC_TRUE;
114  }
115  }
116  }
117  return NC_NOERR;
118 }
119 
136 int
137 rec_detach_scales(NC_GRP_INFO_T *grp, int dimid, hid_t dimscaleid)
138 {
139  NC_VAR_INFO_T *var;
140  NC_GRP_INFO_T *child_grp;
141  int d, i;
142  int retval;
143 
144  assert(grp && grp->name && dimid >= 0 && dimscaleid >= 0);
145  LOG((3, "%s: grp->name %s", __func__, grp->name));
146 
147  /* If there are any child groups, detach dimscale there, if needed. */
148  for (child_grp = grp->children; child_grp; child_grp = child_grp->l.next)
149  if ((retval = rec_detach_scales(child_grp, dimid, dimscaleid)))
150  return retval;
151 
152  /* Find any vars that use this dimension id. */
153  for (i=0; i < grp->vars.nelems; i++)
154  {
155  var = grp->vars.value[i];
156  if (!var) continue;
157  for (d = 0; d < var->ndims; d++)
158  if (var->dimids[d] == dimid && !var->dimscale)
159  {
160  LOG((2, "%s: detaching scale for dimid %d to var %s",
161  __func__, var->dimids[d], var->name));
162  if (var->created)
163  if (var->dimscale_attached && var->dimscale_attached[d])
164  {
165  if (H5DSdetach_scale(var->hdf_datasetid, dimscaleid, d) < 0)
166  return NC_EHDFERR;
167  var->dimscale_attached[d] = NC_FALSE;
168  }
169  }
170  }
171  return NC_NOERR;
172 }
173 
185 int
186 nc4_open_var_grp2(NC_GRP_INFO_T *grp, int varid, hid_t *dataset)
187 {
188  NC_VAR_INFO_T *var;
189 
190  /* Find the requested varid. */
191  if (varid < 0 || varid >= grp->vars.nelems)
192  return NC_ENOTVAR;
193  var = grp->vars.value[varid];
194  if (!var) return NC_ENOTVAR;
195  assert(var->varid == varid);
196 
197  /* Open this dataset if necessary. */
198  if (!var->hdf_datasetid)
199  if ((var->hdf_datasetid = H5Dopen2(grp->hdf_grpid, var->name,
200  H5P_DEFAULT)) < 0)
201  return NC_ENOTVAR;
202 
203  *dataset = var->hdf_datasetid;
204 
205  return NC_NOERR;
206 }
207 
219 int
220 nc4_get_default_fill_value(const NC_TYPE_INFO_T *type_info, void *fill_value)
221 {
222  switch (type_info->nc_typeid)
223  {
224  case NC_CHAR:
225  *(char *)fill_value = NC_FILL_CHAR;
226  break;
227 
228  case NC_STRING:
229  *(char **)fill_value = strdup(NC_FILL_STRING);
230  break;
231 
232  case NC_BYTE:
233  *(signed char *)fill_value = NC_FILL_BYTE;
234  break;
235 
236  case NC_SHORT:
237  *(short *)fill_value = NC_FILL_SHORT;
238  break;
239 
240  case NC_INT:
241  *(int *)fill_value = NC_FILL_INT;
242  break;
243 
244  case NC_UBYTE:
245  *(unsigned char *)fill_value = NC_FILL_UBYTE;
246  break;
247 
248  case NC_USHORT:
249  *(unsigned short *)fill_value = NC_FILL_USHORT;
250  break;
251 
252  case NC_UINT:
253  *(unsigned int *)fill_value = NC_FILL_UINT;
254  break;
255 
256  case NC_INT64:
257  *(long long *)fill_value = NC_FILL_INT64;
258  break;
259 
260  case NC_UINT64:
261  *(unsigned long long *)fill_value = NC_FILL_UINT64;
262  break;
263 
264  case NC_FLOAT:
265  *(float *)fill_value = NC_FILL_FLOAT;
266  break;
267 
268  case NC_DOUBLE:
269  *(double *)fill_value = NC_FILL_DOUBLE;
270  break;
271 
272  default:
273  return NC_EINVAL;
274  }
275 
276  return NC_NOERR;
277 }
278 
290 static int
291 get_fill_value(NC_HDF5_FILE_INFO_T *h5, NC_VAR_INFO_T *var, void **fillp)
292 {
293  size_t size;
294  int retval;
295 
296  /* Find out how much space we need for this type's fill value. */
297  if (var->type_info->nc_type_class == NC_VLEN)
298  size = sizeof(nc_vlen_t);
299  else if (var->type_info->nc_type_class == NC_STRING)
300  size = sizeof(char *);
301  else
302  {
303  if ((retval = nc4_get_typelen_mem(h5, var->type_info->nc_typeid, 0, &size)))
304  return retval;
305  }
306  assert(size);
307 
308  /* Allocate the space. */
309  if (!((*fillp) = calloc(1, size)))
310  return NC_ENOMEM;
311 
312  /* If the user has set a fill_value for this var, use, otherwise
313  * find the default fill value. */
314  if (var->fill_value)
315  {
316  LOG((4, "Found a fill value for var %s", var->name));
317  if (var->type_info->nc_type_class == NC_VLEN)
318  {
319  nc_vlen_t *in_vlen = (nc_vlen_t *)(var->fill_value), *fv_vlen = (nc_vlen_t *)(*fillp);
320 
321  fv_vlen->len = in_vlen->len;
322  if (!(fv_vlen->p = malloc(size * in_vlen->len)))
323  {
324  free(*fillp);
325  *fillp = NULL;
326  return NC_ENOMEM;
327  }
328  memcpy(fv_vlen->p, in_vlen->p, in_vlen->len * size);
329  }
330  else if (var->type_info->nc_type_class == NC_STRING)
331  {
332  if (*(char **)var->fill_value)
333  if (!(**(char ***)fillp = strdup(*(char **)var->fill_value)))
334  {
335  free(*fillp);
336  *fillp = NULL;
337  return NC_ENOMEM;
338  }
339  }
340  else
341  memcpy((*fillp), var->fill_value, size);
342  }
343  else
344  {
345  if (nc4_get_default_fill_value(var->type_info, *fillp))
346  {
347  /* Note: release memory, but don't return error on failure */
348  free(*fillp);
349  *fillp = NULL;
350  }
351  }
352 
353  return NC_NOERR;
354 }
355 
373 int
374 nc4_get_hdf_typeid(NC_HDF5_FILE_INFO_T *h5, nc_type xtype,
375  hid_t *hdf_typeid, int endianness)
376 {
377  NC_TYPE_INFO_T *type;
378  hid_t typeid = 0;
379  int retval = NC_NOERR;
380 
381  assert(hdf_typeid && h5);
382 
383  *hdf_typeid = -1;
384 
385  /* Determine an appropriate HDF5 datatype */
386  if (xtype == NC_NAT)
387  /* NAT = 'Not A Type' (c.f. NaN) */
388  return NC_EBADTYPE;
389  else if (xtype == NC_CHAR || xtype == NC_STRING)
390  {
391  /* NC_CHAR & NC_STRING types create a new HDF5 datatype */
392  if (xtype == NC_CHAR)
393  {
394  if ((typeid = H5Tcopy(H5T_C_S1)) < 0)
395  return NC_EHDFERR;
396  if (H5Tset_strpad(typeid, H5T_STR_NULLTERM) < 0)
397  BAIL(NC_EVARMETA);
398  if(H5Tset_cset(typeid, H5T_CSET_ASCII) < 0)
399  BAIL(NC_EVARMETA);
400 
401  /* Take ownership of the newly created HDF5 datatype */
402  *hdf_typeid = typeid;
403  typeid = 0;
404  }
405  else
406  {
407  if ((typeid = H5Tcopy(H5T_C_S1)) < 0)
408  return NC_EHDFERR;
409  if (H5Tset_size(typeid, H5T_VARIABLE) < 0)
410  BAIL(NC_EVARMETA);
411  if(H5Tset_cset(typeid, H5T_CSET_UTF8) < 0)
412  BAIL(NC_EVARMETA);
413 
414  /* Take ownership of the newly created HDF5 datatype */
415  *hdf_typeid = typeid;
416  typeid = 0;
417  }
418  }
419  else
420  {
421  /* All other types use an existing HDF5 datatype */
422  switch (xtype)
423  {
424  case NC_BYTE: /* signed 1 byte integer */
425  if (endianness == NC_ENDIAN_LITTLE)
426  typeid = H5T_STD_I8LE;
427  else if (endianness == NC_ENDIAN_BIG)
428  typeid = H5T_STD_I8BE;
429  else
430  typeid = H5T_NATIVE_SCHAR;
431  break;
432 
433  case NC_SHORT: /* signed 2 byte integer */
434  if (endianness == NC_ENDIAN_LITTLE)
435  typeid = H5T_STD_I16LE;
436  else if (endianness == NC_ENDIAN_BIG)
437  typeid = H5T_STD_I16BE;
438  else
439  typeid = H5T_NATIVE_SHORT;
440  break;
441 
442  case NC_INT:
443  if (endianness == NC_ENDIAN_LITTLE)
444  typeid = H5T_STD_I32LE;
445  else if (endianness == NC_ENDIAN_BIG)
446  typeid = H5T_STD_I32BE;
447  else
448  typeid = H5T_NATIVE_INT;
449  break;
450 
451  case NC_UBYTE:
452  if (endianness == NC_ENDIAN_LITTLE)
453  typeid = H5T_STD_U8LE;
454  else if (endianness == NC_ENDIAN_BIG)
455  typeid = H5T_STD_U8BE;
456  else
457  typeid = H5T_NATIVE_UCHAR;
458  break;
459 
460  case NC_USHORT:
461  if (endianness == NC_ENDIAN_LITTLE)
462  typeid = H5T_STD_U16LE;
463  else if (endianness == NC_ENDIAN_BIG)
464  typeid = H5T_STD_U16BE;
465  else
466  typeid = H5T_NATIVE_USHORT;
467  break;
468 
469  case NC_UINT:
470  if (endianness == NC_ENDIAN_LITTLE)
471  typeid = H5T_STD_U32LE;
472  else if (endianness == NC_ENDIAN_BIG)
473  typeid = H5T_STD_U32BE;
474  else
475  typeid = H5T_NATIVE_UINT;
476  break;
477 
478  case NC_INT64:
479  if (endianness == NC_ENDIAN_LITTLE)
480  typeid = H5T_STD_I64LE;
481  else if (endianness == NC_ENDIAN_BIG)
482  typeid = H5T_STD_I64BE;
483  else
484  typeid = H5T_NATIVE_LLONG;
485  break;
486 
487  case NC_UINT64:
488  if (endianness == NC_ENDIAN_LITTLE)
489  typeid = H5T_STD_U64LE;
490  else if (endianness == NC_ENDIAN_BIG)
491  typeid = H5T_STD_U64BE;
492  else
493  typeid = H5T_NATIVE_ULLONG;
494  break;
495 
496  case NC_FLOAT:
497  if (endianness == NC_ENDIAN_LITTLE)
498  typeid = H5T_IEEE_F32LE;
499  else if (endianness == NC_ENDIAN_BIG)
500  typeid = H5T_IEEE_F32BE;
501  else
502  typeid = H5T_NATIVE_FLOAT;
503  break;
504 
505  case NC_DOUBLE:
506  if (endianness == NC_ENDIAN_LITTLE)
507  typeid = H5T_IEEE_F64LE;
508  else if (endianness == NC_ENDIAN_BIG)
509  typeid = H5T_IEEE_F64BE;
510  else
511  typeid = H5T_NATIVE_DOUBLE;
512  break;
513 
514  default:
515  /* Maybe this is a user defined type? */
516  if (nc4_find_type(h5, xtype, &type))
517  return NC_EBADTYPE;
518  if (!type)
519  return NC_EBADTYPE;
520  typeid = type->hdf_typeid;
521  break;
522  }
523  assert(typeid);
524 
525  /* Copy the HDF5 datatype, so the function operates uniformly */
526  if ((*hdf_typeid = H5Tcopy(typeid)) < 0)
527  return NC_EHDFERR;
528  typeid = 0;
529  }
530  assert(*hdf_typeid != -1);
531 
532 exit:
533  if (typeid > 0 && H5Tclose(typeid) < 0)
534  BAIL2(NC_EHDFERR);
535  return retval;
536 }
537 
550 static int
551 check_for_vara(nc_type *mem_nc_type, NC_VAR_INFO_T *var, NC_HDF5_FILE_INFO_T *h5)
552 {
553  int retval;
554 
555  /* If mem_nc_type is NC_NAT, it means we want to use the file type
556  * as the mem type as well. */
557  assert(mem_nc_type);
558  if (*mem_nc_type == NC_NAT)
559  *mem_nc_type = var->type_info->nc_typeid;
560  assert(*mem_nc_type);
561 
562  /* No NC_CHAR conversions, you pervert! */
563  if (var->type_info->nc_typeid != *mem_nc_type &&
564  (var->type_info->nc_typeid == NC_CHAR || *mem_nc_type == NC_CHAR))
565  return NC_ECHAR;
566 
567  /* If we're in define mode, we can't read or write data. */
568  if (h5->flags & NC_INDEF)
569  {
570  if (h5->cmode & NC_CLASSIC_MODEL)
571  return NC_EINDEFINE;
572  if ((retval = nc4_enddef_netcdf4_file(h5)))
573  return retval;
574  }
575 
576  return NC_NOERR;
577 }
578 
579 #ifdef LOGGING
580 
583 static void
584 log_dim_info(NC_VAR_INFO_T *var, hsize_t *fdims, hsize_t *fmaxdims,
585  hsize_t *start, hsize_t *count)
586 {
587  int d2;
588 
589  /* Print some debugging info... */
590  LOG((4, "%s: var name %s ndims %d", __func__, var->name, var->ndims));
591  LOG((4, "File space, and requested:"));
592  for (d2 = 0; d2 < var->ndims; d2++)
593  {
594  LOG((4, "fdims[%d]=%Ld fmaxdims[%d]=%Ld", d2, fdims[d2], d2,
595  fmaxdims[d2]));
596  LOG((4, "start[%d]=%Ld count[%d]=%Ld", d2, start[d2], d2, count[d2]));
597  }
598 }
599 #endif /* LOGGING */
600 
601 #ifdef USE_PARALLEL4
602 
613 static int
614 set_par_access(NC_HDF5_FILE_INFO_T *h5, NC_VAR_INFO_T *var, hid_t xfer_plistid)
615 {
616  /* If netcdf is built with parallel I/O, then parallel access can
617  * be used, and, if this file was opened or created for parallel
618  * access, we need to set the transfer mode. */
619  if (h5->parallel)
620  {
621  H5FD_mpio_xfer_t hdf5_xfer_mode;
622 
623  /* Decide on collective or independent. */
624  hdf5_xfer_mode = (var->parallel_access != NC_INDEPENDENT) ?
625  H5FD_MPIO_COLLECTIVE : H5FD_MPIO_INDEPENDENT;
626 
627  /* Set the mode in the transfer property list. */
628  if (H5Pset_dxpl_mpio(xfer_plistid, hdf5_xfer_mode) < 0)
629  return NC_EPARINIT;
630 
631  LOG((4, "%s: %d H5FD_MPIO_COLLECTIVE: %d H5FD_MPIO_INDEPENDENT: %d",
632  __func__, (int)hdf5_xfer_mode, H5FD_MPIO_COLLECTIVE, H5FD_MPIO_INDEPENDENT));
633  }
634  return NC_NOERR;
635 }
636 #endif
637 
662 int
663 nc4_put_vara(NC *nc, int ncid, int varid, const size_t *startp,
664  const size_t *countp, nc_type mem_nc_type, int is_long, void *data)
665 {
666  NC_GRP_INFO_T *grp;
667  NC_HDF5_FILE_INFO_T *h5;
668  NC_VAR_INFO_T *var;
669  NC_DIM_INFO_T *dim;
670  hid_t file_spaceid = 0, mem_spaceid = 0, xfer_plistid = 0;
671  long long unsigned xtend_size[NC_MAX_VAR_DIMS];
672  hsize_t fdims[NC_MAX_VAR_DIMS], fmaxdims[NC_MAX_VAR_DIMS];
673  hsize_t start[NC_MAX_VAR_DIMS], count[NC_MAX_VAR_DIMS];
674  char *name_to_use;
675  int need_to_extend = 0;
676 #ifdef USE_PARALLEL4
677  int extend_possible = 0;
678 #endif
679  int retval = NC_NOERR, range_error = 0, i, d2;
680  void *bufr = NULL;
681 #ifndef HDF5_CONVERT
682  int need_to_convert = 0;
683  size_t len = 1;
684 #endif
685 #ifdef HDF5_CONVERT
686  hid_t mem_typeid = 0;
687 #endif
688 
689  /* Find our metadata for this file, group, and var. */
690  assert(nc);
691  if ((retval = nc4_find_g_var_nc(nc, ncid, varid, &grp, &var)))
692  return retval;
693  h5 = NC4_DATA(nc);
694  assert(grp && h5 && var && var->name);
695 
696  LOG((3, "%s: var->name %s mem_nc_type %d is_long %d",
697  __func__, var->name, mem_nc_type, is_long));
698 
699  /* Check some stuff about the type and the file. If the file must
700  * be switched from define mode, it happens here. */
701  if ((retval = check_for_vara(&mem_nc_type, var, h5)))
702  return retval;
703 
704  /* Convert from size_t and ptrdiff_t to hssize_t, and hsize_t. */
705  for (i = 0; i < var->ndims; i++)
706  {
707  start[i] = startp[i];
708  count[i] = countp[i];
709  }
710 
711  /* Open this dataset if necessary, also checking for a weird case:
712  * a non-coordinate (and non-scalar) variable that has the same
713  * name as a dimension. */
714  if (var->hdf5_name && strlen(var->hdf5_name) >= strlen(NON_COORD_PREPEND) &&
715  strncmp(var->hdf5_name, NON_COORD_PREPEND, strlen(NON_COORD_PREPEND)) == 0 &&
716  var->ndims)
717  name_to_use = var->hdf5_name;
718  else
719  name_to_use = var->name;
720  if (!var->hdf_datasetid)
721  if ((var->hdf_datasetid = H5Dopen2(grp->hdf_grpid, name_to_use, H5P_DEFAULT)) < 0)
722  return NC_ENOTVAR;
723 
724  /* Get file space of data. */
725  if ((file_spaceid = H5Dget_space(var->hdf_datasetid)) < 0)
726  BAIL(NC_EHDFERR);
727 
728  /* Check to ensure the user selection is
729  * valid. H5Sget_simple_extent_dims gets the sizes of all the dims
730  * and put them in fdims. */
731  if (H5Sget_simple_extent_dims(file_spaceid, fdims, fmaxdims) < 0)
732  BAIL(NC_EHDFERR);
733 
734 #ifdef LOGGING
735  log_dim_info(var, fdims, fmaxdims, start, count);
736 #endif
737 
738  /* Check dimension bounds. Remember that unlimited dimensions can
739  * put data beyond their current length. */
740  for (d2 = 0; d2 < var->ndims; d2++)
741  {
742  dim = var->dim[d2];
743  assert(dim && dim->dimid == var->dimids[d2]);
744  if (!dim->unlimited)
745  {
746 #ifdef RELAX_COORD_BOUND
747  if (start[d2] > (hssize_t)fdims[d2] ||
748  (start[d2] == (hssize_t)fdims[d2] && count[d2] > 0))
749 #else
750  if (start[d2] >= (hssize_t)fdims[d2])
751 #endif
752  BAIL_QUIET(NC_EINVALCOORDS);
753  if (start[d2] + count[d2] > fdims[d2])
754  BAIL_QUIET(NC_EEDGE);
755  }
756  }
757 
758  /* Now you would think that no one would be crazy enough to write
759  a scalar dataspace with one of the array function calls, but you
760  would be wrong. So let's check to see if the dataset is
761  scalar. If it is, we won't try to set up a hyperslab. */
762  if (H5Sget_simple_extent_type(file_spaceid) == H5S_SCALAR)
763  {
764  if ((mem_spaceid = H5Screate(H5S_SCALAR)) < 0)
765  BAIL(NC_EHDFERR);
766  }
767  else
768  {
769  if (H5Sselect_hyperslab(file_spaceid, H5S_SELECT_SET, start, NULL,
770  count, NULL) < 0)
771  BAIL(NC_EHDFERR);
772 
773  /* Create a space for the memory, just big enough to hold the slab
774  we want. */
775  if ((mem_spaceid = H5Screate_simple(var->ndims, count, NULL)) < 0)
776  BAIL(NC_EHDFERR);
777  }
778 
779 #ifndef HDF5_CONVERT
780  /* Are we going to convert any data? (No converting of compound or
781  * opaque types.) */
782  if ((mem_nc_type != var->type_info->nc_typeid || (var->type_info->nc_typeid == NC_INT && is_long)) &&
783  mem_nc_type != NC_COMPOUND && mem_nc_type != NC_OPAQUE)
784  {
785  size_t file_type_size;
786 
787  /* We must convert - allocate a buffer. */
788  need_to_convert++;
789  if (var->ndims)
790  for (d2=0; d2<var->ndims; d2++)
791  len *= countp[d2];
792  LOG((4, "converting data for var %s type=%d len=%d", var->name,
793  var->type_info->nc_typeid, len));
794 
795  /* Later on, we will need to know the size of this type in the
796  * file. */
797  assert(var->type_info->size);
798  file_type_size = var->type_info->size;
799 
800  /* If we're reading, we need bufr to have enough memory to store
801  * the data in the file. If we're writing, we need bufr to be
802  * big enough to hold all the data in the file's type. */
803  if(len > 0)
804  if (!(bufr = malloc(len * file_type_size)))
805  BAIL(NC_ENOMEM);
806  }
807  else
808 #endif /* ifndef HDF5_CONVERT */
809  bufr = data;
810 
811 #ifdef HDF5_CONVERT
812  /* Get the HDF type of the data in memory. */
813  if ((retval = nc4_get_hdf_typeid(h5, mem_nc_type, &mem_typeid,
814  var->type_info->endianness)))
815  BAIL(retval);
816 #endif
817 
818  /* Create the data transfer property list. */
819  if ((xfer_plistid = H5Pcreate(H5P_DATASET_XFER)) < 0)
820  BAIL(NC_EHDFERR);
821 
822  /* Apply the callback function which will detect range
823  * errors. Which one to call depends on the length of the
824  * destination buffer type. */
825 #ifdef HDF5_CONVERT
826  if (H5Pset_type_conv_cb(xfer_plistid, except_func, &range_error) < 0)
827  BAIL(NC_EHDFERR);
828 #endif
829 
830 #ifdef USE_PARALLEL4
831  /* Set up parallel I/O, if needed. */
832  if ((retval = set_par_access(h5, var, xfer_plistid)))
833  BAIL(retval);
834 #endif
835 
836  /* Read/write this hyperslab into memory. */
837  /* Does the dataset have to be extended? If it's already
838  extended to the required size, it will do no harm to reextend
839  it to that size. */
840  if (var->ndims)
841  {
842  for (d2 = 0; d2 < var->ndims; d2++)
843  {
844  dim = var->dim[d2];
845  assert(dim && dim->dimid == var->dimids[d2]);
846  if (dim->unlimited)
847  {
848 #ifdef USE_PARALLEL4
849  extend_possible = 1;
850 #endif
851  if (start[d2] + count[d2] > fdims[d2])
852  {
853  xtend_size[d2] = (long long unsigned)(start[d2] + count[d2]);
854  need_to_extend++;
855  }
856  else
857  xtend_size[d2] = (long long unsigned)fdims[d2];
858 
859  if (start[d2] + count[d2] > dim->len)
860  {
861  dim->len = start[d2] + count[d2];
862  dim->extended = NC_TRUE;
863  }
864  }
865  else
866  {
867  xtend_size[d2] = (long long unsigned)dim->len;
868  }
869  }
870 
871 #ifdef USE_PARALLEL4
872  /* Check if anyone wants to extend */
873  if (extend_possible && h5->parallel && NC_COLLECTIVE == var->parallel_access)
874  {
875  /* Form consensus opinion among all processes about whether to perform
876  * collective I/O
877  */
878  if(MPI_SUCCESS != MPI_Allreduce(MPI_IN_PLACE, &need_to_extend, 1, MPI_INT, MPI_BOR, h5->comm))
879  BAIL(NC_EMPI);
880  }
881 #endif /* USE_PARALLEL4 */
882 
883  /* If we need to extend it, we also need a new file_spaceid
884  to reflect the new size of the space. */
885  if (need_to_extend)
886  {
887  LOG((4, "extending dataset"));
888 #ifdef USE_PARALLEL4
889  if (h5->parallel)
890  {
891  if(NC_COLLECTIVE != var->parallel_access)
892  BAIL(NC_ECANTEXTEND);
893 
894  /* Reach consensus about dimension sizes to extend to */
895  if(MPI_SUCCESS != MPI_Allreduce(MPI_IN_PLACE, xtend_size, var->ndims, MPI_UNSIGNED_LONG_LONG, MPI_MAX, h5->comm))
896  BAIL(NC_EMPI);
897  }
898 #endif /* USE_PARALLEL4 */
899  /* Convert xtend_size back to hsize_t for use with H5Dset_extent */
900  for (d2 = 0; d2 < var->ndims; d2++)
901  fdims[d2] = (hsize_t)xtend_size[d2];
902 
903  if (H5Dset_extent(var->hdf_datasetid, fdims) < 0)
904  BAIL(NC_EHDFERR);
905  if (file_spaceid > 0 && H5Sclose(file_spaceid) < 0)
906  BAIL2(NC_EHDFERR);
907  if ((file_spaceid = H5Dget_space(var->hdf_datasetid)) < 0)
908  BAIL(NC_EHDFERR);
909  if (H5Sselect_hyperslab(file_spaceid, H5S_SELECT_SET,
910  start, NULL, count, NULL) < 0)
911  BAIL(NC_EHDFERR);
912  }
913  }
914 
915 #ifndef HDF5_CONVERT
916  /* Do we need to convert the data? */
917  if (need_to_convert)
918  {
919  if ((retval = nc4_convert_type(data, bufr, mem_nc_type, var->type_info->nc_typeid,
920  len, &range_error, var->fill_value,
921  (h5->cmode & NC_CLASSIC_MODEL), is_long, 0)))
922  BAIL(retval);
923  }
924 #endif
925 
926  /* Write the data. At last! */
927  LOG((4, "about to H5Dwrite datasetid 0x%x mem_spaceid 0x%x "
928  "file_spaceid 0x%x", var->hdf_datasetid, mem_spaceid, file_spaceid));
929  if (H5Dwrite(var->hdf_datasetid, var->type_info->hdf_typeid,
930  mem_spaceid, file_spaceid, xfer_plistid, bufr) < 0)
931  BAIL(NC_EHDFERR);
932 
933  /* Remember that we have written to this var so that Fill Value
934  * can't be set for it. */
935  if (!var->written_to)
936  var->written_to = NC_TRUE;
937 
938  /* For strict netcdf-3 rules, ignore erange errors between UBYTE
939  * and BYTE types. */
940  if ((h5->cmode & NC_CLASSIC_MODEL) &&
941  (var->type_info->nc_typeid == NC_UBYTE || var->type_info->nc_typeid == NC_BYTE) &&
942  (mem_nc_type == NC_UBYTE || mem_nc_type == NC_BYTE) &&
943  range_error)
944  range_error = 0;
945 
946 exit:
947 #ifdef HDF5_CONVERT
948  if (mem_typeid > 0 && H5Tclose(mem_typeid) < 0)
949  BAIL2(NC_EHDFERR);
950 #endif
951  if (file_spaceid > 0 && H5Sclose(file_spaceid) < 0)
952  BAIL2(NC_EHDFERR);
953  if (mem_spaceid > 0 && H5Sclose(mem_spaceid) < 0)
954  BAIL2(NC_EHDFERR);
955  if (xfer_plistid && (H5Pclose(xfer_plistid) < 0))
956  BAIL2(NC_EPARINIT);
957 #ifndef HDF5_CONVERT
958  if (need_to_convert && bufr) free(bufr);
959 #endif
960 
961  /* If there was an error return it, otherwise return any potential
962  range error value. If none, return NC_NOERR as usual.*/
963  if (retval)
964  return retval;
965  if (range_error)
966  return NC_ERANGE;
967  return NC_NOERR;
968 }
969 
995 int
996 nc4_get_vara(NC *nc, int ncid, int varid, const size_t *startp,
997  const size_t *countp, nc_type mem_nc_type, int is_long, void *data)
998 {
999  NC_GRP_INFO_T *grp;
1000  NC_HDF5_FILE_INFO_T *h5;
1001  NC_VAR_INFO_T *var;
1002  NC_DIM_INFO_T *dim;
1003  hid_t file_spaceid = 0, mem_spaceid = 0;
1004  hid_t xfer_plistid = 0;
1005  size_t file_type_size;
1006  hsize_t *xtend_size = NULL, count[NC_MAX_VAR_DIMS];
1007  hsize_t fdims[NC_MAX_VAR_DIMS], fmaxdims[NC_MAX_VAR_DIMS];
1008  hsize_t start[NC_MAX_VAR_DIMS];
1009  char *name_to_use;
1010  void *fillvalue = NULL;
1011  int no_read = 0, provide_fill = 0;
1012  int fill_value_size[NC_MAX_VAR_DIMS];
1013  int scalar = 0, retval = NC_NOERR, range_error = 0, i, d2;
1014  void *bufr = NULL;
1015 #ifdef HDF5_CONVERT
1016  hid_t mem_typeid = 0;
1017 #endif
1018 #ifndef HDF5_CONVERT
1019  int need_to_convert = 0;
1020  size_t len = 1;
1021 #endif
1022 
1023  /* Find our metadata for this file, group, and var. */
1024  assert(nc);
1025  if ((retval = nc4_find_g_var_nc(nc, ncid, varid, &grp, &var)))
1026  return retval;
1027  h5 = NC4_DATA(nc);
1028  assert(grp && h5 && var && var->name);
1029 
1030  LOG((3, "%s: var->name %s mem_nc_type %d is_long %d",
1031  __func__, var->name, mem_nc_type, is_long));
1032 
1033  /* Check some stuff about the type and the file. */
1034  if ((retval = check_for_vara(&mem_nc_type, var, h5)))
1035  return retval;
1036 
1037  /* Convert from size_t and ptrdiff_t to hssize_t, and hsize_t. */
1038  for (i = 0; i < var->ndims; i++)
1039  {
1040  start[i] = startp[i];
1041  count[i] = countp[i];
1042  }
1043 
1044  /* Open this dataset if necessary, also checking for a weird case:
1045  * a non-coordinate (and non-scalar) variable that has the same
1046  * name as a dimension. */
1047  if (var->hdf5_name && strlen(var->hdf5_name) >= strlen(NON_COORD_PREPEND) &&
1048  strncmp(var->hdf5_name, NON_COORD_PREPEND, strlen(NON_COORD_PREPEND)) == 0 &&
1049  var->ndims)
1050  name_to_use = var->hdf5_name;
1051  else
1052  name_to_use = var->name;
1053  if (!var->hdf_datasetid)
1054  if ((var->hdf_datasetid = H5Dopen2(grp->hdf_grpid, name_to_use, H5P_DEFAULT)) < 0)
1055  return NC_ENOTVAR;
1056 
1057  /* Get file space of data. */
1058  if ((file_spaceid = H5Dget_space(var->hdf_datasetid)) < 0)
1059  BAIL(NC_EHDFERR);
1060 
1061  /* Check to ensure the user selection is
1062  * valid. H5Sget_simple_extent_dims gets the sizes of all the dims
1063  * and put them in fdims. */
1064  if (H5Sget_simple_extent_dims(file_spaceid, fdims, fmaxdims) < 0)
1065  BAIL(NC_EHDFERR);
1066 
1067 #ifdef LOGGING
1068  log_dim_info(var, fdims, fmaxdims, start, count);
1069 #endif
1070 
1071  /* Check dimension bounds. Remember that unlimited dimensions can
1072  * put data beyond their current length. */
1073  for (d2 = 0; d2 < var->ndims; d2++) {
1074  dim = var->dim[d2];
1075  assert(dim && dim->dimid == var->dimids[d2]);
1076  if (dim->unlimited)
1077  {
1078  size_t ulen;
1079 
1080  /* We can't go beyond the largest current extent of
1081  the unlimited dim. */
1082  if ((retval = NC4_inq_dim(ncid, dim->dimid, NULL, &ulen)))
1083  BAIL(retval);
1084 
1085  /* Check for out of bound requests. */
1086 #ifdef RELAX_COORD_BOUND
1087  if (start[d2] > (hssize_t)ulen ||
1088  (start[d2] == (hssize_t)ulen && count[d2] > 0))
1089 #else
1090  if (start[d2] >= (hssize_t)ulen && ulen > 0)
1091 #endif
1092  BAIL_QUIET(NC_EINVALCOORDS);
1093  if (start[d2] + count[d2] > ulen)
1094  BAIL_QUIET(NC_EEDGE);
1095 
1096  /* Things get a little tricky here. If we're getting
1097  a GET request beyond the end of this var's
1098  current length in an unlimited dimension, we'll
1099  later need to return the fill value for the
1100  variable. */
1101  if (start[d2] >= (hssize_t)fdims[d2])
1102  fill_value_size[d2] = count[d2];
1103  else if (start[d2] + count[d2] > fdims[d2])
1104  fill_value_size[d2] = count[d2] - (fdims[d2] - start[d2]);
1105  else
1106  fill_value_size[d2] = 0;
1107  count[d2] -= fill_value_size[d2];
1108  if (fill_value_size[d2])
1109  provide_fill++;
1110  }
1111  else
1112  {
1113  /* Check for out of bound requests. */
1114 #ifdef RELAX_COORD_BOUND
1115  if (start[d2] > (hssize_t)fdims[d2] ||
1116  (start[d2] == (hssize_t)fdims[d2] && count[d2] > 0))
1117 #else
1118  if (start[d2] >= (hssize_t)fdims[d2])
1119 #endif
1120  BAIL_QUIET(NC_EINVALCOORDS);
1121  if (start[d2] + count[d2] > fdims[d2])
1122  BAIL_QUIET(NC_EEDGE);
1123 
1124  /* Set the fill value boundary */
1125  fill_value_size[d2] = count[d2];
1126  }
1127  }
1128 
1129  /* A little quirk: if any of the count values are zero, don't
1130  * read. */
1131  for (d2 = 0; d2 < var->ndims; d2++)
1132  if (count[d2] == 0)
1133  no_read++;
1134 
1135  /* Later on, we will need to know the size of this type in the
1136  * file. */
1137  assert(var->type_info->size);
1138  file_type_size = var->type_info->size;
1139 
1140  if (!no_read)
1141  {
1142  /* Now you would think that no one would be crazy enough to write
1143  a scalar dataspace with one of the array function calls, but you
1144  would be wrong. So let's check to see if the dataset is
1145  scalar. If it is, we won't try to set up a hyperslab. */
1146  if (H5Sget_simple_extent_type(file_spaceid) == H5S_SCALAR)
1147  {
1148  if ((mem_spaceid = H5Screate(H5S_SCALAR)) < 0)
1149  BAIL(NC_EHDFERR);
1150  scalar++;
1151  }
1152  else
1153  {
1154  if (H5Sselect_hyperslab(file_spaceid, H5S_SELECT_SET,
1155  start, NULL, count, NULL) < 0)
1156  BAIL(NC_EHDFERR);
1157  /* Create a space for the memory, just big enough to hold the slab
1158  we want. */
1159  if ((mem_spaceid = H5Screate_simple(var->ndims, count, NULL)) < 0)
1160  BAIL(NC_EHDFERR);
1161  }
1162 
1163  /* Fix bug when reading HDF5 files with variable of type
1164  * fixed-length string. We need to make it look like a
1165  * variable-length string, because that's all netCDF-4 data
1166  * model supports, lacking anonymous dimensions. So
1167  * variable-length strings are in allocated memory that user has
1168  * to free, which we allocate here. */
1169  if(var->type_info->nc_type_class == NC_STRING &&
1170  H5Tget_size(var->type_info->hdf_typeid) > 1 &&
1171  !H5Tis_variable_str(var->type_info->hdf_typeid)) {
1172  hsize_t fstring_len;
1173 
1174  if ((fstring_len = H5Tget_size(var->type_info->hdf_typeid)) == 0)
1175  BAIL(NC_EHDFERR);
1176  if (!(*(char **)data = malloc(1 + fstring_len)))
1177  BAIL(NC_ENOMEM);
1178  bufr = *(char **)data;
1179  }
1180 
1181 #ifndef HDF5_CONVERT
1182  /* Are we going to convert any data? (No converting of compound or
1183  * opaque types.) */
1184  if ((mem_nc_type != var->type_info->nc_typeid || (var->type_info->nc_typeid == NC_INT && is_long)) &&
1185  mem_nc_type != NC_COMPOUND && mem_nc_type != NC_OPAQUE)
1186  {
1187  /* We must convert - allocate a buffer. */
1188  need_to_convert++;
1189  if (var->ndims)
1190  for (d2 = 0; d2 < var->ndims; d2++)
1191  len *= countp[d2];
1192  LOG((4, "converting data for var %s type=%d len=%d", var->name,
1193  var->type_info->nc_typeid, len));
1194 
1195  /* If we're reading, we need bufr to have enough memory to store
1196  * the data in the file. If we're writing, we need bufr to be
1197  * big enough to hold all the data in the file's type. */
1198  if(len > 0)
1199  if (!(bufr = malloc(len * file_type_size)))
1200  BAIL(NC_ENOMEM);
1201  }
1202  else
1203 #endif /* ifndef HDF5_CONVERT */
1204  if(!bufr)
1205  bufr = data;
1206 
1207  /* Get the HDF type of the data in memory. */
1208 #ifdef HDF5_CONVERT
1209  if ((retval = nc4_get_hdf_typeid(h5, mem_nc_type, &mem_typeid,
1210  var->type_info->endianness)))
1211  BAIL(retval);
1212 #endif
1213 
1214  /* Create the data transfer property list. */
1215  if ((xfer_plistid = H5Pcreate(H5P_DATASET_XFER)) < 0)
1216  BAIL(NC_EHDFERR);
1217 
1218 #ifdef HDF5_CONVERT
1219  /* Apply the callback function which will detect range
1220  * errors. Which one to call depends on the length of the
1221  * destination buffer type. */
1222  if (H5Pset_type_conv_cb(xfer_plistid, except_func, &range_error) < 0)
1223  BAIL(NC_EHDFERR);
1224 #endif
1225 
1226 #ifdef USE_PARALLEL4
1227  /* Set up parallel I/O, if needed. */
1228  if ((retval = set_par_access(h5, var, xfer_plistid)))
1229  BAIL(retval);
1230 #endif
1231 
1232  /* Read this hyperslab into memory. */
1233  LOG((5, "About to H5Dread some data..."));
1234  if (H5Dread(var->hdf_datasetid, var->type_info->native_hdf_typeid,
1235  mem_spaceid, file_spaceid, xfer_plistid, bufr) < 0)
1236  BAIL(NC_EHDFERR);
1237 
1238 #ifndef HDF5_CONVERT
1239  /* Eventually the block below will go away. Right now it's
1240  needed to support conversions between int/float, and range
1241  checking converted data in the netcdf way. These features are
1242  being added to HDF5 at the HDF5 World Hall of Coding right
1243  now, by a staff of thousands of programming gnomes. */
1244  if (need_to_convert)
1245  {
1246  if ((retval = nc4_convert_type(bufr, data, var->type_info->nc_typeid, mem_nc_type,
1247  len, &range_error, var->fill_value,
1248  (h5->cmode & NC_CLASSIC_MODEL), 0, is_long)))
1249  BAIL(retval);
1250 
1251  /* For strict netcdf-3 rules, ignore erange errors between UBYTE
1252  * and BYTE types. */
1253  if ((h5->cmode & NC_CLASSIC_MODEL) &&
1254  (var->type_info->nc_typeid == NC_UBYTE || var->type_info->nc_typeid == NC_BYTE) &&
1255  (mem_nc_type == NC_UBYTE || mem_nc_type == NC_BYTE) &&
1256  range_error)
1257  range_error = 0;
1258  }
1259 #endif
1260 
1261  /* For strict netcdf-3 rules, ignore erange errors between UBYTE
1262  * and BYTE types. */
1263  if ((h5->cmode & NC_CLASSIC_MODEL) &&
1264  (var->type_info->nc_typeid == NC_UBYTE || var->type_info->nc_typeid == NC_BYTE) &&
1265  (mem_nc_type == NC_UBYTE || mem_nc_type == NC_BYTE) &&
1266  range_error)
1267  range_error = 0;
1268 
1269  } /* endif ! no_read */
1270 
1271  else {
1272 #ifdef USE_PARALLEL4 /* Start block contributed by HDF group. */
1273  /* For collective IO read, some processes may not have any element for reading.
1274  Collective requires all processes to participate, so we use H5Sselect_none
1275  for these processes. */
1276  if(var->parallel_access == NC_COLLECTIVE) {
1277 
1278  /* Create the data transfer property list. */
1279  if ((xfer_plistid = H5Pcreate(H5P_DATASET_XFER)) < 0)
1280  BAIL(NC_EHDFERR);
1281 
1282  if ((retval = set_par_access(h5, var, xfer_plistid)))
1283  BAIL(retval);
1284 
1285  if (H5Sselect_none(file_spaceid)<0)
1286  BAIL(NC_EHDFERR);
1287 
1288  /* Since no element will be selected, we just get the memory space the same as the file space.
1289  */
1290  if((mem_spaceid = H5Dget_space(var->hdf_datasetid))<0)
1291  BAIL(NC_EHDFERR);
1292  if (H5Sselect_none(mem_spaceid)<0)
1293  BAIL(NC_EHDFERR);
1294 
1295  /* Read this hyperslab into memory. */
1296  LOG((5, "About to H5Dread some data..."));
1297  if (H5Dread(var->hdf_datasetid, var->type_info->native_hdf_typeid,
1298  mem_spaceid, file_spaceid, xfer_plistid, bufr) < 0)
1299  BAIL(NC_EHDFERR);
1300  }
1301 #endif /* End ifdef USE_PARALLEL4 */
1302  }
1303  /* Now we need to fake up any further data that was asked for,
1304  using the fill values instead. First skip past the data we
1305  just read, if any. */
1306  if (!scalar && provide_fill)
1307  {
1308  void *filldata;
1309  size_t real_data_size = 0;
1310  size_t fill_len;
1311 
1312  /* Skip past the real data we've already read. */
1313  if (!no_read)
1314  for (real_data_size = file_type_size, d2 = 0; d2 < var->ndims; d2++)
1315  real_data_size *= (count[d2] - start[d2]);
1316 
1317  /* Get the fill value from the HDF5 variable. Memory will be
1318  * allocated. */
1319  if (get_fill_value(h5, var, &fillvalue) < 0)
1320  BAIL(NC_EHDFERR);
1321 
1322  /* How many fill values do we need? */
1323  for (fill_len = 1, d2 = 0; d2 < var->ndims; d2++)
1324  fill_len *= (fill_value_size[d2] ? fill_value_size[d2] : 1);
1325 
1326  /* Copy the fill value into the rest of the data buffer. */
1327  filldata = (char *)data + real_data_size;
1328  for (i = 0; i < fill_len; i++)
1329  {
1330 
1331  if (var->type_info->nc_type_class == NC_STRING)
1332  {
1333  if (*(char **)fillvalue)
1334  {
1335  if (!(*(char **)filldata = strdup(*(char **)fillvalue)))
1336  BAIL(NC_ENOMEM);
1337  }
1338  else
1339  *(char **)filldata = NULL;
1340  }
1341  else if(var->type_info->nc_type_class == NC_VLEN) {
1342  if(fillvalue) {
1343  memcpy(filldata,fillvalue,file_type_size);
1344  } else {
1345  *(char **)filldata = NULL;
1346  }
1347  } else
1348  memcpy(filldata, fillvalue, file_type_size);
1349  filldata = (char *)filldata + file_type_size;
1350  }
1351  }
1352 
1353 exit:
1354 #ifdef HDF5_CONVERT
1355  if (mem_typeid > 0 && H5Tclose(mem_typeid) < 0)
1356  BAIL2(NC_EHDFERR);
1357 #endif
1358  if (file_spaceid > 0)
1359  {
1360  if (H5Sclose(file_spaceid) < 0)
1361  BAIL2(NC_EHDFERR);
1362  }
1363  if (mem_spaceid > 0)
1364  {
1365  if (H5Sclose(mem_spaceid) < 0)
1366  BAIL2(NC_EHDFERR);
1367  }
1368  if (xfer_plistid > 0)
1369  {
1370  if (H5Pclose(xfer_plistid) < 0)
1371  BAIL2(NC_EHDFERR);
1372  }
1373 #ifndef HDF5_CONVERT
1374  if (need_to_convert && bufr != NULL)
1375  free(bufr);
1376 #endif
1377  if (xtend_size)
1378  free(xtend_size);
1379  if (fillvalue)
1380  {
1381  if (var->type_info->nc_type_class == NC_VLEN)
1382  nc_free_vlen((nc_vlen_t *)fillvalue);
1383  else if (var->type_info->nc_type_class == NC_STRING && *(char **)fillvalue)
1384  free(*(char **)fillvalue);
1385  free(fillvalue);
1386  }
1387 
1388  /* If there was an error return it, otherwise return any potential
1389  range error value. If none, return NC_NOERR as usual.*/
1390  if (retval)
1391  return retval;
1392  if (range_error)
1393  return NC_ERANGE;
1394  return NC_NOERR;
1395 }
1396 
1411 static int
1412 put_att_grpa(NC_GRP_INFO_T *grp, int varid, NC_ATT_INFO_T *att)
1413 {
1414  hid_t datasetid = 0, locid;
1415  hid_t attid = 0, spaceid = 0, file_typeid = 0;
1416  hsize_t dims[1]; /* netcdf attributes always 1-D. */
1417  htri_t attr_exists;
1418  int retval = NC_NOERR;
1419  void *data;
1420  int phoney_data = 99;
1421 
1422  assert(att->name);
1423  LOG((3, "%s: varid %d att->attnum %d att->name %s att->nc_typeid %d att->len %d",
1424  __func__, varid, att->attnum, att->name,
1425  att->nc_typeid, att->len));
1426 
1427  /* If the file is read-only, return an error. */
1428  if (grp->nc4_info->no_write)
1429  BAIL(NC_EPERM);
1430 
1431  /* Get the hid to attach the attribute to, or read it from. */
1432  if (varid == NC_GLOBAL)
1433  locid = grp->hdf_grpid;
1434  else
1435  {
1436  if ((retval = nc4_open_var_grp2(grp, varid, &datasetid)))
1437  BAIL(retval);
1438  locid = datasetid;
1439  }
1440 
1441  /* Delete the att if it exists already. */
1442  if ((attr_exists = H5Aexists(locid, att->name)) < 0)
1443  BAIL(NC_EHDFERR);
1444  if (attr_exists)
1445  {
1446  if (H5Adelete(locid, att->name) < 0)
1447  BAIL(NC_EHDFERR);
1448  }
1449 
1450  /* Get the length ready, and find the HDF type we'll be
1451  * writing. */
1452  dims[0] = att->len;
1453  if ((retval = nc4_get_hdf_typeid(grp->nc4_info, att->nc_typeid,
1454  &file_typeid, 0)))
1455  BAIL(retval);
1456 
1457  /* Even if the length is zero, HDF5 won't let me write with a
1458  * NULL pointer. So if the length of the att is zero, point to
1459  * some phoney data (which won't be written anyway.)*/
1460  if (!dims[0])
1461  data = &phoney_data;
1462  else if (att->data)
1463  data = att->data;
1464  else if (att->stdata)
1465  data = att->stdata;
1466  else
1467  data = att->vldata;
1468 
1469  /* NC_CHAR types require some extra work. The space ID is set to
1470  * scalar, and the type is told how long the string is. If it's
1471  * really zero length, set the size to 1. (The fact that it's
1472  * really zero will be marked by the NULL dataspace, but HDF5
1473  * doesn't allow me to set the size of the type to zero.)*/
1474  if (att->nc_typeid == NC_CHAR)
1475  {
1476  size_t string_size = dims[0];
1477  if (!string_size)
1478  {
1479  string_size = 1;
1480  if ((spaceid = H5Screate(H5S_NULL)) < 0)
1481  BAIL(NC_EATTMETA);
1482  }
1483  else
1484  {
1485  if ((spaceid = H5Screate(H5S_SCALAR)) < 0)
1486  BAIL(NC_EATTMETA);
1487  }
1488  if (H5Tset_size(file_typeid, string_size) < 0)
1489  BAIL(NC_EATTMETA);
1490  if (H5Tset_strpad(file_typeid, H5T_STR_NULLTERM) < 0)
1491  BAIL(NC_EATTMETA);
1492  }
1493  else
1494  {
1495  if (!att->len)
1496  {
1497  if ((spaceid = H5Screate(H5S_NULL)) < 0)
1498  BAIL(NC_EATTMETA);
1499  }
1500  else
1501  {
1502  if ((spaceid = H5Screate_simple(1, dims, NULL)) < 0)
1503  BAIL(NC_EATTMETA);
1504  }
1505  }
1506  if ((attid = H5Acreate(locid, att->name, file_typeid, spaceid,
1507  H5P_DEFAULT)) < 0)
1508  BAIL(NC_EATTMETA);
1509 
1510  /* Write the values, (even if length is zero). */
1511  if (H5Awrite(attid, file_typeid, data) < 0)
1512  BAIL(NC_EATTMETA);
1513 
1514 exit:
1515  if (file_typeid && H5Tclose(file_typeid))
1516  BAIL2(NC_EHDFERR);
1517  if (attid > 0 && H5Aclose(attid) < 0)
1518  BAIL2(NC_EHDFERR);
1519  if (spaceid > 0 && H5Sclose(spaceid) < 0)
1520  BAIL2(NC_EHDFERR);
1521  return retval;
1522 }
1523 
1535 static int
1536 write_attlist(NC_ATT_INFO_T *attlist, int varid, NC_GRP_INFO_T *grp)
1537 {
1538  NC_ATT_INFO_T *att;
1539  int retval;
1540 
1541  for (att = attlist; att; att = att->l.next)
1542  {
1543  if (att->dirty)
1544  {
1545  LOG((4, "%s: writing att %s to varid %d", __func__, att->name, varid));
1546  if ((retval = put_att_grpa(grp, varid, att)))
1547  return retval;
1548  att->dirty = NC_FALSE;
1549  att->created = NC_TRUE;
1550  }
1551  }
1552  return NC_NOERR;
1553 }
1554 
1555 
1574 static int
1575 write_coord_dimids(NC_VAR_INFO_T *var)
1576 {
1577  hsize_t coords_len[1];
1578  hid_t c_spaceid = -1, c_attid = -1;
1579  int ret = 0;
1580 
1581  /* Write our attribute. */
1582  coords_len[0] = var->ndims;
1583  if ((c_spaceid = H5Screate_simple(1, coords_len, coords_len)) < 0) ret++;
1584  if (!ret && (c_attid = H5Acreate(var->hdf_datasetid, COORDINATES, H5T_NATIVE_INT,
1585  c_spaceid, H5P_DEFAULT)) < 0) ret++;
1586  if (!ret && H5Awrite(c_attid, H5T_NATIVE_INT, var->dimids) < 0) ret++;
1587 
1588  /* Close up shop. */
1589  if (c_spaceid > 0 && H5Sclose(c_spaceid) < 0) ret++;
1590  if (c_attid > 0 && H5Aclose(c_attid) < 0) ret++;
1591  return ret ? NC_EHDFERR : 0;
1592 }
1593 
1604 static int
1605 write_netcdf4_dimid(hid_t datasetid, int dimid)
1606 {
1607  hid_t dimid_spaceid, dimid_attid;
1608  htri_t attr_exists;
1609 
1610  /* Create the space. */
1611  if ((dimid_spaceid = H5Screate(H5S_SCALAR)) < 0)
1612  return NC_EHDFERR;
1613 
1614  /* Does the attribute already exist? If so, don't try to create it. */
1615  if ((attr_exists = H5Aexists(datasetid, NC_DIMID_ATT_NAME)) < 0)
1616  return NC_EHDFERR;
1617  if (attr_exists)
1618  dimid_attid = H5Aopen_by_name(datasetid, ".", NC_DIMID_ATT_NAME,
1619  H5P_DEFAULT, H5P_DEFAULT);
1620  else
1621  /* Create the attribute if needed. */
1622  dimid_attid = H5Acreate(datasetid, NC_DIMID_ATT_NAME,
1623  H5T_NATIVE_INT, dimid_spaceid, H5P_DEFAULT);
1624  if (dimid_attid < 0)
1625  return NC_EHDFERR;
1626 
1627 
1628  /* Write it. */
1629  LOG((4, "%s: writing secret dimid %d", __func__, dimid));
1630  if (H5Awrite(dimid_attid, H5T_NATIVE_INT, &dimid) < 0)
1631  return NC_EHDFERR;
1632 
1633  /* Close stuff*/
1634  if (H5Sclose(dimid_spaceid) < 0)
1635  return NC_EHDFERR;
1636  if (H5Aclose(dimid_attid) < 0)
1637  return NC_EHDFERR;
1638 
1639  return NC_NOERR;
1640 }
1641 
1652 static int
1653 var_create_dataset(NC_GRP_INFO_T *grp, NC_VAR_INFO_T *var, nc_bool_t write_dimid)
1654 {
1655  hid_t plistid = 0, access_plistid = 0, typeid = 0, spaceid = 0;
1656  hsize_t chunksize[H5S_MAX_RANK], dimsize[H5S_MAX_RANK], maxdimsize[H5S_MAX_RANK];
1657  int d;
1658  void *fillp = NULL;
1659  NC_DIM_INFO_T *dim = NULL;
1660  char *name_to_use;
1661  int retval = NC_NOERR;
1662 
1663  LOG((3, "%s:: name %s", __func__, var->name));
1664 
1665  /* Scalar or not, we need a creation property list. */
1666  if ((plistid = H5Pcreate(H5P_DATASET_CREATE)) < 0)
1667  BAIL(NC_EHDFERR);
1668  if ((access_plistid = H5Pcreate(H5P_DATASET_ACCESS)) < 0)
1669  BAIL(NC_EHDFERR);
1670 
1671  /* RJ: this suppose to be FALSE that is defined in H5 private.h as 0 */
1672  if (H5Pset_obj_track_times(plistid,0)<0)
1673  BAIL(NC_EHDFERR);
1674 
1675  /* Find the HDF5 type of the dataset. */
1676  if ((retval = nc4_get_hdf_typeid(grp->nc4_info, var->type_info->nc_typeid, &typeid,
1677  var->type_info->endianness)))
1678  BAIL(retval);
1679 
1680  /* Figure out what fill value to set, if any. */
1681  if (var->no_fill)
1682  {
1683  /* Required to truly turn HDF5 fill values off */
1684  if (H5Pset_fill_time(plistid, H5D_FILL_TIME_NEVER) < 0)
1685  BAIL(NC_EHDFERR);
1686  }
1687  else
1688  {
1689  if ((retval = get_fill_value(grp->nc4_info, var, &fillp)))
1690  BAIL(retval);
1691 
1692  /* If there is a fill value, set it. */
1693  if (fillp)
1694  {
1695  if (var->type_info->nc_type_class == NC_STRING)
1696  {
1697  if (H5Pset_fill_value(plistid, typeid, fillp) < 0)
1698  BAIL(NC_EHDFERR);
1699  }
1700  else
1701  {
1702  /* The fill value set in HDF5 must always be presented as
1703  * a native type, even if the endianness for this dataset
1704  * is non-native. HDF5 will translate the fill value to
1705  * the target endiannesss. */
1706  hid_t fill_typeid = 0;
1707 
1708  if ((retval = nc4_get_hdf_typeid(grp->nc4_info, var->type_info->nc_typeid, &fill_typeid,
1709  NC_ENDIAN_NATIVE)))
1710  BAIL(retval);
1711  if (H5Pset_fill_value(plistid, fill_typeid, fillp) < 0)
1712  {
1713  if (H5Tclose(fill_typeid) < 0)
1714  BAIL(NC_EHDFERR);
1715  BAIL(NC_EHDFERR);
1716  }
1717  if (H5Tclose(fill_typeid) < 0)
1718  BAIL(NC_EHDFERR);
1719  }
1720  }
1721  }
1722 
1723  /* If the user wants to shuffle the data, set that up now. */
1724  if (var->shuffle) {
1725  if (H5Pset_shuffle(plistid) < 0)
1726  BAIL(NC_EHDFERR);
1727  }
1728 
1729  /* If the user wants to deflate the data, set that up now. */
1730  if (var->deflate) {
1731  if (H5Pset_deflate(plistid, var->deflate_level) < 0)
1732  BAIL(NC_EHDFERR);
1733  } else if(var->filterid) {
1734  /* Handle szip case here */
1735  if(var->filterid == H5Z_FILTER_SZIP) {
1736  int options_mask;
1737  int bits_per_pixel;
1738  if(var->nparams != 2)
1739  BAIL(NC_EFILTER);
1740  options_mask = (int)var->params[0];
1741  bits_per_pixel = (int)var->params[1];
1742  if(H5Pset_szip(plistid, options_mask, bits_per_pixel) < 0)
1743  BAIL(NC_EFILTER);
1744  } else {
1745  herr_t code = H5Pset_filter(plistid, var->filterid, H5Z_FLAG_MANDATORY, var->nparams, var->params);
1746  if(code < 0) {
1747  BAIL(NC_EFILTER);
1748  }
1749  }
1750  }
1751 
1752  /* If the user wants to fletcher error correcton, set that up now. */
1753  if (var->fletcher32)
1754  if (H5Pset_fletcher32(plistid) < 0)
1755  BAIL(NC_EHDFERR);
1756 
1757  /* If ndims non-zero, get info for all dimensions. We look up the
1758  dimids and get the len of each dimension. We need this to create
1759  the space for the dataset. In netCDF a dimension length of zero
1760  means an unlimited dimension. */
1761  if (var->ndims)
1762  {
1763  int unlimdim = 0;
1764 
1765  /* Check to see if any unlimited dimensions are used in this var. */
1766  for (d = 0; d < var->ndims; d++) {
1767  dim = var->dim[d];
1768  assert(dim && dim->dimid == var->dimids[d]);
1769  if (dim->unlimited)
1770  unlimdim++;
1771  }
1772 
1773  /* If there are no unlimited dims, and no filters, and the user
1774  * has not specified chunksizes, use contiguous variable for
1775  * better performance. */
1776  if (!var->shuffle && !var->deflate && !var->fletcher32 &&
1777  (var->chunksizes == NULL || !var->chunksizes[0]) && !unlimdim)
1778  var->contiguous = NC_TRUE;
1779 
1780  /* Gather current & maximum dimension sizes, along with chunk sizes */
1781  for (d = 0; d < var->ndims; d++)
1782  {
1783  dim = var->dim[d];
1784  assert(dim && dim->dimid == var->dimids[d]);
1785  dimsize[d] = dim->unlimited ? NC_HDF5_UNLIMITED_DIMSIZE : dim->len;
1786  maxdimsize[d] = dim->unlimited ? H5S_UNLIMITED : (hsize_t)dim->len;
1787  if (!var->contiguous) {
1788  if (var->chunksizes[d])
1789  chunksize[d] = var->chunksizes[d];
1790  else
1791  {
1792  size_t type_size;
1793  if (var->type_info->nc_type_class == NC_STRING)
1794  type_size = sizeof(char *);
1795  else
1796  type_size = var->type_info->size;
1797 
1798  /* Unlimited dim always gets chunksize of 1. */
1799  if (dim->unlimited)
1800  chunksize[d] = 1;
1801  else
1802  chunksize[d] = pow((double)DEFAULT_CHUNK_SIZE/type_size,
1803  1/(double)(var->ndims - unlimdim));
1804 
1805  /* If the chunksize is greater than the dim
1806  * length, make it the dim length. */
1807  if (!dim->unlimited && chunksize[d] > dim->len)
1808  chunksize[d] = dim->len;
1809 
1810  /* Remember the computed chunksize */
1811  var->chunksizes[d] = chunksize[d];
1812  }
1813  }
1814  }
1815 
1816  if (var->contiguous)
1817  {
1818  if (H5Pset_layout(plistid, H5D_CONTIGUOUS) < 0)
1819  BAIL(NC_EHDFERR);
1820  }
1821  else
1822  {
1823  if (H5Pset_chunk(plistid, var->ndims, chunksize) < 0)
1824  BAIL(NC_EHDFERR);
1825  }
1826 
1827  /* Create the dataspace. */
1828  if ((spaceid = H5Screate_simple(var->ndims, dimsize, maxdimsize)) < 0)
1829  BAIL(NC_EHDFERR);
1830  }
1831  else
1832  {
1833  if ((spaceid = H5Screate(H5S_SCALAR)) < 0)
1834  BAIL(NC_EHDFERR);
1835  }
1836 
1837  /* Turn on creation order tracking. */
1838  if (H5Pset_attr_creation_order(plistid, H5P_CRT_ORDER_TRACKED|
1839  H5P_CRT_ORDER_INDEXED) < 0)
1840  BAIL(NC_EHDFERR);
1841 
1842  /* Set per-var chunk cache, for chunked datasets. */
1843  if (!var->contiguous && var->chunk_cache_size)
1844  if (H5Pset_chunk_cache(access_plistid, var->chunk_cache_nelems,
1845  var->chunk_cache_size, var->chunk_cache_preemption) < 0)
1846  BAIL(NC_EHDFERR);
1847 
1848  /* At long last, create the dataset. */
1849  name_to_use = var->hdf5_name ? var->hdf5_name : var->name;
1850  LOG((4, "%s: about to H5Dcreate2 dataset %s of type 0x%x", __func__,
1851  name_to_use, typeid));
1852  if ((var->hdf_datasetid = H5Dcreate2(grp->hdf_grpid, name_to_use, typeid,
1853  spaceid, H5P_DEFAULT, plistid, access_plistid)) < 0)
1854  BAIL(NC_EHDFERR);
1855  var->created = NC_TRUE;
1856  var->is_new_var = NC_FALSE;
1857 
1858  /* If this is a dimscale, mark it as such in the HDF5 file. Also
1859  * find the dimension info and store the dataset id of the dimscale
1860  * dataset. */
1861  if (var->dimscale)
1862  {
1863  if (H5DSset_scale(var->hdf_datasetid, var->name) < 0)
1864  BAIL(NC_EHDFERR);
1865 
1866  /* If this is a multidimensional coordinate variable, write a
1867  * coordinates attribute. */
1868  if (var->ndims > 1)
1869  if ((retval = write_coord_dimids(var)))
1870  BAIL(retval);
1871 
1872  /* If desired, write the netCDF dimid. */
1873  if (write_dimid)
1874  if ((retval = write_netcdf4_dimid(var->hdf_datasetid, var->dimids[0])))
1875  BAIL(retval);
1876  }
1877 
1878 
1879  /* Write attributes for this var. */
1880  if ((retval = write_attlist(var->att, var->varid, grp)))
1881  BAIL(retval);
1882  var->attr_dirty = NC_FALSE;
1883 
1884 exit:
1885  if (typeid > 0 && H5Tclose(typeid) < 0)
1886  BAIL2(NC_EHDFERR);
1887  if (plistid > 0 && H5Pclose(plistid) < 0)
1888  BAIL2(NC_EHDFERR);
1889  if (access_plistid > 0 && H5Pclose(access_plistid) < 0)
1890  BAIL2(NC_EHDFERR);
1891  if (spaceid > 0 && H5Sclose(spaceid) < 0)
1892  BAIL2(NC_EHDFERR);
1893  if (fillp)
1894  {
1895  if (var->type_info->nc_type_class == NC_VLEN)
1896  nc_free_vlen((nc_vlen_t *)fillp);
1897  else if (var->type_info->nc_type_class == NC_STRING && *(char **)fillp)
1898  free(*(char **)fillp);
1899  free(fillp);
1900  }
1901 
1902  return retval;
1903 }
1904 
1915 int
1916 nc4_adjust_var_cache(NC_GRP_INFO_T *grp, NC_VAR_INFO_T * var)
1917 {
1918  size_t chunk_size_bytes = 1;
1919  int d;
1920  int retval;
1921 
1922  /* Nothing to be done. */
1923  if (var->contiguous)
1924  return NC_NOERR;
1925 #ifdef USE_PARALLEL4
1926  return NC_NOERR;
1927 #endif
1928 
1929  /* How many bytes in the chunk? */
1930  for (d = 0; d < var->ndims; d++)
1931  chunk_size_bytes *= var->chunksizes[d];
1932  if (var->type_info->size)
1933  chunk_size_bytes *= var->type_info->size;
1934  else
1935  chunk_size_bytes *= sizeof(char *);
1936 
1937  /* If the chunk cache is too small, and the user has not changed
1938  * the default value of the chunk cache size, then increase the
1939  * size of the cache. */
1940  if (var->chunk_cache_size == CHUNK_CACHE_SIZE)
1941  if (chunk_size_bytes > var->chunk_cache_size)
1942  {
1943  var->chunk_cache_size = chunk_size_bytes * DEFAULT_CHUNKS_IN_CACHE;
1944  if (var->chunk_cache_size > MAX_DEFAULT_CACHE_SIZE)
1945  var->chunk_cache_size = MAX_DEFAULT_CACHE_SIZE;
1946  if ((retval = nc4_reopen_dataset(grp, var)))
1947  return retval;
1948  }
1949 
1950  return NC_NOERR;
1951 }
1952 
1963 static int
1964 commit_type(NC_GRP_INFO_T *grp, NC_TYPE_INFO_T *type)
1965 {
1966  int retval;
1967 
1968  assert(grp && type);
1969 
1970  /* Did we already record this type? */
1971  if (type->committed)
1972  return NC_NOERR;
1973 
1974  /* Is this a compound type? */
1975  if (type->nc_type_class == NC_COMPOUND)
1976  {
1977  NC_FIELD_INFO_T *field;
1978  hid_t hdf_base_typeid, hdf_typeid;
1979 
1980  if ((type->hdf_typeid = H5Tcreate(H5T_COMPOUND, type->size)) < 0)
1981  return NC_EHDFERR;
1982  LOG((4, "creating compound type %s hdf_typeid 0x%x", type->name,
1983  type->hdf_typeid));
1984 
1985  for (field = type->u.c.field; field; field = field->l.next)
1986  {
1987  if ((retval = nc4_get_hdf_typeid(grp->nc4_info, field->nc_typeid,
1988  &hdf_base_typeid, type->endianness)))
1989  return retval;
1990 
1991  /* If this is an array, create a special array type. */
1992  if (field->ndims)
1993  {
1994  int d;
1995  hsize_t dims[NC_MAX_VAR_DIMS];
1996 
1997  for (d = 0; d < field->ndims; d++)
1998  dims[d] = field->dim_size[d];
1999  if ((hdf_typeid = H5Tarray_create(hdf_base_typeid, field->ndims,
2000  dims, NULL)) < 0)
2001  {
2002  if (H5Tclose(hdf_base_typeid) < 0)
2003  return NC_EHDFERR;
2004  return NC_EHDFERR;
2005  }
2006  if (H5Tclose(hdf_base_typeid) < 0)
2007  return NC_EHDFERR;
2008  }
2009  else
2010  hdf_typeid = hdf_base_typeid;
2011  LOG((4, "inserting field %s offset %d hdf_typeid 0x%x", field->name,
2012  field->offset, hdf_typeid));
2013  if (H5Tinsert(type->hdf_typeid, field->name, field->offset,
2014  hdf_typeid) < 0)
2015  return NC_EHDFERR;
2016  if (H5Tclose(hdf_typeid) < 0)
2017  return NC_EHDFERR;
2018  }
2019  }
2020  else if (type->nc_type_class == NC_VLEN)
2021  {
2022  /* Find the HDF typeid of the base type of this vlen. */
2023  if ((retval = nc4_get_hdf_typeid(grp->nc4_info, type->u.v.base_nc_typeid,
2024  &type->u.v.base_hdf_typeid, type->endianness)))
2025  return retval;
2026 
2027  /* Create a vlen type. */
2028  if ((type->hdf_typeid = H5Tvlen_create(type->u.v.base_hdf_typeid)) < 0)
2029  return NC_EHDFERR;
2030  }
2031  else if (type->nc_type_class == NC_OPAQUE)
2032  {
2033  /* Create the opaque type. */
2034  if ((type->hdf_typeid = H5Tcreate(H5T_OPAQUE, type->size)) < 0)
2035  return NC_EHDFERR;
2036  }
2037  else if (type->nc_type_class == NC_ENUM)
2038  {
2039  NC_ENUM_MEMBER_INFO_T *enum_m;
2040 
2041  if (!type->u.e.enum_member)
2042  return NC_EINVAL;
2043 
2044  /* Find the HDF typeid of the base type of this enum. */
2045  if ((retval = nc4_get_hdf_typeid(grp->nc4_info, type->u.e.base_nc_typeid,
2046  &type->u.e.base_hdf_typeid, type->endianness)))
2047  return retval;
2048 
2049  /* Create an enum type. */
2050  if ((type->hdf_typeid = H5Tenum_create(type->u.e.base_hdf_typeid)) < 0)
2051  return NC_EHDFERR;
2052 
2053  /* Add all the members to the HDF5 type. */
2054  for (enum_m = type->u.e.enum_member; enum_m; enum_m = enum_m->l.next)
2055  if (H5Tenum_insert(type->hdf_typeid, enum_m->name, enum_m->value) < 0)
2056  return NC_EHDFERR;
2057  }
2058  else
2059  {
2060  LOG((0, "Unknown class: %d", type->nc_type_class));
2061  return NC_EBADTYPE;
2062  }
2063 
2064  /* Commit the type. */
2065  if (H5Tcommit(grp->hdf_grpid, type->name, type->hdf_typeid) < 0)
2066  return NC_EHDFERR;
2067  type->committed = NC_TRUE;
2068  LOG((4, "just committed type %s, HDF typeid: 0x%x", type->name,
2069  type->hdf_typeid));
2070 
2071  /* Later we will always use the native typeid. In this case, it is
2072  * a copy of the same type pointed to by hdf_typeid, but it's
2073  * easier to maintain a copy. */
2074  if ((type->native_hdf_typeid = H5Tget_native_type(type->hdf_typeid,
2075  H5T_DIR_DEFAULT)) < 0)
2076  return NC_EHDFERR;
2077 
2078  return NC_NOERR;
2079 }
2080 
2091 static int
2092 write_nc3_strict_att(hid_t hdf_grpid)
2093 {
2094  hid_t attid = 0, spaceid = 0;
2095  int one = 1;
2096  int retval = NC_NOERR;
2097  htri_t attr_exists;
2098 
2099  /* If the attribute already exists, call that a success and return
2100  * NC_NOERR. */
2101  if ((attr_exists = H5Aexists(hdf_grpid, NC3_STRICT_ATT_NAME)) < 0)
2102  return NC_EHDFERR;
2103  if (attr_exists)
2104  return NC_NOERR;
2105 
2106  /* Create the attribute to mark this as a file that needs to obey
2107  * strict netcdf-3 rules. */
2108  if ((spaceid = H5Screate(H5S_SCALAR)) < 0)
2109  BAIL(NC_EFILEMETA);
2110  if ((attid = H5Acreate(hdf_grpid, NC3_STRICT_ATT_NAME,
2111  H5T_NATIVE_INT, spaceid, H5P_DEFAULT)) < 0)
2112  BAIL(NC_EFILEMETA);
2113  if (H5Awrite(attid, H5T_NATIVE_INT, &one) < 0)
2114  BAIL(NC_EFILEMETA);
2115 
2116 exit:
2117  if (spaceid > 0 && (H5Sclose(spaceid) < 0))
2118  BAIL2(NC_EFILEMETA);
2119  if (attid > 0 && (H5Aclose(attid) < 0))
2120  BAIL2(NC_EFILEMETA);
2121  return retval;
2122 }
2123 
2132 static int
2133 create_group(NC_GRP_INFO_T *grp)
2134 {
2135  hid_t gcpl_id = 0;
2136  int retval = NC_NOERR;;
2137 
2138  assert(grp);
2139 
2140  /* If this is not the root group, create it in the HDF5 file. */
2141  if (grp->parent)
2142  {
2143  /* Create group, with link_creation_order set in the group
2144  * creation property list. */
2145  if ((gcpl_id = H5Pcreate(H5P_GROUP_CREATE)) < 0)
2146  return NC_EHDFERR;
2147 
2148  /* RJ: this suppose to be FALSE that is defined in H5 private.h as 0 */
2149  if (H5Pset_obj_track_times(gcpl_id,0)<0)
2150  BAIL(NC_EHDFERR);
2151 
2152  if (H5Pset_link_creation_order(gcpl_id, H5P_CRT_ORDER_TRACKED|H5P_CRT_ORDER_INDEXED) < 0)
2153  BAIL(NC_EHDFERR);
2154  if (H5Pset_attr_creation_order(gcpl_id, H5P_CRT_ORDER_TRACKED|H5P_CRT_ORDER_INDEXED) < 0)
2155  BAIL(NC_EHDFERR);
2156  if ((grp->hdf_grpid = H5Gcreate2(grp->parent->hdf_grpid, grp->name,
2157  H5P_DEFAULT, gcpl_id, H5P_DEFAULT)) < 0)
2158  BAIL(NC_EHDFERR);
2159  if (H5Pclose(gcpl_id) < 0)
2160  BAIL(NC_EHDFERR);
2161  }
2162  else
2163  {
2164  /* Since this is the root group, we have to open it. */
2165  if ((grp->hdf_grpid = H5Gopen2(grp->nc4_info->hdfid, "/", H5P_DEFAULT)) < 0)
2166  BAIL(NC_EFILEMETA);
2167  }
2168  return NC_NOERR;
2169 
2170 exit:
2171  if (gcpl_id > 0 && H5Pclose(gcpl_id) < 0)
2172  BAIL2(NC_EHDFERR);
2173  if (grp->hdf_grpid > 0 && H5Gclose(grp->hdf_grpid) < 0)
2174  BAIL2(NC_EHDFERR);
2175  return retval;
2176 }
2177 
2189 static int
2190 attach_dimscales(NC_GRP_INFO_T *grp)
2191 {
2192  NC_VAR_INFO_T *var;
2193  NC_DIM_INFO_T *dim1;
2194  int d, i;
2195  int retval = NC_NOERR;
2196 
2197  /* Attach dimension scales. */
2198  for (i=0; i < grp->vars.nelems; i++)
2199  {
2200  var = grp->vars.value[i];
2201  if (!var) continue;
2202  /* Scales themselves do not attach. But I really wish they
2203  * would. */
2204  if (var->dimscale)
2205  {
2206  /* If this is a multidimensional coordinate variable, it will
2207  * have a special coords attribute (read earlier) with a list
2208  * of the dimensions for this variable. */
2209  }
2210  else /* not a dimscale... */
2211  {
2212  /* Find the scale for each dimension and attach it. */
2213  for (d = 0; d < var->ndims; d++)
2214  {
2215  /* Is there a dimscale for this dimension? */
2216  if (var->dimscale_attached)
2217  {
2218  if (!var->dimscale_attached[d])
2219  {
2220  hid_t dim_datasetid; /* Dataset ID for dimension */
2221  dim1 = var->dim[d];
2222  assert(dim1 && dim1->dimid == var->dimids[d]);
2223 
2224  LOG((2, "%s: attaching scale for dimid %d to var %s",
2225  __func__, var->dimids[d], var->name));
2226 
2227  /* Find dataset ID for dimension */
2228  if (dim1->coord_var)
2229  dim_datasetid = dim1->coord_var->hdf_datasetid;
2230  else
2231  dim_datasetid = dim1->hdf_dimscaleid;
2232  assert(dim_datasetid > 0);
2233  if (H5DSattach_scale(var->hdf_datasetid, dim_datasetid, d) < 0)
2234  BAIL(NC_EHDFERR);
2235  var->dimscale_attached[d] = NC_TRUE;
2236  }
2237 
2238  /* If we didn't find a dimscale to attach, that's a problem! */
2239  if (!var->dimscale_attached[d])
2240  {
2241  LOG((0, "no dimscale found!"));
2242  return NC_EDIMSCALE;
2243  }
2244  }
2245  }
2246  }
2247  }
2248 
2249 exit:
2250  return retval;
2251 }
2252 
2263 static int
2264 var_exists(hid_t grpid, char *name, nc_bool_t *exists)
2265 {
2266  htri_t link_exists;
2267 
2268  /* Reset the boolean */
2269  *exists = NC_FALSE;
2270 
2271  /* Check if the object name exists in the group */
2272  if ((link_exists = H5Lexists(grpid, name, H5P_DEFAULT)) < 0)
2273  return NC_EHDFERR;
2274  if (link_exists)
2275  {
2276  H5G_stat_t statbuf;
2277 
2278  /* Get info about the object */
2279  if (H5Gget_objinfo(grpid, name, 1, &statbuf) < 0)
2280  return NC_EHDFERR;
2281 
2282  if (H5G_DATASET == statbuf.type)
2283  *exists = NC_TRUE;
2284  }
2285 
2286  return NC_NOERR;
2287 }
2288 
2303 int
2304 remove_coord_atts(hid_t hdf_datasetid)
2305 {
2306  htri_t attr_exists;
2307 
2308  /* If the variable dataset has an optional NC_DIMID_ATT_NAME
2309  * attribute, delete it. */
2310  if ((attr_exists = H5Aexists(hdf_datasetid, NC_DIMID_ATT_NAME)) < 0)
2311  return NC_EHDFERR;
2312  if (attr_exists)
2313  {
2314  if (H5Adelete(hdf_datasetid, NC_DIMID_ATT_NAME) < 0)
2315  return NC_EHDFERR;
2316  }
2317 
2318  /* (We could do a better job here and verify that the attributes are
2319  * really dimension scale 'CLASS' & 'NAME' attributes, but that would be
2320  * poking about in the HDF5 DimScale internal data) */
2321  if ((attr_exists = H5Aexists(hdf_datasetid, HDF5_DIMSCALE_CLASS_ATT_NAME)) < 0)
2322  return NC_EHDFERR;
2323  if (attr_exists)
2324  {
2325  if (H5Adelete(hdf_datasetid, HDF5_DIMSCALE_CLASS_ATT_NAME) < 0)
2326  return NC_EHDFERR;
2327  }
2328  if ((attr_exists = H5Aexists(hdf_datasetid, HDF5_DIMSCALE_NAME_ATT_NAME)) < 0)
2329  return NC_EHDFERR;
2330  if (attr_exists)
2331  {
2332  if (H5Adelete(hdf_datasetid, HDF5_DIMSCALE_NAME_ATT_NAME) < 0)
2333  return NC_EHDFERR;
2334  }
2335  return NC_NOERR;
2336 }
2337 
2352 static int
2353 write_var(NC_VAR_INFO_T *var, NC_GRP_INFO_T *grp, nc_bool_t write_dimid)
2354 {
2355  nc_bool_t replace_existing_var = NC_FALSE;
2356  int retval;
2357 
2358  LOG((4, "%s: writing var %s", __func__, var->name));
2359 
2360  /* If the variable has already been created & the fill value changed,
2361  * indicate that the existing variable should be replaced. */
2362  if (var->created && var->fill_val_changed)
2363  {
2364  replace_existing_var = NC_TRUE;
2365  var->fill_val_changed = NC_FALSE;
2366  /* If the variable is going to be replaced,
2367  we need to flag any other attributes associated
2368  with the variable as 'dirty', or else
2369  *only* the fill value attribute will be copied over
2370  and the rest will be lost. See:
2371 
2372  * https://github.com/Unidata/netcdf-c/issues/239 */
2373 
2374  flag_atts_dirty(&var->att);
2375  }
2376 
2377  /* Is this a coordinate var that has already been created in
2378  * the HDF5 file as a dimscale dataset? Check for dims with the
2379  * same name in this group. If there is one, check to see if
2380  * this object exists in the HDF group. */
2381  if (var->became_coord_var)
2382  {
2383  NC_DIM_INFO_T *d1;
2384 
2385  for (d1 = grp->dim; d1; d1 = d1->l.next)
2386  if (!strcmp(d1->name, var->name))
2387  {
2388  nc_bool_t exists;
2389 
2390  if ((retval = var_exists(grp->hdf_grpid, var->name, &exists)))
2391  return retval;
2392  if (exists)
2393  {
2394  /* Indicate that the variable already exists, and should be replaced */
2395  replace_existing_var = NC_TRUE;
2396  flag_atts_dirty(&var->att);
2397  break;
2398  }
2399  }
2400  }
2401 
2402  /* Check dims if the variable will be replaced, so that the dimensions
2403  * will be de-attached and re-attached correctly. */
2404  /* (Note: There's a temptation to merge this loop over the dimensions with
2405  * the prior loop over dimensions, but that blurs the line over the
2406  * purpose of them, so they are currently separate. If performance
2407  * becomes an issue here, it would be possible to merge them. -QAK)
2408  */
2409  if (replace_existing_var)
2410  {
2411  NC_DIM_INFO_T *d1;
2412 
2413  for (d1 = grp->dim; d1; d1 = d1->l.next)
2414  if (!strcmp(d1->name, var->name))
2415  {
2416  nc_bool_t exists;
2417 
2418  if ((retval = var_exists(grp->hdf_grpid, var->name, &exists)))
2419  return retval;
2420  if (exists)
2421  {
2422  hid_t dim_datasetid; /* Dataset ID for dimension */
2423 
2424  /* Find dataset ID for dimension */
2425  if (d1->coord_var)
2426  dim_datasetid = d1->coord_var->hdf_datasetid;
2427  else
2428  dim_datasetid = d1->hdf_dimscaleid;
2429  assert(dim_datasetid > 0);
2430 
2431  /* If we're replacing an existing dimscale dataset, go to
2432  * every var in the file and detach this dimension scale,
2433  * because we have to delete it. */
2434  if ((retval = rec_detach_scales(grp->nc4_info->root_grp,
2435  var->dimids[0], dim_datasetid)))
2436  return retval;
2437  break;
2438  }
2439  }
2440  }
2441 
2442  /* If this is not a dimension scale, do this stuff. */
2443  if (var->was_coord_var && var->dimscale_attached)
2444  {
2445  /* If the variable already exists in the file, Remove any dimension scale
2446  * attributes from it, if they exist. */
2447  /* (The HDF5 Dimension Scale API should really have an API routine
2448  * for making a dataset not a scale. -QAK) */
2449  if (var->created)
2450  {
2451  if ((retval = remove_coord_atts(var->hdf_datasetid)))
2452  BAIL(retval);
2453  }
2454 
2455  if (var->dimscale_attached)
2456  {
2457  int d;
2458 
2459  /* If this is a regular var, detach all its dim scales. */
2460  for (d = 0; d < var->ndims; d++)
2461  if (var->dimscale_attached[d])
2462  {
2463  hid_t dim_datasetid; /* Dataset ID for dimension */
2464  NC_DIM_INFO_T *dim1 = var->dim[d];
2465  assert(dim1 && dim1->dimid == var->dimids[d]);
2466 
2467  /* Find dataset ID for dimension */
2468  if (dim1->coord_var)
2469  dim_datasetid = dim1->coord_var->hdf_datasetid;
2470  else
2471  dim_datasetid = dim1->hdf_dimscaleid;
2472  assert(dim_datasetid > 0);
2473 
2474  if (H5DSdetach_scale(var->hdf_datasetid, dim_datasetid, d) < 0)
2475  BAIL(NC_EHDFERR);
2476  var->dimscale_attached[d] = NC_FALSE;
2477  }
2478  }
2479  }
2480 
2481  /* Delete the HDF5 dataset that is to be replaced. */
2482  if (replace_existing_var)
2483  {
2484  /* Free the HDF5 dataset id. */
2485  if (var->hdf_datasetid && H5Dclose(var->hdf_datasetid) < 0)
2486  BAIL(NC_EHDFERR);
2487  var->hdf_datasetid = 0;
2488 
2489  /* Now delete the variable. */
2490  if (H5Gunlink(grp->hdf_grpid, var->name) < 0)
2491  return NC_EDIMMETA;
2492  }
2493 
2494  /* Create the dataset. */
2495  if (var->is_new_var || replace_existing_var)
2496  {
2497  if ((retval = var_create_dataset(grp, var, write_dimid)))
2498  return retval;
2499  }
2500  else
2501  {
2502  if (write_dimid && var->ndims)
2503  if ((retval = write_netcdf4_dimid(var->hdf_datasetid, var->dimids[0])))
2504  BAIL(retval);
2505  }
2506 
2507  if (replace_existing_var)
2508  {
2509  /* If this is a dimension scale, reattach the scale everywhere it
2510  * is used. (Recall that netCDF dimscales are always 1-D). */
2511  if(var->dimscale)
2512  {
2513  if ((retval = rec_reattach_scales(grp->nc4_info->root_grp,
2514  var->dimids[0], var->hdf_datasetid)))
2515  return retval;
2516  }
2517  /* If it's not a dimension scale, clear the dimscale attached flags,
2518  * so the dimensions are re-attached. */
2519  else
2520  {
2521  if (var->dimscale_attached)
2522  memset(var->dimscale_attached, 0, sizeof(nc_bool_t) * var->ndims);
2523  }
2524  }
2525 
2526  /* Clear coord. var state transition flags */
2527  var->was_coord_var = NC_FALSE;
2528  var->became_coord_var = NC_FALSE;
2529 
2530  /* Now check the attributes for this var. */
2531  if (var->attr_dirty)
2532  {
2533  /* Write attributes for this var. */
2534  if ((retval = write_attlist(var->att, var->varid, grp)))
2535  BAIL(retval);
2536  var->attr_dirty = NC_FALSE;
2537  }
2538 
2539  return NC_NOERR;
2540 exit:
2541  return retval;
2542 }
2543 
2556 static int
2557 write_dim(NC_DIM_INFO_T *dim, NC_GRP_INFO_T *grp, nc_bool_t write_dimid)
2558 {
2559  int retval;
2560  int i;
2561 
2562  /* If there's no dimscale dataset for this dim, create one,
2563  * and mark that it should be hidden from netCDF as a
2564  * variable. (That is, it should appear as a dimension
2565  * without an associated variable.) */
2566  if (0 == dim->hdf_dimscaleid)
2567  {
2568  hid_t spaceid, create_propid;
2569  hsize_t dims[1], max_dims[1], chunk_dims[1] = {1};
2570  char dimscale_wo_var[NC_MAX_NAME];
2571 
2572  LOG((4, "%s: creating dim %s", __func__, dim->name));
2573 
2574  /* Sanity check */
2575  assert(NULL == dim->coord_var);
2576 
2577  /* Create a property list. If this dimension scale is
2578  * unlimited (i.e. it's an unlimited dimension), then set
2579  * up chunking, with a chunksize of 1. */
2580  if ((create_propid = H5Pcreate(H5P_DATASET_CREATE)) < 0)
2581  BAIL(NC_EHDFERR);
2582 
2583  /* RJ: this suppose to be FALSE that is defined in H5 private.h as 0 */
2584  if (H5Pset_obj_track_times(create_propid,0)<0)
2585  BAIL(NC_EHDFERR);
2586 
2587  dims[0] = dim->len;
2588  max_dims[0] = dim->len;
2589  if (dim->unlimited)
2590  {
2591  max_dims[0] = H5S_UNLIMITED;
2592  if (H5Pset_chunk(create_propid, 1, chunk_dims) < 0)
2593  BAIL(NC_EHDFERR);
2594  }
2595 
2596  /* Set up space. */
2597  if ((spaceid = H5Screate_simple(1, dims, max_dims)) < 0)
2598  BAIL(NC_EHDFERR);
2599 
2600  if (H5Pset_attr_creation_order(create_propid, H5P_CRT_ORDER_TRACKED|
2601  H5P_CRT_ORDER_INDEXED) < 0)
2602  BAIL(NC_EHDFERR);
2603 
2604  /* Create the dataset that will be the dimension scale. */
2605  LOG((4, "%s: about to H5Dcreate1 a dimscale dataset %s", __func__, dim->name));
2606  if ((dim->hdf_dimscaleid = H5Dcreate1(grp->hdf_grpid, dim->name, H5T_IEEE_F32BE,
2607  spaceid, create_propid)) < 0)
2608  BAIL(NC_EHDFERR);
2609 
2610  /* Close the spaceid and create_propid. */
2611  if (H5Sclose(spaceid) < 0)
2612  BAIL(NC_EHDFERR);
2613  if (H5Pclose(create_propid) < 0)
2614  BAIL(NC_EHDFERR);
2615 
2616  /* Indicate that this is a scale. Also indicate that not
2617  * be shown to the user as a variable. It is hidden. It is
2618  * a DIM WITHOUT A VARIABLE! */
2619  sprintf(dimscale_wo_var, "%s%10d", DIM_WITHOUT_VARIABLE, (int)dim->len);
2620  if (H5DSset_scale(dim->hdf_dimscaleid, dimscale_wo_var) < 0)
2621  BAIL(NC_EHDFERR);
2622  }
2623 
2624  /* Did we extend an unlimited dimension? */
2625  if (dim->extended)
2626  {
2627  NC_VAR_INFO_T *v1 = NULL;
2628 
2629  assert(dim->unlimited);
2630  /* If this is a dimension without a variable, then update
2631  * the secret length information at the end of the NAME
2632  * attribute. */
2633  for (i=0; i < grp->vars.nelems; i++)
2634  {
2635  if (grp->vars.value[i] && !strcmp(grp->vars.value[i]->name, dim->name))
2636  {
2637  v1 = grp->vars.value[i];
2638  break;
2639  }
2640  }
2641  if (v1)
2642  {
2643  hsize_t *new_size = NULL;
2644  int d1;
2645 
2646  /* Extend the dimension scale dataset to reflect the new
2647  * length of the dimension. */
2648  if (!(new_size = malloc(v1->ndims * sizeof(hsize_t))))
2649  BAIL(NC_ENOMEM);
2650  for (d1 = 0; d1 < v1->ndims; d1++)
2651  {
2652  assert(v1->dim[d1] && v1->dim[d1]->dimid == v1->dimids[d1]);
2653  new_size[d1] = v1->dim[d1]->len;
2654  }
2655  if (H5Dset_extent(v1->hdf_datasetid, new_size) < 0) {
2656  free(new_size);
2657  BAIL(NC_EHDFERR);
2658  }
2659  free(new_size);
2660  }
2661  }
2662 
2663  /* If desired, write the secret dimid. This will be used instead of
2664  * the dimid that the dimension would otherwise receive based on
2665  * creation order. This can be necessary when dims and their
2666  * coordinate variables were created in different order. */
2667  if (write_dimid && dim->hdf_dimscaleid)
2668  if ((retval = write_netcdf4_dimid(dim->hdf_dimscaleid, dim->dimid)))
2669  BAIL(retval);
2670 
2671  return NC_NOERR;
2672 exit:
2673 
2674  return retval;
2675 }
2676 
2693 int
2694 nc4_rec_detect_need_to_preserve_dimids(NC_GRP_INFO_T *grp, nc_bool_t *bad_coord_orderp)
2695 {
2696  NC_VAR_INFO_T *var;
2697  NC_GRP_INFO_T *child_grp;
2698  int last_dimid = -1;
2699  int retval;
2700  int i;
2701 
2702  /* Iterate over variables in this group */
2703  for (i=0; i < grp->vars.nelems; i++)
2704  {
2705  var = grp->vars.value[i];
2706  if (!var) continue;
2707  /* Only matters for dimension scale variables, with non-scalar dimensionality */
2708  if (var->dimscale && var->ndims)
2709  {
2710  /* If the user writes coord vars in a different order then he
2711  * defined their dimensions, then, when the file is reopened, the
2712  * order of the dimids will change to match the order of the coord
2713  * vars. Detect if this is about to happen. */
2714  if (var->dimids[0] < last_dimid)
2715  {
2716  LOG((5, "%s: %s is out of order coord var", __func__, var->name));
2717  *bad_coord_orderp = NC_TRUE;
2718  return NC_NOERR;
2719  }
2720  last_dimid = var->dimids[0];
2721 
2722  /* If there are multidimensional coordinate variables defined, then
2723  * it's also necessary to preserve dimension IDs when the file is
2724  * reopened ... */
2725  if (var->ndims > 1)
2726  {
2727  LOG((5, "%s: %s is multidimensional coord var", __func__, var->name));
2728  *bad_coord_orderp = NC_TRUE;
2729  return NC_NOERR;
2730  }
2731 
2732  /* Did the user define a dimension, end define mode, reenter define
2733  * mode, and then define a coordinate variable for that dimension?
2734  * If so, dimensions will be out of order. */
2735  if (var->is_new_var || var->became_coord_var)
2736  {
2737  LOG((5, "%s: coord var defined after enddef/redef", __func__));
2738  *bad_coord_orderp = NC_TRUE;
2739  return NC_NOERR;
2740  }
2741  }
2742  }
2743 
2744  /* If there are any child groups, check them also for this condition. */
2745  for (child_grp = grp->children; child_grp; child_grp = child_grp->l.next)
2746  if ((retval = nc4_rec_detect_need_to_preserve_dimids(child_grp, bad_coord_orderp)))
2747  return retval;
2748 
2749  return NC_NOERR;
2750 }
2751 
2764 int
2765 nc4_rec_write_metadata(NC_GRP_INFO_T *grp, nc_bool_t bad_coord_order)
2766 {
2767  NC_DIM_INFO_T *dim = NULL;
2768  NC_VAR_INFO_T *var = NULL;
2769  NC_GRP_INFO_T *child_grp = NULL;
2770  int coord_varid = -1;
2771  int var_index = 0;
2772 
2773  int retval;
2774  assert(grp && grp->name && grp->hdf_grpid);
2775  LOG((3, "%s: grp->name %s, bad_coord_order %d", __func__, grp->name, bad_coord_order));
2776 
2777  /* Write global attributes for this group. */
2778  if ((retval = write_attlist(grp->att, NC_GLOBAL, grp)))
2779  return retval;
2780  /* Set the pointers to the beginning of the list of dims & vars in this
2781  * group. */
2782  dim = grp->dim;
2783  if (var_index < grp->vars.nelems)
2784  var = grp->vars.value[var_index];
2785 
2786  /* Because of HDF5 ordering the dims and vars have to be stored in
2787  * this way to ensure that the dims and coordinate vars come out in
2788  * the correct order. */
2789  while (dim || var)
2790  {
2791  nc_bool_t found_coord, wrote_coord;
2792 
2793  /* Write non-coord dims in order, stopping at the first one that
2794  * has an associated coord var. */
2795  for (found_coord = NC_FALSE; dim && !found_coord; dim = dim->l.next)
2796  {
2797  if (!dim->coord_var)
2798  {
2799  if ((retval = write_dim(dim, grp, bad_coord_order)))
2800  return retval;
2801  }
2802  else
2803  {
2804  coord_varid = dim->coord_var->varid;
2805  found_coord = NC_TRUE;
2806  }
2807  }
2808 
2809  /* Write each var. When we get to the coord var we are waiting
2810  * for (if any), then we break after writing it. */
2811  for (wrote_coord = NC_FALSE; var && !wrote_coord; )
2812  {
2813  if ((retval = write_var(var, grp, bad_coord_order)))
2814  return retval;
2815  if (found_coord && var->varid == coord_varid)
2816  wrote_coord = NC_TRUE;
2817  if (++var_index < grp->vars.nelems)
2818  var = grp->vars.value[var_index];
2819  else
2820  var = NULL;
2821  }
2822  } /* end while */
2823 
2824  if ((retval = attach_dimscales(grp)))
2825  return retval;
2826 
2827  /* If there are any child groups, write their metadata. */
2828  for (child_grp = grp->children; child_grp; child_grp = child_grp->l.next)
2829  if ((retval = nc4_rec_write_metadata(child_grp, bad_coord_order)))
2830  return retval;
2831 
2832  return NC_NOERR;
2833 }
2834 
2844 int
2845 nc4_rec_write_groups_types(NC_GRP_INFO_T *grp)
2846 {
2847  NC_GRP_INFO_T *child_grp;
2848  NC_TYPE_INFO_T *type;
2849  int retval;
2850 
2851  assert(grp && grp->name);
2852  LOG((3, "%s: grp->name %s", __func__, grp->name));
2853 
2854  /* Create the group in the HDF5 file if it doesn't exist. */
2855  if (!grp->hdf_grpid)
2856  if ((retval = create_group(grp)))
2857  return retval;
2858 
2859  /* If this is the root group of a file with strict NC3 rules, write
2860  * an attribute. But don't leave the attribute open. */
2861  if (!grp->parent && (grp->nc4_info->cmode & NC_CLASSIC_MODEL))
2862  if ((retval = write_nc3_strict_att(grp->hdf_grpid)))
2863  return retval;
2864 
2865  /* If there are any user-defined types, write them now. */
2866  for (type = grp->type; type; type = type->l.next)
2867  if ((retval = commit_type(grp, type)))
2868  return retval;
2869 
2870  /* If there are any child groups, write their groups and types. */
2871  for (child_grp = grp->children; child_grp; child_grp = child_grp->l.next)
2872  if ((retval = nc4_rec_write_groups_types(child_grp)))
2873  return retval;
2874 
2875  return NC_NOERR;
2876 }
2877 
2905 int
2906 nc4_convert_type(const void *src, void *dest,
2907  const nc_type src_type, const nc_type dest_type,
2908  const size_t len, int *range_error,
2909  const void *fill_value, int strict_nc3, int src_long,
2910  int dest_long)
2911 {
2912  char *cp, *cp1;
2913  float *fp, *fp1;
2914  double *dp, *dp1;
2915  int *ip, *ip1;
2916  signed long *lp, *lp1;
2917  short *sp, *sp1;
2918  signed char *bp, *bp1;
2919  unsigned char *ubp, *ubp1;
2920  unsigned short *usp, *usp1;
2921  unsigned int *uip, *uip1;
2922  long long *lip, *lip1;
2923  unsigned long long *ulip, *ulip1;
2924  size_t count = 0;
2925 
2926  *range_error = 0;
2927  LOG((3, "%s: len %d src_type %d dest_type %d src_long %d dest_long %d",
2928  __func__, len, src_type, dest_type, src_long, dest_long));
2929 
2930  /* OK, this is ugly. If you can think of anything better, I'm open
2931  to suggestions!
2932 
2933  Note that we don't use a default fill value for type
2934  NC_BYTE. This is because Lord Voldemort cast a nofilleramous spell
2935  at Harry Potter, but it bounced off his scar and hit the netcdf-4
2936  code.
2937  */
2938  switch (src_type)
2939  {
2940  case NC_CHAR:
2941  switch (dest_type)
2942  {
2943  case NC_CHAR:
2944  for (cp = (char *)src, cp1 = dest; count < len; count++)
2945  *cp1++ = *cp++;
2946  break;
2947  default:
2948  LOG((0, "%s: Uknown destination type.", __func__));
2949  }
2950  break;
2951 
2952  case NC_BYTE:
2953  switch (dest_type)
2954  {
2955  case NC_BYTE:
2956  for (bp = (signed char *)src, bp1 = dest; count < len; count++)
2957  *bp1++ = *bp++;
2958  break;
2959  case NC_UBYTE:
2960  for (bp = (signed char *)src, ubp = dest; count < len; count++)
2961  {
2962  if (*bp < 0)
2963  (*range_error)++;
2964  *ubp++ = *bp++;
2965  }
2966  break;
2967  case NC_SHORT:
2968  for (bp = (signed char *)src, sp = dest; count < len; count++)
2969  *sp++ = *bp++;
2970  break;
2971  case NC_USHORT:
2972  for (bp = (signed char *)src, usp = dest; count < len; count++)
2973  {
2974  if (*bp < 0)
2975  (*range_error)++;
2976  *usp++ = *bp++;
2977  }
2978  break;
2979  case NC_INT:
2980  if (dest_long)
2981  {
2982  for (bp = (signed char *)src, lp = dest; count < len; count++)
2983  *lp++ = *bp++;
2984  break;
2985  }
2986  else
2987  {
2988  for (bp = (signed char *)src, ip = dest; count < len; count++)
2989  *ip++ = *bp++;
2990  break;
2991  }
2992  case NC_UINT:
2993  for (bp = (signed char *)src, uip = dest; count < len; count++)
2994  {
2995  if (*bp < 0)
2996  (*range_error)++;
2997  *uip++ = *bp++;
2998  }
2999  break;
3000  case NC_INT64:
3001  for (bp = (signed char *)src, lip = dest; count < len; count++)
3002  *lip++ = *bp++;
3003  break;
3004  case NC_UINT64:
3005  for (bp = (signed char *)src, ulip = dest; count < len; count++)
3006  {
3007  if (*bp < 0)
3008  (*range_error)++;
3009  *ulip++ = *bp++;
3010  }
3011  break;
3012  case NC_FLOAT:
3013  for (bp = (signed char *)src, fp = dest; count < len; count++)
3014  *fp++ = *bp++;
3015  break;
3016  case NC_DOUBLE:
3017  for (bp = (signed char *)src, dp = dest; count < len; count++)
3018  *dp++ = *bp++;
3019  break;
3020  default:
3021  LOG((0, "%s: unexpected dest type. src_type %d, dest_type %d",
3022  __func__, src_type, dest_type));
3023  return NC_EBADTYPE;
3024  }
3025  break;
3026 
3027  case NC_UBYTE:
3028  switch (dest_type)
3029  {
3030  case NC_BYTE:
3031  for (ubp = (unsigned char *)src, bp = dest; count < len; count++)
3032  {
3033  if (!strict_nc3 && *ubp > X_SCHAR_MAX)
3034  (*range_error)++;
3035  *bp++ = *ubp++;
3036  }
3037  break;
3038  case NC_SHORT:
3039  for (ubp = (unsigned char *)src, sp = dest; count < len; count++)
3040  *sp++ = *ubp++;
3041  break;
3042  case NC_UBYTE:
3043  for (ubp = (unsigned char *)src, ubp1 = dest; count < len; count++)
3044  *ubp1++ = *ubp++;
3045  break;
3046  case NC_USHORT:
3047  for (ubp = (unsigned char *)src, usp = dest; count < len; count++)
3048  *usp++ = *ubp++;
3049  break;
3050  case NC_INT:
3051  if (dest_long)
3052  {
3053  for (ubp = (unsigned char *)src, lp = dest; count < len; count++)
3054  *lp++ = *ubp++;
3055  break;
3056  }
3057  else
3058  {
3059  for (ubp = (unsigned char *)src, ip = dest; count < len; count++)
3060  *ip++ = *ubp++;
3061  break;
3062  }
3063  case NC_UINT:
3064  for (ubp = (unsigned char *)src, uip = dest; count < len; count++)
3065  *uip++ = *ubp++;
3066  break;
3067  case NC_INT64:
3068  for (ubp = (unsigned char *)src, lip = dest; count < len; count++)
3069  *lip++ = *ubp++;
3070  break;
3071  case NC_UINT64:
3072  for (ubp = (unsigned char *)src, ulip = dest; count < len; count++)
3073  *ulip++ = *ubp++;
3074  break;
3075  case NC_FLOAT:
3076  for (ubp = (unsigned char *)src, fp = dest; count < len; count++)
3077  *fp++ = *ubp++;
3078  break;
3079  case NC_DOUBLE:
3080  for (ubp = (unsigned char *)src, dp = dest; count < len; count++)
3081  *dp++ = *ubp++;
3082  break;
3083  default:
3084  LOG((0, "%s: unexpected dest type. src_type %d, dest_type %d",
3085  __func__, src_type, dest_type));
3086  return NC_EBADTYPE;
3087  }
3088  break;
3089 
3090  case NC_SHORT:
3091  switch (dest_type)
3092  {
3093  case NC_UBYTE:
3094  for (sp = (short *)src, ubp = dest; count < len; count++)
3095  {
3096  if (*sp > X_UCHAR_MAX || *sp < 0)
3097  (*range_error)++;
3098  *ubp++ = *sp++;
3099  }
3100  break;
3101  case NC_BYTE:
3102  for (sp = (short *)src, bp = dest; count < len; count++)
3103  {
3104  if (*sp > X_SCHAR_MAX || *sp < X_SCHAR_MIN)
3105  (*range_error)++;
3106  *bp++ = *sp++;
3107  }
3108  break;
3109  case NC_SHORT:
3110  for (sp = (short *)src, sp1 = dest; count < len; count++)
3111  *sp1++ = *sp++;
3112  break;
3113  case NC_USHORT:
3114  for (sp = (short *)src, usp = dest; count < len; count++)
3115  {
3116  if (*sp < 0)
3117  (*range_error)++;
3118  *usp++ = *sp++;
3119  }
3120  break;
3121  case NC_INT:
3122  if (dest_long)
3123  for (sp = (short *)src, lp = dest; count < len; count++)
3124  *lp++ = *sp++;
3125  else
3126  for (sp = (short *)src, ip = dest; count < len; count++)
3127  *ip++ = *sp++;
3128  break;
3129  case NC_UINT:
3130  for (sp = (short *)src, uip = dest; count < len; count++)
3131  {
3132  if (*sp < 0)
3133  (*range_error)++;
3134  *uip++ = *sp++;
3135  }
3136  break;
3137  case NC_INT64:
3138  for (sp = (short *)src, lip = dest; count < len; count++)
3139  *lip++ = *sp++;
3140  break;
3141  case NC_UINT64:
3142  for (sp = (short *)src, ulip = dest; count < len; count++)
3143  {
3144  if (*sp < 0)
3145  (*range_error)++;
3146  *ulip++ = *sp++;
3147  }
3148  break;
3149  case NC_FLOAT:
3150  for (sp = (short *)src, fp = dest; count < len; count++)
3151  *fp++ = *sp++;
3152  break;
3153  case NC_DOUBLE:
3154  for (sp = (short *)src, dp = dest; count < len; count++)
3155  *dp++ = *sp++;
3156  break;
3157  default:
3158  LOG((0, "%s: unexpected dest type. src_type %d, dest_type %d",
3159  __func__, src_type, dest_type));
3160  return NC_EBADTYPE;
3161  }
3162  break;
3163 
3164  case NC_USHORT:
3165  switch (dest_type)
3166  {
3167  case NC_UBYTE:
3168  for (usp = (unsigned short *)src, ubp = dest; count < len; count++)
3169  {
3170  if (*usp > X_UCHAR_MAX)
3171  (*range_error)++;
3172  *ubp++ = *usp++;
3173  }
3174  break;
3175  case NC_BYTE:
3176  for (usp = (unsigned short *)src, bp = dest; count < len; count++)
3177  {
3178  if (*usp > X_SCHAR_MAX)
3179  (*range_error)++;
3180  *bp++ = *usp++;
3181  }
3182  break;
3183  case NC_SHORT:
3184  for (usp = (unsigned short *)src, sp = dest; count < len; count++)
3185  {
3186  if (*usp > X_SHORT_MAX)
3187  (*range_error)++;
3188  *sp++ = *usp++;
3189  }
3190  break;
3191  case NC_USHORT:
3192  for (usp = (unsigned short *)src, usp1 = dest; count < len; count++)
3193  *usp1++ = *usp++;
3194  break;
3195  case NC_INT:
3196  if (dest_long)
3197  for (usp = (unsigned short *)src, lp = dest; count < len; count++)
3198  *lp++ = *usp++;
3199  else
3200  for (usp = (unsigned short *)src, ip = dest; count < len; count++)
3201  *ip++ = *usp++;
3202  break;
3203  case NC_UINT:
3204  for (usp = (unsigned short *)src, uip = dest; count < len; count++)
3205  *uip++ = *usp++;
3206  break;
3207  case NC_INT64:
3208  for (usp = (unsigned short *)src, lip = dest; count < len; count++)
3209  *lip++ = *usp++;
3210  break;
3211  case NC_UINT64:
3212  for (usp = (unsigned short *)src, ulip = dest; count < len; count++)
3213  *ulip++ = *usp++;
3214  break;
3215  case NC_FLOAT:
3216  for (usp = (unsigned short *)src, fp = dest; count < len; count++)
3217  *fp++ = *usp++;
3218  break;
3219  case NC_DOUBLE:
3220  for (usp = (unsigned short *)src, dp = dest; count < len; count++)
3221  *dp++ = *usp++;
3222  break;
3223  default:
3224  LOG((0, "%s: unexpected dest type. src_type %d, dest_type %d",
3225  __func__, src_type, dest_type));
3226  return NC_EBADTYPE;
3227  }
3228  break;
3229 
3230  case NC_INT:
3231  if (src_long)
3232  {
3233  switch (dest_type)
3234  {
3235  case NC_UBYTE:
3236  for (lp = (long *)src, ubp = dest; count < len; count++)
3237  {
3238  if (*lp > X_UCHAR_MAX || *lp < 0)
3239  (*range_error)++;
3240  *ubp++ = *lp++;
3241  }
3242  break;
3243  case NC_BYTE:
3244  for (lp = (long *)src, bp = dest; count < len; count++)
3245  {
3246  if (*lp > X_SCHAR_MAX || *lp < X_SCHAR_MIN)
3247  (*range_error)++;
3248  *bp++ = *lp++;
3249  }
3250  break;
3251  case NC_SHORT:
3252  for (lp = (long *)src, sp = dest; count < len; count++)
3253  {
3254  if (*lp > X_SHORT_MAX || *lp < X_SHORT_MIN)
3255  (*range_error)++;
3256  *sp++ = *lp++;
3257  }
3258  break;
3259  case NC_USHORT:
3260  for (lp = (long *)src, usp = dest; count < len; count++)
3261  {
3262  if (*lp > X_USHORT_MAX || *lp < 0)
3263  (*range_error)++;
3264  *usp++ = *lp++;
3265  }
3266  break;
3267  case NC_INT: /* src is long */
3268  if (dest_long)
3269  {
3270  for (lp = (long *)src, lp1 = dest; count < len; count++)
3271  {
3272  if (*lp > X_LONG_MAX || *lp < X_LONG_MIN)
3273  (*range_error)++;
3274  *lp1++ = *lp++;
3275  }
3276  }
3277  else /* dest is int */
3278  {
3279  for (lp = (long *)src, ip = dest; count < len; count++)
3280  {
3281  if (*lp > X_INT_MAX || *lp < X_INT_MIN)
3282  (*range_error)++;
3283  *ip++ = *lp++;
3284  }
3285  }
3286  break;
3287  case NC_UINT:
3288  for (lp = (long *)src, uip = dest; count < len; count++)
3289  {
3290  if (*lp > X_UINT_MAX || *lp < 0)
3291  (*range_error)++;
3292  *uip++ = *lp++;
3293  }
3294  break;
3295  case NC_INT64:
3296  for (lp = (long *)src, lip = dest; count < len; count++)
3297  *lip++ = *lp++;
3298  break;
3299  case NC_UINT64:
3300  for (lp = (long *)src, ulip = dest; count < len; count++)
3301  {
3302  if (*lp < 0)
3303  (*range_error)++;
3304  *ulip++ = *lp++;
3305  }
3306  break;
3307  case NC_FLOAT:
3308  for (lp = (long *)src, fp = dest; count < len; count++)
3309  *fp++ = *lp++;
3310  break;
3311  case NC_DOUBLE:
3312  for (lp = (long *)src, dp = dest; count < len; count++)
3313  *dp++ = *lp++;
3314  break;
3315  default:
3316  LOG((0, "%s: unexpected dest type. src_type %d, dest_type %d",
3317  __func__, src_type, dest_type));
3318  return NC_EBADTYPE;
3319  }
3320  }
3321  else
3322  {
3323  switch (dest_type)
3324  {
3325  case NC_UBYTE:
3326  for (ip = (int *)src, ubp = dest; count < len; count++)
3327  {
3328  if (*ip > X_UCHAR_MAX || *ip < 0)
3329  (*range_error)++;
3330  *ubp++ = *ip++;
3331  }
3332  break;
3333  case NC_BYTE:
3334  for (ip = (int *)src, bp = dest; count < len; count++)
3335  {
3336  if (*ip > X_SCHAR_MAX || *ip < X_SCHAR_MIN)
3337  (*range_error)++;
3338  *bp++ = *ip++;
3339  }
3340  break;
3341  case NC_SHORT:
3342  for (ip = (int *)src, sp = dest; count < len; count++)
3343  {
3344  if (*ip > X_SHORT_MAX || *ip < X_SHORT_MIN)
3345  (*range_error)++;
3346  *sp++ = *ip++;
3347  }
3348  break;
3349  case NC_USHORT:
3350  for (ip = (int *)src, usp = dest; count < len; count++)
3351  {
3352  if (*ip > X_USHORT_MAX || *ip < 0)
3353  (*range_error)++;
3354  *usp++ = *ip++;
3355  }
3356  break;
3357  case NC_INT: /* src is int */
3358  if (dest_long)
3359  {
3360  for (ip = (int *)src, lp1 = dest; count < len; count++)
3361  {
3362  if (*ip > X_LONG_MAX || *ip < X_LONG_MIN)
3363  (*range_error)++;
3364  *lp1++ = *ip++;
3365  }
3366  }
3367  else /* dest is int */
3368  {
3369  for (ip = (int *)src, ip1 = dest; count < len; count++)
3370  {
3371  if (*ip > X_INT_MAX || *ip < X_INT_MIN)
3372  (*range_error)++;
3373  *ip1++ = *ip++;
3374  }
3375  }
3376  break;
3377  case NC_UINT:
3378  for (ip = (int *)src, uip = dest; count < len; count++)
3379  {
3380  if (*ip > X_UINT_MAX || *ip < 0)
3381  (*range_error)++;
3382  *uip++ = *ip++;
3383  }
3384  break;
3385  case NC_INT64:
3386  for (ip = (int *)src, lip = dest; count < len; count++)
3387  *lip++ = *ip++;
3388  break;
3389  case NC_UINT64:
3390  for (ip = (int *)src, ulip = dest; count < len; count++)
3391  {
3392  if (*ip < 0)
3393  (*range_error)++;
3394  *ulip++ = *ip++;
3395  }
3396  break;
3397  case NC_FLOAT:
3398  for (ip = (int *)src, fp = dest; count < len; count++)
3399  *fp++ = *ip++;
3400  break;
3401  case NC_DOUBLE:
3402  for (ip = (int *)src, dp = dest; count < len; count++)
3403  *dp++ = *ip++;
3404  break;
3405  default:
3406  LOG((0, "%s: unexpected dest type. src_type %d, dest_type %d",
3407  __func__, src_type, dest_type));
3408  return NC_EBADTYPE;
3409  }
3410  }
3411  break;
3412 
3413  case NC_UINT:
3414  switch (dest_type)
3415  {
3416  case NC_UBYTE:
3417  for (uip = (unsigned int *)src, ubp = dest; count < len; count++)
3418  {
3419  if (*uip > X_UCHAR_MAX)
3420  (*range_error)++;
3421  *ubp++ = *uip++;
3422  }
3423  break;
3424  case NC_BYTE:
3425  for (uip = (unsigned int *)src, bp = dest; count < len; count++)
3426  {
3427  if (*uip > X_SCHAR_MAX)
3428  (*range_error)++;
3429  *bp++ = *uip++;
3430  }
3431  break;
3432  case NC_SHORT:
3433  for (uip = (unsigned int *)src, sp = dest; count < len; count++)
3434  {
3435  if (*uip > X_SHORT_MAX)
3436  (*range_error)++;
3437  *sp++ = *uip++;
3438  }
3439  break;
3440  case NC_USHORT:
3441  for (uip = (unsigned int *)src, usp = dest; count < len; count++)
3442  {
3443  if (*uip > X_USHORT_MAX)
3444  (*range_error)++;
3445  *usp++ = *uip++;
3446  }
3447  break;
3448  case NC_INT:
3449  if (dest_long)
3450  for (uip = (unsigned int *)src, lp = dest; count < len; count++)
3451  {
3452  if (*uip > X_LONG_MAX)
3453  (*range_error)++;
3454  *lp++ = *uip++;
3455  }
3456  else
3457  for (uip = (unsigned int *)src, ip = dest; count < len; count++)
3458  {
3459  if (*uip > X_INT_MAX)
3460  (*range_error)++;
3461  *ip++ = *uip++;
3462  }
3463  break;
3464  case NC_UINT:
3465  for (uip = (unsigned int *)src, uip1 = dest; count < len; count++)
3466  {
3467  if (*uip > X_UINT_MAX)
3468  (*range_error)++;
3469  *uip1++ = *uip++;
3470  }
3471  break;
3472  case NC_INT64:
3473  for (uip = (unsigned int *)src, lip = dest; count < len; count++)
3474  *lip++ = *uip++;
3475  break;
3476  case NC_UINT64:
3477  for (uip = (unsigned int *)src, ulip = dest; count < len; count++)
3478  *ulip++ = *uip++;
3479  break;
3480  case NC_FLOAT:
3481  for (uip = (unsigned int *)src, fp = dest; count < len; count++)
3482  *fp++ = *uip++;
3483  break;
3484  case NC_DOUBLE:
3485  for (uip = (unsigned int *)src, dp = dest; count < len; count++)
3486  *dp++ = *uip++;
3487  break;
3488  default:
3489  LOG((0, "%s: unexpected dest type. src_type %d, dest_type %d",
3490  __func__, src_type, dest_type));
3491  return NC_EBADTYPE;
3492  }
3493  break;
3494 
3495  case NC_INT64:
3496  switch (dest_type)
3497  {
3498  case NC_UBYTE:
3499  for (lip = (long long *)src, ubp = dest; count < len; count++)
3500  {
3501  if (*lip > X_UCHAR_MAX || *lip < 0)
3502  (*range_error)++;
3503  *ubp++ = *lip++;
3504  }
3505  break;
3506  case NC_BYTE:
3507  for (lip = (long long *)src, bp = dest; count < len; count++)
3508  {
3509  if (*lip > X_SCHAR_MAX || *lip < X_SCHAR_MIN)
3510  (*range_error)++;
3511  *bp++ = *lip++;
3512  }
3513  break;
3514  case NC_SHORT:
3515  for (lip = (long long *)src, sp = dest; count < len; count++)
3516  {
3517  if (*lip > X_SHORT_MAX || *lip < X_SHORT_MIN)
3518  (*range_error)++;
3519  *sp++ = *lip++;
3520  }
3521  break;
3522  case NC_USHORT:
3523  for (lip = (long long *)src, usp = dest; count < len; count++)
3524  {
3525  if (*lip > X_USHORT_MAX || *lip < 0)
3526  (*range_error)++;
3527  *usp++ = *lip++;
3528  }
3529  break;
3530  case NC_UINT:
3531  for (lip = (long long *)src, uip = dest; count < len; count++)
3532  {
3533  if (*lip > X_UINT_MAX || *lip < 0)
3534  (*range_error)++;
3535  *uip++ = *lip++;
3536  }
3537  break;
3538  case NC_INT:
3539  if (dest_long)
3540  for (lip = (long long *)src, lp = dest; count < len; count++)
3541  {
3542  if (*lip > X_LONG_MAX || *lip < X_LONG_MIN)
3543  (*range_error)++;
3544  *lp++ = *lip++;
3545  }
3546  else
3547  for (lip = (long long *)src, ip = dest; count < len; count++)
3548  {
3549  if (*lip > X_INT_MAX || *lip < X_INT_MIN)
3550  (*range_error)++;
3551  *ip++ = *lip++;
3552  }
3553  break;
3554  case NC_INT64:
3555  for (lip = (long long *)src, lip1 = dest; count < len; count++)
3556  *lip1++ = *lip++;
3557  break;
3558  case NC_UINT64:
3559  for (lip = (long long *)src, ulip = dest; count < len; count++)
3560  {
3561  if (*lip < 0)
3562  (*range_error)++;
3563  *ulip++ = *lip++;
3564  }
3565  break;
3566  case NC_FLOAT:
3567  for (lip = (long long *)src, fp = dest; count < len; count++)
3568  *fp++ = *lip++;
3569  break;
3570  case NC_DOUBLE:
3571  for (lip = (long long *)src, dp = dest; count < len; count++)
3572  *dp++ = *lip++;
3573  break;
3574  default:
3575  LOG((0, "%s: unexpected dest type. src_type %d, dest_type %d",
3576  __func__, src_type, dest_type));
3577  return NC_EBADTYPE;
3578  }
3579  break;
3580 
3581  case NC_UINT64:
3582  switch (dest_type)
3583  {
3584  case NC_UBYTE:
3585  for (ulip = (unsigned long long *)src, ubp = dest; count < len; count++)
3586  {
3587  if (*ulip > X_UCHAR_MAX)
3588  (*range_error)++;
3589  *ubp++ = *ulip++;
3590  }
3591  break;
3592  case NC_BYTE:
3593  for (ulip = (unsigned long long *)src, bp = dest; count < len; count++)
3594  {
3595  if (*ulip > X_SCHAR_MAX)
3596  (*range_error)++;
3597  *bp++ = *ulip++;
3598  }
3599  break;
3600  case NC_SHORT:
3601  for (ulip = (unsigned long long *)src, sp = dest; count < len; count++)
3602  {
3603  if (*ulip > X_SHORT_MAX)
3604  (*range_error)++;
3605  *sp++ = *ulip++;
3606  }
3607  break;
3608  case NC_USHORT:
3609  for (ulip = (unsigned long long *)src, usp = dest; count < len; count++)
3610  {
3611  if (*ulip > X_USHORT_MAX)
3612  (*range_error)++;
3613  *usp++ = *ulip++;
3614  }
3615  break;
3616  case NC_UINT:
3617  for (ulip = (unsigned long long *)src, uip = dest; count < len; count++)
3618  {
3619  if (*ulip > X_UINT_MAX)
3620  (*range_error)++;
3621  *uip++ = *ulip++;
3622  }
3623  break;
3624  case NC_INT:
3625  if (dest_long)
3626  for (ulip = (unsigned long long *)src, lp = dest; count < len; count++)
3627  {
3628  if (*ulip > X_LONG_MAX)
3629  (*range_error)++;
3630  *lp++ = *ulip++;
3631  }
3632  else
3633  for (ulip = (unsigned long long *)src, ip = dest; count < len; count++)
3634  {
3635  if (*ulip > X_INT_MAX)
3636  (*range_error)++;
3637  *ip++ = *ulip++;
3638  }
3639  break;
3640  case NC_INT64:
3641  for (ulip = (unsigned long long *)src, lip = dest; count < len; count++)
3642  {
3643  if (*ulip > X_INT64_MAX)
3644  (*range_error)++;
3645  *lip++ = *ulip++;
3646  }
3647  break;
3648  case NC_UINT64:
3649  for (ulip = (unsigned long long *)src, ulip1 = dest; count < len; count++)
3650  *ulip1++ = *ulip++;
3651  break;
3652  case NC_FLOAT:
3653  for (ulip = (unsigned long long *)src, fp = dest; count < len; count++)
3654  *fp++ = *ulip++;
3655  break;
3656  case NC_DOUBLE:
3657  for (ulip = (unsigned long long *)src, dp = dest; count < len; count++)
3658  *dp++ = *ulip++;
3659  break;
3660  default:
3661  LOG((0, "%s: unexpected dest type. src_type %d, dest_type %d",
3662  __func__, src_type, dest_type));
3663  return NC_EBADTYPE;
3664  }
3665  break;
3666 
3667  case NC_FLOAT:
3668  switch (dest_type)
3669  {
3670  case NC_UBYTE:
3671  for (fp = (float *)src, ubp = dest; count < len; count++)
3672  {
3673  if (*fp > X_UCHAR_MAX || *fp < 0)
3674  (*range_error)++;
3675  *ubp++ = *fp++;
3676  }
3677  break;
3678  case NC_BYTE:
3679  for (fp = (float *)src, bp = dest; count < len; count++)
3680  {
3681  if (*fp > (double)X_SCHAR_MAX || *fp < (double)X_SCHAR_MIN)
3682  (*range_error)++;
3683  *bp++ = *fp++;
3684  }
3685  break;
3686  case NC_SHORT:
3687  for (fp = (float *)src, sp = dest; count < len; count++)
3688  {
3689  if (*fp > (double)X_SHORT_MAX || *fp < (double)X_SHORT_MIN)
3690  (*range_error)++;
3691  *sp++ = *fp++;
3692  }
3693  break;
3694  case NC_USHORT:
3695  for (fp = (float *)src, usp = dest; count < len; count++)
3696  {
3697  if (*fp > X_USHORT_MAX || *fp < 0)
3698  (*range_error)++;
3699  *usp++ = *fp++;
3700  }
3701  break;
3702  case NC_UINT:
3703  for (fp = (float *)src, uip = dest; count < len; count++)
3704  {
3705  if (*fp > X_UINT_MAX || *fp < 0)
3706  (*range_error)++;
3707  *uip++ = *fp++;
3708  }
3709  break;
3710  case NC_INT:
3711  if (dest_long)
3712  for (fp = (float *)src, lp = dest; count < len; count++)
3713  {
3714  if (*fp > (double)X_LONG_MAX || *fp < (double)X_LONG_MIN)
3715  (*range_error)++;
3716  *lp++ = *fp++;
3717  }
3718  else
3719  for (fp = (float *)src, ip = dest; count < len; count++)
3720  {
3721  if (*fp > (double)X_INT_MAX || *fp < (double)X_INT_MIN)
3722  (*range_error)++;
3723  *ip++ = *fp++;
3724  }
3725  break;
3726  case NC_INT64:
3727  for (fp = (float *)src, lip = dest; count < len; count++)
3728  {
3729  if (*fp > X_INT64_MAX || *fp <X_INT64_MIN)
3730  (*range_error)++;
3731  *lip++ = *fp++;
3732  }
3733  break;
3734  case NC_UINT64:
3735  for (fp = (float *)src, lip = dest; count < len; count++)
3736  {
3737  if (*fp > X_UINT64_MAX || *fp < 0)
3738  (*range_error)++;
3739  *lip++ = *fp++;
3740  }
3741  break;
3742  case NC_FLOAT:
3743  for (fp = (float *)src, fp1 = dest; count < len; count++)
3744  {
3745  /* if (*fp > X_FLOAT_MAX || *fp < X_FLOAT_MIN)
3746  (*range_error)++;*/
3747  *fp1++ = *fp++;
3748  }
3749  break;
3750  case NC_DOUBLE:
3751  for (fp = (float *)src, dp = dest; count < len; count++)
3752  *dp++ = *fp++;
3753  break;
3754  default:
3755  LOG((0, "%s: unexpected dest type. src_type %d, dest_type %d",
3756  __func__, src_type, dest_type));
3757  return NC_EBADTYPE;
3758  }
3759  break;
3760 
3761  case NC_DOUBLE:
3762  switch (dest_type)
3763  {
3764  case NC_UBYTE:
3765  for (dp = (double *)src, ubp = dest; count < len; count++)
3766  {
3767  if (*dp > X_UCHAR_MAX || *dp < 0)
3768  (*range_error)++;
3769  *ubp++ = *dp++;
3770  }
3771  break;
3772  case NC_BYTE:
3773  for (dp = (double *)src, bp = dest; count < len; count++)
3774  {
3775  if (*dp > X_SCHAR_MAX || *dp < X_SCHAR_MIN)
3776  (*range_error)++;
3777  *bp++ = *dp++;
3778  }
3779  break;
3780  case NC_SHORT:
3781  for (dp = (double *)src, sp = dest; count < len; count++)
3782  {
3783  if (*dp > X_SHORT_MAX || *dp < X_SHORT_MIN)
3784  (*range_error)++;
3785  *sp++ = *dp++;
3786  }
3787  break;
3788  case NC_USHORT:
3789  for (dp = (double *)src, usp = dest; count < len; count++)
3790  {
3791  if (*dp > X_USHORT_MAX || *dp < 0)
3792  (*range_error)++;
3793  *usp++ = *dp++;
3794  }
3795  break;
3796  case NC_UINT:
3797  for (dp = (double *)src, uip = dest; count < len; count++)
3798  {
3799  if (*dp > X_UINT_MAX || *dp < 0)
3800  (*range_error)++;
3801  *uip++ = *dp++;
3802  }
3803  break;
3804  case NC_INT:
3805  if (dest_long)
3806  for (dp = (double *)src, lp = dest; count < len; count++)
3807  {
3808  if (*dp > X_LONG_MAX || *dp < X_LONG_MIN)
3809  (*range_error)++;
3810  *lp++ = *dp++;
3811  }
3812  else
3813  for (dp = (double *)src, ip = dest; count < len; count++)
3814  {
3815  if (*dp > X_INT_MAX || *dp < X_INT_MIN)
3816  (*range_error)++;
3817  *ip++ = *dp++;
3818  }
3819  break;
3820  case NC_INT64:
3821  for (dp = (double *)src, lip = dest; count < len; count++)
3822  {
3823  if (*dp > X_INT64_MAX || *dp < X_INT64_MIN)
3824  (*range_error)++;
3825  *lip++ = *dp++;
3826  }
3827  break;
3828  case NC_UINT64:
3829  for (dp = (double *)src, lip = dest; count < len; count++)
3830  {
3831  if (*dp > X_UINT64_MAX || *dp < 0)
3832  (*range_error)++;
3833  *lip++ = *dp++;
3834  }
3835  break;
3836  case NC_FLOAT:
3837  for (dp = (double *)src, fp = dest; count < len; count++)
3838  {
3839  if (*dp > X_FLOAT_MAX || *dp < X_FLOAT_MIN)
3840  (*range_error)++;
3841  *fp++ = *dp++;
3842  }
3843  break;
3844  case NC_DOUBLE:
3845  for (dp = (double *)src, dp1 = dest; count < len; count++)
3846  {
3847  /* if (*dp > X_DOUBLE_MAX || *dp < X_DOUBLE_MIN) */
3848  /* (*range_error)++; */
3849  *dp1++ = *dp++;
3850  }
3851  break;
3852  default:
3853  LOG((0, "%s: unexpected dest type. src_type %d, dest_type %d",
3854  __func__, src_type, dest_type));
3855  return NC_EBADTYPE;
3856  }
3857  break;
3858 
3859  default:
3860  LOG((0, "%s: unexpected src type. src_type %d, dest_type %d",
3861  __func__, src_type, dest_type));
3862  return NC_EBADTYPE;
3863  }
3864  return NC_NOERR;
3865 }
3866 
3879 int
3880 nc4_rec_match_dimscales(NC_GRP_INFO_T *grp)
3881 {
3882  NC_GRP_INFO_T *g;
3883  NC_VAR_INFO_T *var;
3884  NC_DIM_INFO_T *dim;
3885  int retval = NC_NOERR;
3886  int i;
3887 
3888  assert(grp && grp->name);
3889  LOG((4, "%s: grp->name %s", __func__, grp->name));
3890 
3891  /* Perform var dimscale match for child groups. */
3892  for (g = grp->children; g; g = g->l.next)
3893  if ((retval = nc4_rec_match_dimscales(g)))
3894  return retval;
3895 
3896  /* Check all the vars in this group. If they have dimscale info,
3897  * try and find a dimension for them. */
3898  for (i=0; i < grp->vars.nelems; i++)
3899  {
3900  int ndims;
3901  int d;
3902  var = grp->vars.value[i];
3903  if (!var) continue;
3904  /* Check all vars and see if dim[i] != NULL if dimids[i] valid. */
3905  ndims = var->ndims;
3906  for (d = 0; d < ndims; d++)
3907  {
3908  if (var->dim[d] == NULL) {
3909  nc4_find_dim(grp, var->dimids[d], &var->dim[d], NULL);
3910  }
3911  /* assert(var->dim[d] && var->dim[d]->dimid == var->dimids[d]); */
3912  }
3913 
3914  /* Skip dimension scale variables */
3915  if (!var->dimscale)
3916  {
3917  int d;
3918 
3919  /* Are there dimscales for this variable? */
3920  if (var->dimscale_hdf5_objids)
3921  {
3922  for (d = 0; d < var->ndims; d++)
3923  {
3924  nc_bool_t finished = NC_FALSE;
3925 
3926  LOG((5, "%s: var %s has dimscale info...", __func__, var->name));
3927  /* Look at all the dims in this group to see if they
3928  * match. */
3929  for (g = grp; g && !finished; g = g->parent)
3930  {
3931  for (dim = g->dim; dim; dim = dim->l.next)
3932  {
3933  if (var->dimscale_hdf5_objids[d].fileno[0] == dim->hdf5_objid.fileno[0] &&
3934  var->dimscale_hdf5_objids[d].objno[0] == dim->hdf5_objid.objno[0] &&
3935  var->dimscale_hdf5_objids[d].fileno[1] == dim->hdf5_objid.fileno[1] &&
3936  var->dimscale_hdf5_objids[d].objno[1] == dim->hdf5_objid.objno[1])
3937  {
3938  LOG((4, "%s: for dimension %d, found dim %s",
3939  __func__, d, dim->name));
3940  var->dimids[d] = dim->dimid;
3941  var->dim[d] = dim;
3942  finished = NC_TRUE;
3943  break;
3944  }
3945  } /* next dim */
3946  } /* next grp */
3947  LOG((5, "%s: dimid for this dimscale is %d", __func__, var->type_info->nc_typeid));
3948  } /* next var->dim */
3949  }
3950  /* No dimscales for this var! Invent phony dimensions. */
3951  else
3952  {
3953  hid_t spaceid = 0;
3954  hsize_t *h5dimlen = NULL, *h5dimlenmax = NULL;
3955  int dataset_ndims;
3956 
3957  /* Find the space information for this dimension. */
3958  if ((spaceid = H5Dget_space(var->hdf_datasetid)) < 0)
3959  return NC_EHDFERR;
3960 
3961  /* Get the len of each dim in the space. */
3962  if (var->ndims)
3963  {
3964  if (!(h5dimlen = malloc(var->ndims * sizeof(hsize_t))))
3965  return NC_ENOMEM;
3966  if (!(h5dimlenmax = malloc(var->ndims * sizeof(hsize_t))))
3967  {
3968  free(h5dimlen);
3969  return NC_ENOMEM;
3970  }
3971  if ((dataset_ndims = H5Sget_simple_extent_dims(spaceid, h5dimlen,
3972  h5dimlenmax)) < 0) {
3973  free(h5dimlenmax);
3974  free(h5dimlen);
3975  return NC_EHDFERR;
3976  }
3977  if (dataset_ndims != var->ndims) {
3978  free(h5dimlenmax);
3979  free(h5dimlen);
3980  return NC_EHDFERR;
3981  }
3982  }
3983  else
3984  {
3985  /* Make sure it's scalar. */
3986  if (H5Sget_simple_extent_type(spaceid) != H5S_SCALAR)
3987  return NC_EHDFERR;
3988  }
3989 
3990  /* Release the space object. */
3991  if (H5Sclose(spaceid) < 0) {
3992  free(h5dimlen);
3993  free(h5dimlenmax);
3994  return NC_EHDFERR;
3995  }
3996 
3997  /* Create a phony dimension for each dimension in the
3998  * dataset, unless there already is one the correct
3999  * size. */
4000  for (d = 0; d < var->ndims; d++)
4001  {
4002  /* Is there already a phony dimension of the correct size? */
4003  for (dim = grp->dim; dim; dim = dim->l.next)
4004  if ((dim->len == h5dimlen[d]) &&
4005  ((h5dimlenmax[d] == H5S_UNLIMITED && dim->unlimited) ||
4006  (h5dimlenmax[d] != H5S_UNLIMITED && !dim->unlimited)))
4007  break;
4008 
4009  /* Didn't find a phony dim? Then create one. */
4010  if (!dim)
4011  {
4012  char phony_dim_name[NC_MAX_NAME + 1];
4013 
4014  LOG((3, "%s: creating phony dim for var %s", __func__, var->name));
4015  if ((retval = nc4_dim_list_add(&grp->dim, &dim))) {
4016  free(h5dimlenmax);
4017  free(h5dimlen);
4018  return retval;
4019  }
4020  dim->dimid = grp->nc4_info->next_dimid++;
4021  sprintf(phony_dim_name, "phony_dim_%d", dim->dimid);
4022  if (!(dim->name = strdup(phony_dim_name))) {
4023  free(h5dimlenmax);
4024  free(h5dimlen);
4025  return NC_ENOMEM;
4026  }
4027  dim->len = h5dimlen[d];
4028  dim->hash = hash_fast(phony_dim_name, strlen(phony_dim_name));
4029  if (h5dimlenmax[d] == H5S_UNLIMITED)
4030  dim->unlimited = NC_TRUE;
4031  }
4032 
4033  /* The variable must remember the dimid. */
4034  var->dimids[d] = dim->dimid;
4035  var->dim[d] = dim;
4036  } /* next dim */
4037 
4038  /* Free the memory we malloced. */
4039  free(h5dimlen);
4040  free(h5dimlenmax);
4041  }
4042  }
4043  }
4044 
4045  return retval;
4046 }
4047 
4061 int
4062 nc4_get_typelen_mem(NC_HDF5_FILE_INFO_T *h5, nc_type xtype, int is_long,
4063  size_t *len)
4064 {
4065  NC_TYPE_INFO_T *type;
4066  int retval;
4067 
4068  LOG((4, "%s xtype: %d", __func__, xtype));
4069  assert(len);
4070 
4071  /* If this is an atomic type, the answer is easy. */
4072  switch (xtype)
4073  {
4074  case NC_BYTE:
4075  case NC_CHAR:
4076  case NC_UBYTE:
4077  *len = sizeof(char);
4078  return NC_NOERR;
4079  case NC_SHORT:
4080  case NC_USHORT:
4081  *len = sizeof(short);
4082  return NC_NOERR;
4083  case NC_INT:
4084  case NC_UINT:
4085  if (is_long)
4086  *len = sizeof(long);
4087  else
4088  *len = sizeof(int);
4089  return NC_NOERR;
4090  case NC_FLOAT:
4091  *len = sizeof(float);
4092  return NC_NOERR;
4093  case NC_DOUBLE:
4094  *len = sizeof(double);
4095  return NC_NOERR;
4096  case NC_INT64:
4097  case NC_UINT64:
4098  *len = sizeof(long long);
4099  return NC_NOERR;
4100  case NC_STRING:
4101  *len = sizeof(char *);
4102  return NC_NOERR;
4103  }
4104 
4105  /* See if var is compound type. */
4106  if ((retval = nc4_find_type(h5, xtype, &type)))
4107  return retval;
4108 
4109  if (!type)
4110  return NC_EBADTYPE;
4111 
4112  *len = type->size;
4113 
4114  LOG((5, "type->size: %d", type->size));
4115 
4116  return NC_NOERR;
4117 }
4118 
4131 int
4132 nc4_get_typeclass(const NC_HDF5_FILE_INFO_T *h5, nc_type xtype, int *type_class)
4133 {
4134  int retval = NC_NOERR;
4135 
4136  LOG((4, "%s xtype: %d", __func__, xtype));
4137  assert(type_class);
4138 
4139  /* If this is an atomic type, the answer is easy. */
4140  if (xtype <= NC_STRING)
4141  {
4142  switch (xtype)
4143  {
4144  case NC_BYTE:
4145  case NC_UBYTE:
4146  case NC_SHORT:
4147  case NC_USHORT:
4148  case NC_INT:
4149  case NC_UINT:
4150  case NC_INT64:
4151  case NC_UINT64:
4152  /* NC_INT is class used for all integral types */
4153  *type_class = NC_INT;
4154  break;
4155 
4156  case NC_FLOAT:
4157  case NC_DOUBLE:
4158  /* NC_FLOAT is class used for all floating-point types */
4159  *type_class = NC_FLOAT;
4160  break;
4161 
4162  case NC_CHAR:
4163  *type_class = NC_CHAR;
4164  break;
4165 
4166  case NC_STRING:
4167  *type_class = NC_STRING;
4168  break;
4169 
4170  default:
4171  BAIL(NC_EBADTYPE);
4172  }
4173  }
4174  else
4175  {
4176  NC_TYPE_INFO_T *type;
4177 
4178  /* See if it's a used-defined type */
4179  if ((retval = nc4_find_type(h5, xtype, &type)))
4180  BAIL(retval);
4181  if (!type)
4182  BAIL(NC_EBADTYPE);
4183 
4184  *type_class = type->nc_type_class;
4185  }
4186 
4187 exit:
4188  return retval;
4189 }
4190 
4200 void
4201 reportobject(int log, hid_t id, unsigned int type)
4202 {
4203  char name[MAXNAME];
4204  ssize_t len;
4205  const char* typename = NULL;
4206 
4207  len = H5Iget_name(id, name, MAXNAME);
4208  if(len < 0) return;
4209  name[len] = '\0';
4210 
4211  switch (type) {
4212  case H5F_OBJ_FILE: typename = "File"; break;
4213  case H5F_OBJ_DATASET: typename = "Dataset"; break;
4214  case H5F_OBJ_GROUP: typename = "Group"; break;
4215  case H5F_OBJ_DATATYPE: typename = "Datatype"; break;
4216  case H5F_OBJ_ATTR:
4217  typename = "Attribute";
4218  len = H5Aget_name(id, MAXNAME, name);
4219  if(len < 0) len = 0;
4220  name[len] = '\0';
4221  break;
4222  default: typename = "<unknown>"; break;
4223  }
4224  if(log) {
4225 #ifdef LOGGING
4226  LOG((0,"Type = %s(%8" PRId64 ") name='%s'",typename,id,name));
4227 #endif
4228  } else {
4229  fprintf(stderr,"Type = %s(%8" PRId64 ") name='%s'",typename,id,name);
4230  }
4231 }
4232 
4243 static void
4244 reportopenobjectsT(int log, hid_t fid, int ntypes, unsigned int* otypes)
4245 {
4246  int t,i;
4247  ssize_t ocount;
4248  size_t maxobjs = -1;
4249  hid_t* idlist = NULL;
4250 
4251  if(log) {
4252 #ifdef LOGGING
4253  LOG((0,"\nReport: open objects on %" PRId64 "\n",fid));
4254 #endif
4255  } else {
4256  fprintf(stdout,"\nReport: open objects on %" PRId64 "\n",fid);
4257  }
4258  maxobjs = H5Fget_obj_count(fid,H5F_OBJ_ALL);
4259  if(idlist != NULL) free(idlist);
4260  idlist = (hid_t*)malloc(sizeof(hid_t)*maxobjs);
4261  for(t=0;t<ntypes;t++) {
4262  unsigned int ot = otypes[t];
4263  ocount = H5Fget_obj_ids(fid,ot,maxobjs,idlist);
4264  for(i=0;i<ocount;i++) {
4265  hid_t o = idlist[i];
4266  reportobject(log,o,ot);
4267  }
4268  }
4269  if(idlist != NULL) free(idlist);
4270 }
4271 
4280 void
4281 reportopenobjects(int log, hid_t fid)
4282 {
4283  reportopenobjectsT(log, fid,5,OTYPES);
4284 }
4285 
4297 int
4298 NC4_hdf5get_libversion(unsigned* major,unsigned* minor,unsigned* release)
4299 {
4300  if(H5get_libversion(major,minor,release) < 0)
4301  return NC_EHDFERR;
4302  return NC_NOERR;
4303 }
4304 
4315 int
4316 NC4_hdf5get_superblock(struct NC_HDF5_FILE_INFO* h5, int* idp)
4317 {
4318  int stat = NC_NOERR;
4319  unsigned super;
4320  hid_t plist = -1;
4321  if((plist = H5Fget_create_plist(h5->hdfid)) < 0)
4322  {stat = NC_EHDFERR; goto done;}
4323  if(H5Pget_version(plist, &super, NULL, NULL, NULL) < 0)
4324  {stat = NC_EHDFERR; goto done;}
4325  if(idp) *idp = (int)super;
4326 done:
4327  if(plist >= 0) H5Pclose(plist);
4328  return stat;
4329 }
4330 
4331 static int NC4_get_strict_att(NC_HDF5_FILE_INFO_T*);
4332 static int NC4_walk(hid_t, int*);
4333 
4358 int
4359 NC4_isnetcdf4(struct NC_HDF5_FILE_INFO* h5)
4360 {
4361  int stat;
4362  int isnc4 = 0;
4363  int count;
4364 
4365  /* Look for NC3_STRICT_ATT_NAME */
4366  isnc4 = NC4_get_strict_att(h5);
4367  if(isnc4 > 0)
4368  goto done;
4369  /* attribute did not exist */
4370  /* => last resort: walk the HDF5 file looking for markers */
4371  count = 0;
4372  stat = NC4_walk(h5->root_grp->hdf_grpid, &count);
4373  if(stat != NC_NOERR)
4374  isnc4 = 0;
4375  else /* Threshold is at least two matches */
4376  isnc4 = (count >= 2);
4377 
4378 done:
4379  return isnc4;
4380 }
4381 
4390 static int
4391 NC4_get_strict_att(NC_HDF5_FILE_INFO_T* h5)
4392 {
4393  hid_t grp = -1;
4394  hid_t attid = -1;
4395 
4396  /* Get root group */
4397  grp = h5->root_grp->hdf_grpid; /* get root group */
4398  /* Try to extract the NC3_STRICT_ATT_NAME attribute */
4399  attid = H5Aopen_name(grp, NC3_STRICT_ATT_NAME);
4400  H5Aclose(attid);
4401  return attid;
4402 }
4403 
4413 static int
4414 NC4_walk(hid_t gid, int* countp)
4415 {
4416  int ncstat = NC_NOERR;
4417  int i,j,na;
4418  ssize_t len;
4419  hsize_t nobj;
4420  herr_t err;
4421  int otype;
4422  hid_t grpid, dsid;
4423  char name[NC_HDF5_MAX_NAME];
4424 
4425  /* walk group members of interest */
4426  err = H5Gget_num_objs(gid, &nobj);
4427  if(err < 0) return err;
4428 
4429  for(i = 0; i < nobj; i++) {
4430  /* Get name & kind of object in the group */
4431  len = H5Gget_objname_by_idx(gid,(hsize_t)i,name,(size_t)NC_HDF5_MAX_NAME);
4432  if(len < 0) return len;
4433 
4434  otype = H5Gget_objtype_by_idx(gid,(size_t)i);
4435  switch(otype) {
4436  case H5G_GROUP:
4437  grpid = H5Gopen(gid,name);
4438  NC4_walk(grpid,countp);
4439  H5Gclose(grpid);
4440  break;
4441  case H5G_DATASET: /* variables */
4442  /* Check for phony_dim */
4443  if(strcmp(name,"phony_dim")==0)
4444  *countp = *countp + 1;
4445  dsid = H5Dopen(gid,name);
4446  na = H5Aget_num_attrs(dsid);
4447  for(j = 0; j < na; j++) {
4448  hid_t aid = H5Aopen_idx(dsid,(unsigned int) j);
4449  if(aid >= 0) {
4450  const char** p;
4451  ssize_t len = H5Aget_name(aid, NC_HDF5_MAX_NAME, name);
4452  if(len < 0) return len;
4453  /* Is this a netcdf-4 marker attribute */
4454  for(p=NC_RESERVED_VARATT_LIST;*p;p++) {
4455  if(strcmp(name,*p) == 0) {
4456  *countp = *countp + 1;
4457  }
4458  }
4459  }
4460  H5Aclose(aid);
4461  }
4462  H5Dclose(dsid);
4463  break;
4464  default:/* ignore */
4465  break;
4466  }
4467  }
4468  return ncstat;
4469 }
#define NC_FILL_UBYTE
Default fill value.
Definition: netcdf.h:72
#define NC_ENOMEM
Memory allocation (malloc) failure.
Definition: netcdf.h:395
#define NC_CHAR
ISO/ASCII character.
Definition: netcdf.h:35
#define NC_FILL_CHAR
Default fill value.
Definition: netcdf.h:67
#define NC_UBYTE
unsigned 1 byte int
Definition: netcdf.h:41
#define NC_EDIMSCALE
Problem with HDF5 dimscales.
Definition: netcdf.h:450
#define NC_CLASSIC_MODEL
Enforce classic model on netCDF-4.
Definition: netcdf.h:135
#define NC_ERANGE
Math result not representable.
Definition: netcdf.h:394
#define NC_MAX_VAR_DIMS
max per variable dimensions
Definition: netcdf.h:266
#define NC_FILL_UINT
Default fill value.
Definition: netcdf.h:74
#define NC_UINT
unsigned 4-byte int
Definition: netcdf.h:43
#define NC_EHDFERR
Error at HDF5 layer.
Definition: netcdf.h:427
#define NC_OPAQUE
opaque types
Definition: netcdf.h:53
#define NC_EINVALCOORDS
Index exceeds dimension bound.
Definition: netcdf.h:347
#define NC_INT64
signed 8-byte int
Definition: netcdf.h:44
#define NC_STRING
string
Definition: netcdf.h:46
#define NC_DOUBLE
double precision floating point number
Definition: netcdf.h:40
#define NC_FILL_UINT64
Default fill value.
Definition: netcdf.h:76
#define NC_INDEPENDENT
Use with nc_var_par_access() to set parallel access to independent.
Definition: netcdf_par.h:23
#define NC_ECANTEXTEND
Attempt to extend dataset during ind.
Definition: netcdf.h:456
Main header file for the Parallel C API.
int nc_type
The nc_type type is just an int.
Definition: netcdf.h:24
#define NC_FILL_STRING
Default fill value.
Definition: netcdf.h:77
#define NC_COLLECTIVE
Use with nc_var_par_access() to set parallel access to collective.
Definition: netcdf_par.h:25
#define H5Z_FILTER_SZIP
ID of HDF SZIP filter.
Definition: dvarinq.c:15
#define NC_BYTE
signed 1 byte integer
Definition: netcdf.h:34
#define NC_EINDEFINE
Operation not allowed in define mode.
Definition: netcdf.h:340
#define NC_FILL_INT
Default fill value.
Definition: netcdf.h:69
size_t len
Length of VL data (in base type units)
Definition: netcdf.h:668
#define NC_ENDIAN_LITTLE
In HDF5 files you can set the endianness of variables with nc_def_var_endian().
Definition: netcdf.h:276
#define NC_EATTMETA
Problem with attribute metadata.
Definition: netcdf.h:433
#define NC_VLEN
vlen (variable-length) types
Definition: netcdf.h:52
#define NC_EMPI
MPI operation failed.
Definition: netcdf.h:457
#define NC_EFILEMETA
Problem with file metadata.
Definition: netcdf.h:431
#define NC_EFILTER
Filter operation failed.
Definition: netcdf.h:459
#define NC_EBADTYPE
Not a netcdf data type.
Definition: netcdf.h:357
#define NC_EEDGE
Start+count exceeds dimension bound.
Definition: netcdf.h:385
#define NC_EDIMMETA
Problem with dimension metadata.
Definition: netcdf.h:432
#define NC_EINVAL
Invalid Argument.
Definition: netcdf.h:325
#define NC_INT
signed 4 byte integer
Definition: netcdf.h:37
#define NC_ENDIAN_BIG
In HDF5 files you can set the endianness of variables with nc_def_var_endian().
Definition: netcdf.h:277
#define NC_MAX_NAME
Maximum for classic library.
Definition: netcdf.h:265
void * p
Pointer to VL data.
Definition: netcdf.h:669
#define NC_NAT
Not A Type.
Definition: netcdf.h:33
EXTERNL int nc_free_vlen(nc_vlen_t *vl)
Free memory in a VLEN object.
Definition: dvlen.c:31
#define NC_USHORT
unsigned 2-byte int
Definition: netcdf.h:42
#define NC_EPARINIT
Error initializing for parallel access.
Definition: netcdf.h:441
#define NC_FILL_FLOAT
Default fill value.
Definition: netcdf.h:70
#define NC_EVARMETA
Problem with variable metadata.
Definition: netcdf.h:434
This is the type of arrays of vlens.
Definition: netcdf.h:667
#define NC_SHORT
signed 2 byte integer
Definition: netcdf.h:36
#define NC_ENDIAN_NATIVE
In HDF5 files you can set the endianness of variables with nc_def_var_endian().
Definition: netcdf.h:275
#define NC_ENOTVAR
Variable not found.
Definition: netcdf.h:369
#define MAXNAME
Max HDF5 name.
Definition: nc4hdf.c:37
#define NC_EPERM
Write to read only.
Definition: netcdf.h:326
#define NC_NOERR
No Error.
Definition: netcdf.h:315
#define NC_ENUM
enum types
Definition: netcdf.h:54
#define NC_ECHAR
Attempt to convert between text & numbers.
Definition: netcdf.h:376
#define NC_FILL_SHORT
Default fill value.
Definition: netcdf.h:68
#define NC_FILL_DOUBLE
Default fill value.
Definition: netcdf.h:71
#define NC_COMPOUND
compound types
Definition: netcdf.h:55
#define NC_FILL_USHORT
Default fill value.
Definition: netcdf.h:73
#define NC_FILL_BYTE
Default fill value.
Definition: netcdf.h:66
#define NC_GLOBAL
Attribute id to put/get a global attribute.
Definition: netcdf.h:238
#define NC_FLOAT
single precision floating point number
Definition: netcdf.h:39
#define NC_FILL_INT64
Default fill value.
Definition: netcdf.h:75
#define NC_UINT64
unsigned 8-byte int
Definition: netcdf.h:45

Return to the Main Unidata NetCDF page.
Generated on Fri May 11 2018 21:22:25 for NetCDF. NetCDF is a Unidata library.