[gdal-dev] Countour polygons instead of lines

Vincent Mora vincent.mora at oslandia.com
Fri Oct 20 10:33:12 PDT 2017


Hi all,

Thanks a lot for the feedback.

I have also resorted to polygonize (geos-polygonize) and ran into the
classification problem that you mention. There is also the issue of
performance if you have a lot of complex polygons, on a high res raster,
contours would take a couple of seconds, and geos-polygonize enought
time for me to Ctrl+C before the end.

At the time I needed to tile the problem anyway (split heigth classes
into smaller polygons such that the spatial index is usable, which is
not the case if you have huge polygons) so the classification issue was
the only problem I faced in the end.

I've started the work with the "easy" part: making properly oriented and
classified closed rings (see https://trac.osgeo.org/gdal/ticket/6526).

I'll also try the gdal_polygonize approach, but I don't think it's the
same thingl: with raster-classif + gdal_polygonize, if you have 3
classes 1,2,3 polygons from 1 can touch a polygon from 3, with contour
lines there will always be a polygon of class 2 in between.

Or am I missing something.

If someone believes I'm doing something useless and/or stupid (that also
happens), please do tell.

V.



Le 20/10/2017 à 01:15, Roberto Ribeiro a écrit :
> 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
> <mailto: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 <mailto: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
>         <https://trac.osgeo.org/gdal/ticket/3129>
>
>
>         _______________________________________________
>         gdal-dev mailing list
>         gdal-dev at lists.osgeo.org <mailto:gdal-dev at lists.osgeo.org>
>         https://lists.osgeo.org/mailman/listinfo/gdal-dev
>         <https://lists.osgeo.org/mailman/listinfo/gdal-dev>
>
>
>
>     _______________________________________________
>     gdal-dev mailing list
>     gdal-dev at lists.osgeo.org <mailto:gdal-dev at lists.osgeo.org>
>     https://lists.osgeo.org/mailman/listinfo/gdal-dev
>     <https://lists.osgeo.org/mailman/listinfo/gdal-dev>
>
>

-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/gdal-dev/attachments/20171020/dcc754b4/attachment-0001.html>


More information about the gdal-dev mailing list