triangulation fails with large grids of points (solved?)

GillesDuvert commented 4 years ago

Hello,

I'm considering using delaunator for speeding up GDL's TRIANGULATE function, that badly needs it!

I've made the changes in a branch of my fork of GDL It happens that GDL's TRIANGULATE function is frequently used on regular grids of points (ex: pixels positions on a camera). Certainly not the primary intended use of a triangulation, but "it is the use that makes the function" :smile: . A beta-tester reported failures, so I investigated.

I found that weird triangles were produced in a manner similar to issue #8 when the number of points was large, example: 3017x2018 points in a normalized box of size ~1 (that is, each x and y values are divided by the absolute max of the x and y range. This trick is intended to make the (mathematical) triangulation independent of the magnitude of the x and y values, be it kilometers or microns -- a trick necessary when using legacy triangulation programs of the 60's...)

The following image is the triangulation obtained on a regular grid 3017x2018 , each of the 7340942 triangles painted with a random color. Screenshot_delaunator_3017x2018 The movie of the triangles being painted is something to be seen. Everything is fine until the sweeping hull is near the border and then it goes awry.

I have tracked back the culprit to the size of the hash_key. I changed the size of the haskey by adding 3000 entries, and it worked.

I am afraid delaunator is unfortunately not robust to non-square high number of grid-aligned points:

The size m_hash_size = static_cast<std::size_t>(std::ceil(std::sqrt(n))); is too small for large "non-square" problems. At least with the way the hashkey is computed. Of course, speed is related to the hash size so it is important to maintain it small...
Possibly, there is an adverse effect of the dependency of pseudo_angle() to the actual values of dx and dy. It will make the angle comparison uncertain for close points (very small values of dx-dy)

This drove me to try smaller non-square arrays and I found that wrong triangles were possible for rather small values, the smallest I found was 1602x408.

By multiplying the size of the hastable by the golden ratio (after all this seems related to the ratio of the problem's rectangle to the implicit squareness of sqrt()), all the tests I made worked beautifully, with almost no speed loss.

abellgithub commented 4 years ago

Thanks for the report. I'm not terribly surprised by this. At one point I believed there to be a potential problem with using the hash to find the next point, but things worked, so I left it alone. Do you have specific code changes that you're proposing? If expanding the hash fixes the problem, then there's a more general problem. If not having sufficient entries can cause a failure, then the algorithm is broken.

GillesDuvert commented 4 years ago

The code is here for the calling part, and here for the delaunator.cpp changes. Basically I did this, (1.618... being the golden number). It is just a "works for me" without absolutely no proof of adequateness and sufficiency. However, here the beta tester reports success with a 30000^2 problem, so it may be ok. The only other change was to provision for the list of duplicate points, as GDL is the free clone of IDL and IDL permits to return the list of eventual duplicate points.

abellgithub commented 4 years ago

I don't understand why expanding the size of the hash would make things always work or why you chose the golden number. Why were you thinking that this value would be appropriate?

GillesDuvert commented 4 years ago

I'm pretty sure it will not always work, sorry if I seemed to say so. Concerning the golden number, I experimented that the triangulation was always OK for square 'arrays' and also for arrays with ratio 'width'/'height' larger than a value close to the golden number.

abellgithub commented 4 years ago

Can you tell me the X/Y grid values that showed the issue? Do you have a simple test case?

GillesDuvert commented 4 years ago

As I said, a 1602x408 grid was enough. That is, the x values were [0./1601...1601./1601] and y=[0./1601..407/1601], double precision. It gives this: Screenshot_20200829_170500

The problem does not occur (for this small number of points) if the grid values are [0,1,2....1601] etc.

abellgithub commented 10 months ago

Hi. I have no idea why this got dropped. I must have gotten busy. I'd like to recreate to see if I can fully understand. Do you have a good way I can reproduce the graphics?

GillesDuvert commented 10 months ago

Sorry this is a bit old for me, when I wrote this report I was temporarily a specialist of delaunator :smile: but since then, a lot of water flowed under the bridge as we say in France. As for the implementation in GDL, AFAIK we have received no further comments, hopefully it's because there are no flaws, realistically it's because there are not so many people aware that GDL may be a good replacement for their expensive IDL software and give it a try, and, besides, few people in the possibly involved community will need huge triangulations anyway. Some links in this comment are dead, but the changes have been long ago merged into the GDL repository. Concerning the graphics and the tests, it was done using GDL, so i's in the form of a GDL procedure (two, in fact). I can attach the procedures, but you'll need to install GDL...

abellgithub commented 10 months ago

Yes. Much water under the bridge for me too. I haven't even looked at this repo in years. Thanks for the update.

abellgithub / delaunator-cpp

triangulation fails with large grids of points (solved?) #19