Skip to content

Commit 97c9580

Browse files
PeterPetriknyalldawson
authored andcommitted
update MDAL to 0.0.7 (3di support, projections support, bugfixes)
1 parent d17912c commit 97c9580

19 files changed

+1147
-20
lines changed

external/mdal/api/mdal.h

Lines changed: 5 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -83,6 +83,11 @@ MDAL_EXPORT MeshH MDAL_LoadMesh( const char *meshFile );
8383
//! Closes mesh, frees the memory
8484
MDAL_EXPORT void MDAL_CloseMesh( MeshH mesh );
8585

86+
/**
87+
* Returns mesh projection
88+
* not thread-safe and valid only till next call
89+
*/
90+
MDAL_EXPORT const char *MDAL_M_projection( MeshH mesh );
8691
//! Returns vertex count for the mesh
8792
MDAL_EXPORT int MDAL_M_vertexCount( MeshH mesh );
8893
//! Returns vertex X coord for the mesh

external/mdal/cmake/FindNetCDF.cmake

Lines changed: 31 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,31 @@
1+
# - Find NetCDF
2+
# Find the native NetCDF includes and library
3+
# Copyright (c) 2018, Peter Petrik <zilolv at gmail dot com>
4+
#
5+
# This module returns these variables for the rest of the project to use.
6+
#
7+
# NETCDF_FOUND - True if NetCDF found including required interfaces (see below)
8+
# NETCDF_LIBRARY - All netcdf related libraries
9+
# NETCDF_INCLUDE_DIR - All directories to include
10+
11+
IF (NETCDF_INCLUDE_DIR AND NETCDF_LIBRARY)
12+
# Already in cache, be silent
13+
SET (NETCDF_FIND_QUIETLY TRUE)
14+
ENDIF ()
15+
16+
FIND_PACKAGE(PkgConfig)
17+
PKG_CHECK_MODULES(PC_NETCDF QUIET netcdf)
18+
SET(NETCDF_DEFINITIONS ${PC_NETCDF_CFLAGS_OTHER})
19+
20+
FIND_PATH (NETCDF_INCLUDE_DIR netcdf.h
21+
HINTS ${PC_NETCDF_INCLUDEDIR} ${PC_NETCDF_INCLUDE_DIRS} ${NETCDF_PREFIX}/include
22+
PATH_SUFFIXES libnetcdf )
23+
24+
FIND_LIBRARY (NETCDF_LIBRARY
25+
NAMES netcdf libnetcdf
26+
HINTS HINTS ${PC_NETCDF_LIBDIR} ${PC_NETCDF_LIBRARY_DIRS} ${NETCDF_PREFIX}/lib)
27+
28+
INCLUDE (FindPackageHandleStandardArgs)
29+
FIND_PACKAGE_HANDLE_STANDARD_ARGS (NetCDF
30+
DEFAULT_MSG NETCDF_LIBRARY NETCDF_INCLUDE_DIR)
31+

external/mdal/cmake_templates/mdal_config.hpp.in

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -7,6 +7,8 @@
77

88
#cmakedefine HAVE_GDAL
99

10+
#cmakedefine HAVE_NETCDF
11+
1012
#endif // MDAL_CONFIG_HPP
1113

1214

external/mdal/frmts/mdal_3di.cpp

