[GIS] Read RAT (raster attribute table) using GDAL or other python libraries

gdalpythonraster-attribute-table

Referring to this post, is there any update on the reason why GDAL function band.GetDefaultRAT() always returns NoneType? I can't find much documentation.

I've tried to create a RAT for testing (using code from here), nonetheless if I then try to read the RAT using band.GetDefaultRAT() I get NoneType, altough the RAT is surely there (I can see it in as xml and using gdalinfo or ArcCatalog). Is there anything I am missing?

I am running from Jupyter Notebook using conda (4.7.11) environment in Anaconda installed on Windows 10 OS, GDAL 2.3.3, released 2018/12/14

# create a RAT for testing
import os
from osgeo import gdal
import numpy as np
import pandas as pd


fn = 'foo.tif'
ds = gdal.Open(fn, gdal.GA_Update)
rb = ds.GetRasterBand(1)

# Get unique values in the band
u = np.unique(rb.ReadAsArray())

# Create a random field as well
r = np.random.randint(0, 1000, size=u.size)

# Create and populate the RAT
rat = gdal.RasterAttributeTable()
rat.CreateColumn('VALUE', gdal.GFT_Integer, gdal.GFU_Generic)
rat.CreateColumn('RANDOM', gdal.GFT_Integer, gdal.GFU_Generic)
for i in range(u.size):
    rat.SetValueAsInt(i, 0, int(u[i]))
    rat.SetValueAsInt(i, 1, int(r[i]))

# Associate with the band
ds.FlushCache()
rb.SetDefaultRAT(rat)

# Close the dataset and delete variables
ds = None
rat = None
rb = None

# try to read the RAT

ds = gdal.Open(fn)
rb = ds.GetRasterBand(1)
rat=ds.GetRasterBand(1).GetDefaultRAT()
print(type(rat))

RETURNS:

<class 'NoneType'>

Output of gdalinfo -stat foot.tif where RAT is visible


PS D:\Dev\Jupiter\notebooks\data> gdalinfo -stats foo.tif
Driver: GTiff/GeoTIFF
Files: foo.tif
foo.tif.aux.xml
Size is 538, 233
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 = (5.970019531250017,47.775634765625000)
Pixel Size = (0.008335632551115,-0.008350238063305)
Metadata:
AREA_OR_POINT=Area
Image Structure Metadata:
INTERLEAVE=BAND
Corner Coordinates:
Upper Left ( 5.9700195, 47.7756348) ( 5d58'12.07"E, 47d46'32.29"N)
Lower Left ( 5.9700195, 45.8300293) ( 5d58'12.07"E, 45d49'48.11"N)
Upper Right ( 10.4545898, 47.7756348) ( 10d27'16.52"E, 47d46'32.29"N)
Lower Right ( 10.4545898, 45.8300293) ( 10d27'16.52"E, 45d49'48.11"N)
Center ( 8.2123047, 46.8028320) ( 8d12'44.30"E, 46d48'10.20"N)
Band 1 Block=538x15 Type=Byte, ColorInterp=Gray
Min=1.000 Max=5.000
Minimum=1.000, Maximum=5.000, Mean=2.342, StdDev=1.332
NoData Value=0
Metadata:
STATISTICS_MAXIMUM=5
STATISTICS_MEAN=2.3424814343154
STATISTICS_MINIMUM=1
STATISTICS_STDDEV=1.3315515074732
<GDALRasterAttributeTable>
<FieldDefn index="0">
<Name>VALUE</Name>
<Type>0</Type>
<Usage>0</Usage>
</FieldDefn>
<FieldDefn index="1">
<Name>RANDOM</Name>
<Type>0</Type>
<Usage>0</Usage>
</FieldDefn>
<Row index="0">
<F>0</F>
<F>636</F>
</Row>
<Row index="1">
<F>1</F>
<F>205</F>
</Row>
<Row index="2">
<F>2</F>
<F>552</F>
</Row>
<Row index="3">
<F>3</F>
<F>33</F>
</Row>
<Row index="4">
<F>4</F>
<F>308</F>
</Row>
<Row index="5">
<F>5</F>
<F>742</F>
</Row>
</GDALRasterAttributeTable>

Best Answer

You are opening the raster in read-only mode.

rs = gdal.Open(fn)

Open it in update mode instead.

rs = gdal.Open(fn, gdal.GA_Update)

With this change your code works in gdal 2.4. However, it does not work for me with gdal 2.2 and 2.3. This may be a result of this commit which appears to fix GeoTIFF GetDefaultRAT() in gdal 2.4 and 3.0.

Related Question