Field3D
SparseFieldIO.cpp
Go to the documentation of this file.
1//----------------------------------------------------------------------------//
2
3/*
4 * Copyright (c) 2009 Sony Pictures Imageworks Inc
5 *
6 * All rights reserved.
7 *
8 * Redistribution and use in source and binary forms, with or without
9 * modification, are permitted provided that the following conditions
10 * are met:
11 *
12 * Redistributions of source code must retain the above copyright
13 * notice, this list of conditions and the following disclaimer.
14 * Redistributions in binary form must reproduce the above copyright
15 * notice, this list of conditions and the following disclaimer in the
16 * documentation and/or other materials provided with the
17 * distribution. Neither the name of Sony Pictures Imageworks nor the
18 * names of its contributors may be used to endorse or promote
19 * products derived from this software without specific prior written
20 * permission.
21 *
22 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
23 * "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
24 * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS
25 * FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE
26 * COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
27 * INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
28 * (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
29 * SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
30 * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT,
31 * STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
32 * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED
33 * OF THE POSSIBILITY OF SUCH DAMAGE.
34 */
35
36//----------------------------------------------------------------------------//
37
42//----------------------------------------------------------------------------//
43
44#include <boost/intrusive_ptr.hpp>
45
46#include <boost/thread/thread.hpp>
47#include <boost/thread/mutex.hpp>
48
49#include "InitIO.h"
50#include "SparseFieldIO.h"
51#include "Types.h"
52
53//----------------------------------------------------------------------------//
54
55using namespace boost;
56using namespace std;
57
58//----------------------------------------------------------------------------//
59
61
62//----------------------------------------------------------------------------//
63// Field3D namespaces
64//----------------------------------------------------------------------------//
65
66using namespace Exc;
67using namespace Hdf5Util;
68
69//----------------------------------------------------------------------------//
70// Anonymous namespace
71//----------------------------------------------------------------------------//
72
73namespace {
74
75//----------------------------------------------------------------------------//
76
77template <typename Data_T>
78struct ReadThreadingState
79{
80 ReadThreadingState(const OgIGroup &i_location,
82 const size_t i_numVoxels,
83 const size_t i_numBlocks,
84 const size_t i_numOccupiedBlocks,
85 const bool i_isCompressed,
86 const std::vector<size_t> &i_blockIdxToDatasetIdx)
87 : location(i_location),
88 blocks(i_blocks),
89 numVoxels(i_numVoxels),
90 numBlocks(i_numBlocks),
91 numOccupiedBlocks(i_numOccupiedBlocks),
92 isCompressed(i_isCompressed),
93 blockIdxToDatasetIdx(i_blockIdxToDatasetIdx),
94 nextBlockToRead(0)
95 { }
96 // Data members
97 const OgIGroup &location;
99 const size_t numVoxels;
100 const size_t numBlocks;
101 const size_t numOccupiedBlocks;
102 const bool isCompressed;
103 const std::vector<size_t> &blockIdxToDatasetIdx;
104 size_t nextBlockToRead;
105 // Mutexes
106 boost::mutex readMutex;
107};
108
109//----------------------------------------------------------------------------//
110
111template <typename Data_T>
112class ReadBlockOp
113{
114public:
115 ReadBlockOp(ReadThreadingState<Data_T> &state, const size_t threadId)
116 : m_state(state)
117 {
118 // Set up the compression cache
119 const uLong srcLen = m_state.numVoxels * sizeof(Data_T);
121 m_cache.resize(cmpLenBound);
122 // Initialize the reader
123 m_readerPtr.reset(
124 new OgSparseDataReader<Data_T>(m_state.location, m_state.numVoxels,
125 m_state.numOccupiedBlocks,
126 m_state.isCompressed));
127 m_reader = m_readerPtr.get();
128 // Set the thread id
129 m_reader->setThreadId(threadId);
130 }
131 void operator() ()
132 {
133 // Get next block to read
134 size_t blockIdx;
135 {
136 boost::mutex::scoped_lock lock(m_state.readMutex);
137 blockIdx = m_state.nextBlockToRead;
138 m_state.nextBlockToRead++;
139 }
140 // Loop over blocks until we run out
141 while (blockIdx < m_state.numBlocks) {
142 if (m_state.blocks[blockIdx].isAllocated) {
143 const size_t datasetIdx = m_state.blockIdxToDatasetIdx[blockIdx];
144 m_reader->readBlock(datasetIdx, m_state.blocks[blockIdx].data);
145 }
146 // Get next block idx
147 {
148 boost::mutex::scoped_lock lock(m_state.readMutex);
149 blockIdx = m_state.nextBlockToRead;
150 m_state.nextBlockToRead++;
151 }
152 }
153 }
154private:
155 // Data members ---
156 ReadThreadingState<Data_T> &m_state;
157 std::vector<uint8_t> m_cache;
158 boost::shared_ptr<OgSparseDataReader<Data_T> > m_readerPtr;
160};
161
162//----------------------------------------------------------------------------//
163
164template <typename Data_T>
165struct ThreadingState
166{
167 ThreadingState(OgOCDataset<Data_T> &i_data,
169 const size_t i_numVoxels,
170 const size_t i_numBlocks,
171 const std::vector<uint8_t> &i_isAllocated)
172 : data(i_data),
173 blocks(i_blocks),
174 numVoxels(i_numVoxels),
175 numBlocks(i_numBlocks),
176 isAllocated(i_isAllocated),
177 nextBlockToCompress(0),
178 nextBlockToWrite(0)
179 {
180 // Find first in-use block
181 for (size_t i = 0; i < numBlocks; ++i) {
182 if (blocks[i].isAllocated) {
183 nextBlockToCompress = i;
184 nextBlockToWrite = i;
185 return;
186 }
187 }
188 // If we get here, there are no active blocks. Set to numBlocks
189 nextBlockToCompress = numBlocks;
190 nextBlockToWrite = numBlocks;
191 }
192 // Data members
195 const size_t numVoxels;
196 const size_t numBlocks;
197 const std::vector<uint8_t> isAllocated;
198 size_t nextBlockToCompress;
199 size_t nextBlockToWrite;
200 // Mutexes
201 boost::mutex compressMutex;
202};
203
204//----------------------------------------------------------------------------//
205
206template <typename Data_T>
207class WriteBlockOp
208{
209public:
210 WriteBlockOp(ThreadingState<Data_T> &state, const size_t threadId)
211 : m_state(state), m_threadId(threadId)
212 {
213 const uLong srcLen = m_state.numVoxels * sizeof(Data_T);
215 m_cache.resize(cmpLenBound);
216 }
217 void operator() ()
218 {
219 const int level = 1;
220 // Get next block to compress
221 size_t blockIdx;
222 {
223 boost::mutex::scoped_lock lock(m_state.compressMutex);
224 blockIdx = m_state.nextBlockToCompress++;
225 // Step counter to next
226 while (m_state.nextBlockToCompress < m_state.numBlocks) {
227 if (m_state.blocks[m_state.nextBlockToCompress].isAllocated) {
228 break;
229 }
230 m_state.nextBlockToCompress++;
231 }
232 }
233 // Loop over blocks until we run out
234 while (blockIdx < m_state.numBlocks) {
235 if (m_state.blocks[blockIdx].isAllocated) {
236 // Block data as bytes
237 const uint8_t *srcData =
238 reinterpret_cast<const uint8_t *>(m_state.blocks[blockIdx].data);
239 // Length of compressed data is stored here
240 const uLong srcLen = m_state.numVoxels * sizeof(Data_T);
243 // Perform compression
244 const int status = compress2(&m_cache[0], &cmpLen,
246 // Error check
247 if (status != Z_OK) {
248 std::cout << "ERROR: Couldn't compress in SparseFieldIO." << std::endl
249 << " Level: " << level << std::endl
250 << " Status: " << status << std::endl
251 << " srcLen: " << srcLen << std::endl
252 << " cmpLenBound: " << cmpLenBound << std::endl
253 << " cmpLen: " << cmpLen << std::endl;
254 return;
255 }
256 // Wait to write data
257 while (m_state.nextBlockToWrite != blockIdx) {
258 // Spin
259 boost::this_thread::sleep(boost::posix_time::microseconds(1));
260 }
261 // Do the writing
262 m_state.data.addData(cmpLen, &m_cache[0]);
263 // Let next block write
264 m_state.nextBlockToWrite++;
265 while (m_state.nextBlockToWrite < m_state.numBlocks){
266 // Increment to next
267 if (m_state.blocks[m_state.nextBlockToWrite].isAllocated) {
268 break;
269 }
270 m_state.nextBlockToWrite++;
271 }
272 }
273 // Get next block idx
274 {
275 boost::mutex::scoped_lock lock(m_state.compressMutex);
276 blockIdx = m_state.nextBlockToCompress++;
277 // Step counter to next
278 while (m_state.nextBlockToCompress < m_state.numBlocks) {
279 if (m_state.blocks[m_state.nextBlockToCompress].isAllocated) {
280 break;
281 }
282 m_state.nextBlockToCompress++;
283 }
284 }
285 }
286 }
287private:
288 // Data members ---
289 ThreadingState<Data_T> &m_state;
290 std::vector<uint8_t> m_cache;
291 const size_t m_threadId;
292};
293
294//----------------------------------------------------------------------------//
295
296} // Anonymous namespace
297
298//----------------------------------------------------------------------------//
299// Static members
300//----------------------------------------------------------------------------//
301
302const int SparseFieldIO::k_versionNumber(1);
303const std::string SparseFieldIO::k_versionAttrName("version");
304const std::string SparseFieldIO::k_extentsStr("extents");
305const std::string SparseFieldIO::k_extentsMinStr("extents_min");
306const std::string SparseFieldIO::k_extentsMaxStr("extents_max");
307const std::string SparseFieldIO::k_dataWindowStr("data_window");
308const std::string SparseFieldIO::k_dataWindowMinStr("data_window_min");
309const std::string SparseFieldIO::k_dataWindowMaxStr("data_window_max");
310const std::string SparseFieldIO::k_componentsStr("components");
311const std::string SparseFieldIO::k_dataStr("data");
312const std::string SparseFieldIO::k_blockOrderStr("block_order");
313const std::string SparseFieldIO::k_numBlocksStr("num_blocks");
314const std::string SparseFieldIO::k_blockResStr("block_res");
315const std::string SparseFieldIO::k_bitsPerComponentStr("bits_per_component");
316const std::string SparseFieldIO::k_numOccupiedBlocksStr("num_occupied_blocks");
317const std::string SparseFieldIO::k_isCompressed("data_is_compressed");
318
319//----------------------------------------------------------------------------//
320
322SparseFieldIO::read(hid_t layerGroup, const std::string &filename,
323 const std::string &layerPath,
324 DataTypeEnum typeEnum)
325{
326 Box3i extents, dataW;
327 int components;
328 int blockOrder;
329 int numBlocks;
330 V3i blockRes;
331
332 if (layerGroup == -1) {
333 Msg::print(Msg::SevWarning, "Bad layerGroup.");
334 return FieldBase::Ptr();
335 }
336
337 int version;
338 if (!readAttribute(layerGroup, k_versionAttrName, 1, version))
339 throw MissingAttributeException("Couldn't find attribute: " +
341
342 if (version != k_versionNumber)
343 throw UnsupportedVersionException("SparseField version not supported: " +
345
346 if (!readAttribute(layerGroup, k_extentsStr, 6, extents.min.x))
347 throw MissingAttributeException("Couldn't find attribute: " +
349
351 throw MissingAttributeException("Couldn't find attribute: " +
353
354 if (!readAttribute(layerGroup, k_componentsStr, 1, components))
355 throw MissingAttributeException("Couldn't find attribute: " +
357
358 // Read block order
359 if (!readAttribute(layerGroup, k_blockOrderStr, 1, blockOrder))
360 throw MissingAttributeException("Couldn't find attribute: " +
362
363 // Read number of blocks total
364 if (!readAttribute(layerGroup, k_numBlocksStr, 1, numBlocks))
365 throw MissingAttributeException("Couldn't find attribute: " +
367
368 // Read block resolution in each dimension
369 if (!readAttribute(layerGroup, k_blockResStr, 3, blockRes.x))
370 throw MissingAttributeException("Couldn't find attribute: " +
372
373 // ... Check that it matches the # reported by summing the active blocks
374
375 int numCalculatedBlocks = blockRes.x * blockRes.y * blockRes.z;
376 if (numCalculatedBlocks != numBlocks)
377 throw FileIntegrityException("Incorrect block count in SparseFieldIO::read");
378
379 // Call the appropriate read function based on the data type ---
380
381 FieldBase::Ptr result;
382
383 int occupiedBlocks;
384 if (!readAttribute(layerGroup, k_numOccupiedBlocksStr, 1, occupiedBlocks))
385 throw MissingAttributeException("Couldn't find attribute: " +
387
388 // Check the data type ---
389
390 int bits;
392 throw MissingAttributeException("Couldn't find attribute: " +
394
395 bool isHalf = false;
396 bool isFloat = false;
397 bool isDouble = false;
398
399 switch (bits) {
400 case 16:
401 isHalf = true;
402 break;
403 case 64:
404 isDouble = true;
405 break;
406 case 32:
407 default:
408 isFloat = true;
409 }
410
411 // Finally, read the data ---
412
413 if (components == 1) {
414 if (isHalf && typeEnum == DataTypeHalf) {
416 field->setSize(extents, dataW);
417 field->setBlockOrder(blockOrder);
418 readData<half>(layerGroup, numBlocks, filename, layerPath, field);
419 result = field;
420 } else if (isFloat && typeEnum == DataTypeFloat) {
422 field->setSize(extents, dataW);
423 field->setBlockOrder(blockOrder);
424 readData<float>(layerGroup, numBlocks, filename, layerPath, field);
425 result = field;
426 } else if (isDouble && typeEnum == DataTypeDouble) {
428 field->setSize(extents, dataW);
429 field->setBlockOrder(blockOrder);
430 readData<double>(layerGroup, numBlocks, filename, layerPath, field);
431 result = field;
432 }
433 } else if (components == 3) {
434 if (isHalf && typeEnum == DataTypeVecHalf) {
436 field->setSize(extents, dataW);
437 field->setBlockOrder(blockOrder);
438 readData<V3h>(layerGroup, numBlocks, filename, layerPath, field);
439 result = field;
440 } else if (isFloat && typeEnum == DataTypeVecFloat) {
442 field->setSize(extents, dataW);
443 field->setBlockOrder(blockOrder);
444 readData<V3f>(layerGroup, numBlocks, filename, layerPath, field);
445 result = field;
446 } else if (isDouble && typeEnum == DataTypeVecDouble) {
448 field->setSize(extents, dataW);
449 field->setBlockOrder(blockOrder);
450 readData<V3d>(layerGroup, numBlocks, filename, layerPath, field);
451 result = field;
452 }
453 }
454
455 return result;
456}
457
458//----------------------------------------------------------------------------//
459
461SparseFieldIO::read(const OgIGroup &layerGroup, const std::string &filename,
462 const std::string &layerPath, OgDataType typeEnum)
463{
464 Box3i extents, dataW;
465 int blockOrder;
466 int numBlocks;
467 V3i blockRes;
468
469 if (!layerGroup.isValid()) {
470 throw MissingGroupException("Invalid group in SparseFieldIO::read()");
471 }
472
473 // Check version ---
474
476 layerGroup.findAttribute<int>(k_versionAttrName);
477 if (!versionAttr.isValid()) {
478 throw MissingAttributeException("Couldn't find attribute: " +
480 }
481 const int version = versionAttr.value();
482
483 if (version != k_versionNumber) {
484 throw UnsupportedVersionException("SparseField version not supported: " +
486 }
487
488 // Get extents ---
489
491 layerGroup.findAttribute<veci32_t>(k_extentsMinStr);
493 layerGroup.findAttribute<veci32_t>(k_extentsMaxStr);
494 if (!extMinAttr.isValid()) {
495 throw MissingAttributeException("Couldn't find attribute " +
497 }
498 if (!extMaxAttr.isValid()) {
499 throw MissingAttributeException("Couldn't find attribute " +
501 }
502
503 extents.min = extMinAttr.value();
504 extents.max = extMaxAttr.value();
505
506 // Get data window ---
507
509 layerGroup.findAttribute<veci32_t>(k_dataWindowMinStr);
511 layerGroup.findAttribute<veci32_t>(k_dataWindowMaxStr);
512 if (!dwMinAttr.isValid()) {
513 throw MissingAttributeException("Couldn't find attribute " +
515 }
516 if (!dwMaxAttr.isValid()) {
517 throw MissingAttributeException("Couldn't find attribute " +
519 }
520
521 dataW.min = dwMinAttr.value();
522 dataW.max = dwMaxAttr.value();
523
524 // Get num components ---
525
527 layerGroup.findAttribute<uint8_t>(k_componentsStr);
528 if (!numComponentsAttr.isValid()) {
529 throw MissingAttributeException("Couldn't find attribute " +
531 }
532
533 // Read block order ---
534
536 layerGroup.findAttribute<uint8_t>(k_blockOrderStr);
537 if (!blockOrderAttr.isValid()) {
538 throw MissingAttributeException("Couldn't find attribute: " +
540 }
541 blockOrder = blockOrderAttr.value();
542
543 // Read number of blocks total ---
544
546 layerGroup.findAttribute<uint32_t>(k_numBlocksStr);
547 if (!numBlocksAttr.isValid()) {
548 throw MissingAttributeException("Couldn't find attribute: " +
550 }
551 numBlocks = numBlocksAttr.value();
552
553 // Read block resolution in each dimension ---
554
556 layerGroup.findAttribute<veci32_t>(k_blockResStr);
557 if (!blockResAttr.isValid()) {
558 throw MissingAttributeException("Couldn't find attribute: " +
560 }
561 blockRes = blockResAttr.value();
562
563 // ... Check that it matches the # reported by summing the active blocks
564
565 int numCalculatedBlocks = blockRes.x * blockRes.y * blockRes.z;
566 if (numCalculatedBlocks != numBlocks) {
567 throw FileIntegrityException("Incorrect block count in "
568 "SparseFieldIO::read()");
569 }
570
571 // Call the appropriate read function based on the data type ---
572
575 if (!occupiedBlocksAttr.isValid()) {
576 throw MissingAttributeException("Couldn't find attribute: " +
578 }
579
580 // Check if the data is compressed ---
581
583 layerGroup.findAttribute<uint8_t>(k_isCompressed);
584 if (!isCompressedAttr.isValid()) {
585 throw MissingAttributeException("Couldn't find attribute: " +
587 }
588
589 // Finally, read the data ---
590
591 FieldBase::Ptr result;
592
594
595 if (isCompressedAttr.value() == 0) {
596 typeOnDisk = layerGroup.datasetType(k_dataStr);
597 } else {
598 typeOnDisk = layerGroup.compressedDatasetType(k_dataStr);
599 }
600
601 if (typeEnum == typeOnDisk) {
602 if (typeEnum == F3DFloat16) {
603 result = readData<float16_t>(layerGroup, extents, dataW, blockOrder,
604 numBlocks, filename, layerPath);
605 } else if (typeEnum == F3DFloat32) {
606 result = readData<float32_t>(layerGroup, extents, dataW, blockOrder,
607 numBlocks, filename, layerPath);
608 } else if (typeEnum == F3DFloat64) {
609 result = readData<float64_t>(layerGroup, extents, dataW, blockOrder,
610 numBlocks, filename, layerPath);
611 } else if (typeEnum == F3DVec16) {
612 result = readData<vec16_t>(layerGroup, extents, dataW, blockOrder,
613 numBlocks, filename, layerPath);
614 } else if (typeEnum == F3DVec32) {
615 result = readData<vec32_t>(layerGroup, extents, dataW, blockOrder,
616 numBlocks, filename, layerPath);
617 } else if (typeEnum == F3DVec64) {
618 result = readData<vec64_t>(layerGroup, extents, dataW, blockOrder,
619 numBlocks, filename, layerPath);
620 }
621 }
622
623 return result;
624}
625
626//----------------------------------------------------------------------------//
627
628bool
629SparseFieldIO::write(hid_t layerGroup, FieldBase::Ptr field)
630{
631 if (layerGroup == -1) {
632 Msg::print(Msg::SevWarning, "Bad layerGroup.");
633 return false;
634 }
635
636 // Add version attribute
638 1, k_versionNumber)) {
639 Msg::print(Msg::SevWarning, "Error adding version attribute.");
640 return false;
641 }
642
655
656 bool success = true;
657 if (halfField) {
659 } else if (floatField) {
661 } else if (doubleField) {
663 } else if (vecHalfField) {
665 } else if (vecFloatField) {
667 } else if (vecDoubleField) {
669 } else {
670 throw WriteLayerException("SparseFieldIO::write does not support the given "
671 "SparseField template parameter");
672 }
673
674 return success;
675}
676
677//----------------------------------------------------------------------------//
678
679bool
680SparseFieldIO::write(OgOGroup &layerGroup, FieldBase::Ptr field)
681{
682 using namespace Exc;
683
684 // Add version attribute
686
699
700 bool success = true;
701
702 if (floatField) {
704 }
705 else if (halfField) {
707 }
708 else if (doubleField) {
710 }
711 else if (vecFloatField) {
713 }
714 else if (vecHalfField) {
716 }
717 else if (vecDoubleField) {
719 }
720 else {
721 throw WriteLayerException("SparseFieldIO does not support the given "
722 "SparseField template parameter");
723 }
724
725 return success;
726}
727
728//----------------------------------------------------------------------------//
729
730template <class Data_T>
732SparseFieldIO::readData(const OgIGroup &location, const Box3i &extents,
733 const Box3i &dataW, const size_t blockOrder,
734 const size_t numBlocks, const std::string &filename,
735 const std::string &layerPath)
736{
737 using namespace std;
738 using namespace Exc;
739 using namespace Sparse;
740
742 result->setSize(extents, dataW);
743 result->setBlockOrder(blockOrder);
744
746 const int components = FieldTraits<Data_T>::dataDims();
747 const size_t numVoxels = (1 << (result->m_blockOrder * 3));
748 const int valuesPerBlock = (1 << (result->m_blockOrder * 3)) * components;
749
750 // Read the number of occupied blocks ---
751
753 location.findAttribute<uint32_t>(k_numOccupiedBlocksStr);
754 if (!occupiedBlocksAttr.isValid()) {
755 throw MissingAttributeException("Couldn't find attribute: " +
757 }
758 const size_t occupiedBlocks = occupiedBlocksAttr.value();
759
760 // Set up the dynamic read info ---
761
762 if (dynamicLoading) {
763 // Set up the field reference
765 result->addReference(filename, layerPath, valuesPerBlock, numVoxels,
766 occupiedBlocks);
767 }
768
769 // Read the block info data sets ---
770
771 SparseBlock<Data_T> *blocks = result->m_blocks;
772
773 // ... Read the isAllocated array and set up the block mapping array
774 std::vector<size_t> blockIdxToDatasetIdx(numBlocks);
775
776 {
777 // Grab the data
778 vector<uint8_t> isAllocated(numBlocks);
780 location.findDataset<uint8_t>("block_is_allocated_data");
781 if (!isAllocatedData.isValid()) {
782 throw MissingGroupException("Couldn't find block_is_allocated_data: ");
783 }
784 isAllocatedData.getData(0, &isAllocated[0], OGAWA_THREAD);
785 // Allocate the blocks and set up the block mapping array
786 for (size_t i = 0, nextBlockOnDisk = 0; i < numBlocks; ++i) {
787 blocks[i].isAllocated = isAllocated[i];
788 if (!dynamicLoading && isAllocated[i]) {
789 blocks[i].resize(numVoxels);
790 // Update the block mapping array
791 blockIdxToDatasetIdx[i] = nextBlockOnDisk;
793 }
794 }
795 }
796
797 // ... Read the emptyValue array ---
798
799 {
800 // Grab the data
801 vector<Data_T> emptyValue(numBlocks);
803 location.findDataset<Data_T>("block_empty_value_data");
804 if (!emptyValueData.isValid()) {
805 throw MissingGroupException("Couldn't find block_empty_value_data: ");
806 }
807 emptyValueData.getData(0, &emptyValue[0], OGAWA_THREAD);
808 // Fill in the field
809 for (size_t i = 0; i < numBlocks; ++i) {
810 blocks[i].emptyValue = emptyValue[i];
811 }
812 }
813
814 // Read the data ---
815
816 // Check whether data is compressed
818 location.findAttribute<uint8_t>(k_isCompressed);
819 const bool isCompressed = isCompressedAttr.value() != 0;
820
821 if (occupiedBlocks > 0) {
822 if (dynamicLoading) {
823 // Defer loading to the sparse cache
824 result->setupReferenceBlocks();
825 } else {
826 // Threading state
827 ReadThreadingState<Data_T> state(location, blocks, numVoxels, numBlocks,
828 occupiedBlocks, isCompressed,
829 blockIdxToDatasetIdx);
830 // Number of threads
831 const size_t numThreads = numIOThreads();
832 // Launch threads
833 boost::thread_group threads;
834 for (size_t i = 0; i < numThreads; ++i) {
835 threads.create_thread(ReadBlockOp<Data_T>(state, i));
836 }
837 threads.join_all();
838 }
839 }
840
841 return result;
842}
843
844//----------------------------------------------------------------------------//
845// Template implementations
846//----------------------------------------------------------------------------//
847
849template <class Data_T>
850bool SparseFieldIO::writeInternal(hid_t layerGroup,
851 typename SparseField<Data_T>::Ptr field)
852{
853 using namespace std;
854 using namespace Exc;
855 using namespace Hdf5Util;
856 using namespace Sparse;
857
858 Box3i ext(field->extents()), dw(field->dataWindow());
859
860 int components = FieldTraits<Data_T>::dataDims();
861
862 int valuesPerBlock = (1 << (field->m_blockOrder * 3)) * components;
863
864 // Add extents attribute ---
865
866 int extents[6] =
867 { ext.min.x, ext.min.y, ext.min.z, ext.max.x, ext.max.y, ext.max.z };
868
869 if (!writeAttribute(layerGroup, k_extentsStr, 6, extents[0])) {
870 Msg::print(Msg::SevWarning, "Error adding size attribute.");
871 return false;
872 }
873
874 // Add data window attribute ---
875
876 int dataWindow[6] =
877 { dw.min.x, dw.min.y, dw.min.z, dw.max.x, dw.max.y, dw.max.z };
878
879 if (!writeAttribute(layerGroup, k_dataWindowStr, 6, dataWindow[0])) {
880 Msg::print(Msg::SevWarning, "Error adding size attribute.");
881 return false;
882 }
883
884 // Add components attribute ---
885
886 if (!writeAttribute(layerGroup, k_componentsStr, 1, components)) {
887 Msg::print(Msg::SevWarning, "Error adding components attribute.");
888 return false;
889 }
890
891 // Add block order attribute ---
892
893 int blockOrder = field->m_blockOrder;
894
895 if (!writeAttribute(layerGroup, k_blockOrderStr, 1, blockOrder)) {
896 Msg::print(Msg::SevWarning, "Error adding block order attribute.");
897 return false;
898 }
899
900 // Add number of blocks attribute ---
901
902 V3i &blockRes = field->m_blockRes;
903 int numBlocks = blockRes.x * blockRes.y * blockRes.z;
904
905 if (!writeAttribute(layerGroup, k_numBlocksStr, 1, numBlocks)) {
906 Msg::print(Msg::SevWarning, "Error adding number of blocks attribute.");
907 return false;
908 }
909
910 // Add block resolution in each dimension ---
911
912 if (!writeAttribute(layerGroup, k_blockResStr, 3, blockRes.x)) {
913 Msg::print(Msg::SevWarning, "Error adding block res attribute.");
914 return false;
915 }
916
917 // Add the bits per component attribute ---
918
921 Msg::print(Msg::SevWarning, "Error adding bits per component attribute.");
922 return false;
923 }
924
925 // Write the block info data sets ---
926
927 SparseBlock<Data_T> *blocks = field->m_blocks;
928
929 // ... Write the isAllocated array
930 {
931 vector<char> isAllocated(numBlocks);
932 for (int i = 0; i < numBlocks; ++i) {
933 isAllocated[i] = static_cast<char>(blocks[i].isAllocated);
934 }
935 writeSimpleData<char>(layerGroup, "block_is_allocated_data", isAllocated);
936 }
937
938 // ... Write the emptyValue array
939 {
940 vector<Data_T> emptyValue(numBlocks);
941 for (int i = 0; i < numBlocks; ++i) {
942 emptyValue[i] = static_cast<Data_T>(blocks[i].emptyValue);
943 }
944 writeSimpleData<Data_T>(layerGroup, "block_empty_value_data", emptyValue);
945 }
946
947 // Count the number of occupied blocks ---
948 int occupiedBlocks = 0;
949 for (int i = 0; i < numBlocks; ++i) {
950 if (blocks[i].isAllocated) {
951 occupiedBlocks++;
952 }
953 }
954
955 if (!writeAttribute(layerGroup, k_numOccupiedBlocksStr, 1, occupiedBlocks)) {
956 throw WriteAttributeException("Couldn't add attribute " +
958 }
959
960 if (occupiedBlocks > 0) {
961
962 // Make the memory data space
963 hsize_t memDims[1];
964 memDims[0] = valuesPerBlock;
967
968 // Make the file data space
969 hsize_t fileDims[2];
970 fileDims[0] = occupiedBlocks;
971 fileDims[1] = valuesPerBlock;
974
975 // Set up gzip property list
979 chunkSize[0] = 1;
980 chunkSize[1] = valuesPerBlock;
981 if (gzipAvailable) {
983 if (status < 0) {
984 return false;
985 }
987 if (status < 0) {
988 return false;
989 }
990 }
991
992 // Add the data set
995 fileDataSpace.id(),
997 if (dataSet.id() < 0)
998 throw CreateDataSetException("Couldn't create data set in "
999 "SparseFieldIO::writeInternal");
1000
1001 // For each allocated block ---
1002
1003 int nextBlockIdx = 0;
1004 hsize_t offset[2];
1005 hsize_t count[2];
1006 herr_t status;
1007
1008 for (int i = 0; i < numBlocks; ++i) {
1009 if (blocks[i].isAllocated) {
1010 offset[0] = nextBlockIdx; // Index of next block
1011 offset[1] = 0; // Index of first data in block. Always 0
1012 count[0] = 1; // Number of columns to read. Always 1
1013 count[1] = valuesPerBlock; // Number of values in one column
1015 offset, NULL, count, NULL);
1016 if (status < 0) {
1018 "Couldn't select slab " +
1019 boost::lexical_cast<std::string>(nextBlockIdx));
1020 }
1021 Data_T *data = field->m_blocks[i].data;
1023 memDataSpace.id(),
1024 fileDataSpace.id(), H5P_DEFAULT, data);
1025 if (status < 0) {
1027 "Couldn't write slab " +
1028 boost::lexical_cast<std::string>(nextBlockIdx));
1029 }
1030 // Increment nextBlockIdx
1031 nextBlockIdx++;
1032 }
1033 }
1034
1035 } // if occupiedBlocks > 0
1036
1037 return true;
1038
1039}
1040
1041//----------------------------------------------------------------------------//
1042
1043template <class Data_T>
1044bool SparseFieldIO::writeInternal(OgOGroup &layerGroup,
1045 typename SparseField<Data_T>::Ptr field)
1046{
1047 using namespace Exc;
1048 using namespace Sparse;
1049
1050 SparseBlock<Data_T> *blocks = field->m_blocks;
1051
1052 const int components = FieldTraits<Data_T>::dataDims();
1054 const V3i &blockRes = field->m_blockRes;
1055 const size_t numBlocks = blockRes.x * blockRes.y * blockRes.z;
1056 const size_t numVoxels = (1 << (field->m_blockOrder * 3));
1057
1058 const Box3i ext(field->extents()), dw(field->dataWindow());
1059
1060 // Add attributes ---
1061
1064
1067
1069
1071
1073 field->m_blockOrder);
1074
1076
1078
1080
1081 // Write the isAllocated array
1082 std::vector<uint8_t> isAllocated(numBlocks);
1083 for (size_t i = 0; i < numBlocks; ++i) {
1084 isAllocated[i] = static_cast<uint8_t>(blocks[i].isAllocated);
1085 }
1086 OgODataset<uint8_t> isAllocatedData(layerGroup, "block_is_allocated_data");
1087 isAllocatedData.addData(numBlocks, &isAllocated[0]);
1088
1089 // Write the emptyValue array
1090 std::vector<Data_T> emptyValue(numBlocks);
1091 for (size_t i = 0; i < numBlocks; ++i) {
1092 emptyValue[i] = static_cast<Data_T>(blocks[i].emptyValue);
1093 }
1094 OgODataset<Data_T> emptyValueData(layerGroup, "block_empty_value_data");
1095 emptyValueData.addData(numBlocks, &emptyValue[0]);
1096
1097 // Count the number of occupied blocks
1098 int occupiedBlocks = 0;
1099 for (size_t i = 0; i < numBlocks; ++i) {
1100 if (blocks[i].isAllocated) {
1101 occupiedBlocks++;
1102 }
1103 }
1106 occupiedBlocks);
1107
1108 // Add data to file ---
1109
1110 // Create the compressed dataset regardless of whether there are blocks
1111 // to write.
1112 OgOCDataset<Data_T> data(layerGroup, k_dataStr);
1113 // Write data if there is any
1114 if (occupiedBlocks > 0) {
1115 // Threading state
1116 ThreadingState<Data_T> state(data, blocks, numVoxels, numBlocks,
1117 isAllocated);
1118 // Number of threads
1119 const size_t numThreads = numIOThreads();
1120 // Launch threads
1121 boost::thread_group threads;
1122 for (size_t i = 0; i < numThreads; ++i) {
1123 threads.create_thread(WriteBlockOp<Data_T>(state, i));
1124 }
1125 threads.join_all();
1126 }
1127
1128 return true;
1129}
1130
1131//----------------------------------------------------------------------------//
1132
1133template <class Data_T>
1134bool SparseFieldIO::readData(hid_t location,
1135 int numBlocks,
1136 const std::string &filename,
1137 const std::string &layerPath,
1138 typename SparseField<Data_T>::Ptr result)
1139{
1140 using namespace std;
1141 using namespace Exc;
1142 using namespace Hdf5Util;
1143 using namespace Sparse;
1144
1145 int occupiedBlocks;
1146
1148
1149 int components = FieldTraits<Data_T>::dataDims();
1150 int numVoxels = (1 << (result->m_blockOrder * 3));
1151 int valuesPerBlock = numVoxels * components;
1152
1153 // Read the number of occupied blocks ---
1154
1155 if (!readAttribute(location, k_numOccupiedBlocksStr, 1, occupiedBlocks))
1156 throw MissingAttributeException("Couldn't find attribute: " +
1158
1159 // Set up the dynamic read info ---
1160
1161 if (dynamicLoading) {
1162 // Set up the field reference
1163 result->addReference(filename, layerPath,
1164 valuesPerBlock, numVoxels,
1165 occupiedBlocks);
1166 }
1167
1168 // Read the block info data sets ---
1169
1170 SparseBlock<Data_T> *blocks = result->m_blocks;
1171
1172 // ... Read the isAllocated array
1173
1174 {
1175 vector<char> isAllocated(numBlocks);
1176 readSimpleData<char>(location, "block_is_allocated_data", isAllocated);
1177 for (int i = 0; i < numBlocks; ++i) {
1178 blocks[i].isAllocated = isAllocated[i];
1179 if (!dynamicLoading && isAllocated[i]) {
1180 blocks[i].resize(numVoxels);
1181 }
1182 }
1183 }
1184
1185 // ... Read the emptyValue array ---
1186
1187 {
1188 vector<Data_T> emptyValue(numBlocks);
1189 readSimpleData<Data_T>(location, "block_empty_value_data", emptyValue);
1190 for (int i = 0; i < numBlocks; ++i) {
1191 blocks[i].emptyValue = emptyValue[i];
1192 }
1193 }
1194
1195 // Read the data ---
1196
1197 if (occupiedBlocks > 0) {
1198
1199 if (dynamicLoading) {
1200
1201 result->setupReferenceBlocks();
1202
1203 } else {
1204
1205 size_t b = 0, bend = b + numBlocks;
1206
1207 SparseDataReader<Data_T> reader(location, valuesPerBlock, occupiedBlocks);
1208
1209 // We'll read at most 50meg at a time
1210 static const long maxMemPerPass = 50*1024*1024;
1211
1212 for (int nextBlockIdx = 0;;) {
1213
1214 long mem = 0;
1215 std::vector<Data_T*> memoryList;
1216
1217 for (; b != bend && mem < maxMemPerPass; ++b) {
1218 if (blocks[b].isAllocated) {
1219 mem += sizeof(Data_T)*numVoxels;
1220 memoryList.push_back(blocks[b].data);
1221 }
1222 }
1223
1224 // all done.
1225 if (!memoryList.size()) {
1226 break;
1227 }
1228
1229 reader.readBlockList(nextBlockIdx, memoryList);
1230 nextBlockIdx += memoryList.size();
1231 }
1232
1233 }
1234
1235 } // if occupiedBlocks > 0
1236
1237 return true;
1238
1239}
1240
1241//----------------------------------------------------------------------------//
1242
1244
1245//----------------------------------------------------------------------------//
Contains the initIO function.
FIELD3D_API size_t numIOThreads()
Returns the number of I/O threads to use.
Definition InitIO.cpp:92
Imath::V3i V3i
Definition SpiMathLib.h:71
Imath::Box3i Box3i
Definition SpiMathLib.h:77
#define FIELD3D_MTX_T
Definition StdMathLib.h:99
Field3D::V3i veci32_t
Definition Traits.h:94
OgDataType
Enumerates the various uses for Ogawa-level groups.
Definition Traits.h:125
@ F3DFloat32
Definition Traits.h:142
@ F3DFloat64
Definition Traits.h:143
@ F3DVec64
Definition Traits.h:148
@ F3DVec32
Definition Traits.h:147
@ F3DFloat16
Definition Traits.h:141
@ F3DVec16
Definition Traits.h:146
DataTypeEnum
Definition Traits.h:108
@ DataTypeFloat
Definition Traits.h:112
@ DataTypeHalf
Definition Traits.h:109
@ DataTypeVecHalf
Definition Traits.h:114
@ DataTypeVecDouble
Definition Traits.h:116
@ DataTypeVecFloat
Definition Traits.h:115
@ DataTypeDouble
Definition Traits.h:113
Contains typedefs for the commonly used types in Field3D.
boost::intrusive_ptr< FieldBase > Ptr
Definition Field.h:97
const Box3i & extents() const
Returns the extents of the data. This signifies the relevant area that the data exists over....
Definition Field.h:249
const Box3i & dataWindow() const
Returns the data window. Any coordinate inside this window is safe to pass to value() in the Field su...
Definition Field.h:253
static int dataDims()
Definition Traits.h:178
Scoped object - creates a dataset on creation and closes it on destruction.
Definition Hdf5Util.h:264
Scoped object - creates a dataspace on creation and closes it on destruction.
Definition Hdf5Util.h:235
This class gets used by SparseFieldIO and SparseFileManager to read the block data....
This Field subclass stores voxel data in block-allocated arrays.
Block * m_blocks
Array of blocks. Not using std::vector since SparseBlock is noncopyable.
void addReference(const std::string &filename, const std::string &layerPath, int valuesPerBlock, int numVoxels, int occupiedBlocks)
Internal function to create a Reference for the current field, for use in dynamic reading.
boost::intrusive_ptr< SparseField > Ptr
int m_blockOrder
Block order (size = 2^blockOrder)
void setupReferenceBlocks()
Internal function to setup the Reference's block pointers, for use with dynamic reading.
V3i m_blockRes
Block array resolution.
bool doLimitMemUse() const
Returns whether to limit memory usage and do dynamic loading for sparse fields.
static SparseFileManager & singleton()
Returns a reference to the singleton instance.
FIELD3D_API bool readAttribute(hid_t location, const std::string &attrName, std::string &value)
Reads a string attribute.
FIELD3D_API bool writeAttribute(hid_t location, const std::string &attrName, const std::string &value)
Writes a string attribute.
FIELD3D_API bool checkHdf5Gzip()
Checks whether gzip is available in the current hdf5 library.
Definition Hdf5Util.cpp:722
Namespace for Exception objects.
Definition Exception.h:57
Contains utility functions and classes for Hdf5 files.
Definition Hdf5Util.h:86
@ SevWarning
Definition Log.h:68
FIELD3D_API void print(Severity severity, const std::string &message)
Sends the string to the assigned output, prefixing the message with the severity.
Definition Log.cpp:70
Namespace for sparse field specifics.
#define FIELD3D_NAMESPACE_SOURCE_CLOSE
Definition ns.h:60
static int h5bits()
Storage for one individual block of a SparseField.
bool isAllocated
Whether the block is allocated or not.
void resize(int n)
Alloc data.
Data_T emptyValue
The value to use if the block isn't allocated. We allow setting this per block so that we for example...
Data_T * data
Pointer to data. Null if block is unallocated.