Lines changed: 233 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,233 @@
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_3di.hpp"
7+
8+
MDAL::Loader3Di::Loader3Di( const std::string &fileName )
9+
: LoaderCF( fileName )
10+
{
11+
}
12+
13+
MDAL::CFDimensions MDAL::Loader3Di::populateDimensions()
14+
{
15+
CFDimensions dims;
16+
size_t count;
17+
int ncid;
18+
19+
// 2D Mesh
20+
mNcFile.getDimension( "nMesh2D_nodes", &count, &ncid );
21+
dims.setDimension( CFDimensions::Face2D, count, ncid );
22+
23+
mNcFile.getDimension( "nCorner_Nodes", &count, &ncid );
24+
dims.setDimension( CFDimensions::MaxVerticesInFace, count, ncid );
25+
26+
// Vertices count is populated later in populateFacesAndVertices
27+
// it is not known from the array dimensions
28+
29+
// Time
30+
mNcFile.getDimension( "time", &count, &ncid );
31+
dims.setDimension( CFDimensions::Time, count, ncid );
32+
33+
return dims;
34+
}
35+
36+
void MDAL::Loader3Di::populateFacesAndVertices( MDAL::Mesh *mesh )
37+
{
38+
assert( mesh );
39+
size_t faceCount = mDimensions.size( CFDimensions::Face2D );
40+
mesh->faces.resize( faceCount );
41+
size_t verticesInFace = mDimensions.size( CFDimensions::MaxVerticesInFace );
42+
size_t arrsize = faceCount * verticesInFace;
43+
std::map<std::string, size_t> xyToVertex2DId;
44+
45+
// X coordinate
46+
int ncidX = mNcFile.getVarId( "Mesh2DContour_x" );
47+
double fillX = mNcFile.getFillValue( ncidX );
48+
std::vector<double> faceVerticesX( arrsize );
49+
if ( nc_get_var_double( mNcFile.handle(), ncidX, faceVerticesX.data() ) ) throw MDAL_Status::Err_UnknownFormat;
50+
51+
// Y coordinate
52+
int ncidY = mNcFile.getVarId( "Mesh2DContour_y" );
53+
double fillY = mNcFile.getFillValue( ncidY );
54+
std::vector<double> faceVerticesY( arrsize );
55+
if ( nc_get_var_double( mNcFile.handle(), ncidY, faceVerticesY.data() ) ) throw MDAL_Status::Err_UnknownFormat;
56+
57+
// now populate create faces and backtrack which vertices
58+
// are used in multiple faces
59+
for ( size_t faceId = 0; faceId < faceCount; ++faceId )
60+
{
61+
Face face;
62+
63+
for ( size_t faceVertexId = 0; faceVertexId < verticesInFace; ++faceVertexId )
64+
{
65+
size_t arrId = faceId * verticesInFace + faceVertexId;
66+
Vertex vertex;
67+
vertex.x = faceVerticesX[arrId];
68+
vertex.y = faceVerticesY[arrId];
69+
vertex.z = 0;
70+
71+
if ( MDAL::equals( vertex.x, fillX ) || MDAL::equals( vertex.y, fillY ) )
72+
break;
73+
74+
75+
size_t vertexId;
76+
77+
std::string key = std::to_string( vertex.x ) + "," + std::to_string( vertex.y );
78+
const auto it = xyToVertex2DId.find( key );
79+
if ( it == xyToVertex2DId.end() )
80+
{
81+
// new vertex
82+
vertexId = mesh->vertices.size();
83+
xyToVertex2DId[key] = vertexId;
84+
mesh->vertices.push_back( vertex );
85+
}
86+
else
87+
{
88+
// existing vertex
89+
vertexId = it->second;
90+
}
91+
92+
face.push_back( vertexId );
93+
94+
}
95+
96+
mesh->faces[faceId] = face;
97+
}
98+
99+
// Only now we have number of vertices, since we identified vertices that
100+
// are used in multiple faces
101+
mDimensions.setDimension( CFDimensions::Vertex2D, mesh->vertices.size() );
102+
}
103+
104+
void MDAL::Loader3Di::addBedElevation( MDAL::Mesh *mesh )
105+
{
106+
assert( mesh );
107+
if ( mesh->faces.empty() )
108+
return;
109+
110+
size_t faceCount = mesh->faces.size();
111+
112+
// read Z coordinate of 3di computation nodes centers
113+
int ncidZ = mNcFile.getVarId( "Mesh2DFace_zcc" );
114+
double fillZ = mNcFile.getFillValue( ncidZ );
115+
std::vector<double> coordZ( faceCount );
116+
if ( nc_get_var_double( mNcFile.handle(), ncidZ, coordZ.data() ) )
117+
return; //error reading the array
118+
119+
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();
130+
for ( size_t i = 0; i < faceCount; ++i )
131+
{
132+
dataset->values[i].x = MDAL::safeValue( coordZ[i], fillZ );
133+
}
134+
group->datasets.push_back( dataset );
135+
mesh->datasetGroups.push_back( group );
136+
}
137+
138+
std::string MDAL::Loader3Di::getCoordinateSystemVariableName()
139+
{
140+
return "projected_coordinate_system";
141+
}
142+
143+
std::set<std::string> MDAL::Loader3Di::ignoreNetCDFVariables()
144+
{
145+
std::set<std::string> ignore_variables;
146+
147+
ignore_variables.insert( "projected_coordinate_system" );
148+
ignore_variables.insert( "time" );
149+
150+
std::vector<std::string> meshes;
151+
meshes.push_back( "Mesh1D" );
152+
meshes.push_back( "Mesh2D" );
153+
154+
for ( const std::string &mesh : meshes )
155+
{
156+
ignore_variables.insert( mesh );
157+
ignore_variables.insert( mesh + "Node_id" );
158+
ignore_variables.insert( mesh + "Node_type" );
159+
160+
ignore_variables.insert( mesh + "Face_xcc" );
161+
ignore_variables.insert( mesh + "Face_ycc" );
162+
ignore_variables.insert( mesh + "Face_zcc" );
163+
ignore_variables.insert( mesh + "Contour_x" );
164+
ignore_variables.insert( mesh + "Contour_y" );
165+
ignore_variables.insert( mesh + "Face_sumax" );
166+
167+
ignore_variables.insert( mesh + "Line_id" );
168+
ignore_variables.insert( mesh + "Line_xcc" );
169+
ignore_variables.insert( mesh + "Line_ycc" );
170+
ignore_variables.insert( mesh + "Line_zcc" );
171+
ignore_variables.insert( mesh + "Line_type" );
172+
}
173+
174+
return ignore_variables;
175+
}
176+
177+
std::string MDAL::Loader3Di::nameSuffix( MDAL::CFDimensions::Type type )
178+
{
179+
MDAL_UNUSED( type );
180+
return "";
181+
}
182+
183+
void MDAL::Loader3Di::parseNetCDFVariableMetadata( int varid, const std::string &variableName, std::string &name, bool *is_vector, bool *is_x )
184+
{
185+
*is_vector = false;
186+
*is_x = true;
187+
188+
std::string long_name = mNcFile.getAttrStr( "long_name", varid );
189+
if ( long_name.empty() )
190+
{
191+
std::string standard_name = mNcFile.getAttrStr( "standard_name", varid );
192+
if ( standard_name.empty() )
193+
{
194+
name = variableName;
195+
}
196+
else
197+
{
198+
if ( MDAL::contains( standard_name, "_x_" ) )
199+
{
200+
*is_vector = true;
201+
name = MDAL::replace( standard_name, "_x_", "" );
202+
}
203+
else if ( MDAL::contains( standard_name, "_y_" ) )
204+
{
205+
*is_vector = true;
206+
*is_x = false;
207+
name = MDAL::replace( standard_name, "_y_", "" );
208+
}
209+
else
210+
{
211+
name = standard_name;
212+
}
213+
}
214+
}
215+
else
216+
{
217+
if ( MDAL::contains( long_name, " in x direction" ) )
218+
{
219+
*is_vector = true;
220+
name = MDAL::replace( long_name, " in x direction", "" );
221+
}
222+
else if ( MDAL::contains( long_name, " in y direction" ) )
223+
{
224+
*is_vector = true;
225+
*is_x = false;
226+
name = MDAL::replace( long_name, " in y direction", "" );
227+
}
228+
else
229+
{
230+
name = long_name;
231+
}
232+
}
233+
}

