Support for unstructured grid with a vertical dimension

[dennissergeev]

✨ Feature Request

Firstly, thanks very much for creating this library, I'm excited to use it in my research!

Secondly, would it be possible to generalise geovista.Transform.from_unstructured() for it to take in arrays with more than 2 dimensions, i.e. vertical levels (and possibly time)?

Motivation

I am trying to visualise LFRic output (happy to eventually contribute to the gallery by the way), and while the from_unstructured() method works great for 2D arrays, I would like to plot something w.r.t. model height. Currently there's no obvious way to do this as that function accepts only longitudes and latitudes, throwing an error if I pass a full 3D array as data.

This is probably related to the "isosurfaces support" to-do item.

[bjlittle]

Hey @dennissergeev,

Stoked to hear that you're excited about geovista... I am too!

Prior to geovista, I was playing with animating model time-series and also providing a slider so the user can manually control the render of spatial data with another single dimension (such as time). To be honest, that worked reasonably well. Here's a proof-of-concept I'm currently working on for geovista... it's showcasing an unstructured ocean model with variable resolution quad-faces:

Apologies for the poor quality GIF, it doesn't really do the colormap justice... but I think you get the idea of the direction I'm headed in with that kind of capability.

So time-series support is definitely in the pipeline, perhaps along with ipywidget integration for jupyter... but the overriding goal here is for geovista to easily integrate into the existing scientific python ecosystem e.g., scripting, notebooks, https://panel.holoviz.org/ et al... and this should be a level playing field for the community, with regards to whether you're coming to geovista from iris, xarray, cupy, numpy or just where ever i.e., geovista is specifically agnostic of iris.

There are some fundamentals that I'm keen to bank first for geovista such as CRS, planar projection and texture mapping support, but later isosurfaces, isolines, streamplots, winds visualization, mesh warping for ocean bathymetry and surface altitudes, web services tiling etc are all in play 👍

With regards to the vertical I've not quite thought that far ahead to be honest (geovista is just born, so I'm focusing on a solid core first and building outwards). Nevertheless, what kinda viz are you particularly interested in?

I was mulling over geodesic ribbon trajectory slices through unstructured meshes to show vertical levels in the planar, but I've not ventured too far into this space other than daydreaming about it... however I am motivated by concrete use cases though 😉

So if you can be more specific and give me some examples that would really help 👍

[dennissergeev]

Thanks for the speedy response @bjlittle!

Here's a proof-of-concept I'm currently working on for geovista...

Looks awesome! I would love to use this once you're happy with the prototype.

the overriding goal here is for geovista to easily integrate into the existing scientific python ecosystem

I agree that's the way to go. I'm mostly using iris because I mostly use the UM in my research (and now starting to use LFRic).

There are some fundamentals that I'm keen to bank first for geovista

Fair enough. I can try to help with some "good first-time issues"...

Nevertheless, what kinda viz are you particularly interested in?

I can of course just plot some slices of the vertical coordinate, but eventually I would like to create plots with the vertical dimension, i.e. planet's radius, such as vertical cross-sections, isosurfaces, wind vectors. As an example, here's what I made with PyVista a while ago: https://dennissergeev.github.io/exoconvection-apj-2020/

If this kind of visualisation can be combined with the time animation as in your gif above, that would be even more awesome.

[bjlittle]

@dennissergeev Awesome, I've come across that visualization before, super cool, congrats!

Okay, I can see what you're aiming for, and that all makes sense, thanks. I'm guessing your Trappist-1e 3D viz could also be projected e.g., say the sphere to Plate Carree, but still have the 3D vertical convection layers on top, thus taking advantage of a 2D planar but viewed in 3D, if you get my meaning 🤔

I think that your original 3D viz is perhaps currently possible in ParaView. I know that @tinyendian (Wolfgang Hayek @ NIWA) might be able to comment or point you to what's possible there too - if you're interested?

[dennissergeev]

Thanks @bjlittle!

I'm guessing your Trappist-1e 3D viz could also be projected e.g., say the sphere to Plate Carree,

