Skip to content

Commit d43f637

Browse files
committedDec 7, 2018
[MDAL] update MDAL to 0.1.0 (new API)
1 parent d6701f2 commit d43f637

23 files changed

+1928
-699
lines changed
 

‎external/mdal/api/mdal.h

Lines changed: 105 additions & 33 deletions
Original file line numberDiff line numberDiff line change
@@ -63,8 +63,9 @@ enum MDAL_Status
6363
Warn_NodeNotUnique
6464
};
6565

66-
//! Mesh
6766
typedef void *MeshH;
67+
typedef void *MeshVertexIteratorH;
68+
typedef void *MeshFaceIteratorH;
6869
typedef void *DatasetGroupH;
6970
typedef void *DatasetH;
7071

@@ -78,7 +79,10 @@ MDAL_EXPORT MDAL_Status MDAL_LastStatus();
7879
/// MESH
7980
///////////////////////////////////////////////////////////////////////////////////////
8081

81-
//! Loads mesh file. On error see MDAL_LastStatus for error type This effectively loads whole mesh in-memory
82+
/**
83+
* Loads mesh file. On error see MDAL_LastStatus for error type
84+
* This may effectively load whole mesh in-memory for some providers
85+
*/
8286
MDAL_EXPORT MeshH MDAL_LoadMesh( const char *meshFile );
8387
//! Closes mesh, frees the memory
8488
MDAL_EXPORT void MDAL_CloseMesh( MeshH mesh );
@@ -88,20 +92,18 @@ MDAL_EXPORT void MDAL_CloseMesh( MeshH mesh );
8892
* not thread-safe and valid only till next call
8993
*/
9094
MDAL_EXPORT const char *MDAL_M_projection( MeshH mesh );
95+
96+
/**
97+
* Returns mesh extent in native projection
98+
* Returns NaN on error
99+
*/
100+
MDAL_EXPORT void MDAL_M_extent( MeshH mesh, double *minX, double *maxX, double *minY, double *maxY );
91101
//! Returns vertex count for the mesh
92102
MDAL_EXPORT int MDAL_M_vertexCount( MeshH mesh );
93-
//! Returns vertex X coord for the mesh
94-
MDAL_EXPORT double MDAL_M_vertexXCoordinatesAt( MeshH mesh, int index );
95-
//! Returns vertex Y coord for the mesh
96-
MDAL_EXPORT double MDAL_M_vertexYCoordinatesAt( MeshH mesh, int index );
97-
//! Returns vertex Z coord for the mesh
98-
MDAL_EXPORT double MDAL_M_vertexZCoordinatesAt( MeshH mesh, int index );
99103
//! Returns face count for the mesh
100104
MDAL_EXPORT int MDAL_M_faceCount( MeshH mesh );
101-
//! Returns number of vertices face consist of, e.g. 3 for triangle
102-
MDAL_EXPORT int MDAL_M_faceVerticesCountAt( MeshH mesh, int index );
103-
//! Returns vertex index for face
104-
MDAL_EXPORT int MDAL_M_faceVerticesIndexAt( MeshH mesh, int face_index, int vertex_index );
105+
//! Returns maximum number of vertices face can consist of, e.g. 4 for regular quad mesh
106+
MDAL_EXPORT int MDAL_M_faceVerticesMaximumCount( MeshH mesh );
105107

106108
/**
107109
* Loads dataset file. On error see MDAL_LastStatus for error type.
@@ -110,17 +112,74 @@ MDAL_EXPORT int MDAL_M_faceVerticesIndexAt( MeshH mesh, int face_index, int vert
110112
* can be freed manually with MDAL_CloseDataset if needed
111113
*/
112114
MDAL_EXPORT void MDAL_M_LoadDatasets( MeshH mesh, const char *datasetFile );
113-
114115
//! Returns dataset groups count
115116
MDAL_EXPORT int MDAL_M_datasetGroupCount( MeshH mesh );
116-
117117
//! Returns dataset group handle
118118
MDAL_EXPORT DatasetGroupH MDAL_M_datasetGroup( MeshH mesh, int index );
119119

120+
///////////////////////////////////////////////////////////////////////////////////////
121+
/// MESH VERTICES
122+
///////////////////////////////////////////////////////////////////////////////////////
123+
124+
/**
125+
* Returns iterator to the mesh vertices
126+
* For some formats this may effectively load all vertices in-memory until iterator is closed
127+
*/
128+
MDAL_EXPORT MeshVertexIteratorH MDAL_M_vertexIterator( MeshH mesh );
129+
130+
/**
131+
* Returns vertices from iterator for the mesh
132+
* \param iterator mesh data iterator
133+
* \param verticesCount maximum number or vertices to be written to buffer
134+
* \param coordinates must be allocated to 3* verticesCount items to store x1, y1, z1, ..., xN, yN, zN coordinates
135+
* \returns number of vertices written in the buffer
136+
*/
137+
MDAL_EXPORT int MDAL_VI_next( MeshVertexIteratorH iterator, int verticesCount, double *coordinates );
138+
139+
//! Closes mesh data iterator, frees the memory
140+
MDAL_EXPORT void MDAL_VI_close( MeshVertexIteratorH iterator );
141+
142+
///////////////////////////////////////////////////////////////////////////////////////
143+
/// MESH FACES
144+
///////////////////////////////////////////////////////////////////////////////////////
145+
146+
/**
147+
* Returns iterator to the mesh faces
148+
* For some formats this may effectively load all faces in-memory until iterator is closed
149+
*/
150+
MDAL_EXPORT MeshFaceIteratorH MDAL_M_faceIterator( MeshH mesh );
151+
152+
/**
153+
* Returns next faces from iterator for the mesh
154+
*
155+
* Reading stops when vertexIndicesBuffer capacity is full / faceOffsetsBuffer
156+
* capacity is full / end of faces is reached, whatever comes first
157+
*
158+
* \param iterator mesh data iterator
159+
* \param faceOffsetsBufferLen size of faceOffsetsBuffer, minimum 1
160+
* \param faceOffsetsBuffer allocated array to store face offset in vertexIndicesBuffer for given face.
161+
* To find number of vertices of face i, calculate faceOffsetsBuffer[i] - faceOffsetsBuffer[i-1]
162+
* \param vertexIndicesBufferLen size of vertexIndicesBuffer, minimum is MDAL_M_faceVerticesMaximumCount()
163+
* \param vertexIndicesBuffer writes vertex indexes for faces
164+
* faceOffsetsBuffer[i-1] is index where the vertices for face i begins,
165+
* \returns number of faces written in the buffer
166+
*/
167+
MDAL_EXPORT int MDAL_FI_next( MeshFaceIteratorH iterator,
168+
int faceOffsetsBufferLen,
169+
int *faceOffsetsBuffer,
170+
int vertexIndicesBufferLen,
171+
int *vertexIndicesBuffer );
172+
173+
//! Closes mesh data iterator, frees the memory
174+
MDAL_EXPORT void MDAL_FI_close( MeshFaceIteratorH iterator );
175+
120176
///////////////////////////////////////////////////////////////////////////////////////
121177
/// DATASET GROUPS
122178
///////////////////////////////////////////////////////////////////////////////////////
123179

180+
//! Returns dataset parent mesh
181+
MDAL_EXPORT MeshH MDAL_G_mesh( DatasetGroupH group );
182+
124183
//! Returns dataset count in group
125184
MDAL_EXPORT int MDAL_G_datasetCount( DatasetGroupH group );
126185

@@ -154,6 +213,12 @@ MDAL_EXPORT bool MDAL_G_hasScalarData( DatasetGroupH group );
154213
//! Whether dataset is on vertices
155214
MDAL_EXPORT bool MDAL_G_isOnVertices( DatasetGroupH group );
156215

216+
/**
217+
* Returns the min and max values of the group
218+
* Returns NaN on error
219+
*/
220+
MDAL_EXPORT void MDAL_G_minimumMaximum( DatasetGroupH group, double *min, double *max );
221+
157222
///////////////////////////////////////////////////////////////////////////////////////
158223
/// DATASETS
159224
///////////////////////////////////////////////////////////////////////////////////////
@@ -167,32 +232,39 @@ MDAL_EXPORT double MDAL_D_time( DatasetH dataset );
167232
//! Returns number of values
168233
MDAL_EXPORT int MDAL_D_valueCount( DatasetH dataset );
169234

170-
/**
171-
* Returns scalar value associated with the index from the dataset
172-
* for nodata return numeric_limits<double>:quiet_NaN
173-
*/
174-
MDAL_EXPORT double MDAL_D_value( DatasetH dataset, int valueIndex );
235+
//! Returns whether dataset is valid
236+
MDAL_EXPORT bool MDAL_D_isValid( DatasetH dataset );
175237

176-
/**
177-
* Returns X value associated with the index from the vector dataset
178-
* for nodata return numeric_limits<double>:quiet_NaN
179-
*/
180-
MDAL_EXPORT double MDAL_D_valueX( DatasetH dataset, int valueIndex );
238+
//! Data type to be returned by MDAL_D_data
239+
enum MDAL_DataType
240+
{
241+
SCALAR_DOUBLE, //!< Double value for scalar datasets
242+
VECTOR_2D_DOUBLE, //!< Double, double value for vector datasets
243+
ACTIVE_INTEGER //!< Integer, active flag for dataset faces. Some formats support switching off the element for particular timestep.
244+
};
181245

182246
/**
183-
* Returns Y value associated with the index from the vector dataset
184-
* for nodata return numeric_limits<double>:quiet_NaN
247+
* Populates buffer with values from the dataset
248+
* for nodata, returned is numeric_limits<double>::quiet_NaN
249+
*
250+
* \param dataset handle to dataset
251+
* \param indexStart index of face/vertex to start reading of values to the buffer
252+
* \param count number of values to be written to the buffer
253+
* \param dataType type of values to be written to the buffer
254+
* \param buffer output array to be populated with the values. must be already allocated
255+
* For SCALAR_DOUBLE, the minimum size must be valuesCount * size_of(double)
256+
* For VECTOR_2D_DOUBLE, the minimum size must be valuesCount * 2 * size_of(double).
257+
* Values are returned as x1, y1, x2, y2, ..., xN, yN
258+
* For ACTIVE_INTEGER, the minimum size must be valuesCount * size_of(int)
259+
* \returns number of values written to buffer. If return value != count requested, see MDAL_LastStatus() for error type
185260
*/
186-
MDAL_EXPORT double MDAL_D_valueY( DatasetH dataset, int valueIndex );
261+
MDAL_EXPORT int MDAL_D_data( DatasetH dataset, int indexStart, int count, MDAL_DataType dataType, void *buffer );
187262

188263
/**
189-
* Whether element is active - should be taken into account
190-
* Some formats support switching off the element for particular timestep
264+
* Returns the min and max values of the dataset
265+
* Returns NaN on error
191266
*/
192-
MDAL_EXPORT bool MDAL_D_active( DatasetH dataset, int faceIndex );
193-
194-
//! Returns whether dataset is valid
195-
MDAL_EXPORT bool MDAL_D_isValid( DatasetH dataset );
267+
MDAL_EXPORT void MDAL_D_minimumMaximum( DatasetH dataset, double *min, double *max );
196268

197269
#ifdef __cplusplus
198270
}

‎external/mdal/frmts/mdal_2dm.cpp

Lines changed: 57 additions & 48 deletions
Original file line numberDiff line numberDiff line change
@@ -17,6 +17,46 @@
1717
#include "mdal.h"
1818
#include "mdal_utils.hpp"
1919

20+
MDAL::Mesh2dm::Mesh2dm( size_t verticesCount,
21+
size_t facesCount,
22+
size_t faceVerticesMaximumCount,
23+
MDAL::BBox extent,
24+
const std::string &uri,
25+
const std::map<size_t, size_t> vertexIDtoIndex )
26+
: MemoryMesh( verticesCount, facesCount, faceVerticesMaximumCount, extent, uri )
27+
, mVertexIDtoIndex( vertexIDtoIndex )
28+
{
29+
}
30+
31+
MDAL::Mesh2dm::~Mesh2dm() = default;
32+
33+
bool _parse_vertex_id_gaps( std::map<size_t, size_t> &vertexIDtoIndex, size_t vertexIndex, size_t vertexID, MDAL_Status *status )
34+
{
35+
if ( vertexIndex == vertexID )
36+
return false;
37+
38+
std::map<size_t, size_t>::iterator search = vertexIDtoIndex.find( vertexID );
39+
if ( search != vertexIDtoIndex.end() )
40+
{
41+
if ( status ) *status = MDAL_Status::Warn_ElementNotUnique;
42+
return true;
43+
}
44+
45+
vertexIDtoIndex[vertexID] = vertexIndex;
46+
return false;
47+
}
48+
49+
size_t MDAL::Mesh2dm::vertexIndex( size_t vertexID ) const
50+
{
51+
auto ni2i = mVertexIDtoIndex.find( vertexID );
52+
if ( ni2i != mVertexIDtoIndex.end() )
53+
{
54+
return ni2i->second; // convert from ID to index
55+
}
56+
return vertexID;
57+
}
58+
59+
2060
MDAL::Loader2dm::Loader2dm( const std::string &meshFile ):
2161
mMeshFile( meshFile )
2262
{
@@ -71,7 +111,6 @@ std::unique_ptr<MDAL::Mesh> MDAL::Loader2dm::load( MDAL_Status *status )
71111

72112
size_t faceIndex = 0;
73113
size_t vertexIndex = 0;
74-
std::map<size_t, size_t> faceIDtoIndex;
75114
std::map<size_t, size_t> vertexIDtoIndex;
76115

77116
while ( std::getline( in, line ) )
@@ -81,20 +120,11 @@ std::unique_ptr<MDAL::Mesh> MDAL::Loader2dm::load( MDAL_Status *status )
81120
chunks = split( line, " ", SplitBehaviour::SkipEmptyParts );
82121
assert( faceIndex < faceCount );
83122

84-
size_t elemID = toSizeT( chunks[1] );
85-
86-
std::map<size_t, size_t>::iterator search = faceIDtoIndex.find( elemID );
87-
if ( search != faceIDtoIndex.end() )
88-
{
89-
if ( status ) *status = MDAL_Status::Warn_ElementNotUnique;
90-
continue;
91-
}
92-
faceIDtoIndex[elemID] = faceIndex;
93123
Face &face = faces[faceIndex];
94124
face.resize( 4 );
95125
// Right now we just store node IDs here - we will convert them to node indices afterwards
96126
for ( size_t i = 0; i < 4; ++i )
97-
face[i] = toSizeT( chunks[i + 2] );
127+
face[i] = toSizeT( chunks[i + 2] ) - 1; // 2dm is numbered from 1
98128

99129
faceIndex++;
100130
}
@@ -103,21 +133,12 @@ std::unique_ptr<MDAL::Mesh> MDAL::Loader2dm::load( MDAL_Status *status )
103133
chunks = split( line, " ", SplitBehaviour::SkipEmptyParts );
104134
assert( faceIndex < faceCount );
105135

106-
size_t elemID = toSizeT( chunks[1] );
107-
108-
std::map<size_t, size_t>::iterator search = faceIDtoIndex.find( elemID );
109-
if ( search != faceIDtoIndex.end() )
110-
{
111-
if ( status ) *status = MDAL_Status::Warn_ElementNotUnique;
112-
continue;
113-
}
114-
faceIDtoIndex[elemID] = faceIndex;
115136
Face &face = faces[faceIndex];
116137
face.resize( 3 );
117138
// Right now we just store node IDs here - we will convert them to node indices afterwards
118139
for ( size_t i = 0; i < 3; ++i )
119140
{
120-
face[i] = toSizeT( chunks[i + 2] );
141+
face[i] = toSizeT( chunks[i + 2] ) - 1; // 2dm is numbered from 1
121142
}
122143

123144
faceIndex++;
@@ -132,31 +153,16 @@ std::unique_ptr<MDAL::Mesh> MDAL::Loader2dm::load( MDAL_Status *status )
132153
chunks = split( line, " ", SplitBehaviour::SkipEmptyParts );
133154
assert( faceIndex < faceCount );
134155

135-
size_t elemID = toSizeT( chunks[1] );
136-
137-
std::map<size_t, size_t>::iterator search = faceIDtoIndex.find( elemID );
138-
if ( search != faceIDtoIndex.end() )
139-
{
140-
if ( status ) *status = MDAL_Status::Warn_ElementNotUnique;
141-
continue;
142-
}
143-
faceIDtoIndex[elemID] = faceIndex;
156+
//size_t elemID = toSizeT( chunks[1] );
144157
assert( false ); //TODO mark element as unusable
145158

