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