Skip to content

Commit f742044

Browse files
PeterPetrikPeterPetrik
committedApr 18, 2019
update MDAL 0.3.1
·
ltr-3_40final-3_8_0
1 parent fd492ab commit f742044

File tree

8 files changed

+120
-20
lines changed

8 files changed

+120
-20
lines changed
 

‎external/mdal/frmts/mdal_2dm.cpp

Lines changed: 31 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -13,6 +13,7 @@
1313
#include <map>
1414
#include <cassert>
1515
#include <limits>
16+
#include <algorithm>
1617

1718
#include "mdal_2dm.hpp"
1819
#include "mdal.h"
@@ -64,6 +65,19 @@ size_t MDAL::Mesh2dm::vertexIndex( size_t vertexID ) const
6465
return vertexID;
6566
}
6667

68+
size_t MDAL::Mesh2dm::maximumVertexId() const
69+
{
70+
size_t maxIndex = verticesCount() - 1;
71+
if ( mVertexIDtoIndex.empty() )
72+
return maxIndex;
73+
else
74+
{
75+
// std::map is sorted!
76+
size_t maxID = mVertexIDtoIndex.rbegin()->first;
77+
return std::max( maxIndex, maxID );
78+
}
79+
}
80+
6781
MDAL::Driver2dm::Driver2dm():
6882
Driver( DRIVER_NAME,
6983
"2DM Mesh File",
@@ -146,6 +160,7 @@ std::unique_ptr<MDAL::Mesh> MDAL::Driver2dm::load( const std::string &meshFile,
146160
size_t faceIndex = 0;
147161
size_t vertexIndex = 0;
148162
std::map<size_t, size_t> vertexIDtoIndex;
163+
size_t lastVertexID = 0;
149164

150165
while ( std::getline( in, line ) )
151166
{
@@ -205,7 +220,22 @@ std::unique_ptr<MDAL::Mesh> MDAL::Driver2dm::load( const std::string &meshFile,
205220
else if ( startsWith( line, "ND" ) )
206221
{
207222
chunks = split( line, ' ' );
208-
size_t nodeID = toSizeT( chunks[1] ) - 1; // 2dm is numbered from 1
223+
size_t nodeID = toSizeT( chunks[1] );
224+
225+
if ( nodeID != 0 )
226+
{
227+
// specification of 2DM states that ID should be positive integer numbered from 1
228+
// but it seems some formats do not respect that
229+
if ( ( lastVertexID != 0 ) && ( nodeID <= lastVertexID ) )
230+
{
231+
// the algorithm requires that the file has NDs orderer by index
232+
if ( status ) *status = MDAL_Status::Err_InvalidData;
233+
return nullptr;
234+
}
235+
lastVertexID = nodeID;
236+
}
237+
nodeID -= 1; // 2dm is numbered from 1
238+
209239
_parse_vertex_id_gaps( vertexIDtoIndex, vertexIndex, nodeID, status );
210240
assert( vertexIndex < vertexCount );
211241
Vertex &vertex = vertices[vertexIndex];

‎external/mdal/frmts/mdal_2dm.hpp

Lines changed: 19 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -29,21 +29,26 @@ namespace MDAL
2929
~Mesh2dm() override;
3030

3131

32-
//! Some formats supports gaps in the vertex indexing, but we return continuos array from MDAL
33-
//! for most of the formats this returns
32+
//! HYDRO_AS-2D supports gaps in the vertex indexing,
33+
//! but we use continuos array of vertices in MDAL
3434
//! \param vertexID internal index/ID of the vertex that native format uses
35-
//! \returns index of the vertex in the continuous array of vertices we returned by readVertices()
35+
//! \returns index of the vertex in the continuous array of vertices we returned by readVertices().
36+
//! For invalid vertexID it is returned index that is out of vertices array bounds.
3637
virtual size_t vertexIndex( size_t vertexID ) const;
3738

39+
//! Returns maximum vertex ID.
40+
//! For meshes without gaps in vertex indexing, it is vertex count - 1
41+
virtual size_t maximumVertexId() const;
42+
3843
private:
39-
// 2dm supports "gaps" in the mesh indexing
40-
// Store only the indices that have different index and ID
41-
// https://github.com/lutraconsulting/MDAL/issues/51
44+
//! 2dm supports "gaps" in the mesh indexing
45+
//! Store only the indices that have different index and ID
46+
//! https://github.com/lutraconsulting/MDAL/issues/51
4247
std::map<size_t, size_t> mVertexIDtoIndex;
4348
};
4449

4550
/**
46-
* 2DM format specification used in TUFLOW and BASEMENET solvers
51+
* 2DM format specification used in TUFLOW, HYDRO_AS-2D and BASEMENET solvers
4752
* Text file format representing mesh vertices (ND) and faces (E**)
4853
* ND id x y z
4954
* Supports lines, triangles and polygons up to 9 vertices (implemented triangles and quads)
@@ -60,6 +65,13 @@ namespace MDAL
6065
*
6166
* Note that some 2dm formats do have some extra columns after mat_id column with
6267
* data with unknown origin/name (e.g. tests/data/2dm/regular_grid.2dm)
68+
*
69+
* HYDRO_AS-2D also allows gaps in vertex indexing. In this case we support only files
70+
* where the vertices are sorted by ID in the source file (limitation of the implementation)
71+
*
72+
* Vertex/Face IDs should be indexed from 1. We support indexing from 0 for datasets in xmdf format,
73+
* but not for ascii dat format (since the loop is from 0 to maximumVertexId() which is in this case
74+
* numberical_limits<size_t>::max(); (limitation of the implementation)
6375
*/
6476
class Driver2dm: public Driver
6577
{

‎external/mdal/frmts/mdal_ascii_dat.cpp

Lines changed: 30 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -13,6 +13,7 @@
1313
#include <map>
1414
#include <cassert>
1515
#include <memory>
16+
#include <limits>
1617

1718
#include "mdal_ascii_dat.hpp"
1819
#include "mdal.h"
@@ -103,7 +104,8 @@ void MDAL::DriverAsciiDat::loadOldFormat( std::ifstream &in,
103104
if ( cardType == "ND" && items.size() >= 2 )
104105
{
105106
size_t fileNodeCount = toSizeT( items[1] );
106-
if ( mesh->verticesCount() != fileNodeCount )
107+
size_t meshIdCount = maximumId( mesh ) + 1;
108+
if ( meshIdCount != fileNodeCount )
107109
EXIT_WITH_ERROR( MDAL_Status::Err_IncompatibleMesh );
108110
}
109111
else if ( cardType == "SCALAR" || cardType == "VECTOR" )
@@ -167,7 +169,8 @@ void MDAL::DriverAsciiDat::loadNewFormat( std::ifstream &in,
167169
if ( cardType == "ND" && items.size() >= 2 )
168170
{
169171
size_t fileNodeCount = toSizeT( items[1] );
170-
if ( mesh->verticesCount() != fileNodeCount )
172+
size_t meshIdCount = maximumId( mesh ) + 1;
173+
if ( meshIdCount != fileNodeCount )
171174
EXIT_WITH_ERROR( MDAL_Status::Err_IncompatibleMesh );
172175
}
173176
else if ( cardType == "NC" && items.size() >= 2 )
@@ -247,6 +250,15 @@ void MDAL::DriverAsciiDat::loadNewFormat( std::ifstream &in,
247250
}
248251
}
249252

253+
size_t MDAL::DriverAsciiDat::maximumId( const MDAL::Mesh *mesh ) const
254+
{
255+
const Mesh2dm *m2dm = dynamic_cast<const Mesh2dm *>( mesh );
256+
if ( m2dm )
257+
return m2dm->maximumVertexId();
258+
else
259+
return mesh->verticesCount() - 1;
260+
}
261+
250262
/**
251263
* The DAT format contains "datasets" and each dataset has N-outputs. One output
252264
* represents data for all vertices/faces for one timestep
@@ -265,6 +277,14 @@ void MDAL::DriverAsciiDat::load( const std::string &datFile, MDAL::Mesh *mesh, M
265277
return;
266278
}
267279

280+
size_t mID = maximumId( mesh );
281+
if ( mID == std::numeric_limits<size_t>::max() )
282+
{
283+
// This happens when mesh is 2DM and vertices are numbered from 0
284+
if ( status ) *status = MDAL_Status::Err_IncompatibleMesh;
285+
return;
286+
}
287+
268288
std::ifstream in( mDatFile, std::ifstream::in );
269289
std::string line;
270290
if ( !std::getline( in, line ) )
@@ -316,17 +336,22 @@ void MDAL::DriverAsciiDat::readVertexTimestep(
316336

317337
const Mesh2dm *m2dm = dynamic_cast<const Mesh2dm *>( mesh );
318338
double *values = dataset->values();
319-
for ( size_t i = 0; i < mesh->verticesCount(); ++i )
339+
size_t meshIdCount = maximumId( mesh ) + 1; // these are native format indexes (IDs). For formats without gaps it equals vertex array index
340+
size_t vertexCount = mesh->verticesCount();
341+
342+
for ( size_t id = 0; id < meshIdCount; ++id )
320343
{
321344
std::string line;
322345
std::getline( stream, line );
323346
std::vector<std::string> tsItems = split( line, ' ' );
324347

325348
size_t index;
326349
if ( m2dm )
327-
index = m2dm->vertexIndex( i );
350+
index = m2dm->vertexIndex( id ); //this index may be out of values array
328351
else
329-
index = i;
352+
index = id;
353+
354+
if ( index >= vertexCount ) continue;
330355

331356
if ( isVector )
332357
{

‎external/mdal/frmts/mdal_ascii_dat.hpp

Lines changed: 9 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -35,6 +35,9 @@ namespace MDAL
3535
* such dataset, the dataset name contains "_els_" substring
3636
* (e.g. depth_els_1.dat)
3737
*
38+
* HYDRO_AS-2D solver can have mesh that has numbering gaps, but
39+
* speficies values for even missing indexes in dataset file
40+
*
3841
* In one file, there is always one dataset group stored.
3942
*
4043
* Sometime the "older" datasets may have some part of the
@@ -60,6 +63,12 @@ namespace MDAL
6063
void loadOldFormat( std::ifstream &in, Mesh *mesh, MDAL_Status *status ) const;
6164
void loadNewFormat( std::ifstream &in, Mesh *mesh, MDAL_Status *status ) const;
6265

66+
//! Gets maximum (native) index.
67+
//! For meshes without indexing gap it is vertexCount - 1
68+
//! For some HYDRO_AS-2D meshes with indexing gaps, it returns
69+
//! maximum native index of the vertex in defined in the mesh
70+
size_t maximumId( const Mesh *mesh ) const;
71+
6372
void readVertexTimestep(
6473
const Mesh *mesh,
6574
std::shared_ptr<DatasetGroup> group,

‎external/mdal/frmts/mdal_gdal.cpp

Lines changed: 26 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -317,7 +317,28 @@ void MDAL::DriverGdal::addDataToOutput( GDALRasterBandH raster_band, std::shared
317317
{
318318
assert( raster_band );
319319

320-
double nodata = GDALGetRasterNoDataValue( raster_band, nullptr );
320+
// nodata
321+
int pbSuccess;
322+
double nodata = GDALGetRasterNoDataValue( raster_band, &pbSuccess );
323+
if ( pbSuccess == 0 ) nodata = std::numeric_limits<double>::quiet_NaN();
324+
bool hasNoData = !std::isnan( nodata );
325+
326+
// offset and scale
327+
double offset = 0.0;
328+
double scale = GDALGetRasterScale( raster_band, &pbSuccess );
329+
if ( ( pbSuccess == 0 ) || MDAL::equals( scale, 0.0 ) || std::isnan( scale ) )
330+
{
331+
scale = 1.0;
332+
}
333+
else
334+
{
335+
offset = GDALGetRasterOffset( raster_band, &pbSuccess );
336+
if ( ( pbSuccess == 0 ) || std::isnan( offset ) )
337+
{
338+
offset = 0.0;
339+
}
340+
}
341+
321342
unsigned int mXSize = meshGDALDataset()->mXSize;
322343
unsigned int mYSize = meshGDALDataset()->mYSize;
323344

@@ -349,8 +370,11 @@ void MDAL::DriverGdal::addDataToOutput( GDALRasterBandH raster_band, std::shared
349370
{
350371
unsigned int idx = x + mXSize * y;
351372
double val = mPafScanline[x];
352-
if ( !MDAL::equals( val, nodata ) )
373+
if ( !hasNoData || !MDAL::equals( val, nodata ) )
353374
{
375+
// Apply scale and offset
376+
val = val * scale + offset;
377+
354378
// values is prepolulated with NODATA values, so store only legal values
355379
if ( is_vector )
356380
{

‎external/mdal/frmts/mdal_xdmf.cpp

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -153,7 +153,7 @@ size_t MDAL::XdmfFunctionDataset::scalarData( size_t indexStart, size_t count, d
153153
assert( mType != FunctionType::Join );
154154

155155
if ( mType == FunctionType::Subtract )
156-
return substractFunction( indexStart, count, buffer );
156+
return subtractFunction( indexStart, count, buffer );
157157

158158
if ( mType == FunctionType::Flow )
159159
return flowFunction( indexStart, count, buffer );
@@ -176,7 +176,7 @@ size_t MDAL::XdmfFunctionDataset::activeData( size_t indexStart, size_t count, i
176176
return count;
177177
}
178178

179-
size_t MDAL::XdmfFunctionDataset::substractFunction( size_t indexStart, size_t count, double *buffer )
179+
size_t MDAL::XdmfFunctionDataset::subtractFunction( size_t indexStart, size_t count, double *buffer )
180180
{
181181
std::vector<double> buf( 2 * count, std::numeric_limits<double>::quiet_NaN() );
182182
size_t copyVals = extractRawData( indexStart, count, 2, buf );

‎external/mdal/frmts/mdal_xdmf.hpp

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -86,7 +86,7 @@ namespace MDAL
8686
* Currently we do not use any fancy bison/flex based
8787
* expression parsing, just supporting few types of
8888
* most common function types:
89-
* - substraction (A-B)
89+
* - subtraction (A-B)
9090
* - join ( [A, B] vector)
9191
* - magnitude
9292
*
@@ -126,7 +126,7 @@ namespace MDAL
126126
size_t activeData( size_t indexStart, size_t count, int *buffer ) override;
127127

128128
private:
129-
size_t substractFunction( size_t indexStart, size_t count, double *buffer );
129+
size_t subtractFunction( size_t indexStart, size_t count, double *buffer );
130130
size_t flowFunction( size_t indexStart, size_t count, double *buffer );
131131
size_t joinFunction( size_t indexStart, size_t count, double *buffer );
132132
size_t extractRawData( size_t indexStart, size_t count, size_t nDatasets, std::vector<double> &buf );

‎external/mdal/mdal.cpp

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -22,7 +22,7 @@ static MDAL_Status sLastStatus;
2222

2323
const char *MDAL_Version()
2424
{
25-
return "0.3.0";
25+
return "0.3.1";
2626
}
2727

2828
MDAL_Status MDAL_LastStatus()

0 commit comments

Comments
 (0)
Please sign in to comment.