PostGIS

Spatial and Geographic objects for PostgreSQL

• Home
• Download
• Documentation
• Development
• Support
• 

Selecting only pixels of particular range of values with ST_Reclass

Author: Regina Obe  2016/03/13  in tips ( newbie, raster )

This raster question comes up quite a bit on PostGIS mailing lists and stack overflow and the best answer often involves the often forgotten ST_Reclass function that has existed since PostGIS 2.0.
People often resort to the much slower though more flexible ST_MapAlgebra or dumping out their rasters as Pixel valued polygons they then filter with WHERE val > 90, where ST_Reclass does the same thing but orders of magnitude faster.

The question goes something like this

I have this raster measuring contaminent levels of ammonia in the pixel values in band 1, and I want to rank contamination by 0 low, 1 medium, 2 high. Then I want to find the area of this contamination or do some other crazy geometric thing to it.

The basic strategy is to reduce your raster into a simpler one that contains 0s, 1s, 2s where 0 ends up being marked as nodata. Then you mark off 0 (or whatever number you choose) as nodata. In the end you end up with a fairly simple raster that is easy to vectorize or keep as raster and do stats on.

So for example here we use ST_Reclass to reclassify all pixel values >= 0 and <= 90 as 0, >90 < 100 as 1, and 100 to 1000 as 2. The last argument makes value of 0 represent nodata.

SELECT ST_Reclass(rast, 1, 
    '[0-90]:0,(90-100):1,[100-1000):2', 
        '4BUI', 0) AS rast
    FROM  sometable
    WHERE filename = '123.tif';

Will give you a new set of rasters that have pixel values of 1,2 and nodata.

Now the next part of the question goes, how do you do geometric operations on this new set, like for example computing the centroid in this case of each toxic level area

WITH cgeoms AS ( SELECT ST_DumpAsPolygons( 
    ST_Reclass(rast, 1, 
    '[0-90]:0,(90-100):1,[100-1000):2', 
        '4BUI', 0), 1
    ) AS gval
    FROM  sometable
    WHERE filename = '123.tif' )
SELECT ST_Centroid( 
    ST_Union( (gval).geom ) ) AS geom,
    (gval).val
FROM cgeoms
GROUP BY (gval).val;

There are also many stats you can do with functions such as ST_Histogram and ST_ValueCount so you do not have to resort to vectorizing your raster to do basic stats on it.

More Posts

Why can't I do CREATE EXTENSION postgis_sfcgal

PostGIS 2.2.0 came out this month, and the SFCGAL extension that offers advanced 3D and volumetric support, in addition to some extended 2D functions like ST_ApproximateMedialAxis became a standard PostgreSQL extension and seems to be a fairly popular extension.

I’ve seen several reports on GIS Stack Exchange of people trying to install PostGIS SFCGAL and getting things such as

ERROR: could not open extension control file /Applications/Postgres.app/Contents/Versions/9.5/share/postgresql/extension/postgis_sfcgal.control

as in this question: SFCGAL in PostGIS problem.

There are two main causes for this:

Read more...

---

Enabling PostGIS raster drivers and out-of-db raster support

I’ve been asked this question in some shape or form at least 3 times, mostly from people puzzled why they get this error. The last iteration went something like this:

I can’t use ST_AsPNG when doing something like

SELECT ST_AsPNG(rast) 
    FROM  sometable;

Gives error: Warning: pg_query(): Query failed: ERROR: rt_raster_to_gdal: Could not load the output GDAL driver.

Read more...

---