Interpolations.jl icon indicating copy to clipboard operation
Interpolations.jl copied to clipboard

Interpolation on a multi-dimensional curvilinear grid

Open briochemc opened this issue 3 years ago • 1 comments

Maybe it's already doable with Interpolations.jl but I could not figure it out. I have a curvilinear grid of the globe with some data and I want to interpolate it (a fairly common task in Earth sciences). So, for the knots of the grid, I have two m×n arrays, LAT and LON. For the data I have another m×n array V.

I would like to be able to do something like

itp = interpolate((LAT, LON), V)

Here is a non-working minimal example to test (with a plot to illustrate the type of grid):

lat = -10:10
lon = -20:20
# The curvilinear grid:
LAT = [x + 10cos(y/40) for x in lat, y in lon]
LON = [y + 0.01x^3 for x in lat, y in lon]
# Just a plot to illustrate:
using Plots
lonrows = vcat(([x; NaN] for x in eachrow(LON))...)
loncols = vcat(([x; NaN] for x in eachcol(LON))...)
latrows = vcat(([x; NaN] for x in eachrow(LAT))...)
latcols = vcat(([x; NaN] for x in eachcol(LAT))...)
plot([lonrows; loncols], [latrows; latcols], lab="curvilinear grid")
# The data on that grid:
V = rand(size(LON)...)
# Now I would like to interpolate:
using Interpolations
itp = interpolate(LON, LAT, V) # throws no method matching interpolate(::Matrix, ::Matrix, ::Matrix)
itp = interpolate((LON, LAT), V) # throws no method matching interpolate(::Tuple{Matrix, Matrix}, ::Matrix)

The curvilinear grid of this example:

Screen Shot 2021-06-18 at 11 28 57 am

But maybe this is for another package than Interpolations.jl?

briochemc avatar Jun 18 '21 01:06 briochemc

The closest this package comes is Gridded.

itp = interpolate((lat,lon), V, Gridded(Linear()))

What I'm trying to think about is if this could be solved by mapping from the curvilinear grid to the flat grid.

I apologize for the response taking so long.

mkitti avatar Aug 08 '21 17:08 mkitti