Hello, I'm trying to save a byte raster to a geopackage. I ran into a problem that the geopackage raster driver does not behave as expected when working with bytes and color tables. I have written a python script that explains the issue (attached):
The script creates 4 rasters in similar ways but prints different outputs: GTIFF(no palette): Minimum=0, Maximum=9, Color Table=False, band count=1 GTIFF(palette): Minimum=0, Maximum=9, Color Table=True, band count=1 GPKG(no palette): Minimum=0, Maximum=9, Color Table=False, band count=4 GPKG(palette): Minimum=0, Maximum=255, Color Table=False, band count=4 The problem is that even though we set a color table in the geopackage when writing, it reads as having no color table when reading. Worse yet, the original raster values are lost when written with a color table. Is there any workaround or fix for this? We want to be able to keep non-image byte data with a palette to geopackage. Regards, Ben
import gdal import numpy as np template_array = np.array([[0, 1, 2, 3, 4, 5, 6, 7, 8, 9]]) wgs84_geo_wkt = 'GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AUTHORITY["EPSG","4326"]]' color_table = gdal.ColorTable() for i in range(10): r = int((9 - i) * 255 / 9) color_table.SetColorEntry(i, (r, 127, 127, 255)) def write(driver, path, palette): ds = gdal.GetDriverByName(driver).Create(path, *template_array.shape[::-1], 1, gdal.GDT_Byte) ds.SetGeoTransform((30, 0.01, 0.0, 30, 0.0, -0.01)) ds.SetProjection(wgs84_geo_wkt) band = ds.GetRasterBand(1) band.WriteArray(template_array) if palette: band.SetColorTable(color_table) band.FlushCache() del band del ds def read(name, palette, path): ds = gdal.OpenEx(path, gdal.OF_READONLY) band = ds.GetRasterBand(1) stats = band.GetStatistics(True, True) ct = band.GetColorTable() band_count = ds.RasterCount print( f"{name}({'palette' if palette else 'no palette'}): Minimum={stats[0]:.0f}, Maximum={stats[1]:.0f}," f" Color Table={ct is not None}, band count={band_count}") for name, path, pal in [ ('GTIFF', r'D:\temp\demo_tif.tif', False), ('GTIFF', r'D:\temp\demo_tif_palette.tif', True), ('GPKG', r'D:\temp\demo_gpkg.gpkg', False), ('GPKG', r'D:\temp\demo_gpkg_palette.gpkg', True), ]: write(name, path, pal) read(name, pal, path)
_______________________________________________ gdal-dev mailing list gdal-dev@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/gdal-dev