[raster] raster2pgsql -t auto tile size follows GeoTIFF Natural Block size

Commonly found in Cloud Optimized GeoTIFFs.

Closes https://github.com/postgis/postgis/pull/408
Closes #4413

git-svn-id: http://svn.osgeo.org/postgis/trunk@17470 b70326c6-7e19-0410-871a-916f4a2858ee

diff --git a/NEWS b/NEWS
index 2c9b3cb6f413963c6eefc4d61c2dd340534e8844..31ffb614ac2f608280ea2d1c6649f488845e0a57 100644 (file)
--- a/NEWS
+++ b/NEWS
@@ -153,6 +153,9 @@ PostGIS 3.0.0
   - #4198, Add ST_ConstrainedDelaunayTriangles SFCGAL function (Darafei
            Praliaskouski)
   - #4412, Support ingesting rasters with NODATA=NaN (Darafei Praliaskouski)
+  - #4413, Raster tile size follows GeoTIFF block size on raster2pgsql -t auto
+           (Darafei Praliaskouski)
+
 
 
 PostGIS 2.5.0

diff --git a/raster/loader/raster2pgsql.c b/raster/loader/raster2pgsql.c
index 797196c0fc49953d2cb2119b1b06635da589b69c..4d57db6407418d8669e8f8b6d509e1f1853628b6 100644 (file)
--- a/raster/loader/raster2pgsql.c
+++ b/raster/loader/raster2pgsql.c
@@ -450,60 +450,35 @@ usage() {
        ));
 }
 
