NetCDF 4.8.0
Loading...
Searching...
No Matches
nc4hdf.c
Go to the documentation of this file.
1/* Copyright 2018, University Corporation for Atmospheric
2 * Research. See the COPYRIGHT file for copying and redistribution
3 * conditions. */
17#include "config.h"
18#include "netcdf.h"
19#include "nc4internal.h"
20#include "ncdispatch.h"
21#include "hdf5internal.h"
22#include "hdf5err.h" /* For BAIL2 */
23#include "hdf5debug.h"
24#include <math.h>
25
26#ifdef HAVE_INTTYPES_H
27#define __STDC_FORMAT_MACROS
28#include <inttypes.h>
29#endif
30
31#define NC_HDF5_MAX_NAME 1024
41static int
42flag_atts_dirty(NCindex *attlist) {
43
44 NC_ATT_INFO_T *att = NULL;
45 int i;
46
47 if(attlist == NULL) {
48 return NC_NOERR;
49 }
50
51 for(i=0;i<ncindexsize(attlist);i++) {
52 att = (NC_ATT_INFO_T*)ncindexith(attlist,i);
53 if(att == NULL) continue;
54 att->dirty = NC_TRUE;
55 }
56
57 return NC_NOERR;
58}
59
76int
77rec_reattach_scales(NC_GRP_INFO_T *grp, int dimid, hid_t dimscaleid)
78{
79 NC_VAR_INFO_T *var;
80 NC_GRP_INFO_T *child_grp;
81 int d, i;
82 int retval;
83
84 assert(grp && grp->hdr.name && dimid >= 0 && dimscaleid >= 0);
85 LOG((3, "%s: grp->hdr.name %s", __func__, grp->hdr.name));
86
87 /* If there are any child groups, attach dimscale there, if needed. */
88 for (i = 0; i < ncindexsize(grp->children); i++)
89 {
90 child_grp = (NC_GRP_INFO_T*)ncindexith(grp->children, i);
91 assert(child_grp);
92 if ((retval = rec_reattach_scales(child_grp, dimid, dimscaleid)))
93 return retval;
94 }
95
96 /* Find any vars that use this dimension id. */
97 for (i = 0; i < ncindexsize(grp->vars); i++)
98 {
99 NC_HDF5_VAR_INFO_T *hdf5_var;
100
101 var = (NC_VAR_INFO_T*)ncindexith(grp->vars,i);
102 assert(var && var->format_var_info);
103 hdf5_var = (NC_HDF5_VAR_INFO_T *)var->format_var_info;
104
105 hdf5_var = (NC_HDF5_VAR_INFO_T*)var->format_var_info;
106 assert(hdf5_var != NULL);
107 for (d = 0; d < var->ndims; d++)
108 {
109 if (var->dimids[d] == dimid && !hdf5_var->dimscale)
110 {
111 LOG((2, "%s: attaching scale for dimid %d to var %s",
112 __func__, var->dimids[d], var->hdr.name));
113 if (var->created)
114 {
115 if (H5DSattach_scale(hdf5_var->hdf_datasetid,
116 dimscaleid, d) < 0)
117 return NC_EHDFERR;
118 hdf5_var->dimscale_attached[d] = NC_TRUE;
119 }
120 }
121 }
122 }
123 return NC_NOERR;
124}
125
142int
143rec_detach_scales(NC_GRP_INFO_T *grp, int dimid, hid_t dimscaleid)
144{
145 NC_VAR_INFO_T *var;
146 NC_GRP_INFO_T *child_grp;
147 int d, i;
148 int retval;
149
150 assert(grp && grp->hdr.name && dimid >= 0 && dimscaleid >= 0);
151 LOG((3, "%s: grp->hdr.name %s", __func__, grp->hdr.name));
152
153 /* If there are any child groups, detach dimscale there, if needed. */
154 for(i=0;i<ncindexsize(grp->children);i++) {
155 child_grp = (NC_GRP_INFO_T*)ncindexith(grp->children,i);
156 if(child_grp == NULL) continue;
157 if ((retval = rec_detach_scales(child_grp, dimid, dimscaleid)))
158 return retval;
159 }
160
161 /* Find any vars that use this dimension id. */
162 for (i = 0; i < ncindexsize(grp->vars); i++)
163 {
164 NC_HDF5_VAR_INFO_T *hdf5_var;
165 var = (NC_VAR_INFO_T*)ncindexith(grp->vars, i);
166 assert(var && var->format_var_info);
167 hdf5_var = (NC_HDF5_VAR_INFO_T *)var->format_var_info;
168
169 for (d = 0; d < var->ndims; d++)
170 {
171 if (var->dimids[d] == dimid && !hdf5_var->dimscale)
172 {
173 LOG((2, "%s: detaching scale for dimid %d to var %s",
174 __func__, var->dimids[d], var->hdr.name));
175 if (var->created)
176 {
177 if (hdf5_var->dimscale_attached && hdf5_var->dimscale_attached[d])
178 {
179 if (H5DSdetach_scale(hdf5_var->hdf_datasetid,
180 dimscaleid, d) < 0)
181 return NC_EHDFERR;
182 hdf5_var->dimscale_attached[d] = NC_FALSE;
183 }
184 }
185 }
186 }
187 }
188 return NC_NOERR;
189}
190
202int
203nc4_open_var_grp2(NC_GRP_INFO_T *grp, int varid, hid_t *dataset)
204{
205 NC_VAR_INFO_T *var;
206 NC_HDF5_VAR_INFO_T *hdf5_var;
207
208 assert(grp && grp->format_grp_info && dataset);
209
210 /* Find the requested varid. */
211 if (!(var = (NC_VAR_INFO_T *)ncindexith(grp->vars, varid)))
212 return NC_ENOTVAR;
213 assert(var && var->hdr.id == varid && var->format_var_info);
214 hdf5_var = (NC_HDF5_VAR_INFO_T *)var->format_var_info;
215
216 /* Open this dataset if necessary. */
217 if (!hdf5_var->hdf_datasetid)
218 {
219 NC_HDF5_GRP_INFO_T *hdf5_grp;
220 hdf5_grp = (NC_HDF5_GRP_INFO_T *)grp->format_grp_info;
221
222 if ((hdf5_var->hdf_datasetid = H5Dopen2(hdf5_grp->hdf_grpid,
223 var->hdr.name, H5P_DEFAULT)) < 0)
224 return NC_ENOTVAR;
225 }
226
227 *dataset = hdf5_var->hdf_datasetid;
228
229 return NC_NOERR;
230}
231
249int
250nc4_get_hdf_typeid(NC_FILE_INFO_T *h5, nc_type xtype,
251 hid_t *hdf_typeid, int endianness)
252{
253 NC_TYPE_INFO_T *type;
254 hid_t typeid = 0;
255 int retval = NC_NOERR;
256
257 assert(hdf_typeid && h5);
258
259 *hdf_typeid = -1;
260
261 /* Determine an appropriate HDF5 datatype */
262 if (xtype == NC_NAT)
263 return NC_EBADTYPE;
264 else if (xtype == NC_CHAR || xtype == NC_STRING)
265 {
266 /* NC_CHAR & NC_STRING types create a new HDF5 datatype */
267 if (xtype == NC_CHAR)
268 {
269 if ((typeid = H5Tcopy(H5T_C_S1)) < 0)
270 return NC_EHDFERR;
271 if (H5Tset_strpad(typeid, H5T_STR_NULLTERM) < 0)
272 BAIL(NC_EVARMETA);
273 if(H5Tset_cset(typeid, H5T_CSET_ASCII) < 0)
274 BAIL(NC_EVARMETA);
275
276 /* Take ownership of the newly created HDF5 datatype */
277 *hdf_typeid = typeid;
278 typeid = 0;
279 }
280 else
281 {
282 if ((typeid = H5Tcopy(H5T_C_S1)) < 0)
283 return NC_EHDFERR;
284 if (H5Tset_size(typeid, H5T_VARIABLE) < 0)
285 BAIL(NC_EVARMETA);
286 if(H5Tset_cset(typeid, H5T_CSET_UTF8) < 0)
287 BAIL(NC_EVARMETA);
288
289 /* Take ownership of the newly created HDF5 datatype */
290 *hdf_typeid = typeid;
291 typeid = 0;
292 }
293 }
294 else
295 {
296 /* All other types use an existing HDF5 datatype */
297 switch (xtype)
298 {
299 case NC_BYTE: /* signed 1 byte integer */
300 if (endianness == NC_ENDIAN_LITTLE)
301 typeid = H5T_STD_I8LE;
302 else if (endianness == NC_ENDIAN_BIG)
303 typeid = H5T_STD_I8BE;
304 else
305 typeid = H5T_NATIVE_SCHAR;
306 break;
307
308 case NC_SHORT: /* signed 2 byte integer */
309 if (endianness == NC_ENDIAN_LITTLE)
310 typeid = H5T_STD_I16LE;
311 else if (endianness == NC_ENDIAN_BIG)
312 typeid = H5T_STD_I16BE;
313 else
314 typeid = H5T_NATIVE_SHORT;
315 break;
316
317 case NC_INT:
318 if (endianness == NC_ENDIAN_LITTLE)
319 typeid = H5T_STD_I32LE;
320 else if (endianness == NC_ENDIAN_BIG)
321 typeid = H5T_STD_I32BE;
322 else
323 typeid = H5T_NATIVE_INT;
324 break;
325
326 case NC_UBYTE:
327 if (endianness == NC_ENDIAN_LITTLE)
328 typeid = H5T_STD_U8LE;
329 else if (endianness == NC_ENDIAN_BIG)
330 typeid = H5T_STD_U8BE;
331 else
332 typeid = H5T_NATIVE_UCHAR;
333 break;
334
335 case NC_USHORT:
336 if (endianness == NC_ENDIAN_LITTLE)
337 typeid = H5T_STD_U16LE;
338 else if (endianness == NC_ENDIAN_BIG)
339 typeid = H5T_STD_U16BE;
340 else
341 typeid = H5T_NATIVE_USHORT;
342 break;
343
344 case NC_UINT:
345 if (endianness == NC_ENDIAN_LITTLE)
346 typeid = H5T_STD_U32LE;
347 else if (endianness == NC_ENDIAN_BIG)
348 typeid = H5T_STD_U32BE;
349 else
350 typeid = H5T_NATIVE_UINT;
351 break;
352
353 case NC_INT64:
354 if (endianness == NC_ENDIAN_LITTLE)
355 typeid = H5T_STD_I64LE;
356 else if (endianness == NC_ENDIAN_BIG)
357 typeid = H5T_STD_I64BE;
358 else
359 typeid = H5T_NATIVE_LLONG;
360 break;
361
362 case NC_UINT64:
363 if (endianness == NC_ENDIAN_LITTLE)
364 typeid = H5T_STD_U64LE;
365 else if (endianness == NC_ENDIAN_BIG)
366 typeid = H5T_STD_U64BE;
367 else
368 typeid = H5T_NATIVE_ULLONG;
369 break;
370
371 case NC_FLOAT:
372 if (endianness == NC_ENDIAN_LITTLE)
373 typeid = H5T_IEEE_F32LE;
374 else if (endianness == NC_ENDIAN_BIG)
375 typeid = H5T_IEEE_F32BE;
376 else
377 typeid = H5T_NATIVE_FLOAT;
378 break;
379
380 case NC_DOUBLE:
381 if (endianness == NC_ENDIAN_LITTLE)
382 typeid = H5T_IEEE_F64LE;
383 else if (endianness == NC_ENDIAN_BIG)
384 typeid = H5T_IEEE_F64BE;
385 else
386 typeid = H5T_NATIVE_DOUBLE;
387 break;
388
389 default:
390 /* Maybe this is a user defined type? */
391 if (nc4_find_type(h5, xtype, &type))
392 return NC_EBADTYPE;
393 if (!type)
394 return NC_EBADTYPE;
395 typeid = ((NC_HDF5_TYPE_INFO_T *)type->format_type_info)->hdf_typeid;
396 break;
397 }
398 assert(typeid);
399
400 /* Copy the HDF5 datatype, so the function operates uniformly */
401 if ((*hdf_typeid = H5Tcopy(typeid)) < 0)
402 return NC_EHDFERR;
403 typeid = 0;
404 }
405 assert(*hdf_typeid != -1);
406
407exit:
408 if (typeid > 0 && H5Tclose(typeid) < 0)
409 BAIL2(NC_EHDFERR);
410 return retval;
411}
412
427static int
428put_att_grpa(NC_GRP_INFO_T *grp, int varid, NC_ATT_INFO_T *att)
429{
430 NC_HDF5_GRP_INFO_T *hdf5_grp;
431 hid_t datasetid = 0, locid;
432 hid_t attid = 0, spaceid = 0, file_typeid = 0;
433 hid_t existing_att_typeid = 0, existing_attid = 0, existing_spaceid = 0;
434 hsize_t dims[1]; /* netcdf attributes always 1-D. */
435 htri_t attr_exists;
436 void *data;
437 int phoney_data = 99;
438 int retval = NC_NOERR;
439
440 assert(att->hdr.name && grp && grp->format_grp_info);
441 LOG((3, "%s: varid %d att->hdr.id %d att->hdr.name %s att->nc_typeid %d "
442 "att->len %d", __func__, varid, att->hdr.id, att->hdr.name,
443 att->nc_typeid, att->len));
444
445 /* Get HDF5-specific group info. */
446 hdf5_grp = (NC_HDF5_GRP_INFO_T *)grp->format_grp_info;
447
448 /* If the file is read-only, return an error. */
449 if (grp->nc4_info->no_write)
450 BAIL(NC_EPERM);
451
452 /* Get the hid to attach the attribute to, or read it from. */
453 if (varid == NC_GLOBAL)
454 locid = hdf5_grp->hdf_grpid;
455 else
456 {
457 if ((retval = nc4_open_var_grp2(grp, varid, &datasetid)))
458 BAIL(retval);
459 locid = datasetid;
460 }
461
462 /* Get the length ready, and find the HDF type we'll be
463 * writing. */
464 dims[0] = att->len;
465 if ((retval = nc4_get_hdf_typeid(grp->nc4_info, att->nc_typeid,
466 &file_typeid, 0)))
467 BAIL(retval);
468
469 /* Even if the length is zero, HDF5 won't let me write with a
470 * NULL pointer. So if the length of the att is zero, point to
471 * some phoney data (which won't be written anyway.)*/
472 if (!dims[0])
473 data = &phoney_data;
474 else if (att->data)
475 data = att->data;
476 else if (att->stdata)
477 data = att->stdata;
478 else
479 data = att->vldata;
480
481 /* NC_CHAR types require some extra work. The space ID is set to
482 * scalar, and the type is told how long the string is. If it's
483 * really zero length, set the size to 1. (The fact that it's
484 * really zero will be marked by the NULL dataspace, but HDF5
485 * doesn't allow me to set the size of the type to zero.)*/
486 if (att->nc_typeid == NC_CHAR)
487 {
488 size_t string_size = dims[0];
489 if (!string_size)
490 {
491 string_size = 1;
492 if ((spaceid = H5Screate(H5S_NULL)) < 0)
493 BAIL(NC_EATTMETA);
494 }
495 else
496 {
497 if ((spaceid = H5Screate(H5S_SCALAR)) < 0)
498 BAIL(NC_EATTMETA);
499 }
500 if (H5Tset_size(file_typeid, string_size) < 0)
501 BAIL(NC_EATTMETA);
502 if (H5Tset_strpad(file_typeid, H5T_STR_NULLTERM) < 0)
503 BAIL(NC_EATTMETA);
504 }
505 else
506 {
507 if (!att->len)
508 {
509 if ((spaceid = H5Screate(H5S_NULL)) < 0)
510 BAIL(NC_EATTMETA);
511 }
512 else
513 {
514 if ((spaceid = H5Screate_simple(1, dims, NULL)) < 0)
515 BAIL(NC_EATTMETA);
516 }
517 }
518
519 /* Does the att exists already? */
520 if ((attr_exists = H5Aexists(locid, att->hdr.name)) < 0)
521 BAIL(NC_EHDFERR);
522 if (attr_exists)
523 {
524 hssize_t npoints;
525
526 /* Open the attribute. */
527 if ((existing_attid = H5Aopen(locid, att->hdr.name, H5P_DEFAULT)) < 0)
528 BAIL(NC_EATTMETA);
529
530 /* Find the type of the existing attribute. */
531 if ((existing_att_typeid = H5Aget_type(existing_attid)) < 0)
532 BAIL(NC_EATTMETA);
533
534 /* How big is the attribute? */
535 if ((existing_spaceid = H5Aget_space(existing_attid)) < 0)
536 BAIL(NC_EATTMETA);
537 if ((npoints = H5Sget_simple_extent_npoints(existing_spaceid)) < 0)
538 BAIL(NC_EATTMETA);
539
540 /* For text attributes the size is specified in the datatype
541 and it is enough to compare types using H5Tequal(). */
542 if (!H5Tequal(file_typeid, existing_att_typeid) ||
543 (att->nc_typeid != NC_CHAR && npoints != att->len))
544 {
545 /* The attribute exists but we cannot re-use it. */
546
547 /* Delete the attribute. */
548 if (H5Adelete(locid, att->hdr.name) < 0)
549 BAIL(NC_EHDFERR);
550
551 /* Re-create the attribute with the type and length
552 reflecting the new value (or values). */
553 if ((attid = H5Acreate(locid, att->hdr.name, file_typeid, spaceid,
554 H5P_DEFAULT)) < 0)
555 BAIL(NC_EATTMETA);
556
557 /* Write the values, (even if length is zero). */
558 if (H5Awrite(attid, file_typeid, data) < 0)
559 BAIL(NC_EATTMETA);
560 }
561 else
562 {
563 /* The attribute exists and we can re-use it. */
564
565 /* Write the values, re-using the existing attribute. */
566 if (H5Awrite(existing_attid, file_typeid, data) < 0)
567 BAIL(NC_EATTMETA);
568 }
569 }
570 else
571 {
572 /* The attribute does not exist yet. */
573
574 /* Create the attribute. */
575 if ((attid = H5Acreate(locid, att->hdr.name, file_typeid, spaceid,
576 H5P_DEFAULT)) < 0)
577 BAIL(NC_EATTMETA);
578
579 /* Write the values, (even if length is zero). */
580 if (H5Awrite(attid, file_typeid, data) < 0)
581 BAIL(NC_EATTMETA);
582 }
583
584exit:
585 if (file_typeid && H5Tclose(file_typeid))
586 BAIL2(NC_EHDFERR);
587 if (attid > 0 && H5Aclose(attid) < 0)
588 BAIL2(NC_EHDFERR);
589 if (existing_att_typeid && H5Tclose(existing_att_typeid))
590 BAIL2(NC_EHDFERR);
591 if (existing_attid > 0 && H5Aclose(existing_attid) < 0)
592 BAIL2(NC_EHDFERR);
593 if (spaceid > 0 && H5Sclose(spaceid) < 0)
594 BAIL2(NC_EHDFERR);
595 if (existing_spaceid > 0 && H5Sclose(existing_spaceid) < 0)
596 BAIL2(NC_EHDFERR);
597 return retval;
598}
599
611static int
612write_attlist(NCindex *attlist, int varid, NC_GRP_INFO_T *grp)
613{
614 NC_ATT_INFO_T *att;
615 int retval;
616 int i;
617
618 for(i = 0; i < ncindexsize(attlist); i++)
619 {
620 att = (NC_ATT_INFO_T *)ncindexith(attlist, i);
621 assert(att);
622 if (att->dirty)
623 {
624 LOG((4, "%s: writing att %s to varid %d", __func__, att->hdr.name, varid));
625 if ((retval = put_att_grpa(grp, varid, att)))
626 return retval;
627 att->dirty = NC_FALSE;
628 att->created = NC_TRUE;
629 }
630 }
631 return NC_NOERR;
632}
633
647static int
648write_coord_dimids(NC_VAR_INFO_T *var)
649{
650 NC_HDF5_VAR_INFO_T *hdf5_var;
651 hsize_t coords_len[1];
652 hid_t c_spaceid = -1, c_attid = -1;
653 int retval = NC_NOERR;
654
655 assert(var && var->format_var_info);
656
657 /* Get HDF5-specific var info. */
658 hdf5_var = (NC_HDF5_VAR_INFO_T *)var->format_var_info;
659
660 /* Set up space for attribute. */
661 coords_len[0] = var->ndims;
662 if ((c_spaceid = H5Screate_simple(1, coords_len, coords_len)) < 0)
663 BAIL(NC_EHDFERR);
664
665 /* Create the attribute. */
666 if ((c_attid = H5Acreate(hdf5_var->hdf_datasetid, COORDINATES,
667 H5T_NATIVE_INT, c_spaceid, H5P_DEFAULT)) < 0)
668 BAIL(NC_EHDFERR);
669
670 /* Write our attribute. */
671 if (H5Awrite(c_attid, H5T_NATIVE_INT, var->dimids) < 0)
672 BAIL(NC_EHDFERR);
673
674exit:
675 if (c_spaceid >= 0 && H5Sclose(c_spaceid) < 0)
676 BAIL2(NC_EHDFERR);
677 if (c_attid >= 0 && H5Aclose(c_attid) < 0)
678 BAIL2(NC_EHDFERR);
679 return retval;
680}
681
692static int
693write_netcdf4_dimid(hid_t datasetid, int dimid)
694{
695 hid_t dimid_spaceid = -1, dimid_attid = -1;
696 htri_t attr_exists;
697 int retval = NC_NOERR;
698
699 /* Create the space. */
700 if ((dimid_spaceid = H5Screate(H5S_SCALAR)) < 0)
701 BAIL(NC_EHDFERR);
702
703 /* Does the attribute already exist? If so, don't try to create it. */
704 if ((attr_exists = H5Aexists(datasetid, NC_DIMID_ATT_NAME)) < 0)
705 BAIL(NC_EHDFERR);
706 if (attr_exists)
707 dimid_attid = H5Aopen_by_name(datasetid, ".", NC_DIMID_ATT_NAME,
708 H5P_DEFAULT, H5P_DEFAULT);
709 else
710 /* Create the attribute if needed. */
711 dimid_attid = H5Acreate(datasetid, NC_DIMID_ATT_NAME,
712 H5T_NATIVE_INT, dimid_spaceid, H5P_DEFAULT);
713 if (dimid_attid < 0)
714 BAIL(NC_EHDFERR);
715
716
717 /* Write it. */
718 LOG((4, "%s: writing secret dimid %d", __func__, dimid));
719 if (H5Awrite(dimid_attid, H5T_NATIVE_INT, &dimid) < 0)
720 BAIL(NC_EHDFERR);
721
722exit:
723 /* Close stuff*/
724 if (dimid_spaceid >= 0 && H5Sclose(dimid_spaceid) < 0)
725 BAIL2(NC_EHDFERR);
726 if (dimid_attid >= 0 && H5Aclose(dimid_attid) < 0)
727 BAIL2(NC_EHDFERR);
728
729 return retval;
730}
731
746static int
747var_create_dataset(NC_GRP_INFO_T *grp, NC_VAR_INFO_T *var, nc_bool_t write_dimid)
748{
749 NC_HDF5_GRP_INFO_T *hdf5_grp;
750 NC_HDF5_VAR_INFO_T *hdf5_var;
751 hid_t plistid = 0, access_plistid = 0, typeid = 0, spaceid = 0;
752 hsize_t chunksize[H5S_MAX_RANK], dimsize[H5S_MAX_RANK], maxdimsize[H5S_MAX_RANK];
753 int d;
754 void *fillp = NULL;
755 NC_DIM_INFO_T *dim = NULL;
756 char *name_to_use;
757 int retval;
758 unsigned int* params = NULL;
759
760 assert(grp && grp->format_grp_info && var && var->format_var_info);
761
762 LOG((3, "%s:: name %s", __func__, var->hdr.name));
763
764 /* Get HDF5-specific group and var info. */
765 hdf5_grp = (NC_HDF5_GRP_INFO_T *)grp->format_grp_info;
766 hdf5_var = (NC_HDF5_VAR_INFO_T *)var->format_var_info;
767
768 /* Scalar or not, we need a creation property list. */
769 if ((plistid = H5Pcreate(H5P_DATASET_CREATE)) < 0)
770 BAIL(NC_EHDFERR);
771 if ((access_plistid = H5Pcreate(H5P_DATASET_ACCESS)) < 0)
772 BAIL(NC_EHDFERR);
773
774 /* Turn off object tracking times in HDF5. */
775 if (H5Pset_obj_track_times(plistid, 0) < 0)
776 BAIL(NC_EHDFERR);
777
778 /* Find the HDF5 type of the dataset. */
779 if ((retval = nc4_get_hdf_typeid(grp->nc4_info, var->type_info->hdr.id, &typeid,
780 var->type_info->endianness)))
781 BAIL(retval);
782
783 /* Figure out what fill value to set, if any. */
784 if (var->no_fill)
785 {
786 /* Required to truly turn HDF5 fill values off */
787 if (H5Pset_fill_time(plistid, H5D_FILL_TIME_NEVER) < 0)
788 BAIL(NC_EHDFERR);
789 }
790 else
791 {
792 if ((retval = nc4_get_fill_value(grp->nc4_info, var, &fillp)))
793 BAIL(retval);
794
795 /* If there is a fill value, set it. */
796 if (fillp)
797 {
798 if (var->type_info->nc_type_class == NC_STRING)
799 {
800 if (H5Pset_fill_value(plistid, typeid, fillp) < 0)
801 BAIL(NC_EHDFERR);
802 }
803 else
804 {
805 /* The fill value set in HDF5 must always be presented as
806 * a native type, even if the endianness for this dataset
807 * is non-native. HDF5 will translate the fill value to
808 * the target endiannesss. */
809 hid_t fill_typeid = 0;
810
811 if ((retval = nc4_get_hdf_typeid(grp->nc4_info, var->type_info->hdr.id, &fill_typeid,
813 BAIL(retval);
814 if (H5Pset_fill_value(plistid, fill_typeid, fillp) < 0)
815 {
816 if (H5Tclose(fill_typeid) < 0)
817 BAIL(NC_EHDFERR);
818 BAIL(NC_EHDFERR);
819 }
820 if (H5Tclose(fill_typeid) < 0)
821 BAIL(NC_EHDFERR);
822 }
823 }
824 }
825
826 /* If the user wants to fletcher error correction, set that up now. */
827 /* Since it is a checksum of sorts, flatcher is always applied first */
828 if (var->fletcher32)
829 if (H5Pset_fletcher32(plistid) < 0)
830 BAIL(NC_EHDFERR);
831
832 /* If the user wants to shuffle the data, set that up now. */
833 if (var->shuffle) {
834 if (H5Pset_shuffle(plistid) < 0)
835 BAIL(NC_EHDFERR);
836 }
837
838 /* If the user wants to compress the data, using either zlib
839 * (a.k.a deflate) or szip, or another filter, set that up now.
840 * Szip and zip can be turned on
841 * either directly with nc_def_var_szip/deflate(), or using
842 * nc_def_var_filter(). If the user
843 * has specified a filter, it will be applied here. */
844 if(var->filters != NULL) {
845 int j;
846 NClist* filters = (NClist*)var->filters;
847 for(j=0;j<nclistlength(filters);j++) {
848 struct NC_HDF5_Filter* fi = (struct NC_HDF5_Filter*)nclistget(filters,j);
849 {
850 if(fi->filterid == H5Z_FILTER_DEFLATE) {/* Handle zip case here */
851 unsigned level;
852 if(fi->nparams != 1)
853 BAIL(NC_EFILTER);
854 level = (int)fi->params[0];
855 if(H5Pset_deflate(plistid, level) < 0)
856 BAIL(NC_EFILTER);
857 } else if(fi->filterid == H5Z_FILTER_SZIP) {/* Handle szip case here */
858 int options_mask;
859 int bits_per_pixel;
860 if(fi->nparams != 2)
861 BAIL(NC_EFILTER);
862 options_mask = (int)fi->params[0];
863 bits_per_pixel = (int)fi->params[1];
864 if(H5Pset_szip(plistid, options_mask, bits_per_pixel) < 0)
865 BAIL(NC_EFILTER);
866 } else {
867 herr_t code = H5Pset_filter(plistid, fi->filterid,
868#if 0
869 H5Z_FLAG_MANDATORY,
870#else
871 H5Z_FLAG_OPTIONAL,
872#endif
873 fi->nparams, fi->params);
874 if(code < 0)
875 BAIL(NC_EFILTER);
876 }
877 }
878 }
879 }
880
881 /* If ndims non-zero, get info for all dimensions. We look up the
882 dimids and get the len of each dimension. We need this to create
883 the space for the dataset. In netCDF a dimension length of zero
884 means an unlimited dimension. */
885 if (var->ndims)
886 {
887 int unlimdim = 0;
888
889 /* Check to see if any unlimited dimensions are used in this var. */
890 for (d = 0; d < var->ndims; d++) {
891 dim = var->dim[d];
892 assert(dim && dim->hdr.id == var->dimids[d]);
893 if (dim->unlimited)
894 unlimdim++;
895 }
896
897 /* If there are no unlimited dims, and no filters, and the user
898 * has not specified chunksizes, use contiguous variable for
899 * better performance. */
900 if (!var->shuffle && !var->fletcher32 && nclistlength((NClist*)var->filters) == 0 &&
901 (var->chunksizes == NULL || !var->chunksizes[0]) && !unlimdim)
902 var->storage = NC_CONTIGUOUS;
903
904 /* Gather current & maximum dimension sizes, along with chunk
905 * sizes. */
906 for (d = 0; d < var->ndims; d++)
907 {
908 dim = var->dim[d];
909 assert(dim && dim->hdr.id == var->dimids[d]);
910 dimsize[d] = dim->unlimited ? NC_HDF5_UNLIMITED_DIMSIZE : dim->len;
911 maxdimsize[d] = dim->unlimited ? H5S_UNLIMITED : (hsize_t)dim->len;
912 if (var->storage == NC_CHUNKED)
913 {
914 if (var->chunksizes[d])
915 chunksize[d] = var->chunksizes[d];
916 else
917 {
918 size_t type_size;
919 if (var->type_info->nc_type_class == NC_STRING)
920 type_size = sizeof(char *);
921 else
922 type_size = var->type_info->size;
923
924 /* Unlimited dim always gets chunksize of 1. */
925 if (dim->unlimited)
926 chunksize[d] = 1;
927 else
928 chunksize[d] = pow((double)DEFAULT_CHUNK_SIZE/type_size,
929 1/(double)(var->ndims - unlimdim));
930
931 /* If the chunksize is greater than the dim
932 * length, make it the dim length. */
933 if (!dim->unlimited && chunksize[d] > dim->len)
934 chunksize[d] = dim->len;
935
936 /* Remember the computed chunksize */
937 var->chunksizes[d] = chunksize[d];
938 }
939 }
940 }
941
942 /* Create the dataspace. */
943 if ((spaceid = H5Screate_simple(var->ndims, dimsize, maxdimsize)) < 0)
944 BAIL(NC_EHDFERR);
945 }
946 else
947 {
948 if ((spaceid = H5Screate(H5S_SCALAR)) < 0)
949 BAIL(NC_EHDFERR);
950 }
951
952 /* Set the var storage to contiguous, compact, or chunked. Don't
953 * try to set chunking for scalar vars, they will default to
954 * contiguous if not set to compact. */
955 if (var->storage == NC_CONTIGUOUS)
956 {
957 if (H5Pset_layout(plistid, H5D_CONTIGUOUS) < 0)
958 BAIL(NC_EHDFERR);
959 }
960 else if (var->storage == NC_COMPACT)
961 {
962 if (H5Pset_layout(plistid, H5D_COMPACT) < 0)
963 BAIL(NC_EHDFERR);
964 }
965 else if (var->ndims)
966 {
967 if (H5Pset_chunk(plistid, var->ndims, chunksize) < 0)
968 BAIL(NC_EHDFERR);
969 }
970
971 /* Turn on creation order tracking. */
972 if (H5Pset_attr_creation_order(plistid, H5P_CRT_ORDER_TRACKED|
973 H5P_CRT_ORDER_INDEXED) < 0)
974 BAIL(NC_EHDFERR);
975
976 /* Set per-var chunk cache, for chunked datasets. */
977 if (var->storage == NC_CHUNKED && var->chunk_cache_size)
978 if (H5Pset_chunk_cache(access_plistid, var->chunk_cache_nelems,
979 var->chunk_cache_size, var->chunk_cache_preemption) < 0)
980 BAIL(NC_EHDFERR);
981
982 /* At long last, create the dataset. */
983 name_to_use = var->alt_name ? var->alt_name : var->hdr.name;
984 LOG((4, "%s: about to H5Dcreate2 dataset %s of type 0x%x", __func__,
985 name_to_use, typeid));
986 if ((hdf5_var->hdf_datasetid = H5Dcreate2(hdf5_grp->hdf_grpid, name_to_use, typeid,
987 spaceid, H5P_DEFAULT, plistid, access_plistid)) < 0)
988 BAIL(NC_EHDFERR);
989 var->created = NC_TRUE;
990 var->is_new_var = NC_FALSE;
991
992 /* Always write the hidden coordinates attribute, which lists the
993 * dimids of this var. When present, this speeds opens. When no
994 * present, dimscale matching is used. */
995 if (var->ndims > 1)
996 if ((retval = write_coord_dimids(var)))
997 BAIL(retval);
998
999 /* If this is a dimscale, mark it as such in the HDF5 file. Also
1000 * find the dimension info and store the dataset id of the dimscale
1001 * dataset. */
1002 if (hdf5_var->dimscale)
1003 {
1004 if (H5DSset_scale(hdf5_var->hdf_datasetid, var->hdr.name) < 0)
1005 BAIL(NC_EHDFERR);
1006
1007 /* If this is a multidimensional coordinate variable, write a
1008 * coordinates attribute. */
1009 /* if (var->ndims > 1) */
1010 /* if ((retval = write_coord_dimids(var))) */
1011 /* BAIL(retval); */
1012
1013 /* If desired, write the netCDF dimid. */
1014 if (write_dimid)
1015 if ((retval = write_netcdf4_dimid(hdf5_var->hdf_datasetid, var->dimids[0])))
1016 BAIL(retval);
1017 }
1018
1019 /* Write attributes for this var. */
1020 if ((retval = write_attlist(var->att, var->hdr.id, grp)))
1021 BAIL(retval);
1022 var->attr_dirty = NC_FALSE;
1023
1024exit:
1025 nullfree(params);
1026 if (typeid > 0 && H5Tclose(typeid) < 0)
1027 BAIL2(NC_EHDFERR);
1028 if (plistid > 0 && H5Pclose(plistid) < 0)
1029 BAIL2(NC_EHDFERR);
1030 if (access_plistid > 0 && H5Pclose(access_plistid) < 0)
1031 BAIL2(NC_EHDFERR);
1032 if (spaceid > 0 && H5Sclose(spaceid) < 0)
1033 BAIL2(NC_EHDFERR);
1034 if (fillp)
1035 {
1036 if (var->type_info->nc_type_class == NC_VLEN)
1037 nc_free_vlen((nc_vlen_t *)fillp);
1038 else if (var->type_info->nc_type_class == NC_STRING && *(char **)fillp)
1039 free(*(char **)fillp);
1040 free(fillp);
1041 }
1042
1043 return retval;
1044}
1045
1059int
1060nc4_adjust_var_cache(NC_GRP_INFO_T *grp, NC_VAR_INFO_T *var)
1061{
1062 size_t chunk_size_bytes = 1;
1063 int d;
1064 int retval;
1065
1066 /* Nothing to be done for contiguous or compact data. */
1067 if (var->storage != NC_CHUNKED)
1068 return NC_NOERR;
1069
1070#ifdef USE_PARALLEL4
1071 /* Don't set cache for files using parallel I/O. */
1072 if (grp->nc4_info->parallel)
1073 return NC_NOERR;
1074#endif
1075
1076 /* How many bytes in the chunk? */
1077 for (d = 0; d < var->ndims; d++)
1078 chunk_size_bytes *= var->chunksizes[d];
1079 if (var->type_info->size)
1080 chunk_size_bytes *= var->type_info->size;
1081 else
1082 chunk_size_bytes *= sizeof(char *);
1083
1084 /* If the chunk cache is too small, and the user has not changed
1085 * the default value of the chunk cache size, then increase the
1086 * size of the cache. */
1087 if (var->chunk_cache_size == CHUNK_CACHE_SIZE)
1088 if (chunk_size_bytes > var->chunk_cache_size)
1089 {
1090 var->chunk_cache_size = chunk_size_bytes * DEFAULT_CHUNKS_IN_CACHE;
1091 if (var->chunk_cache_size > MAX_DEFAULT_CACHE_SIZE)
1092 var->chunk_cache_size = MAX_DEFAULT_CACHE_SIZE;
1093 if ((retval = nc4_reopen_dataset(grp, var)))
1094 return retval;
1095 }
1096
1097 return NC_NOERR;
1098}
1099
1115static int
1116commit_type(NC_GRP_INFO_T *grp, NC_TYPE_INFO_T *type)
1117{
1118 NC_HDF5_GRP_INFO_T *hdf5_grp;
1119 NC_HDF5_TYPE_INFO_T *hdf5_type;
1120 hid_t base_hdf_typeid;
1121 int retval;
1122
1123 assert(grp && grp->format_grp_info && type && type->format_type_info);
1124
1125 /* Get HDF5-specific group and type info. */
1126 hdf5_grp = (NC_HDF5_GRP_INFO_T *)grp->format_grp_info;
1127 hdf5_type = (NC_HDF5_TYPE_INFO_T *)type->format_type_info;
1128
1129 /* Did we already record this type? */
1130 if (type->committed)
1131 return NC_NOERR;
1132
1133 /* Is this a compound type? */
1134 if (type->nc_type_class == NC_COMPOUND)
1135 {
1136 NC_FIELD_INFO_T *field;
1137 hid_t hdf_base_typeid, hdf_typeid;
1138 int i;
1139
1140 if ((hdf5_type->hdf_typeid = H5Tcreate(H5T_COMPOUND, type->size)) < 0)
1141 return NC_EHDFERR;
1142 LOG((4, "creating compound type %s hdf_typeid 0x%x", type->hdr.name,
1143 hdf5_type->hdf_typeid));
1144
1145 for(i=0;i<nclistlength(type->u.c.field);i++)
1146 {
1147 field = (NC_FIELD_INFO_T *)nclistget(type->u.c.field, i);
1148 assert(field);
1149 if ((retval = nc4_get_hdf_typeid(grp->nc4_info, field->nc_typeid,
1150 &hdf_base_typeid, type->endianness)))
1151 return retval;
1152
1153 /* If this is an array, create a special array type. */
1154 if (field->ndims)
1155 {
1156 int d;
1157 hsize_t dims[NC_MAX_VAR_DIMS];
1158
1159 for (d = 0; d < field->ndims; d++)
1160 dims[d] = field->dim_size[d];
1161 if ((hdf_typeid = H5Tarray_create(hdf_base_typeid, field->ndims,
1162 dims, NULL)) < 0)
1163 {
1164 if (H5Tclose(hdf_base_typeid) < 0)
1165 return NC_EHDFERR;
1166 return NC_EHDFERR;
1167 }
1168 if (H5Tclose(hdf_base_typeid) < 0)
1169 return NC_EHDFERR;
1170 }
1171 else
1172 hdf_typeid = hdf_base_typeid;
1173 LOG((4, "inserting field %s offset %d hdf_typeid 0x%x", field->hdr.name,
1174 field->offset, hdf_typeid));
1175 if (H5Tinsert(hdf5_type->hdf_typeid, field->hdr.name, field->offset,
1176 hdf_typeid) < 0)
1177 return NC_EHDFERR;
1178 if (H5Tclose(hdf_typeid) < 0)
1179 return NC_EHDFERR;
1180 }
1181 }
1182 else if (type->nc_type_class == NC_VLEN)
1183 {
1184 /* Find the HDF typeid of the base type of this vlen. */
1185 if ((retval = nc4_get_hdf_typeid(grp->nc4_info, type->u.v.base_nc_typeid,
1186 &base_hdf_typeid, type->endianness)))
1187 return retval;
1188
1189 /* Create a vlen type. */
1190 if ((hdf5_type->hdf_typeid = H5Tvlen_create(base_hdf_typeid)) < 0)
1191 return NC_EHDFERR;
1192 }
1193 else if (type->nc_type_class == NC_OPAQUE)
1194 {
1195 /* Create the opaque type. */
1196 if ((hdf5_type->hdf_typeid = H5Tcreate(H5T_OPAQUE, type->size)) < 0)
1197 return NC_EHDFERR;
1198 }
1199 else if (type->nc_type_class == NC_ENUM)
1200 {
1201 NC_ENUM_MEMBER_INFO_T *enum_m;
1202 int i;
1203
1204 if (nclistlength(type->u.e.enum_member) == 0)
1205 return NC_EINVAL;
1206
1207 /* Find the HDF typeid of the base type of this enum. */
1208 if ((retval = nc4_get_hdf_typeid(grp->nc4_info, type->u.e.base_nc_typeid,
1209 &base_hdf_typeid, type->endianness)))
1210 return retval;
1211
1212 /* Create an enum type. */
1213 if ((hdf5_type->hdf_typeid = H5Tenum_create(base_hdf_typeid)) < 0)
1214 return NC_EHDFERR;
1215
1216 /* Add all the members to the HDF5 type. */
1217 for(i=0;i<nclistlength(type->u.e.enum_member);i++) {
1218 enum_m = (NC_ENUM_MEMBER_INFO_T*)nclistget(type->u.e.enum_member,i);
1219 if (H5Tenum_insert(hdf5_type->hdf_typeid, enum_m->name, enum_m->value) < 0)
1220 return NC_EHDFERR;
1221 }
1222 }
1223 else
1224 {
1225 LOG((0, "Unknown class: %d", type->nc_type_class));
1226 return NC_EBADTYPE;
1227 }
1228
1229 /* Commit the type. */
1230 if (H5Tcommit(hdf5_grp->hdf_grpid, type->hdr.name, hdf5_type->hdf_typeid) < 0)
1231 return NC_EHDFERR;
1232 type->committed = NC_TRUE;
1233 LOG((4, "just committed type %s, HDF typeid: 0x%x", type->hdr.name,
1234 hdf5_type->hdf_typeid));
1235
1236 /* Later we will always use the native typeid. In this case, it is
1237 * a copy of the same type pointed to by hdf_typeid, but it's
1238 * easier to maintain a copy. */
1239 if ((hdf5_type->native_hdf_typeid = H5Tget_native_type(hdf5_type->hdf_typeid,
1240 H5T_DIR_DEFAULT)) < 0)
1241 return NC_EHDFERR;
1242
1243 return NC_NOERR;
1244}
1245
1256static int
1257write_nc3_strict_att(hid_t hdf_grpid)
1258{
1259 hid_t attid = 0, spaceid = 0;
1260 int one = 1;
1261 int retval = NC_NOERR;
1262 htri_t attr_exists;
1263
1264 /* If the attribute already exists, call that a success and return
1265 * NC_NOERR. */
1266 if ((attr_exists = H5Aexists(hdf_grpid, NC3_STRICT_ATT_NAME)) < 0)
1267 return NC_EHDFERR;
1268 if (attr_exists)
1269 return NC_NOERR;
1270
1271 /* Create the attribute to mark this as a file that needs to obey
1272 * strict netcdf-3 rules. */
1273 if ((spaceid = H5Screate(H5S_SCALAR)) < 0)
1274 BAIL(NC_EFILEMETA);
1275 if ((attid = H5Acreate(hdf_grpid, NC3_STRICT_ATT_NAME,
1276 H5T_NATIVE_INT, spaceid, H5P_DEFAULT)) < 0)
1277 BAIL(NC_EFILEMETA);
1278 if (H5Awrite(attid, H5T_NATIVE_INT, &one) < 0)
1279 BAIL(NC_EFILEMETA);
1280
1281exit:
1282 if (spaceid > 0 && (H5Sclose(spaceid) < 0))
1283 BAIL2(NC_EFILEMETA);
1284 if (attid > 0 && (H5Aclose(attid) < 0))
1285 BAIL2(NC_EFILEMETA);
1286 return retval;
1287}
1288
1301static int
1302create_group(NC_GRP_INFO_T *grp)
1303{
1304 NC_HDF5_GRP_INFO_T *hdf5_grp, *parent_hdf5_grp;
1305 hid_t gcpl_id = -1;
1306 int retval = NC_NOERR;;
1307
1308 assert(grp && grp->format_grp_info && grp->parent &&
1309 grp->parent->format_grp_info);
1310
1311 /* Get HDF5 specific group info for group and parent. */
1312 hdf5_grp = (NC_HDF5_GRP_INFO_T *)grp->format_grp_info;
1313 parent_hdf5_grp = (NC_HDF5_GRP_INFO_T *)grp->parent->format_grp_info;
1314 assert(parent_hdf5_grp->hdf_grpid);
1315
1316 /* Create group, with link_creation_order set in the group
1317 * creation property list. */
1318 if ((gcpl_id = H5Pcreate(H5P_GROUP_CREATE)) < 0)
1319 BAIL(NC_EHDFERR);
1320
1321 /* Set track_times to be FALSE. */
1322 if (H5Pset_obj_track_times(gcpl_id, 0) < 0)
1323 BAIL(NC_EHDFERR);
1324
1325 /* Tell HDF5 to keep track of objects in creation order. */
1326 if (H5Pset_link_creation_order(gcpl_id, H5P_CRT_ORDER_TRACKED|H5P_CRT_ORDER_INDEXED) < 0)
1327 BAIL(NC_EHDFERR);
1328
1329 /* Tell HDF5 to keep track of attributes in creation order. */
1330 if (H5Pset_attr_creation_order(gcpl_id, H5P_CRT_ORDER_TRACKED|H5P_CRT_ORDER_INDEXED) < 0)
1331 BAIL(NC_EHDFERR);
1332
1333 /* Create the group. */
1334 if ((hdf5_grp->hdf_grpid = H5Gcreate2(parent_hdf5_grp->hdf_grpid, grp->hdr.name,
1335 H5P_DEFAULT, gcpl_id, H5P_DEFAULT)) < 0)
1336 BAIL(NC_EHDFERR);
1337
1338exit:
1339 if (gcpl_id > -1 && H5Pclose(gcpl_id) < 0)
1340 BAIL2(NC_EHDFERR);
1341 if (retval)
1342 if (hdf5_grp->hdf_grpid > 0 && H5Gclose(hdf5_grp->hdf_grpid) < 0)
1343 BAIL2(NC_EHDFERR);
1344 return retval;
1345}
1346
1358static int
1359attach_dimscales(NC_GRP_INFO_T *grp)
1360{
1361 NC_VAR_INFO_T *var;
1362 NC_HDF5_VAR_INFO_T *hdf5_var;
1363 int d, v;
1364
1365 /* Attach dimension scales. */
1366 for (v = 0; v < ncindexsize(grp->vars); v++)
1367 {
1368 /* Get pointer to var and HDF5-specific var info. */
1369 var = (NC_VAR_INFO_T *)ncindexith(grp->vars, v);
1370 assert(var && var->format_var_info);
1371 hdf5_var = (NC_HDF5_VAR_INFO_T *)var->format_var_info;
1372
1373 /* Scales themselves do not attach. But I really wish they
1374 * would. */
1375 if (hdf5_var->dimscale)
1376 continue;
1377
1378 /* Find the scale for each dimension, if any, and attach it. */
1379 for (d = 0; d < var->ndims; d++)
1380 {
1381 /* Is there a dimscale for this dimension? */
1382 if (hdf5_var->dimscale_attached)
1383 {
1384 if (!hdf5_var->dimscale_attached[d])
1385 {
1386 hid_t dsid; /* Dataset ID for dimension */
1387 assert(var->dim[d] && var->dim[d]->hdr.id == var->dimids[d] &&
1388 var->dim[d]->format_dim_info);
1389
1390 LOG((2, "%s: attaching scale for dimid %d to var %s",
1391 __func__, var->dimids[d], var->hdr.name));
1392
1393 /* Find dataset ID for dimension */
1394 if (var->dim[d]->coord_var)
1395 dsid = ((NC_HDF5_VAR_INFO_T *)(var->dim[d]->coord_var->format_var_info))->hdf_datasetid;
1396 else
1397 dsid = ((NC_HDF5_DIM_INFO_T *)var->dim[d]->format_dim_info)->hdf_dimscaleid;
1398 assert(dsid > 0);
1399
1400 /* Attach the scale. */
1401 if (H5DSattach_scale(hdf5_var->hdf_datasetid, dsid, d) < 0)
1402 return NC_EHDFERR;
1403 hdf5_var->dimscale_attached[d] = NC_TRUE;
1404 }
1405 }
1406 }
1407 }
1408
1409 return NC_NOERR;
1410}
1411
1422static int
1423var_exists(hid_t grpid, char *name, nc_bool_t *exists)
1424{
1425 htri_t link_exists;
1426
1427 /* Reset the boolean */
1428 *exists = NC_FALSE;
1429
1430 /* Check if the object name exists in the group */
1431 if ((link_exists = H5Lexists(grpid, name, H5P_DEFAULT)) < 0)
1432 return NC_EHDFERR;
1433 if (link_exists)
1434 {
1435 H5G_stat_t statbuf;
1436
1437 /* Get info about the object */
1438 if (H5Gget_objinfo(grpid, name, 1, &statbuf) < 0)
1439 return NC_EHDFERR;
1440
1441 if (H5G_DATASET == statbuf.type)
1442 *exists = NC_TRUE;
1443 }
1444
1445 return NC_NOERR;
1446}
1447
1463static int
1464remove_coord_atts(hid_t hdf_datasetid)
1465{
1466 htri_t attr_exists;
1467
1468 /* If the variable dataset has an optional NC_DIMID_ATT_NAME
1469 * attribute, delete it. */
1470 if ((attr_exists = H5Aexists(hdf_datasetid, NC_DIMID_ATT_NAME)) < 0)
1471 return NC_EHDFERR;
1472 if (attr_exists)
1473 {
1474 if (H5Adelete(hdf_datasetid, NC_DIMID_ATT_NAME) < 0)
1475 return NC_EHDFERR;
1476 }
1477
1478 /* Remove the dimension scale 'CLASS' & 'NAME' attributes. */
1479 if ((attr_exists = H5Aexists(hdf_datasetid,
1480 HDF5_DIMSCALE_CLASS_ATT_NAME)) < 0)
1481 return NC_EHDFERR;
1482 if (attr_exists)
1483 {
1484 if (H5Adelete(hdf_datasetid, HDF5_DIMSCALE_CLASS_ATT_NAME) < 0)
1485 return NC_EHDFERR;
1486 }
1487 if ((attr_exists = H5Aexists(hdf_datasetid,
1488 HDF5_DIMSCALE_NAME_ATT_NAME)) < 0)
1489 return NC_EHDFERR;
1490 if (attr_exists)
1491 {
1492 if (H5Adelete(hdf_datasetid, HDF5_DIMSCALE_NAME_ATT_NAME) < 0)
1493 return NC_EHDFERR;
1494 }
1495 return NC_NOERR;
1496}
1497
1512static int
1513write_var(NC_VAR_INFO_T *var, NC_GRP_INFO_T *grp, nc_bool_t write_dimid)
1514{
1515 NC_HDF5_GRP_INFO_T *hdf5_grp;
1516 NC_HDF5_VAR_INFO_T *hdf5_var;
1517 nc_bool_t replace_existing_var = NC_FALSE;
1518 int retval;
1519
1520 assert(var && var->format_var_info && grp && grp->format_grp_info);
1521
1522 LOG((4, "%s: writing var %s", __func__, var->hdr.name));
1523
1524 /* Get HDF5-specific group and var info. */
1525 hdf5_grp = (NC_HDF5_GRP_INFO_T *)grp->format_grp_info;
1526 hdf5_var = (NC_HDF5_VAR_INFO_T *)var->format_var_info;
1527
1528 /* If the variable has already been created & the fill value changed,
1529 * indicate that the existing variable should be replaced. */
1530 if (var->created && var->fill_val_changed)
1531 {
1532 replace_existing_var = NC_TRUE;
1533 var->fill_val_changed = NC_FALSE;
1534 /* If the variable is going to be replaced, we need to flag any
1535 other attributes associated with the variable as 'dirty', or
1536 else *only* the fill value attribute will be copied over and
1537 the rest will be lost. See
1538 https://github.com/Unidata/netcdf-c/issues/239 */
1539 flag_atts_dirty(var->att);
1540 }
1541
1542 /* Is this a coordinate var that has already been created in
1543 * the HDF5 file as a dimscale dataset? Check for dims with the
1544 * same name in this group. If there is one, check to see if
1545 * this object exists in the HDF group. */
1546 if (var->became_coord_var)
1547 {
1548 if ((NC_DIM_INFO_T *)ncindexlookup(grp->dim, var->hdr.name))
1549 {
1550 nc_bool_t exists;
1551
1552 if ((retval = var_exists(hdf5_grp->hdf_grpid, var->hdr.name, &exists)))
1553 return retval;
1554 if (exists)
1555 {
1556 /* Indicate that the variable already exists, and should
1557 * be replaced. */
1558 replace_existing_var = NC_TRUE;
1559 flag_atts_dirty(var->att);
1560 }
1561 }
1562 }
1563
1564 /* Check dims if the variable will be replaced, so that the
1565 * dimensions will be de-attached and re-attached correctly. */
1566 if (replace_existing_var)
1567 {
1568 NC_DIM_INFO_T *d1;
1569
1570 /* Is there a dim with this var's name? */
1571 if ((d1 = (NC_DIM_INFO_T *)ncindexlookup(grp->dim, var->hdr.name)))
1572 {
1573 nc_bool_t exists;
1574 assert(d1->format_dim_info && d1->hdr.name);
1575
1576 if ((retval = var_exists(hdf5_grp->hdf_grpid, var->hdr.name, &exists)))
1577 return retval;
1578 if (exists)
1579 {
1580 hid_t dsid;
1581
1582 /* Find dataset ID for dimension */
1583 if (d1->coord_var)
1584 dsid = ((NC_HDF5_VAR_INFO_T *)d1->coord_var->format_var_info)->hdf_datasetid;
1585 else
1586 dsid = ((NC_HDF5_DIM_INFO_T *)d1->format_dim_info)->hdf_dimscaleid;
1587 assert(dsid > 0);
1588
1589 /* If we're replacing an existing dimscale dataset, go to
1590 * every var in the file and detach this dimension scale,
1591 * because we have to delete it. */
1592 if ((retval = rec_detach_scales(grp->nc4_info->root_grp,
1593 var->dimids[0], dsid)))
1594 return retval;
1595 }
1596 }
1597 }
1598
1599 /* If this is not a dimension scale, remove any attached scales,
1600 * and delete dimscale attributes from the var. */
1601 if (var->was_coord_var && hdf5_var->dimscale_attached)
1602 {
1603 int d;
1604
1605 /* If the variable already exists in the file, Remove any dimension scale
1606 * attributes from it, if they exist. */
1607 if (var->created)
1608 if ((retval = remove_coord_atts(hdf5_var->hdf_datasetid)))
1609 return retval;
1610
1611 /* If this is a regular var, detach all its dim scales. */
1612 for (d = 0; d < var->ndims; d++)
1613 {
1614 if (hdf5_var->dimscale_attached[d])
1615 {
1616 hid_t dsid; /* Dataset ID for dimension */
1617 assert(var->dim[d] && var->dim[d]->hdr.id == var->dimids[d] &&
1618 var->dim[d]->format_dim_info);
1619
1620 /* Find dataset ID for dimension */
1621 if (var->dim[d]->coord_var)
1622 dsid = ((NC_HDF5_VAR_INFO_T *)var->dim[d]->coord_var->format_var_info)->hdf_datasetid;
1623 else
1624 dsid = ((NC_HDF5_DIM_INFO_T *)var->dim[d]->format_dim_info)->hdf_dimscaleid;
1625 assert(dsid > 0);
1626
1627 /* Detach this dim scale. */
1628 if (H5DSdetach_scale(hdf5_var->hdf_datasetid, dsid, d) < 0)
1629 return NC_EHDFERR;
1630 hdf5_var->dimscale_attached[d] = NC_FALSE;
1631 }
1632 }
1633 }
1634
1635 /* Delete the HDF5 dataset that is to be replaced. */
1636 if (replace_existing_var)
1637 {
1638 /* Free the HDF5 dataset id. */
1639 if (hdf5_var->hdf_datasetid && H5Dclose(hdf5_var->hdf_datasetid) < 0)
1640 return NC_EHDFERR;
1641 hdf5_var->hdf_datasetid = 0;
1642
1643 /* Now delete the variable. */
1644 if (H5Gunlink(hdf5_grp->hdf_grpid, var->hdr.name) < 0)
1645 return NC_EDIMMETA;
1646 }
1647
1648 /* Create the dataset. */
1649 if (var->is_new_var || replace_existing_var)
1650 {
1651 if ((retval = var_create_dataset(grp, var, write_dimid)))
1652 return retval;
1653 }
1654 else
1655 {
1656 if (write_dimid && var->ndims)
1657 if ((retval = write_netcdf4_dimid(hdf5_var->hdf_datasetid,
1658 var->dimids[0])))
1659 return retval;
1660 }
1661
1662 if (replace_existing_var)
1663 {
1664 /* If this is a dimension scale, reattach the scale everywhere it
1665 * is used. (Recall that netCDF dimscales are always 1-D). */
1666 if(hdf5_var->dimscale)
1667 {
1668 if ((retval = rec_reattach_scales(grp->nc4_info->root_grp,
1669 var->dimids[0], hdf5_var->hdf_datasetid)))
1670 return retval;
1671 }
1672 /* If it's not a dimension scale, clear the dimscale attached flags,
1673 * so the dimensions are re-attached. */
1674 else
1675 {
1676 if (hdf5_var->dimscale_attached)
1677 memset(hdf5_var->dimscale_attached, 0, sizeof(nc_bool_t) * var->ndims);
1678 }
1679 }
1680
1681 /* Clear coord. var state transition flags */
1682 var->was_coord_var = NC_FALSE;
1683 var->became_coord_var = NC_FALSE;
1684
1685 /* Now check the attributes for this var. */
1686 if (var->attr_dirty)
1687 {
1688 /* Write attributes for this var. */
1689 if ((retval = write_attlist(var->att, var->hdr.id, grp)))
1690 return retval;
1691 var->attr_dirty = NC_FALSE;
1692 }
1693
1694 return NC_NOERR;
1695}
1696
1710int
1711nc4_create_dim_wo_var(NC_DIM_INFO_T *dim)
1712{
1713 NC_HDF5_DIM_INFO_T *hdf5_dim;
1714 NC_HDF5_GRP_INFO_T *hdf5_grp;
1715 hid_t spaceid = -1, create_propid = -1;
1716 hsize_t dims[1], max_dims[1], chunk_dims[1] = {1};
1717 char dimscale_wo_var[NC_MAX_NAME];
1718 int retval = NC_NOERR;
1719
1720 LOG((4, "%s: creating dim %s", __func__, dim->hdr.name));
1721
1722 /* Sanity check */
1723 assert(!dim->coord_var);
1724
1725 /* Get HDF5-specific dim and group info. */
1726 hdf5_grp = (NC_HDF5_GRP_INFO_T *)dim->container->format_grp_info;
1727 hdf5_dim = (NC_HDF5_DIM_INFO_T *)dim->format_dim_info;
1728
1729 /* Create a property list. */
1730 if ((create_propid = H5Pcreate(H5P_DATASET_CREATE)) < 0)
1731 BAIL(NC_EHDFERR);
1732
1733 /* Turn off recording of times associated with this object. */
1734 if (H5Pset_obj_track_times(create_propid, 0) < 0)
1735 BAIL(NC_EHDFERR);
1736
1737 /* Set size of dataset to size of dimension. */
1738 dims[0] = dim->len;
1739 max_dims[0] = dim->len;
1740
1741 /* If this dimension scale is unlimited (i.e. it's an unlimited
1742 * dimension), then set up chunking, with a chunksize of 1. */
1743 if (dim->unlimited)
1744 {
1745 max_dims[0] = H5S_UNLIMITED;
1746 if (H5Pset_chunk(create_propid, 1, chunk_dims) < 0)
1747 BAIL(NC_EHDFERR);
1748 }
1749
1750 /* Set up space. */
1751 if ((spaceid = H5Screate_simple(1, dims, max_dims)) < 0)
1752 BAIL(NC_EHDFERR);
1753
1754 /* Turn on creation-order tracking. */
1755 if (H5Pset_attr_creation_order(create_propid, H5P_CRT_ORDER_TRACKED|
1756 H5P_CRT_ORDER_INDEXED) < 0)
1757 BAIL(NC_EHDFERR);
1758
1759 /* Create the dataset that will be the dimension scale. */
1760 LOG((4, "%s: about to H5Dcreate1 a dimscale dataset %s", __func__,
1761 dim->hdr.name));
1762 if ((hdf5_dim->hdf_dimscaleid = H5Dcreate2(hdf5_grp->hdf_grpid, dim->hdr.name,
1763 H5T_IEEE_F32BE, spaceid,
1764 H5P_DEFAULT, create_propid,
1765 H5P_DEFAULT)) < 0)
1766 BAIL(NC_EHDFERR);
1767
1768 /* Indicate that this is a scale. Also indicate that not
1769 * be shown to the user as a variable. It is hidden. It is
1770 * a DIM WITHOUT A VARIABLE! */
1771 sprintf(dimscale_wo_var, "%s%10d", DIM_WITHOUT_VARIABLE, (int)dim->len);
1772 if (H5DSset_scale(hdf5_dim->hdf_dimscaleid, dimscale_wo_var) < 0)
1773 BAIL(NC_EHDFERR);
1774
1775 /* Since this dimension was created out of order, we cannot rely on
1776 * it getting the correct dimid on file open. We must assign it
1777 * explicitly. */
1778 if ((retval = write_netcdf4_dimid(hdf5_dim->hdf_dimscaleid, dim->hdr.id)))
1779 BAIL(retval);
1780
1781exit:
1782 if (spaceid > 0 && H5Sclose(spaceid) < 0)
1783 BAIL2(NC_EHDFERR);
1784 if (create_propid > 0 && H5Pclose(create_propid) < 0)
1785 BAIL2(NC_EHDFERR);
1786 return retval;
1787}
1788
1801static int
1802write_dim(NC_DIM_INFO_T *dim, NC_GRP_INFO_T *grp, nc_bool_t write_dimid)
1803{
1804 NC_HDF5_DIM_INFO_T *hdf5_dim;
1805 int retval = NC_NOERR;
1806
1807 assert(dim && dim->format_dim_info && grp && grp->format_grp_info);
1808
1809 /* Get HDF5-specific dim and group info. */
1810 hdf5_dim = (NC_HDF5_DIM_INFO_T *)dim->format_dim_info;
1811
1812 /* If there's no dimscale dataset for this dim, create one,
1813 * and mark that it should be hidden from netCDF as a
1814 * variable. (That is, it should appear as a dimension
1815 * without an associated variable.) */
1816 if (!hdf5_dim->hdf_dimscaleid)
1817 if ((retval = nc4_create_dim_wo_var(dim)))
1818 BAIL(retval);
1819
1820 /* Did we extend an unlimited dimension? */
1821 if (dim->extended)
1822 {
1823 NC_VAR_INFO_T *v1 = NULL;
1824
1825 assert(dim->unlimited);
1826
1827 /* If this is a dimension with an associated coordinate var,
1828 * then update the length of that coord var. */
1829 v1 = dim->coord_var;
1830 if (v1)
1831 {
1832 NC_HDF5_VAR_INFO_T *hdf5_v1;
1833 hsize_t *new_size;
1834 int d1;
1835
1836 hdf5_v1 = (NC_HDF5_VAR_INFO_T *)v1->format_var_info;
1837
1838 /* Extend the dimension scale dataset to reflect the new
1839 * length of the dimension. */
1840 if (!(new_size = malloc(v1->ndims * sizeof(hsize_t))))
1841 BAIL(NC_ENOMEM);
1842 for (d1 = 0; d1 < v1->ndims; d1++)
1843 {
1844 assert(v1->dim[d1] && v1->dim[d1]->hdr.id == v1->dimids[d1]);
1845 new_size[d1] = v1->dim[d1]->len;
1846 }
1847 if (H5Dset_extent(hdf5_v1->hdf_datasetid, new_size) < 0)
1848 BAIL(NC_EHDFERR);
1849 free(new_size);
1850 }
1851 }
1852
1853 /* If desired, write the secret dimid. This will be used instead of
1854 * the dimid that the dimension would otherwise receive based on
1855 * creation order. This can be necessary when dims and their
1856 * coordinate variables were created in different order. */
1857 if (write_dimid && hdf5_dim->hdf_dimscaleid)
1858 if ((retval = write_netcdf4_dimid(hdf5_dim->hdf_dimscaleid, dim->hdr.id)))
1859 BAIL(retval);
1860
1861exit:
1862
1863 return retval;
1864}
1865
1878int
1879nc4_rec_write_metadata(NC_GRP_INFO_T *grp, nc_bool_t bad_coord_order)
1880{
1881 NC_DIM_INFO_T *dim = NULL;
1882 NC_VAR_INFO_T *var = NULL;
1883 NC_GRP_INFO_T *child_grp = NULL;
1884 int coord_varid = -1;
1885 int var_index = 0;
1886 int dim_index = 0;
1887 int retval;
1888 int i;
1889
1890 assert(grp && grp->hdr.name &&
1891 ((NC_HDF5_GRP_INFO_T *)(grp->format_grp_info))->hdf_grpid);
1892 LOG((3, "%s: grp->hdr.name %s, bad_coord_order %d", __func__, grp->hdr.name,
1893 bad_coord_order));
1894
1895 /* Write global attributes for this group. */
1896 if ((retval = write_attlist(grp->att, NC_GLOBAL, grp)))
1897 return retval;
1898
1899 /* Set the pointers to the beginning of the list of dims & vars in this
1900 * group. */
1901 dim = (NC_DIM_INFO_T *)ncindexith(grp->dim, dim_index);
1902 var = (NC_VAR_INFO_T *)ncindexith(grp->vars, var_index);
1903
1904 /* Because of HDF5 ordering the dims and vars have to be stored in
1905 * this way to ensure that the dims and coordinate vars come out in
1906 * the correct order. */
1907 while (dim || var)
1908 {
1909 nc_bool_t found_coord, wrote_coord;
1910
1911 /* Write non-coord dims in order, stopping at the first one that
1912 * has an associated coord var. */
1913 for (found_coord = NC_FALSE; dim && !found_coord; )
1914 {
1915 if (!dim->coord_var)
1916 {
1917 if ((retval = write_dim(dim, grp, bad_coord_order)))
1918 return retval;
1919 }
1920 else
1921 {
1922 coord_varid = dim->coord_var->hdr.id;
1923 found_coord = NC_TRUE;
1924 }
1925 dim = (NC_DIM_INFO_T *)ncindexith(grp->dim, ++dim_index);
1926 }
1927
1928 /* Write each var. When we get to the coord var we are waiting
1929 * for (if any), then we break after writing it. */
1930 for (wrote_coord = NC_FALSE; var && !wrote_coord; )
1931 {
1932 if ((retval = write_var(var, grp, bad_coord_order)))
1933 return retval;
1934 if (found_coord && var->hdr.id == coord_varid)
1935 wrote_coord = NC_TRUE;
1936 var = (NC_VAR_INFO_T *)ncindexith(grp->vars, ++var_index);
1937 }
1938 } /* end while */
1939
1940 /* Attach dimscales to vars in this group. */
1941 if ((retval = attach_dimscales(grp)))
1942 return retval;
1943
1944 /* If there are any child groups, write their metadata. */
1945 for (i = 0; i < ncindexsize(grp->children); i++)
1946 {
1947 child_grp = (NC_GRP_INFO_T *)ncindexith(grp->children, i);
1948 assert(child_grp);
1949 if ((retval = nc4_rec_write_metadata(child_grp, bad_coord_order)))
1950 return retval;
1951 }
1952 return NC_NOERR;
1953}
1954
1964int
1966{
1967 NC_GRP_INFO_T *child_grp;
1968 NC_HDF5_GRP_INFO_T *hdf5_grp;
1969 NC_TYPE_INFO_T *type;
1970 int retval;
1971 int i;
1972
1973 assert(grp && grp->hdr.name && grp->format_grp_info);
1974 LOG((3, "%s: grp->hdr.name %s", __func__, grp->hdr.name));
1975
1976 /* Get HDF5-specific group info. */
1977 hdf5_grp = (NC_HDF5_GRP_INFO_T *)grp->format_grp_info;
1978
1979 /* Create the group in the HDF5 file if it doesn't exist. */
1980 if (!hdf5_grp->hdf_grpid)
1981 if ((retval = create_group(grp)))
1982 return retval;
1983
1984 /* If this is the root group of a file with strict NC3 rules, write
1985 * an attribute. But don't leave the attribute open. */
1986 if (!grp->parent && (grp->nc4_info->cmode & NC_CLASSIC_MODEL))
1987 if ((retval = write_nc3_strict_att(hdf5_grp->hdf_grpid)))
1988 return retval;
1989
1990 /* If there are any user-defined types, write them now. */
1991 for(i=0;i<ncindexsize(grp->type);i++) {
1992 type = (NC_TYPE_INFO_T *)ncindexith(grp->type, i);
1993 assert(type);
1994 if ((retval = commit_type(grp, type)))
1995 return retval;
1996 }
1997
1998 /* If there are any child groups, write their groups and types. */
1999 for(i=0;i<ncindexsize(grp->children);i++) {
2000 if((child_grp = (NC_GRP_INFO_T*)ncindexith(grp->children,i)) == NULL) continue;
2001 if ((retval = nc4_rec_write_groups_types(child_grp)))
2002 return retval;
2003 }
2004 return NC_NOERR;
2005}
2006
2020int
2021nc4_rec_match_dimscales(NC_GRP_INFO_T *grp)
2022{
2023 NC_GRP_INFO_T *g;
2024 NC_VAR_INFO_T *var;
2025 NC_DIM_INFO_T *dim;
2026 int retval = NC_NOERR;
2027 int i;
2028
2029 assert(grp && grp->hdr.name);
2030 LOG((4, "%s: grp->hdr.name %s", __func__, grp->hdr.name));
2031
2032 /* Perform var dimscale match for child groups. */
2033 for (i = 0; i < ncindexsize(grp->children); i++)
2034 {
2035 g = (NC_GRP_INFO_T*)ncindexith(grp->children, i);
2036 assert(g);
2037 if ((retval = nc4_rec_match_dimscales(g)))
2038 return retval;
2039 }
2040
2041 /* Check all the vars in this group. If they have dimscale info,
2042 * try and find a dimension for them. */
2043 for (i = 0; i < ncindexsize(grp->vars); i++)
2044 {
2045 NC_HDF5_VAR_INFO_T *hdf5_var;
2046 int ndims;
2047 int d;
2048
2049 /* Get pointer to var and to the HDF5-specific var info. */
2050 var = (NC_VAR_INFO_T*)ncindexith(grp->vars, i);
2051 assert(var && var->format_var_info);
2052 hdf5_var = (NC_HDF5_VAR_INFO_T *)var->format_var_info;
2053
2054 /* Check all vars and see if dim[i] != NULL if dimids[i] valid. */
2055 /* This loop is very odd. Under normal circumstances, var->dimid[d] is zero
2056 (from the initial calloc) which is a legitimate dimid. The code does not
2057 distinquish this case from the dimscale case where the id might actually
2058 be defined.
2059 The original nc4_find_dim searched up the group tree looking for the given
2060 dimid in one of the dim lists associated with each ancestor group.
2061 I changed nc4_fnd_dim to use the dimid directly using h5->alldims.
2062 However, here that is incorrect because it will find the dimid 0 always
2063 (if any dimensions were defined). Except that when dimscale dimids have
2064 been defined, one or more of the values in var->dimids will have a
2065 legitimate value.
2066 The solution I choose is to modify nc4_var_list_add to initialize dimids to
2067 illegal values (-1). This is another example of the problems with dimscales.
2068 */
2069 ndims = var->ndims;
2070 for (d = 0; d < ndims; d++)
2071 {
2072 if (var->dim[d] == NULL) {
2073 nc4_find_dim(grp, var->dimids[d], &var->dim[d], NULL);
2074 }
2075 /* assert(var->dim[d] && var->dim[d]->hdr.id == var->dimids[d]); */
2076 }
2077
2078 /* Skip dimension scale variables */
2079 if (!hdf5_var->dimscale)
2080 {
2081 int d;
2082 int j;
2083
2084 /* Are there dimscales for this variable? */
2085 if (hdf5_var->dimscale_hdf5_objids)
2086 {
2087 for (d = 0; d < var->ndims; d++)
2088 {
2089 nc_bool_t finished = NC_FALSE;
2090 LOG((5, "%s: var %s has dimscale info...", __func__, var->hdr.name));
2091
2092 /* Check this and parent groups. */
2093 for (g = grp; g && !finished; g = g->parent)
2094 {
2095 /* Check all dims in this group. */
2096 for (j = 0; j < ncindexsize(g->dim); j++)
2097 {
2098 /* Get the HDF5 specific dim info. */
2099 NC_HDF5_DIM_INFO_T *hdf5_dim;
2100 dim = (NC_DIM_INFO_T *)ncindexith(g->dim, j);
2101 assert(dim && dim->format_dim_info);
2102 hdf5_dim = (NC_HDF5_DIM_INFO_T *)dim->format_dim_info;
2103
2104 /* Check for exact match of fileno/objid arrays
2105 * to find identical objects in HDF5 file. */
2106#if H5_VERSION_GE(1,12,0)
2107 int token_cmp;
2108 if (H5Otoken_cmp(hdf5_var->hdf_datasetid, &hdf5_var->dimscale_hdf5_objids[d].token, &hdf5_dim->hdf5_objid.token, &token_cmp) < 0)
2109 return NC_EHDFERR;
2110
2111 if (hdf5_var->dimscale_hdf5_objids[d].fileno == hdf5_dim->hdf5_objid.fileno &&
2112 token_cmp == 0)
2113#else
2114 if (hdf5_var->dimscale_hdf5_objids[d].fileno[0] == hdf5_dim->hdf5_objid.fileno[0] &&
2115 hdf5_var->dimscale_hdf5_objids[d].objno[0] == hdf5_dim->hdf5_objid.objno[0] &&
2116 hdf5_var->dimscale_hdf5_objids[d].fileno[1] == hdf5_dim->hdf5_objid.fileno[1] &&
2117 hdf5_var->dimscale_hdf5_objids[d].objno[1] == hdf5_dim->hdf5_objid.objno[1])
2118#endif
2119 {
2120 LOG((4, "%s: for dimension %d, found dim %s", __func__,
2121 d, dim->hdr.name));
2122 var->dimids[d] = dim->hdr.id;
2123 var->dim[d] = dim;
2124 finished = NC_TRUE;
2125 break;
2126 }
2127 } /* next dim */
2128 } /* next grp */
2129 LOG((5, "%s: dimid for this dimscale is %d", __func__,
2130 var->type_info->hdr.id));
2131 } /* next var->dim */
2132 }
2133 /* No dimscales for this var! Invent phony dimensions. */
2134 else
2135 {
2136 hid_t spaceid = 0;
2137 hsize_t *h5dimlen = NULL, *h5dimlenmax = NULL;
2138 int dataset_ndims;
2139
2140 /* Find the space information for this dimension. */
2141 if ((spaceid = H5Dget_space(hdf5_var->hdf_datasetid)) < 0)
2142 return NC_EHDFERR;
2143
2144 /* Get the len of each dim in the space. */
2145 if (var->ndims)
2146 {
2147 if (!(h5dimlen = malloc(var->ndims * sizeof(hsize_t))))
2148 return NC_ENOMEM;
2149 if (!(h5dimlenmax = malloc(var->ndims * sizeof(hsize_t))))
2150 {
2151 free(h5dimlen);
2152 return NC_ENOMEM;
2153 }
2154 if ((dataset_ndims = H5Sget_simple_extent_dims(spaceid, h5dimlen,
2155 h5dimlenmax)) < 0) {
2156 free(h5dimlenmax);
2157 free(h5dimlen);
2158 return NC_EHDFERR;
2159 }
2160 if (dataset_ndims != var->ndims) {
2161 free(h5dimlenmax);
2162 free(h5dimlen);
2163 return NC_EHDFERR;
2164 }
2165 }
2166 else
2167 {
2168 /* Make sure it's scalar. */
2169 if (H5Sget_simple_extent_type(spaceid) != H5S_SCALAR)
2170 return NC_EHDFERR;
2171 }
2172
2173 /* Release the space object. */
2174 if (H5Sclose(spaceid) < 0) {
2175 free(h5dimlen);
2176 free(h5dimlenmax);
2177 return NC_EHDFERR;
2178 }
2179
2180 /* Create a phony dimension for each dimension in the
2181 * dataset, unless there already is one the correct
2182 * size. */
2183 for (d = 0; d < var->ndims; d++)
2184 {
2185 int k;
2186 int match;
2187 /* Is there already a phony dimension of the correct size? */
2188 for(match=-1,k=0;k<ncindexsize(grp->dim);k++) {
2189 if((dim = (NC_DIM_INFO_T*)ncindexith(grp->dim,k)) == NULL) continue;
2190 if ((dim->len == h5dimlen[d]) &&
2191 ((h5dimlenmax[d] == H5S_UNLIMITED && dim->unlimited) ||
2192 (h5dimlenmax[d] != H5S_UNLIMITED && !dim->unlimited)))
2193 {match = k; break;}
2194 }
2195
2196 /* Didn't find a phony dim? Then create one. */
2197 if (match < 0)
2198 {
2199 char phony_dim_name[NC_MAX_NAME + 1];
2200 sprintf(phony_dim_name, "phony_dim_%d", grp->nc4_info->next_dimid);
2201 LOG((3, "%s: creating phony dim for var %s", __func__, var->hdr.name));
2202 if ((retval = nc4_dim_list_add(grp, phony_dim_name, h5dimlen[d], -1, &dim)))
2203 {
2204 free(h5dimlenmax);
2205 free(h5dimlen);
2206 return retval;
2207 }
2208 /* Create struct for HDF5-specific dim info. */
2209 if (!(dim->format_dim_info = calloc(1, sizeof(NC_HDF5_DIM_INFO_T))))
2210 return NC_ENOMEM;
2211 if (h5dimlenmax[d] == H5S_UNLIMITED)
2212 dim->unlimited = NC_TRUE;
2213 }
2214
2215 /* The variable must remember the dimid. */
2216 var->dimids[d] = dim->hdr.id;
2217 var->dim[d] = dim;
2218 } /* next dim */
2219
2220 /* Free the memory we malloced. */
2221 free(h5dimlen);
2222 free(h5dimlenmax);
2223 }
2224 }
2225 }
2226
2227 return retval;
2228}
2229
2241void
2242reportobject(int uselog, hid_t id, unsigned int type)
2243{
2244 char name[NC_HDF5_MAX_NAME];
2245 ssize_t len;
2246 const char* typename = NULL;
2247 long long printid = (long long)id;
2248
2249 len = H5Iget_name(id, name, NC_HDF5_MAX_NAME);
2250 if(len < 0) return;
2251 name[len] = '\0';
2252
2253 switch (type) {
2254 case H5F_OBJ_FILE: typename = "File"; break;
2255 case H5F_OBJ_DATASET: typename = "Dataset"; break;
2256 case H5F_OBJ_GROUP: typename = "Group"; break;
2257 case H5F_OBJ_DATATYPE: typename = "Datatype"; break;
2258 case H5F_OBJ_ATTR:
2259 typename = "Attribute";
2260 len = H5Aget_name(id, NC_HDF5_MAX_NAME, name);
2261 if(len < 0) len = 0;
2262 name[len] = '\0';
2263 break;
2264 default: typename = "<unknown>"; break;
2265 }
2266#ifdef LOGGING
2267 if(uselog) {
2268 LOG((0,"Type = %s(%lld) name='%s'",typename,printid,name));
2269 } else
2270#endif
2271 {
2272 fprintf(stderr,"Type = %s(%lld) name='%s'",typename,printid,name);
2273 }
2274
2275}
2276
2287static void
2288reportopenobjectsT(int uselog, hid_t fid, int ntypes, unsigned int* otypes)
2289{
2290 int t,i;
2291 ssize_t ocount;
2292 size_t maxobjs = -1;
2293 hid_t* idlist = NULL;
2294
2295 /* Always report somehow */
2296#ifdef LOGGING
2297 if(uselog)
2298 LOG((0,"\nReport: open objects on %lld",(long long)fid));
2299 else
2300#endif
2301 fprintf(stdout,"\nReport: open objects on %lld\n",(long long)fid);
2302 maxobjs = H5Fget_obj_count(fid,H5F_OBJ_ALL);
2303 if(idlist != NULL) free(idlist);
2304 idlist = (hid_t*)malloc(sizeof(hid_t)*maxobjs);
2305 for(t=0;t<ntypes;t++) {
2306 unsigned int ot = otypes[t];
2307 ocount = H5Fget_obj_ids(fid,ot,maxobjs,idlist);
2308 for(i=0;i<ocount;i++) {
2309 hid_t o = idlist[i];
2310 reportobject(uselog,o,ot);
2311 }
2312 }
2313 if(idlist != NULL) free(idlist);
2314}
2315
2324void
2325reportopenobjects(int uselog, hid_t fid)
2326{
2327 unsigned int OTYPES[5] = {H5F_OBJ_FILE, H5F_OBJ_DATASET, H5F_OBJ_GROUP,
2328 H5F_OBJ_DATATYPE, H5F_OBJ_ATTR};
2329
2330 reportopenobjectsT(uselog, fid ,5, OTYPES);
2331}
2332
2340void
2341showopenobjects5(NC_FILE_INFO_T* h5)
2342{
2343 NC_HDF5_FILE_INFO_T *hdf5_info;
2344
2345 assert(h5 && h5->format_file_info);
2346 hdf5_info = (NC_HDF5_FILE_INFO_T *)h5->format_file_info;
2347
2348 fprintf(stderr,"===== begin showopenobjects =====\n");
2349 reportopenobjects(0,hdf5_info->hdfid);
2350 fprintf(stderr,"===== end showopenobjects =====\n");
2351 fflush(stderr);
2352}
2353
2362void
2364{
2365 NC_FILE_INFO_T* h5 = NULL;
2366
2367 /* Find our metadata for this file. */
2368 if (nc4_find_nc_grp_h5(ncid, NULL, NULL, &h5) != NC_NOERR)
2369 fprintf(stderr,"failed\n");
2370 else
2371 showopenobjects5(h5);
2372 fflush(stderr);
2373}
2374
2386int
2387NC4_hdf5get_libversion(unsigned* major,unsigned* minor,unsigned* release)
2388{
2389 if(H5get_libversion(major,minor,release) < 0)
2390 return NC_EHDFERR;
2391 return NC_NOERR;
2392}
2393
2404int
2405NC4_hdf5get_superblock(struct NC_FILE_INFO* h5, int* idp)
2406{
2407 NC_HDF5_FILE_INFO_T *hdf5_info;
2408 int stat = NC_NOERR;
2409 unsigned super;
2410 hid_t plist = -1;
2411
2412 assert(h5 && h5->format_file_info);
2413 hdf5_info = (NC_HDF5_FILE_INFO_T *)h5->format_file_info;
2414
2415 if((plist = H5Fget_create_plist(hdf5_info->hdfid)) < 0)
2416 {stat = NC_EHDFERR; goto done;}
2417 if(H5Pget_version(plist, &super, NULL, NULL, NULL) < 0)
2418 {stat = NC_EHDFERR; goto done;}
2419 if(idp) *idp = (int)super;
2420done:
2421 if(plist >= 0) H5Pclose(plist);
2422 return stat;
2423}
2424
2425static int NC4_strict_att_exists(NC_FILE_INFO_T*);
2426static int NC4_walk(hid_t, int*);
2427
2453int
2454NC4_isnetcdf4(struct NC_FILE_INFO* h5)
2455{
2456 int stat;
2457 int isnc4 = 0;
2458 int exists;
2459 int count;
2460
2461 /* Look for NC3_STRICT_ATT_NAME */
2462 exists = NC4_strict_att_exists(h5);
2463 if(exists)
2464 goto done;
2465 /* attribute did not exist */
2466 /* => last resort: walk the HDF5 file looking for markers */
2467 count = 0;
2468 stat = NC4_walk(((NC_HDF5_GRP_INFO_T *)(h5->root_grp->format_grp_info))->hdf_grpid,
2469 &count);
2470 if(stat != NC_NOERR)
2471 isnc4 = 0;
2472 else /* Threshold is at least two matches */
2473 isnc4 = (count >= 2);
2474
2475done:
2476 return isnc4;
2477}
2478
2487static int
2488NC4_strict_att_exists(NC_FILE_INFO_T *h5)
2489{
2490 hid_t grpid = -1;
2491 htri_t attr_exists;
2492
2493 /* Get root group ID. */
2494 grpid = ((NC_HDF5_GRP_INFO_T *)(h5->root_grp->format_grp_info))->hdf_grpid;
2495
2496 /* See if the NC3_STRICT_ATT_NAME attribute exists */
2497 if ((attr_exists = H5Aexists(grpid, NC3_STRICT_ATT_NAME)) < 0)
2498 return 1;
2499 return (attr_exists?1:0);
2500}
2501
2511static int
2512NC4_walk(hid_t gid, int* countp)
2513{
2514 int ncstat = NC_NOERR;
2515 int i,j,na;
2516 ssize_t len;
2517 hsize_t nobj;
2518 herr_t err;
2519 int otype;
2520 hid_t grpid, dsid;
2521 char name[NC_HDF5_MAX_NAME];
2522
2523 /* walk group members of interest */
2524 err = H5Gget_num_objs(gid, &nobj);
2525 if(err < 0) return err;
2526
2527 for(i = 0; i < nobj; i++) {
2528 /* Get name & kind of object in the group */
2529 len = H5Gget_objname_by_idx(gid,(hsize_t)i,name,(size_t)NC_HDF5_MAX_NAME);
2530 if(len < 0) return len;
2531
2532 otype = H5Gget_objtype_by_idx(gid,(size_t)i);
2533 switch(otype) {
2534 case H5G_GROUP:
2535 grpid = H5Gopen(gid,name);
2536 NC4_walk(grpid,countp);
2537 H5Gclose(grpid);
2538 break;
2539 case H5G_DATASET: /* variables */
2540 /* Check for phony_dim */
2541 if(strcmp(name,"phony_dim")==0)
2542 *countp = *countp + 1;
2543 dsid = H5Dopen(gid,name);
2544 na = H5Aget_num_attrs(dsid);
2545 for(j = 0; j < na; j++) {
2546 hid_t aid = H5Aopen_idx(dsid,(unsigned int) j);
2547 if(aid >= 0) {
2548 const NC_reservedatt* ra;
2549 ssize_t len = H5Aget_name(aid, NC_HDF5_MAX_NAME, name);
2550 if(len < 0) return len;
2551 /* Is this a netcdf-4 marker attribute */
2552 /* Is this a netcdf-4 marker attribute */
2553 ra = NC_findreserved(name);
2554 if(ra != NULL)
2555 *countp = *countp + 1;
2556 }
2557 H5Aclose(aid);
2558 }
2559 H5Dclose(dsid);
2560 break;
2561 default:/* ignore */
2562 break;
2563 }
2564 }
2565 return ncstat;
2566}
2567
EXTERNL int nc_free_vlen(nc_vlen_t *vl)
Free memory in a VLEN object.
Definition dvlen.c:41
int nc4_adjust_var_cache(NC_GRP_INFO_T *grp, NC_VAR_INFO_T *var)
Definition nc4hdf.c:1060
int nc4_rec_write_groups_types(NC_GRP_INFO_T *grp)
Definition nc4hdf.c:1965
static int var_exists(hid_t grpid, char *name, nc_bool_t *exists)
Definition nc4hdf.c:1423
void showopenobjects(int ncid)
Definition nc4hdf.c:2363
static int write_var(NC_VAR_INFO_T *var, NC_GRP_INFO_T *grp, nc_bool_t write_dimid)
Definition nc4hdf.c:1513
static int put_att_grpa(NC_GRP_INFO_T *grp, int varid, NC_ATT_INFO_T *att)
Definition nc4hdf.c:428
static int remove_coord_atts(hid_t hdf_datasetid)
Definition nc4hdf.c:1464
static int write_coord_dimids(NC_VAR_INFO_T *var)
Definition nc4hdf.c:648
static void reportopenobjectsT(int uselog, hid_t fid, int ntypes, unsigned int *otypes)
Definition nc4hdf.c:2288
static int write_netcdf4_dimid(hid_t datasetid, int dimid)
Definition nc4hdf.c:693
static int write_attlist(NCindex *attlist, int varid, NC_GRP_INFO_T *grp)
Definition nc4hdf.c:612
int nc4_get_hdf_typeid(NC_FILE_INFO_T *h5, nc_type xtype, hid_t *hdf_typeid, int endianness)
Definition nc4hdf.c:250
static int commit_type(NC_GRP_INFO_T *grp, NC_TYPE_INFO_T *type)
Definition nc4hdf.c:1116
static int NC4_walk(hid_t, int *)
Definition nc4hdf.c:2512
int rec_detach_scales(NC_GRP_INFO_T *grp, int dimid, hid_t dimscaleid)
Definition nc4hdf.c:143
void reportopenobjects(int uselog, hid_t fid)
Definition nc4hdf.c:2325
#define NC_HDF5_MAX_NAME
Definition nc4hdf.c:31
int rec_reattach_scales(NC_GRP_INFO_T *grp, int dimid, hid_t dimscaleid)
Definition nc4hdf.c:77
void showopenobjects5(NC_FILE_INFO_T *h5)
Definition nc4hdf.c:2341
static int write_nc3_strict_att(hid_t hdf_grpid)
Definition nc4hdf.c:1257
int nc4_rec_write_metadata(NC_GRP_INFO_T *grp, nc_bool_t bad_coord_order)
Definition nc4hdf.c:1879
static int var_create_dataset(NC_GRP_INFO_T *grp, NC_VAR_INFO_T *var, nc_bool_t write_dimid)
Definition nc4hdf.c:747
static int write_dim(NC_DIM_INFO_T *dim, NC_GRP_INFO_T *grp, nc_bool_t write_dimid)
Definition nc4hdf.c:1802
int NC4_hdf5get_libversion(unsigned *major, unsigned *minor, unsigned *release)
Definition nc4hdf.c:2387
static int create_group(NC_GRP_INFO_T *grp)
Definition nc4hdf.c:1302
int nc4_rec_match_dimscales(NC_GRP_INFO_T *grp)
Definition nc4hdf.c:2021
int nc4_create_dim_wo_var(NC_DIM_INFO_T *dim)
Definition nc4hdf.c:1711
static int flag_atts_dirty(NCindex *attlist)
Definition nc4hdf.c:42
static int NC4_strict_att_exists(NC_FILE_INFO_T *)
Definition nc4hdf.c:2488
int NC4_isnetcdf4(struct NC_FILE_INFO *h5)
Definition nc4hdf.c:2454
int NC4_hdf5get_superblock(struct NC_FILE_INFO *h5, int *idp)
Definition nc4hdf.c:2405
void reportobject(int uselog, hid_t id, unsigned int type)
Definition nc4hdf.c:2242
static int attach_dimscales(NC_GRP_INFO_T *grp)
Definition nc4hdf.c:1359
int nc4_open_var_grp2(NC_GRP_INFO_T *grp, int varid, hid_t *dataset)
Definition nc4hdf.c:203
const NC_reservedatt * NC_findreserved(const char *name)
int nc4_dim_list_add(NC_GRP_INFO_T *grp, const char *name, size_t len, int assignedid, NC_DIM_INFO_T **dim)
int nc4_find_type(const NC_FILE_INFO_T *h5, nc_type typeid, NC_TYPE_INFO_T **type)
int nc4_find_dim(NC_GRP_INFO_T *grp, int dimid, NC_DIM_INFO_T **dim, NC_GRP_INFO_T **dim_grp)
int nc4_find_nc_grp_h5(int ncid, NC **nc, NC_GRP_INFO_T **grp, NC_FILE_INFO_T **h5)
int nc4_get_fill_value(NC_FILE_INFO_T *h5, NC_VAR_INFO_T *var, void **fillp)
Definition nc4var.c:1264
Main header file for the C API.
#define NC_EBADTYPE
Not a netcdf data type.
Definition netcdf.h:372
#define NC_UINT
unsigned 4-byte int
Definition netcdf.h:44
#define NC_ENDIAN_NATIVE
In HDF5 files you can set the endianness of variables with nc_def_var_endian().
Definition netcdf.h:289
#define NC_EFILTER
Filter operation failed.
Definition netcdf.h:475
#define NC_INT
signed 4 byte integer
Definition netcdf.h:38
#define NC_ENDIAN_BIG
In HDF5 files you can set the endianness of variables with nc_def_var_endian().
Definition netcdf.h:291
#define NC_MAX_VAR_DIMS
max per variable dimensions
Definition netcdf.h:277
#define DIM_WITHOUT_VARIABLE
Definition netcdf.h:493
#define NC_BYTE
signed 1 byte integer
Definition netcdf.h:35
#define NC_EPERM
Write to read only.
Definition netcdf.h:341
#define NC_VLEN
vlen (variable-length) types
Definition netcdf.h:53
#define NC_NAT
Not A Type.
Definition netcdf.h:34
#define NC_DOUBLE
double precision floating point number
Definition netcdf.h:41
#define NC_UBYTE
unsigned 1 byte int
Definition netcdf.h:42
#define NC_FLOAT
single precision floating point number
Definition netcdf.h:40
#define NC_ENOMEM
Memory allocation (malloc) failure.
Definition netcdf.h:410
#define NC_COMPOUND
compound types
Definition netcdf.h:56
#define NC_CHUNKED
In HDF5 files you can set storage for each variable to be either contiguous or chunked,...
Definition netcdf.h:298
#define NC_SHORT
signed 2 byte integer
Definition netcdf.h:37
#define NC_ENUM
enum types
Definition netcdf.h:55
#define NC_INT64
signed 8-byte int
Definition netcdf.h:45
#define NC_EATTMETA
Problem with attribute metadata.
Definition netcdf.h:449
#define NC_EHDFERR
Error at HDF5 layer.
Definition netcdf.h:443
#define NC_CONTIGUOUS
In HDF5 files you can set storage for each variable to be either contiguous or chunked,...
Definition netcdf.h:299
#define NC_GLOBAL
Attribute id to put/get a global attribute.
Definition netcdf.h:249
#define NC_UINT64
unsigned 8-byte int
Definition netcdf.h:46
#define NC_COMPACT
In HDF5 files you can set storage for each variable to be either contiguous or chunked,...
Definition netcdf.h:300
#define NC_EFILEMETA
Problem with file metadata.
Definition netcdf.h:447
#define NC_ENOTVAR
Variable not found.
Definition netcdf.h:384
#define NC_EINVAL
Invalid Argument.
Definition netcdf.h:340
#define NC_CLASSIC_MODEL
Enforce classic model on netCDF-4.
Definition netcdf.h:139
#define NC_ENDIAN_LITTLE
In HDF5 files you can set the endianness of variables with nc_def_var_endian().
Definition netcdf.h:290
#define NC_MAX_NAME
Maximum for classic library.
Definition netcdf.h:276
#define NC_NOERR
No Error.
Definition netcdf.h:330
#define NC_EDIMMETA
Problem with dimension metadata.
Definition netcdf.h:448
#define NC_USHORT
unsigned 2-byte int
Definition netcdf.h:43
#define NC_OPAQUE
opaque types
Definition netcdf.h:54
#define NC_STRING
string
Definition netcdf.h:47
#define NC_CHAR
ISO/ASCII character.
Definition netcdf.h:36
#define NC_EVARMETA
Problem with variable metadata.
Definition netcdf.h:450
int nc_type
The nc_type type is just an int.
Definition netcdf.h:25
This is the type of arrays of vlens.
Definition netcdf.h:698