lutraconsulting / MDAL

Mesh Data Abstraction Library
http://www.mdal.xyz/
MIT License
160 stars 50 forks source link

Render U and V face based scalar variables as vector in QGIS 3.18 Mesh layer on an unstructured triangular mesh. #355

Closed james-morrison-mowi closed 3 years ago

james-morrison-mowi commented 3 years ago

I'm using QGIS 3.18.0 to read in data from NetCDF files from a custom hydrodynamic model which uses a triangular unstructured grid.

I do postprocessing on the NetCDF files to ensure the mesh file can be viewed in QGIS. The mesh displays correctly and the U and V variables appear as Scalar variables which are applied to the cell centers (or face_nodes) of the mesh. QGIS doesn't seem to be able to parse them as a vector in the properties of the mesh.

When I look at gridded NetCDF files from other sources I see u and v being parsed as a vector, Are U and V variables on faces able to be viewed as vectors?

Do I have to apply them to the face_x, face_y coordinates? I am using xarray to postprocess the NetCDF file and based the initial process on this example which doesn't include a vector variable.

I have pasted the output of ncdump below from one output.

netcdf main_grid20 {
dimensions:
        time = 1897 ;
        node = 35986 ;
        face = 63301 ;
        nmax_face = 3 ;
        edge = 99384 ;
        nmax_edge = 2 ;
variables:
        float time(time) ;
                time:_FillValue = NaNf ;
        int64 mesh2d ;
                mesh2d:cf_role = "mesh_topology" ;
                mesh2d:long_name = "Topology data of 2D mesh" ;
                mesh2d:topology_dimension = 2LL ;
                mesh2d:node_coordinates = "node_x node_y" ;
                mesh2d:face_node_connectivity = "face_nodes" ;
                mesh2d:edge_node_connectivity = "edge_nodes" ;
        double node_x(node) ;
                node_x:_FillValue = NaN ;
        double node_y(node) ;
                node_y:_FillValue = NaN ;
        double face_x(face) ;
                face_x:_FillValue = NaN ;
        double face_y(face) ;
                face_y:_FillValue = NaN ;
        double face_nodes(face, nmax_face) ;
                face_nodes:_FillValue = -1. ;
                face_nodes:cf_role = "face_node_connectivity" ;
                face_nodes:long_name = "Vertex nodes of mesh faces (counterclockwise)" ;
                face_nodes:start_index = 1LL ;
                face_nodes:coordinates = "face_x face_y" ;
        double edge_x(edge) ;
                edge_x:_FillValue = NaN ;
        double edge_y(edge) ;
                edge_y:_FillValue = NaN ;
        int64 edge_nodes(edge, nmax_edge) ;
                edge_nodes:_FillValue = -1LL ;
                edge_nodes:cf_role = "edge_node_connectivity" ;
                edge_nodes:long_name = "Start and end nodes of mesh edges" ;
                edge_nodes:start_index = 1LL ;
                edge_nodes:coordinates = "edge_y edge_x" ;
        double u10(time, face) ;
                u10:_FillValue = -32767. ;
                u10:long_name = "10 metre U wind component" ;
                u10:units = "meters s-1" ;
                u10:missing_value = -32767LL ;
                u10:coordinates = "face_x face_y" ;
        double v10(time, face) ;
                v10:_FillValue = -32767. ;
                v10:long_name = "10 metre V wind component" ;
                v10:units = "meters s-1" ;
                v10:missing_value = -32767LL ;
                v10:coordinates = "face_x face_y" ;
        int64 data_on_edges(edge) ;
                data_on_edges:mesh = "mesh2d" ;
                data_on_edges:coordinates = "edge_y edge_x" ;

// global attributes:
                :Conventions = "CF-1.8 UGRID-1.0" ;
                :coordinates = "node_x node_y" ;
PeterPetrik commented 3 years ago

Hi,

the logic for UGRID files is in the https://github.com/lutraconsulting/MDAL/blob/0afc30c187a67ebee7cf89c7ef02d3fdc722f53e/mdal/frmts/mdal_ugrid.cpp#L551 file on lines 551-588.

you need to change long name to one of the supported variants, e.g. u component of 10 metre wind or x-component of 10 metre wind.

Kind regards, Peter

james-morrison-mowi commented 3 years ago

Thanks, I can get it to appear now as a Vector Type within the Source tab and it is possible to display the Vector data as a contour but the Show Vectors tab is still greyed out.

Update I see the Vector tab can be enabled by clicking the arrow on the Sources tab, but the vectors still do not appear.

I suspect it relates to the variable being assigned to the face as opposed directly to face_x, face_y?

I'm not sure how to restructure the dataset to represent vectors visually. The ncdump output for the related variables are shown below:

        double u(time, face) ;
                u:_FillValue = NaN ;
                u:mesh = "mesh2d" ;
                u:long_name = "u component of velocity" ;
                u:units = "meters s-1" ;
                u:coordinates = "face_x face_y" ;
        double v(time, face) ;
                v:_FillValue = NaN ;
                v:mesh = "mesh2d" ;
                v:long_name = "v component of velocity" ;
                v:units = "meters s-1" ;
                v:coordinates = "face_x face_y" ;
PeterPetrik commented 3 years ago

it should work with face vectors too. this looks like you need to tweak the vector arrow settings. try fixed length option or so.

james-morrison-mowi commented 3 years ago

I tried fixed length, I also see an issue with the Contours if I untick and then re-tick show contours it won't redisplay the contours (the mesh itself displays ok).

The mesh itself is in British National Grid coordinates although I haven't written the crs to the NetCDF file.

If I go to the Symbology tab and click the arrow icon next to velocity the arrows will appear but if I click on either the contour or arrow symbol after greying it the layer won't reappear (the mesh does). This happens on Windows and the nightly build on Linux from a few days ago.

PeterPetrik commented 3 years ago

best to create gif and upload here. Also maybe consider buying a SLA with us (lutraconsulting.co.uk/support) and we can do screensharing and fix your bugs as part of SLA?