Yes, certainly. Visualising it both in a planar projection or as a 3D sphere would be great! I guess that's why you're envisaging CRS support - to have an API similar to cartopy?

I think that you're original 3D viz is perhaps currently possible in ParaView.

That would be useful to know, so any tips from Wolfgang would be appreciated.
At the same time, this is not urgent and I'm happy to wait/contribute for geovista to mature enough to support this.

P.S. Great logo by the way!

[tinyendian]

Hi @dennissergeev and @bjlittle, visualisations such as Trappist-1e (which looks great btw!!) can definitely be done with ParaView, including time-dependence of the model data, spherical projection, moving cameras, interactive web publication, ...

I haven't used PyVista before, but, looking at the examples, it seems to use a similar philosophy as ParaView, providing a simplifying layer on top of VTK. ParaView comes with a few perks, such as a powerful GUI, and a client-server mode, which is useful for very large models that require an HPC backend to handle. Its Python scripting capabilities are comprehensive and integrate nicely with the GUI, but the price to pay is that the scripts look somewhat obscure compared to PyVista, the latter seems a bit more Pythonic.

Very happy to help, if you'd like to give ParaView a try (it can read CF-netCDF and LFRic output), or if you have any questions.

Btw, I worked in the exoplanet group in Exeter 10 years ago as a postdoc, great to see the work that you guys are doing there! And say hi to Nathan for me please 🙂.

[dennissergeev]

hi @tinyendian,

Thanks! Yes, the exoplanet group is growing stronger, so expect more cool 3D visualisations of exoplanet simulations 😄!

I actually used ParaView a little bit in the past, though only as a GUI. Now I prefer using PyVista because of its pythonic API and integration with Jupyter, but in the meantime if I wanted to give ParaView another go and visualise LFRic output, where should I start?

[tinyendian]

Hi @dennissergeev, the easiest way to start with ParaView is a conda installation of the LFRic reader plugin, which will install the latest ParaView GUI as a dependency,

conda install -c conda-forge lfric_reader

If you have a ParaView build on your system already (along with the required build tools and dependencies), you can also easily build the reader yourself, see the instructions on https://github.com/niwa/lfric_reader.

Note that Jupyter integration was recently added to ParaView, too, https://blog.kitware.com/paraview-jupyter-notebook/, but I haven't tried this out yet.

With regards to Python scripting, I agree that PyVista coding looks a lot nicer and straightforward. ParaView really has its strengths when it comes to interactive visualisation (have a look at, e.g., linked cameras for comparing two models), and handling very large datasets (hundreds of millions of cells) efficiently. Maybe there will be ways to integrate Iris/PyVista and ParaView a bit in the future, e.g., via automatically generated ParaView scripts (ParaView can be fully controlled by Python scripting).

Looking forward to seeing more exoplanet simulations 😍!

[dennissergeev]

Thanks very much @tinyendian! I'll give it a go.

[github-actions]

In order to maintain a backlog of relevant issues, we automatically label them as stale after 180 days of inactivity.

If this issue is still important to you, then please comment on this issue and the stale label will be removed.

Otherwise this issue will be automatically closed in 28 days time.

[github-actions]

This stale issue has been automatically closed due to a lack of community activity.

If you still care about this issue, then please either:

• Re-open this issue, if you have sufficient permissions, or
• Add a comment pinging @bjlittle who will re-open on your behalf.

[bjlittle]

Hey @dennissergeev,

Just wanted to let you know that I want to get my ducks-in-a-row to schedule time for supporting this issue. The first step in this journey will be to extend the capability of the geovista.bridge to ingest mesh data with a vertical dimension.

I'm guessing that you've got a tonne of archived data that you'd be keen for geovista to visualise out-the-box... is it possible for you to share it with me somehow?

Cheers 😄

[dennissergeev]

Hey @bjlittle!

This is exciting to hear! Happy to help test / develop things, just ping me.
And I'm happy to supply some sample data, of course. Here's an output file from one of my recent simulations:
https://we.tl/t-gt1W6B4mkl (about 100 Mb)

