# Gradient vector of a matrix of geographical extent

I have a geographical matrix `M` of dimensions `nLat=181` and `nLon=221` , with `lats = 30:0.25:90` and `lons=-30:0.25:25`.
I want to compute the gradient vector of `M` w.r.t to ` lats , lons`.
I have tried `ForwardDiff , Iterations , ScatteredIterations , Flux ` and my problem is that they all accept continuous functions as input.
I have tried the method interpolating , or this one, but to no avail.

Here is my current code

``````x = lats
y = lons
points = (x,y)
samples = M
itp = Interpolations.interpolate(points, samples, Gridded(Linear()))
``````

I get the following error

``````DimensionMismatch("arrays could not be broadcast to a common size; got a dimension with lengths 181 and 221")
``````

You could adapt the 1D solution in the post linked to 2D as follows:

``````using Interpolations

f(x,y) = x^2 + y^2    # function used just to generate data

x = y = -2:0.02:2
z = f.(x, y')         #  z is the numeric data at coordinates (x, y)

itp = interpolate((x,y), z, Gridded(Linear()));
``````

If the data matrix `z` is noisy, we could fit smoothing 2D splines before estimating the gradients. An example using Dierckx:

``````using Interpolations, Dierckx

f(x,y) = x^2 + y^2
x, y = -2:0.02:2, -3:0.02:3
nx, ny = length(x), length(y)
z = f.(x, y') + rand(nx, ny)            # noisy data matrix over grid defined by (x,y)

s0 = 0.01*nx*ny*hypot(extrema(z)...)    # test different smoothing factors
spl = Spline2D(x, y, z; kx=3, ky=3, s=s0)
zsmooth = evalgrid(spl, x, y)
itp = interpolate((x,y), zsmooth, Gridded(Linear()));
Thank you very much, that works 