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 ESMValGroup/ESMValCore#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.