The mesh is C48 and it has 38 vertical levels. Currently, the output netCDF file does not store the level height values (only level numbers), but in this simulation I used the um_L38_29t_9s_40km level set with the following eta spacing:

     (/                                                                                                                                 &   
        0.0000000_r_def,  0.0005095_r_def,  0.0020380_r_def,  0.0045854_r_def, &
        0.0081519_r_def,  0.0127373_r_def,  0.0183417_r_def,  0.0249651_r_def, &
        0.0326074_r_def,  0.0412688_r_def,  0.0509491_r_def,  0.0616485_r_def, &
        0.0733668_r_def,  0.0861040_r_def,  0.0998603_r_def,  0.1146356_r_def, &
        0.1304298_r_def,  0.1472430_r_def,  0.1650752_r_def,  0.1839264_r_def, &
        0.2037966_r_def,  0.2246857_r_def,  0.2465938_r_def,  0.2695209_r_def, &
        0.2934670_r_def,  0.3184321_r_def,  0.3444162_r_def,  0.3714396_r_def, &
        0.3998142_r_def,  0.4298913_r_def,  0.4620737_r_def,  0.4968308_r_def, &
        0.5347160_r_def,  0.5763897_r_def,  0.6230643_r_def,  0.6772068_r_def, &
        0.7443435_r_def,  0.8383348_r_def,  1.0000000_r_def /)

Let me know if you need more info!

[bjlittle]

@dennissergeev Starting to play in this space, so I just thought that I'd give you a minor update...

I've been working on a use case with an oceanographer at the Met Office, where they have 50 depth levels (bathymetry) of ORCA ocean model data (on tri-polar curvilinear grid), which has been pre-processed through a canny edge filter to detect potential temperature gradients. Initially, they want to visualise the data as a point cloud, which is an ideal first use case.

I was able to easily extend geovista to create the geo-referenced point cloud, which kinda looks like this:

• I've rendered the points with their levels values
• Used the cmocean deep colormap, which is super lush
• Added a texture map base layer with opacity (to confirm geo-refrencing)
• Added 10m Natural Earth coastlines (to confirm geo-refrencing)

votemper.mp4

Next steps are to:

• make it easy to reproject the point cloud (say to Plate Carree), which might make interactive interrogation a wee bit more intuitive
• extend to support non-polydata, and integrate into the geovista.bridge API

In hindsight, I should have rendered this example with culling=back to make the things slightly less confusing, but this is just a PoC at the moment.

Baby steps.

[bjlittle]

votemper-coastlines.mp4

[bjlittle]

votemper-coastlines-widget.mp4

[dennissergeev]

Sorry for my slow response, I was on holiday :)

This looks amazing, many thanks for sharing! Is this functionality available in the latest release? (No rush, just wondering.)

And since you've shown that point cloud rendering is possible, I guess extending this to vectors/arrows to show the flow velocity in 3D should be quite straightforward, right?

[bjlittle]

Hey @dennissergeev,

Hope you had an awesome holiday.

This functionality will be available in the forthcoming 0.2.0 release.

Rendering point clouds is a nice entry point into this space, so 3D grids and vectors should follow, yeah.

There are a few other things that I also need to tidy-up, so I'm going to release little and often, if I can. But I'll certainly ping you when I make more progress 👍

[dennissergeev]

Excellent, thank you @bjlittle!

[dennissergeev]

Hi @bjlittle! Great to see the steady pace of geovista releases 😄
Does the latest release include this functionality or is it still work in progress for a future release?
Cheers!

[bjlittle]

Hey @dennissergeev,

Lovely to hear from you again, as always. Great question 😄

Okay, so the 0.2.0 point-cloud support is the first step in that direction, but ultimately it's not what you actually want.

My focus at the moment is banking stable core functionality with the view to presenting a lightning talk at SciPy 2023 in July, where at that point geovista finally steps out from the shadows. I'm really interested to see the initial reaction from the scientific community to geovista. I think it's going to be a wee bit of a watershed moment, which I'm equally nervous and excited about...

Anyways, in terms of what's roughly coming up next and when:

