Skip to content

Commit 435b594

Browse files
PeterPetrikPeterPetrik
committedMay 7, 2019
update to MDAL 0.3.2
·
ltr-3_40final-3_8_0
1 parent e213bde commit 435b594

20 files changed

+1208
-106
lines changed
 

‎external/mdal/frmts/mdal_3di.cpp

Lines changed: 5 additions & 11 deletions
Original file line numberDiff line numberDiff line change
@@ -20,24 +20,24 @@ MDAL::Driver3Di *MDAL::Driver3Di::create()
2020
return new Driver3Di();
2121
}
2222

23-
MDAL::CFDimensions MDAL::Driver3Di::populateDimensions( const NetCDFFile &ncFile )
23+
MDAL::CFDimensions MDAL::Driver3Di::populateDimensions( )
2424
{
2525
CFDimensions dims;
2626
size_t count;
2727
int ncid;
2828

2929
// 2D Mesh
30-
ncFile.getDimension( "nMesh2D_nodes", &count, &ncid );
30+
mNcFile.getDimension( "nMesh2D_nodes", &count, &ncid );
3131
dims.setDimension( CFDimensions::Face2D, count, ncid );
3232

33-
ncFile.getDimension( "nCorner_Nodes", &count, &ncid );
33+
mNcFile.getDimension( "nCorner_Nodes", &count, &ncid );
3434
dims.setDimension( CFDimensions::MaxVerticesInFace, count, ncid );
3535

3636
// Vertices count is populated later in populateFacesAndVertices
3737
// it is not known from the array dimensions
3838

3939
// Time
40-
ncFile.getDimension( "time", &count, &ncid );
40+
mNcFile.getDimension( "time", &count, &ncid );
4141
dims.setDimension( CFDimensions::Time, count, ncid );
4242

4343
return dims;
@@ -111,7 +111,7 @@ void MDAL::Driver3Di::populateFacesAndVertices( Vertices &vertices, Faces &faces
111111
mDimensions.setDimension( CFDimensions::Vertex2D, vertices.size() );
112112
}
113113

114-
void MDAL::Driver3Di::addBedElevation( MDAL::Mesh *mesh )
114+
void MDAL::Driver3Di::addBedElevation( MemoryMesh *mesh )
115115
{
116116
assert( mesh );
117117
if ( 0 == mesh->facesCount() )
@@ -189,12 +189,6 @@ std::set<std::string> MDAL::Driver3Di::ignoreNetCDFVariables()
189189
return ignore_variables;
190190
}
191191

192-
std::string MDAL::Driver3Di::nameSuffix( MDAL::CFDimensions::Type type )
193-
{
194-
MDAL_UNUSED( type );
195-
return "";
196-
}
197-
198192
void MDAL::Driver3Di::parseNetCDFVariableMetadata( int varid, const std::string &variableName, std::string &name, bool *is_vector, bool *is_x )
199193
{
200194
*is_vector = false;

‎external/mdal/frmts/mdal_3di.hpp

Lines changed: 2 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -44,12 +44,11 @@ namespace MDAL
4444
Driver3Di *create() override;
4545

4646
private:
47-
CFDimensions populateDimensions( const NetCDFFile &ncFile ) override;
47+
CFDimensions populateDimensions( ) override;
4848
void populateFacesAndVertices( Vertices &vertices, Faces &faces ) override;
49-
void addBedElevation( Mesh *mesh ) override;
49+
void addBedElevation( MemoryMesh *mesh ) override;
5050
std::string getCoordinateSystemVariableName() override;
5151
std::set<std::string> ignoreNetCDFVariables() override;
52-
std::string nameSuffix( CFDimensions::Type type ) override;
5352
void parseNetCDFVariableMetadata( int varid, const std::string &variableName,
5453
std::string &name, bool *is_vector, bool *is_x ) override;
5554

‎external/mdal/frmts/mdal_cf.cpp

Lines changed: 8 additions & 10 deletions
Original file line numberDiff line numberDiff line change
@@ -79,16 +79,13 @@ MDAL::cfdataset_info_map MDAL::DriverCF::parseDatasetGroupInfo()
7979
continue;
8080

8181
size_t arr_size = mDimensions.size( mDimensions.type( dimid ) ) * nTimesteps;
82-
std::string suffix = nameSuffix( mDimensions.type( dimid ) );
8382

8483
// Get name, if it is vector and if it is x or y
8584
std::string name;
8685
bool is_vector = true;
8786
bool is_x = false;
8887

8988
parseNetCDFVariableMetadata( varid, variable_name, name, &is_vector, &is_x );
90-
if ( !suffix.empty() )
91-
name = name + " (" + suffix + ")";
9289

9390
// Add it to the map
9491
auto it = dsinfo_map.find( name );
@@ -136,12 +133,12 @@ static void populate_vals( bool is_vector, double *vals, size_t i,
136133
{
137134
if ( is_vector )
138135
{
139-
vals[2 * i] = MDAL::safeValue( vals_x[idx], fill_val_x );
140-
vals[2 * i + 1] = MDAL::safeValue( vals_y[idx], fill_val_y );
136+
vals[2 * i] = MDAL::safeValue( vals_x[idx], fill_val_x );
137+
vals[2 * i + 1] = MDAL::safeValue( vals_y[idx], fill_val_y );
141138
}
142139
else
143140
{
144-
vals[i] = MDAL::safeValue( vals_x[idx], fill_val_x );
141+
vals[i] = MDAL::safeValue( vals_x[idx], fill_val_x );
145142
}
146143
}
147144

@@ -190,7 +187,7 @@ void MDAL::DriverCF::addDatasetGroups( MDAL::Mesh *mesh, const std::vector<doubl
190187

191188
// read X data
192189
double fill_val_x = mNcFile.getFillValue( dsi.ncid_x );
193-
std::vector<double> vals_x( dsi.arr_size );
190+
std::vector<double> vals_x( dsi.arr_size, std::numeric_limits<double>::quiet_NaN() );
194191
if ( nc_get_var_double( mNcFile.handle(), dsi.ncid_x, vals_x.data() ) ) CF_THROW_ERR;
195192

196193
// read Y data if vector
@@ -199,7 +196,7 @@ void MDAL::DriverCF::addDatasetGroups( MDAL::Mesh *mesh, const std::vector<doubl
199196
if ( dsi.is_vector )
200197
{
201198
fill_val_y = mNcFile.getFillValue( dsi.ncid_y );
202-
vals_y.resize( dsi.arr_size );
199+
vals_y.resize( dsi.arr_size, std::numeric_limits<double>::quiet_NaN() );
203200
if ( nc_get_var_double( mNcFile.handle(), dsi.ncid_y, vals_y.data() ) ) CF_THROW_ERR;
204201
}
205202

@@ -254,7 +251,8 @@ bool MDAL::DriverCF::canRead( const std::string &uri )
254251
{
255252
NetCDFFile ncFile;
256253
ncFile.openFile( uri );
257-
populateDimensions( ncFile );
254+
mNcFile = ncFile;
255+
populateDimensions( );
258256
}
259257
catch ( MDAL_Status )
260258
{
@@ -318,7 +316,7 @@ std::unique_ptr< MDAL::Mesh > MDAL::DriverCF::load( const std::string &fileName,
318316
mNcFile.openFile( mFileName );
319317

320318
// Parse dimensions
321-
mDimensions = populateDimensions( mNcFile );
319+
mDimensions = populateDimensions( );
322320

323321
// Create mMesh
324322
Faces faces;

‎external/mdal/frmts/mdal_cf.hpp

Lines changed: 3 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -63,6 +63,7 @@ namespace MDAL
6363

6464
//! NetCDF Climate and Forecast (CF) Metadata Conventions
6565
//! http://cfconventions.org
66+
//! and http://ugrid-conventions.github.io/ugrid-conventions/
6667
class DriverCF: public Driver
6768
{
6869
public:
@@ -74,12 +75,11 @@ namespace MDAL
7475
std::unique_ptr< Mesh > load( const std::string &fileName, MDAL_Status *status ) override;
7576

7677
protected:
77-
virtual CFDimensions populateDimensions( const NetCDFFile &ncFile ) = 0;
78+
virtual CFDimensions populateDimensions( ) = 0;
7879
virtual void populateFacesAndVertices( Vertices &vertices, Faces &faces ) = 0;
79-
virtual void addBedElevation( MDAL::Mesh *mesh ) = 0;
80+
virtual void addBedElevation( MDAL::MemoryMesh *mesh ) = 0;
8081
virtual std::string getCoordinateSystemVariableName() = 0;
8182
virtual std::set<std::string> ignoreNetCDFVariables() = 0;
82-
virtual std::string nameSuffix( CFDimensions::Type type ) = 0;
8383
virtual void parseNetCDFVariableMetadata( int varid, const std::string &variableName,
8484
std::string &name, bool *is_vector, bool *is_x ) = 0;
8585

‎external/mdal/frmts/mdal_flo2d.cpp

Lines changed: 6 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -109,6 +109,7 @@ void MDAL::DriverFlo2D::parseCADPTSFile( const std::string &datFileName, std::ve
109109
// CADPTS.DAT - COORDINATES OF CELL CENTERS (ELEM NUM, X, Y)
110110
while ( std::getline( cadptsStream, line ) )
111111
{
112+
line = MDAL::rtrim( line );
112113
std::vector<std::string> lineParts = MDAL::split( line, ' ' );
113114
if ( lineParts.size() != 3 )
114115
{
@@ -140,6 +141,7 @@ void MDAL::DriverFlo2D::parseFPLAINFile( std::vector<double> &elevations,
140141

141142
while ( std::getline( fplainStream, line ) )
142143
{
144+
line = MDAL::rtrim( line );
143145
std::vector<std::string> lineParts = MDAL::split( line, ' ' );
144146
if ( lineParts.size() != 7 )
145147
{
@@ -220,6 +222,7 @@ void MDAL::DriverFlo2D::parseTIMDEPFile( const std::string &datFileName, const s
220222

221223
while ( std::getline( inStream, line ) )
222224
{
225+
line = MDAL::rtrim( line );
223226
std::vector<std::string> lineParts = MDAL::split( line, ' ' );
224227
if ( lineParts.size() == 1 )
225228
{
@@ -301,6 +304,7 @@ void MDAL::DriverFlo2D::parseDEPTHFile( const std::string &datFileName, const st
301304
// DEPTH.OUT - COORDINATES (ELEM NUM, X, Y, MAX DEPTH)
302305
while ( std::getline( depthStream, line ) )
303306
{
307+
line = MDAL::rtrim( line );
304308
if ( vertex_idx == nVertices ) throw MDAL_Status::Err_IncompatibleMesh;
305309

306310
std::vector<std::string> lineParts = MDAL::split( line, ' ' );
@@ -348,6 +352,7 @@ void MDAL::DriverFlo2D::parseVELFPVELOCFile( const std::string &datFileName )
348352
{
349353
if ( vertex_idx == nVertices ) throw MDAL_Status::Err_IncompatibleMesh;
350354

355+
line = MDAL::rtrim( line );
351356
std::vector<std::string> lineParts = MDAL::split( line, ' ' );
352357
if ( lineParts.size() != 4 )
353358
{
@@ -378,6 +383,7 @@ void MDAL::DriverFlo2D::parseVELFPVELOCFile( const std::string &datFileName )
378383
{
379384
if ( vertex_idx == nVertices ) throw MDAL_Status::Err_IncompatibleMesh;
380385

386+
line = MDAL::rtrim( line );
381387
std::vector<std::string> lineParts = MDAL::split( line, ' ' );
382388
if ( lineParts.size() != 4 )
383389
{

‎external/mdal/frmts/mdal_gdal.cpp

Lines changed: 1 addition & 42 deletions
Original file line numberDiff line numberDiff line change
@@ -396,47 +396,6 @@ void MDAL::DriverGdal::addDataToOutput( GDALRasterBandH raster_band, std::shared
396396
}
397397
}
398398

399-
void MDAL::DriverGdal::activateFaces( std::shared_ptr<MemoryDataset> tos )
400-
{
401-
// only for data on vertices
402-
if ( !tos->group()->isOnVertices() )
403-
return;
404-
405-
bool isScalar = tos->group()->isScalar();
406-
407-
// Activate only Faces that do all Vertex's outputs with some data
408-
int *active = tos->active();
409-
const double *values = tos->constValues();
410-
411-
for ( unsigned int idx = 0; idx < meshGDALDataset()->mNVolumes; ++idx )
412-
{
413-
Face elem = mMesh->faces.at( idx );
414-
for ( size_t i = 0; i < 4; ++i )
415-
{
416-
const size_t vertexIndex = elem[i];
417-
if ( isScalar )
418-
{
419-
double val = values[vertexIndex];
420-
if ( std::isnan( val ) )
421-
{
422-
active[idx] = 0; //NOT ACTIVE
423-
break;
424-
}
425-
}
426-
else
427-
{
428-
double x = values[2 * vertexIndex];
429-
double y = values[2 * vertexIndex + 1];
430-
if ( std::isnan( x ) || std::isnan( y ) )
431-
{
432-
active[idx] = 0; //NOT ACTIVE
433-
break;
434-
}
435-
}
436-
}
437-
}
438-
}
439-
440399
void MDAL::DriverGdal::addDatasetGroups()
441400
{
442401
// Add dataset to mMesh
@@ -465,7 +424,7 @@ void MDAL::DriverGdal::addDatasetGroups()
465424
{
466425
addDataToOutput( raster_bands[i], dataset, is_vector, i == 0 );
467426
}
468-
activateFaces( dataset );
427+
MDAL::activateFaces( mMesh.get(), dataset );
469428
dataset->setStatistics( MDAL::calculateStatistics( dataset ) );
470429
group->datasets.push_back( dataset );
471430
}

‎external/mdal/frmts/mdal_gdal.hpp

Lines changed: 0 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -90,7 +90,6 @@ namespace MDAL
9090
metadata_hash parseMetadata( GDALMajorObjectH gdalBand, const char *pszDomain = nullptr );
9191
void addDataToOutput( GDALRasterBandH raster_band, std::shared_ptr<MemoryDataset> tos, bool is_vector, bool is_x );
9292
bool addSrcProj();
93-
void activateFaces( std::shared_ptr<MemoryDataset> tos );
9493
void addDatasetGroups();
9594
void createMesh();
9695
void parseRasterBands( const GdalDataset *cfGDALDataset );

‎external/mdal/frmts/mdal_netcdf.cpp

Lines changed: 10 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -103,21 +103,21 @@ std::string NetCDFFile::getAttrStr( const std::string &name, const std::string &
103103
return getAttrStr( attr_name, arr_id );
104104
}
105105

106-
std::string NetCDFFile::getAttrStr( const std::string &name, int varid ) const
106+
std::string NetCDFFile::getAttrStr( const std::string &attr_name, int varid ) const
107107
{
108108
assert( mNcid != 0 );
109109

110110
size_t attlen = 0;
111111

112-
if ( nc_inq_attlen( mNcid, varid, name.c_str(), &attlen ) )
112+
if ( nc_inq_attlen( mNcid, varid, attr_name.c_str(), &attlen ) )
113113
{
114114
// attribute is missing
115115
return std::string();
116116
}
117117

118118
char *string_attr = static_cast<char *>( malloc( attlen + 1 ) );
119119

120-
if ( nc_get_att_text( mNcid, varid, name.c_str(), string_attr ) ) throw MDAL_Status::Err_UnknownFormat;
120+
if ( nc_get_att_text( mNcid, varid, attr_name.c_str(), string_attr ) ) throw MDAL_Status::Err_UnknownFormat;
121121
string_attr[attlen] = '\0';
122122

123123
std::string res( string_attr );
@@ -139,6 +139,7 @@ double NetCDFFile::getAttrDouble( int varid, const std::string &attr_name ) cons
139139
return res;
140140
}
141141

142+
142143
int NetCDFFile::getVarId( const std::string &name )
143144
{
144145
int ncid_val;
@@ -153,3 +154,9 @@ void NetCDFFile::getDimension( const std::string &name, size_t *val, int *ncid_v
153154
if ( nc_inq_dimid( mNcid, name.c_str(), ncid_val ) != NC_NOERR ) throw MDAL_Status::Err_UnknownFormat;
154155
if ( nc_inq_dimlen( mNcid, *ncid_val, val ) != NC_NOERR ) throw MDAL_Status::Err_UnknownFormat;
155156
}
157+
158+
bool NetCDFFile::hasDimension( const std::string &name ) const
159+
{
160+
int ncid_val;
161+
return nc_inq_dimid( mNcid, name.c_str(), &ncid_val ) == NC_NOERR;
162+
}

‎external/mdal/frmts/mdal_netcdf.hpp

Lines changed: 9 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -28,11 +28,19 @@ class NetCDFFile
2828

2929
int getAttrInt( const std::string &name, const std::string &attr_name ) const;
3030
double getAttrDouble( int varid, const std::string &attr_name ) const;
31+
/**
32+
* Get string attribute
33+
* \param name name of the variable
34+
* \param attr_name name of the attribute of the variable
35+
* \returns empty string if attribute is missing, else attribute value
36+
*/
3137
std::string getAttrStr( const std::string &name, const std::string &attr_name ) const;
32-
std::string getAttrStr( const std::string &name, int varid ) const;
38+
std::string getAttrStr( const std::string &attr_name, int varid ) const;
3339
double getFillValue( int varid ) const;
3440
int getVarId( const std::string &name );
3541
void getDimension( const std::string &name, size_t *val, int *ncid_val ) const;
42+
bool hasDimension( const std::string &name ) const;
43+
3644
private:
3745
int mNcid; // C handle to the file
3846
};

‎external/mdal/frmts/mdal_selafin.cpp

Lines changed: 606 additions & 0 deletions
Large diffs are not rendered by default.

‎external/mdal/frmts/mdal_selafin.hpp

Lines changed: 99 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,99 @@
1+
/*
2+
MDAL - Mesh Data Abstraction Library (MIT License)
3+
Copyright (C) 2019 Peter Petrik (zilolv at gmail dot com)
4+
*/
5+
6+
#ifndef MDAL_SELAFIN_HPP
7+
#define MDAL_SELAFIN_HPP
8+
9+
#include <string>
10+
#include <memory>
11+
#include <map>
12+
#include <iostream>
13+
#include <fstream>
14+
15+
#include "mdal_data_model.hpp"
16+
#include "mdal_memory_data_model.hpp"
17+
#include "mdal.h"
18+
#include "mdal_driver.hpp"
19+
20+
namespace MDAL
21+
{
22+
class SerafinStreamReader
23+
{
24+
public:
25+
SerafinStreamReader();
26+
void initialize( const std::string &fileName );
27+
28+
std::string read_string( size_t len );
29+
std::vector<double> read_double_arr( size_t len );
30+
std::vector<int> read_int_arr( size_t len );
31+
std::vector<size_t> read_size_t_arr( size_t len );
32+
33+
double read_double( );
34+
int read_int( );
35+
size_t read_sizet( );
36+
37+
size_t remainingBytes();
38+
private:
39+
void ignore_array_length( );
40+
std::string read_string_without_length( size_t len );
41+
void ignore( int len );
42+
bool getStreamPrecision();
43+
44+
std::string mFileName;
45+
bool mStreamInFloatPrecision = true;
46+
bool mIsNativeLittleEndian = true;
47+
long long mFileSize = -1;
48+
std::ifstream mIn;
49+
};
50+
51+
/**
52+
* Serafin format (also called Selafin)
53+
*
54+
* Binary format for triangular mesh with datasets defined on vertices
55+
* http://www.opentelemac.org/downloads/Archive/v6p0/telemac2d_user_manual_v6p0.pdf Appendix 3
56+
* https://www.gdal.org/drv_selafin.html
57+
*/
58+
class DriverSelafin: public Driver
59+
{
60+
public:
61+
DriverSelafin();
62+
~DriverSelafin() override;
63+
DriverSelafin *create() override;
64+
65+
bool canRead( const std::string &uri ) override;
66+
std::unique_ptr< Mesh > load( const std::string &meshFile, MDAL_Status *status ) override;
67+
68+
private:
69+
typedef std::map<double, std::vector<double> > timestep_map; //TIME (sorted), nodeVal
70+
71+
void createMesh( double xOrigin,
72+
double yOrigin,
73+
size_t nElems,
74+
size_t nPoints,
75+
size_t nPointsPerElem,
76+
std::vector<size_t> &ikle,
77+
std::vector<double> &x,
78+
std::vector<double> &y );
79+
void addData( const std::vector<std::string> &var_names, const std::vector<timestep_map> &data, size_t nPoints );
80+
void parseFile( std::vector<std::string> &var_names,
81+
double *xOrigin,
82+
double *yOrigin,
83+
size_t *nElem,
84+
size_t *nPoint,
85+
size_t *nPointsPerElem,
86+
std::vector<size_t> &ikle,
87+
std::vector<double> &x,
88+
std::vector<double> &y,
89+
std::vector<timestep_map> &data );
90+
91+
bool getStreamPrecision( std::ifstream &in );
92+
93+
std::unique_ptr< MDAL::MemoryMesh > mMesh;
94+
std::string mFileName;
95+
SerafinStreamReader mReader;
96+
};
97+
98+
} // namespace MDAL
99+
#endif //MDAL_SELAFIN_HPP

‎external/mdal/frmts/mdal_ugrid.cpp

Lines changed: 323 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,323 @@
1+
/*
2+
MDAL - Mesh Data Abstraction Library (MIT License)
3+
Copyright (C) 2019 Peter Petrik (zilolv at gmail dot com)
4+
*/
5+
6+
#include "mdal_ugrid.hpp"
7+
#include "mdal_utils.hpp"
8+
9+
#include <netcdf.h>
10+
#include <assert.h>
11+
12+
MDAL::DriverUgrid::DriverUgrid()
13+
: DriverCF(
14+
"Ugrid",
15+
"UGRID Results",
16+
"*.nc" )
17+
{
18+
}
19+
20+
MDAL::DriverUgrid *MDAL::DriverUgrid::create()
21+
{
22+
return new DriverUgrid();
23+
}
24+
25+
std::string MDAL::DriverUgrid::findMeshName( int dimension, bool optional ) const
26+
{
27+
const std::vector<std::string> variables = mNcFile.readArrNames();
28+
for ( const std::string &varName : variables )
29+
{
30+
const std::string cf_role = mNcFile.getAttrStr( varName, "cf_role" );
31+
if ( cf_role == "mesh_topology" )
32+
{
33+
int topology_dimension = mNcFile.getAttrInt( varName, "topology_dimension" );
34+
if ( topology_dimension == dimension )
35+
{
36+
return varName;
37+
}
38+
}
39+
}
40+
if ( optional )
41+
return "";
42+
else
43+
throw MDAL_Status::Err_UnknownFormat;
44+
}
45+
46+
std::string MDAL::DriverUgrid::nodeZVariableName() const
47+
{
48+
// looks like mesh attributes does not have node_z array name
49+
// reference
50+
return mMesh2dName + "_node_z";
51+
}
52+
53+
MDAL::CFDimensions MDAL::DriverUgrid::populateDimensions( )
54+
{
55+
CFDimensions dims;
56+
size_t count;
57+
int ncid;
58+
59+
mMesh1dName = findMeshName( 1, true ); // optional, may not be present
60+
mMesh2dName = findMeshName( 2, false ); // force
61+
62+
// 2D Mesh
63+
64+
// node_dimension is usually something like nMesh2D_node
65+
// number of nodes/vertices in the mesh
66+
const std::string mesh2dNode = mNcFile.getAttrStr( mMesh2dName, "node_dimension" );
67+
mNcFile.getDimension( mesh2dNode, &count, &ncid );
68+
dims.setDimension( CFDimensions::Vertex2D, count, ncid );
69+
70+
// face_dimension is usually something like nMesh2D_face
71+
// number of faces in the mesh
72+
const std::string mesh2dFace = mNcFile.getAttrStr( mMesh2dName, "face_dimension" );
73+
mNcFile.getDimension( mesh2dFace, &count, &ncid );
74+
dims.setDimension( CFDimensions::Face2D, count, ncid );
75+
76+
// edge_dimension is usually something like nMesh2D_edge
77+
// number of edges in the mesh
78+
const std::string mesh2dEdge = mNcFile.getAttrStr( mMesh2dName, "edge_dimension" );
79+
mNcFile.getDimension( mesh2dEdge, &count, &ncid );
80+
dims.setDimension( CFDimensions::Face2DEdge, count, ncid );
81+
82+
// max_face_nodes_dimension is usually something like max_nMesh2D_face_nodes
83+
// maximum number of vertices in faces
84+
const std::string mesh2dMaxNodesInFace = mNcFile.getAttrStr( mMesh2dName, "max_face_nodes_dimension" );
85+
mNcFile.getDimension( mesh2dMaxNodesInFace, &count, &ncid );
86+
dims.setDimension( CFDimensions::MaxVerticesInFace, count, ncid );
87+
88+
// Time
89+
mNcFile.getDimension( "time", &count, &ncid );
90+
dims.setDimension( CFDimensions::Time, count, ncid );
91+
92+
return dims;
93+
}
94+
95+
void MDAL::DriverUgrid::populateFacesAndVertices( Vertices &vertices, Faces &faces )
96+
{
97+
populateVertices( vertices );
98+
populateFaces( faces );
99+
}
100+
101+
void MDAL::DriverUgrid::populateVertices( MDAL::Vertices &vertices )
102+
{
103+
assert( vertices.empty() );
104+
size_t vertexCount = mDimensions.size( CFDimensions::Vertex2D );
105+
vertices.resize( vertexCount );
106+
Vertex *vertexPtr = vertices.data();
107+
108+
// Parse 2D Mesh
109+
// node_coordinates should be something like Mesh2D_node_x Mesh2D_node_y
110+
std::string verticesXName, verticesYName;
111+
parse2VariablesFromAttribute( mMesh2dName, "node_coordinates", verticesXName, verticesYName, false );
112+
const std::vector<double> vertices2D_x = mNcFile.readDoubleArr( verticesXName, vertexCount );
113+
const std::vector<double> vertices2D_y = mNcFile.readDoubleArr( verticesYName, vertexCount );
114+
115+
std::vector<double> vertices2D_z;
116+
if ( mNcFile.hasArr( nodeZVariableName() ) )
117+
{
118+
vertices2D_z = mNcFile.readDoubleArr( nodeZVariableName(), vertexCount );
119+
}
120+
121+
for ( size_t i = 0; i < vertexCount; ++i, ++vertexPtr )
122+
{
123+
vertexPtr->x = vertices2D_x[i];
124+
vertexPtr->y = vertices2D_y[i];
125+
if ( !vertices2D_z.empty() )
126+
vertexPtr->z = vertices2D_z[i];
127+
}
128+
}
129+
130+
void MDAL::DriverUgrid::populateFaces( MDAL::Faces &faces )
131+
{
132+
assert( faces.empty() );
133+
size_t faceCount = mDimensions.size( CFDimensions::Face2D );
134+
faces.resize( faceCount );
135+
136+
// Parse 2D Mesh
137+
// face_node_connectivity is usually something like Mesh2D_face_nodes
138+
const std::string mesh2dFaceNodeConnectivity = mNcFile.getAttrStr( mMesh2dName, "face_node_connectivity" );
139+
140+
size_t verticesInFace = mDimensions.size( CFDimensions::MaxVerticesInFace );
141+
int fill_val = mNcFile.getAttrInt( mesh2dFaceNodeConnectivity, "_FillValue" );
142+
int start_index = mNcFile.getAttrInt( mesh2dFaceNodeConnectivity, "start_index" );
143+
std::vector<int> face_nodes_conn = mNcFile.readIntArr( mesh2dFaceNodeConnectivity, faceCount * verticesInFace );
144+
145+
for ( size_t i = 0; i < faceCount; ++i )
146+
{
147+
size_t nVertices = verticesInFace;
148+
std::vector<size_t> idxs;
149+
150+
for ( size_t j = 0; j < verticesInFace; ++j )
151+
{
152+
size_t idx = verticesInFace * i + j;
153+
int val = face_nodes_conn[idx];
154+
155+
if ( fill_val == val )
156+
{
157+
// found fill val
158+
nVertices = j;
159+
assert( nVertices > 1 );
160+
break;
161+
}
162+
else
163+
{
164+
idxs.push_back( static_cast<size_t>( val - start_index ) );
165+
}
166+
}
167+
faces[i] = idxs;
168+
}
169+
170+
}
171+
172+
void MDAL::DriverUgrid::addBedElevation( MDAL::MemoryMesh *mesh )
173+
{
174+
if ( mNcFile.hasArr( nodeZVariableName() ) ) MDAL::addBedElevationDatasetGroup( mesh, mesh->vertices );
175+
}
176+
177+
std::string MDAL::DriverUgrid::getCoordinateSystemVariableName()
178+
{
179+
std::string coordinate_system_variable;
180+
181+
// first try to get the coordinate system variable from grid definition
182+
if ( mNcFile.hasArr( nodeZVariableName() ) )
183+
{
184+
coordinate_system_variable = mNcFile.getAttrStr( nodeZVariableName(), "grid_mapping" );
185+
}
186+
187+
// if automatic discovery fails, try to check some hardcoded common variables that store projection
188+
if ( coordinate_system_variable.empty() )
189+
{
190+
if ( mNcFile.hasArr( "projected_coordinate_system" ) )
191+
coordinate_system_variable = "projected_coordinate_system";
192+
else if ( mNcFile.hasArr( "wgs84" ) )
193+
coordinate_system_variable = "wgs84";
194+
}
195+
196+
// return, may be empty
197+
return coordinate_system_variable;
198+
}
199+
200+
std::set<std::string> MDAL::DriverUgrid::ignoreNetCDFVariables()
201+
{
202+
std::set<std::string> ignore_variables;
203+
204+
ignore_variables.insert( "projected_coordinate_system" );
205+
ignore_variables.insert( "time" );
206+
ignore_variables.insert( "timestep" );
207+
208+
std::vector<std::string> meshes;
209+
if ( mNcFile.hasArr( mMesh1dName ) )
210+
meshes.push_back( mMesh1dName );
211+
meshes.push_back( mMesh2dName );
212+
213+
for ( const std::string &mesh : meshes )
214+
{
215+
std::string xName, yName;
216+
ignore_variables.insert( mesh );
217+
parse2VariablesFromAttribute( mesh, "node_coordinates", xName, yName, true );
218+
ignore_variables.insert( xName );
219+
ignore_variables.insert( yName );
220+
ignore_variables.insert( mNcFile.getAttrStr( mesh, "edge_node_connectivity" ) );
221+
parse2VariablesFromAttribute( mesh, "edge_coordinates", xName, yName, true );
222+
if ( !xName.empty() )
223+
{
224+
ignore_variables.insert( xName );
225+
ignore_variables.insert( mNcFile.getAttrStr( xName, "bounds" ) );
226+
}
227+
if ( !yName.empty() )
228+
{
229+
ignore_variables.insert( yName );
230+
ignore_variables.insert( mNcFile.getAttrStr( yName, "bounds" ) );
231+
}
232+
ignore_variables.insert( mNcFile.getAttrStr( mesh, "face_node_connectivity" ) );
233+
parse2VariablesFromAttribute( mesh, "face_coordinates", xName, yName, true );
234+
if ( !xName.empty() )
235+
{
236+
ignore_variables.insert( xName );
237+
ignore_variables.insert( mNcFile.getAttrStr( xName, "bounds" ) );
238+
}
239+
if ( !yName.empty() )
240+
{
241+
ignore_variables.insert( yName );
242+
ignore_variables.insert( mNcFile.getAttrStr( yName, "bounds" ) );
243+
}
244+
ignore_variables.insert( mNcFile.getAttrStr( mesh, "edge_face_connectivity" ) );
245+
}
246+
247+
return ignore_variables;
248+
}
249+
250+
void MDAL::DriverUgrid::parseNetCDFVariableMetadata( int varid, const std::string &variableName, std::string &name, bool *is_vector, bool *is_x )
251+
{
252+
*is_vector = false;
253+
*is_x = true;
254+
255+
std::string long_name = mNcFile.getAttrStr( "long_name", varid );
256+
if ( long_name.empty() )
257+
{
258+
std::string standard_name = mNcFile.getAttrStr( "standard_name", varid );
259+
if ( standard_name.empty() )
260+
{
261+
name = variableName;
262+
}
263+
else
264+
{
265+
if ( MDAL::contains( standard_name, "_x_" ) )
266+
{
267+
*is_vector = true;
268+
name = MDAL::replace( standard_name, "_x_", "" );
269+
}
270+
else if ( MDAL::contains( standard_name, "_y_" ) )
271+
{
272+
*is_vector = true;
273+
*is_x = false;
274+
name = MDAL::replace( standard_name, "_y_", "" );
275+
}
276+
else
277+
{
278+
name = standard_name;
279+
}
280+
}
281+
}
282+
else
283+
{
284+
if ( MDAL::contains( long_name, ", x-component" ) )
285+
{
286+
*is_vector = true;
287+
name = MDAL::replace( long_name, ", x-component", "" );
288+
}
289+
else if ( MDAL::contains( long_name, ", y-component" ) )
290+
{
291+
*is_vector = true;
292+
*is_x = false;
293+
name = MDAL::replace( long_name, ", y-component", "" );
294+
}
295+
else
296+
{
297+
name = long_name;
298+
}
299+
}
300+
}
301+
302+
void MDAL::DriverUgrid::parse2VariablesFromAttribute( const std::string &name, const std::string &attr_name,
303+
std::string &var1, std::string &var2, bool optional ) const
304+
{
305+
const std::string mesh2dNodeCoordinates = mNcFile.getAttrStr( name, attr_name );
306+
const std::vector<std::string> chunks = MDAL::split( mesh2dNodeCoordinates, ' ' );
307+
308+
if ( chunks.size() != 2 )
309+
{
310+
if ( optional )
311+
{
312+
var1 = "";
313+
var2 = "";
314+
}
315+
else
316+
throw MDAL_Status::Err_UnknownFormat;
317+
}
318+
else
319+
{
320+
var1 = chunks[0];
321+
var2 = chunks[1];
322+
}
323+
}

‎external/mdal/frmts/mdal_ugrid.hpp

Lines changed: 52 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,52 @@
1+
/*
2+
MDAL - Mesh Data Abstraction Library (MIT License)
3+
Copyright (C) 2019 Peter Petrik (zilolv at gmail dot com)
4+
*/
5+
6+
#ifndef MDAL_UGRID_HPP
7+
#define MDAL_UGRID_HPP
8+
9+
#include <map>
10+
#include <string>
11+
#include <stddef.h>
12+
13+
#include "mdal_cf.hpp"
14+
#include "mdal_driver.hpp"
15+
16+
namespace MDAL
17+
{
18+
/**
19+
* Driver of UGRID file format.
20+
*
21+
* The result UGRID NetCDF file is strictly based on CF-conventions 1.6
22+
*/
23+
class DriverUgrid: public DriverCF
24+
{
25+
public:
26+
DriverUgrid();
27+
~DriverUgrid() override = default;
28+
DriverUgrid *create() override;
29+
30+
private:
31+
CFDimensions populateDimensions( ) override;
32+
void populateFacesAndVertices( Vertices &vertices, Faces &faces ) override;
33+
void populateVertices( Vertices &vertices );
34+
void populateFaces( Faces &faces );
35+
void addBedElevation( MemoryMesh *mesh ) override;
36+
std::string getCoordinateSystemVariableName() override;
37+
std::set<std::string> ignoreNetCDFVariables() override;
38+
void parseNetCDFVariableMetadata( int varid, const std::string &variableName,
39+
std::string &name, bool *is_vector, bool *is_x ) override;
40+
41+
void parse2VariablesFromAttribute( const std::string &name, const std::string &attr_name,
42+
std::string &var1, std::string &var2,
43+
bool optional ) const;
44+
std::string findMeshName( int dimension, bool optional ) const;
45+
std::string mMesh2dName;
46+
std::string mMesh1dName;
47+
std::string nodeZVariableName() const;
48+
};
49+
50+
} // namespace MDAL
51+
52+
#endif // MDAL_UGRID_HPP

‎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.1";
25+
return "0.3.2";
2626
}
2727

2828
MDAL_Status MDAL_LastStatus()

‎external/mdal/mdal_data_model.cpp

Lines changed: 10 additions & 20 deletions
Original file line numberDiff line numberDiff line change
@@ -224,6 +224,16 @@ std::string MDAL::Mesh::driverName() const
224224

225225
MDAL::Mesh::~Mesh() = default;
226226

227+
std::shared_ptr<MDAL::DatasetGroup> MDAL::Mesh::group( const std::string &name )
228+
{
229+
for ( auto grp : datasetGroups )
230+
{
231+
if ( grp->name() == name )
232+
return grp;
233+
}
234+
return std::shared_ptr<MDAL::DatasetGroup>();
235+
}
236+
227237
void MDAL::Mesh::setSourceCrs( const std::string &str )
228238
{
229239
mCrs = MDAL::trim( str );
@@ -239,26 +249,6 @@ void MDAL::Mesh::setSourceCrsFromEPSG( int code )
239249
setSourceCrs( std::string( "EPSG:" ) + std::to_string( code ) );
240250
}
241251

242-
void MDAL::Mesh::setExtent( const BBox &extent )
243-
{
244-
mExtent = extent;
245-
}
246-
247-
void MDAL::Mesh::setFaceVerticesMaximumCount( size_t faceVerticesMaximumCount )
248-
{
249-
mFaceVerticesMaximumCount = faceVerticesMaximumCount;
250-
}
251-
252-
void MDAL::Mesh::setFacesCount( size_t facesCount )
253-
{
254-
mFacesCount = facesCount;
255-
}
256-
257-
void MDAL::Mesh::setVerticesCount( size_t verticesCount )
258-
{
259-
mVerticesCount = verticesCount;
260-
}
261-
262252
size_t MDAL::Mesh::verticesCount() const
263253
{
264254
return mVerticesCount;

‎external/mdal/mdal_data_model.hpp

Lines changed: 3 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -164,16 +164,14 @@ namespace MDAL
164164
void setSourceCrsFromWKT( const std::string &wkt );
165165
void setSourceCrsFromEPSG( int code );
166166

167-
void setVerticesCount( size_t verticesCount );
168-
void setFacesCount( size_t facesCount );
169-
void setFaceVerticesMaximumCount( size_t faceVerticesMaximumCount );
170-
void setExtent( const BBox &extent );
171-
172167
virtual std::unique_ptr<MDAL::MeshVertexIterator> readVertices() = 0;
173168
virtual std::unique_ptr<MDAL::MeshFaceIterator> readFaces() = 0;
174169

175170
DatasetGroups datasetGroups;
176171

172+
//! Find a dataset group by name
173+
std::shared_ptr<DatasetGroup> group( const std::string &name );
174+
177175
size_t verticesCount() const;
178176
size_t facesCount() const;
179177
std::string uri() const;

‎external/mdal/mdal_driver_manager.cpp

Lines changed: 4 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -8,6 +8,7 @@
88
#include "frmts/mdal_2dm.hpp"
99
#include "frmts/mdal_ascii_dat.hpp"
1010
#include "frmts/mdal_binary_dat.hpp"
11+
#include "frmts/mdal_selafin.hpp"
1112
#include "mdal_utils.hpp"
1213

1314
#ifdef HAVE_HDF5
@@ -21,6 +22,7 @@
2122
#endif
2223

2324
#ifdef HAVE_NETCDF
25+
#include "frmts/mdal_ugrid.hpp"
2426
#include "frmts/mdal_3di.hpp"
2527
#include "frmts/mdal_sww.hpp"
2628
#endif
@@ -121,6 +123,7 @@ MDAL::DriverManager::DriverManager()
121123
{
122124
// MESH DRIVERS
123125
mDrivers.push_back( std::make_shared<MDAL::Driver2dm>() );
126+
mDrivers.push_back( std::make_shared<MDAL::DriverSelafin>() );
124127

125128
#ifdef HAVE_HDF5
126129
mDrivers.push_back( std::make_shared<MDAL::DriverFlo2D>() );
@@ -130,6 +133,7 @@ MDAL::DriverManager::DriverManager()
130133
#ifdef HAVE_NETCDF
131134
mDrivers.push_back( std::make_shared<MDAL::Driver3Di>() );
132135
mDrivers.push_back( std::make_shared<MDAL::DriverSWW>() );
136+
mDrivers.push_back( std::make_shared<MDAL::DriverUgrid>() );
133137
#endif
134138

135139
#if defined HAVE_GDAL && defined HAVE_NETCDF

‎external/mdal/mdal_utils.cpp

Lines changed: 57 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -323,14 +323,16 @@ bool MDAL::equals( double val1, double val2, double eps )
323323

324324
double MDAL::safeValue( double val, double nodata, double eps )
325325
{
326+
if ( std::isnan( val ) )
327+
return val;
328+
329+
if ( std::isnan( nodata ) )
330+
return val;
331+
326332
if ( equals( val, nodata, eps ) )
327-
{
328333
return std::numeric_limits<double>::quiet_NaN();
329-
}
330-
else
331-
{
332-
return val;
333-
}
334+
335+
return val;
334336
}
335337

336338
double MDAL::parseTimeUnits( const std::string &units )
@@ -541,3 +543,52 @@ void MDAL::addFaceScalarDatasetGroup( MDAL::Mesh *mesh,
541543
group->setStatistics( MDAL::calculateStatistics( group ) );
542544
mesh->datasetGroups.push_back( group );
543545
}
546+
547+
void MDAL::activateFaces( MDAL::MemoryMesh *mesh, std::shared_ptr<MemoryDataset> dataset )
548+
{
549+
// only for data on vertices
550+
if ( !dataset->group()->isOnVertices() )
551+
return;
552+
553+
bool isScalar = dataset->group()->isScalar();
554+
555+
// Activate only Faces that do all Vertex's outputs with some data
556+
int *active = dataset->active();
557+
const double *values = dataset->constValues();
558+
const size_t nFaces = mesh->facesCount();
559+
560+
for ( unsigned int idx = 0; idx < nFaces; ++idx )
561+
{
562+
Face elem = mesh->faces.at( idx );
563+
for ( size_t i = 0; i < elem.size(); ++i )
564+
{
565+
const size_t vertexIndex = elem[i];
566+
if ( isScalar )
567+
{
568+
double val = values[vertexIndex];
569+
if ( std::isnan( val ) )
570+
{
571+
active[idx] = 0; //NOT ACTIVE
572+
break;
573+
}
574+
}
575+
else
576+
{
577+
double x = values[2 * vertexIndex];
578+
double y = values[2 * vertexIndex + 1];
579+
if ( std::isnan( x ) || std::isnan( y ) )
580+
{
581+
active[idx] = 0; //NOT ACTIVE
582+
break;
583+
}
584+
}
585+
}
586+
}
587+
}
588+
589+
bool MDAL::isNativeLittleEndian()
590+
{
591+
// https://stackoverflow.com/a/4181991/2838364
592+
int n = 1;
593+
return ( *( char * )&n == 1 );
594+
}

‎external/mdal/mdal_utils.hpp

Lines changed: 5 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -24,6 +24,9 @@
2424

2525
namespace MDAL
2626
{
27+
// endianness
28+
bool isNativeLittleEndian();
29+
2730
// numbers
2831
bool equals( double val1, double val2, double eps = std::numeric_limits<double>::epsilon() );
2932

@@ -112,6 +115,8 @@ namespace MDAL
112115
void addBedElevationDatasetGroup( MDAL::Mesh *mesh, const Vertices &vertices );
113116
//! Adds face scalar dataset group to mesh
114117
void addFaceScalarDatasetGroup( MDAL::Mesh *mesh, const std::vector<double> &values, const std::string &name );
118+
//! Loop through all faces and activate those which has all 4 values on vertices valid
119+
void activateFaces( MDAL::MemoryMesh *mesh, std::shared_ptr<MemoryDataset> dataset );
115120

116121
} // namespace MDAL
117122
#endif //MDAL_UTILS_HPP

‎src/providers/mdal/CMakeLists.txt

Lines changed: 4 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -42,6 +42,7 @@ IF (WITH_INTERNAL_MDAL)
4242
${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_2dm.cpp
4343
${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_ascii_dat.cpp
4444
${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_binary_dat.cpp
45+
${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_selafin.cpp
4546
)
4647

4748
SET(MDAL_LIB_HDRS
@@ -54,6 +55,7 @@ IF (WITH_INTERNAL_MDAL)
5455
${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_2dm.hpp
5556
${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_ascii_dat.hpp
5657
${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_binary_dat.hpp
58+
${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_selafin.hpp
5759
)
5860

5961
IF(HDF5_FOUND)
@@ -91,12 +93,14 @@ IF (WITH_INTERNAL_MDAL)
9193
${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_3di.cpp
9294
${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_netcdf.cpp
9395
${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_sww.cpp
96+
${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_ugrid.cpp
9497
)
9598
SET(MDAL_LIB_HDRS ${MDAL_LIB_HDRS}
9699
${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_cf.hpp
97100
${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_3di.hpp
98101
${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_netcdf.hpp
99102
${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_sww.hpp
103+
${CMAKE_SOURCE_DIR}/external/mdal/frmts/mdal_ugrid.hpp
100104
)
101105
SET (HAVE_NETCDF TRUE)
102106
ENDIF(NETCDF_FOUND)

0 commit comments

Comments
 (0)
Please sign in to comment.