Re: [gdal-dev] Check validity of geometries before writing them into vector tiles?

[Blake Thompson]

All,

Geometry validity with polygons is really hard. The only solution at Mapbox
we found to fix this properly for Vector Tiles that we create was to modify
an existing algorithm to solve all the potential problems of geometry
validity. The result of this is https://github.com/mapbox/wagyu/ which we
use in our Vector Tile creation scripts now. Unfortunately it doesn't have
great performance when using extremely high numbers of points in polygons.
However, for vector tiles this typically should not be needed if
simplification is done prior to VT creation. If there is supported needed
around this problem I am more then willing to lend my expertise.

Thanks,

Blake Thompson

On Wed, Feb 14, 2018 at 10:53 AM, Rahkonen Jukka (MML) <
jukka.rahko...@maanmittauslaitos.fi> wrote:

> Even Rouault wrote:
>
> > Perhaps you could play with the SIMPLIFICATION and
> SIMPLIFICATION_MAX_ZOOM options ?
> Sure, as I wrote I admit that my test did not make much sense, but trying
> things before reading the manual sometimes reveals something interesting.
>
> > Perhaps you should also use an already simplified layer for the lowest
> zoom level (see the CONF option)
> For sure yes, coordinate space of 4096x4096 is far too small for this
> dataset.
>
> > Are you sure you get polygons with less than 4 points ? Normally they
> should be discarded.
> Quite sure yes by looking at what the ST_IsValidReason from the SQLite
> dialect prints
> ST_IsValidReason(geometry) (String) = Invalid: Toxic Geometry ... too few
> points
>
> There are other variants of invalid geometries and actually a few valid as
> well. Here is one example with an invalid component
>
>   MULTIPOLYGON (((516384 6815744,516384 6815744,532768 6815744,516384
> 6815744)),
> ((516384 6815744,516384 6815744)),((516384 6815744,516384 6815744)))
>
> My zero tile is here http://www.latuviitta.org/downloads/0.pbf
>
> -Jukka Rahkonen-
> ___
> gdal-dev mailing list
> gdal-dev@lists.osgeo.org
> https://lists.osgeo.org/mailman/listinfo/gdal-dev
___
gdal-dev mailing list
gdal-dev@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/gdal-dev

Re: [gdal-dev] ZSTD compression for GeoTIFF

[Helmut Kudrnovsky]

Even Rouault-2 wrote
> Hi,
> 
> I've just pushed in libtiff upstream, and its internal copy in GDAL,
> support for zstandard/
> zstd compression [1], and the appropriate changes in the GDAL GeoTIFF
> driver.
> This requires building libtiff (or GDAL with its internal libtiff copy)
> against libzstd >= 1.0
> My findings are that single-core compression is at least twice faster with
> zstd than with 
> deflate, for equivalent or better compression rate. Decompression is a bit
> faster, but not 
> that significantly.
> 
> So a compelling use case for zstd is potentially in complex processing
> chains where you 
> have many intermediate internal products. At that point of adoption, I'd
> not recommend it 
> for consumer-facing products.
> 
> Even
> 
> [1] https://github.com/facebook/zstd

 fyi

GRASS trunk has also already implemented ztsd as a raster compression
option.(1)

And it's available for testing already in the OSGeo4W-winGRASS trunk daily
builds.


(1) https://grass.osgeo.org/grass75/manuals/r.compress.html



-
best regards
Helmut
--
Sent from: http://osgeo-org.1560.x6.nabble.com/GDAL-Dev-f3742093.html
___
gdal-dev mailing list
gdal-dev@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/gdal-dev

Re: [gdal-dev] Creating Cloud Optimized GeoTIFFs

[Kurt Schwehr]

Added the delete and validate works now.  Thanks!

On Wed, Feb 14, 2018 at 11:28 AM, Even Rouault 
wrote:

> On mercredi 14 février 2018 11:00:01 CET Kurt Schwehr wrote:
>
> > (working from pretty close to head...)
>
> >
>
> > A follow up question on how to create cogeo files...
>
> >
>
> > I appear to need to do additional work to set BLOCK_OFFSET_?_?. How do I
>
> > go about this?
>
> >
>
> > I'm trying to create a cogeo without gdal_translate in java by using
> create
>
> > and createcopy like this:
>
> >
>
> > https://gist.github.com/schwehr/39680bec7fd8e3e3f840122ea3bafc65
>
>
>
> I guess the issue is the one mentinoned in
>
> http://gdal.org/java/org/gdal/gdal/Driver.html#CreateCopy-
> java.lang.String-org.gdal.gdal.Dataset-int-java.util.Vector-org.gdal.gdal.
> ProgressCallback-
>
>
>
> """
>
> At the end of dataset manipulation, the delete() method must be called on
> the returned dataset otherwise data might not be properly flushed to the
> disk.
>
> """
>
>
>
> --
>
> Spatialys - Geospatial professional services
>
> http://www.spatialys.com
>



-- 
--
http://schwehr.org
___
gdal-dev mailing list
gdal-dev@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/gdal-dev

Re: [gdal-dev] Creating Cloud Optimized GeoTIFFs

[Even Rouault]

On mercredi 14 février 2018 11:00:01 CET Kurt Schwehr wrote:
> (working from pretty close to head...)
> 
> A follow up question on how to create cogeo files...
> 
> I appear to need to do additional work to set BLOCK_OFFSET_?_?.  How do I
> go about this?
> 
> I'm trying to create a cogeo without gdal_translate in java by using create
> and createcopy like this:
> 
> https://gist.github.com/schwehr/39680bec7fd8e3e3f840122ea3bafc65

I guess the issue is the one mentinoned in
http://gdal.org/java/org/gdal/gdal/Driver.html#CreateCopy-java.lang.String-org.gdal.gdal.Dataset-int-java.util.Vector-org.gdal.gdal.ProgressCallback-

"""
At the end of dataset manipulation, the delete() method *must* be called on the 
returned dataset otherwise data might not be properly flushed to the disk.
"""

-- 
Spatialys - Geospatial professional services
http://www.spatialys.com
___
gdal-dev mailing list
gdal-dev@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/gdal-dev

Re: [gdal-dev] Creating Cloud Optimized GeoTIFFs

[Kurt Schwehr]

(working from pretty close to head...)

A follow up question on how to create cogeo files...

I appear to need to do additional work to set BLOCK_OFFSET_?_?.  How do I
go about this?

I'm trying to create a cogeo without gdal_translate in java by using create
and createcopy like this:

https://gist.github.com/schwehr/39680bec7fd8e3e3f840122ea3bafc65

However,  I get this failure when checking the result:

validate_cloud_optimized_geotiff.py cogeo.tif
Traceback (most recent call last):
  File "validate_cloud_optimized_geotiff.py", line 144, in validate
data_offset = int(main_band.GetMetadataItem('BLOCK_OFFSET_0_0', 'TIFF'))
TypeError: int() argument must be a string, a bytes-like object or a
number, not 'NoneType

If I use gdal_translate on the result form the java, it validates.

gdal_translate -co COMPRESS=LZW -co TILED=YES -co COPY_SRC_OVERVIEWS=YES
cogeo.tif cogeo2.tif
validate_cloud_optimized_geotiff.py cogeo2.tif
cogeo2.tif is a valid cloud optimized GeoTIFF

gdalinfo -mdd all -listmdd -mm cogeo.tif > cogeo.tif.info
gdalinfo -mdd all -listmdd -mm cogeo2.tif > cogeo2.tif.info

diff -u cogeo{,2}.tif.info
--- cogeo.tif.info 2018-02-14 10:54:27.050423325 -0800
+++ cogeo2.tif.info 2018-02-14 10:54:41.862227514 -0800
@@ -1,5 +1,5 @@
 Driver: GTiff/GeoTIFF