external/mdal/frmts/mdal_3di.hpp

Lines changed: 62 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,62 @@
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_3DI_HPP
7+
#define MDAL_3DI_HPP
8+
9+
#include <map>
10+
#include <string>
11+
#include <stddef.h>
12+
13+
#include "mdal_cf.hpp"
14+
15+
namespace MDAL
16+
{
17+
18+
/**
19+
* Loader of 3Di file format.
20+
*
21+
* The result 3Di NetCDF file is based on CF-conventions with some additions.
22+
* It is unstructured grid with data stored in NetCDF/HDF5 file format.
23+
* A division is made between a 1D and 2D which can be distinguished through
24+
* the prefixes “MESH2D” and “MESH1D”. For both meshes the information is present
25+
* in coordinate, id, and type variables.
26+
*
27+
* A version of the data scheme is not present yet.
28+
*
29+
* The 2D Mesh consists of calculation Nodes that represents centers of Faces.
30+
* There is no concept of Vertices in the file. The vertices that forms a face
31+
* are specified by X,Y coordinates in the "Face Contours" arrays. The "lines"
32+
* represent the face's edges and are again specified by X,Y coordinate of the
33+
* line center. Data is specified on calculation nodes (i.e. dataset defined on faces)
34+
* and on lines (i.e. dataset defined on edges - not implemented yet)
35+
*
36+
* The 1D Mesh is present too, but not parsed yet.
37+
*/
38+
class Loader3Di: public LoaderCF
39+
{
40+
public:
41+
Loader3Di( const std::string &fileName );
42+
~Loader3Di() override = default;
43+
44+
private:
45+
CFDimensions populateDimensions() override;
46+
void populateFacesAndVertices( MDAL::Mesh *mesh ) override;
47+
void addBedElevation( MDAL::Mesh *mesh ) override;
48+
std::string getCoordinateSystemVariableName() override;
49+
std::set<std::string> ignoreNetCDFVariables() override;
50+
std::string nameSuffix( CFDimensions::Type type ) override;
51+
void parseNetCDFVariableMetadata( int varid, const std::string &variableName,
52+
std::string &name, bool *is_vector, bool *is_x ) override;
53+
54+
//! Returns number of vertices
55+
size_t parse2DMesh();
56+
57+
void addBedElevationDatasetOnFaces();
58+
};
59+
60+
} // namespace MDAL
61+
62+
#endif // MDAL_3DI_HPP