-static void calc_tile_size(
-       int dimX, int dimY,
-       int *tileX, int *tileY
-) {
-       int i = 0;
-       int j = 0;
-       int min = 30;
-       int max = 100;
-
-       int d = 0;
-       double r = 0;
-       /*int _d = 0;*/
-       double _r = -1;
-       int _i = 0;
-
-       /* j = 0, X */
-       for (j = 0; j < 2; j++) {
-               _i = 0;
-               /*_d = 0;*/
-               _r = -1;
-
-               if (j < 1 && dimX <= max) {
-                       *tileX = dimX;
-                       continue;
-               }
-               else if (dimY <= max) {
-                       *tileY = dimY;
-                       continue;
-               }
-
-               for (i = max; i >= min; i--) {
-                       if (j < 1) {
-                               d = dimX / i;
-                               r = (double) dimX / (double) i;
-
-                       }
-                       else {
-                               d = dimY / i;
-                               r = (double) dimY / (double) i;
-                       }
-                       r = r - (double) d;
-
-                       if (FLT_EQ(_r, -1.0) || (r < _r) || FLT_EQ(r, _r))
+static void
+calc_tile_size(uint32_t dimX, uint32_t dimY, int *tileX, int *tileY)
+{
+       uint32_t min_tile_size = 30;
+       uint32_t max_tile_size = 300;
+
+       for (uint8_t current_dimension = 0; current_dimension <= 1; current_dimension++)
+       {
+               uint32_t img_size = (current_dimension == 0) ? dimX : dimY;
+               uint32_t best_gap = max_tile_size;
+               uint32_t best_size = img_size;
+
+               if (img_size > max_tile_size)
+               {
+                       for (uint32_t tile_size = max_tile_size; tile_size >= min_tile_size; tile_size--)
                        {
-                               /*_d = d;*/
-                               _r = r;
-                               _i = i;
+                               uint32_t gap = img_size % tile_size;
+                               if (gap < best_gap)
+                               {
+                                       best_gap = gap;
+                                       best_size = tile_size;
+                               }
                        }
                }
 
-               if (j < 1)
-                       *tileX = _i;
+               if (current_dimension == 0)
+                       *tileX = best_size;
                else
-                       *tileY = _i;
+                       *tileY = best_size;
        }
 }
 
@@ -1569,6 +1544,8 @@ convert_raster(int idx, RTLOADERCFG *config, RASTERINFO *info, STRINGBUFFER *til
        int _tile_size[2] = {0, 0};
        int xtile = 0;
        int ytile = 0;
+       int naturalx = 1;
+       int naturaly = 1;
        double gt[6] = {0.};
        const char* pszProjectionRef = NULL;
        int tilesize = 0;
@@ -1708,37 +1685,7 @@ convert_raster(int idx, RTLOADERCFG *config, RASTERINFO *info, STRINGBUFFER *til
        info->dim[0] = GDALGetRasterXSize(hdsSrc);
        info->dim[1] = GDALGetRasterYSize(hdsSrc);
 
-       /* tile size is "auto" */
-       if (
-               config->tile_size[0] == -1 &&
-               config->tile_size[1] == -1
-       ) {
-               calc_tile_size(
-                       info->dim[0], info->dim[1],
-                       &(config->tile_size[0]), &(config->tile_size[1])
-               );
-
-               rtinfo(_("Using computed tile size: %dx%d"), config->tile_size[0], config->tile_size[1]);
-       }
-
-       /* decide on tile size */
-       if (!config->tile_size[0])
-               info->tile_size[0] = info->dim[0];
-       else
-               info->tile_size[0] = config->tile_size[0];
-       if (!config->tile_size[1])
-               info->tile_size[1] = info->dim[1];
-       else
-               info->tile_size[1] = config->tile_size[1];
-
-       /* number of tiles */
-       if ((uint32_t)info->tile_size[0] != info->dim[0])
-               ntiles[0] = (info->dim[0] + info->tile_size[0]  - 1) / info->tile_size[0];
-       if ((uint32_t)info->tile_size[1] != info->dim[1])
-               ntiles[1] = (info->dim[1] + info->tile_size[1]  - 1) / info->tile_size[1];
-
-       /* estimate size of 1 tile */
-       tilesize = info->tile_size[0] * info->tile_size[1];
+       tilesize = 0;
 
        /* go through bands for attributes */
        for (i = 0; i < info->nband_count; i++) {
@@ -1753,6 +1700,7 @@ convert_raster(int idx, RTLOADERCFG *config, RASTERINFO *info, STRINGBUFFER *til
                        GDALClose(hdsSrc);
                        return 0;
                }
+               GDALGetBlockSize(hbandSrc, &naturalx, &naturaly);
 
                /* convert data type to that of postgis raster */
                info->bandtype[i] = rt_util_gdal_datatype_to_pixtype(info->gdalbandtype[i]);
@@ -1769,10 +1717,40 @@ convert_raster(int idx, RTLOADERCFG *config, RASTERINFO *info, STRINGBUFFER *til
                                info->nodataval[i] = 0;
                }
 
-               /* update estimated size of 1 tile */
-               tilesize *= rt_pixtype_size(info->bandtype[i]);
+               /* update estimated byte size of 1 pixel */
+               tilesize += rt_pixtype_size(info->bandtype[i]);
+       }
+
+       /* tile size is "auto" */
+       if (config->tile_size[0] == -1 && config->tile_size[1] == -1)
+       {
+               calc_tile_size((naturalx > 1) ? (uint32_t)naturalx : info->dim[0],
+                              (naturaly > 1) ? (uint32_t)naturaly : info->dim[1],
+                              &(config->tile_size[0]),
+                              &(config->tile_size[1]));
+
+               rtinfo(_("Using computed tile size: %dx%d"), config->tile_size[0], config->tile_size[1]);
        }
 
+       /* decide on tile size */
+       if (!config->tile_size[0])
+               info->tile_size[0] = info->dim[0];
+       else
+               info->tile_size[0] = config->tile_size[0];
+       if (!config->tile_size[1])
+               info->tile_size[1] = info->dim[1];
+       else
+               info->tile_size[1] = config->tile_size[1];
+
+       /* number of tiles */
+       if ((uint32_t)info->tile_size[0] != info->dim[0])
+               ntiles[0] = (info->dim[0] + info->tile_size[0] - 1) / info->tile_size[0];
+       if ((uint32_t)info->tile_size[1] != info->dim[1])
+               ntiles[1] = (info->dim[1] + info->tile_size[1] - 1) / info->tile_size[1];
+
+       /* estimate size of 1 tile */
+       tilesize *= info->tile_size[0] * info->tile_size[1];
+
        /* roughly estimate size of one tile and all bands */
        tilesize *= 1.1;
        if (tilesize > MAXTILESIZE)

diff --git a/raster/test/regress/Makefile.in b/raster/test/regress/Makefile.in
index 402fb19d95967833295480b874e81b392760b15e..44a006c932f80d38f271356838fb3b7a3600f00f 100644 (file)
--- a/raster/test/regress/Makefile.in
+++ b/raster/test/regress/Makefile.in
@@ -154,7 +154,8 @@ TEST_LOADER = \
        loader/BasicOutDB \
        loader/Tiled10x10 \
        loader/Tiled10x10Copy \
-       loader/Tiled8x8
+       loader/Tiled8x8 \
+       loader/TiledAuto
 
 TESTS = $(TEST_FIRST) \
        $(TEST_METADATA) $(TEST_IO) $(TEST_BASIC_FUNC) \

diff --git a/raster/test/regress/loader/TiledAuto-post.pl b/raster/test/regress/loader/TiledAuto-post.pl
new file mode 100755 (executable)
index 0000000..97227db
--- /dev/null
+++ b/raster/test/regress/loader/TiledAuto-post.pl
@@ -0,0 +1 @@
+unlink "loader/TiledAuto.tif";

diff --git a/raster/test/regress/loader/TiledAuto-pre.pl b/raster/test/regress/loader/TiledAuto-pre.pl
new file mode 100755 (executable)
index 0000000..6d05f9c
--- /dev/null
+++ b/raster/test/regress/loader/TiledAuto-pre.pl
@@ -0,0 +1 @@
+link "loader/testraster2.tif", "loader/TiledAuto.tif";

diff --git a/raster/test/regress/loader/TiledAuto.opts b/raster/test/regress/loader/TiledAuto.opts
new file mode 100644 (file)
index 0000000..96fbc85
--- /dev/null
+++ b/raster/test/regress/loader/TiledAuto.opts
@@ -0,0 +1 @@
+-t auto -C

diff --git a/raster/test/regress/loader/TiledAuto.select.expected b/raster/test/regress/loader/TiledAuto.select.expected
new file mode 100644 (file)
index 0000000..5379c47
--- /dev/null
+++ b/raster/test/regress/loader/TiledAuto.select.expected
@@ -0,0 +1 @@
+448

diff --git a/raster/test/regress/loader/TiledAuto.select.sql b/raster/test/regress/loader/TiledAuto.select.sql
new file mode 100644 (file)
index 0000000..81d8317
--- /dev/null
+++ b/raster/test/regress/loader/TiledAuto.select.sql
@@ -0,0 +1 @@
+select count(*) from loadedrast;

diff --git a/raster/test/regress/loader/testraster2.tif b/raster/test/regress/loader/testraster2.tif
new file mode 100644 (file)
index 0000000..685de7d
Binary files /dev/null and b/raster/test/regress/loader/testraster2.tif differ