-Files: cogeo.tif
+Files: cogeo2.tif
 Size is 743, 372
 Coordinate System is:
 GEOGCS["WGS 84",
@@ -19,8 +19,8 @@
 Metadata:
   AREA_OR_POINT=Area
 Metadata (DERIVED_SUBDATASETS):
-  DERIVED_SUBDATASET_1_NAME=DERIVED_SUBDATASET:LOGAMPLITUDE:cogeo.tif
-  DERIVED_SUBDATASET_1_DESC=log10 of amplitude of input bands from
cogeo.tif
+  DERIVED_SUBDATASET_1_NAME=DERIVED_SUBDATASET:LOGAMPLITUDE:cogeo2.tif
+  DERIVED_SUBDATASET_1_DESC=log10 of amplitude of input bands from
cogeo2.tif
 Image Structure Metadata:
   COMPRESSION=LZW
   INTERLEAVE=BAND

gdalinfo -mdd all -listmdd -mm cogeo.tif
Driver: GTiff/GeoTIFF
Files: cogeo.tif
Size is 743, 372
Coordinate System is:
GEOGCS["WGS 84",
DATUM["WGS_1984",
SPHEROID["WGS 84",6378137,298.257223563,
AUTHORITY["EPSG","7030"]],
AUTHORITY["EPSG","6326"]],
PRIMEM["Greenwich",0],
UNIT["degree",0.0174532925199433],
AUTHORITY["EPSG","4326"]]
Origin = (116.349975262217256,39.975075059082918)
Pixel Size = (0.000134747292618,-0.000134747292618)
Metadata domains:
  IMAGE_STRUCTURE
  (default)
  DERIVED_SUBDATASETS
Metadata:
  AREA_OR_POINT=Area
Metadata (DERIVED_SUBDATASETS):
  DERIVED_SUBDATASET_1_NAME=DERIVED_SUBDATASET:LOGAMPLITUDE:cogeo.tif
  DERIVED_SUBDATASET_1_DESC=log10 of amplitude of input bands from cogeo.tif
Image Structure Metadata:
  COMPRESSION=LZW
  INTERLEAVE=BAND
Corner Coordinates:
Upper Left  ( 116.3499753,  39.9750751) (116d20'59.91"E, 39d58'30.27"N)
Lower Left  ( 116.3499753,  39.9249491) (116d20'59.91"E, 39d55'29.82"N)
Upper Right ( 116.4500925,  39.9750751) (116d27' 0.33"E, 39d58'30.27"N)
Lower Right ( 116.4500925,  39.9249491) (116d27' 0.33"E, 39d55'29.82"N)
Center  ( 116.4000339,  39.9500121) (116d24' 0.12"E, 39d57' 0.04"N)
Band 1 Block=256x256 Type=Float32, ColorInterp=Gray
  Description = B5
Computed Min/Max=0.000,0.000
  Overviews: 372x186, 186x93, 93x47, 47x24, 24x12


On Sun, Dec 10, 2017 at 5:34 PM, Kurt Schwehr  wrote:

> Thanks Even for the feedback.  There have be a few offline discussions
> going on and Even added some notes to the document on Trac:
>
> https://trac.osgeo.org/gdal/wiki/CloudOptimizedGeoTIFF?
> action=diff=11
> 
>
> This stems from me switching Earth Engine from LZW to DEFLATE when
> exporting GeoTIFFs (and I added tiling).  We have a report from a user that
> ENVI 5.4 (and 5.1) Classic can't read the resulting images but QGIS &
> ArcGIS can.  I'm reverting exports back to LZW compression.
>
> On Wed, Nov 15, 2017 at 6:02 AM, Even Rouault 
> wrote:
>
>> On mardi 14 novembre 2017 14:20:58 CET Kurt Schwehr wrote:
>>
>> > Hi Even,
>>
>> >
>>
>> > I have some follow up questions on Cloud Optimized GeoTIFFs:
>>
>>
>>
>> The main constraint of C.O.G is that all the IFD definitions are at the
>> beginning of the file, to avoid seeking at various points in it. Other
>> parameters are pretty much unspecified.
>>
>>
>>
>> >
>>
>> > * Is there a preferred/better INTERLEAVE if there is more than one band?
>>
>>
>>
>> Depends on access patterns. If as soon as you process one pixel you need
>> to access the value for all bands, then INTERLEAVE=PIXEL is better, and it
>> will result in smaller sizes of StripOffsets/TileOffsets and
>> StripByteCount/TileByteCount arrays
>>
>>
>>
>> > * Is there a preferred tile blocksize? You have 512 in your examples.
>> Are
>>
>> > there any major trade offs between using 128, 256, 512, or 1024 for x
>> and y
>>
>> > block sizes?
>>
>>
>>
>> Too small blocksizes will result in larger ...Offsets and ...ByteCount
>> arrays.
>>
>>
>>
>> 

Re: [gdal-dev] Check validity of geometries before writing them into vector tiles?

[Rahkonen Jukka (MML)]

Even Rouault wrote:

> Perhaps you could play with the SIMPLIFICATION and SIMPLIFICATION_MAX_ZOOM 
> options ?
Sure, as I wrote I admit that my test did not make much sense, but trying 
things before reading the manual sometimes reveals something interesting.  

> Perhaps you should also use an already simplified layer for the lowest zoom 
> level (see the CONF option)
For sure yes, coordinate space of 4096x4096 is far too small for this dataset.

> Are you sure you get polygons with less than 4 points ? Normally they should 
> be discarded.
Quite sure yes by looking at what the ST_IsValidReason from the SQLite dialect 
prints 
ST_IsValidReason(geometry) (String) = Invalid: Toxic Geometry ... too few points

There are other variants of invalid geometries and actually a few valid as 
well. Here is one example with an invalid component

  MULTIPOLYGON (((516384 6815744,516384 6815744,532768 6815744,516384 6815744)),
((516384 6815744,516384 6815744)),((516384 6815744,516384 6815744)))

My zero tile is here http://www.latuviitta.org/downloads/0.pbf

-Jukka Rahkonen-
___
gdal-dev mailing list
gdal-dev@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/gdal-dev

Re: [gdal-dev] Check validity of geometries before writing them into vector tiles?

[Even Rouault]

On mercredi 14 février 2018 15:10:57 CET Rahkonen Jukka (MML) wrote:
> Hi,
> 
> I admit that my test was somewhat lunatic but I took some random dataset
> from Finland with 68101 polygons and converted data into MVT with default
> settings which means that minzoom was 0.  As a result 12196 of the source
> polygons were written into the 0-level protobuf tile (in EPSG:3067 gridset)
> and none of the polygons is valid. Most polygons have too few points and
> those which have enough points have self-intersections.
> 
> Perhaps there should be some sort of geometry validator in the writer chain?

There is some geometry validation, but it is not white In my tests, with 
the 
ne_10m_admin_1_states_provinces dataset, when I initially implemented very 
strict 
geometry validation, a number of polygons were completely dropped, so I ended 
up 
implementing a more relaxed logic: for an outer ring, it after conversion to 
integer 
coordinates, it is not valid (with GEOSIsValid() testing), then do : buffer( 
tolerance) followed 
by bufffer(-tolerance) followed by simplifypreservetopology(tolerance) (where 
tolerance is 2 
* tile_dim_in_crs_unit / EXTENT) followed by a new round of integer coordinates 
conversion. 
If that's still not valid, keep it.
For an inner ring, drop it if when included in the outer ring, the resulting 
polygon is not valid

Perhaps you could play with the SIMPLIFICATION and SIMPLIFICATION_MAX_ZOOM 
options ?

Perhaps you should also use an already simplified layer for the lowest zoom 
level (see the 
CONF option)

Are you sure you get polygons with less than 4 points ? Normally they should be 
discarded.

Even

-- 
Spatialys - Geospatial professional services
http://www.spatialys.com
___
gdal-dev mailing list
gdal-dev@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/gdal-dev

[gdal-dev] ZSTD compression for GeoTIFF

[Even Rouault]

Hi,

I've just pushed in libtiff upstream, and its internal copy in GDAL, support 
for zstandard/
zstd compression [1], and the appropriate changes in the GDAL GeoTIFF driver.
This requires building libtiff (or GDAL with its internal libtiff copy) against 
libzstd >= 1.0
My findings are that single-core compression is at least twice faster with zstd 
than with 
deflate, for equivalent or better compression rate. Decompression is a bit 
faster, but not 
that significantly.

So a compelling use case for zstd is potentially in complex processing chains 
where you 
have many intermediate internal products. At that point of adoption, I'd not 
recommend it 
for consumer-facing products.

Even

[1] https://github.com/facebook/zstd

-- 
Spatialys - Geospatial professional services
http://www.spatialys.com
___
gdal-dev mailing list
gdal-dev@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/gdal-dev

[gdal-dev] Check validity of geometries before writing them into vector tiles?

[Rahkonen Jukka (MML)]

Hi,

I admit that my test was somewhat lunatic but I took some random dataset from 
Finland with 68101 polygons and converted data into MVT with default settings 
which means that minzoom was 0.  As a result 12196 of the source polygons were 
written into the 0-level protobuf tile (in EPSG:3067 gridset) and none of the 
polygons is valid. Most polygons have too few points and those which have 
enough points have self-intersections.

Perhaps there should be some sort of geometry validator in the writer chain?

-Jukka Rahkonen-
___
gdal-dev mailing list
gdal-dev@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/gdal-dev