On 15 November 2012 21:25, Chris Barker <chris.bar...@noaa.gov> wrote:
> On Wed, Nov 14, 2012 at 1:50 AM, Ian Thomas <ianthoma...@gmail.com> wrote:
>
> > I think the code used to determine which triangle contains a certain
> point
> > should be factored out into its own TriFinder class,
>
> +1 -- this is a generally useful feature. In fact, it would be nice if
> a lot of this were in a pacakge that deals with triangular meshes,
> apart from MPL altogether (a scikit maybe?)
>
I hadn't considered any wider interest in it. From a selfish matplotlib
point of view, even if it was in say a scikit, we wouldn't want to add an
external dependency to the scikit so we would still need a copy of the
source code within matplotlib anyway, just to do our interpolation for
plotting purposes. The situation would be just like the delaunay code
which exists in a scikit but we include it in the mpl source code to avoid
the external dependency.
Once it is in mpl it is available for other projects to use if they wish,
subject to the BSD-style license. There would be a question of who is
responsible for maintaining it as I don't particularly want to spread
myself thinly on other open source projects when there are a lot of
additions to mpl that I would like to do. There would be the danger for
the other project of the reverse situation of our delaunay example, where
we have code that needs improvement but it is essentially unmaintained,
causing problems for users and embarrassment for developers.
> I have a C++ TriFinder class
> > that I could modify to work within matplotlib, and it is O(log N) so
> should
> > be faster than your version for typical use cases.
>
> What algorithm does this use? Is the code open source and/or availabel
> for other projects?
>
> I'm working on a package for working with unstructured grids in
> general, and also have a use for "what triangle is this point in" code
> for other purposes -- and I havne't found a fast, robust code for this
> yet.
>
It is code I have written recently based on the trapezoidal map algorithm
from the book 'Computational Geometry: Algorithms and Applications' by M.
de Berg et al. It currently exists only on my main linux box, but it will
be open source within mpl for others to use as mentioned above (subject to
acceptance by the mpl developers of course). It is an excellent book that
I wholeheartedly recommend to anyone with a passing interest in
computational geometry. The algorithm itself looks painful initially (as
do many of the algorithms in the book) but it reduces to a surprisingly
small amount of code.
For well-formed triangulations (i.e. no duplicate points, no zero area
triangles and no overlapping triangles) the algorithm is 'sufficiently'
robust. It is not completely robust in the sense of Shewchuk's robust
predicates, but is designed quite cleverly (or luckily) to avoid many of
the pitfalls that occur in naive point-in-triangle tests. For example, a
naive implementation of a point-in-triangle test for a point on or very
near the common edge of two adjacent triangles might return true for both
triangles (which is usually fine), or false for both (which is
catastrophic). The trapezoidal map avoids these problems by reducing the
test to a single point-line test which is therefore guaranteed to return
just one of the two triangles. It is possible to think of a situation that
causes the algorithm to fail at a triangle which has such a small area that
the slopes of two of the sides are within the machine precision of each
other, but it is also possible to use the triangle connectivity to check
for and resolve this problem. I haven't done so yet and need to consider
the tradeoff of effort required vs potential gain.
Someone who required guaranteed robustness could use Shewchuk's point-line
test with the algorithm. I would not do this with mpl however, as the
correctness of the point-line test depends a lot on computer architecture
and compiler flags. This would get out of date quickly and would need keen
monitoring by someone with knowledge of and interest in the area, plus
access to a lot of different computer architectures and OSes. I fail on
all of those counts!
> >> particularly as only a few days ago I committed to writing a triangular
> grid
> >> interpolator for quad grids
>
> what is a triangular interpolator for quad grids? sounds useful, too.
That was poor English from me. I meant interpolating from a triangular
grid *to* a quad grid, typically to make use of the wide range of quad grid
plotting functions like contour, pcolor, etc.
Ian
------------------------------------------------------------------------------
Monitor your physical, virtual and cloud infrastructure from a single
web console. Get in-depth insight into apps, servers, databases, vmware,
SAP, cloud infrastructure, etc. Download 30-day Free Trial.
Pricing starts from $795 for 25 servers or applications!
http://p.sf.net/sfu/zoho_dev2dev_nov
_______________________________________________
Matplotlib-devel mailing list
Matplotlib-devel@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/matplotlib-devel
