Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

write_mdim() stars object to view in QGIS as mesh layer #704

Open
jvandens opened this issue Aug 24, 2024 · 2 comments
Open

write_mdim() stars object to view in QGIS as mesh layer #704

jvandens opened this issue Aug 24, 2024 · 2 comments

Comments

@jvandens
Copy link

Consider the star object below, is there a way to save this using write_mdim() or other method that is compatible with viewing in QGIS as a Mesh Layer? https://docs.qgis.org/3.34/en/docs/user_manual/working_with_mesh/mesh.html#

For reference, this star object contains a polygon geometry stacked in z (depth) dimension and time dimension with 2 attributes (dataset groups in Mesh layer terminology)

image

star = structure(list(`F Coli` = structure(c(1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1), .Dim = c(geometry = 5L, 
z = 10L, times = 8L)), Entero = structure(c(1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1), .Dim = c(geometry = 5L, 
z = 10L, times = 8L))), dimensions = structure(list(geometry = structure(list(
    from = 1L, to = 5L, offset = NA_real_, delta = NA_real_, 
    refsys = structure(list(input = "NAD83 / UTM zone 18N", wkt = "PROJCRS[\"NAD83 / UTM zone 18N\",\n    BASEGEOGCRS[\"NAD83\",\n        DATUM[\"North American Datum 1983\",\n            ELLIPSOID[\"GRS 1980\",6378137,298.257222101,\n                LENGTHUNIT[\"metre\",1]]],\n        PRIMEM[\"Greenwich\",0,\n            ANGLEUNIT[\"degree\",0.0174532925199433]],\n        ID[\"EPSG\",4269]],\n    CONVERSION[\"UTM zone 18N\",\n        METHOD[\"Transverse Mercator\",\n            ID[\"EPSG\",9807]],\n        PARAMETER[\"Latitude of natural origin\",0,\n            ANGLEUNIT[\"Degree\",0.0174532925199433],\n            ID[\"EPSG\",8801]],\n        PARAMETER[\"Longitude of natural origin\",-75,\n            ANGLEUNIT[\"Degree\",0.0174532925199433],\n            ID[\"EPSG\",8802]],\n        PARAMETER[\"Scale factor at natural origin\",0.9996,\n            SCALEUNIT[\"unity\",1],\n            ID[\"EPSG\",8805]],\n        PARAMETER[\"False easting\",500000,\n            LENGTHUNIT[\"metre\",1],\n            ID[\"EPSG\",8806]],\n        PARAMETER[\"False northing\",0,\n            LENGTHUNIT[\"metre\",1],\n            ID[\"EPSG\",8807]]],\n    CS[Cartesian,2],\n        AXIS[\"(E)\",east,\n            ORDER[1],\n            LENGTHUNIT[\"metre\",1]],\n        AXIS[\"(N)\",north,\n            ORDER[2],\n            LENGTHUNIT[\"metre\",1]],\n    ID[\"EPSG\",26918]]"), class = "crs"), 
    point = FALSE, values = structure(list(structure(list(structure(c(498756.910095215, 
    510957.220092773, 522344.630126953, 512631.440124512, 498756.910095215, 
    4285731.00012207, 4304042.50012207, 4298020.00012207, 4279347.50012207, 
    4285731.00012207), .Dim = c(5L, 2L))), class = c("XY", "POLYGON", 
    "sfg")), structure(list(structure(c(512631.440124512, 522344.630126953, 
    533818.380126953, 525469.630126953, 512631.440124512, 4279347.50012207, 
    4298020.00012207, 4292652.50012207, 4273846.50012207, 4279347.50012207
    ), .Dim = c(5L, 2L))), class = c("XY", "POLYGON", "sfg")), 
        structure(list(structure(c(525469.630126953, 533818.380126953, 
        545675.130126953, 538335.190124512, 525469.630126953, 
        4273846.50012207, 4292652.50012207, 4287713.50012207, 
        4268440.50012207, 4273846.50012207), .Dim = c(5L, 2L))), class = c("XY", 
        "POLYGON", "sfg")), structure(list(structure(c(538335.190124512, 
        545675.130126953, 558535.25012207, 552430.810119629, 
        538335.190124512, 4268440.50012207, 4287713.50012207, 
        4283077.50012207, 4262513.00012207, 4268440.50012207), .Dim = c(5L, 
        2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
            structure(c(552430.810119629, 558535.25012207, 572934.75012207, 
            568734.440124512, 552430.810119629, 4262513.00012207, 
            4283077.50012207, 4279116.00012207, 4256515.00012207, 
            4262513.00012207), .Dim = c(5L, 2L))), class = c("XY", 
        "POLYGON", "sfg"))), class = c("sfc_POLYGON", "sfc"), precision = 0, bbox = structure(c(xmin = 498756.910095215, 
    ymin = 4256515.00012207, xmax = 572934.75012207, ymax = 4304042.50012207
    ), class = "bbox"), crs = structure(list(input = "NAD83 / UTM zone 18N", 
        wkt = "PROJCRS[\"NAD83 / UTM zone 18N\",\n    BASEGEOGCRS[\"NAD83\",\n        DATUM[\"North American Datum 1983\",\n            ELLIPSOID[\"GRS 1980\",6378137,298.257222101,\n                LENGTHUNIT[\"metre\",1]]],\n        PRIMEM[\"Greenwich\",0,\n            ANGLEUNIT[\"degree\",0.0174532925199433]],\n        ID[\"EPSG\",4269]],\n    CONVERSION[\"UTM zone 18N\",\n        METHOD[\"Transverse Mercator\",\n            ID[\"EPSG\",9807]],\n        PARAMETER[\"Latitude of natural origin\",0,\n            ANGLEUNIT[\"Degree\",0.0174532925199433],\n            ID[\"EPSG\",8801]],\n        PARAMETER[\"Longitude of natural origin\",-75,\n            ANGLEUNIT[\"Degree\",0.0174532925199433],\n            ID[\"EPSG\",8802]],\n        PARAMETER[\"Scale factor at natural origin\",0.9996,\n            SCALEUNIT[\"unity\",1],\n            ID[\"EPSG\",8805]],\n        PARAMETER[\"False easting\",500000,\n            LENGTHUNIT[\"metre\",1],\n            ID[\"EPSG\",8806]],\n        PARAMETER[\"False northing\",0,\n            LENGTHUNIT[\"metre\",1],\n            ID[\"EPSG\",8807]]],\n    CS[Cartesian,2],\n        AXIS[\"(E)\",east,\n            ORDER[1],\n            LENGTHUNIT[\"metre\",1]],\n        AXIS[\"(N)\",north,\n            ORDER[2],\n            LENGTHUNIT[\"metre\",1]],\n    ID[\"EPSG\",26918]]"), class = "crs"), n_empty = 0L)), class = "dimension"), 
    z = structure(list(from = 1L, to = 10L, offset = 1, delta = 1, 
        refsys = NA_character_, point = FALSE, values = NULL), class = "dimension"), 
    times = structure(list(from = 1L, to = 8L, offset = NA_real_, 
        delta = NA_real_, refsys = NA_character_, point = FALSE, 
        values = structure(list(start = c(0, 0.020833333954215, 
        0.0625, 0.10416666418314, 0.145833343267441, 0.1875, 
        0.22916667163372, 0.270833313465118), end = c(0.020833333954215, 
        0.0625, 0.10416666418314, 0.145833343267441, 0.1875, 
        0.22916667163372, 0.270833313465118, 0.3125)), class = "intervals")), class = "dimension")), raster = structure(list(
    affine = c(0, 0), dimensions = c(NA_character_, NA_character_
    ), curvilinear = FALSE, blocksizes = NULL), class = "stars_raster"), class = "dimensions"), class = "stars")
@edzer
Copy link
Member

edzer commented Aug 24, 2024

write_mdim() will write out the geometries to e.g. a netCDF file with discrete geometries (polygons); I don't know how the mesh stuff in qgis handles those, but it seems they could be quads. Any ideas or suggestions, @mdsumner ?

@jvandens
Copy link
Author

reprex_model.nc.zip

Thanks @edzer ,if it helps @mdsumner , attached is the above star saved with write_mdim(). When loaded into QGIS, it seems to recognize some of the data that is in there, but is confused by the geometry and other dimensions (see first 2 screenshots)

For reference, the mesh are polygon quads, see last screenshot for what they should look like.

There should be 5 faces, 12 vertices repeated over 10 depths/elevations and 8 timesteps. With 2 datasets (entero and Fcoli) at each face.

image

image

image

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Assignees
No one assigned
Labels
None yet
Projects
None yet
Milestone
No milestone
Development

No branches or pull requests
2 participants
@edzer @jvandens