<div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr">Hi,<br></div><br><div class="gmail_quote"><div dir="ltr" class="gmail_attr">On Mon, Mar 11, 2019 at 1:21 AM Nicklas Avén <<a href="mailto:nicklas.aven@jordogskog.no">nicklas.aven@jordogskog.no</a>> wrote:<br></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">Hello<br>
<br>
<br>
I have been working on a delauney-triangulation code.<br>
<br>
Link to the code in the bottom<br>
<br>
I need it booth for packaging maps to TilelessMap, but I also need it to <br>
be able to implement editing in the app.<br>
<br>
For packaging I will use it in PostGIS.<br>
<br>
<br>
The results so far looks so promising so I might suggest that it maybe <br>
could go into PostGIS core. Now it is in an extension.<br></blockquote><div><br></div><div>PostGIS is GPL2, your code is GPL3. Please change before I look at it so I can still do something for PostGIS in Delaunay :)</div><div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">All core logic comes from this paper<br>
<br>
<a href="https://www.newcastle.edu.au/__data/assets/pdf_file/0019/22519/23_A-fast-algortithm-for-generating-constrained-Delaunay-triangulations.pdf" rel="noreferrer" target="_blank">https://www.newcastle.edu.au/__data/assets/pdf_file/0019/22519/23_A-fast-algortithm-for-generating-constrained-Delaunay-triangulations.pdf</a><br>
<br>
<br>
What I need is constrained triangulation, so I get only triangles in the <br>
polygons, not in between.<br></blockquote><div><br></div><div>This is cool. There's a hidden function in SFCGAL that does the CDT: <a href="https://trac.osgeo.org/postgis/ticket/4198">https://trac.osgeo.org/postgis/ticket/4198</a> - but in the end it's the same algo as Tesselate?</div><div><br></div><div>CDT is needed for properly generating isochrones: <a href="https://www.patreon.com/posts/isochrones-are-20933638">https://www.patreon.com/posts/isochrones-are-20933638</a> - so I would love to see CDT in core replacing ST_DelaunayTriangles and helping to iron the issues!</div><div><br></div><div>CDT will also help replacing ST_ConcaveHull with a lot faster thing. </div><div><br></div><div><br></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">
From what I know PostGIS then only have ST_Tesselate, which is horribly <br>
slow and very, very picky with what it gets. Since it is more picky than <br>
ST_Isvalid it is very difficult to handle.<br>
<br>
<br>
ST_DelaunayTriangles is quite fast, but cannot constrain (if I haven't <br>
missed something).<br>
<br>
<br>
Anyway, what I need is also the possibility to get the result as indexes <br>
into the original point arrays. That is what the GPU wants. I haven't <br>
implemented that yet, just because I couldn't wrap my head around how to <br>
create arrays in PostgreSQL.<br></blockquote><div><br></div><div>I believe this is also useful for compressed versions of TINs, that should be stored as half edges.</div><div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">I also have a plan to suggest an extension to the twkb-format so the <br>
triangel indexes can be a part of a special geometry-type.<br></blockquote><div><br></div><div>TRIANGLE is already part of geometry. What if it was just supported in serialization as is?</div><div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">
<br>
Now in TilelessMap I encode the triangle indexes as 3D pointsand keep <br>
them in a separate column. But with everything in a twkb geometry we <br>
have a format ready for pushing into a gpu. Just decode the points and <br>
tri-index into float arrays and serve to the GPU and it will be happy.<br>
<br>
That is what I do in Tilelessmap.<br>
<br>
<br>
<br>
Back to the triangulation code.<br>
<br>
It seems pretty fast.<br>
<br>
This is from the largest polygon I have got with 249667 vertex points in <br>
one single polygon (st_dump gives 1 row):<br>
<br>
With ST_DelaunayTriangles unconstrained<br>
<br>
test=# EXPLAIN ANALYZE select ST_DelaunayTriangles(geom) geom from <br>
my_riks WHERE ID = 14681;<br>
QUERY PLAN<br>
-----------------------------------------------------------------------------------------------------------------------------<br>
Index Scan using my_riks_pkey on my_riks (cost=0.29..8.30 rows=1 <br>
width=32) (actual time=2381.912..2386.171 rows=1 loops=1)<br>
Index Cond: (id = 14681)<br>
Planning time: 0.056 ms<br>
Execution time: 2386.185 ms<br>
(4 rows)<br>
<br>
Time: 2386.537 ms (00:02.387)<br>
<br>
<br>
With the new code constrained<br>
<br>
test=# EXPLAIN ANALYZE select triangulate(geom) geom from my_riks WHERE <br>
ID = 14681;<br>
QUERY PLAN<br>
---------------------------------------------------------------------------------------------------------------------------<br>
Index Scan using my_riks_pkey on my_riks (cost=0.29..8.30 rows=1 <br>
width=32) (actual time=522.360..527.080 rows=1 loops=1)<br>
Index Cond: (id = 14681)<br>
Planning time: 0.036 ms<br>
Execution time: 527.095 ms<br>
(4 rows)<br>
<br>
Time: 527.471 ms<br>
<br>
Getting unconstrained from the new code is about the same speed, which <br>
is a little strange since the job to constrain is all done after the <br>
first round of Delauney-triangulation<br>
<br>
test=# EXPLAIN ANALYZE select triangulate(geom,0) geom from my_riks <br>
WHERE ID = 14681;<br>
QUERY PLAN<br>
---------------------------------------------------------------------------------------------------------------------------<br>
Index Scan using my_riks_pkey on my_riks (cost=0.29..8.30 rows=1 <br>
width=32) (actual time=516.255..522.467 rows=1 loops=1)<br>
Index Cond: (id = 14681)<br>
Planning time: 0.039 ms<br>
Execution time: 522.482 ms<br>
(4 rows)<br>
<br>
Time: 522.793 ms<br>
<br>
For smaller geometries the new code is about twice as fast as <br>
ST_DelaunayTriangles<br>
<br>
St_Tessalate cannot be compared on larger polygons, it takes minutes.<br>
<br>
<br>
The new code now handles some invalidities. There are also possibilities <br>
to handle more.<br></blockquote><div><br></div><div>What do you use for inCircle predicate?</div><div><br></div><div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">
<br>
Right now the code doesn't always stop gracefully if you find a polygon <br>
it cannot handle.<br>
<br>
I have even had the code unstoppable, so I had to shut down the machine.<br>
<br>
So, be warned if yo want to test the code, it is a lot of error handling <br>
to do.<br></blockquote><div><br></div><div>Would be cool if it can be made a part of PostGIS. Then we can utilize the coverage infrastructure making sure we can generate tests.</div><div><br></div><div>Also, it is a single-geometry-input thing? Then we can make a fuzzer for it and let Fuzzie look for the crashers. </div><div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">
<br>
But if you can spend a reboot at worst, I would be glad if you can test <br>
it. I have tested on a few million polygons, which seems to work great.<br>
<br>
The code also needs a lot of commenting to be easy to follow (and clean <br>
up and refactoring). I will do that especially if there is interest in <br>
the code.<br>
<br>
To build as extension run make in folder postgis_extension<br>
<br>
The code is here:<br>
<br>
<a href="https://gitlab.com/nicklasaven/triangulation" rel="noreferrer" target="_blank">https://gitlab.com/nicklasaven/triangulation</a><br>
<br>
<br>
Thanks<br>
<br>
<br>
Nicklas Avén<br>
<br>
_______________________________________________<br>
postgis-devel mailing list<br>
<a href="mailto:postgis-devel@lists.osgeo.org" target="_blank">postgis-devel@lists.osgeo.org</a><br>
<a href="https://lists.osgeo.org/mailman/listinfo/postgis-devel" rel="noreferrer" target="_blank">https://lists.osgeo.org/mailman/listinfo/postgis-devel</a></blockquote></div><br clear="all"><div><br></div>-- <br><div dir="ltr" class="gmail_signature"><div dir="ltr"><div><div>Darafei Praliaskouski</div><div>Support me: <a href="http://patreon.com/komzpa" target="_blank">http://patreon.com/komzpa</a></div></div></div></div></div></div></div>