Unstructured Mesh Layers

[wonder-sk]

QGIS Enhancement: Unstructured Mesh Layer

Date 2018/03/28

Author Martin Dobias (@wonder-sk), Peter Petrik (@PeterPetrik)

Contact wonder dot sk at gmail dot com, zilolv at gmail dot com

Maintainer @PeterPetrik, @wonder-sk

Version QGIS 3.2 or 3.4

Summary

This QEP proposes implementation of a new map layer type: in addition to raster and vector layers,
there will be mesh layer supporting display of data on structured or unstructured meshes (grids).

The idea is to introduce foundations to QGIS so that it is possible to create visualizations like this one (done with Crayfish plugin):

What is a mesh?

In our context, a mesh is a collection of vertices, edges and faces in 2D or 3D space:

• vertices - XY(Z) points (in the layer's coordinate reference system)
• edges - connect pairs of vertices
• faces - sets of edges forming a closed shape - typically triangles or quadrilaterals (quads), rarely polygons with higher number of vertices

Mesh gives us information about the spatial structure. In addition to the mesh we have datasets that assign a value to every vertex. For example, ice cap thickness at particular moment of time. A single file may contain multiple datasets - typically multiple quantities (e.g. water depth, water flow) that may be varying in time (time being represented in discrete timesteps, so each quantity may have N arrays, one for each timestep). Datasets do not have to vary in time (e.g. maximum water depth over the whole simulation).

Here is an example of a triangular mesh with numbered vertices:

The following table gives an idea what information is stored in datasets. Table columns represent indices of mesh vertices, each row represents one dataset. The first two are scalar datasets, the latter two are datasets with 2D vectors.

	1	2	3	...	13
Water depth at time=0s	5	5	5	...	2
Water depth at time=60s	6	5	3	...	4
Water flow at time=0s	[1,2]	[2,2]	[3,2]	...	[1,2]
Water flow at time=60s	[3,2]	[3,2]	[2,2]	...	[4,2]

In some cases datasets assign values to faces or edges instead of assigning values to vertices.

We can visualize the data by assigning colors to values (similarly to how it is done with "Singleband pseudocolor" raster rendering) and interpolating data between vertices according to the mesh topology. It is common that some quantities are 2D vectors rather than being simple scalar values (e.g. wind direction). For such quantities it is very desired to display arrows indicating vector direction.

When would I use it?

Most often this kind of representation is used when preparing data for simulation software or when viewing results of physical simulations, typically for meteorology, oceanography, hydrological or hydraulic models. All computation in such software is done on meshes, with values (physical quantities) usually stored in vertices (less commonly in edges or faces). Results usually comprise of various quantities (e.g. wind speed, water depth) which may be also time-varying (e.g. calculated water flow estimates in 5 minute intervals).

Some of the modelling software packages are free and open source (e.g. AnuGA, EPANET), there are some freeware packages (e.g. Basement, HEC-RAS) as well as many commercial pieces of software.

Do we really need this in QGIS? Maybe raster and vector layers are good enough?

It will be very useful to have native support for this data representation. Mesh-based data cannot be properly represented either as vector or raster layers. Vector data in QGIS use "simple features" representation which looses topological relationships between vertices, edges and faces, moreover we would need potentially large number of attributes for features. Efficient access to data is also a concern - dataset sizes can easily get to gigabytes. Reusing existing vector data support would be therefore impractical and complicated - just like as if we wanted to use ordinary vector layers for point clouds: in theory it is possible, but the nature of the data is quite different.

Raster layers are not an answer either. Sometimes users take mesh-based data and export those to raster grid files, but this always comes at the expense of loss of information and it only makes sense for visualization of results, it does not allow manipulation of the raw data. Rasters in fact can be thought of as a special case of meshes - having all vertices in a grid with regular spacing and all faces being quads.

Mesh layers also come useful when dealing with weather data (historical or predicted). Although such data are commonly available as raster grids (with some formats such as GRIB being supported by GDAL), interpolation and rendering of vector data adds extra value when visualizing such datasets.

To summarize, native support for mesh data opens QGIS to a whole new set of use cases.

Isn't this is already handled by plugins?

Only to some degree. There is Crayfish plugin that deals with mesh data, but it some issues when it comes to distribution of it: the data access and rendering are written in C++ for efficiency which means that it needs users to either manually compile the c++ library from source code or to download pre-built binaries. Once support for mesh layers is available in QGIS 3, we would like to update Crayfish to become Python-only plugin making use of mesh layers, still providing support for various mesh formats and adding extra GUI functionality.

There is also Meshlayer plugin from Oslandia which is written in Python and uses OpenGL for rendering, however it has various limitations (only supports triangular meshes, no support for element-centered data, ...)

Proposed Solution

A mesh layer (QgsMeshLayer) consists of the following main components:

• data provider - just like how providers are used for vector and raster layers, this is an abstraction for access to mesh structure and to data on top of the mesh.

• renderer - takes care of rendering the mesh layer to map images according to the configuration specified by the user

• cache - since fetching data from data provider may be slow, a cache is used to keep some data in memory for fast access

Data provider

There may be multiple data providers available for different file formats, but one layer will always point to a single data provider (QgsMeshDataProvider)

The main responsibilities of a data provider are:

• return mesh structure (list of nodes, edges and faces)
• return metadata about datasets (similar to getting a list of raster bands and their metadata)
• return a window of data from a dataset

In the initial implementation we assume that data providers will be read-only. In the future this may change and there may be calls to modify mesh structure or add/remove/modify datasets.

Renderer

Mesh layer renderer (QgsMeshLayerRenderer) is responsible for drawing map of given map extent.
It is derived from QgsMapLayerRenderer and it gets called from QGIS map rendering engine.

Renderer renders data from the active dataset - users will be able to choose which dataset is the active one.

Rendering consists of three steps and the user should be able to configure which steps will be used:

• rendering of contours. This interpolates values from the active dataset and assigns a color to each pixel covered by the mesh. For vector datasets, magnitudes of vectors are used as the values.

  

• rendering of vectors. This draws arrow symbols with appropriate scaling for data from the active dataset. Obviously this is only possible for vector datasets.

  

• rendering of mesh. The active dataset is not used here - only mesh structure gets rendered, optionally with labels.

  

The rendering algorithm consists of the following steps:

1. Given map extent, use spatial index to figure out which elements of the mesh are needed for rendering
2. Using the list of elements, figure out which data blocks will be needed from the active dataset.
3. Try to fetch the required data blocks from the cache. If the blocks are not available in the cache, get them from the data provider and temporarily store in the cache.
4. Walk elements of the mesh in map extent and render one by one (first contours, then vectors, then mesh).

Cache

In order to make sure that the data access is fast enough, we will need to keep a cache. It will consist of the following parts:

1. Mesh structure. This data structure will be accessed very often and it is relatively small so it can be cached as a whole.

2. Data blocks. Access to the data from data providers may be relatively slow and very often the map renderer will need access to the same or similar subset of data over time. At the same time the total amount of data stored in data blocks may be large (in the order of gigabytes), so this part of the cache will have to be limited with a maximum cache size and data blocks would be replaced with the least recently used (LRU) strategy.

3. "Derived" triangular mesh. Used to speed up and simplify rendering, as well as providing mesh for rendering in 3D map view. More details below.

4. Spatial index. For map rendering and layer identification we need fast spatial queries ("which mesh elements are in this rectangle").

Derived mesh

Reasons why we need derived triangular mesh:

• mesh structure in general may contain triangles, quadrilaterals or even polygons with higher number of vertices. Interpolation on shapes other than triangles is more expensive/difficult and for rendering in 3D we need to have a triangular mesh anyway.
• when rendering 1D elements (edges), we need to turn them into quads according to their thickness (and quads are turned into a pair of triangles).
• when map CRS is different from the mesh layer CRS, we need to do on-the-fly reprojection. Keeping derived mesh in the map CRS saves us from reprojecting on every rendering.
• level of detail - when rendering a large area, we can do simplification and use fewer triangles for rendering, for smaller area we may want to split elements such as quads into more than two triangles for better rendering quality.
  

Derived mesh is kept for particular map settings parameters (especially scale and CRS). Whenever the map settings are changed, the derived mesh needs to be invalidated. The assumption is that often these map settings parameters stay the same over multiple rendering calls.

GUI

Loading of mesh layers should be possible via the integrated browser panel and/or through the unified add layer dialog.

Once a mesh layer is loaded, it will have a map layer properties dialog similar to the ones for vector/raster layers, with renderer configuration similar to the GUI mockups above.

Supported data formats

Out of the box, we would like to support few generic file formats that are useful for wide range of uses and do not introduce any additional dependencies, e.g. GRIB format (for weather data). Support for other formats may be added by plugins.

Affected files

Most of the development will be done in new files (e.g. src/core/mesh/qgsmeshlayer*.cpp).

We will need to add a new layer type in QgsMapLayer and add support for it in various bits of QGIS code.

Q/A

Would it be possible to implement support for any custom format?

Yes, by implementation of subclass of QgsMeshDataProvider and registering it in the list of providers available. This could be done in Python plugins.

How will you handle time component?

This QEP intentionally tries to stay away from dealing with time component as that is quite a large topic itself. We assume that plugins may introduce convenience functionality for handling of the temporal nature of mesh layers (e.g. easy moving between time steps, animation, export of videos). In the future we would like to look into having a separate QEP to cover time component within QGIS for all types of map layers.

Would it be possible to visualize mesh data in 3D?

Yes, probably not from the very beginning, but it should be possible at some point.

Would it be possible to load very large files with results?

Yes, data provider does not store any data internally in memory. The mesh layer cache will gather data from provider when needed. It will be configurable (e.g. limit for memory usage) so the performance should be good even for large files.

Can I carry out Processing on my mesh layer?

We may provide some basic Processing modules for conversion of mesh layers to raster/vector layers. For advanced processing (e.g. mesh calculator, time manipulation) you can use Crayfish or develop a python plugin.

Is it the end of Crayfish?

No, at least not in the short term. Crayfish will just get rid of the C++ part (and become pure Python plugin) and serve as a time manager with some additional functionality (plots, animations). But if we manage to implement support for time management in QGIS and integrate other parts of it as well, maybe there will not be much need for it anymore...

Backwards Compatibility

Not applicable

Votes

(required)

[vpicavet]

Hi, nice QEP, very glad to see this coming into QGIS !

Could you elaborate on the limitations you found in the MeshLayer Plugin ? We would like to be sure that all cases handled by this plugin would be handled as well by native support.

Do you plan to use OpenGL rendering as well for mesh layers ? It proved to be really efficient for the plugin.

@vmora will certainly react to this. We are really interested in this area of work and have various use cases. We would for sure be able to take part in the development, should funding be assured.

[anvpires]

Big +1

[rouault]

Just wondering if the low-level part of this QEP (file format handling, spatial indexing, caching) wouldn't make sense to belong to some sort of a stand-alone MDAL (Mesh Data Abstraction Library), that would join the greater family of GDAL, PDAL, etc...

[wonder-sk]

Hi @vpicavet

MeshLayer plugin limitations to my best knowledge from reading the code (I may be wrong):

• data input - only supports triangles
  • often one needs to deal with quads, sometimes even elements with higher number of vertices
  • some meshes use 1D elements (edges)
• rendering
  • it is not possible render vector data (arrows)
  • mesh structure cannot be rendered (nodes/elements)
  • continuous rendering only supports interpolation between min/max value (not multiple colors for different values)
• data access
  • does not implement caching of data - providers either need to load all data to memory or need to implement caching individually
  • does not have spatial index - needs linear search over all triangles to answer a query "what value is at this coordinate"

Regarding OpenGL rendering for 2D maps, I am a bit skeptical. My breakdown of pros/cons would be:

Pros:

• faster than rendering on CPU

Cons:

• code is more difficult to read and maintain (writing raw OpenGL code and pixel shaders is hard!)
• only solves one piece of rendering demands (rendering of vector arrows would be more difficult with OpenGL)

If rendering of a mesh layer on CPU takes 50ms and rendering with OpenGL is done in 1ms, then I would probably argue that the disadvantages are greater than the performance benefit...

[wonder-sk]

@rouault Yes that may as well happen at some point. I would prefer to do a gradual approach with implementation directly in QGIS first and then if there is enough interest, move the low-level code and drivers to a separate library...

MDAL would be actually a great name as it starts with my initials :-D The next thing would be to figure out whether to pronounce it "em-dal" or "moodle" !

[rouault]

I would prefer to do a gradual approach with implementation directly in QGIS first and then if there is enough interest, move the low-level code and drivers to a separate library...

Sure. You'll have to consider the licensing of those parts from the beginning though. Not wanting to start a flame war about licensing, but a GPL'ed MDAL could limit its audience.

[nyalldawson]

@rouault

I guess ideally a standalone library would avoid the qt dependency too...

[cbertelli]

I think the vision of Martin is a great value for this development!
Anyway, when I read the enhancement proposal I had mixed feelings. I agree with @rouault. Being a bad guy, I may say that QGIS should not try to solve all problems by itself, that taking advantage from libraries like GDAL, PDAL, Numpy and from the continuous development of GRASS is is a good idea. Cooperation is a great advantage of Open Source.
On the other hand, @wonder-sk criticises the foundation of QGIS on GDAL simple features, which is a strong, efficient, but very constraining foundation. My point of view is database-centric, while @wonder-sk is process-centric, anyway it's a very good idea to widen the data types that QGIS can handle, but it could be an opportunity to think about the future, sharing the idea of a topology based model beside simple features, looking forward a stronger integration with OSM and other database models. Why not joining forces with GRASS in speculating a future for Open GIS?

[mdsumner]

This is the first clear manifesto against the confines of traditional GIS that I've seen in the public domain, thank you! I hope this gets traction beyond the scope of QGIS, it's a foundation needed by many. (I definitely support @rouault point about that, QGIS is awesome but it's not the whole world - it is a wide audience though, so long-term it's the right place to get support for truly general tools).

[wonder-sk]

@cbertelli

Being a bad guy, I may say that QGIS should not try to solve all problems by itself, that taking advantage from libraries like GDAL, PDAL, Numpy and from the continuous development of GRASS is is a good idea.

Absolutely agreed. However as of now (to my knowledge) there is no such open source library to support I/O, rendering and editing of this kind of data.

@wonder-sk criticises the foundation of QGIS on GDAL simple features [...]

Don't get me wrong - I have nothing against vector simple features! That representation is great for vector layers, it is just not the best representation to use for meshes.

[mdsumner]

I believe this is a meta-topic of the highest order, it is criticized in the small in different contexts, has pockets in various places and only by a large-enough unification, a truly "not-just-simple features" with enough key players will it get traction soon.

This discussion is relevant: https://lists.osgeo.org/pipermail/gdal-dev/2017-October/047464.html and there are many others. If anyone wants to pursue this topic beyond this context I'm keen ;) My small experiment is here: https://github.com/hypertidy/silicate and I'll be exploring this proposal in detail, it's the most interesting big picture topic for spatial I've seen in open source.

[cbertelli]

@wonder-sk. Unfortunately MDAL cannot be called Moodle because Moodle is another thing (see moodle.org, a well known and dependable OS e-learning project) and Em-dal is not so appealing ;.), but @rouault's idea of an abstraction library is a good idea.

