BELFEM 0.9.3
Berkeley Lab Finite Element Framework
Loading...
Searching...
No Matches
hdf5_tools.hpp
Go to the documentation of this file.
1/*
2 * BELFEM -- The Berkeley Lab Finite Element Framework
3 * Copyright (c) 2026, The Regents of the University of California,
4 * through Lawrence Berkeley National Laboratory (subject to receipt of any required
5 * approvals from the U.S. Dept. of Energy). All rights reserved.
6 *
7 * Developers: Christian Messe, Gregory Giard
8 *
9 * See the top-level LICENSE file for the complete license and disclaimer.
10 */
11
12#ifndef BELFEM_HDF5_TOOLS_HPP
13#define BELFEM_HDF5_TOOLS_HPP
14
15#include <cstring>
16
17#include "stringtools.hpp"
18#include "assert.hpp"
19#include "hdf5_types.hpp"
20#include "cl_Vector.hpp"
21#include "cl_Matrix.hpp"
22#include "cl_Cell.hpp"
23
24namespace belfem
25{
26 namespace hdf5
27 {
28//------------------------------------------------------------------------------
29
30#ifdef BELFEM_HDF5
38 struct DatasetHandles
39 {
40 hid_t mDataSet = H5I_INVALID_HID ;
41 hid_t mDataType = H5I_INVALID_HID ;
42 hid_t mDataSpace = H5I_INVALID_HID ;
43
44 DatasetHandles() = default ;
45
46 // a copy or a move would close the same identifiers twice
47 DatasetHandles( const DatasetHandles & ) = delete ;
48 DatasetHandles & operator=( const DatasetHandles & ) = delete ;
49 DatasetHandles( DatasetHandles && ) = delete ;
50 DatasetHandles & operator=( DatasetHandles && ) = delete ;
51
52 ~DatasetHandles()
53 {
54 if( mDataSpace >= 0 )
55 {
56 H5Sclose( mDataSpace );
57 }
58 if( mDataType >= 0 )
59 {
60 H5Tclose( mDataType );
61 }
62 if( mDataSet >= 0 )
63 {
64 H5Dclose( mDataSet );
65 }
66 }
67 };
68#endif
69
70//------------------------------------------------------------------------------
71
76 inline bool
77 dataset_exists( hid_t aLoc, const std::string & aLabel )
78 {
79#ifdef BELFEM_HDF5
80 hid_t tDataSet = 0;
81 return H5Lexists( aLoc, aLabel.c_str(), tDataSet );
82#else
83 return false;
84#endif
85 }
86
87//------------------------------------------------------------------------------
88
93 inline bool
94 group_exists( hid_t aLoc, const std::string & aLabel )
95 {
96#ifdef BELFEM_HDF5
97 hid_t tGroup = 0;
98 return H5Lexists( aLoc, aLabel.c_str(), tGroup );
99#else
100 return false;
101#endif
102 }
103
104//------------------------------------------------------------------------------
105
106// The whole function is guarded, not just its body: its SIGNATURE names
107// H5L_info2_t, which hdf5_types.hpp does not stub. With the body alone guarded
108// this header stopped compiling under -DUSE_HDF5=OFF, and since cl_HDF5.hpp
109// includes it and src/io always compiles cl_HDF5.cpp, that broke the whole
110// tree rather than only the HDF5 callers
111#ifdef BELFEM_HDF5
112
113// collect_names() and get_groups() below are written against the 1.12 link and
114// object interfaces ( H5L_info2_t, H5O_info2_t, H5Literate2,
115// H5Oget_info_by_name3 ). belfem_find_package carries no version test, so the
116// requirement is stated here, where it is created -- and this header is
117// included across the whole tree, so on an older HDF5 the failure would
118// otherwise be a wall of template errors in every module rather than one line
119#if ! H5_VERSION_GE( 1, 12, 0 )
120#error "BELFEM requires HDF5 1.12 or newer ( hdf5_tools.hpp uses H5Literate2 and H5Oget_info_by_name3 )."
121#endif
122
123 inline herr_t
124 collect_names(
125 hid_t aLoc,
126 const char * aLabel,
127 const H5L_info2_t * aInfo,
128 void * aData )
129 {
130 auto * tLabels = static_cast<std::vector<std::string> *>( aData );
131
132 // GROUPS only. H5Literate2 visits every link at the location,
133 // datasets and named datatypes included, and a caller that opens
134 // each returned name as a group fails on the first one that is not.
135 // A soft or external link is skipped as well rather than followed:
136 // resolving it can leave the location, which is not what listing a
137 // location means
138 if( aInfo != nullptr && aInfo->type != H5L_TYPE_HARD )
139 {
140 return 0;
141 }
142
143 H5O_info2_t tInfo ;
144
145 // a failed query is an ERROR, not an entry to skip. Returning 0
146 // here would let a permission failure or a corrupt object produce a
147 // short list that the caller cannot tell from a complete one. The
148 // callback contract is: negative aborts the iteration and becomes
149 // H5Literate2's return value, which get_groups() checks
150 if( H5Oget_info_by_name3( aLoc, aLabel, & tInfo,
151 H5O_INFO_BASIC, H5P_DEFAULT ) < 0 )
152 {
153 return -1;
154 }
155
156 // not a group is not an error -- it is the thing this filter exists
157 // to drop, so the iteration continues
158 if( tInfo.type != H5O_TYPE_GROUP )
159 {
160 return 0;
161 }
162
163 tLabels->emplace_back( aLabel );
164
165 return 0;
166 }
167
168#endif
169
173 inline Cell< std::string >
175 {
176#ifdef BELFEM_HDF5
177 Cell< std::string > tGroups ;
178
179 hsize_t tIndex = 0 ;
180
181 // H5_ITER_INC, not H5_ITER_NATIVE : native order is documented as
182 // whatever the library finds fastest, and callers that treat the
183 // first entry as a reference ( db2exo picks the grid of the first
184 // table it sees ) would then depend on it
185 const herr_t tStatus = H5Literate2(
186 aLoc,
187 H5_INDEX_NAME,
188 H5_ITER_INC,
189 & tIndex,
190 collect_names,
191 & tGroups.vector_data() );
192
193 // runs once per file, so the always-active tier is affordable and
194 // the message survives a release build. Silently returning a short
195 // list here would be indistinguishable from a location that really
196 // holds fewer groups
197 BELFEM_ERROR( tStatus >= 0,
198 "H5Literate2 failed while listing the groups of an HDF5 location." );
199
200 return tGroups ;
201
202#else
203 return {};
204#endif
205 }
206
207//------------------------------------------------------------------------------
208
213 template < typename T >
214 bool
216 {
217#ifdef BELFEM_HDF5
218 return H5Tget_size( datatype<T>() ) == sizeof( T );
219#else
220 return false;
221#endif
222 }
223
224//------------------------------------------------------------------------------
225
246 template < typename T >
247 void
248 check_read_datatype( const hid_t aFileType, const std::string & aLabel )
249 {
250#ifdef BELFEM_HDF5
251 BELFEM_ERROR( H5Tget_class( aFileType ) == H5Tget_class( datatype<T>() ),
252 "datatype class mismatch reading '%s': the file holds class %i, but %s is class %i.\n"
253 "HDF5 can convert between widths, not between classes - this file does not hold what the reader expects.",
254 aLabel.c_str(),
255 ( int ) H5Tget_class( aFileType ),
256 datatype_string<T>().c_str(),
257 ( int ) H5Tget_class( datatype<T>() ) );
258#endif
259 }
260
261//------------------------------------------------------------------------------
262
272 template < typename T >
273 void
275 hid_t & aLoc,
276 const std::string & aLabel,
277 const T & aValue,
278 herr_t & aStatus
279 )
280 {
281#ifdef BELFEM_HDF5
282 BELFEM_ERROR( ! hdf5::dataset_exists( aLoc, aLabel ),
283 "Dataset %s of type %s does already exist.",
284 aLabel.c_str(),
285 datatype_string<T>().c_str() );
286
288 "Error in datatype size of type %s.",
289 datatype_string<T>().c_str() );
290
291 // FILE type: fixed width and little endian, so the on-disk
292 // layout does not depend on this build's native type widths
293 hid_t tDataType = H5Tcopy( filetype<T>() );
294
295 hsize_t tDims[ 1 ] = { 1 };
296
297 hid_t tDataSpace
298 = H5Screate_simple( 1, tDims, nullptr );
299
300 hid_t tDataSet = H5Dcreate(
301 aLoc,
302 aLabel.c_str(),
303 tDataType,
304 tDataSpace,
305 H5P_DEFAULT,
306 H5P_DEFAULT,
307 H5P_DEFAULT );
308
309 aStatus = H5Dwrite(
310 tDataSet,
311 datatype< T >(), // MEMORY type: native, not the file type
312 H5S_ALL,
313 H5S_ALL,
314 H5P_DEFAULT,
315 & aValue );
316
317 H5Sclose( tDataSpace );
318 H5Tclose( tDataType );
319 H5Dclose( tDataSet );
320
321 BELFEM_ASSERT( aStatus == 0,
322 "Something went wrong while trying to store scalar %s of type %s",
323 aLabel.c_str(),
324 datatype_string<T>().c_str());
325#endif
326 }
327//------------------------------------------------------------------------------
328
329 template < typename T >
330 void
332 hid_t & aLoc,
333 const std::string & aLabel,
334 T & aValue,
335 herr_t & aStatus
336 )
337 {
338#ifdef BELFEM_HDF5
339 BELFEM_ERROR( hdf5::dataset_exists( aLoc, aLabel ),
340 "Dataset %s of type %s does not exist.",
341 aLabel.c_str(),
342 datatype_string<T>().c_str() );
343
345 "Error in datatype size of type %s.",
346 datatype_string<T>().c_str() );
347
348 DatasetHandles tHandles ;
349
350 tHandles.mDataSet = H5Dopen1( aLoc, aLabel.c_str() );
351
352 tHandles.mDataType = H5Dget_type( tHandles.mDataSet );
353
354 // the file type must be convertible into T; the read
355 // below then uses datatype<T>() so HDF5 does the conversion
356 check_read_datatype< T >( tHandles.mDataType, aLabel );
357
358 tHandles.mDataSpace = H5Dget_space( tHandles.mDataSet );
359
360 aStatus = H5Dread(
361 tHandles.mDataSet,
362 datatype< T >(), // MEMORY type: let HDF5 convert
363 H5S_ALL,
364 H5S_ALL,
365 H5P_DEFAULT,
366 &aValue );
367
368 BELFEM_ASSERT( aStatus == 0,
369 "Something went wrong while trying to load scalar %s of type %s",
370 aLabel.c_str(),
371 datatype_string<T>().c_str() );
372
373#endif
374 }
375
376//------------------------------------------------------------------------------
377
378 void
380 hid_t & aLoc,
381 const std::string & aLabel,
382 const bool & aValue,
383 herr_t & aStatus );
384
385//------------------------------------------------------------------------------
386
387 void
389 hid_t & aLoc,
390 const std::string & aLabel,
391 bool & aValue,
392 herr_t & aStatus );
393
394//------------------------------------------------------------------------------
395
396 void
398 hid_t & aLoc,
399 const std::string & aLabel,
400 const std::string & aValue,
401 herr_t & aStatus );
402
403//------------------------------------------------------------------------------
404
405 void
407 hid_t & aLoc,
408 const std::string & aLabel,
409 std::string & aValue,
410 herr_t & aStatus );
411
412//------------------------------------------------------------------------------
413 template < typename T >
414 void
416 hid_t & aLoc,
417 const std::string & aLabel,
418 const T * aData,
419 const hsize_t & aLength,
420 herr_t & aStatus )
421 {
422#ifdef BELFEM_HDF5
423 BELFEM_ERROR( ! hdf5::dataset_exists( aLoc, aLabel ),
424 "Dataset %s of type Vector<%s> does already exist.",
425 aLabel.c_str(),
426 datatype_string<T>().c_str() );
427
429 "Error in datatype size of type %s.",
430 datatype_string<T>().c_str() );
431
432 // FILE type: fixed width and little endian, so the on-disk
433 // layout does not depend on this build's native type widths
434 hid_t tDataType = H5Tcopy( filetype<T>() );
435
436 hsize_t tDims[ 1 ];
437 tDims[ 0 ] = aLength;
438
439 hid_t tDataSpace
440 = H5Screate_simple( 1, tDims, nullptr );
441
442 hid_t tDataSet = H5Dcreate(
443 aLoc,
444 aLabel.c_str(),
445 tDataType,
446 tDataSpace,
447 H5P_DEFAULT,
448 H5P_DEFAULT,
449 H5P_DEFAULT );
450
451 if( aLength > 0 )
452 {
453 aStatus = H5Dwrite(
454 tDataSet,
455 datatype< T >(), // MEMORY type: native, not the file type
456 H5S_ALL,
457 H5S_ALL,
458 H5P_DEFAULT,
459 &aData[ 0 ] );
460 }
461
462 H5Sclose( tDataSpace );
463 H5Tclose( tDataType );
464 H5Dclose( tDataSet );
465
466 BELFEM_ASSERT( aStatus == 0,
467 "Something went wrong while trying to store vector %s of type %s",
468 aLabel.c_str(),
469 datatype_string<T>().c_str() );
470#endif
471 }
472
473//------------------------------------------------------------------------------
474
475 hsize_t
477 hid_t & aLoc,
478 const std::string & aLabel );
479
480//------------------------------------------------------------------------------
481
482 template < typename T >
483 void
485 hid_t & aLoc,
486 const std::string & aLabel,
487 T * aData,
488 const hsize_t aLength,
489 herr_t & aStatus )
490 {
491#ifdef BELFEM_HDF5
492
493 BELFEM_ERROR( hdf5::dataset_exists( aLoc, aLabel ),
494 "Dataset %s of type Vector<%s> does not exist.",
495 aLabel.c_str(),
496 datatype_string<T>().c_str() );
497
499 "Error in datatype size of type %s.",
500 datatype_string<T>().c_str() );
501
502 DatasetHandles tHandles ;
503
504 tHandles.mDataSet = H5Dopen1( aLoc, aLabel.c_str() );
505
506 tHandles.mDataType = H5Dget_type( tHandles.mDataSet );
507
508 // the file type must be convertible into T; the read
509 // below then uses datatype<T>() so HDF5 does the conversion
510 check_read_datatype< T >( tHandles.mDataType, aLabel );
511
512 tHandles.mDataSpace = H5Dget_space( tHandles.mDataSet );
513
514 hsize_t tDims[ 1 ];
515
516 // the extent buffer holds one entry; a dataset of another rank is
517 // refused before the extents are read into it
518 const int tRank = H5Sget_simple_extent_ndims( tHandles.mDataSpace );
519
520 BELFEM_ERROR( tRank == 1, "Dataset %s is not a 1-D array (rank %i)",
521 aLabel.c_str(), tRank );
522
523 aStatus = H5Sget_simple_extent_dims( tHandles.mDataSpace, tDims, nullptr );
524
525 BELFEM_ERROR( tDims[ 0 ] == aLength, "Lengths do not match: is %lu, expect %lu.",
526 ( long unsigned int ) tDims[ 0 ], ( long unsigned int ) aLength );
527
528 if( aLength > 0 )
529 {
530 aStatus = H5Dread(
531 tHandles.mDataSet,
532 datatype< T >(), // MEMORY type: let HDF5 convert
533 H5S_ALL,
534 H5S_ALL,
535 H5P_DEFAULT,
536 &aData[ 0 ] );
537 }
538 else
539 {
540 aStatus = 0;
541 }
542
543 BELFEM_ASSERT( aStatus == 0,
544 "Something went wrong while trying to load vector %s of type %s",
545 aLabel.c_str(),
546 datatype_string<T>().c_str() );
547#endif
548 }
549//------------------------------------------------------------------------------
550
551 template < typename T >
552 void
554 hid_t & aLoc,
555 const std::string & aLabel,
556 const Vector< T > & aVector,
557 herr_t & aStatus )
558 {
559#ifdef BELFEM_HDF5
560
561 BELFEM_ERROR( ! hdf5::dataset_exists( aLoc, aLabel ),
562 "Dataset %s of type Vector<%s> does already exist.",
563 aLabel.c_str(),
564 datatype_string<T>().c_str() );
565
567 "Error in datatype size of type %s.",
568 datatype_string<T>().c_str() );
569
570 // FILE type: fixed width and little endian, so the on-disk
571 // layout does not depend on this build's native type widths
572 hid_t tDataType = H5Tcopy( filetype<T>() );
573 hid_t tDataSet = 0;
574
575 hsize_t tLength = aVector.length();
576 hsize_t tDims[ 1 ];
577 tDims[ 0 ] = tLength;
578
579 hid_t tDataSpace
580 = H5Screate_simple( 1, tDims, nullptr );
581
582 tDataSet = H5Dcreate(
583 aLoc,
584 aLabel.c_str(),
585 tDataType,
586 tDataSpace,
587 H5P_DEFAULT,
588 H5P_DEFAULT,
589 H5P_DEFAULT );
590
591 if( tLength > 0 )
592 {
593 T* tData = ( T* ) malloc( tLength * sizeof( T ) );
594
595 std::memcpy(
596 tData,
597 aVector.data(),
598 tLength * sizeof( T ) );
599
600 aStatus = H5Dwrite(
601 tDataSet,
602 datatype< T >(), // MEMORY type: native, not the file type
603 H5S_ALL,
604 H5S_ALL,
605 H5P_DEFAULT,
606 & tData[ 0 ]);
607
608 free( tData );
609 }
610
611 H5Sclose( tDataSpace );
612 H5Tclose( tDataType );
613 H5Dclose( tDataSet );
614
615 BELFEM_ASSERT( aStatus == 0,
616 "Something went wrong while trying to store vector %s of type %s",
617 aLabel.c_str(),
618 datatype_string<T>().c_str() );
619#endif
620 }
621
622//------------------------------------------------------------------------------
623
624 template < typename T >
625 void
627 hid_t & aLoc,
628 const std::string & aLabel,
629 Vector< T > & aVector,
630 herr_t & aStatus )
631 {
632#ifdef BELFEM_HDF5
633 BELFEM_ERROR( hdf5::dataset_exists( aLoc, aLabel ),
634 "Dataset %s of type Vector<%s> does not exist.",
635 aLabel.c_str(),
636 datatype_string<T>().c_str() );
637
639 "Error in datatype size of type %s.",
640 datatype_string<T>().c_str() );
641
642 DatasetHandles tHandles ;
643
644 tHandles.mDataSet = H5Dopen1( aLoc, aLabel.c_str() );
645
646 tHandles.mDataType = H5Dget_type( tHandles.mDataSet );
647
648 // the file type must be convertible into T; the read
649 // below then uses datatype<T>() so HDF5 does the conversion
650 check_read_datatype< T >( tHandles.mDataType, aLabel );
651
652 tHandles.mDataSpace = H5Dget_space( tHandles.mDataSet );
653
654 hsize_t tDims[ 1 ];
655
656 // the extent buffer holds one entry; a dataset of another rank is
657 // refused before the extents are read into it
658 const int tRank = H5Sget_simple_extent_ndims( tHandles.mDataSpace );
659
660 BELFEM_ERROR( tRank == 1, "Dataset %s is not a 1-D array (rank %i)",
661 aLabel.c_str(), tRank );
662
663 aStatus = H5Sget_simple_extent_dims( tHandles.mDataSpace, tDims, nullptr );
664
665 hsize_t tLength = tDims[ 0 ];
666
667 aVector.set_size( tLength );
668
669 if( tLength > 0 )
670 {
671 T* tData = ( T* ) malloc( tLength * sizeof( T ) );
672
673 aStatus = H5Dread(
674 tHandles.mDataSet,
675 datatype< T >(), // MEMORY type: let HDF5 convert
676 H5S_ALL,
677 H5S_ALL,
678 H5P_DEFAULT,
679 &tData[ 0 ] );
680
681 std::memcpy(
682 aVector.data(),
683 tData,
684 tLength * sizeof( T ) );
685
686 free( tData );
687 }
688 else
689 {
690 aStatus = 0;
691 }
692
693 BELFEM_ASSERT( aStatus == 0,
694 "Something went wrong while trying to load vector %s of type %s",
695 aLabel.c_str(),
696 datatype_string<T>().c_str() );
697#endif
698 }
699//------------------------------------------------------------------------------
700
701 template < typename T >
702 void
704 hid_t & aLoc,
705 const std::string & aLabel,
706 const Matrix< T > & aMatrix,
707 herr_t & aStatus,
708 const bool aTranspose = false )
709 {
710#ifdef BELFEM_HDF5
711 BELFEM_ERROR( ! hdf5::dataset_exists( aLoc, aLabel ),
712 "Dataset %s of type Matrix<%s> does already exist.",
713 aLabel.c_str(),
714 datatype_string<T>().c_str() );
715
717 "Error in datatype size of type %s.",
718 datatype_string<T>().c_str() );
719
720 // FILE type: fixed width and little endian, so the on-disk
721 // layout does not depend on this build's native type widths
722 hid_t tDataType = H5Tcopy( filetype<T>() );
723 hid_t tDataSet = 0;
724
725 hsize_t tDims[ 2 ];
726
727 if ( aTranspose )
728 {
729 tDims[ 0 ] = aMatrix.n_cols();
730 tDims[ 1 ] = aMatrix.n_rows();
731
732 }
733 else
734 {
735 tDims[ 0 ] = aMatrix.n_rows();
736 tDims[ 1 ] = aMatrix.n_cols();
737 }
738
739 hid_t tDataSpace
740 = H5Screate_simple( 2, tDims, nullptr );
741
742 tDataSet = H5Dcreate(
743 aLoc,
744 aLabel.c_str(),
745 tDataType,
746 tDataSpace,
747 H5P_DEFAULT,
748 H5P_DEFAULT,
749 H5P_DEFAULT );
750
751 if( tDims[ 0 ]*tDims[ 1 ] > 0 )
752 {
753 T** tData = ( T** ) malloc( tDims[ 0 ] * sizeof( T * ) );
754 tData[ 0 ] = ( T* ) malloc( tDims[ 0 ]*tDims[ 1 ] * sizeof( T ) );
755
756 for( hsize_t i=0; i<tDims[ 0 ]; ++i )
757 {
758 tData[ i ] = tData[ 0 ]+ i*tDims[ 1 ];
759 }
760
761 if ( aTranspose )
762 {
763 for ( hsize_t i = 0; i < tDims[ 0 ]; ++i )
764 {
765 for ( hsize_t j = 0; j < tDims[ 1 ]; ++j )
766 {
767 tData[ i ][ j ] = aMatrix( j, i );
768 }
769 }
770 }
771 else
772 {
773 for ( hsize_t i = 0; i < tDims[ 0 ]; ++i )
774 {
775 for ( hsize_t j = 0; j < tDims[ 1 ]; ++j )
776 {
777 tData[ i ][ j ] = aMatrix( i, j );
778 }
779 }
780 }
781 aStatus = H5Dwrite(
782 tDataSet,
783 datatype< T >(), // MEMORY type: native, not the file type
784 H5S_ALL,
785 H5S_ALL,
786 H5P_DEFAULT,
787 & tData[ 0 ][ 0 ]);
788
789 free( tData[ 0 ] );
790 free( tData );
791 }
792
793 H5Sclose( tDataSpace );
794 H5Tclose( tDataType );
795 H5Dclose( tDataSet );
796
797 BELFEM_ASSERT( aStatus == 0,
798 "Something went wrong while trying to store matrix %s of type %s",
799 aLabel.c_str(),
800 datatype_string<T>().c_str() );
801#endif
802 }
803
804//------------------------------------------------------------------------------
805
806 template < typename T >
807 void
809 hid_t & aLoc,
810 const std::string & aLabel,
811 Matrix< T > & aMatrix,
812 herr_t & aStatus,
813 const bool aTranspose = false )
814 {
815#ifdef BELFEM_HDF5
816
817 BELFEM_ERROR( hdf5::dataset_exists( aLoc, aLabel ),
818 "Dataset %s of type Matrix<%s> does not exist.",
819 aLabel.c_str(),
820 datatype_string<T>().c_str() );
821
823 "Error in datatype size of type %s.",
824 datatype_string<T>().c_str() );
825
826 DatasetHandles tHandles ;
827
828 tHandles.mDataSet = H5Dopen1( aLoc, aLabel.c_str() );
829
830 tHandles.mDataType = H5Dget_type( tHandles.mDataSet );
831
832 // the file type must be convertible into T; the read
833 // below then uses datatype<T>() so HDF5 does the conversion
834 check_read_datatype< T >( tHandles.mDataType, aLabel );
835
836 tHandles.mDataSpace = H5Dget_space( tHandles.mDataSet );
837
838 hsize_t tDims[ 2 ];
839
840 aStatus = H5Sget_simple_extent_dims( tHandles.mDataSpace, tDims, nullptr );
841
842 // allocate memory for output (transpose swaps the on-disk dimensions)
843 if ( aTranspose )
844 {
845 aMatrix.set_size( tDims[ 1 ], tDims[ 0 ] );
846 }
847 else
848 {
849 aMatrix.set_size( tDims[ 0 ], tDims[ 1 ] );
850 }
851
852 if( tDims[ 0 ]*tDims[ 1 ] > 0 )
853 {
854 T** tData = ( T** )
855 malloc( tDims[ 0 ]*sizeof( T* ) );
856
857 tData[ 0 ] = ( T* )
858 malloc( tDims[ 0 ]* tDims[ 1 ] * sizeof( T ) );
859
860 for( hsize_t i=1; i<tDims[ 0 ]; ++i )
861 {
862 tData[ i ] = tData[ 0 ]+ i*tDims[ 1 ];
863 }
864
865 aStatus = H5Dread(
866 tHandles.mDataSet,
867 datatype< T >(), // MEMORY type: let HDF5 convert
868 H5S_ALL,
869 H5S_ALL,
870 H5P_DEFAULT,
871 &tData[ 0 ][ 0 ] );
872
873 if ( aTranspose )
874 {
875 // on-disk (tDims[0] x tDims[1]) -> memory (tDims[1] x tDims[0])
876 for ( hsize_t i = 0; i < tDims[ 0 ]; ++i )
877 {
878 for ( hsize_t j = 0; j < tDims[ 1 ]; ++j )
879 {
880 aMatrix( j, i ) = tData[ i ][ j ];
881 }
882 }
883 }
884 else
885 {
886 for ( hsize_t j = 0; j < tDims[ 1 ]; ++j )
887 {
888 for ( hsize_t i = 0; i < tDims[ 0 ]; ++i )
889 {
890 aMatrix( i, j ) = tData[ i ][ j ];
891 }
892 }
893 }
894 free( tData[ 0 ] );
895 free( tData );
896 }
897 else if( aStatus == 2 )
898 {
899 aStatus = 0;
900 }
901
902 BELFEM_ASSERT( aStatus == 0,
903 "Something went wrong while trying to load matrix %s of type %s",
904 aLabel.c_str(),
905 datatype_string<T>().c_str() );
906#endif
907 }
908
909 inline void
911 const std::string & aLabel,
912 const Cell< string > & aCell,
913 herr_t & aStatus )
914 {
915#ifdef BELFEM_HDF5
916
917 DatasetHandles tHandles ;
918
919 hsize_t tNumStrings = aCell.size();
920 hsize_t tDim[ 1 ] = { tNumStrings };
921 tHandles.mDataSpace = H5Screate_simple( 1, tDim, NULL );
922
923 tHandles.mDataType = H5Tcopy( H5T_C_S1 );
924 H5Tset_size( tHandles.mDataType, H5T_VARIABLE );
925
926 tHandles.mDataSet = H5Dcreate2(
927 aLoc,
928 aLabel.c_str(),
929 tHandles.mDataType,
930 tHandles.mDataSpace,
931 H5P_DEFAULT,
932 H5P_DEFAULT,
933 H5P_DEFAULT );
934
935 BELFEM_ERROR( tHandles.mDataSet > -1, "Error creating dataset" );
936
937 Cell< const char * > tBuffer( tNumStrings, nullptr );
938 for ( hsize_t i=0; i<tNumStrings; ++i )
939 {
940 tBuffer( i ) = aCell( i ).c_str();
941 }
942
943 aStatus = H5Dwrite(
944 tHandles.mDataSet,
945 tHandles.mDataType,
946 H5S_ALL,
947 H5S_ALL,
948 H5P_DEFAULT,
949 tBuffer.data() );
950
951 BELFEM_ERROR( aStatus ==0, "Error writing data" );
952#endif
953 }
954
955 inline void
957 const std::string & aLabel,
958 Cell< string > & aCell,
959 herr_t & aStatus )
960 {
961#ifdef BELFEM_HDF5
962 aCell.clear();
963
964 BELFEM_ERROR( hdf5::dataset_exists( aLoc, aLabel ),
965 "Dataset %s does not exist at expected location.", aLabel.c_str() );
966
967 DatasetHandles tHandles ;
968
969 tHandles.mDataSet = H5Dopen1( aLoc, aLabel.c_str() );
970
971 tHandles.mDataSpace = H5Dget_space( tHandles.mDataSet );
972
973 hsize_t tDim[ 1 ] ;
974
975 // the extent buffer holds one entry; a dataset of another rank is
976 // refused before the extents are read into it
977 const int tRank = H5Sget_simple_extent_ndims( tHandles.mDataSpace );
978
979 BELFEM_ERROR( tRank == 1, "Dataset %s is not a 1-D array (rank %i)",
980 aLabel.c_str(), tRank );
981
982 BELFEM_ERROR( H5Sget_simple_extent_dims( tHandles.mDataSpace, tDim, nullptr ) > -1 ,
983 "Error getting dataspace dimensions" );
984
985 hsize_t tNumStrings = tDim[ 0 ];
986
987 // This must match the datatype used in writing -- INCLUDING the
988 // character set: H5T_C_S1 defaults to ASCII, but h5py writes
989 // vlen strings as UTF-8, and HDF5 refuses the UTF-8 -> ASCII
990 // conversion path, failing the H5Dread below with no better
991 // diagnostic than a nonzero status. Copy the cset from the
992 // dataset itself so BELFEM-written ASCII files behave as before
993 // and Python-written UTF-8 files become readable.
994 tHandles.mDataType = H5Tcopy( H5T_C_S1 );
995 H5Tset_size( tHandles.mDataType, H5T_VARIABLE );
996
997 hid_t tFileType = H5Dget_type( tHandles.mDataSet );
998 H5Tset_cset( tHandles.mDataType, H5Tget_cset( tFileType ) );
999 H5Tclose( tFileType );
1000
1001 // HDF5 allocates the strings during the read. The null slots let
1002 // the reclaim below free exactly what a failed read left allocated,
1003 // so it runs on both paths, before the refusal and before the owner
1004 // closes the type and the dataspace it needs. An empty Cell may
1005 // hand out a null data(), which H5Dread accepts at zero elements
1006 // and H5Dvlen_reclaim refuses, hence the size guard.
1007 Cell< char * > tBuffer( tNumStrings, nullptr );
1008
1009 const herr_t tReadStatus = H5Dread( tHandles.mDataSet, tHandles.mDataType,
1010 H5S_ALL, H5S_ALL, H5P_DEFAULT, tBuffer.data() );
1011
1012 if( tReadStatus == 0 )
1013 {
1014 aCell.set_size( tNumStrings, "" );
1015
1016 for ( hsize_t i=0; i<tNumStrings; ++i )
1017 {
1018 aCell( i ) = std::string( tBuffer( i ) );
1019 }
1020 }
1021
1022 if( tNumStrings > 0 )
1023 {
1024 H5Dvlen_reclaim( tHandles.mDataType, tHandles.mDataSpace, H5P_DEFAULT, tBuffer.data() );
1025 }
1026
1027 BELFEM_ERROR( tReadStatus == 0, "Failed to read strings from file." );
1028#endif
1029 }
1030//------------------------------------------------------------------------------
1031
1032 inline hid_t
1033 open_group( const std::string & aLabel, hid_t aParent )
1034 {
1035#ifdef BELFEM_HDF5
1036 if ( hdf5::group_exists( aParent, aLabel ) )
1037 {
1038
1039 return H5Gopen2(
1040 aParent,
1041 aLabel.c_str(),
1042 H5P_DEFAULT );
1043 }
1044 else
1045 {
1046 BELFEM_ERROR( false,
1047 "Group %s does not exist.",
1048 aLabel.c_str() );
1049
1050 return -1;
1051 }
1052#else
1053 return 0 ;
1054#endif
1055 }
1056
1057//------------------------------------------------------------------------------
1058
1059 inline herr_t
1061 {
1062#ifdef BELFEM_HDF5
1063 return H5Gclose( aGroup );
1064#else
1065 return 0 ;
1066#endif
1067 }
1068
1069//------------------------------------------------------------------------------
1070
1071 string
1072 create_tree( const Cell< string > & aTree, const string & aLabel );
1073
1074//------------------------------------------------------------------------------
1075
1076//------------------------------------------------------------------------------
1077 } /* namespace hdf5 */
1078} /* namespace belfem */
1079
1080#endif //BELFEM_HDF5_TOOLS_HPP
#define BELFEM_ERROR(aCheck,...)
Definition assert.hpp:264
#define BELFEM_ASSERT(aCheck,...)
Definition assert.hpp:244
A bounds-checked wrapper around std::vector.
Definition cl_Cell.hpp:45
void set_size(const size_t aSize)
New elements are value-initialized ( zero for arithmetic T, nullptr for pointers ).
Definition cl_Cell.hpp:199
size_t size() const
Definition cl_Cell.hpp:189
void clear()
Empties the Cell.
Definition cl_Cell.hpp:249
std::vector< T > & vector_data()
The underlying vector.
Definition cl_Cell.hpp:146
T * data()
Returns a borrowed pointer into contiguous storage, possibly nullptr when the Cell is empty ( test si...
Definition cl_Cell.hpp:129
Dense column-major matrix.
Definition cl_BZ_Matrix.hpp:28
void set_size(const size_t aNumRows, const size_t aNumCols)
Definition cl_AR_Matrix.hpp:164
size_t n_rows() const
Definition cl_AR_Matrix.hpp:183
size_t n_cols() const
Definition cl_AR_Matrix.hpp:191
Column vector.
Definition cl_BZ_Vector.hpp:41
T * data()
Returns a borrowed raw buffer under the layout contract stated in cl_Vector.hpp; may be nullptr when ...
Definition cl_AR_Vector.hpp:163
void set_size(const size_t aNumRows)
change the size of the vector
Definition cl_AR_Vector.hpp:215
size_t length() const
get the length of the vector
Definition cl_AR_Vector.hpp:235
Definition cl_HDF5_Dataset.hpp:21
void save_vector_to_file(hid_t &aLoc, const std::string &aLabel, const Vector< T > &aVector, herr_t &aStatus)
Definition hdf5_tools.hpp:553
bool group_exists(hid_t aLoc, const std::string &aLabel)
test if a group exists ( tests for a link of that name; the link type is not checked )
Definition hdf5_tools.hpp:94
void save_string_to_file(hid_t &aFileID, const std::string &aLabel, const std::string &aValue, herr_t &aStatus)
Definition hdf5_tools.cpp:121
void load_scalar_from_file(hid_t &aLoc, const std::string &aLabel, T &aValue, herr_t &aStatus)
Definition hdf5_tools.hpp:331
void save_array_to_file(hid_t &aLoc, const std::string &aLabel, const T *aData, const hsize_t &aLength, herr_t &aStatus)
Definition hdf5_tools.hpp:415
void load_bool_from_file(hid_t &aFileID, const std::string &aLabel, bool &aValue, herr_t &aStatus)
Definition hdf5_tools.cpp:75
void save_strings_to_file(hid_t &aLoc, const std::string &aLabel, const Cell< string > &aCell, herr_t &aStatus)
Definition hdf5_tools.hpp:910
void save_scalar_to_file(hid_t &aLoc, const std::string &aLabel, const T &aValue, herr_t &aStatus)
saves a scalar value to a file file must be open
Definition hdf5_tools.hpp:274
hid_t open_group(const std::string &aLabel, hid_t aParent)
Definition hdf5_tools.hpp:1033
void load_strings_from_file(hid_t &aLoc, const std::string &aLabel, Cell< string > &aCell, herr_t &aStatus)
Definition hdf5_tools.hpp:956
void load_matrix_from_file(hid_t &aLoc, const std::string &aLabel, Matrix< T > &aMatrix, herr_t &aStatus, const bool aTranspose=false)
Definition hdf5_tools.hpp:808
void save_matrix_to_file(hid_t &aLoc, const std::string &aLabel, const Matrix< T > &aMatrix, herr_t &aStatus, const bool aTranspose=false)
Definition hdf5_tools.hpp:703
bool check_datatye_size()
this function returns true of both the HDF5 datatype and the passed datatype have the same size
Definition hdf5_tools.hpp:215
hid_t datatype()
Definition hdf5_types.hpp:36
Cell< std::string > get_groups(hid_t aLoc)
list groups that exist in location
Definition hdf5_tools.hpp:174
hid_t filetype()
Definition hdf5_types.hpp:56
hsize_t get_array_size(hid_t &aFileID, const std::string &aLabel)
Definition hdf5_tools.cpp:223
bool dataset_exists(hid_t aLoc, const std::string &aLabel)
test if a dataset exists ( tests for a link of that name; the link type is not checked )
Definition hdf5_tools.hpp:77
string create_tree(const Cell< string > &aTree, const string &aLabel)
Definition hdf5_tools.cpp:258
void load_vector_from_file(hid_t &aLoc, const std::string &aLabel, Vector< T > &aVector, herr_t &aStatus)
Definition hdf5_tools.hpp:626
void check_read_datatype(const hid_t aFileType, const std::string &aLabel)
Guard for every read: the FILE datatype must be convertible into T.
Definition hdf5_tools.hpp:248
void load_string_from_file(hid_t &aFileID, const std::string &aLabel, std::string &aValue, herr_t &aStatus)
Definition hdf5_tools.cpp:178
void load_array_from_file(hid_t &aLoc, const std::string &aLabel, T *aData, const hsize_t aLength, herr_t &aStatus)
Definition hdf5_tools.hpp:484
void save_bool_to_file(hid_t &aFileID, const std::string &aLabel, const bool &aValue, herr_t &aStatus)
Definition hdf5_tools.cpp:21
herr_t close_group(hid_t aGroup)
Definition hdf5_tools.hpp:1060
USER GUIDES:
Definition cl_Capacitor.cpp:16
std::string datatype_string()
Definition stringtools.hpp:64
int hsize_t
Definition hdf5_types.hpp:22
int hid_t
Definition hdf5_types.hpp:20
int herr_t
Definition hdf5_types.hpp:21