external/mdal/frmts/mdal_ascii_dat.cpp

Lines changed: 4 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -69,7 +69,7 @@ void MDAL::LoaderAsciiDat::load( MDAL::Mesh *mesh, MDAL_Status *status )
6969
oldFormat = true;
7070
isVector = ( line == "VECTOR" );
7171

72-
group.reset( new DatasetGroup() );
72+
group = std::make_shared< DatasetGroup >();
7373
group->uri = mDatFile;
7474
group->setName( name );
7575
group->isScalar = !isVector;
@@ -118,7 +118,7 @@ void MDAL::LoaderAsciiDat::load( MDAL::Mesh *mesh, MDAL_Status *status )
118118
}
119119
isVector = cardType == "BEGVEC";
120120

121-
group.reset( new DatasetGroup() );
121+
group = std::make_shared< DatasetGroup >();
122122
group->uri = mDatFile;
123123
group->setName( name );
124124
group->isScalar = !isVector;
@@ -197,7 +197,7 @@ void MDAL::LoaderAsciiDat::readVertexTimestep(
197197
size_t vertexCount = mesh->vertices.size();
198198
size_t faceCount = mesh->faces.size();
199199

200-
std::shared_ptr<MDAL::Dataset> dataset( new MDAL::Dataset );
200+
std::shared_ptr<MDAL::Dataset> dataset = std::make_shared< MDAL::Dataset >();
201201
dataset->time = t / 3600.; // TODO read TIMEUNITS
202202
dataset->values.resize( vertexCount );
203203
dataset->active.resize( faceCount );
@@ -267,7 +267,7 @@ void MDAL::LoaderAsciiDat::readFaceTimestep(
267267

268268
size_t faceCount = mesh->faces.size();
269269

270-
std::shared_ptr<MDAL::Dataset> dataset( new MDAL::Dataset );
270+
std::shared_ptr<MDAL::Dataset> dataset = std::make_shared< MDAL::Dataset >();
271271
dataset->time = t / 3600.;
272272
dataset->values.resize( faceCount );
273273
dataset->parent = group.get();

0 commit comments

Comments
 (0)