ST_Intersection — Returns a raster or a set of geometry-pixelvalue pairs representing the shared portion of two rasters or the geometrical intersection of a vectorization of the raster and a geometry.
setof geomval ST_Intersection(geometry geom, raster rast, integer band_num=1);
setof geomval ST_Intersection(raster rast, geometry geom);
setof geomval ST_Intersection(raster rast, integer band, geometry geomin);
raster ST_Intersection(raster rast1, raster rast2, double precision[] nodataval);
raster ST_Intersection(raster rast1, raster rast2, text returnband, double precision[] nodataval);
raster ST_Intersection(raster rast1, integer band1, raster rast2, integer band2, double precision[] nodataval);
raster ST_Intersection(raster rast1, integer band1, raster rast2, integer band2, text returnband, double precision[] nodataval);
Returns a raster or a set of geometry-pixelvalue pairs representing the shared portion of two rasters or the geometrical intersection of a vectorization of the raster and a geometry.
The first three variants, returning a setof geomval, works in vector space. The raster is first vectorized (using ST_DumpAsPolygons) into a set of geomval rows and those rows are then intersected with the geometry using the ST_Intersection (geometry, geometry) PostGIS function. Geometries intersecting only with a nodata value area of a raster returns an empty geometry. They are normally excluded from the results by the proper usage of ST_Intersects in the WHERE clause.
You can access the geometry and the value parts of the resulting set of geomval by surrounding them with parenthesis and adding '.geom' or '.val' at the end of the expression. e.g. (ST_Intersection(rast, geom)).geom
The other variants, returning a raster, works in raster space. They are using the two rasters version of ST_MapAlgebraExpr to perform the intersection.
The extent of the resulting raster corresponds to the geometrical intersection of the two raster extents. The resulting raster includes 'BAND1', 'BAND2' or 'BOTH' bands, following what is passed as the returnband parameter. Nodata value areas present in any band results in nodata value areas in every bands of the result. In other words, any pixel intersecting with a nodata value pixel becomes a nodata value pixel in the result.
Rasters resulting from ST_Intersection must have a nodata value assigned for areas not intersecting. You can define or replace the nodata value for any resulting band by providing a nodataval[] array of one or two nodata values depending if you request 'BAND1', 'BAND2' or 'BOTH' bands. The first value in the array replace the nodata value in the first band and the second value replace the nodata value in the second band. If one input band do not have a nodata value defined and none are provided as an array, one is chosen using the ST_MinPossibleValue function. All variant accepting an array of nodata value can also accept a single value which will be assigned to each requested band.
In all variants, if no band number is specified band 1 is assumed. If you need an intersection between a raster and geometry that returns a raster, refer to ST_Clip.
|
|
|
For more control over the result extent or NODATA handling, use the two-raster variant of ST_MapAlgebraExpr. To compute a raster-space intersection with a geometry, use ST_Clip. Use this function together with ST_Intersects and an index on the raster or geometry column when filtering table data. |
Enhanced: 2.0.0 - Intersection in the raster space was introduced. In earlier pre-2.0.0 versions, only intersection performed in vector space were supported.
SELECT
foo.rid,
foo.gid,
GeometryType((foo.geomval).geom) AS geom_type,
(foo.geomval).val AS value
FROM (
SELECT
A.rid,
g.gid,
ST_Intersection(A.rast, g.geom) As geomval
FROM dummy_rast AS A
CROSS JOIN (
VALUES
(1, ST_Point(3427928, 5793243.85) ),
(2, 'LINESTRING(3427927.85 5793243.75,3427927.8 5793243.75,3427927.8 5793243.8)'::geometry),
(3, 'LINESTRING(1 2,3 4)'::geometry)
) As g(gid, geom)
WHERE A.rid = 2
) AS foo;
rid | gid | geom_type | value -----+-----+----------------------+------- 2 | 1 | POINT | 249 2 | 1 | POINT | 253 2 | 2 | POINT | 254 2 | 2 | POINT | 251 2 | 2 | POINT | 253 2 | 2 | LINESTRING | 252 2 | 2 | MULTILINESTRING | 250 2 | 3 | GEOMETRYCOLLECTION |
The same line example rendered in a common coordinate frame. The raster footprint and complete input line remain visible behind the returned intersection pieces.
WITH source AS (
SELECT
rast,
ST_GeomFromText(
'LINESTRING(3427927.85 5793243.75,'
'3427927.8 5793243.75,3427927.8 5793243.8)'
) AS input_geometry
FROM dummy_rast
WHERE rid = 2
), intersection AS (
SELECT
rast,
input_geometry,
ST_Collect(gv.geom) AS result_geometry
FROM source
CROSS JOIN LATERAL ST_Intersection(rast, input_geometry) AS gv
GROUP BY rast, input_geometry
), layers AS (
SELECT
'raster' AS title,
ST_Translate(ST_Envelope(rast), -3427927.75, -5793243.75) AS geom
FROM intersection
UNION ALL
SELECT
'input',
ST_Translate(input_geometry, -3427927.75, -5793243.75)
FROM intersection
UNION ALL
SELECT
'intersection',
ST_Translate(result_geometry, -3427927.75, -5793243.75)
FROM intersection
)
SELECT title, geom AS geom
FROM layers
ORDER BY CASE title
WHEN 'raster' THEN 1
WHEN 'input' THEN 2
ELSE 3
END;
title | geom --------------+------------------------------------------------------------- raster | POLYGON((0 0,0 0.25,0.25 0.25,0.25 0,0 0)) input | LINESTRING(0.1 0,0.05 0,0.05 0.05) intersection | GEOMETRYCOLLECTION
geomval, ST_Intersects, ST_MapAlgebraExpr, ST_Clip, ST_AsText