I would like to get some advice on the most efficient way to return a list of unique values of a discrete-valued raster using Python and GDAL.
I had thought that the most obvious way would be to examine the raster's attribute table, but if I do band.GetDefaultRAT()
on the band of a raster dataset that contains an attribute table (the table is visible in ArcCatalog, anyway), the result is always None
:
>>> rat = band.GetDefaultRAT()
>>> rat == None
True
In that case, I end up having to scan through each cell of the raster and build a list of unique values manually. Is this the only way to do it?
Or is there a way to build an attribute table with Python and GDAL, then query it for a list of unique values?
-
Which version of GDAL are you using btw?R.K.– R.K.2012年09月13日 07:29:48 +00:00Commented Sep 13, 2012 at 7:29
-
stupid question but what is RAT?nickves– nickves2012年09月13日 09:28:33 +00:00Commented Sep 13, 2012 at 9:28
-
RAT stands for Raster Attribute Table.Markus M.– Markus M.2016年05月03日 02:59:43 +00:00Commented May 3, 2016 at 2:59
1 Answer 1
If I understood correctly, you can use np.unique function from numpy lib:
from osgeo import gdal
import numpy as np
ds = gdal.Open("myimg.ext")
band = ds.GetRasterBand(1)
array = np.array(band.ReadAsArray())
values = np.unique(array)
or you can one-shot it:
values = np.unique(np.array(ds.GetRasterBand(1).ReadAsArray()))
-
Why wrap
band.ReadAsArray()
in anp.array
call? Doesn't it already return a numpy array?jpmc26– jpmc262018年07月23日 18:08:50 +00:00Commented Jul 23, 2018 at 18:08 -
1Yes it does. Wrapping it in a np.array has no performance drawbacks, as it is not a copy but it uses the same memory address and helps your IDE identify the object so you'll have auto-completion enabled.nickves– nickves2018年07月25日 10:32:55 +00:00Commented Jul 25, 2018 at 10:32