If the input data is a sparse grid instead of a function evaluated on a dense grid, the nice solution above using NaNs may leak:

Here below another approach using Dierckx.jl to interpolate a bounding polygon (and respective values) which is then merged with the input data and fed to PyPlot.jl’s `tricontour()`

.

The main advantages of this workflow is that irregular input data can be consumed and “arbitrary” polygons can be defined.

*PS: Someone versed in PyPlot, may step-in with further advice.*

```
using Dierckx, PyPlot
# Define some input data inside a triangle
f(x,y) = cos(2π*x) * sin(2π*y)
Np = 2_000; x = zeros(Np); y = zeros(Np)
for i in 1:Np
x[i] = rand(0:0.01:1); y[i] = rand(0:0.01:1-x[i])
end
z = f.(x,y)
# define bounding polygon vertices for contouring (here a triangle)
polygon = [0. 0.; 1. 0.; 0. 1.; 0. 0.] # closed triangle (4 points)
N = size(polygon,1)
# Interpolate polygon border using linear parametric splines, and respective values using cubic
t = collect(1:N)
spl1 = ParametricSpline(t, polygon', k=1) # linear splines
tfine = LinRange(0, N, 1000) # polygon boundary samples = 1000
xypoly1 = evaluate(spl1,tfine) # 2 x 1000 matrix
spl2 = Spline2D(x, y, z; kx=3, ky=3, s=1e-4)
zpoly2 = evaluate(spl2, xypoly1[1,:], xypoly1[2,:]) # cubic splines
# merge input data with polygon boundary
xnew = [x; xypoly1[1,:]]
ynew = [y; xypoly1[2,:]]
znew = [z; zpoly2]
xyz = [xnew ynew znew]
xyz = unique(xyz, dims= 1) # keep only unique points
# plot the data support and filled contours
clf()
PyPlot.scatter(x,y,s=0.1)
plt.plot(xypoly1[1,:],xypoly1[2,:], c="red", lw=1)
clf()
PyPlot.tricontourf(xyz[:,1], xyz[:,2], xyz[:,3])
PyPlot.colorbar()
PyPlot.tricontour(xyz[:,1], xyz[:,2], xyz[:,3], colors="black", linewidths=0.5)
```