[Matplotlib-users] Matplotlib with invalid triangulations

[Ian Thomas]

Hello Pat,

The solution is to use the function Triangulation.set_mask() to mask out
the zero-area triangles.  The masked-out triangles will be ignored in
subsequent calls to LinearTriInterpolator, tricontourf, etc.  For example:

# +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-
import matplotlib.tri as mtri
import numpy as np

# manually construct an invalid triangulation having a zero area element
x = np.array([0.0, 1.0, 1.0, 1.0, 2.0])
y = np.array([1.0, 0.0, 2.0, 1.0, 1.0])
z = np.zeros(5)

triangles = np.array( [[0, 1, 2], [1, 3, 2], [1, 4, 2], [0, 4, 1]])

# create a Matplotlib Triangulation
triang = mtri.Triangulation(x,y,triangles)

# ---------- start of new code ----------
xy = np.dstack((triang.x[triang.triangles], triang.y[triang.triangles]))  #
shape (ntri,3,2)
twice_area = np.cross(xy[:,1,:] - xy[:,0,:], xy[:,2,:] - xy[:,0,:])  #
shape (ntri)
mask = twice_area < 1e-10  # shape (ntri)

if np.any(mask):
    triang.set_mask(mask)
# ---------- end of new code ----------

# to perform the linear interpolation
interpolator = mtri.LinearTriInterpolator(triang, z)
m_z = interpolator(1.0, 1.0)
# +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-

Note that I have used a small positive number to test the triangle areas
against rather than zero.  This is to avoid problems with rounding errors.
You may need to alter this threshold.

Ian


On 19 November 2016 at 12:24, Pat Prodanovic <pprodano at gmail.com> wrote:

> Hello,
>
> I am using JR Shewchuk's Triangle to create triangulations for use in
> floodplain modeling. I am using a city's digital terrain model with
> hundreds of thousands of breaklines that constrain where triangles can form
> in the triangulations (streets, rivers, etc). Triangle does this very
> efficiently.
>
> Sometimes the input topology I am using has bad inputs which makes
> Triangle create zero area elements. When I import these triangulations to
> Matplotlib I get the error that such triangulations are invalid (when using
> the LinearTriInterpolator() method). I understand the trapezoid map
> algorithm implemented requires only valid triangulations. So far, so good.
>
> The option of cleaning the input topology before using Matplotlib exists,
> but it is really cumbersome. Rather than topology cleaning, am I am able to
> somehow throw out the zero area elements from the call to
> LinearTriInterpolator() method, and still use the triangulation in
> Matplotlib? My other alternative is to use something other than trapezoidal
> map algorithm, but this is just not computationally efficient.
>
> I've reproduced the following example that illustrates the problem in a
> small code snippet. Any suggestions?
>
> Thanks,
>
> Pat Prodanovic
>
> # +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-
> import matplotlib.tri as mtri
> import numpy as np
>
> # manually construct an invalid triangulation having a zero area element
> x = np.array([0.0, 1.0, 1.0, 1.0, 2.0])
> y = np.array([1.0, 0.0, 2.0, 1.0, 1.0])
> z = np.zeros(5)
>
> triangles = np.array( [[0, 1, 2], [1, 3, 2], [1, 4, 2], [0, 4, 1]])
>
> # create a Matplotlib Triangulation
> triang = mtri.Triangulation(x,y,triangles)
>
> # to perform the linear interpolation
> interpolator = mtri.LinearTriInterpolator(triang, z)
> m_z = interpolator(1.0, 1.0)
> # +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-
>
> _______________________________________________
> Matplotlib-users mailing list
> Matplotlib-users at python.org
> https://mail.python.org/mailman/listinfo/matplotlib-users
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mail.python.org/pipermail/matplotlib-users/attachments/20161121/e966b5b7/attachment.html>