146159
faceIndex++;
147160
}
148161
else if ( startsWith( line, "ND" ) )
149162
{
150163
chunks = split( line, " ", SplitBehaviour::SkipEmptyParts );
151-
size_t nodeID = toSizeT( chunks[1] );
152-
153-
std::map<size_t, size_t>::iterator search = vertexIDtoIndex.find( nodeID );
154-
if ( search != vertexIDtoIndex.end() )
155-
{
156-
if ( status ) *status = MDAL_Status::Warn_NodeNotUnique;
157-
continue;
158-
}
159-
vertexIDtoIndex[nodeID] = vertexIndex;
164+
size_t nodeID = toSizeT( chunks[1] ) - 1; // 2dm is numbered from 1
165+
_parse_vertex_id_gaps( vertexIDtoIndex, vertexIndex, nodeID, status );
160166
assert( vertexIndex < vertexCount );
161167
Vertex &vertex = vertices[vertexIndex];
162168
vertex.x = toDouble( chunks[2] );
@@ -178,25 +184,28 @@ std::unique_ptr<MDAL::Mesh> MDAL::Loader2dm::load( MDAL_Status *status )
178184
{
179185
face[nd] = ni2i->second; // convert from ID to index
180186
}
181-
else
187+
else if ( vertices.size() < nodeID )
182188
{
183-
assert( false ); //TODO mark element as unusable
184-
185189
if ( status ) *status = MDAL_Status::Warn_ElementWithInvalidNode;
186190
}
187191
}
188-
189192
//TODO check validity of the face
190193
//check that we have distinct nodes
191194
}
192195

193-
std::unique_ptr< Mesh > mesh( new Mesh );
194-
mesh->uri = mMeshFile;
196+
std::unique_ptr< MemoryMesh > mesh(
197+
new Mesh2dm(
198+
vertices.size(),
199+
faces.size(),
200+
4, //maximum quads
201+
computeExtent( vertices ),
202+
mMeshFile,
203+
vertexIDtoIndex
204+
)
205+
);
195206
mesh->faces = faces;
196207
mesh->vertices = vertices;
197-
mesh->faceIDtoIndex = faceIDtoIndex;
198-
mesh->vertexIDtoIndex = vertexIDtoIndex;
199-
mesh->addBedElevationDataset();
208+
mesh->addBedElevationDataset( vertices, faces );
200209

201210
return mesh;
202211
}

‎external/mdal/frmts/mdal_2dm.hpp

Lines changed: 26 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -10,10 +10,36 @@
1010
#include <memory>
1111

1212
#include "mdal_data_model.hpp"
13+
#include "mdal_memory_data_model.hpp"
1314
#include "mdal.h"
1415

1516
namespace MDAL
1617
{
18+
class Mesh2dm: public MemoryMesh
19+
{
20+
public:
21+
Mesh2dm( size_t verticesCount,
22+
size_t facesCount,
23+
size_t faceVerticesMaximumCount,
24+
BBox extent,
25+
const std::string &uri,
26+
const std::map<size_t, size_t> vertexIDtoIndex
27+
);
28+
~Mesh2dm() override;
29+
30+
31+
//! Some formats supports gaps in the vertex indexing, but we return continuos array from MDAL
32+
//! for most of the formats this returns
33+
//! \param vertexID internal index/ID of the vertex that native format uses
34+
//! \returns index of the vertex in the continuous array of vertices we returned by readVertices()
35+
virtual size_t vertexIndex( size_t vertexID ) const;
36+
37+
private:
38+
// 2dm supports "gaps" in the mesh indexing
39+
// Store only the indices that have different index and ID
40+
// https://github.com/lutraconsulting/MDAL/issues/51
41+
std::map<size_t, size_t> mVertexIDtoIndex;
42+
};
1743

1844
class Loader2dm
1945
{

‎external/mdal/frmts/mdal_3di.cpp

Lines changed: 24 additions & 20 deletions
Original file line numberDiff line numberDiff line change
@@ -33,11 +33,11 @@ MDAL::CFDimensions MDAL::Loader3Di::populateDimensions()
3333
return dims;
3434
}
3535

36-
void MDAL::Loader3Di::populateFacesAndVertices( MDAL::Mesh *mesh )
36+
void MDAL::Loader3Di::populateFacesAndVertices( Vertices &vertices, Faces &faces )
3737
{
38-
assert( mesh );
38+
assert( vertices.empty() );
3939
size_t faceCount = mDimensions.size( CFDimensions::Face2D );
40-
mesh->faces.resize( faceCount );
40+
faces.resize( faceCount );
4141
size_t verticesInFace = mDimensions.size( CFDimensions::MaxVerticesInFace );
4242
size_t arrsize = faceCount * verticesInFace;
4343
std::map<std::string, size_t> xyToVertex2DId;
@@ -79,9 +79,9 @@ void MDAL::Loader3Di::populateFacesAndVertices( MDAL::Mesh *mesh )
7979
if ( it == xyToVertex2DId.end() )
8080
{
8181
// new vertex
82-
vertexId = mesh->vertices.size();
82+
vertexId = vertices.size();
8383
xyToVertex2DId[key] = vertexId;
84-
mesh->vertices.push_back( vertex );
84+
vertices.push_back( vertex );
8585
}
8686
else
8787
{
@@ -93,21 +93,21 @@ void MDAL::Loader3Di::populateFacesAndVertices( MDAL::Mesh *mesh )
9393

9494
}
9595

96-
mesh->faces[faceId] = face;
96+
faces[faceId] = face;
9797
}
9898

9999
// Only now we have number of vertices, since we identified vertices that
100100
// are used in multiple faces
101-
mDimensions.setDimension( CFDimensions::Vertex2D, mesh->vertices.size() );
101+
mDimensions.setDimension( CFDimensions::Vertex2D, vertices.size() );
102102
}
103103

104104
void MDAL::Loader3Di::addBedElevation( MDAL::Mesh *mesh )
105105
{
106106
assert( mesh );
107-
if ( mesh->faces.empty() )
107+
if ( 0 == mesh->facesCount() )
108108
return;
109109

110-
size_t faceCount = mesh->faces.size();
110+
size_t faceCount = mesh->facesCount();
111111

112112
// read Z coordinate of 3di computation nodes centers
113113
int ncidZ = mNcFile.getVarId( "Mesh2DFace_zcc" );
@@ -117,20 +117,24 @@ void MDAL::Loader3Di::addBedElevation( MDAL::Mesh *mesh )
117117
return; //error reading the array
118118

119119

120-
std::shared_ptr<DatasetGroup> group = std::make_shared< DatasetGroup >();
121-
group->isOnVertices = false;
122-
group->isScalar = true;
123-
group->setName( "Bed Elevation" );
124-
group->uri = mesh->uri;
125-
std::shared_ptr<MDAL::Dataset> dataset = std::make_shared< Dataset >();
126-
dataset->time = 0.0;
127-
dataset->values.resize( faceCount );
128-
dataset->active.resize( faceCount );
129-
dataset->parent = group.get();
120+
std::shared_ptr<DatasetGroup> group = std::make_shared< DatasetGroup >(
121+
mesh,
122+
mesh->uri(),
123+
"Bed Elevation"
124+
);
125+
126+
group->setIsOnVertices( false );
127+
group->setIsScalar( true );
128+
129+
std::shared_ptr<MDAL::MemoryDataset> dataset = std::make_shared< MemoryDataset >( group.get() );
130+
dataset->setTime( 0.0 );
131+
double *values = dataset->values();
130132
for ( size_t i = 0; i < faceCount; ++i )
131133
{
132-
dataset->values[i].x = MDAL::safeValue( coordZ[i], fillZ );
134+
values[i] = MDAL::safeValue( coordZ[i], fillZ );
133135
}
136+
dataset->setStatistics( MDAL::calculateStatistics( dataset ) );
137+
group->setStatistics( MDAL::calculateStatistics( group ) );
134138
group->datasets.push_back( dataset );
135139
mesh->datasetGroups.push_back( group );
136140
}

‎external/mdal/frmts/mdal_3di.hpp

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -43,8 +43,8 @@ namespace MDAL
4343

4444
private:
4545
CFDimensions populateDimensions() override;
46-
void populateFacesAndVertices( MDAL::Mesh *mesh ) override;
47-
void addBedElevation( MDAL::Mesh *mesh ) override;
46+
void populateFacesAndVertices( Vertices &vertices, Faces &faces ) override;
47+
void addBedElevation( Mesh *mesh ) override;
4848
std::string getCoordinateSystemVariableName() override;
4949
std::set<std::string> ignoreNetCDFVariables() override;
5050
std::string nameSuffix( CFDimensions::Type type ) override;

‎external/mdal/frmts/mdal_ascii_dat.cpp

Lines changed: 49 additions & 44 deletions
Original file line numberDiff line numberDiff line change
@@ -17,6 +17,7 @@
1717
#include "mdal_ascii_dat.hpp"
1818
#include "mdal.h"
1919
#include "mdal_utils.hpp"
20+
#include "mdal_2dm.hpp"
2021

2122
#include <math.h>
2223

@@ -69,10 +70,12 @@ void MDAL::LoaderAsciiDat::load( MDAL::Mesh *mesh, MDAL_Status *status )
6970
oldFormat = true;
7071
isVector = ( line == "VECTOR" );
7172

72-
group = std::make_shared< DatasetGroup >();
73-
group->uri = mDatFile;
74-
group->setName( name );
75-
group->isScalar = !isVector;
73+
group = std::make_shared< DatasetGroup >(
74+
mesh,
75+
mDatFile,
76+
name
77+
);
78+
group->setIsScalar( !isVector );
7679
}
7780
else
7881
EXIT_WITH_ERROR( MDAL_Status::Err_UnknownFormat );
@@ -83,10 +86,14 @@ void MDAL::LoaderAsciiDat::load( MDAL::Mesh *mesh, MDAL_Status *status )
8386
faceCentered = true;
8487

8588
if ( group )
86-
group->isOnVertices = !faceCentered;
89+
group->setIsOnVertices( !faceCentered );
8790

8891
while ( std::getline( in, line ) )
8992
{
93+
// Replace tabs by spaces,
94+
// since basement v.2.8 uses tabs instead of spaces (e.g. 'TS 0\t0.0')
95+
line = replace( line, "\t", " " );
96+
9097
std::vector<std::string> items = split( line, " ", SplitBehaviour::SkipEmptyParts );
9198
if ( items.size() < 1 )
9299
continue; // empty line?? let's skip it
@@ -95,13 +102,13 @@ void MDAL::LoaderAsciiDat::load( MDAL::Mesh *mesh, MDAL_Status *status )
95102
if ( cardType == "ND" && items.size() >= 2 )
96103
{
97104
size_t fileNodeCount = toSizeT( items[1] );
98-
if ( mesh->vertexIDtoIndex.size() != fileNodeCount )
105+
if ( mesh->verticesCount() != fileNodeCount )
99106
EXIT_WITH_ERROR( MDAL_Status::Err_IncompatibleMesh );
100107
}
101108
else if ( !oldFormat && cardType == "NC" && items.size() >= 2 )
102109
{
103110
size_t fileElemCount = toSizeT( items[1] );
104-
if ( mesh->faceIDtoIndex.size() != fileElemCount )
111+
if ( mesh->facesCount() != fileElemCount )
105112
EXIT_WITH_ERROR( MDAL_Status::Err_IncompatibleMesh );
106113
}
107114
else if ( !oldFormat && cardType == "OBJTYPE" )
@@ -118,11 +125,13 @@ void MDAL::LoaderAsciiDat::load( MDAL::Mesh *mesh, MDAL_Status *status )
118125
}
119126
isVector = cardType == "BEGVEC";
120127

121-
group = std::make_shared< DatasetGroup >();
122-
group->uri = mDatFile;
123-
group->setName( name );
124-
group->isScalar = !isVector;
125-
group->isOnVertices = !faceCentered;
128+
group = std::make_shared< DatasetGroup >(
129+
mesh,
130+
mDatFile,
131+
name
132+
);
133+
group->setIsScalar( !isVector );
134+
group->setIsOnVertices( !faceCentered );
126135
}
127136
else if ( !oldFormat && cardType == "ENDDS" )
128137
{
@@ -131,6 +140,7 @@ void MDAL::LoaderAsciiDat::load( MDAL::Mesh *mesh, MDAL_Status *status )
131140
debug( "ENDDS card for no active dataset!" );
132141
EXIT_WITH_ERROR( MDAL_Status::Err_UnknownFormat );
133142
}
143+
group->setStatistics( MDAL::calculateStatistics( group ) );
134144
mesh->datasetGroups.push_back( group );
135145
group.reset();
136146
}
@@ -179,6 +189,7 @@ void MDAL::LoaderAsciiDat::load( MDAL::Mesh *mesh, MDAL_Status *status )
179189
if ( !group || group->datasets.size() == 0 )
180190
EXIT_WITH_ERROR( MDAL_Status::Err_UnknownFormat );
181191

