PostGIS
Toggle Dark/Light/Auto mode Toggle Dark/Light/Auto mode Toggle Dark/Light/Auto mode Back to homepage

Getting distinct pixel values and pixel value counts of a raster

PostGIS raster has so so many functions and probably at least 10 ways of doing something some much much slower than others. Suppose you have a raster, or you have a raster area of interest – say elevation raster for example, and you want to know the distinct pixel values in the area. The temptation is to reach for ST_Value function in raster, but there is a much much more efficient function to use, and that is the ST_ValueCount function.

ST_ValueCount function is one of many statistical raster functions available with PostGIS 2.0+. It is a set returning function that returns 2 values for each row: a pixel value (value), and a count of pixels (count) in the raster that have that value. It also has variants that allow you to filter for certain pixel values.

This tip was prompted by the question on stackexchange How can I extract all distinct values from a PostGIS Raster?

Examples

Get a distinct pixel value list for band 1 of all raster tiles

Specifying the band number is optional and defaults to 1 if not specified. I like to specify it just cause I work a lot with multi-band rasters. In this case I just have digital elevation data table called dem.

SELECT DISTINCT (pvc).value
 FROM (SELECT ST_ValueCount(dem.rast,1) As pvc
   FROM dem) As f
 ORDER BY (pvc).value;

Now if you are using PostgreSQL 9.3+, you can make this even shorter by using the LATERAL clause (since LATERAL keyword is optional in most cases, you can skip saying LATERAL and write it like this)

SELECT DISTINCT (pvc).value
 FROM dem, ST_ValueCount(dem.rast,1) As pvc
 ORDER BY (pvc).value;

Get a pixel value list for band 1 and total pixels for an area of interest

Now if you have a huge coverage, chances are you only care about a particular area, not the 3 million tiles you have. You might want to also know how many times the value appears just for your area of interest So you really want to combine your arsenal with ST_Clip and ST_Intersects. Here I am using the 9.3 LATERAL short-hand

SELECT  (pvc).value, SUM((pvc).count) As tot_pix
 FROM dem
  INNER JOIN
ST_GeomFromText('LINESTRING(-87.627 41.8819,
 -87.629 41.8830)'
 , 4326) As geom
  ON ST_Intersects(dem.rast, geom),
    ST_ValueCount(ST_Clip(dem.rast,geom),1) As pvc
  GROUP BY (pvc).value
 ORDER BY (pvc).value;