# Contour plot for non-rectangular domain

**URL:** <https://discourse.julialang.org/t/contour-plot-for-non-rectangular-domain/57381>\
**Category:** General Usage\
**Tags:** plotting, visualization, gmt\
**Created:** [March 17, 2021, 1:41pm UTC](https://discourse.julialang.org/t/contour-plot-for-non-rectangular-domain/57381 "2021-03-17T13:41:20Z")\
**Posts on this page:** 7\
**Page:** 1

<div class="post-metadata">

**Author:** ![Gus\_Hart](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/gus_hart/32/6987_2.png) [@Gus\_Hart](https://discourse.julialang.org/u/Gus_Hart)\
**Post date:** [March 17, 2021, 1:41pm UTC](https://discourse.julialang.org/t/contour-plot-for-non-rectangular-domain/57381/1 "2021-03-17T13:41:20Z")

</div>

There was a previous thread on this, but the use case was too different from mine to be informative.

I have a function defined over a triangle and would like to plot it as a contour plot but only in the triangular domain. Is there some way to do this? It’s not hard in mathematica (but seems not to translate to more “normal” plotting software).

For a concrete example to aim at, consider f(x,y) = cos(2pi x)\*sin(2pi y) over the triangle with vertices (0,0), (1,0), (0,1). I’d like the plot to be _blank_ outside the triangle.

---

<div class="post-metadata">

**Author:** ![jacobadenbaum](https://avatars.discourse-cdn.com/v4/letter/j/5daacb/32.png) [@jacobadenbaum](https://discourse.julialang.org/u/jacobadenbaum)\
**Post date:** [March 17, 2021, 2:38pm UTC](https://discourse.julialang.org/t/contour-plot-for-non-rectangular-domain/57381/2 "2021-03-17T14:38:48Z")

</div>

`Plots.jl` ignores `NaN` values, so while this may be a bit hacky, you can do something like this:

```julia
using Plots
f(x,y) = x < 1 - y ? cos(2 * pi * x) * sin(2 * pi * y) : NaN
grid = LinRange(0,1, 1000)
contourf(grid, grid, f)

```

![test](https://global.discourse-cdn.com/julialang/original/3X/0/d/0d49ccf6da091625190cce7a97665bb1040d566d.png)  
You can obviously generalize this to arbitrary domains so long as it is easy for you to check whether a point lies inside it.

[Edit: I misread the domain you were asking about in the initial response]

---

<div class="post-metadata">

**Author:** ![Gus\_Hart](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/gus_hart/32/6987_2.png) [@Gus\_Hart](https://discourse.julialang.org/u/Gus_Hart)\
**Post date:** [March 17, 2021, 5:31pm UTC](https://discourse.julialang.org/t/contour-plot-for-non-rectangular-domain/57381/3 "2021-03-17T17:31:32Z")

</div>

I like that. It works. If there was some way to define a plot region, that would be cleaner, but this totally does the trick. Thank you.

---

<div class="post-metadata">

**Author:** ![rafael.guerra](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/rafael.guerra/32/216610_2.png) [@rafael.guerra](https://discourse.julialang.org/u/rafael.guerra)\
**Post date:** [March 17, 2021, 11:17pm UTC](https://discourse.julialang.org/t/contour-plot-for-non-rectangular-domain/57381/4 "2021-03-17T23:17:53Z")

</div>

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:  
 ![Plots_contourf_problem_of_NaNs](https://global.discourse-cdn.com/julialang/original/3X/1/3/13069ddbaee62c9bab5f4ca035afdd7850f91612.png)

Here below another approach using [Dierckx.jl](https://github.com/kbarbary/Dierckx.jl/blob/master/README.md) to interpolate a bounding polygon (and respective values) which is then merged with the input data and fed to [PyPlot.jl](https://github.com/JuliaPy/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._

```julia
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)

```

![PyPlot_input_data_and_polygon](https://global.discourse-cdn.com/julialang/original/3X/4/4/44c6b2c583163f00f5fb6af2aa1bcee7e3f1381b.png)  
 ![PyPlot_contourf_inside_polygon](https://global.discourse-cdn.com/julialang/original/3X/6/1/61114c9db52324e1aa953051b4c71863c89c4860.png)

---

<div class="post-metadata">

**Author:** ![joa-quim](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/joa-quim/32/227_2.png) [@joa-quim](https://discourse.julialang.org/u/joa-quim)\
**Post date:** [March 18, 2021, 3:34am UTC](https://discourse.julialang.org/t/contour-plot-for-non-rectangular-domain/57381/5 "2021-03-18T03:34:02Z")

</div>

It was already possible to do this with GMT but it was a bit cumbersome, now with master version one can

```julia
f(x,y) = cos(2 * pi * x) * sin(2 * pi * y);
x = linspace(0,1,100);
G = mat2grid(f, x, x); # because countourf cannot yet ingest functions

# Clip inside a triangle for simplicity but the shape could be anything and any number of shapes
contourf(G, clip=[0.1 0.1; 0.5 0.9; 0.9 0.1], fmt=:png, show=true)

```

 ![GMTjl_tmp](https://global.discourse-cdn.com/julialang/original/3X/4/9/49064afba57037a34904223bac42a1aa46f460bd.png)

---

<div class="post-metadata">

**Author:** ![Gus\_Hart](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/gus_hart/32/6987_2.png) [@Gus\_Hart](https://discourse.julialang.org/u/Gus_Hart)\
**Post date:** [March 18, 2021, 1:23pm UTC](https://discourse.julialang.org/t/contour-plot-for-non-rectangular-domain/57381/6 "2021-03-18T13:23:27Z")

</div>

What is GMT?

---

<div class="post-metadata">

**Author:** ![joa-quim](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/joa-quim/32/227_2.png) [@joa-quim](https://discourse.julialang.org/u/joa-quim)\
**Post date:** [March 18, 2021, 1:25pm UTC](https://discourse.julialang.org/t/contour-plot-for-non-rectangular-domain/57381/7 "2021-03-18T13:25:24Z")

</div>

[GMT.jl](https://github.com/GenericMappingTools/GMT.jl) the best mapping package 🙂
