Extract region from a cube with 2D latitude/longitude coordinates

[bouweandela]

✨ Feature Request

Motivation

When working with ocean model data, it would be really convenient if iris.cube.Cube.intersection supported 2D latitude and longitude coordinates. Currently

import iris

cube = iris.load_cube("http://esgf1.dkrz.de/thredds/dodsC/cmip5/cmip5/output1/MPI-M/MPI-ESM-LR/historical/mon/ocean/Omon/r1i1p1/v20120625/tos/tos_Omon_MPI-ESM-LR_historical_r1i1p1_185001-200512.nc")

cube.intersection(
    longitude=(10, 20),
    latitude=(10, 20),
    ignore_bounds=True,
)

throws a

CoordinateMultiDimError: Multi-dimensional coordinate not supported: 'longitude'

These coordinates are 2D because they are close enough to a rectangular grid to make this a meaningful data storage format. Therefore it makes sense to be able to select a rectangular cube slice for this type of data. Any points that are outside the region but inside the rectangular slice could e.g. be masked.

I recently implemented something like this in https://github.com/ESMValGroup/ESMValCore/pull/245, but it would be much nicer if this could be supported in iris.

[trexfeathers]

Thanks @bouweandela! We'd love a PR from you if this is something you'd like to see 🙂

[bouweandela]

Thanks! I'll try to find some time to do this soon.

[trexfeathers]

@bouweandela is this still important to you? We're unlikely to dedicate time to it at the UK Met Office:

• Quite niche
• Quite difficult
• You've demonstrated there are other ways of achieving this

[bouweandela]

Yes, still interesting, I'll try to have a look at this soon. I would have expected this is an interesting feature for the met office too since you've recently added support for meshes? Or is there some other way to extract a region from data that is on non rectangular grid?

[trexfeathers]

I would have expected this is an interesting feature for the met office too since you've recently added support for meshes? Or is there some other way to extract a region from data that is on non rectangular grid?

Ah, this is indeed an important question! But the arrival of truly unstructured meshes has also complicated the question...

• Is the region a rectangle (bounded by lon + lat lines), or bounded by great circles on all sides?
• What to do with cells that are partially in the region?
  • Exclude?
  • Include?
  • Re-shape?

To avoid making assumptions wherever possible, we've been exploring using GeoVista/PyVista for regional extraction (example), since it should be able to support any required case. But it's early days yet.

[ESadek-MO]

We've not forgotten about this, this is a candidate for a future sprint on https://github.com/SciTools/iris/labels/Feature%3A%20GeoVista integration.

[github-actions[bot]]

In order to maintain a backlog of relevant issues, we automatically label them as stale after 500 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.