[postgis-users] Wondering about points3d into contours

[Stephen Woodbridge]

On 9/26/2011 10:48 AM, Sandro Santilli wrote:
> On Sat, Sep 24, 2011 at 11:32:24PM -0400, Stephen Woodbridge wrote:
>> Hi all,
>>
>> I wrote code years ago in C that would take a cloud of 3D points
>> insert them into a Delaunay triangular surface, then slice the
>> surface triangles with a plane(s) into intersecting edges and then
>> composite rings from the edges. It seems like we are very close to
>> being able to do this in PostGIS. I'm not sure if we have the
>> Delaunay triangularization yet.
>>
>> I'm thinking of something like:
>>
>> select b.z, a.* from
>>    st_buildarea(
>>      st_collect(
>>        intersection(
>>          delaunay_triangles(
>>            select the_geom from points),
>>            st_translate(
>>              st_expand(select extents(the_geom) from points, 0.1),
>>              0.0, 0.0, b.z
>>            )
>>          )
>>        )
>>      )
>>    ) as a,
>>    (select 1.0 * generate_series(2, 10, 2) as z) as b
>>   group by b.z;
>>
>> So I'm thinking of something like this is the use case, where this
>> would create contour rings for Z at 2, 4, 6, 8, and 10 based on the
>> triangulated surface created from the points. My SQL is probably
>> broken, but it is only to present the idea.
>>
>> So any ideas on how to do this now?
>> Or if Delaunay triangularization will get added to PostGIS or GEOS?
>
> Delaunay triangularization is in JTS already.
> Should be ported to GEOS first (3.4.0?) and then
> PostGIS (2.1?) could use it.
>
> There's been some interest for the first step but none of that ended
> up turning into actual fundings. It'll surely be helpful to have use
> cases like yours in the C-API interface design phase.
>
> --strk;
>
>    ()   Free GIS&  Flash consultant/developer
>    /\   http://strk.keybit.net/services.html

Hi strk,

Well, I wish I was in a position for offer funding, but I'm not although 
I do occasionally discuss that with clients.

If there is a ticket that you would like me to update with use cases I 
would be happy to do that. This use case is very broad in its 
applicability to real world cases. Think of the points as XY as in 
location of a sensor or some other point of measurement and the Z as the 
measurement value. Creating a 3D triangulated surface and cutting it 
with z-planes creates contours of equal value results. Create a vertical 
plane(s) and you are generating elevation profiles. and if you have 2d 
path like a bicycle route, a GPS track, etc and extrude that into a 
vertical panel the intersects the surface, you intersection becomes the 
measurement profile of the path. If the surface is elevation, then you 
have an elevation profile, if the measurement was chemical 
concentrations in a water body and the track is fish, then the profile 
because an exposure profile for the fish.

So some use cases of this are:

1. driving distance - compute the costs to the node using Dijkstra and 
the XY is the node location and Z is the cost.
2. weather air-presure - XY is sensor location Z is the bars, rings are 
isobars
3. XY location, Z=altitude, height contours one of the more obvious cases

And I think that being able to build a pipeline of functions like I 
offered above would be extremely useful. in all of these use cases.

While I do not need it specifically at the moment, anyone trying to 
render 3D landscapes might want to be able to extract the triangles or 
might want to be able to apply hillshading to the the triangulated 
surface and create a postGIS raster as a result. For the hillshading 
task I believe they would like to get the average normal at the each 
triangle node, this is needed for both Gouraud and Phong shading 
algorithms, being able to extract that might be of interest in the future.

Thank you for the feedback and update.

Best regards,
   -Steve