[PeterPetrik]

We can still pronounce it Methylenedioxyallylamphetamine

[vmora]

@wonder-sk first of all: Great QEP!

I've been working mostly on simulation related projects for the past 4 years. I resorted to a lot of "ad hoc" solutions. Meshlayer tries to be a bit more general (it's used in several projects) but it's just a plugin, and a rather small one at that:: a thousand lines of python. Just a word about the limitations that I believe it has not::

• data input - only supports triangles: Not exactly, but if you have quads, it's the meshdataprovider role to convert them into pairs of triangles.
• continuous rendering only supports interpolation between min/max values: there is also a classified symbology

Limitations really it has: does not implement caching of data and does not have spatial index. Those limitation are structural (i.e. I'm not planing to change that) : the full mesh is sent to the graphic card, only values are fetched dynamically. Usually simulation meshes fit in memory, so not a big deal for my use case. But for a general GIS layer, this is not acceptable.

Concerning OpenGL 2D rendering:

• code is more difficult to read and maintain (writing raw OpenGL code and pixel shaders is hard!):
  well, c++ is not exactly trivial to write either, but you get used to it, same hold for shaders
• only solves one piece of rendering demands (rendering of vector arrows would be more difficult with OpenGL), yes QT is not there to help with that, but it can be done in OpenGL

If rendering of a mesh layer on CPU takes 50ms and rendering with OpenGL is done in 1ms, then I would probably argue that the disadvantages are greater than the performance benefit...

Here is the key issue for me: if the CPU can do that in less than 20ms (50Hz) fair enough, if not, then you need something else. Animation is the most natural way to understand time variation, and IMO a must-have feature for post-processing simulation results.

gradual approach with implementation directly in QGIS first

MDAL is certainly appealing. I don't have a strong opinion on that, but...

Mesh-based data cannot be properly represented either as vector or raster layers

Indeed, but it would be great to have the mesh layer being loaded as a point layer, a polygon layer and a raster layer (with infinite resolutions). The arrow representation would benefit raster too.

So MDAL could do just that: take a mesh format and expose 3 pseudo-datasources.

[wonder-sk]

@vmora Many thanks for your comments.

Regarding meshlayer plugin:

• triangle-only support: my idea is that data providers should not need to worry about conversion to triangles - conversion to triangles should be handled by some intermediate layer as the use of triangles is an internal implementation detail. It is also more complicated than just splitting a quad into two triangles: there should be bilinear interpolation within a quad which can be approximated by sub-diving quads using extra vertices. This issue is more pronounced in elements with more vertices (e.g. hexagons). 1D elements may be represented as two triangles, but their vertex positions should be adjusted depending on map scale to look like they have constant width.
• rendering with classified symbology: if I understand correctly, you either use "classified" symbology and assign intervals to a discrete color, but currently you can't define that 1=red, 2=white, 4=blue and values in between are interpolated (e.g. 1.5 would be light red). I know it is not too difficult to add though...

In terms of memory consumption, my assumption in the QEP is that mesh structure always fits into memory, but datasets do not have to.

Regarding OpenGL rendering:

• I agree with what you say that c++ is not trivial either, on the other hand, vast majority of QGIS code is c++, so there is an expectation that devs already know it, while OpenGL programming is quite niche (and debugging shaders is a pain). Just recently I had a chat with some GRASS devs about Nviz tool - nobody really wants to touch it because the OpenGL code is simply a different world...
• rendering speed: another issue is that people rarely get just plain mesh layer rendered - there are other background layers, effects, labeling, compositing. So even if you use render caching, the other bits of the rendering pipeline may be too slow for a smooth 30/60 fps animation. My take on the animation would be to first prepare full animation frames of the map (with larger delays between frames), cache them and then replay them with desired fast animation. When we add time/animation support to QGIS core (including animation of vectors and rasters), then we will have no other option anyway (unless we would turn also vector/raster rendering into OpenGL).

Loading mesh layer as point/polygon/raster layer: for this I would like to use Processing toolbox, where users could take a mesh layers as input and convert them to whatever representation they need (in fact, Crayfish already supports such algorithm: export nodes, export elements, export data as raster, export contour isolines).

[vmora]

the other bits of the rendering pipeline may be too slow for a smooth 30/60 fps animation

It works for MeshLayer, the other layers are cached (as long as you don't want pan/zoom + animate).

unless we would turn also vector/raster rendering into OpenGL

It wouldn't help if the bottleneck is data fetch, which occurs with every frame.

Loading mesh layer as point/polygon/raster ...

Well, I'm not sure that the extend of what I was thinking was properly expressed in my short sentence.

What the mesh represents is a continuous scalar or vector field (node values + elements interpolation), and I think raster really do the job to represent fields. The only real difference is that you can zoom indefinitely on the mesh representation, but bottom line you are not interested in a resolution that is higher than your viewport (or a multiple of that for eg. antialiasing).

Classified colors is something we already have for rasters. On the other hand, the arrow field representation would be really nice to have for raster layers (not just mesh layers), same holds for nicer classified and continuous colors.

The basic idea is to implement a rasterizer in MDAL (that has more or less the same API as GDAL except for the max resolution). That rasterizer gives you the interpolated field values (one band for scalar, 2/3 for vector fields) at the desired resolution. Then the only thing that need improvements in QGIS is the raster symbology.

So what I propose, with 3 views on the same datasource doesn't involve brand new concepts in GIS. Of course we won't benefit from the new concepts, eg. processing won't see the "graph" aspect of the mesh (a 4th view on the same thing)... but do we truly need the new concepts ?

[wonder-sk]

@vmora Thanks for elaborating on that - now I recall someone from Oslandia mentioned such idea to me longer time ago (either Hugo or VincentP or you?) and they said they would like to have mesh rasterizer driver in GDAL to be able to reuse QGIS raster renderers. Did you pursue that idea any further or were there some obstacles that were difficult to overcome?

Having such mesh rasterizer in QGIS would be a good solution to avoid duplication of some rendering code. And having vector field renderer shared with rasters would be nice too (so that raster grid datasets would not need to be opened as mesh layers).

I will think about this approach a bit more to see if there are any complications coming from it. At the same time, having the possibility to directly access the mesh structure (graph) and raw datasets could be still useful (e.g. for mesh calculator or for mesh structure/data editing).

[rduivenvoorde]

@wonder-sk Thanks for this, it is (to me) very interesting and I learn from your thorough write ups.

I sometimes get grib files, and either open them with QGIS/GDAL, or your Crayfish plugin. I wonder if I understand correctly that the first QgsMeshDataProvider implementation will probably be GDAL, followed by a 'native' QGIS one, being your work in Crayfish?

As grib and other meteo formats often contain TimeRanges, I wonder if the QgsMeshLayer or QgsMeshDataProvider will have machinery to handle those. Actually I would favour native handling of time-ranges in all layers/providers :-)

[haubourg]

@rduivenvoorde side note: I'm still on QGIS server OGC certification and we don't satisfy WMS time support. Having time support there would be cool (isn't it @anitagraser ?)

[wonder-sk]

@rduivenvoorde The idea is that GDAL-based data provider would be included to support GRIB and possibly other data sources. Crayfish would introduce another provider (or a whole bunch of providers - one for each format).

Regarding time steps vs time ranges, at this point the mesh layers / data providers do not really care about that - they work with a bunch of datasets that may or may not contain some time information in metadata. I hope one day there will be a nice QEP for time dimension handling and that's where time support should get formalized.

[vmora]

@wonder-sk I've pieces of informations about raster animation by changing band (netcd from wrf weather forecast) and rasterizing mesh fields with opengl and I whanted to combine that with the "python vrt" support that is available in gdal thanks to @rouault. But I did not find the time yet.

[anitagraser]

@haubourg sure WMS-T support in QGIS Server would be great! I'm using GeoServer for WMS-T now.

[vpicavet]

I will think about this approach a bit more to see if there are any complications coming from it. At the same time, having the possibility to directly access the mesh structure (graph) and raw datasets could be still useful (e.g. for mesh calculator or for mesh structure/data editing).

As said by Vincent the "MDAL" library could offer a raster view of the Mesh for display ( with OpenGL process to render as raster), as well as vector views for points ( middle of cells or nodes) and edges. This would allow for mesh calculator or mesh structure editing too.

What may be needed is a graph structure support in QGIS : this would also be useful to support topology editing, like the PostGIS one. But it could be added later as another view on the mesh layer.

Did you find time to think about such an approach ?

[wonder-sk]

Just earlier today we have had a discussion and decided to go ahead with MDAL library. The initial goals are to provide data structures for representation of mesh data and support for drivers - following the model of GDAL/OGR. It is meant to be a C++ library, MIT licensed.

In the first stage I would consider MDAL to be just a simple (read-only) data access library, without caching, indexing or rasterization support - all those things would be initially handled inside QGIS. Using Vincent's suggestion, for scalar values we would rasterize mesh data and then use ordinary raster pseudo-color renderer. For vector values we want to create a generic vector field renderer that may be reused for vector, raster or mesh data.

Once the basics are implemented inside MDAL, we can start looking into further topics:

• providing raster view (with OpenGL or without it)
• providing vector views (vertices, edges, faces)
• caching, spatial indexing
• write support (i.e. editing of mesh structure and/or associated data)
• mesh calculator
• Python API
• command line tools (e.g. ogr2ogr for mesh data)

@PeterPetrik has just created MDAL repository https://github.com/lutraconsulting/MDAL

[vmora]

@rouault MDAL should be compatible with gdal to enjoy gdal raster goodies, it'd be nice to have compatible raster and vector API, wouldn't it ?

And if we go for an MIT licence.... why not make MDAL part of gdal ?

[rouault]

@@> why not make MDAL part of gdal ?

On a pure technical point of view, there's the QT dependency question ( lutraconsulting/MDAL#1 ) that would be a no-go for GDAL (or that would make the MDAL part necessarily an optionally compiled component of GDAL, which would in particular mean that a lot of packaging wouldn't compile MDAL functionality)

This issue left aside, that might make sense. We have already added recently the GNM (Geographical Network Model) abstraction in GDAL, so why not the mesh one ?

There are pros & cons to being integrated or standalone:

• pros:
  • re-use of all GDAL "infrastructure" (CI, website, etc..)
  • existing packaging
  • if functional bridges between MDAL and GDAL, then easy to use GDAL ones. And if GDAL is a required dependency of MDAL, that would make even more sense
  • GDAL being a required dependency of QGIS, then MDAL would automatically be available.
• cons:
  • perhaps a feeling of being less free for creativity due to GDAL being perceived as "heavy weight" (but as QGIS was envisionned first as being the host, that wouldn't really hold as an argument),
  • GDAL's build system can be a bit painful (not sexy GNUMakefile and nmake makefile.vc)
  • Perhaps GDAL's yearly schedule for feature versions would be too slow ?

[vpicavet]

@rouault @wonder-sk It would be better indeed to drop QT dependency, and incorporate this library inside the GDAL framework. Pros are definitely highers than cons.

And you forgot the best "pro" argument : if MDAL is part of GDAL, we would favor @rouault 's contributions, and that would make a real difference :-)

[ahuarte47]

Hi, nice QEP, a big +1 for MDAL!

What about to develop the Mesh basic struct with a streamable style?

If the Mesh provides a stream or iterator of triangles or faces, instead of arrays of vertices or faces, MDAL could provide a huge set of features to render. Something similar to QgsVectorFeatureIterator for normal vector layers.

QGIS could draw huge sets of TIN's from LiDAR repositories, etc. A similar concept to this.

[rapidlasso]

Thanks for the pointer to our "streaming mesh" work. I also like to advocate against - yet another - two array based mesh representation but suggest to consider one that exists of a stream of 3 events: Vertex, Polygon, Finalization. A more format specific reference to such a "Streaming Mesh" is this one here. It exlains the concept of storing (or exposing an API to) a mesh that allows true streaming access. Such an access has two components: (1) vertices and triangles/polygons can appear interleaved in the sense that vertices just need to be loaded/introduced when the the first triangle/polygon references them and not all of them upfront (such as in many other indexed formats) and (2) there is "finalization" information that tells us when vertices are no longer needed (and can be safely deallocated). Every existing indexed format implicitly fits this model and simply has the worst-case "stream width" which refers to the maximum number of vertices that are active (introduced and not finalized) as the mesh streams by. For the standard indexed mesh format such as PLY all vertices are loaded/introduced first and are finalized only after the last triangle/polygon was read. More details (abstract and paper) are here: Streaming Meshes. Any actual implementation will require only tiny API changes compared to a non-streaming approach. But the practical benefits can be huge as meshes grow in size.

[PeterPetrik]

@ahuarte47 @rapidlasso, thank you for reference to streaming mesh. The support of large datasets and meshes is definitely the goal. Going to look into it, but is there anywhere public header of streaming meshes library available?

[rapidlasso]

For your purpose I could imagine that operating in increments of single vertices and triangles (as done in the sample Streaming Mesh format used by LAStools BLAST engine) might be too granular. But here an as small as possible example. This is our SMA (Streaming Mesh Ascii) format which is also a valid OBJ file. This is a simple tetrahedron.

# nverts 4
# nfaces 4
#
v -1 -1 -1
v 1 1 -1
v 1 -1 1
f 2 1 3
v -1 1 1
f 4 3 1
f -4 2 4
f -3 -2 -1 

The interleaving is easily seen. After the first three vertices we already can specify the first face. Then comes the fourths and last vertex. Hence now the next three (here the final three) faces can be specified as all needed vertices have been introduced in the stream.

Now notice the negative reference '-4' in the third face. This is the last time that the first vertex is referenced so we finalize it. In this moment (once the triangle was rendered or processed or etc) we could potentially already deallocate it's memory. That is what allows streaming. The knowledge when we can destroy things. Like in a video stream. Once we decompressed and displayed a frame and completed using it for predictive decompression of future frames we deallocate it. Otherwise the whole "Finding Nemo" movie would be in memory at the end of watching it ... (-:

We use negative indices to "hide" finalization tags because it's actually supported in the OBJ format specification. Those are rarely used in practice it seems but a complete OBJ reader should implement them. They are relative indices that reference the vertex in the array at the position computed from the total number of introduced vertices plus (hence minus) this (negative) value. So when -4 appears the vertex count is at 4 and 4 - 4 = 0 so that vertex at position 0 in the array gets finalized. Using a hash instead of an array to map from indices to vertices allows to deallocate finalized vertices.

The last face finalizes the other three vertices explicitly, although the end of the file would do that anyhow.

Your API could implement introduction of vertices and faces in batches rather than individuals. For small meshes where streaming does not matter those batches are simply all vertices first, then all triangles, and nothing changes, cause initially most users won't have big meshes anyways. But integrating the concepts from the start - some variety of the above - can really turn out to be a slick move in the long run.

I could share a larger example if needed.

[mustafaucr]

First of all I really congratulate you this is very big development! I am working on a plugin hydrology, hydrodynamic modelling in QGIS, as well. I wonder that. How can we prepare a mesh layer? What are the inputs for that layer? As I can see there is no attribute table. QGIS only show the layer data, but I dont know how I can create a mesh layer.

[saberraz]

You need to add your format to mdal:
https://github.com/lutraconsulting/MDAL

There are already examples of hydrodynamic model outputs, for example 3Di:
https://github.com/lutraconsulting/MDAL/tree/master/mdal/frmts

Currently, a copy of medal sits within QGIS source:
https://github.com/qgis/QGIS/tree/master/external/mdal

Ping us if you need further information or help.

[mustafaucr]

Thank for quick response @saberraz . However I use python for developping is there any python examples? or else I can't find?

[saberraz]

Once the format is available from mdal, you should be able to interact with it similar to other layers in Python. See the latest Crayfish plugin for example.

[PeterPetrik]

http://www.mdal.xyz

[mustafaucr]

I stil can not unserstand how I can create a mesh layer. For instance creating a vector layer is described in here:https://docs.qgis.org/testing/en/docs/pyqgis_developer_cookbook/vector.html#add-features. I want to an example like that. For instance I want to learn how time series initialize and features add?

[saberraz]

Once your format is supported in MDAL, you can add a layer:

QgsMeshLayer (“path/to/file.grb”, “my layer”, “mdal”)

Here is more from QGIS documentation.
https://qgis.org/api/classQgsMeshLayer.html

[mustafaucr]

Well it is okey. However it uses *.grb file. How can I create like a file I cannot see the attribute table like shape file because when I want to create a shape file I add the fields whatever I want because I know which fields I need. Four instance a mesh layer use time series, wind information, temperature information etc. How can I initialize these fields to mash layer.

[wonder-sk]

@mustafaucr what file format would you like to use? Currently QGIS supports read-only access to existing mesh data - you cannot create new files.

Also please note that mesh layers do not have attribute table like vector layers do.

[mustafaucr]

Well it is not possible to create a mesh layer file in QGIS. What about python? It is very new so I guess there is no library for creating any mesh layer file formats.

[wonder-sk]

You can create input data for mesh layer even in a text editor, for example:

• mesh structure defined in .2dm files: https://www.xmswiki.com/wiki/SMS:2D_Mesh_Files_*.2dm
• mesh data defined in .dat files: https://www.xmswiki.com/wiki/SMS:ASCII_Dataset_Files_*.dat

[JosefWN]

Seems awesome! Will it support NetCDF on OS X too? It doesn't seem to right now, using a clean install of QGis 3.4 and the mesh layer (or the raser layer I should add). I'm getting: "Raster layer provided is not valid" and "file.nc is not a valid or recognized data source". EDIT: File is valid though, I can open it just fine in Panoply for example.

Also, in the "Add layer" menu, the mesh option does not seem to be available, only going via the data source manager?

[PeterPetrik]

Hi @Puffton, NetCDF files are supported, but maybe you have some file that has unrecognized syntax. Which installer do you use (https://lutraconsulting.github.io/qgis-mac-packager/ or QGIS official one?)
Maybe best to ask these questions on qgis-developer list, since this page is dedicated to enhancement proposals. Thanks

[JosefWN]

I used the official QGIS installer.