• the main deliverable for 0.3.0 will contain support for projecting lines and coastlines, including the ability to add a configurable graticule (i.e., labelled parallels and meridians) to the mesh to aid visual geolocation
• 0.4.0 will be all about bootstrapping the docs. For me, the docs are the shop window to geovista, so having some decent, quality docs is paramount prior to SciPy.

Both 0.3.0 and 0.4.0 are going to follow on the coattails of 0.2.0 and land prior to SciPy.

If there's still time before SciPy, I'd love to steamline and expose a slicker projection API and also seed patterns for geovista.crs projection classes e.g., address #140

Realistically, post SciPy when the dust settles I'm super keen to finally address this issue.

In my head, that timeline and schedule makes sense... is that okay with you?

[dennissergeev]

Yep, totally okay with me! Good luck at SciPy, I'm confident people will be super excited about geovista - I know my colleagues are!

And do let me know if you want me to test/review anything.

[jbusecke]

Hey everyone, I just came here to
a) thank you for this amazing package. I just made some videos for a project site visit and they look stunning!
b) add a gigantic +1 to supporting 3D (lon, lat, depth) visuzlization. My usual usecases revolve around isosurfaces. This video was done with pyvista-xarray and my first time for me that 3D visualization of my usual ocean/climate model data 'just worked'. Having this sort of plot (with bathymetry) on a 3D globe would be monumental for my work.
Further out, having a way to call this from pyvista-xarray somehow would be the holy grail for the researchers/students in my project but I understand that the basics need to work first 🤗

If there is anything I could do to help you out here, please let me know. Unfortunately I am not well versed in pyvista (and the underlying stack) yet, but could certainly test drive, provide data samples etc.

Again big thanks for this awesome package.

[dennissergeev]

Amazing videos @jbusecke!

[bjlittle]

Thanks so much @jbusecke. Your feedback means a lot and your animations are truly lush 🤩

geovista is really just starting on it's journey, and I'm so glad you found us. Sounds like your use case definitely bolsters the demand for volume rendering for sure 👍

After hanging out with @banesullivan at SciPy 2023, we're keen to join the dots between pyvista-xarray and geovista to open the door for the xarray community to easily render and explore their Earth sciences data using pyvista. I agree with you; the potential is huge.

I also had the pleasure of crossing paths with @mgrover1, who similarly wants to render volumes for radar data... so it's inevitable that this is going to happen, which I'm pretty darn excited about too.

Game on.

[mgrover1]

I would be interested in a sprint on this soon 😄

[tkoyama010]

@jbusecke
It's a great video. Could you post this on the PyVista forum to introduce it to more people? If you can't take the time to do so I will introduce that you created this video and post it.
https://github.com/pyvista/pyvista/discussions/new?category=show-and-tell

[bjlittle]

@dennissergeev @jbusecke @mgrover1

