[MOAB-dev] r1737 - in MOAB/trunk: . mhdf/include mhdf/src test/h5file
kraftche at mcs.anl.gov
kraftche at mcs.anl.gov
Wed Apr 2 19:24:41 CDT 2008
Author: kraftche
Date: 2008-04-02 19:24:41 -0500 (Wed, 02 Apr 2008)
New Revision: 1737
Modified:
MOAB/trunk/ReadHDF5.cpp
MOAB/trunk/WriteHDF5.cpp
MOAB/trunk/WriteHDF5.hpp
MOAB/trunk/mhdf/include/mhdf.h
MOAB/trunk/mhdf/src/tags.c
MOAB/trunk/test/h5file/h5varlen.cpp
Log:
o Do variable-length tag I/O in chunks
- reduces HDF5 read and write of files with many variable-length
tags by an a factor of 20 and 10 respectively.
o Modify mhdf API to pass back data table length for variable-length tags.
o Add test for variable length tags containing handles
o Add test for very large variable-length tag values (i.e. 1M doubles).
Modified: MOAB/trunk/ReadHDF5.cpp
===================================================================
--- MOAB/trunk/ReadHDF5.cpp 2008-04-02 19:11:26 UTC (rev 1736)
+++ MOAB/trunk/ReadHDF5.cpp 2008-04-03 00:24:41 UTC (rev 1737)
@@ -1495,7 +1495,7 @@
mhdf_Status status;
std::string name;
MBErrorCode rval;
- long num_values;
+ long num_values, data_size;
hid_t data[3];
MBTagType mbtype;
assert ((hdf_read_type == 0) || (H5Tget_size(hdf_read_type) == read_size));
@@ -1508,12 +1508,14 @@
if (MB_SUCCESS != rval)
return rval;
- mhdf_openSparseTagData( filePtr, name.c_str(), &num_values, data, &status );
+ mhdf_openSparseTagData( filePtr, name.c_str(), &num_values, &data_size, data, &status );
if (mhdf_isError( &status ) )
{
readUtil->report_error( mhdf_message( &status ) );
return MB_FAILURE;
}
+ // fixed-length tag
+ assert( num_values == data_size );
// Split buffer into two portions: one for handles and one for data.
@@ -1609,6 +1611,9 @@
return MB_SUCCESS;
}
+#define assert_range( ARRAY, BYTES ) \
+ assert( (char*)(ARRAY) >= dataBuffer && ((char*)(ARRAY)) + (BYTES) <= dataBuffer + bufferSize )
+
MBErrorCode ReadHDF5::read_var_len_tag( MBTag tag_handle,
hid_t hdf_read_type,
bool is_handle_type )
@@ -1616,7 +1621,7 @@
mhdf_Status status;
std::string name;
MBErrorCode rval;
- long num_values;
+ long num_values, num_data;
hid_t data[3];
MBTagType mbtype;
// hdf_read_type is NULL (zero) for opaque tag data.
@@ -1632,19 +1637,47 @@
if (MB_SUCCESS != rval)
return rval;
- mhdf_openSparseTagData( filePtr, name.c_str(), &num_values, data, &status );
+ mhdf_openSparseTagData( filePtr, name.c_str(), &num_values, &num_data, data, &status );
if (mhdf_isError( &status ) )
{
readUtil->report_error( mhdf_message( &status ) );
return MB_FAILURE;
}
+
+ // Subdivide buffer into 5 chunks:
+ // Be careful of order so alignment is valid
+ // 1) pointer array (input for MOAB)
+ // 2) offset array (read from file)
+ // 3) entity handles
+ // 4) tag sizes (calculated from file data)
+ // 5) tag data
+ const size_t avg_data_size = (num_data * elem_size) / num_values;
+ const size_t per_ent_size = sizeof(void*) + sizeof(long) + sizeof(MBEntityHandle) + sizeof(int);
+ long num_ent = bufferSize / (per_ent_size + 2*avg_data_size);
+ if (num_ent == 0)
+ num_ent = bufferSize / 4 / sizeof(void*);
+ const size_t data_buffer_size = (bufferSize - num_ent * per_ent_size) / elem_size;
+ const void** const pointer_buffer = reinterpret_cast<const void**>(dataBuffer);
+ long* const end_idx_buffer = reinterpret_cast<long*>(pointer_buffer + num_ent);
+ char* handle_buffer_start = reinterpret_cast<char*>(end_idx_buffer + num_ent);
+ if (((size_t)handle_buffer_start) % sizeof(MBEntityHandle))
+ handle_buffer_start += sizeof(MBEntityHandle) - ((size_t)handle_buffer_start);
+ MBEntityHandle* const handle_buffer = reinterpret_cast<MBEntityHandle*>(handle_buffer_start);
+ int* const size_buffer = reinterpret_cast<int*>(handle_buffer + num_ent);
+ char* const data_buffer = reinterpret_cast<char*>(size_buffer + num_ent);
- // Process each entity individually
- MBEntityHandle id;
- long offset, prev_offset = 0;
- for (long i = 0; i < num_values; ++i) {
- // read entity ID
- mhdf_readSparseTagEntities( data[0], i, 1, handleType, &id, &status );
+
+ // do num_ent blocks of entities
+ long remaining = num_values;
+ long offset = 0;
+ long data_offset = 0;
+ while (remaining) {
+ const long count = remaining < num_ent ? remaining : num_ent;
+ remaining -= count;
+
+ // read entity IDs
+ assert_range( handle_buffer, count * sizeof(MBEntityHandle) );
+ mhdf_readSparseTagEntities( data[0], offset, count, handleType, handle_buffer, &status );
if (mhdf_isError( &status ))
{
readUtil->report_error( mhdf_message( &status ) );
@@ -1654,7 +1687,7 @@
return MB_FAILURE;
}
// convert entity ID to MBEntityHandle
- rval = convert_id_to_handle( &id, 1 );
+ rval = convert_id_to_handle( handle_buffer, count );
if (MB_SUCCESS != rval)
{
mhdf_closeData( filePtr, data[0], &status );
@@ -1663,7 +1696,8 @@
return rval;
}
// read end index of tag value
- mhdf_readSparseTagIndices( data[2], i, 1, H5T_NATIVE_LONG, &offset, &status );
+ assert_range( end_idx_buffer, count * sizeof(long) );
+ mhdf_readSparseTagIndices( data[2], offset, count, H5T_NATIVE_LONG, end_idx_buffer, &status );
if (mhdf_isError( &status ))
{
readUtil->report_error( mhdf_message( &status ) );
@@ -1672,44 +1706,111 @@
mhdf_closeData( filePtr, data[2], &status );
return MB_FAILURE;
}
- // calculate length of tag value
- ++offset;
- int count = (int)(offset - prev_offset);
- assert( count * elem_size <= bufferSize );
- mhdf_readSparseTagValues( data[1], prev_offset, count, hdf_read_type, dataBuffer, &status );
- if (mhdf_isError( &status ))
- {
- readUtil->report_error( mhdf_message( &status ) );
- mhdf_closeData( filePtr, data[0], &status );
- mhdf_closeData( filePtr, data[1], &status );
- mhdf_closeData( filePtr, data[2], &status );
- return MB_FAILURE;
- }
- // for handle-type tags, convert file IDs to MBEntityHandles in tag data
- if (is_handle_type)
- {
- rval = convert_id_to_handle( (MBEntityHandle*)dataBuffer, count );
- if (MB_SUCCESS != rval)
- {
- mhdf_closeData( filePtr, data[0], &status );
- mhdf_closeData( filePtr, data[1], &status );
- mhdf_closeData( filePtr, data[2], &status );
- return rval;
+ // update entity offset
+ offset += count;
+
+ char* pointer = data_buffer;
+ long finished = 0;
+ while (finished < count) {
+
+ // iterate over entities until we have enough to fill the data buffer
+ long i = finished;
+ long prev_end_idx = data_offset - 1;
+ for (; i < count && end_idx_buffer[i] - data_offset < (long)data_buffer_size; ++i) {
+ // calculate tag size for entity
+ const int size = (end_idx_buffer[i] - prev_end_idx);
+ if (size < 0 || (long)size != (end_idx_buffer[i] - prev_end_idx)) {
+ readUtil->report_error( "Invalid end index in variable length tag data for tag: \"%s\"", name.c_str() );
+ mhdf_closeData( filePtr, data[0], &status );
+ mhdf_closeData( filePtr, data[1], &status );
+ mhdf_closeData( filePtr, data[2], &status );
+ return MB_FAILURE;
+ }
+ assert_range( size_buffer + i - finished, sizeof(int) );
+ size_buffer[i-finished] = size * elem_size;
+ assert_range( pointer_buffer + i - finished, sizeof(void*) );
+ pointer_buffer[i-finished] = pointer;
+ assert_range( pointer_buffer[i-finished], size_buffer[i - finished] );
+ pointer += size_buffer[i-finished];
+ prev_end_idx = end_idx_buffer[i];
}
+
+ const size_t num_val = prev_end_idx - data_offset + 1;
+ if (num_val) {
+ // read data
+ assert( num_val <= data_buffer_size );
+ assert_range( data_buffer, num_val * elem_size );
+ mhdf_readSparseTagValues( data[1], data_offset, num_val, hdf_read_type, data_buffer, &status );
+ if (mhdf_isError( &status )) {
+ readUtil->report_error( mhdf_message( &status ) );
+ mhdf_closeData( filePtr, data[0], &status );
+ mhdf_closeData( filePtr, data[1], &status );
+ mhdf_closeData( filePtr, data[2], &status );
+ return MB_FAILURE;
+ }
+ data_offset += num_val;
+
+ if (is_handle_type) {
+ rval = convert_id_to_handle( reinterpret_cast<MBEntityHandle*>(data_buffer), num_val );
+ if (MB_SUCCESS != rval) {
+ mhdf_closeData( filePtr, data[0], &status );
+ mhdf_closeData( filePtr, data[1], &status );
+ mhdf_closeData( filePtr, data[2], &status );
+ return rval;
+ }
+ }
+ // put data in tags
+ rval = iFace->tag_set_data( tag_handle, handle_buffer + finished, i - finished, pointer_buffer, size_buffer );
+ if (MB_SUCCESS != rval) {
+ mhdf_closeData( filePtr, data[0], &status );
+ mhdf_closeData( filePtr, data[1], &status );
+ mhdf_closeData( filePtr, data[2], &status );
+ return rval;
+ }
+
+ finished = i;
+ }
+ else if (finished < i) { // zero-length tag values???
+ finished = i;
+ }
+ // if tag value to big for buffer
+ else {
+ const int size = (end_idx_buffer[i] - prev_end_idx);
+ std::vector<char*> tmp_buffer( size * elem_size );
+ mhdf_readSparseTagValues( data[1], data_offset, size, hdf_read_type, &tmp_buffer[0], &status );
+ if (mhdf_isError( &status )) {
+ readUtil->report_error( mhdf_message( &status ) );
+ mhdf_closeData( filePtr, data[0], &status );
+ mhdf_closeData( filePtr, data[1], &status );
+ mhdf_closeData( filePtr, data[2], &status );
+ return MB_FAILURE;
+ }
+ data_offset += size;
+
+ if (is_handle_type) {
+ rval = convert_id_to_handle( reinterpret_cast<MBEntityHandle*>(&tmp_buffer[0]), size );
+ if (MB_SUCCESS != rval) {
+ mhdf_closeData( filePtr, data[0], &status );
+ mhdf_closeData( filePtr, data[1], &status );
+ mhdf_closeData( filePtr, data[2], &status );
+ return rval;
+ }
+ }
+
+ pointer_buffer[0] = &tmp_buffer[0];
+ size_buffer[0] = size * elem_size;
+ rval = iFace->tag_set_data( tag_handle, handle_buffer + i, 1, pointer_buffer, size_buffer );
+ if (MB_SUCCESS != rval) {
+ mhdf_closeData( filePtr, data[0], &status );
+ mhdf_closeData( filePtr, data[1], &status );
+ mhdf_closeData( filePtr, data[2], &status );
+ return rval;
+ }
+
+ prev_end_idx = end_idx_buffer[i];
+ ++finished;
+ }
}
- // set tag data
- const void* ptrarr[1] = { dataBuffer };
- int bytes = count * elem_size;
- rval = iFace->tag_set_data( tag_handle, &id, 1, ptrarr, &bytes );
- if (MB_SUCCESS != rval)
- {
- mhdf_closeData( filePtr, data[0], &status );
- mhdf_closeData( filePtr, data[1], &status );
- mhdf_closeData( filePtr, data[2], &status );
- return rval;
- }
- // iterate
- prev_offset = offset;
}
mhdf_closeData( filePtr, data[0], &status );
Modified: MOAB/trunk/WriteHDF5.cpp
===================================================================
--- MOAB/trunk/WriteHDF5.cpp 2008-04-02 19:11:26 UTC (rev 1736)
+++ MOAB/trunk/WriteHDF5.cpp 2008-04-03 00:24:41 UTC (rev 1737)
@@ -192,18 +192,24 @@
return (A); \
} while(false)
-bool WriteHDF5::convert_handle_tag( MBEntityHandle* data, size_t count )
+bool WriteHDF5::convert_handle_tag( const MBEntityHandle* source,
+ MBEntityHandle* dest, size_t count ) const
{
- assert( sizeof(MBEntityHandle) == sizeof(id_t) );
bool some_valid = false;
for (size_t i = 0; i < count; ++i) {
- data[i] = idMap.find( data[i] );
- if (data[i])
+ dest[i] = idMap.find( source[i] );
+ if (dest[i])
some_valid = true;
}
return some_valid;
}
+bool WriteHDF5::convert_handle_tag( MBEntityHandle* data, size_t count ) const
+{
+ assert( sizeof(MBEntityHandle) == sizeof(id_t) );
+ return convert_handle_tag( data, data, count );
+}
+
MBErrorCode WriteHDF5::assign_ids( const MBRange& entities, id_t id )
{
MBRange::const_pair_iterator pi;
@@ -1498,7 +1504,7 @@
int mb_size;
MBTagType mb_type;
MBDataType mb_data_type;
- long table_size;
+ long table_size, data_size;
hid_t value_type = 0;
//get tag properties from moab
@@ -1518,10 +1524,13 @@
mhdf_openSparseTagData( filePtr,
name.c_str(),
&table_size,
+ &data_size,
tables,
&status);
CHK_MHDF_ERR_0(status);
assert( tag_data.range.size() + tag_data.offset <= (unsigned long)table_size );
+ // fixed-length tag
+ assert( table_size == data_size );
// Write IDs for tagged entities
rval = write_sparse_ids( tag_data, tables[0] );
@@ -1619,6 +1628,7 @@
hid_t tables[3];
std::string name;
long table_size;
+ long data_table_size;
//get tag properties from moab
if (MB_SUCCESS != iFace->tag_get_name( tag_data.tag_id, name ))
@@ -1643,12 +1653,13 @@
if (mb_size != MB_VARIABLE_LENGTH)
return MB_FAILURE;
-DEBUGOUT((std::string("Tag: ") + name + "\n").c_str());
+DEBUGOUT((std::string("Var Len Tag: ") + name + "\n").c_str());
//open tables to write info
mhdf_openSparseTagData( filePtr,
name.c_str(),
&table_size,
+ &data_table_size,
tables,
&status);
CHK_MHDF_ERR_0(status);
@@ -1659,37 +1670,107 @@
CHK_MB_ERR_2( rval, tables, status );
mhdf_closeData( filePtr, tables[0], &status );
CHK_MHDF_ERR_2(status, tables + 1);
+
+
+ // Split buffer into four chunks ordered such that there are no
+ // alignment issues:
+ // 1) tag data pointer buffer (data from MOAB)
+ // 2) tag offset buffer (to be written)
+ // 3) tag size buffer (data from MOAB)
+ // 4) concatenated tag data buffer (to be written)
+ const size_t quarter_buffer_size = bufferSize / 4;
+ const size_t num_entities = quarter_buffer_size / sizeof(void*);
+ assert( num_entities > 0 );
+ const void** const pointer_buffer = reinterpret_cast<const void**>(dataBuffer);
+ long* const offset_buffer = reinterpret_cast<long*>(pointer_buffer + num_entities);
+ int* const size_buffer = reinterpret_cast<int*>(offset_buffer + num_entities);
+ char* const data_buffer = reinterpret_cast<char*>(size_buffer + num_entities);
+ assert( data_buffer < bufferSize + dataBuffer );
+ const size_t data_buffer_size = dataBuffer + bufferSize - data_buffer;
- // Write each tag value separately
- size_t indx_offset = tag_data.offset;
- size_t data_offset = tag_data.varDataOffset;
- int size;
- unsigned long idx;
- const void* ptr;
- for (MBRange::const_iterator i = tag_data.range.begin();
- i!= tag_data.range.end(); ++i, ++indx_offset) {
-
- rval = iFace->tag_get_data( tag_data.tag_id, &*i, 1, &ptr, &size );
+ // offsets into tables
+ long offset_offset = tag_data.offset; // offset at which to write indices
+ long data_offset = tag_data.varDataOffset; // offset at which to write data buffer
+ long offset = tag_data.varDataOffset; // used to calculate indices
+
+ // iterate in blocks of num_entities entities
+ size_t bytes = 0; // occupied size of data buffer
+ size_t remaining = tag_data.range.size();
+ MBRange::const_iterator i = tag_data.range.begin();
+ while (remaining) {
+ const size_t count = remaining < num_entities ? remaining : num_entities;
+ remaining -= count;
+
+ // get subset of entity handles
+ MBRange::const_iterator e = i; e += count;
+ MBRange subrange; subrange.merge( i, e );
+ i = e;
+
+ // get pointers and sizes for entities from MOAB
+ rval = iFace->tag_get_data( tag_data.tag_id, subrange, pointer_buffer, size_buffer );
CHK_MB_ERR_2(rval, tables + 1, status);
- // Convert MBEntityHandles to file ids
- if (mb_data_type == MB_TYPE_HANDLE) {
- assert( size < bufferSize );
- memcpy( dataBuffer, ptr, size );
- convert_handle_tag( reinterpret_cast<MBEntityHandle*>(dataBuffer),
- size / sizeof(MBEntityHandle) );
- ptr = dataBuffer;
+ // calculate end indices from sizes, and process tag data
+ for (size_t j = 0; j < count; ++j) {
+ const size_t size = size_buffer[j];
+ offset += size / type_size;
+ offset_buffer[j] = offset - 1;
+
+ // if space in data buffer, add current tag value and continue
+ assert(size_buffer[j] >= 0);
+ const void* ptr = pointer_buffer[j];
+
+ // flush buffer if need more room
+ if (bytes + size > data_buffer_size) {
+ // write out tag data buffer
+ if (bytes) { // bytes might be zero if tag value is larger than buffer
+ mhdf_writeSparseTagValues( tables[1], data_offset, bytes / type_size, hdf_type, data_buffer, &status );
+ CHK_MHDF_ERR_2(status, tables + 1);
+ data_offset += bytes / type_size;
+ bytes = 0;
+ }
+ }
+
+ // special case: if tag data is larger than buffer write it w/out buffering
+ if (size > data_buffer_size) {
+ if (mb_data_type == MB_TYPE_HANDLE) {
+ std::vector<MBEntityHandle> tmp_storage(size/sizeof(MBEntityHandle));
+ convert_handle_tag( reinterpret_cast<const MBEntityHandle*>(ptr),
+ &tmp_storage[0], tmp_storage.size() );
+ ptr = &tmp_storage[0];
+ }
+ mhdf_writeSparseTagValues( tables[1], data_offset, size / type_size, hdf_type, ptr, &status );
+ CHK_MHDF_ERR_2(status, tables + 1);
+ data_offset += size / type_size;
+ }
+ // otherwise copy data into buffer to be written during a later ieration
+ else {
+ if (mb_data_type == MB_TYPE_HANDLE)
+ convert_handle_tag( reinterpret_cast<const MBEntityHandle*>(ptr),
+ reinterpret_cast<MBEntityHandle*>(data_buffer + bytes),
+ size/sizeof(MBEntityHandle) );
+ else
+ memcpy( data_buffer + bytes, pointer_buffer[j], size );
+ bytes += size;
+ }
}
- mhdf_writeSparseTagValues( tables[1], data_offset, size/type_size, hdf_type, ptr, &status );
+ // write offsets
+ mhdf_writeSparseTagIndices( tables[2], offset_offset, count, H5T_NATIVE_LONG, offset_buffer, &status );
CHK_MHDF_ERR_2(status, tables + 1);
-
- idx = data_offset + size/type_size - 1;
- mhdf_writeSparseTagIndices( tables[2], indx_offset, 1, H5T_NATIVE_ULONG, &idx, &status );
+ offset_offset += count;
+ }
+ assert( offset_offset == tag_data.offset + tag_data.range.size() );
+
+ // flush data buffer
+ if (bytes) {
+ // write out tag data buffer
+ mhdf_writeSparseTagValues( tables[1], data_offset, bytes / type_size, hdf_type, data_buffer, &status );
CHK_MHDF_ERR_2(status, tables + 1);
- data_offset += size/type_size;
+ data_offset += bytes / type_size;
}
-
+ assert( offset == data_offset );
+
mhdf_closeData( filePtr, tables[1], &status );
CHK_MHDF_ERR_1(status, tables[2]);
mhdf_closeData( filePtr, tables[2], &status );
Modified: MOAB/trunk/WriteHDF5.hpp
===================================================================
--- MOAB/trunk/WriteHDF5.hpp 2008-04-02 19:11:26 UTC (rev 1736)
+++ MOAB/trunk/WriteHDF5.hpp 2008-04-03 00:24:41 UTC (rev 1737)
@@ -311,7 +311,10 @@
* the file or at least one of the handles is NULL (zero). false
* otherwise
*/
- bool convert_handle_tag( MBEntityHandle* data, size_t count );
+ bool convert_handle_tag( MBEntityHandle* data, size_t count ) const;
+ bool convert_handle_tag( const MBEntityHandle* source,
+ MBEntityHandle* dest,
+ size_t count ) const;
/** Get IDs of adjacent entities.
*
Modified: MOAB/trunk/mhdf/include/mhdf.h
===================================================================
--- MOAB/trunk/mhdf/include/mhdf.h 2008-04-02 19:11:26 UTC (rev 1736)
+++ MOAB/trunk/mhdf/include/mhdf.h 2008-04-03 00:24:41 UTC (rev 1737)
@@ -1842,11 +1842,17 @@
* Open the file objects containing all sparse data for a given tag in. The
* sparse data is stored in a pair of objects. The first is a vector of
* global IDs. The second is a vector of tag values for each entity specified
- * in the list of global IDs.
+ * in the list of global IDs. For variable-length tags, a third table
+ * containing end offsets for each tag value is returned in the third
+ * position of the output hid_t array (see mhdf_readSparseTagIndices.)
*
*\param file_handle The file.
*\param tag_name The tag.
- *\param num_values_out The number of entities for which tag values are stored.
+ *\param num_entity_out The number of entities for which there are tag values.
+ *\param num_values_out The number of data values. For fixed-length tags,
+ * this is the same as num_entity_out. For
+ * variable-length tags, it is the total number of
+ * values in the data table.
*\param entities_and_values_out The handles to the pair of file objects.
* The first is the vector of global IDs. The second
* is the list of corresponding tag values. The third
@@ -1858,6 +1864,7 @@
void
mhdf_openSparseTagData( mhdf_FileHandle file_handle,
const char* tag_name,
+ long* num_entity_out,
long* num_values_out,
hid_t entities_and_values_out[3],
mhdf_Status* status );
Modified: MOAB/trunk/mhdf/src/tags.c
===================================================================
--- MOAB/trunk/mhdf/src/tags.c 2008-04-02 19:11:26 UTC (rev 1736)
+++ MOAB/trunk/mhdf/src/tags.c 2008-04-03 00:24:41 UTC (rev 1737)
@@ -1405,7 +1405,7 @@
H5Gclose( elem_id );
*num_values_out = (long)size;
- if (data_id > 0)
+ if (data_id >= 0)
mhdf_setOkay( status );
API_END_H( 1 );
@@ -1601,12 +1601,13 @@
void
mhdf_openSparseTagData( mhdf_FileHandle file_handle,
const char* tag_name,
+ long* num_entity_out,
long* num_values_out,
hid_t handles_out[3],
mhdf_Status* status )
{
hid_t tag_id, index_id, data_id, offset_id = -1;
- hsize_t size1, size2;
+ hsize_t num_ent, data_size, num_data;
int rval;
unsigned idx;
API_BEGIN;
@@ -1614,14 +1615,14 @@
tag_id = get_tag( file_handle, tag_name, status );
if (tag_id < 0) return ;
- index_id = mhdf_open_table( tag_id, SPARSE_ENTITY_NAME, 1, &size1, status );
+ index_id = mhdf_open_table( tag_id, SPARSE_ENTITY_NAME, 1, &num_ent, status );
if (index_id < 0)
{
H5Gclose( tag_id );
return ;
}
- data_id = mhdf_open_table( tag_id, SPARSE_VALUES_NAME, 1, &size2, status );
+ data_id = mhdf_open_table( tag_id, SPARSE_VALUES_NAME, 1, &data_size, status );
if (data_id < 0)
{
H5Gclose( tag_id );
@@ -1638,12 +1639,9 @@
return ;
}
- /* If variable length...
- * Replace if data table (size2) with size of offset table. For variable-length
- * tags, the index and offset tables should be the same size, but the size of
- * the data table could be anything. */
+ /* If variable length... */
if (rval) {
- offset_id = mhdf_open_table( tag_id, TAG_VAR_INDICES, 1, &size2, status );
+ offset_id = mhdf_open_table( tag_id, TAG_VAR_INDICES, 1, &num_data, status );
if (offset_id < 0) {
H5Gclose( tag_id );
H5Dclose( index_id );
@@ -1651,9 +1649,13 @@
return ;
}
}
+ /* Otherwise the number of values is the same as the size of the data table */
+ else {
+ num_data = data_size;
+ }
H5Gclose( tag_id );
- if (size1 != size2)
+ if (num_ent != num_data)
{
mhdf_setFail( status, "Data length mismatch for sparse tag data -- invalid file.");
if (offset_id >= 0)
@@ -1662,7 +1664,9 @@
H5Dclose( data_id );
return ;
}
- *num_values_out = (long)size1;
+ *num_entity_out = (long)num_ent;
+ if (num_values_out)
+ *num_values_out = (long)data_size;
handles_out[0] = index_id;
handles_out[1] = data_id;
Modified: MOAB/trunk/test/h5file/h5varlen.cpp
===================================================================
--- MOAB/trunk/test/h5file/h5varlen.cpp 2008-04-02 19:11:26 UTC (rev 1736)
+++ MOAB/trunk/test/h5file/h5varlen.cpp 2008-04-03 00:24:41 UTC (rev 1737)
@@ -11,6 +11,8 @@
void test_var_length_data_big();
+void test_var_length_big_data();
+
void test_var_length_opaque();
void test_var_length_mesh_data();
@@ -21,6 +23,8 @@
void test_var_length_default_opaque();
+void test_var_length_handle_tag();
+
void create_mesh( MBInterface& mb );
void create_big_mesh( MBInterface& mb );
@@ -49,12 +53,14 @@
int err_count = 0;
err_count += RUN_TEST( test_var_length_no_data );
err_count += RUN_TEST( test_var_length_data );
- err_count += RUN_TEST( test_var_length_data_big );
+ err_count += RUN_TEST( test_var_length_big_data );
err_count += RUN_TEST( test_var_length_opaque );
err_count += RUN_TEST( test_var_length_mesh_data );
err_count += RUN_TEST( test_var_length_default_data );
err_count += RUN_TEST( test_var_length_mesh_opaque );
err_count += RUN_TEST( test_var_length_default_opaque );
+ err_count += RUN_TEST( test_var_length_handle_tag );
+ err_count += RUN_TEST( test_var_length_data_big );
return err_count;
}
@@ -162,11 +168,91 @@
}
+void calculate_big_value( MBInterface& moab, MBEntityHandle vert, size_t size, double* data )
+{
+ // Make values like Fibonacci numbers, except use X and Y coords
+ // rather than 0 and 1 as first two values.
+
+ CHECK( size >= 3 );
+ MBErrorCode rval = moab.get_coords( &vert, 1, data );
+ CHECK_ERR(rval);
+
+ for (size_t j = 2; j < size; ++j)
+ data[j] = data[j-2] + data[j-1];
+ CHECK_ERR(rval);
+}
+
+
+void test_var_length_big_data()
+{
+ MBErrorCode rval;
+ MBCore moab1, moab2;
+ MBInterface &mb1 = moab1, &mb2 = moab2;
+ MBTag tag;
+
+ create_mesh( mb1 );
+ rval = mb1.tag_create_variable_length( "test_tag", MB_TAG_SPARSE, MB_TYPE_DOUBLE, tag );
+ CHECK_ERR( rval );
+
+ // choose 3 vertices upon which to set data
+ MBRange range;
+ rval = mb1.get_entities_by_type( 0, MBVERTEX, range );
+ CHECK_ERR(rval);
+ MBEntityHandle verts[3] = { range.front(),
+ *(range.begin() += range.size()/3),
+ *(range.begin() += 2*range.size()/3) };
+
+ // set 1-millon value tag data on three vertices
+ std::vector<double> data(1000000);
+ for (int i = 0; i < 3; ++i) {
+ calculate_big_value( mb1, verts[i], data.size(), &data[0] );
+ const void* ptr = &data[0];
+ const int size = data.size() * sizeof(double);
+ rval = mb1.tag_set_data( tag, verts + i, 1, &ptr, &size );
+ CHECK_ERR(rval);
+ }
+
+ read_write( "test_var_length_big_data.h5m", mb1, mb2 );
+ compare_tags( "test_tag", mb1, mb2 );
+
+ // check 3 tagged vertices
+ rval = mb2.tag_get_handle( "test_tag", tag );
+ CHECK_ERR(rval);
+ range.clear();
+ rval = mb2.get_entities_by_type_and_tag( 0, MBVERTEX, &tag, 0, 1, range, MBInterface::UNION );
+ CHECK_ERR(rval);
+ CHECK_EQUAL( range.size(), (MBEntityHandle)3 );
+
+ // check tag values
+ for (MBRange::const_iterator i = range.begin(); i != range.end(); ++i) {
+ // calculate expected value
+ const MBEntityHandle h = *i;
+ calculate_big_value( mb2, h, data.size(), &data[0] );
+
+ // get actual value
+ const void* ptr;
+ int size;
+ rval = mb2.tag_get_data( tag, &h, 1, &ptr, &size );
+ CHECK_ERR(rval);
+ CHECK_EQUAL( data.size() * sizeof(double), (size_t)size );
+
+ // compare values
+ const double* act_data = reinterpret_cast<const double*>(ptr);
+ int wrong_count = 0;
+ for (size_t j = 0; j < data.size(); ++j)
+ if (act_data[j] != data[j])
+ ++wrong_count;
+ CHECK_EQUAL( 0, wrong_count );
+ }
+}
+
+
+
void test_var_length_opaque()
{
MBCore moab;
create_mesh( moab );
- test_var_length_data_common( "tst_var_length_data_opaque.h5m", moab, true );
+ test_var_length_data_common( "test_var_length_opaque.h5m", moab, true );
}
void test_global_value_common( bool mesh_value )
@@ -360,6 +446,93 @@
test_global_opaque_common( false );
}
+
+void test_var_length_handle_tag()
+{
+ MBErrorCode rval;
+ MBCore moab1, moab2;
+ MBInterface &mb1 = moab1, &mb2 = moab2;
+ MBTag tag;
+ MBRange::const_iterator i;
+
+ create_mesh( mb1 );
+ rval = mb1.tag_create_variable_length( "test_tag", MB_TAG_SPARSE, MB_TYPE_HANDLE, tag );
+ CHECK_ERR( rval );
+
+ // Get all entities
+ MBRange range;
+ rval = mb1.get_entities_by_handle( 0, range );
+ CHECK_ERR(rval);
+
+ // For each entity, if it is a vertex store its own handle
+ // in its tag. Otherwise store the element connectivity list
+ // in the tag. Skip entity sets.
+ MBEntityHandle num_tagged_entities = 0;
+ for (i = range.begin(); i != range.end(); ++i) {
+ MBEntityHandle h = *i;
+ MBEntityType type = mb1.type_from_handle( h );
+ if (type == MBVERTEX) {
+ const int size = sizeof(MBEntityHandle);
+ const void* ptr = &h;
+ rval = mb1.tag_set_data( tag, &h, 1, &ptr, &size );
+ CHECK_ERR(rval);
+ ++num_tagged_entities;
+ }
+ else if (type != MBENTITYSET) {
+ int size = 0;
+ const MBEntityHandle* conn = 0;
+ rval = mb1.get_connectivity( h, conn, size );
+ CHECK_ERR(rval);
+ size *= sizeof(MBEntityHandle);
+ const void* ptr = conn;
+ rval = mb1.tag_set_data( tag, &h, 1, &ptr, &size );
+ CHECK_ERR(rval);
+ ++num_tagged_entities;
+ }
+ }
+
+ read_write( "test_var_length_handle_tag.h5m", mb1, mb2 );
+ compare_tags( "test_tag", mb1, mb2 );
+
+ // check number of tagged entities
+ rval = mb2.tag_get_handle( "test_tag", tag );
+ CHECK_ERR(rval);
+ range.clear();
+ for (MBEntityType t = MBVERTEX; t != MBENTITYSET; ++t) {
+ rval = mb2.get_entities_by_type_and_tag( 0, t, &tag, 0, 1, range, MBInterface::UNION );
+ CHECK_ERR(rval);
+ }
+ CHECK_EQUAL( num_tagged_entities, range.size() );
+
+ // check tag values
+ for (i = range.begin(); i != range.end(); ++i) {
+ MBEntityHandle h = *i;
+
+ const void* ptr;
+ int size;
+ rval = mb2.tag_get_data( tag, &h, 1, &ptr, &size );
+ CHECK_ERR(rval);
+
+ CHECK_EQUAL( (size_t)0, (size_t)size % sizeof(MBEntityHandle) );
+ size /= sizeof(MBEntityHandle);
+ const MBEntityHandle* handles = reinterpret_cast<const MBEntityHandle*>(ptr);
+
+ if (mb2.type_from_handle(h) == MBVERTEX) {
+ CHECK_EQUAL( 1, size );
+ CHECK_EQUAL( h, *handles );
+ }
+ else {
+ int len;
+ const MBEntityHandle* conn;
+ rval = mb2.get_connectivity( h, conn, len );
+ CHECK_ERR(rval);
+ CHECK_EQUAL( len, size );
+ for (int j = 0; j < len; ++j)
+ CHECK_EQUAL( conn[j], handles[j] );
+ }
+ }
+}
+
void create_structured_quad_mesh( MBInterface& mb, int x, int y )
{
MBErrorCode rval;
@@ -394,7 +567,7 @@
void create_big_mesh( MBInterface& mb )
{
- create_structured_quad_mesh( mb, 100, 100 );
+ create_structured_quad_mesh( mb, 1000, 700 );
}
void compare_tags( const char* name, MBInterface& mb1, MBInterface& mb2 )
@@ -451,4 +624,3 @@
remove( filename );
CHECK_ERR(rval);
}
-
More information about the moab-dev
mailing list