192+
group->setStatistics( MDAL::calculateStatistics( group ) );
182193
mesh->datasetGroups.push_back( group );
183194
group.reset();
184195
}
@@ -194,65 +205,62 @@ void MDAL::LoaderAsciiDat::readVertexTimestep(
194205
{
195206
assert( group );
196207

197-
size_t vertexCount = mesh->vertices.size();
198-
size_t faceCount = mesh->faces.size();
208+
size_t vertexCount = mesh->verticesCount();
209+
size_t faceCount = mesh->facesCount();
199210

200-
std::shared_ptr<MDAL::Dataset> dataset = std::make_shared< MDAL::Dataset >();
201-
dataset->time = t / 3600.; // TODO read TIMEUNITS
202-
dataset->values.resize( vertexCount );
203-
dataset->active.resize( faceCount );
204-
dataset->parent = group.get();
211+
std::shared_ptr<MDAL::MemoryDataset> dataset = std::make_shared< MDAL::MemoryDataset >( group.get() );
212+
dataset->setTime( t / 3600. ); // TODO read TIMEUNITS
205213

214+
int *active = dataset->active();
206215
// only for new format
207216
for ( size_t i = 0; i < faceCount; ++i )
208217
{
209218
if ( hasStatus )
210219
{
211220
std::string line;
212221
std::getline( stream, line );
213-
dataset->active[i] = toBool( line );
222+
active[i] = toBool( line );
214223
}
215-
else
216-
dataset->active[i] = true;
217224
}
218225

219-
for ( size_t i = 0; i < mesh->vertexIDtoIndex.size(); ++i )
226+
const Mesh2dm *m2dm = dynamic_cast<const Mesh2dm *>( mesh );
227+
double *values = dataset->values();
228+
for ( size_t i = 0; i < mesh->verticesCount(); ++i )
220229
{
221230
std::string line;
222231
std::getline( stream, line );
223232
std::vector<std::string> tsItems = split( line, " ", SplitBehaviour::SkipEmptyParts );
224233

225-
auto idx = mesh->vertexIDtoIndex.find( i + 1 ); // ID in 2dm are numbered from 1
226-
if ( idx == mesh->vertexIDtoIndex.end() )
227-
continue; // node ID that does not exist in the mesh
228-
229-
size_t index = idx->second;
234+
size_t index;
235+
if ( m2dm )
236+
index = m2dm->vertexIndex( i );
237+
else
238+
index = i;
230239

231240
if ( isVector )
232241
{
233242
if ( tsItems.size() >= 2 ) // BASEMENT files with vectors have 3 columns
234243
{
235-
dataset->values[index].x = toDouble( tsItems[0] );
236-
dataset->values[index].y = toDouble( tsItems[1] );
244+
values[2 * index] = toDouble( tsItems[0] );
245+
values[2 * index + 1] = toDouble( tsItems[1] );
237246
}
238247
else
239248
{
240249
debug( "invalid timestep line" );
241-
dataset->values[index].noData = true;
242250
}
243251
}
244252
else
245253
{
246254
if ( tsItems.size() >= 1 )
247-
dataset->values[index].x = toDouble( tsItems[0] );
255+
values[index] = toDouble( tsItems[0] );
248256
else
249257
{
250258
debug( "invalid timestep line" );
251-
dataset->values[index].noData = true;
252259
}
253260
}
254261
}
255262

263+
dataset->setStatistics( MDAL::calculateStatistics( dataset ) );
256264
group->datasets.push_back( dataset );
257265
}
258266

@@ -265,13 +273,11 @@ void MDAL::LoaderAsciiDat::readFaceTimestep(
265273
{
266274
assert( group );
267275

268-
size_t faceCount = mesh->faces.size();
269-
270-
std::shared_ptr<MDAL::Dataset> dataset = std::make_shared< MDAL::Dataset >();
271-
dataset->time = t / 3600.;
272-
dataset->values.resize( faceCount );
273-
dataset->parent = group.get();
276+
size_t faceCount = mesh->facesCount();
274277

278+
std::shared_ptr<MDAL::MemoryDataset> dataset = std::make_shared< MDAL::MemoryDataset >( group.get() );
279+
dataset->setTime( t / 3600. );
280+
double *values = dataset->values();
275281
// TODO: hasStatus
276282
for ( size_t index = 0; index < faceCount; ++index )
277283
{
@@ -283,26 +289,25 @@ void MDAL::LoaderAsciiDat::readFaceTimestep(
283289
{
284290
if ( tsItems.size() >= 2 ) // BASEMENT files with vectors have 3 columns
285291
{
286-
dataset->values[index].x = toDouble( tsItems[0] );
287-
dataset->values[index].y = toDouble( tsItems[1] );
292+
values[2 * index] = toDouble( tsItems[0] );
293+
values[2 * index + 1] = toDouble( tsItems[1] );
288294
}
289295
else
290296
{
291297
debug( "invalid timestep line" );
292-
dataset->values[index].noData = true;
293298
}
294299
}
295300
else
296301
{
297302
if ( tsItems.size() >= 1 )
298-
dataset->values[index].x = toDouble( tsItems[0] );
303+
values[index] = toDouble( tsItems[0] );
299304
else
300305
{
301306
debug( "invalid timestep line" );
302-
dataset->values[index].noData = true;
303307
}
304308
}
305309
}
306310

311+
dataset->setStatistics( MDAL::calculateStatistics( dataset ) );
307312
group->datasets.push_back( dataset );
308313
}

‎external/mdal/frmts/mdal_binary_dat.cpp

Lines changed: 38 additions & 27 deletions
Original file line numberDiff line numberDiff line change
@@ -103,8 +103,8 @@ void MDAL::LoaderBinaryDat::load( MDAL::Mesh *mesh, MDAL_Status *status )
103103
if ( !in )
104104
EXIT_WITH_ERROR( MDAL_Status::Err_FileNotFound ); // Couldn't open the file
105105

106-
size_t vertexCount = mesh->vertices.size();
107-
size_t elemCount = mesh->faces.size();
106+
size_t vertexCount = mesh->verticesCount();
107+
size_t elemCount = mesh->facesCount();
108108

109109
int card = 0;
110110
int version;
@@ -125,15 +125,19 @@ void MDAL::LoaderBinaryDat::load( MDAL::Mesh *mesh, MDAL_Status *status )
125125
if ( version != CT_VERSION ) // Version should be 3000
126126
EXIT_WITH_ERROR( MDAL_Status::Err_UnknownFormat );
127127

128-
std::shared_ptr<DatasetGroup> group = std::make_shared< DatasetGroup >(); // DAT datasets
129-
group->uri = mDatFile;
130-
group->isOnVertices = true;
128+
std::shared_ptr<DatasetGroup> group = std::make_shared< DatasetGroup >(
129+
mesh,
130+
mDatFile
131+
); // DAT datasets
132+
group->setIsOnVertices( true );
131133

132134
// in TUFLOW results there could be also a special timestep (99999) with maximums
133135
// we will put it into a separate dataset
134-
std::shared_ptr<DatasetGroup> groupMax = std::make_shared< DatasetGroup >();
135-
groupMax->uri = mDatFile;
136-
groupMax->isOnVertices = true;
136+
std::shared_ptr<DatasetGroup> groupMax = std::make_shared< DatasetGroup >(
137+
mesh,
138+
mDatFile
139+
);
140+
groupMax->setIsOnVertices( true );
137141

138142
while ( card != CT_ENDDS )
139143
{
@@ -166,13 +170,13 @@ void MDAL::LoaderBinaryDat::load( MDAL::Mesh *mesh, MDAL_Status *status )
166170
break;
167171

168172
case CT_BEGSCL:
169-
group->isScalar = true;
170-
groupMax->isScalar = true;
173+
group->setIsScalar( true );
174+
groupMax->setIsScalar( true );
171175
break;
172176

173177
case CT_BEGVEC:
174-
group->isScalar = false;
175-
groupMax->isScalar = false;
178+
group->setIsScalar( false );
179+
groupMax->setIsScalar( false );
176180
break;
177181

178182
case CT_VECTYPE:
@@ -231,8 +235,14 @@ void MDAL::LoaderBinaryDat::load( MDAL::Mesh *mesh, MDAL_Status *status )
231235

232236
if ( !group || group->datasets.size() == 0 )
233237
EXIT_WITH_ERROR( MDAL_Status::Err_UnknownFormat );
234-
238+
group->setStatistics( MDAL::calculateStatistics( group ) );
235239
mesh->datasetGroups.push_back( group );
240+
241+
if ( groupMax && groupMax->datasets.size() > 0 )
242+
{
243+
groupMax->setStatistics( MDAL::calculateStatistics( groupMax ) );
244+
mesh->datasetGroups.push_back( groupMax );
245+
}
236246
}
237247

238248
bool MDAL::LoaderBinaryDat::readVertexTimestep( const MDAL::Mesh *mesh,
@@ -243,17 +253,15 @@ bool MDAL::LoaderBinaryDat::readVertexTimestep( const MDAL::Mesh *mesh,
243253
int sflg,
244254
std::ifstream &in )
245255
{
246-
assert( group && groupMax && ( group->isScalar == groupMax->isScalar ) );
247-
bool isScalar = group->isScalar;
256+
assert( group && groupMax && ( group->isScalar() == groupMax->isScalar() ) );
257+
bool isScalar = group->isScalar();
248258

249-
size_t vertexCount = mesh->vertices.size();
250-
size_t faceCount = mesh->faces.size();
259+
size_t vertexCount = mesh->verticesCount();
260+
size_t faceCount = mesh->facesCount();
251261

252-
std::shared_ptr<MDAL::Dataset> dataset = std::make_shared< MDAL::Dataset >();
253-
dataset->values.resize( vertexCount );
254-
dataset->active.resize( faceCount );
255-
dataset->parent = group.get();
262+
std::shared_ptr<MDAL::MemoryDataset> dataset = std::make_shared< MDAL::MemoryDataset >( group.get() );
256263

264+
int *activeFlags = dataset->active();
257265
bool active = true;
258266
for ( size_t i = 0; i < faceCount; ++i )
259267
{
@@ -263,9 +271,10 @@ bool MDAL::LoaderBinaryDat::readVertexTimestep( const MDAL::Mesh *mesh,
263271
return true; //error
264272

265273
}
266-
dataset->active[i] = active;
274+
activeFlags[i] = active;
267275
}
268276

277+
double *values = dataset->values();
269278
for ( size_t i = 0; i < vertexCount; ++i )
270279
{
271280
if ( !isScalar )
@@ -277,8 +286,8 @@ bool MDAL::LoaderBinaryDat::readVertexTimestep( const MDAL::Mesh *mesh,
277286
if ( read( in, reinterpret_cast< char * >( &y ), 4 ) )
278287
return true; //error
279288

280-
dataset->values[i].x = static_cast< double >( x );
281-
dataset->values[i].y = static_cast< double >( y );
289+
values[2 * i] = static_cast< double >( x );
290+
values[2 * i + 1] = static_cast< double >( y );
282291
}
283292
else
284293
{
@@ -287,18 +296,20 @@ bool MDAL::LoaderBinaryDat::readVertexTimestep( const MDAL::Mesh *mesh,
287296
if ( read( in, reinterpret_cast< char * >( &scalar ), 4 ) )
288297
return true; //error
289298

290-
dataset->values[i].x = static_cast< double >( scalar );
299+
values[i] = static_cast< double >( scalar );
291300
}
292301
}
293302

294303
if ( MDAL::equals( time, 99999.0 ) ) // Special TUFLOW dataset with maximus
295304
{
296-
dataset->time = time;
305+
dataset->setTime( time );
306+
dataset->setStatistics( MDAL::calculateStatistics( dataset ) );
297307
groupMax->datasets.push_back( dataset );
298308
}
299309
else
300310
{
301-
dataset->time = time; // TODO read TIMEUNITS
311+
dataset->setTime( time ); // TODO read TIMEUNITS
312+
dataset->setStatistics( MDAL::calculateStatistics( dataset ) );
302313
group->datasets.push_back( dataset );
303314
}
304315
return false; //OK

‎external/mdal/frmts/mdal_cf.cpp

Lines changed: 39 additions & 36 deletions
Original file line numberDiff line numberDiff line change
@@ -132,49 +132,36 @@ MDAL::cfdataset_info_map MDAL::LoaderCF::parseDatasetGroupInfo()
132132
return dsinfo_map;
133133
}
134134

135-
static void populate_vals( bool is_vector, std::vector<MDAL::Value> &vals, size_t i,
135+
static void populate_vals( bool is_vector, double *vals, size_t i,
136136
const std::vector<double> &vals_x, const std::vector<double> &vals_y,
137137
size_t idx, double fill_val_x, double fill_val_y )
138138
{
139-
140-
vals[i].x = MDAL::safeValue( vals_x[idx], fill_val_x );
141139
if ( is_vector )
142140
{
143-
vals[i].y = MDAL::safeValue( vals_y[idx], fill_val_y );
141+
vals[2 * i] = MDAL::safeValue( vals_x[idx], fill_val_x );
142+
vals[2 * i + 1] = MDAL::safeValue( vals_y[idx], fill_val_y );
144143
}
145-
}
146-
147-
static void populate_nodata( std::vector<MDAL::Value> &vals, size_t from_i, size_t to_i )
148-
{
149-
for ( size_t i = from_i; i < to_i; ++i )
144+
else
150145
{
151-
vals[i].noData = true;
152-
vals[i].x = std::numeric_limits<double>::quiet_NaN();
153-
vals[i].y = std::numeric_limits<double>::quiet_NaN();
146+
vals[i] = MDAL::safeValue( vals_x[idx], fill_val_x );
154147
}
155148
}
156149

157-
158-
std::shared_ptr<MDAL::Dataset> MDAL::LoaderCF::createFace2DDataset( size_t ts, const MDAL::CFDatasetGroupInfo &dsi,
150+
std::shared_ptr<MDAL::Dataset> MDAL::LoaderCF::createFace2DDataset( std::shared_ptr<DatasetGroup> group, size_t ts, const MDAL::CFDatasetGroupInfo &dsi,
159151
const std::vector<double> &vals_x, const std::vector<double> &vals_y,
160152
double fill_val_x, double fill_val_y )
161153
{
162154
assert( dsi.outputType == CFDimensions::Face2D );
163155
size_t nFaces2D = mDimensions.size( CFDimensions::Face2D );
164156
size_t nLine1D = mDimensions.size( CFDimensions::Line1D );
165157

166-
std::shared_ptr<MDAL::Dataset> dataset = std::make_shared<MDAL::Dataset>();
167-
dataset->values.resize( mDimensions.faceCount() );
168-
169-
populate_nodata( dataset->values,
170-
0,
171-
nLine1D );
158+
std::shared_ptr<MDAL::MemoryDataset> dataset = std::make_shared<MDAL::MemoryDataset>( group.get() );
172159

173160
for ( size_t i = 0; i < nFaces2D; ++i )
174161
{
175162
size_t idx = ts * nFaces2D + i;
176163
populate_vals( dsi.is_vector,
177-
dataset->values,
164+
dataset->values(),
178165
nLine1D + i,
179166
vals_x,
180167
vals_y,
@@ -194,10 +181,13 @@ void MDAL::LoaderCF::addDatasetGroups( MDAL::Mesh *mesh, const std::vector<doubl
194181
{
195182
const CFDatasetGroupInfo dsi = it.second;
196183
// Create a dataset group
197-
std::shared_ptr<MDAL::DatasetGroup> group = std::make_shared<MDAL::DatasetGroup>();
198-
group->uri = mFileName;
199-
group->setName( dsi.name );
200-
group->isScalar = !dsi.is_vector;
184+
std::shared_ptr<MDAL::DatasetGroup> group = std::make_shared<MDAL::DatasetGroup>(
185+
mesh,
186+
mFileName,
187+
dsi.name
188+
);
189+
group->setIsScalar( !dsi.is_vector );
190+
group->setIsOnVertices( false );
201191

202192
// read X data
203193
double fill_val_x = mNcFile.getFillValue( dsi.ncid_x );
@@ -221,18 +211,20 @@ void MDAL::LoaderCF::addDatasetGroups( MDAL::Mesh *mesh, const std::vector<doubl
221211

222212
if ( dsi.outputType == CFDimensions::Face2D )
223213
{
224-
group->isOnVertices = false;
225-
dataset = createFace2DDataset( ts, dsi, vals_x, vals_y, fill_val_x, fill_val_y );
214+
dataset = createFace2DDataset( group, ts, dsi, vals_x, vals_y, fill_val_x, fill_val_y );
226215
}
227216

228-
dataset->parent = group.get();
229-
dataset->time = time;
217+
dataset->setTime( time );
218+
dataset->setStatistics( MDAL::calculateStatistics( dataset ) );
230219
group->datasets.push_back( dataset );
231220
}
232221

233222
// Add to mesh
234223
if ( !group->datasets.empty() )
224+
{
225+
group->setStatistics( MDAL::calculateStatistics( group ) );
235226
mesh->datasetGroups.push_back( group );
227+
}
236228
}
237229
}
238230

@@ -299,9 +291,6 @@ std::unique_ptr< MDAL::Mesh > MDAL::LoaderCF::load( MDAL_Status *status )
299291
{
300292
if ( status ) *status = MDAL_Status::None;
301293

302-
std::unique_ptr< MDAL::Mesh > mesh( new MDAL::Mesh );
303-
mesh->uri = mFileName;
304-
305294
//Dimensions dims;
306295
std::vector<double> times;
307296

@@ -314,7 +303,21 @@ std::unique_ptr< MDAL::Mesh > MDAL::LoaderCF::load( MDAL_Status *status )
314303
mDimensions = populateDimensions();
315304

316305
// Create mMesh
317-
populateFacesAndVertices( mesh.get() );
306+
Faces faces;
307+
Vertices vertices;
308+
populateFacesAndVertices( vertices, faces );
309+
std::unique_ptr< MemoryMesh > mesh(
310+
new MemoryMesh(
311+
vertices.size(),
312+
faces.size(),
313+
mDimensions.MaxVerticesInFace,
314+
computeExtent( vertices ),
315+
mFileName
316+
)
317+
);
318+
mesh->faces = faces;
319+
mesh->vertices = vertices;
320+
318321
addBedElevation( mesh.get() );
319322
setProjection( mesh.get() );
320323

@@ -326,14 +329,14 @@ std::unique_ptr< MDAL::Mesh > MDAL::LoaderCF::load( MDAL_Status *status )
326329

327330
// Create datasets
328331
addDatasetGroups( mesh.get(), times, dsinfo_map );
332+
333+
return mesh;
329334
}
330335
catch ( MDAL_Status error )
331336
{
332337
if ( status ) *status = ( error );
333-
if ( mesh ) mesh.reset();
338+
return std::unique_ptr<Mesh>();
334339
}
335-
336-
return mesh;
337340
}
338341

339342
//////////////////////////////////////////////////////////////////////////////////////

‎external/mdal/frmts/mdal_cf.hpp

Lines changed: 8 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -76,7 +76,7 @@ namespace MDAL
7676

7777
protected:
7878
virtual CFDimensions populateDimensions() = 0;
79-
virtual void populateFacesAndVertices( MDAL::Mesh *mesh ) = 0;
79+
virtual void populateFacesAndVertices( Vertices &vertices, Faces &faces ) = 0;
8080
virtual void addBedElevation( MDAL::Mesh *mesh ) = 0;
8181
virtual std::string getCoordinateSystemVariableName() = 0;
8282
virtual std::set<std::string> ignoreNetCDFVariables() = 0;
@@ -87,10 +87,13 @@ namespace MDAL
8787
void setProjection( MDAL::Mesh *m );
8888
cfdataset_info_map parseDatasetGroupInfo();
8989
void parseTime( std::vector<double> &times );
90-
std::shared_ptr<MDAL::Dataset> createFace2DDataset( size_t ts, const MDAL::CFDatasetGroupInfo &dsi,
91-
const std::vector<double> &vals_x,
92-
const std::vector<double> &vals_y,
93-
double fill_val_x, double fill_val_y );
90+
std::shared_ptr<MDAL::Dataset> createFace2DDataset(
91+
std::shared_ptr<MDAL::DatasetGroup> group,
92+
size_t ts,
93+
const MDAL::CFDatasetGroupInfo &dsi,
94+
const std::vector<double> &vals_x,
95+
const std::vector<double> &vals_y,
96+
double fill_val_x, double fill_val_y );
9497

9598
void addDatasetGroups( Mesh *mesh,
9699
const std::vector<double> &times,

‎external/mdal/frmts/mdal_gdal.cpp

Lines changed: 70 additions & 46 deletions
Original file line numberDiff line numberDiff line change
@@ -258,14 +258,16 @@ void MDAL::LoaderGdal::parseRasterBands( const MDAL::GdalDataset *cfGDALDataset
258258
}
259259
}
260260

261-
void MDAL::LoaderGdal::addDataToOutput( GDALRasterBandH raster_band, std::shared_ptr<Dataset> tos, bool is_vector, bool is_x )
261+
void MDAL::LoaderGdal::addDataToOutput( GDALRasterBandH raster_band, std::shared_ptr<MemoryDataset> tos, bool is_vector, bool is_x )
262262
{
263263
assert( raster_band );
264264

265265
double nodata = GDALGetRasterNoDataValue( raster_band, nullptr );
266266
unsigned int mXSize = meshGDALDataset()->mXSize;
267267
unsigned int mYSize = meshGDALDataset()->mYSize;
268268

269+
double *values = tos->values();
270+
269271
for ( unsigned int y = 0; y < mYSize; ++y )
270272
{
271273
// buffering per-line
@@ -292,53 +294,65 @@ void MDAL::LoaderGdal::addDataToOutput( GDALRasterBandH raster_band, std::shared
292294
{
293295
unsigned int idx = x + mXSize * y;
294296
double val = mPafScanline[x];
295-
bool noData = false;
296-
if ( MDAL::equals( val, nodata ) )
297-
{
298-
// store all nodata value as this hardcoded number
299-
val = MDAL_NODATA;
300-
noData = true;
301-
}
302-
303-
if ( is_vector )
297+
if ( !MDAL::equals( val, nodata ) )
304298
{
305-
if ( is_x )
299+
// values is prepolulated with NODATA values, so store only legal values
300+
if ( is_vector )
306301
{
307-
tos->values[idx].x = val;
308-
tos->values[idx].noData = noData;
302+
if ( is_x )
303+
{
304+
values[2 * idx] = val;
305+
}
306+
else
307+
{
308+
values[2 * idx + 1] = val;
309+
}
309310
}
310311
else
311312
{
312-
tos->values[idx].y = val;
313-
tos->values[idx].noData = noData;
313+
values[idx] = val;
314314
}
315315
}
316-
else
317-
{
318-
tos->values[idx].x = val;
319-
tos->values[idx].noData = noData;
320-
}
321316
}
322317
}
323318
}
324319

325-
void MDAL::LoaderGdal::activateFaces( std::shared_ptr<Dataset> tos )
320+
void MDAL::LoaderGdal::activateFaces( std::shared_ptr<MemoryDataset> tos )
326321
{
322+
// only for data on vertices
323+
if ( !tos->group()->isOnVertices() )
324+
return;
325+
326+
bool isScalar = tos->group()->isScalar();
327+
327328
// Activate only Faces that do all Vertex's outputs with some data
329+
int *active = tos->active();
330+
const double *values = tos->constValues();
331+
328332
for ( unsigned int idx = 0; idx < meshGDALDataset()->mNVolumes; ++idx )
329333
{
330334
Face elem = mMesh->faces.at( idx );
331-
332-
if ( tos->values[elem[0]].noData ||
333-
tos->values[elem[1]].noData ||
334-
tos->values[elem[2]].noData ||
335-
tos->values[elem[3]].noData )
335+
for ( size_t i = 0; i < 4; ++i )
336336
{
337-
tos->active[idx] = 0; //NOT ACTIVE
338-
}
339-
else
340-
{
341-
tos->active[idx] = 1; //ACTIVE
337+
if ( isScalar )
338+
{
339+
double val = values[elem[i]];
340+
if ( std::isnan( val ) )
341+
{
342+
active[idx] = 0; //NOT ACTIVE
343+
break;
344+
}
345+
}
346+
else
347+
{
348+
double x = values[elem[2 * i]];
349+
double y = values[elem[2 * i + 1]];
350+
if ( std::isnan( x ) || std::isnan( y ) )
351+
{
352+
active[idx] = 0; //NOT ACTIVE
353+
break;
354+
}
355+
}
342356
}
343357
}
344358
}
@@ -348,32 +362,35 @@ void MDAL::LoaderGdal::addDatasetGroups()
348362
// Add dataset to mMesh
349363
for ( data_hash::const_iterator band = mBands.begin(); band != mBands.end(); band++ )
350364
{
351-
std::shared_ptr<DatasetGroup> group = std::make_shared< DatasetGroup >();
352-
group->uri = mFileName;
353-
group->setName( band->first );
354-
group->isOnVertices = true;
365+
if ( band->second.empty() )
366+
continue;
367+
368+
std::shared_ptr<DatasetGroup> group = std::make_shared< DatasetGroup >(
369+
mMesh.get(),
370+
mFileName,
371+
band->first
372+
);
373+
group->setIsOnVertices( true );
374+
bool is_vector = ( band->second.begin()->second.size() > 1 );
375+
group->setIsScalar( !is_vector );
355376

356377
for ( timestep_map::const_iterator time_step = band->second.begin(); time_step != band->second.end(); time_step++ )
357378
{
358379
std::vector<GDALRasterBandH> raster_bands = time_step->second;
359-
bool is_vector = ( raster_bands.size() > 1 );
360-
361-
std::shared_ptr<MDAL::Dataset> dataset = std::make_shared< MDAL::Dataset >();
362-
group->isScalar = !is_vector;
363-
364-
dataset->time = time_step->first;
365-
dataset->values.resize( meshGDALDataset()->mNPoints );
366-
dataset->active.resize( meshGDALDataset()->mNVolumes );
367-
dataset->parent = group.get();
380+
std::shared_ptr<MDAL::MemoryDataset> dataset = std::make_shared< MDAL::MemoryDataset >( group.get() );
368381

382+
dataset->setTime( time_step->first );
369383
for ( std::vector<GDALRasterBandH>::size_type i = 0; i < raster_bands.size(); ++i )
370384
{
371385
addDataToOutput( raster_bands[i], dataset, is_vector, i == 0 );
372386
}
373387
activateFaces( dataset );
374-
388+
dataset->setStatistics( MDAL::calculateStatistics( dataset ) );
375389
group->datasets.push_back( dataset );
376390
}
391+
392+
// TODO use GDALComputeRasterMinMax
393+
group->setStatistics( MDAL::calculateStatistics( group ) );
377394
mMesh->datasetGroups.push_back( group );
378395
}
379396
}
@@ -386,7 +403,14 @@ void MDAL::LoaderGdal::createMesh()
386403
Faces faces( meshGDALDataset()->mNVolumes );
387404
initFaces( vertices, faces, is_longitude_shifted );
388405

389-
mMesh.reset( new Mesh() );
406+
mMesh.reset( new MemoryMesh(
407+
vertices.size(),
408+
faces.size(),
409+
4, //maximum quads
410+
computeExtent( vertices ),
411+
mFileName
412+
)
413+
);
390414
mMesh->vertices = vertices;
391415
mMesh->faces = faces;
392416
bool proj_added = addSrcProj();

‎external/mdal/frmts/mdal_gdal.hpp

Lines changed: 3 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -82,17 +82,17 @@ namespace MDAL
8282
bool meshes_equals( const GdalDataset *ds1, const GdalDataset *ds2 ) const;
8383

8484
metadata_hash parseMetadata( GDALMajorObjectH gdalBand, const char *pszDomain = nullptr );
85-
void addDataToOutput( GDALRasterBandH raster_band, std::shared_ptr<Dataset> tos, bool is_vector, bool is_x );
85+
void addDataToOutput( GDALRasterBandH raster_band, std::shared_ptr<MemoryDataset> tos, bool is_vector, bool is_x );
8686
bool addSrcProj();
87-
void activateFaces( std::shared_ptr<Dataset> tos );
87+
void activateFaces( std::shared_ptr<MemoryDataset> tos );
8888
void addDatasetGroups();
8989
void createMesh();
9090
void parseRasterBands( const GdalDataset *cfGDALDataset );
9191

9292
const std::string mFileName;
9393
const std::string mDriverName; /* GDAL driver name */
9494
double *mPafScanline; /* temporary buffer for reading one raster line */
95-
std::unique_ptr< Mesh > mMesh;
95+
std::unique_ptr< MemoryMesh > mMesh;
9696
gdal_datasets_vector gdal_datasets;
9797
data_hash mBands; /* raster bands GDAL handle */
9898
};

‎external/mdal/frmts/mdal_hdf5.cpp

Lines changed: 241 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,241 @@
1+
/*
2+
MDAL - Mesh Data Abstraction Library (MIT License)
3+
Copyright (C) 2018 Lutra Consulting Limited
4+
*/
5+
6+
#include "mdal_hdf5.hpp"
7+
8+
9+
HdfFile::HdfFile( const std::string &path )
10+
: d( std::make_shared< Handle >( H5Fopen( path.c_str(), H5F_ACC_RDONLY, H5P_DEFAULT ) ) )
11+
{
12+
}
13+
14+
bool HdfFile::isValid() const { return d->id >= 0; }
15+
16+
hid_t HdfFile::id() const { return d->id; }
17+
18+
HdfGroup::HdfGroup( hid_t file, const std::string &path )
19+
: d( std::make_shared< Handle >( H5Gopen( file, path.c_str() ) ) )
20+
{}
21+
22+
bool HdfGroup::isValid() const { return d->id >= 0; }
23+
24+
hid_t HdfGroup::id() const { return d->id; }
25+
26+
hid_t HdfGroup::file_id() const { return H5Iget_file_id( d->id ); }
27+
28+
std::string HdfGroup::name() const
29+
{
30+
char name[HDF_MAX_NAME];
31+
H5Iget_name( d->id, name, HDF_MAX_NAME );
32+
return std::string( name );
33+
}
34+
35+
std::vector<std::string> HdfGroup::groups() const { return objects( H5G_GROUP ); }
36+
37+
std::vector<std::string> HdfGroup::datasets() const { return objects( H5G_DATASET ); }
38+
39+
std::vector<std::string> HdfGroup::objects() const { return objects( H5G_UNKNOWN ); }
40+
41+
std::string HdfGroup::childPath( const std::string &childName ) const { return name() + "/" + childName; }
42+
43+
std::vector<std::string> HdfGroup::objects( H5G_obj_t type ) const
44+
{
45+
std::vector<std::string> lst;
46+
47+
hsize_t nobj;
48+
H5Gget_num_objs( d->id, &nobj );
49+
for ( hsize_t i = 0; i < nobj; ++i )
50+
{
51+
if ( type == H5G_UNKNOWN || H5Gget_objtype_by_idx( d->id, i ) == type )
52+
{
53+
char name[HDF_MAX_NAME];
54+
H5Gget_objname_by_idx( d->id, i, name, ( size_t )HDF_MAX_NAME );
55+
lst.push_back( std::string( name ) );
56+
}
57+
}
58+
return lst;
59+
}
60+
61+
HdfAttribute::HdfAttribute( hid_t obj_id, const std::string &attr_name )
62+
: d( std::make_shared< Handle >( H5Aopen( obj_id, attr_name.c_str(), H5P_DEFAULT ) ) )
63+
{}
64+
65+
bool HdfAttribute::isValid() const { return d->id >= 0; }
66+
67+
hid_t HdfAttribute::id() const { return d->id; }
68+
69+
std::string HdfAttribute::readString() const
70+
{
71+
char name[HDF_MAX_NAME];
72+
hid_t datatype = H5Tcopy( H5T_C_S1 );
73+
H5Tset_size( datatype, HDF_MAX_NAME );
74+
herr_t status = H5Aread( d->id, datatype, name );
75+
if ( status < 0 )
76+
{
77+
//MDAL::debug("Failed to read data!");
78+
return std::string();
79+
}
80+
H5Tclose( datatype );
81+
return std::string( name );
82+
}
83+
84+
HdfDataset::HdfDataset( hid_t file, const std::string &path )
85+
: d( std::make_shared< Handle >( H5Dopen2( file, path.c_str(), H5P_DEFAULT ) ) )
86+
{}
87+
88+
bool HdfDataset::isValid() const { return d->id >= 0; }
89+
90+
hid_t HdfDataset::id() const { return d->id; }
91+
92+
std::vector<hsize_t> HdfDataset::dims() const
93+
{
94+
hid_t sid = H5Dget_space( d->id );
95+
std::vector<hsize_t> d( H5Sget_simple_extent_ndims( sid ) );
96+
H5Sget_simple_extent_dims( sid, d.data(), NULL );
97+
H5Sclose( sid );
98+
return d;
99+
}
100+
101+
hsize_t HdfDataset::elementCount() const
102+
{
103+
hsize_t count = 1;
104+
for ( hsize_t dsize : dims() )
105+
count *= dsize;
106+
return count;
107+
}
108+
109+
H5T_class_t HdfDataset::type() const
110+
{
111+
hid_t tid = H5Dget_type( d->id );
112+
H5T_class_t t_class = H5Tget_class( tid );
113+
H5Tclose( tid );
114+
return t_class;
115+
}
116+
117+
std::vector<uchar> HdfDataset::readArrayUint8( const std::vector<hsize_t> offsets, const std::vector<hsize_t> counts ) const { return readArray<uchar>( H5T_NATIVE_UINT8, offsets, counts ); }
118+
119+
std::vector<float> HdfDataset::readArray( const std::vector<hsize_t> offsets, const std::vector<hsize_t> counts ) const { return readArray<float>( H5T_NATIVE_FLOAT, offsets, counts ); }
120+
121+
std::vector<double> HdfDataset::readArrayDouble( const std::vector<hsize_t> offsets, const std::vector<hsize_t> counts ) const { return readArray<double>( H5T_NATIVE_DOUBLE, offsets, counts ); }
122+
123+
std::vector<int> HdfDataset::readArrayInt( const std::vector<hsize_t> offsets, const std::vector<hsize_t> counts ) const { return readArray<int>( H5T_NATIVE_INT, offsets, counts ); }
124+
125+
std::vector<uchar> HdfDataset::readArrayUint8() const { return readArray<uchar>( H5T_NATIVE_UINT8 ); }
126+
127+
std::vector<float> HdfDataset::readArray() const { return readArray<float>( H5T_NATIVE_FLOAT ); }
128+
129+
std::vector<double> HdfDataset::readArrayDouble() const { return readArray<double>( H5T_NATIVE_DOUBLE ); }
130+
131+
std::vector<int> HdfDataset::readArrayInt() const { return readArray<int>( H5T_NATIVE_INT ); }
132+
133+
std::vector<std::string> HdfDataset::readArrayString() const
134+
{
135+
std::vector<std::string> ret;
136+
137+
hid_t datatype = H5Tcopy( H5T_C_S1 );
138+
H5Tset_size( datatype, HDF_MAX_NAME );
139+
140+
std::vector<HdfString> arr = readArray<HdfString>( datatype );
141+
142+
H5Tclose( datatype );
143+
144+
for ( const HdfString &str : arr )
145+
{
146+
std::string dat = std::string( str.data );
147+
ret.push_back( MDAL::trim( dat ) );
148+
}
149+
150+
return ret;
151+
}
152+
153+
float HdfDataset::readFloat() const
154+
{
155+
if ( elementCount() != 1 )
156+
{
157+
MDAL::debug( "Not scalar!" );
158+
return 0;
159+
}
160+
161+
float value;
162+
herr_t status = H5Dread( d->id, H5T_NATIVE_FLOAT, H5S_ALL, H5S_ALL, H5P_DEFAULT, &value );
163+
if ( status < 0 )
164+
{
165+
MDAL::debug( "Failed to read data!" );
166+
return 0;
167+
}
168+
return value;
169+
}
170+
171+
std::string HdfDataset::readString() const
172+
{
173+
if ( elementCount() != 1 )
174+
{
175+
MDAL::debug( "Not scalar!" );
176+
return std::string();
177+
}
178+
179+
char name[HDF_MAX_NAME];
180+
hid_t datatype = H5Tcopy( H5T_C_S1 );
181+
H5Tset_size( datatype, HDF_MAX_NAME );
182+
herr_t status = H5Dread( d->id, datatype, H5S_ALL, H5S_ALL, H5P_DEFAULT, name );
183+
if ( status < 0 )
184+
{
185+
MDAL::debug( "Failed to read data!" );
186+
return std::string();
187+
}
188+
H5Tclose( datatype );
189+
return std::string( name );
190+
}
191+
192+
193+
HdfDataspace::HdfDataspace( const std::vector<hsize_t> &dims )
194+
: d( std::make_shared< Handle >( H5Screate_simple(
195+
dims.size(),
196+
dims.data(),
197+
dims.data()
198+
)
199+
)
200+
)
201+
{
202+
}
203+
204+
HdfDataspace::HdfDataspace( hid_t dataset )
205+
: d( std::make_shared< Handle >( H5Dget_space( dataset ) ) )
206+
{
207+
}
208+
209+
void HdfDataspace::selectHyperslab( hsize_t start, hsize_t count )
210+
{
211+
// this function works only for 1D arrays
212+
assert( H5Sget_simple_extent_ndims( d->id ) == 1 );
213+
214+
herr_t status = H5Sselect_hyperslab( d->id, H5S_SELECT_SET, &start, NULL, &count, NULL );
215+
if ( status < 0 )
216+
{
217+
MDAL::debug( "Failed to select 1D hyperslab!" );
218+
}
219+
}
220+
221+
void HdfDataspace::selectHyperslab( const std::vector<hsize_t> offsets,
222+
const std::vector<hsize_t> counts )
223+
{
224+
assert( H5Sget_simple_extent_ndims( d->id ) == static_cast<int>( offsets.size() ) );
225+
assert( offsets.size() == counts.size() );
226+
227+
herr_t status = H5Sselect_hyperslab( d->id,
228+
H5S_SELECT_SET,
229+
offsets.data(),
230+
NULL,
231+
counts.data(),
232+
NULL );
233+
if ( status < 0 )
234+
{
235+
MDAL::debug( "Failed to select 1D hyperslab!" );
236+
}
237+
}
238+
239+
bool HdfDataspace::isValid() const { return d->id >= 0; }
240+
241+
hid_t HdfDataspace::id() const { return d->id; }

‎external/mdal/frmts/mdal_hdf5.hpp

Lines changed: 91 additions & 141 deletions
Original file line numberDiff line numberDiff line change
@@ -6,6 +6,7 @@
66
#ifndef MDAL_HDF5_HPP
77
#define MDAL_HDF5_HPP
88

9+
910
/** A simple C++ wrapper around HDF5 library API */
1011

1112
// for compatibility (older hdf5 version in Travis)
@@ -18,6 +19,8 @@ typedef unsigned char uchar;
1819
#include <memory>
1920
#include <vector>
2021
#include <string>
22+
#include <numeric>
23+
2124
#include <assert.h>
2225
#include "stdlib.h"
2326

@@ -31,7 +34,8 @@ template <int TYPE> inline void hdfClose( hid_t id ) { MDAL_UNUSED( id ); assert
3134
template <> inline void hdfClose<H5I_FILE>( hid_t id ) { H5Fclose( id ); }
3235
template <> inline void hdfClose<H5I_GROUP>( hid_t id ) { H5Gclose( id ); }
3336
template <> inline void hdfClose<H5I_DATASET>( hid_t id ) { H5Dclose( id ); }
34-
template <> inline void hdfClose<H5I_ATTR>( hid_t id ) { H5Dclose( id ); }
37+
template <> inline void hdfClose<H5I_ATTR>( hid_t id ) { H5Aclose( id ); }
38+
template <> inline void hdfClose<H5I_DATASPACE>( hid_t id ) { H5Sclose( id ); }
3539

3640
template <int TYPE>
3741
class HdfH
@@ -47,25 +51,24 @@ class HdfH
4751
class HdfGroup;
4852
class HdfDataset;
4953
class HdfAttribute;
54+
class HdfDataspace;
5055

5156
class HdfFile
5257
{
5358
public:
5459
typedef HdfH<H5I_FILE> Handle;
5560

56-
HdfFile( const std::string &path )
57-
: d( std::make_shared< Handle >( H5Fopen( path.c_str(), H5F_ACC_RDONLY, H5P_DEFAULT ) ) )
58-
{
59-
}
61+
HdfFile( const std::string &path );
6062

61-
bool isValid() const { return d->id >= 0; }
62-
hid_t id() const { return d->id; }
63+
bool isValid() const;
64+
hid_t id() const;
6365

6466
inline std::vector<std::string> groups() const;
6567

6668
inline HdfGroup group( const std::string &path ) const;
6769
inline HdfDataset dataset( const std::string &path ) const;
6870
inline HdfAttribute attribute( const std::string &attr_name ) const;
71+
inline bool pathExists( const std::string &path ) const;
6972

7073
protected:
7174
std::shared_ptr<Handle> d;
@@ -76,49 +79,27 @@ class HdfGroup
7679
public:
7780
typedef HdfH<H5I_GROUP> Handle;
7881

79-
HdfGroup( hid_t file, const std::string &path )
80-
: d( std::make_shared< Handle >( H5Gopen( file, path.c_str() ) ) )
81-
{}
82+
HdfGroup( hid_t file, const std::string &path );
8283

83-
bool isValid() const { return d->id >= 0; }
84-
hid_t id() const { return d->id; }
85-
hid_t file_id() const { return H5Iget_file_id( d->id ); }
84+
bool isValid() const;
85+
hid_t id() const;
86+
hid_t file_id() const;
8687

87-
std::string name() const
88-
{
89-
char name[HDF_MAX_NAME];
90-
H5Iget_name( d->id, name, HDF_MAX_NAME );
91-
return std::string( name );
92-
}
88+
std::string name() const;
9389

94-
std::vector<std::string> groups() const { return objects( H5G_GROUP ); }
95-
std::vector<std::string> datasets() const { return objects( H5G_DATASET ); }
96-
std::vector<std::string> objects() const { return objects( H5G_UNKNOWN ); }
90+
std::vector<std::string> groups() const;
91+
std::vector<std::string> datasets() const;
92+
std::vector<std::string> objects() const;
9793

98-
std::string childPath( const std::string &childName ) const { return name() + "/" + childName; }
94+
std::string childPath( const std::string &childName ) const;
9995

10096
inline HdfGroup group( const std::string &groupName ) const;
10197
inline HdfDataset dataset( const std::string &dsName ) const;
10298
inline HdfAttribute attribute( const std::string &attr_name ) const;
99+
inline bool pathExists( const std::string &path ) const;
103100

104101
protected:
105-
std::vector<std::string> objects( H5G_obj_t type ) const
106-
{
107-
std::vector<std::string> lst;
108-
109-
hsize_t nobj;
110-
H5Gget_num_objs( d->id, &nobj );
111-
for ( hsize_t i = 0; i < nobj; ++i )
112-
{
113-
if ( type == H5G_UNKNOWN || H5Gget_objtype_by_idx( d->id, i ) == type )
114-
{
115-
char name[HDF_MAX_NAME];
116-
H5Gget_objname_by_idx( d->id, i, name, ( size_t )HDF_MAX_NAME );
117-
lst.push_back( std::string( name ) );
118-
}
119-
}
120-
return lst;
121-
}
102+
std::vector<std::string> objects( H5G_obj_t type ) const;
122103

123104
protected:
124105
std::shared_ptr<Handle> d;
@@ -130,95 +111,70 @@ class HdfAttribute
130111
public:
131112
typedef HdfH<H5I_ATTR> Handle;
132113

133-
HdfAttribute( hid_t obj_id, const std::string &attr_name )
134-
: d( std::make_shared< Handle >( H5Aopen( obj_id, attr_name.c_str(), H5P_DEFAULT ) ) )
135-
{}
114+
HdfAttribute( hid_t obj_id, const std::string &attr_name );
136115

137-
bool isValid() const { return d->id >= 0; }
138-
hid_t id() const { return d->id; }
116+
bool isValid() const;
117+
hid_t id() const;
139118

140-
std::string readString() const
141-
{
142-
char name[HDF_MAX_NAME];
143-
hid_t datatype = H5Tcopy( H5T_C_S1 );
144-
H5Tset_size( datatype, HDF_MAX_NAME );
145-
herr_t status = H5Aread( d->id, datatype, name );
146-
if ( status < 0 )
147-
{
148-
//MDAL::debug("Failed to read data!");
149-
return std::string();
150-
}
151-
H5Tclose( datatype );
152-
return std::string( name );
153-
}
119+
std::string readString() const;
154120
protected:
155121
std::shared_ptr<Handle> d;
156122
};
157123

158-
class HdfDataset
124+
class HdfDataspace
159125
{
160126
public:
161-
typedef HdfH<H5I_DATASET> Handle;
162-
163-
HdfDataset( hid_t file, const std::string &path )
164-
: d( std::make_shared< Handle >( H5Dopen2( file, path.c_str(), H5P_DEFAULT ) ) )
165-
{}
166-
167-
bool isValid() const { return d->id >= 0; }
168-
hid_t id() const { return d->id; }
169-
170-
std::vector<hsize_t> dims() const
171-
{
172-
hid_t sid = H5Dget_space( d->id );
173-
std::vector<hsize_t> d( H5Sget_simple_extent_ndims( sid ) );
174-
H5Sget_simple_extent_dims( sid, d.data(), NULL );
175-
H5Sclose( sid );
176-
return d;
177-
}
178-
179-
hsize_t elementCount() const
180-
{
181-
hsize_t count = 1;
182-
for ( hsize_t dsize : dims() )
183-
count *= dsize;
184-
return count;
185-
}
186-
187-
H5T_class_t type() const
188-
{
189-
hid_t tid = H5Dget_type( d->id );
190-
H5T_class_t t_class = H5Tget_class( tid );
191-
H5Tclose( tid );
192-
return t_class;
193-
}
127+
typedef HdfH<H5I_DATASPACE> Handle;
128+
//! memory dataspace for simple N-D array
129+
HdfDataspace( const std::vector<hsize_t> &dims );
130+
//! dataspace of the dataset
131+
HdfDataspace( hid_t dataset );
132+
//! select from 1D array
133+
void selectHyperslab( hsize_t start, hsize_t count );
134+
//! select from N-D array
135+
void selectHyperslab( const std::vector<hsize_t> offsets,
136+
const std::vector<hsize_t> counts );
137+
138+
bool isValid() const;
139+
hid_t id() const;
194140

195-
std::vector<uchar> readArrayUint8() const { return readArray<uchar>( H5T_NATIVE_UINT8 ); }
196-
197-
std::vector<float> readArray() const { return readArray<float>( H5T_NATIVE_FLOAT ); }
141+
protected:
142+
std::shared_ptr<Handle> d;
143+
};
198144

199-
std::vector<double> readArrayDouble() const { return readArray<double>( H5T_NATIVE_DOUBLE ); }
145+
class HdfDataset
146+
{
147+
public:
148+
typedef HdfH<H5I_DATASET> Handle;
200149

201-
std::vector<int> readArrayInt() const { return readArray<int>( H5T_NATIVE_INT ); }
150+
HdfDataset( hid_t file, const std::string &path );
202151

203-
std::vector<std::string> readArrayString() const
204-
{
205-
std::vector<std::string> ret;
152+
bool isValid() const;
153+
hid_t id() const;
206154

207-
hid_t datatype = H5Tcopy( H5T_C_S1 );
208-
H5Tset_size( datatype, HDF_MAX_NAME );
155+
std::vector<hsize_t> dims() const;
209156

210-
std::vector<HdfString> arr = readArray<HdfString>( datatype );
157+
hsize_t elementCount() const;
211158

212-
H5Tclose( datatype );
159+
H5T_class_t type() const;
213160

214-
for ( const HdfString &str : arr )
215-
{
216-
std::string dat = std::string( str.data );
217-
ret.push_back( MDAL::trim( dat ) );
218-
}
161+
//! Reads full array into vector
162+
//! Array can have any number of dimenstions
163+
//! and it is fully read into 1D vector
164+
std::vector<uchar> readArrayUint8() const;
165+
std::vector<float> readArray() const;
166+
std::vector<double> readArrayDouble() const;
167+
std::vector<int> readArrayInt() const;
168+
std::vector<std::string> readArrayString() const;
219169

220-
return ret;
221-
}
170+
//! Reads part of the N-D array into vector,
171+
//! for each dimension specified by offset and count
172+
//! size of offsets and counts must be same as rank (number of dims) of dataset
173+
//! the results array is 1D
174+
std::vector<uchar> readArrayUint8( const std::vector<hsize_t> offsets, const std::vector<hsize_t> counts ) const;
175+
std::vector<float> readArray( const std::vector<hsize_t> offsets, const std::vector<hsize_t> counts ) const;
176+
std::vector<double> readArrayDouble( const std::vector<hsize_t> offsets, const std::vector<hsize_t> counts ) const;
177+
std::vector<int> readArrayInt( const std::vector<hsize_t> offsets, const std::vector<hsize_t> counts ) const;
222178

223179

224180
template <typename T> std::vector<T> readArray( hid_t mem_type_id ) const
@@ -234,45 +190,35 @@ class HdfDataset
234190
return data;
235191
}
236192

237-
float readFloat() const
193+
template <typename T> std::vector<T> readArray( hid_t mem_type_id,
194+
const std::vector<hsize_t> offsets,
195+
const std::vector<hsize_t> counts ) const
238196
{
239-
if ( elementCount() != 1 )
240-
{
241-
MDAL::debug( "Not scalar!" );
242-
return 0;
243-
}
197+
HdfDataspace dataspace( d->id );
198+
dataspace.selectHyperslab( offsets, counts );
244199

245-
float value;
246-
herr_t status = H5Dread( d->id, H5T_NATIVE_FLOAT, H5S_ALL, H5S_ALL, H5P_DEFAULT, &value );
247-
if ( status < 0 )
248-
{
249-
MDAL::debug( "Failed to read data!" );
250-
return 0;
251-
}
252-
return value;
253-
}
200+
hsize_t totalItems = 1;
201+
for ( auto it = counts.begin(); it != counts.end(); ++it )
202+
totalItems *= *it;
254203

255-
std::string readString() const
256-
{
257-
if ( elementCount() != 1 )
258-
{
259-
MDAL::debug( "Not scalar!" );
260-
return std::string();
261-
}
204+
std::vector<hsize_t> dims = {totalItems};
205+
HdfDataspace memspace( dims );
206+
memspace.selectHyperslab( 0, totalItems );
262207

263-
char name[HDF_MAX_NAME];
264-
hid_t datatype = H5Tcopy( H5T_C_S1 );
265-
H5Tset_size( datatype, HDF_MAX_NAME );
266-
herr_t status = H5Dread( d->id, datatype, H5S_ALL, H5S_ALL, H5P_DEFAULT, name );
208+
std::vector<T> data( totalItems );
209+
herr_t status = H5Dread( d->id, mem_type_id, memspace.id(), dataspace.id(), H5P_DEFAULT, data.data() );
267210
if ( status < 0 )
268211
{
269212
MDAL::debug( "Failed to read data!" );
270-
return std::string();
213+
return std::vector<T>();
271214
}
272-
H5Tclose( datatype );
273-
return std::string( name );
215+
return data;
274216
}
275217

218+
float readFloat() const;
219+
220+
std::string readString() const;
221+
276222
protected:
277223
std::shared_ptr<Handle> d;
278224
};
@@ -291,4 +237,8 @@ inline HdfAttribute HdfFile::attribute( const std::string &attr_name ) const { r
291237

292238
inline HdfAttribute HdfGroup::attribute( const std::string &attr_name ) const { return HdfAttribute( d->id, attr_name ); }
293239

240+
inline bool HdfFile::pathExists( const std::string &path ) const { return H5Lexists( d->id, path.c_str(), H5P_DEFAULT ) > 0; }
241+
242+
inline bool HdfGroup::pathExists( const std::string &path ) const { return H5Lexists( d->id, path.c_str(), H5P_DEFAULT ) > 0; }
243+
294244
#endif // MDAL_HDF5_HPP

‎external/mdal/frmts/mdal_xmdf.cpp

Lines changed: 141 additions & 57 deletions
Original file line numberDiff line numberDiff line change
@@ -12,14 +12,86 @@
1212
#include <string>
1313
#include <vector>
1414
#include <memory>
15+
#include <algorithm>
1516

17+
MDAL::XmdfDataset::~XmdfDataset() = default;
18+
19+
MDAL::XmdfDataset::XmdfDataset( DatasetGroup *grp, const HdfDataset &valuesDs, const HdfDataset &activeDs, hsize_t timeIndex )
20+
: Dataset( grp )
21+
, mHdf5DatasetValues( valuesDs )
22+
, mHdf5DatasetActive( activeDs )
23+
, mTimeIndex( timeIndex )
24+
{
25+
}
26+
27+
const HdfDataset &MDAL::XmdfDataset::dsValues() const
28+
{
29+
return mHdf5DatasetValues;
30+
}
31+
32+
const HdfDataset &MDAL::XmdfDataset::dsActive() const
33+
{
34+
return mHdf5DatasetActive;
35+
}
36+
37+
hsize_t MDAL::XmdfDataset::timeIndex() const
38+
{
39+
return mTimeIndex;
40+
}
41+
42+
43+
size_t MDAL::XmdfDataset::scalarData( size_t indexStart, size_t count, double *buffer )
44+
{
45+
assert( group()->isScalar() ); //checked in C API interface
46+
std::vector<hsize_t> offsets = {timeIndex(), indexStart};
47+
std::vector<hsize_t> counts = {1, count};
48+
std::vector<float> values = dsValues().readArray( offsets, counts );
49+
const float *input = values.data();
50+
for ( size_t j = 0; j < count; ++j )
51+
{
52+
buffer[j] = double( input[j] );
53+
}
54+
return count;
55+
}
56+
57+
size_t MDAL::XmdfDataset::vectorData( size_t indexStart, size_t count, double *buffer )
58+
{
59+
assert( !group()->isScalar() ); //checked in C API interface
60+
std::vector<hsize_t> offsets = {timeIndex(), indexStart, 0};
61+
std::vector<hsize_t> counts = {1, count, 2};
62+
std::vector<float> values = dsValues().readArray( offsets, counts );
63+
const float *input = values.data();
64+
for ( size_t j = 0; j < count; ++j )
65+
{
66+
buffer[2 * j] = double( input[2 * j] );
67+
buffer[2 * j + 1] = double( input[2 * j + 1] );
68+
}
69+
70+
return count;
71+
}
72+
73+
size_t MDAL::XmdfDataset::activeData( size_t indexStart, size_t count, int *buffer )
74+
{
75+
std::vector<hsize_t> offsets = {timeIndex(), indexStart};
76+
std::vector<hsize_t> counts = {1, count};
77+
std::vector<uchar> active = dsActive().readArrayUint8( offsets, counts );
78+
const uchar *input = active.data();
79+
for ( size_t j = 0; j < count; ++j )
80+
{
81+
buffer[j] = bool( input[ j ] );
82+
}
83+
return count;
84+
}
85+
86+
///////////////////////////////////////////////////////////////////////////////////////
1687

1788
MDAL::LoaderXmdf::LoaderXmdf( const std::string &datFile )
1889
: mDatFile( datFile )
1990
{}
2091

2192
void MDAL::LoaderXmdf::load( MDAL::Mesh *mesh, MDAL_Status *status )
2293
{
94+
mMesh = mesh;
2395
if ( status ) *status = MDAL_Status::None;
2496

2597
HdfFile file( mDatFile );
@@ -38,8 +110,8 @@ void MDAL::LoaderXmdf::load( MDAL::Mesh *mesh, MDAL_Status *status )
38110

39111
// TODO: check version?
40112

41-
size_t vertexCount = mesh->vertices.size();
42-
size_t faceCount = mesh->faces.size();
113+
size_t vertexCount = mesh->verticesCount();
114+
size_t faceCount = mesh->facesCount();
43115
std::vector<std::string> rootGroups = file.groups();
44116
if ( rootGroups.size() != 1 )
45117
{
@@ -53,31 +125,40 @@ void MDAL::LoaderXmdf::load( MDAL::Mesh *mesh, MDAL_Status *status )
53125

54126
DatasetGroups groups; // DAT outputs data
55127

56-
HdfGroup gTemporal = gMesh.group( "Temporal" );
57-
if ( gTemporal.isValid() )
128+
if ( gMesh.pathExists( "Temporal" ) )
58129
{
59-
addDatasetGroupsFromXmdfGroup( groups, gTemporal, vertexCount, faceCount );
130+
HdfGroup gTemporal = gMesh.group( "Temporal" );
131+
if ( gTemporal.isValid() )
132+
{
133+
addDatasetGroupsFromXmdfGroup( groups, gTemporal, vertexCount, faceCount );
134+
}
60135
}
61136

62-
HdfGroup gMaximums = gMesh.group( "Maximums" );
63-
if ( gMaximums.isValid() )
137+
if ( gMesh.pathExists( "Temporal" ) )
64138
{
65-
for ( const std::string &name : gMaximums.groups() )
139+
HdfGroup gMaximums = gMesh.group( "Maximums" );
140+
if ( gMaximums.isValid() )
66141
{
67-
HdfGroup g = gMaximums.group( name );
68-
std::shared_ptr<MDAL::DatasetGroup> maxGroup = readXmdfGroupAsDatasetGroup( g, name + "/Maximums", vertexCount, faceCount );
69-
if ( maxGroup->datasets.size() != 1 )
70-
MDAL::debug( "Maximum dataset should have just one timestep!" );
71-
else
72-
groups.push_back( maxGroup );
142+
for ( const std::string &name : gMaximums.groups() )
143+
{
144+
HdfGroup g = gMaximums.group( name );
145+
std::shared_ptr<MDAL::DatasetGroup> maxGroup = readXmdfGroupAsDatasetGroup( g, name + "/Maximums", vertexCount, faceCount );
146+
if ( !maxGroup || maxGroup->datasets.size() != 1 )
147+
MDAL::debug( "Maximum dataset should have just one timestep!" );
148+
else
149+
groups.push_back( maxGroup );
150+
}
73151
}
74152
}
75153

76154
// res_to_res.exe (TUFLOW utiity tool)
77-
HdfGroup gDifference = gMesh.group( "Difference" );
78-
if ( gDifference.isValid() )
155+
if ( gMesh.pathExists( "Difference" ) )
79156
{
80-
addDatasetGroupsFromXmdfGroup( groups, gDifference, vertexCount, faceCount );
157+
HdfGroup gDifference = gMesh.group( "Difference" );
158+
if ( gDifference.isValid() )
159+
{
160+
addDatasetGroupsFromXmdfGroup( groups, gDifference, vertexCount, faceCount );
161+
}
81162
}
82163

83164
mesh->datasetGroups.insert(
@@ -93,18 +174,21 @@ void MDAL::LoaderXmdf::addDatasetGroupsFromXmdfGroup( DatasetGroups &groups, con
93174
{
94175
HdfGroup g = rootGroup.group( name );
95176
std::shared_ptr<DatasetGroup> ds = readXmdfGroupAsDatasetGroup( g, name, vertexCount, faceCount );
96-
groups.push_back( ds );
177+
if ( ds && ds->datasets.size() > 0 )
178+
groups.push_back( ds );
97179
}
98180
}
99181

100182
std::shared_ptr<MDAL::DatasetGroup> MDAL::LoaderXmdf::readXmdfGroupAsDatasetGroup(
101183
const HdfGroup &rootGroup, const std::string &name, size_t vertexCount, size_t faceCount )
102184
{
103-
std::shared_ptr<DatasetGroup> group( new DatasetGroup() );
185+
std::shared_ptr<DatasetGroup> group;
104186
std::vector<std::string> gDataNames = rootGroup.datasets();
105187
if ( !MDAL::contains( gDataNames, "Times" ) ||
106188
!MDAL::contains( gDataNames, "Values" ) ||
107-
!MDAL::contains( gDataNames, "Active" ) )
189+
!MDAL::contains( gDataNames, "Active" ) ||
190+
!MDAL::contains( gDataNames, "Mins" ) ||
191+
!MDAL::contains( gDataNames, "Maxs" ) )
108192
{
109193
MDAL::debug( "ignoring dataset " + name + " - not having required arrays" );
110194
return group;
@@ -113,19 +197,31 @@ std::shared_ptr<MDAL::DatasetGroup> MDAL::LoaderXmdf::readXmdfGroupAsDatasetGrou
113197
HdfDataset dsTimes = rootGroup.dataset( "Times" );
114198
HdfDataset dsValues = rootGroup.dataset( "Values" );
115199
HdfDataset dsActive = rootGroup.dataset( "Active" );
200+
HdfDataset dsMins = rootGroup.dataset( "Mins" );
201+
HdfDataset dsMaxs = rootGroup.dataset( "Maxs" );
116202

117203
std::vector<hsize_t> dimTimes = dsTimes.dims();
118204
std::vector<hsize_t> dimValues = dsValues.dims();
119205
std::vector<hsize_t> dimActive = dsActive.dims();
206+
std::vector<hsize_t> dimMins = dsMins.dims();
207+
std::vector<hsize_t> dimMaxs = dsMaxs.dims();
120208

121-
if ( dimTimes.size() != 1 || ( dimValues.size() != 2 && dimValues.size() != 3 ) || dimActive.size() != 2 )
209+
if ( dimTimes.size() != 1 ||
210+
( dimValues.size() != 2 && dimValues.size() != 3 ) ||
211+
dimActive.size() != 2 ||
212+
dimMins.size() != 1 ||
213+
dimMaxs.size() != 1
214+
)
122215
{
123216
MDAL::debug( "ignoring dataset " + name + " - arrays not having correct dimension counts" );
124217
return group;
125218
}
126219
hsize_t nTimeSteps = dimTimes[0];
127220

128-
if ( dimValues[0] != nTimeSteps || dimActive[0] != nTimeSteps )
221+
if ( dimValues[0] != nTimeSteps ||
222+
dimActive[0] != nTimeSteps ||
223+
dimMins[0] != nTimeSteps ||
224+
dimMaxs[0] != nTimeSteps )
129225
{
130226
MDAL::debug( "ignoring dataset " + name + " - arrays not having correct dimension sizes" );
131227
return group;
@@ -139,44 +235,32 @@ std::shared_ptr<MDAL::DatasetGroup> MDAL::LoaderXmdf::readXmdfGroupAsDatasetGrou
139235
bool isVector = dimValues.size() == 3;
140236

141237
std::vector<float> times = dsTimes.readArray();
142-
std::vector<float> values = dsValues.readArray();
143-
std::vector<uchar> active = dsActive.readArrayUint8();
144238

145-
group->setName( name );
146-
group->isScalar = !isVector;
147-
group->isOnVertices = true;
148-
group->uri = mDatFile;
149-
for ( hsize_t i = 0; i < nTimeSteps; ++i )
150-
{
151-
std::shared_ptr<Dataset> dataset( new Dataset() );
152-
dataset->values.resize( vertexCount );
153-
dataset->active.resize( faceCount );
154-
dataset->parent = group.get();
155-
dataset->time = double( times[i] );
239+
// all fine, set group and return
240+
group = std::make_shared<MDAL::DatasetGroup>(
241+
mMesh,
242+
mDatFile,
243+
name
244+
);
245+
group->setIsScalar( !isVector );
246+
group->setIsOnVertices( true );
156247

157-
if ( isVector )
158-
{
159-
const float *input = values.data() + 2 * i * vertexCount;
160-
for ( size_t j = 0; j < vertexCount; ++j )
161-
{
162-
dataset->values[j].x = double( input[2 * j] );
163-
dataset->values[j].y = double( input[2 * j + 1] );
164-
}
165-
}
166-
else
167-
{
168-
const float *input = values.data() + i * vertexCount;
169-
for ( size_t j = 0; j < vertexCount; ++j )
170-
{
171-
dataset->values[j].x = double( input[j] );
172-
}
173-
}
248+
// lazy loading of min and max of the dataset group
249+
std::vector<float> mins = dsMins.readArray();
250+
std::vector<float> maxs = dsMaxs.readArray();
251+
Statistics stats;
252+
stats.minimum = static_cast<double>( *std::min_element( mins.begin(), mins.end() ) );
253+
stats.maximum = static_cast<double>( *std::max_element( maxs.begin(), maxs.end() ) );
254+
group->setStatistics( stats );
174255

175-
const uchar *input = active.data() + i * faceCount;
176-
for ( size_t j = 0; j < faceCount; ++j )
177-
{
178-
dataset->active[j] = input[j];
179-
}
256+
for ( hsize_t i = 0; i < nTimeSteps; ++i )
257+
{
258+
std::shared_ptr<XmdfDataset> dataset( new XmdfDataset( group.get(), dsValues, dsActive, i ) );
259+
dataset->setTime( double( times[i] ) );
260+
Statistics stats;
261+
stats.minimum = static_cast<double>( mins[i] );
262+
stats.maximum = static_cast<double>( maxs[i] );
263+
dataset->setStatistics( stats );
180264
group->datasets.push_back( dataset );
181265
}
182266

‎external/mdal/frmts/mdal_xmdf.hpp

Lines changed: 35 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -20,13 +20,48 @@
2020
namespace MDAL
2121
{
2222

23+
/**
24+
* The XmdfDataset reads the data directly from HDF5 file
25+
* by usage of hyperslabs retrieval
26+
*
27+
* basically all (timesteps) data for one particular dataset groups
28+
* are stored in single
29+
* 3D arrays (time, x, y) for vector datasets
30+
* 2D arrays (time, x) for scalar datasets
31+
* 2D arrays (time, active) for active flags
32+
*/
33+
class XmdfDataset: public Dataset
34+
{
35+
public:
36+
XmdfDataset( DatasetGroup *grp,
37+
const HdfDataset &valuesDs,
38+
const HdfDataset &activeDs,
39+
hsize_t timeIndex );
40+
~XmdfDataset() override;
41+
42+
size_t scalarData( size_t indexStart, size_t count, double *buffer ) override;
43+
size_t vectorData( size_t indexStart, size_t count, double *buffer ) override;
44+
size_t activeData( size_t indexStart, size_t count, int *buffer ) override;
45+
46+
const HdfDataset &dsValues() const;
47+
const HdfDataset &dsActive() const;
48+
hsize_t timeIndex() const;
49+
50+
private:
51+
HdfDataset mHdf5DatasetValues;
52+
HdfDataset mHdf5DatasetActive;
53+
// index or row where the data for this timestep begins
54+
hsize_t mTimeIndex;
55+
};
56+
2357
class LoaderXmdf
2458
{
2559
public:
2660
LoaderXmdf( const std::string &datFile );
2761
void load( Mesh *mesh, MDAL_Status *status );
2862

2963
private:
64+
MDAL::Mesh *mMesh = nullptr;
3065
std::string mDatFile;
3166
std::shared_ptr<MDAL::DatasetGroup> readXmdfGroupAsDatasetGroup(
3267
const HdfGroup &rootGroup,

‎external/mdal/mdal.cpp

Lines changed: 222 additions & 157 deletions
Large diffs are not rendered by default.

‎external/mdal/mdal_data_model.cpp

Lines changed: 178 additions & 29 deletions
Original file line numberDiff line numberDiff line change
@@ -5,25 +5,87 @@
55

66
#include "mdal_data_model.hpp"
77
#include <assert.h>
8+
#include <math.h>
89
#include <algorithm>
910
#include "mdal_utils.hpp"
1011

11-
bool MDAL::Dataset::isActive( size_t faceIndex )
12+
MDAL::Dataset::~Dataset() = default;
13+
14+
MDAL::Dataset::Dataset( MDAL::DatasetGroup *parent )
15+
: mParent( parent )
16+
{
17+
assert( mParent );
18+
}
19+
20+
size_t MDAL::Dataset::valuesCount() const
1221
{
13-
assert( parent );
14-
if ( parent->isOnVertices )
22+
if ( group()->isOnVertices() )
1523
{
16-
if ( active.size() > faceIndex )
17-
return active[faceIndex];
18-
else
19-
return false;
24+
return mesh()->verticesCount();
2025
}
2126
else
2227
{
23-
return true;
28+
return mesh()->facesCount();
2429
}
2530
}
2631

32+
MDAL::Statistics MDAL::Dataset::statistics() const
33+
{
34+
return mStatistics;
35+
}
36+
37+
void MDAL::Dataset::setStatistics( const MDAL::Statistics &statistics )
38+
{
39+
mStatistics = statistics;
40+
}
41+
42+
MDAL::DatasetGroup *MDAL::Dataset::group() const
43+
{
44+
return mParent;
45+
}
46+
47+
MDAL::Mesh *MDAL::Dataset::mesh() const
48+
{
49+
return mParent->mesh();
50+
}
51+
52+
double MDAL::Dataset::time() const
53+
{
54+
return mTime;
55+
}
56+
57+
void MDAL::Dataset::setTime( double time )
58+
{
59+
mTime = time;
60+
}
61+
62+
bool MDAL::Dataset::isValid() const
63+
{
64+
return mIsValid;
65+
}
66+
67+
void MDAL::Dataset::setIsValid( bool isValid )
68+
{
69+
mIsValid = isValid;
70+
}
71+
72+
MDAL::DatasetGroup::DatasetGroup( MDAL::Mesh *parent,
73+
const std::string &uri,
74+
const std::string &name )
75+
: mParent( parent )
76+
, mUri( uri )
77+
{
78+
assert( mParent );
79+
setName( name );
80+
}
81+
82+
MDAL::DatasetGroup::DatasetGroup( MDAL::Mesh *parent, const std::string &uri )
83+
: mParent( parent )
84+
, mUri( uri )
85+
{
86+
assert( mParent );
87+
}
88+
2789
std::string MDAL::DatasetGroup::getMetadata( const std::string &key )
2890
{
2991
for ( auto &pair : metadata )
@@ -61,9 +123,66 @@ void MDAL::DatasetGroup::setName( const std::string &name )
61123
setMetadata( "name", name );
62124
}
63125

126+
std::string MDAL::DatasetGroup::uri() const
127+
{
128+
return mUri;
129+
}
130+
131+
MDAL::Statistics MDAL::DatasetGroup::statistics() const
132+
{
133+
return mStatistics;
134+
}
135+
136+
void MDAL::DatasetGroup::setStatistics( const Statistics &statistics )
137+
{
138+
mStatistics = statistics;
139+
}
140+
141+
MDAL::Mesh *MDAL::DatasetGroup::mesh() const
142+
{
143+
return mParent;
144+
}
145+
146+
bool MDAL::DatasetGroup::isOnVertices() const
147+
{
148+
return mIsOnVertices;
149+
}
150+
151+
void MDAL::DatasetGroup::setIsOnVertices( bool isOnVertices )
152+
{
153+
// datasets are initialized (e.g. values array, active array) based
154+
// on this property. Do not allow to modify later on.
155+
assert( datasets.empty() );
156+
mIsOnVertices = isOnVertices;
157+
}
158+
159+
bool MDAL::DatasetGroup::isScalar() const
160+
{
161+
return mIsScalar;
162+
}
163+
164+
void MDAL::DatasetGroup::setIsScalar( bool isScalar )
165+
{
166+
// datasets are initialized (e.g. values array, active array) based
167+
// on this property. Do not allow to modify later on.
168+
assert( datasets.empty() );
169+
mIsScalar = isScalar;
170+
}
171+
172+
MDAL::Mesh::Mesh( size_t verticesCount, size_t facesCount, size_t faceVerticesMaximumCount, MDAL::BBox extent, const std::string &uri )
173+
: mVerticesCount( verticesCount )
174+
, mFacesCount( facesCount )
175+
, mFaceVerticesMaximumCount( faceVerticesMaximumCount )
176+
, mExtent( extent )
177+
, mUri( uri )
178+
{
179+
}
180+
181+
MDAL::Mesh::~Mesh() = default;
182+
64183
void MDAL::Mesh::setSourceCrs( const std::string &str )
65184
{
66-
crs = MDAL::trim( str );
185+
mCrs = MDAL::trim( str );
67186
}
68187

69188
void MDAL::Mesh::setSourceCrsFromWKT( const std::string &wkt )
@@ -76,26 +195,56 @@ void MDAL::Mesh::setSourceCrsFromEPSG( int code )
76195
setSourceCrs( std::string( "EPSG:" ) + std::to_string( code ) );
77196
}
78197

79-
void MDAL::Mesh::addBedElevationDataset()
198+
void MDAL::Mesh::setExtent( const BBox &extent )
80199
{
81-
if ( faces.empty() )
82-
return;
200+
mExtent = extent;
201+
}
83202

84-
std::shared_ptr<DatasetGroup> group = std::make_shared< DatasetGroup >();
85-
group->isOnVertices = true;
86-
group->isScalar = true;
87-
group->setName( "Bed Elevation" );
88-
group->uri = uri;
89-
std::shared_ptr<MDAL::Dataset> dataset = std::make_shared< Dataset >();
90-
dataset->time = 0.0;
91-
dataset->values.resize( vertices.size() );
92-
dataset->active.resize( faces.size() );
93-
dataset->parent = group.get();
94-
std::fill( dataset->active.begin(), dataset->active.end(), 1 );
95-
for ( size_t i = 0; i < vertices.size(); ++i )
96-
{
97-
dataset->values[i].x = vertices[i].z;
98-
}
99-
group->datasets.push_back( dataset );
100-
datasetGroups.push_back( group );
203+
void MDAL::Mesh::setFaceVerticesMaximumCount( size_t faceVerticesMaximumCount )
204+
{
205+
mFaceVerticesMaximumCount = faceVerticesMaximumCount;
206+
}
207+
208+
void MDAL::Mesh::setFacesCount( size_t facesCount )
209+
{
210+
mFacesCount = facesCount;
101211
}
212+
213+
void MDAL::Mesh::setVerticesCount( size_t verticesCount )
214+
{
215+
mVerticesCount = verticesCount;
216+
}
217+
218+
size_t MDAL::Mesh::verticesCount() const
219+
{
220+
return mVerticesCount;
221+
}
222+
223+
size_t MDAL::Mesh::facesCount() const
224+
{
225+
return mFacesCount;
226+
}
227+
228+
std::string MDAL::Mesh::uri() const
229+
{
230+
return mUri;
231+
}
232+
233+
MDAL::BBox MDAL::Mesh::extent() const
234+
{
235+
return mExtent;
236+
}
237+
238+
std::string MDAL::Mesh::crs() const
239+
{
240+
return mCrs;
241+
}
242+
243+
size_t MDAL::Mesh::faceVerticesMaximumCount() const
244+
{
245+
return mFaceVerticesMaximumCount;
246+
}
247+
248+
MDAL::MeshVertexIterator::~MeshVertexIterator() = default;
249+
250+
MDAL::MeshFaceIterator::~MeshFaceIterator() = default;

‎external/mdal/mdal_data_model.hpp

Lines changed: 109 additions & 50 deletions
Original file line numberDiff line numberDiff line change
@@ -3,18 +3,20 @@
33
Copyright (C) 2018 Peter Petrik (zilolv at gmail dot com)
44
*/
55

6-
#ifndef MDAL_DEFINES_HPP
7-
#define MDAL_DEFINES_HPP
6+
#ifndef MDAL_DATA_MODEL_HPP
7+
#define MDAL_DATA_MODEL_HPP
88

99
#include <stddef.h>
1010
#include <vector>
1111
#include <memory>
1212
#include <map>
1313
#include <string>
14+
#include "mdal.h"
1415

1516
namespace MDAL
1617
{
1718
class DatasetGroup;
19+
class Mesh;
1820

1921
struct BBox
2022
{
@@ -27,89 +29,146 @@ namespace MDAL
2729
double maxY;
2830
};
2931

30-
31-
typedef struct
32-
{
33-
double x;
34-
double y;
35-
double z; // Bed elevation
36-
} Vertex;
37-
38-
typedef std::vector<size_t> Face;
39-
40-
typedef std::vector<Vertex> Vertices;
41-
typedef std::vector<Face> Faces;
42-
4332
typedef struct
4433
{
45-
double x;
46-
double y;
47-
48-
bool noData = false;
49-
} Value; //Dataset Value
34+
double minimum = std::numeric_limits<double>::quiet_NaN();
35+
double maximum = std::numeric_limits<double>::quiet_NaN();
36+
} Statistics;
5037

5138
typedef std::vector< std::pair< std::string, std::string > > Metadata;
5239

5340
class Dataset
5441
{
5542
public:
56-
double time;
43+
Dataset( DatasetGroup *parent );
44+
virtual ~Dataset();
45+
46+
size_t valuesCount() const;
47+
virtual size_t scalarData( size_t indexStart, size_t count, double *buffer ) = 0;
48+
virtual size_t vectorData( size_t indexStart, size_t count, double *buffer ) = 0;
49+
virtual size_t activeData( size_t indexStart, size_t count, int *buffer ) = 0;
5750

58-
/**
59-
* size - face count if !isOnVertices
60-
* size - vertex count if isOnVertices
61-
*/
62-
std::vector<Value> values;
63-
std::vector<bool> active; // size - face count. Whether the output for this is active...
51+
Statistics statistics() const;
52+
void setStatistics( const Statistics &statistics );
6453

65-
bool isValid = true;
66-
DatasetGroup *parent = nullptr;
54+
bool isValid() const;
55+
void setIsValid( bool isValid );
6756

68-
bool isActive( size_t faceIndex );
57+
DatasetGroup *group() const;
58+
Mesh *mesh() const;
59+
60+
double time() const;
61+
void setTime( double time );
62+
63+
private:
64+
double mTime = std::numeric_limits<double>::quiet_NaN();
65+
bool mIsValid = true;
66+
DatasetGroup *mParent = nullptr;
67+
Statistics mStatistics;
6968
};
7069

7170
typedef std::vector<std::shared_ptr<Dataset>> Datasets;
7271

7372
class DatasetGroup
7473
{
7574
public:
76-
std::string getMetadata( const std::string &key );
75+
DatasetGroup( Mesh *parent,
76+
const std::string &uri );
77+
78+
DatasetGroup( Mesh *parent,
79+
const std::string &uri,
80+
const std::string &name );
7781

82+
std::string getMetadata( const std::string &key );
7883
void setMetadata( const std::string &key, const std::string &val );
7984

8085
std::string name();
8186
void setName( const std::string &name );
8287

8388
Metadata metadata;
84-
85-
bool isScalar = true;
86-
bool isOnVertices = true;
8789
Datasets datasets;
88-
std::string uri; // file/uri from where it came
90+
91+
bool isScalar() const;
92+
void setIsScalar( bool isScalar );
93+
94+
bool isOnVertices() const;
95+
void setIsOnVertices( bool isOnVertices );
96+
97+
std::string uri() const;
98+
99+
Statistics statistics() const;
100+
void setStatistics( const Statistics &statistics );
101+
102+
Mesh *mesh() const;
103+
104+
private:
105+
Mesh *mParent = nullptr;
106+
bool mIsScalar = true;
107+
bool mIsOnVertices = true;
108+
std::string mUri; // file/uri from where it came
109+
Statistics mStatistics;
89110
};
90111

91112
typedef std::vector<std::shared_ptr<DatasetGroup>> DatasetGroups;
92113

93-
struct Mesh
114+
class MeshVertexIterator
94115
{
95-
std::string uri; // file/uri from where it came
96-
std::string crs;
97-
98-
Vertices vertices;
99-
std::map<size_t, size_t> vertexIDtoIndex; // only for 2DM and DAT files
100-
101-
Faces faces;
102-
std::map<size_t, size_t> faceIDtoIndex; // only for 2DM and DAT files
116+
public:
117+
virtual ~MeshVertexIterator();
103118

104-
DatasetGroups datasetGroups;
119+
virtual size_t next( size_t vertexCount, double *coordinates ) = 0;
120+
};
105121

106-
void setSourceCrs( const std::string &str );
107-
void setSourceCrsFromWKT( const std::string &wkt );
108-
void setSourceCrsFromEPSG( int code );
122+
class MeshFaceIterator
123+
{
124+
public:
125+
virtual ~MeshFaceIterator();
109126

110-
void addBedElevationDataset();
127+
virtual size_t next( size_t faceOffsetsBufferLen,
128+
int *faceOffsetsBuffer,
129+
size_t vertexIndicesBufferLen,
130+
int *vertexIndicesBuffer ) = 0;
111131
};
112132

133+
class Mesh
134+
{
135+
public:
136+
Mesh( size_t verticesCount,
137+
size_t facesCount,
138+
size_t faceVerticesMaximumCount,
139+
BBox extent,
140+
const std::string &uri );
141+
virtual ~Mesh();
142+
143+
void setSourceCrs( const std::string &str );
144+
void setSourceCrsFromWKT( const std::string &wkt );
145+
void setSourceCrsFromEPSG( int code );
146+
147+
void setVerticesCount( size_t verticesCount );
148+
void setFacesCount( size_t facesCount );
149+
void setFaceVerticesMaximumCount( size_t faceVerticesMaximumCount );
150+
void setExtent( const BBox &extent );
151+
152+
virtual std::unique_ptr<MDAL::MeshVertexIterator> readVertices() = 0;
153+
virtual std::unique_ptr<MDAL::MeshFaceIterator> readFaces() = 0;
154+
155+
DatasetGroups datasetGroups;
156+
157+
size_t verticesCount() const;
158+
size_t facesCount() const;
159+
std::string uri() const;
160+
BBox extent() const;
161+
std::string crs() const;
162+
size_t faceVerticesMaximumCount() const;
163+
164+
private:
165+
size_t mVerticesCount = 0;
166+
size_t mFacesCount = 0;
167+
size_t mFaceVerticesMaximumCount = 0; //typically 3 or 4, sometimes up to 9
168+
BBox mExtent;
169+
const std::string mUri; // file/uri from where it came
170+
std::string mCrs;
171+
};
113172
} // namespace MDAL
114-
#endif //MDAL_DEFINES_HPP
173+
#endif //MDAL_DATA_MODEL_HPP
115174

Lines changed: 207 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,207 @@
1+
/*
2+
MDAL - Mesh Data Abstraction Library (MIT License)
3+
Copyright (C) 2018 Peter Petrik (zilolv at gmail dot com)
4+
*/
5+
6+
#include "mdal_memory_data_model.hpp"
7+
#include <assert.h>
8+
#include <math.h>
9+
#include <cstring>
10+
#include <algorithm>
11+
#include <iterator>
12+
#include "mdal_utils.hpp"
13+
14+
MDAL::MemoryDataset::MemoryDataset( MDAL::DatasetGroup *grp )
15+
: Dataset( grp )
16+
, mValues( group()->isScalar() ? valuesCount() : 2 * valuesCount(),
17+
std::numeric_limits<double>::quiet_NaN() )
18+
, mActive( group()->isOnVertices() ? mesh()->facesCount() : 0,
19+
1 )
20+
{
21+
}
22+
23+
MDAL::MemoryDataset::~MemoryDataset() = default;
24+
25+
int *MDAL::MemoryDataset::active()
26+
{
27+
return mActive.data();
28+
}
29+
30+
double *MDAL::MemoryDataset::values()
31+
{
32+
return mValues.data();
33+
}
34+
35+
const int *MDAL::MemoryDataset::constActive() const
36+
{
37+
return mActive.data();
38+
}
39+
40+
const double *MDAL::MemoryDataset::constValues() const
41+
{
42+
return mValues.data();
43+
}
44+
45+
size_t MDAL::MemoryDataset::activeData( size_t indexStart, size_t count, int *buffer )
46+
{
47+
if ( group()->isOnVertices() )
48+
{
49+
size_t nValues = mActive.size();
50+
51+
if ( ( count < 1 ) || ( indexStart >= nValues ) )
52+
return 0;
53+
54+
size_t copyValues = std::min( nValues - indexStart, count );
55+
memcpy( buffer, constActive() + indexStart, copyValues * sizeof( int ) );
56+
return copyValues;
57+
}
58+
else
59+
{
60+
memset( buffer, true, count * sizeof( int ) );
61+
}
62+
63+
return count;
64+
}
65+
66+
size_t MDAL::MemoryDataset::scalarData( size_t indexStart, size_t count, double *buffer )
67+
{
68+
assert( group()->isScalar() ); //checked in C API interface
69+
size_t nValues = valuesCount();
70+
assert( mValues.size() == nValues );
71+
72+
if ( ( count < 1 ) || ( indexStart >= nValues ) )
73+
return 0;
74+
75+
size_t copyValues = std::min( nValues - indexStart, count );
76+
memcpy( buffer, constValues() + indexStart, copyValues * sizeof( double ) );
77+
return copyValues;
78+
}
79+
80+
size_t MDAL::MemoryDataset::vectorData( size_t indexStart, size_t count, double *buffer )
81+
{
82+
assert( !group()->isScalar() ); //checked in C API interface
83+
size_t nValues = valuesCount();
84+
assert( mValues.size() == nValues * 2 );
85+
86+
if ( ( count < 1 ) || ( indexStart >= nValues ) )
87+
return 0;
88+
89+
size_t copyValues = std::min( nValues - indexStart, count );
90+
memcpy( buffer, constValues() + 2 * indexStart, 2 * copyValues * sizeof( double ) );
91+
return copyValues;
92+
}
93+
94+
MDAL::MemoryMesh::MemoryMesh( size_t verticesCount, size_t facesCount, size_t faceVerticesMaximumCount, MDAL::BBox extent, const std::string &uri )
95+
: MDAL::Mesh( verticesCount, facesCount, faceVerticesMaximumCount, extent, uri )
96+
{
97+
}
98+
99+
std::unique_ptr<MDAL::MeshVertexIterator> MDAL::MemoryMesh::readVertices()
100+
{
101+
std::unique_ptr<MDAL::MemoryMeshVertexIterator> it( new MemoryMeshVertexIterator( this ) );
102+
return it;
103+
}
104+
105+
std::unique_ptr<MDAL::MeshFaceIterator> MDAL::MemoryMesh::readFaces()
106+
{
107+
std::unique_ptr<MDAL::MemoryMeshFaceIterator> it( new MemoryMeshFaceIterator( this ) );
108+
return it;
109+
}
110+
111+
void MDAL::MemoryMesh::addBedElevationDataset( const MDAL::Vertices &vertices, const MDAL::Faces &faces )
112+
{
113+
MDAL::addBedElevationDatasetGroup( this, vertices, faces );
114+
}
115+
116+
MDAL::MemoryMesh::~MemoryMesh() = default;
117+
118+
MDAL::MemoryMeshVertexIterator::MemoryMeshVertexIterator( const MDAL::MemoryMesh *mesh )
119+
: mMemoryMesh( mesh )
120+
{
121+
122+
}
123+
124+
MDAL::MemoryMeshVertexIterator::~MemoryMeshVertexIterator() = default;
125+
126+
size_t MDAL::MemoryMeshVertexIterator::next( size_t vertexCount, double *coordinates )
127+
{
128+
assert( mMemoryMesh );
129+
assert( coordinates );
130+
131+
size_t maxVertices = mMemoryMesh->verticesCount();
132+
133+
if ( vertexCount > maxVertices )
134+
return 0;
135+
136+
if ( mLastVertexIndex >= maxVertices )
137+
return 0;
138+
139+
size_t i = 0;
140+
141+
while ( true )
142+
{
143+
if ( mLastVertexIndex + i >= maxVertices )
144+
break;
145+
146+
if ( i >= vertexCount )
147+
break;
148+
149+
const Vertex v = mMemoryMesh->vertices[mLastVertexIndex + i];
150+
coordinates[3 * i] = v.x;
151+
coordinates[3 * i + 1] = v.y;
152+
coordinates[3 * i + 2] = v.z;
153+
154+
++i;
155+
}
156+
157+
mLastVertexIndex += i;
158+
return i;
159+
}
160+
161+
MDAL::MemoryMeshFaceIterator::MemoryMeshFaceIterator( const MDAL::MemoryMesh *mesh )
162+
: mMemoryMesh( mesh )
163+
{
164+
}
165+
166+
MDAL::MemoryMeshFaceIterator::~MemoryMeshFaceIterator() = default;
167+
168+
size_t MDAL::MemoryMeshFaceIterator::next(
169+
size_t faceOffsetsBufferLen, int *faceOffsetsBuffer,
170+
size_t vertexIndicesBufferLen, int *vertexIndicesBuffer )
171+
{
172+
assert( mMemoryMesh );
173+
assert( faceOffsetsBuffer );
174+
assert( vertexIndicesBuffer );
175+
176+
size_t maxFaces = mMemoryMesh->facesCount();
177+
size_t faceVerticesMaximumCount = mMemoryMesh->faceVerticesMaximumCount();
178+
size_t vertexIndex = 0;
179+
size_t faceIndex = 0;
180+
181+
while ( true )
182+
{
183+
if ( vertexIndex + faceVerticesMaximumCount > vertexIndicesBufferLen )
184+
break;
185+
186+
if ( faceIndex >= faceOffsetsBufferLen )
187+
break;
188+
189+
if ( mLastFaceIndex + faceIndex >= maxFaces )
190+
break;
191+
192+
const Face f = mMemoryMesh->faces[mLastFaceIndex + faceIndex];
193+
for ( size_t faceVertexIndex = 0; faceVertexIndex < f.size(); ++faceVertexIndex )
194+
{
195+
assert( vertexIndex < vertexIndicesBufferLen );
196+
vertexIndicesBuffer[vertexIndex] = static_cast<int>( f[faceVertexIndex] );
197+
++vertexIndex;
198+
}
199+
200+
assert( faceIndex < faceOffsetsBufferLen );
201+
faceOffsetsBuffer[faceIndex] = static_cast<int>( vertexIndex );
202+
++faceIndex;
203+
}
204+
205+
mLastFaceIndex += faceIndex;
206+
return faceIndex;
207+
}
Lines changed: 125 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,125 @@
1+
/*
2+
MDAL - Mesh Data Abstraction Library (MIT License)
3+
Copyright (C) 2018 Peter Petrik (zilolv at gmail dot com)
4+
*/
5+
6+
#ifndef MDAL_MEMORY_DATA_MODEL_HPP
7+
#define MDAL_MEMORY_DATA_MODEL_HPP
8+
9+
#include <stddef.h>
10+
#include <vector>
11+
#include <memory>
12+
#include <map>
13+
#include <string>
14+
#include "mdal.h"
15+
#include "mdal_data_model.hpp"
16+
17+
namespace MDAL
18+
{
19+
typedef struct
20+
{
21+
double x;
22+
double y;
23+
double z; // Bed elevation
24+
} Vertex;
25+
26+
typedef std::vector<size_t> Face;
27+
typedef std::vector<Vertex> Vertices;
28+
typedef std::vector<Face> Faces;
29+
30+
/**
31+
* The MemoryDataset stores all the data in the memory
32+
*/
33+
class MemoryDataset: public Dataset
34+
{
35+
public:
36+
MemoryDataset( DatasetGroup *grp );
37+
~MemoryDataset() override;
38+
39+
size_t scalarData( size_t indexStart, size_t count, double *buffer ) override;
40+
size_t vectorData( size_t indexStart, size_t count, double *buffer ) override;
41+
size_t activeData( size_t indexStart, size_t count, int *buffer ) override;
42+
43+
/**
44+
* valid pointer only for for dataset defined on vertices
45+
*/
46+
int *active();
47+
double *values();
48+
49+
const int *constActive() const;
50+
const double *constValues() const;
51+
52+
private:
53+
/**
54+
* Stores vector2d/scalar data for dataset in form
55+
* scalars: x1, x2, x3, ..., xN
56+
* vector2D: x1, y1, x2, y2, x3, y3, .... , xN, yN
57+
*
58+
* all values are initialized to std::numerical_limits<double>::quiet_NaN (==NODATA)
59+
*
60+
* size:
61+
* - face count if isOnFaces & isScalar
62+
* - vertex count if isOnVertices & isScalar
63+
* - face count * 2 if isOnFaces & isVector
64+
* - vertex count * 2 if isOnVertices & isVector
65+
*/
66+
std::vector<double> mValues;
67+
/**
68+
* Active flag, whether the face is active or not (disabled)
69+
* Only make sense for dataset defined on vertices with size == face count
70+
* For dataset defined on faces, this is empty vector
71+
*
72+
* Values are initialized by default to 1 (active)
73+
*/
74+
std::vector<int> mActive;
75+
};
76+
77+
class MemoryMesh: public Mesh
78+
{
79+
public:
80+
MemoryMesh( size_t verticesCount,
81+
size_t facesCount,
82+
size_t faceVerticesMaximumCount,
83+
BBox extent,
84+
const std::string &uri );
85+
~MemoryMesh() override;
86+
87+
std::unique_ptr<MDAL::MeshVertexIterator> readVertices() override;
88+
std::unique_ptr<MDAL::MeshFaceIterator> readFaces() override;
89+
90+
void addBedElevationDataset( const Vertices &vertices, const Faces &faces );
91+
92+
Vertices vertices;
93+
Faces faces;
94+
};
95+
96+
class MemoryMeshVertexIterator: public MeshVertexIterator
97+
{
98+
public:
99+
MemoryMeshVertexIterator( const MemoryMesh *mesh );
100+
~MemoryMeshVertexIterator() override;
101+
102+
size_t next( size_t vertexCount, double *coordinates ) override;
103+
104+
const MemoryMesh *mMemoryMesh;
105+
size_t mLastVertexIndex = 0;
106+
107+
};
108+
109+
class MemoryMeshFaceIterator: public MeshFaceIterator
110+
{
111+
public:
112+
MemoryMeshFaceIterator( const MemoryMesh *mesh );
113+
~MemoryMeshFaceIterator() override;
114+
115+
size_t next( size_t faceOffsetsBufferLen,
116+
int *faceOffsetsBuffer,
117+
size_t vertexIndicesBufferLen,
118+
int *vertexIndicesBuffer ) override;
119+
120+
const MemoryMesh *mMemoryMesh;
121+
size_t mLastFaceIndex = 0;
122+
123+
};
124+
} // namespace MDAL
125+
#endif //MDAL_MEMORY_DATA_MODEL_HPP

‎external/mdal/mdal_utils.cpp

Lines changed: 140 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -10,6 +10,7 @@
1010
#include <algorithm>
1111
#include <sstream>
1212
#include <math.h>
13+
#include <assert.h>
1314

1415
bool MDAL::fileExists( const std::string &filename )
1516
{
@@ -254,3 +255,142 @@ double MDAL::parseTimeUnits( const std::string &units )
254255

255256
return divBy;
256257
}
258+
259+
MDAL::Statistics _calculateStatistics( const std::vector<double> &values, size_t count, bool isVector )
260+
{
261+
MDAL::Statistics ret;
262+
263+
double min = std::numeric_limits<double>::quiet_NaN();
264+
double max = std::numeric_limits<double>::quiet_NaN();
265+
bool firstIteration = true;
266+
267+
for ( size_t i = 0; i < count; ++i )
268+
{
269+
double magnitude;
270+
if ( isVector )
271+
{
272+
double x = values[2 * i];
273+
double y = values[2 * i + 1];
274+
if ( isnan( x ) || isnan( y ) )
275+
continue;
276+
magnitude = sqrt( x * x + y * y );
277+
}
278+
else
279+
{
280+
double x = values[i];
281+
if ( isnan( x ) )
282+
continue;
283+
magnitude = x;
284+
}
285+
286+
if ( firstIteration )
287+
{
288+
firstIteration = false;
289+
min = magnitude;
290+
max = magnitude;
291+
}
292+
else
293+
{
294+
if ( magnitude < min )
295+
{
296+
min = magnitude;
297+
}
298+
if ( magnitude > max )
299+
{
300+
max = magnitude;
301+
}
302+
}
303+
}
304+
305+
ret.minimum = min;
306+
ret.maximum = max;
307+
return ret;
308+
}
309+
310+
MDAL::Statistics MDAL::calculateStatistics( std::shared_ptr<DatasetGroup> grp )
311+
{
312+
Statistics ret;
313+
if ( !grp )
314+
return ret;
315+
316+
for ( std::shared_ptr<Dataset> ds : grp->datasets )
317+
{
318+
MDAL::Statistics dsStats = ds->statistics();
319+
combineStatistics( ret, dsStats );
320+
}
321+
return ret;
322+
}
323+
324+
MDAL::Statistics MDAL::calculateStatistics( std::shared_ptr<Dataset> dataset )
325+
{
326+
Statistics ret;
327+
if ( !dataset )
328+
return ret;
329+
330+
bool isVector = !dataset->group()->isScalar();
331+
size_t bufLen = 2000;
332+
std::vector<double> buffer( isVector ? bufLen * 2 : bufLen );
333+
334+
size_t i = 0;
335+
while ( i < dataset->valuesCount() )
336+
{
337+
size_t valsRead;
338+
if ( isVector )
339+
{
340+
valsRead = dataset->vectorData( i, bufLen, buffer.data() );
341+
}
342+
else
343+
{
344+
valsRead = dataset->scalarData( i, bufLen, buffer.data() );
345+
}
346+
MDAL::Statistics dsStats = _calculateStatistics( buffer, valsRead, isVector );
347+
combineStatistics( ret, dsStats );
348+
i += valsRead;
349+
}
350+
351+
return ret;
352+
}
353+
354+
void MDAL::combineStatistics( MDAL::Statistics &main, const MDAL::Statistics &other )
355+
{
356+
if ( isnan( main.minimum ) ||
357+
( !isnan( other.minimum ) && ( main.minimum > other.minimum ) ) )
358+
{
359+
main.minimum = other.minimum;
360+
}
361+
362+
if ( isnan( main.maximum ) ||
363+
( !isnan( other.maximum ) && ( main.maximum < other.maximum ) ) )
364+
{
365+
main.maximum = other.maximum;
366+
}
367+
}
368+
369+
void MDAL::addBedElevationDatasetGroup( MDAL::Mesh *mesh, const Vertices &vertices, const Faces &faces )
370+
{
371+
if ( !mesh )
372+
return;
373+
374+
if ( 0 == mesh->facesCount() )
375+
return;
376+
377+
std::shared_ptr<DatasetGroup> group = std::make_shared< DatasetGroup >(
378+
mesh,
379+
mesh->uri(),
380+
"Bed Elevation"
381+
);
382+
group->setIsOnVertices( true );
383+
group->setIsScalar( true );
384+
385+
std::shared_ptr<MDAL::MemoryDataset> dataset = std::make_shared< MemoryDataset >( group.get() );
386+
dataset->setTime( 0.0 );
387+
double *vals = dataset->values();
388+
for ( size_t i = 0; i < vertices.size(); ++i )
389+
{
390+
vals[i] = vertices[i].z;
391+
}
392+
dataset->setStatistics( MDAL::calculateStatistics( dataset ) );
393+
group->datasets.push_back( dataset );
394+
group->setStatistics( MDAL::calculateStatistics( group ) );
395+
mesh->datasetGroups.push_back( group );
396+
}

‎external/mdal/mdal_utils.hpp

Lines changed: 15 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -12,9 +12,11 @@
1212
#include <limits>
1313

1414
#include "mdal_data_model.hpp"
15+
#include "mdal_memory_data_model.hpp"
1516

1617
// avoid unused variable warnings
1718
#define MDAL_UNUSED(x) (void)x;
19+
#define MDAL_NAN std::numeric_limits<double>::quiet_NaN()
1820

1921
namespace MDAL
2022
{
@@ -88,9 +90,21 @@ namespace MDAL
8890
BBox computeExtent( const Vertices &vertices );
8991

9092
// time
91-
9293
//! Returns a delimiter to get time in hours
9394
double parseTimeUnits( const std::string &units );
9495

96+
// statistics
97+
void combineStatistics( Statistics &main, const Statistics &other );
98+
99+
//! Calculates statistics for dataset group
100+
Statistics calculateStatistics( std::shared_ptr<DatasetGroup> grp );
101+
102+
//! Calculates statistics for dataset
103+
Statistics calculateStatistics( std::shared_ptr<Dataset> dataset );
104+
105+
// mesh & datasets
106+
//! Add bed elevatiom dataset group to mesh
107+
void addBedElevationDatasetGroup( MDAL::Mesh *mesh, const Vertices &vertices, const Faces &faces );
108+
95109
} // namespace MDAL
96110
#endif //MDAL_UTILS_HPP

‎src/providers/mdal/CMakeLists.txt

Lines changed: 3 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -36,6 +36,7 @@ IF (WITH_INTERNAL_MDAL)
3636
${CMAKE_SOURCE_DIR}/external/mdal/mdal_utils.cpp
3737
${CMAKE_SOURCE_DIR}/external/mdal/mdal_loader.cpp
3838
${CMAKE_SOURCE_DIR}/external/mdal/mdal_data_model.cpp
39+
${CMAKE_SOURCE_DIR}/external/mdal/mdal_memory_data_model.cpp
3940
${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_2dm.cpp
4041
${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_ascii_dat.cpp
4142
${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_binary_dat.cpp
@@ -46,13 +47,15 @@ IF (WITH_INTERNAL_MDAL)
4647
${CMAKE_SOURCE_DIR}/external/mdal/mdal_utils.hpp
4748
${CMAKE_SOURCE_DIR}/external/mdal/mdal_loader.hpp
4849
${CMAKE_SOURCE_DIR}/external/mdal/mdal_data_model.hpp
50+
${CMAKE_SOURCE_DIR}/external/mdal/mdal_memory_data_model.hpp
4951
${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_2dm.hpp
5052
${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_ascii_dat.hpp
5153
${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_binary_dat.hpp
5254
)
5355

5456
IF(HDF5_FOUND)
5557
SET(MDAL_LIB_SRCS ${MDAL_LIB_SRCS}
58+
${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_hdf5.cpp
5659
${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_xmdf.cpp
5760
)
5861
SET(MDAL_LIB_HDRS ${MDAL_LIB_HDRS}

0 commit comments

Comments
 (0)
Please sign in to comment.