ST_MapAlgebraExpr — 2バンド版: 二つの入力バンドに対する妥当な PostgreSQL 代数演算で形成された、指定したピクセルタイプとなる 1バンドラスターを生成します。バンドを指定しない場合には、どちらも 1番と仮定します。結果ラスターは、一つ目のラスターのアラインメント (スケール、スキュー、ピクセル角位置) にあわされます。範囲は"extenttype"引数で定義されます。取りうる"extenttype"の値は INTERSECTION, UNION, FIRST, SECOND です。
raster ST_MapAlgebraExpr(raster rast1, raster rast2, text expression, text pixeltype=same_as_rast1_band, text extenttype=INTERSECTION, text nodata1expr=NULL, text nodata2expr=NULL, double precision nodatanodataval=NULL);
raster ST_MapAlgebraExpr(raster rast1, integer band1, raster rast2, integer band2, text expression, text pixeltype=same_as_rast1_band, text extenttype=INTERSECTION, text nodata1expr=NULL, text nodata2expr=NULL, double precision nodatanodataval=NULL);
|
|
|
ST_MapAlgebraExpr は 2.1.0で非推奨になりました。代わりにST_MapAlgebra (expression version) を使います。 |
expressionで定義された妥当な二つのバンドへの PostgreSQL 代数演算を入力ラスター rast1, (rast2) に適用して、一つのバンドを持つラスターを生成します。band1, band2が指定されない場合には、1番バンドと仮定します。新しいラスターは、一つ目のラスターと同じアラインメント (スケール、スキュー、ピクセル隅) を持ちます。新しいラスターは、extenttype引数で定義される範囲になります。
expression
二つのラスターと PostgreSQL 定義済み関数/演算子を含む PostgreSQL 代数式です。関数と演算子は、二つのピクセルがインタセクトするピクセルの値を定めます。たとえば(([rast1] + [rast2])/2.0)::integer といったふうになります。
pixeltype
出力ラスターのピクセルタイプです。必ずST_BandPixelTypeに挙げられたものの一つになるか、省略されるか、NULL に設定されます。引数として渡されないか NULL が渡された場合には、一つ目のラスターのピクセルタイプになります。
extenttype
新しいラスターの範囲を制御します。
INTERSECTION - 新しいラスターの範囲は二つのラスターのインタセクトした領域です。これがデフォルトです。
UNION - 新しいラスターの範囲は二つのラスターの結合です。
FIRST - 新しいラスターの範囲は一つ目のラスターと同じです。
SECOND - 新しいラスターの範囲は二つ目のラスターと同じです。
nodata1expr
rast1が NODATA 値で、特に rast2ピクセルに値がある時に、rast2だけを返すか返すべき値を定義する定数を含む代数式です。
nodata2expr
rast2が NODATA 値で、特に rast2ピクセルに値がある時に、rast1だけを返すか返すべき値を定義する定数を含む代数式です。
nodatanodataval
rast1と rast2のピクセルの両方が NOADTA 値になる場合に返すべき定数です。
pixeltypeが渡された場合には、新しいラスターは、指定されたピクセルタイプのバンドを持ちます。pixeltype として NULL が渡されたりピクセルタイプを指定しない場合には、新しいラスターはrast1と同じピクセルタイプになります。
数式の中で使える語は、元バンドのピクセル値を参照する [rast1.val], [rast2.val]、1始まりの列/行インデックスを参照する[rast1.x], [rast1.y]などです。
Availability: 2.0.0
2バンドのインターセクトしている部分と結合。
元のラスターから 1バンドラスターを生成します。元のラスターバンドの値について 2で割った余りが入ります。
ラスターのいかした集合を生成します。
メートル単位マサチューセッツ州平面でボストン周囲にイかした形状を挿入して、インターセクト部分と結合部分とを示します。
DROP TABLE IF EXISTS fun_shapes;
CREATE TABLE fun_shapes(rid serial PRIMARY KEY, fun_name text, rast raster);
INSERT INTO fun_shapes(fun_name, rast)
VALUES ('ref', ST_AsRaster(ST_MakeEnvelope(235229, 899970, 237229, 901930, 26986), 200, 200, '8BUI', 0, 0));
INSERT INTO fun_shapes(fun_name, rast)
WITH ref(rast) AS (SELECT rast FROM fun_shapes WHERE fun_name = 'ref' )
SELECT 'area' AS fun_name, ST_AsRaster(ST_Buffer(ST_SetSRID(ST_Point(236229, 900930), 26986), 1000),
ref.rast, '8BUI', 10, 0) As rast
FROM ref
UNION ALL
SELECT 'rand bubbles',
ST_AsRaster(
(
SELECT ST_Collect(geom)
FROM (
SELECT ST_Buffer(
ST_SetSRID(ST_Point(236229 + i*random()*100, 900930 + j*random()*100), 26986),
random()*20
) AS geom
FROM generate_series(1, 10) AS i,
generate_series(1, 10) AS j
) AS foo
),
ref.rast, '8BUI', 200, 0
)
FROM ref;
SELECT ST_MapAlgebraExpr(area.rast, bub.rast, '[rast2.val]', '8BUI', 'INTERSECTION', '[rast2.val]', '[rast1.val]') As interrast,
ST_MapAlgebraExpr(area.rast, bub.rast, '[rast2.val]', '8BUI', 'UNION', '[rast2.val]', '[rast1.val]') As unionrast
FROM
(SELECT rast FROM fun_shapes WHERE fun_name = 'area') As area
CROSS JOIN
(SELECT rast FROM fun_shapes WHERE fun_name = 'rand bubbles') As bub
二つの決定的かつ自己完結型のラスターを使た同じインターセクト部分と結合部分の計算。
WITH canvas AS (
SELECT ST_AddBand(
ST_MakeEmptyRaster(80, 80, 0, 80, 1, -1, 0, 0, 0),
1, '8BUI', 40, 0
) AS rast
), inputs AS (
SELECT
rast AS raster_1,
ST_AsRaster(
ST_Buffer(ST_Point(40, 40), 28),
rast,
'8BUI', 220, 0
) AS raster_2
FROM canvas
), variants AS (
SELECT 'raster 1' AS title, raster_1 AS rendered FROM inputs
UNION ALL
SELECT 'raster 2', raster_2 FROM inputs
UNION ALL
SELECT 'intersection', ST_MapAlgebraExpr(
raster_1, raster_2,
'[rast2.val]', '8BUI', 'INTERSECTION',
'[rast2.val]', '[rast1.val]'
) FROM inputs
UNION ALL
SELECT 'union', ST_MapAlgebraExpr(
raster_1, raster_2,
'[rast2.val]', '8BUI', 'UNION',
'[rast2.val]', '[rast1.val]'
) FROM inputs
)
SELECT title, ST_AsPNG(rendered) AS image
FROM variants
ORDER BY CASE title
WHEN 'raster 1' THEN 1
WHEN 'raster 2' THEN 2
WHEN 'intersection' THEN 3
ELSE 4
END;
title | image --------------+--------------------------- raster 1 | PNG image, 80 x 80 pixels raster 2 | PNG image, 56 x 56 pixels intersection | PNG image, 56 x 56 pixels union | PNG image, 80 x 80 pixels
別のバンドとして 1単位 1ピクセルのキャンバスの上にオーバーレイしているラスター。ST_AsPNG を使って画像を描画しています。単一バンドのラスターは全ての灰色に見えます。
WITH mygeoms AS (
SELECT 2 AS bnum, ST_Buffer(ST_Point(1, 5), 10) AS geom
UNION ALL
SELECT 3,
ST_Buffer(
ST_GeomFromText('LINESTRING(50 50,150 150,150 50)'),
10,
'join=bevel'
)
UNION ALL
SELECT 1,
ST_Buffer(
ST_GeomFromText('LINESTRING(60 50,150 150,150 50)'),
5,
'join=bevel'
)
), canvas AS (
SELECT ST_AddBand(
ST_MakeEmptyRaster(
200,
200,
ST_XMin(e)::integer,
ST_YMax(e)::integer,
1, -1, 0, 0),
'8BUI'::text,
0
) AS rast
FROM (
SELECT ST_Extent(geom) AS e
FROM mygeoms
) AS bounds
), rbands AS (
SELECT ARRAY(
SELECT ST_MapAlgebraExpr(
canvas.rast,
ST_AsRaster(m.geom, canvas.rast, '8BUI', 100),
'[rast2.val]', '8BUI', 'FIRST',
'[rast2.val]', '[rast1.val]'
)
FROM mygeoms AS m
CROSS JOIN canvas
ORDER BY m.bnum
) AS rasts
)
SELECT title, ST_AsPNG(rendered) AS image
FROM rbands
CROSS JOIN LATERAL (VALUES
('band 1', rasts[1]),
('band 2', rasts[2]),
('band 3', rasts[3]),
('RGB', ST_AddBand(ST_AddBand(rasts[1], rasts[2]), rasts[3]))
) AS variants(title, rendered);
title | image --------+----------------------------- band 1 | PNG image, 200 x 200 pixels band 2 | PNG image, 200 x 200 pixels band 3 | PNG image, 200 x 200 pixels RGB | PNG image, 200 x 200 pixels
空中写真上に、選択した区画の 2メートル幅の境界線を重ねます。クエリーでは最初に区画ジオメトリーを結合して、元ラスターの破片を領域で切り抜き、小さい破片を結合し、最後に、バンドの先頭 2バンドと、地図代数で描いた第 3バンドを追加します。
新しい 3バンドのラスターを生成します。最初の二つは切り抜いたバンドで、三つ目のバンドはジオメトリーから生成します。クエリーは Windows 用 PostGIS で 3.6秒かかりました。
WITH pr AS (
SELECT ST_Clip(rast, ST_Expand(geom, 50)) AS rast, g.geom
FROM aerials.o_2_boston AS r
INNER JOIN (
SELECT ST_Union(ST_Transform(geom, 26986)) AS geom
FROM landparcels
WHERE pid IN ('0303890000', '0303900000')
) AS g
ON ST_Intersects(rast::geometry, ST_Expand(g.geom, 50))
), prunion AS (
SELECT ST_AddBand(
NULL,
ARRAY[
ST_Union(rast, 1),
ST_Union(rast, 2),
ST_Union(rast, 3)
]
) AS clipped,
geom
FROM pr
GROUP BY geom
)
SELECT ST_AddBand(
ST_Band(clipped, ARRAY[1, 2]),
ST_MapAlgebraExpr(
ST_Band(clipped, 3),
ST_AsRaster(ST_Buffer(ST_Boundary(geom), 2), clipped, '8BUI', 250),
'[rast2.val]',
'8BUI',
'FIRST',
'[rast2.val]',
'[rast1.val]') ) As rast
FROM prunion;
区画境界による上書き。
青線が選択した区画の境界です。 |
上の生成された 4個パネルのキャンバスの例では、外部空中写真や区画テーブルを求めることなしに同じバンド重ね合わせ技術を示しています。