# Making a surface from scattered data

I would to transform scattered data into a surface.

Google told me that is I could do it in Matlab[1,2] or Python  using Delaunay triangulation. But I cannot find a similar example in Julia with VoronoiDelaunay.jl and its related pages.

Does anybody has an example ?

Thank you You can use a Thin Plate Spline / PolyharmonicSpline to do it.

Here is an old gist that has an Julia v0.5ish implementation and example and here is an updated version that works with Julia v1.0

1 Like

just for anyone that needs an example with Julia v1.x:

``````ENV["MPLBACKEND"]="tkagg" # if you need
using PyPlot
using LinearAlgebra

include("polyharmonic_spline.jl")
# contain all code from https://github.com/lstagner/OrbitTomography.jl/blob/master/src/polyharmonic.jl
# and `interpolate()`definition from https://gist.github.com/lstagner/04a05b120e0be7de9915

x,y = randn(500),randn(500)
z = exp.(-(x.^2 .+ y.^2))
S2 = PolyharmonicSpline(2,[x y],z)

n=20
xgrid = ones(n)*range(-3,stop=3,length=n)'
ygrid = range(-3,stop=3,length=n)*ones(n)'

xx = reshape(xgrid,n*n)
yy = reshape(ygrid,n*n)

zz = interpolate(S2,xx,yy)
zgrid = reshape(zz,n,n);

plot_surface(xgrid,ygrid,zgrid,alpha=0.5)
scatter3D(x,y,z,color="r")
show()
`````` 1 Like

With the updated version you don’t need the `interpolate` method. The constructor for the spline returns a interpolating “function”. You can simplify the code to be

``````using PyPlot
using LinearAlgebra

include("polyharmonic_spline.jl")
# contain all code from https://github.com/lstagner/OrbitTomography.jl/blob/master/src/polyharmonic.jl

x,y = randn(500),randn(500)
z = exp.(-(x.^2 .+ y.^2))
S2 = PolyharmonicSpline(2,[x y],z)

n=20
xgrid = ones(n)*range(-3,stop=3,length=n)'
ygrid = range(-3,stop=3,length=n)*ones(n)'

plot_surface(xgrid,ygrid,zgrid,alpha=0.5)
scatter3D(x,y,z,color="r")
show()
``````
2 Likes

I’m not sure why we need spline interpolation to do this. Using the default GR backend, seems like `surface` naturally interpolates via Delaunay triangulation. The following code

``````using Plots

x = randn(1000)
y = randn(1000)
z = @. sin(x)*sin(y)
surface(x,y,z,camera=(0,90))
``````

gives something like