Quick question ... when we get to the point where we have a proof-of-concept for the volume rendering in a feature branch (I'm going to figure out when I can schedule this work), would you guys be happy to kick the tyres to test/play and feedback?

Prior to that it would be ideal to get my hands on some representative data that you wanted to render. Nothing high resolution, if possible, just some data that I can quickly render when developing/testing. Cheers 🍻

Thumbs-up if you're in 👍

[jbusecke]

@jbusecke

It's a great video. Could you post this on the PyVista forum to introduce it to more people? If you can't take the time to do so I will introduce that you created this video and post it.

https://github.com/pyvista/pyvista/discussions/new?category=show-and-tell

Sure thing will do that next week.

[jbusecke]

Prior to that it would be ideal to get my hands on some representative data that you wanted to render. Nothing high resolution, if possible, just some data that I can quickly render when developing/testing. Cheers 🍻

@bjlittle I am thinking mostly about CMIP6 climate model output (this should be very representative for most climate/ocean model output. If you can run this snippet

from xmip.utils import google_cmip_col
from xmip.preprocessing import combined_preprocessing
col = google_cmip_col()
cat = col.search(
    variable_id='so',
    table_id='Omon',
    source_id='CanESM5',
    experiment_id='historical',
    member_id='r1i1p1f1'
)
ds = cat.to_dataset_dict(preprocess=combined_preprocessing)['CMIP.CCCma.CanESM5.historical.Omon.gn']
ds = ds.isel(time=0).squeeze().drop([co for co in ds.coords if co not in ['lon', 'lat']])
ds.to_netcdf('geovista_cmip_demo.nc')

You could get a ~20 MB sample of 3D ocean salinity. Let me know if this works (requires xMIP) otherwise Ill see how to send you this data.

@tkoyama010 pyvista/pyvista#4937

[github-actions]

In order to maintain a backlog of relevant issues, we automatically label them as stale after 180 days of inactivity.

If this issue is still important to you, then please comment on this issue and the stale label will be removed.

Otherwise this issue will be automatically closed in 28 days time.

[tkoyama010]

This issue is still important.

[bjlittle]

@all-contributors please add @dennissergeev for ideas

[allcontributors]

@bjlittle

I've put up a pull request to add @dennissergeev! 🎉

[bjlittle]

@all-contributors please add @mgrover1 for ideas and promotion

[allcontributors]

@bjlittle

I've put up a pull request to add @mgrover1! 🎉

[edmundhenley-mo]

Quick question ... when we get to the point where we have a proof-of-concept for the volume rendering in a feature branch (I'm going to figure out when I can schedule this work), would you guys be happy to kick the tyres to test/play and feedback?

@bjlittle - just found this one - subscribed immediately!
TLDR: once you've got your proof-of-concept done, shout iff you want to extend onwards-and-upwards for some space weather use-cases! Happy to help with some data wrangling & tyre kicking if so - no worries if you've got hands full with use-cases from lower down!

Click for context, and potential good 1st target = WACCM-X

Context

This sort of volume / isosurface rendering could be very interesting for space weather use-cases in Earth's ionosphere, thermosphere and below (where we have more-easily supportable geographic lat/lon coords, as well as more bonkers ones too!). These regions can be highly dynamic, with multi-scale spatiotemporal effects which ~require seeing ~4D representations to understand full picture - risk losing insights if you collapse data straight to 2D or lower.
Geovista's interactivity - ability to move round from dayside to nightside, or look down on the poles - could be super-useful for enabling workflows which start with exploration of the full-fat datasets.

• i.e. basically exactly the same points @dennissergeev and Cristina made in lightning talks on Geovista!

There's some very nice videos here illustrating this sort of thing - look at all of those ripples as the solar storm hits on 17 Mar 2015 - the 4D synoptic view here invaluable - a lower-dimensional representation couldn't cut it!

• Thermospheric temperature
• Ionospheric electron number density - volumetric equivalent of the total electron content in iris examples
• More details on what's being plotted in Figure 1 of Palmroth+ 2021

From dim recollection I believe @eelcodoornbos generated these with Blender (?). In any case, using output from NCAR's WACCM-X model.

WACCM-X: a good first target?

Iff this is of interest, WACCM-X might be a good target for dipping geovista's toes in the space weather water (if you've not already?):

• it's a state-of-the-art community model, coupling into lower atmosphere (Liu+ 2017)
• well-behaved: nice netCDF output, on sane geographic grid, vertical on pressure levels - and data volumes not too big by terrestrial standards, so good for exploration (characteristics in slide 5 of WACCM-X tutorial
• so good place to start!
• Note NCAR's TIE-GCM another good contender: similar characteristics, and very widely used

Later can see if Geovista's ugrid support might be an easy way around trickier things like geomagnetic coordinates in other models, ...

There's WACCM-X output available from the NASA Community Coordinated Modeling Center (CCMC). E.g. for this "St Patrick's day storm", there's run "WACCMX-Weimer-01_2015-03-TP-01_102523_IT_1", and from there:

• Data download
• Vis API

Happy to put some effort into giving this a whirl if of interest!

[bjlittle]

@all-contributors please add @edmundhenley-mo for ideas and userTesting

[allcontributors]

@bjlittle

I've put up a pull request to add @edmundhenley-mo! 🎉

[bjlittle]

@edmundhenley-mo Great! 🚀

Good news! I've finally managed to secure a deployment into a science team to work on this ... it's the only way I'm going to get the time to do it!

I'll come knocking at your door afterwards for sure! 👍