[gdal-dev] Countour polygons instead of lines

[Roberto Ribeiro]

I too have run into this problem in the past, though I believe, if you are
to use gdal_contour you have access to the underlying DEM/TIN/LAS, and if
you have that it's much easier just converting it directly to polygons
(using the technique Jamie posted, or similar ones). It could be a
gdal_contour flag if it makes sense to group it this way (contour lines and
hypsometry maps are closely related, after all), but it makes more sense to
me that it be its own routine.

Turning contour lines into contour polygons is surprisingly difficult, I
found. Closing the polygons is the easy part, but knowing which polygon to
assign which value is considerably trickier, because a polygon ill have
part of its boundary representing one value, and part representing another.
If you have a natively closed polygon with no inner ring, you can start
from there and assign to the polygon its outer boundary's values, then
proceed to the smallest polygon that completely encloses the previous one,
and so on. But it's quite possible you'll have no natively closed loops in
your AoI, most noticeably when doing high resolution work. You can pick a
line at random and do an R-tree search for the closest value, but that can
be pretty heavy if you have a lot of lines (and if you have break lines it
completely throws a wrench at your process). Thinking about it now, it
might work to assign the value during polygon creation time
, but that would mean writing a custom polygonize function. But it might
work.

I've always been curious about it (given that I tried doing it once, but
resorted to the raster approach), does anyone know of an algorithm for this?

2017-10-19 15:57 GMT-02:00 Jamie Adams <jaadfoo at gmail.com>:

> From a user perspective, this can be done pretty easily using gdal_calc.py
> & gdal_polygonize.py. I've used this technique in the past several times.
> This doesn't help with gdal_contour of course, but is a viable way of doing
> the same thing. Maybe worth considering as a new python tool.
>
> Using a random SRTM tile in bash:
>
> for iso in 800 1000 1200 1400 1600
> do
> gdal_calc.py --NoDataValue=0 -A N22E005.tif --calc="A>${iso}" --outfile
> ${iso}.tif
> gdal_polygonize.py -f "ESRI Shapefile" ${iso}.tif ${iso}.shp
> done
>
> On Thu, Oct 19, 2017 at 5:47 AM, Vincent Mora <vincent.mora at oslandia.com>
> wrote:
>
>> Hi Gregers,
>>
>> Thanks for your reply.
>>
>> The question was not "how to" yet, but "shall I" ;o)
>>
>> I think I did something similar to the "ends matching" in a python
>> script that creates those polygons (in a particular context that also
>> requires tiling).
>>
>> What I'd like to know, is if there is a reason why it's not already
>> available in gdal, since there is the need for that and apparently algos
>> to reach the goal.
>>
>> > You cannot always make polygons with them unless you are on a perfect
>> island
>>
>> I'm not sure I understand the that. The polygons I have in mind are
>> polygons with holes (multipolygons for a given class actually), is it
>> related to what you're saying ?
>>
>>
>> V.
>>
>>
>> Le 19/10/2017 à 14:23, Hans Gregers Hedegaard Petersen a écrit :
>> > Hi Vincent,
>> >
>> >> I know this is old stuff (references below), but making polygons
>> instead
>> >> of lines would be a great option for gdal_contour IMO.
>> >>
>> >> It could be also another program included in gdal/app (if it is
>> already,
>> >> I can't find it).
>> >>
>> >> What do you think, shall I add that ? If yes, first or second option ?
>> > I had a similar problem years ago, and after making consistent
>> > orientation [1] it is simply a matter of running through the lines and
>> > matching starting points with endpoints in the linestrings.  This is a
>> > rather fast operation on a spatially indexed dataset.
>> >
>> > You cannot always make polygons with them unless you are on a perfect
>> > island - for instance I processed Denmark in tiles (over 43000 tiles
>> > of 1x1km) and then "glued" the linestrings together afterwards. Speed
>> > was essential, so I first glued all linestrings in a tile (which would
>> > make a large number of closed linestrings), then between tiles (but
>> > only for those that were not yet closed).
>> >
>> >
>> > Cheers,
>> >
>> > Gregers
>> >
>> >
>> > [1]: https://trac.osgeo.org/gdal/ticket/3129
>>
>>
>> _______________________________________________
>> gdal-dev mailing list
>> gdal-dev at lists.osgeo.org
>> https://lists.osgeo.org/mailman/listinfo/gdal-dev
>>
>
>
> _______________________________________________
> gdal-dev mailing list
> gdal-dev at lists.osgeo.org
> https://lists.osgeo.org/mailman/listinfo/gdal-dev
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/gdal-dev/attachments/20171019/a226bd96/attachment-0001.html>