Я использую PostGIS2.0 для пересечения растров и полигонов. Мне трудно понять, какую операцию мне следует использовать, и какой это самый быстрый способ выполнить. Моя проблема заключается в следующем:
- У меня есть многоугольник и растр
- Я хочу найти все пиксели, которые попадают в многоугольник, и получить сумму значений пикселей
- И (обновленная проблема): я получаю массивные значения для некоторых пикселей, которых нет в исходном растре, когда я выполняю запрос
Мне трудно понять, следует ли мне использовать ST_Intersects()
или ST_Intersection()
. Я также не знаю, как лучше всего суммировать мои пиксели. Вот первый подход, который я попробовал (# 1):
SELECT
r.rast
FROM
raster as r,
polygon as p
WHERE
ST_Intersects(r.rast, p.geom)
Это возвращает список rast
значений, с которыми я не уверен, что делать. Я попытался вычислить итоговую статистику, используя, ST_SummaryStats()
но я не уверен, является ли это взвешенной суммой всех пикселей, которые лежат в пределах многоугольника.
SELECT
(result).count,
(result).sum
FROM (
SELECT
ST_SummaryStats(r.rast) As result
FROM
raster As r,
polygon As p
WHERE
ST_Intersects(r.rast, p.geom)
) As tmp
Другой подход, который я попробовал (# 2), использует ST_Intersection()
:
SELECT
(gv).geom,
(gv).val
FROM
(
SELECT
ST_Intersection(r.rast, p.geom) AS gv
FROM
raster as r,
polygon as p
WHERE
ST_Intersects(r.rast, p.geom)
) as foo;
Это возвращает список геометрий, которые я анализирую дальше, но я предполагаю, что это менее эффективно.
Мне непонятно, какой порядок действий также самый быстрый. Должен ли я всегда выбирать raster, polygon
или polygon, raster
, или конвертировать полигон в растр, чтобы он был raster, raster
?
РЕДАКТИРОВАТЬ: Я обновил подход № 2 с некоторыми подробностями из R.K.
ссылки.
Используя подход № 2, я заметил следующую ошибку в результатах, которая является частью причины, по которой я не понял вывод. Вот изображение моего исходного растра и контур многоугольника, который используется для его пересечения, наложенный сверху:
И вот результат пересечения с использованием PostGIS:
Проблема с результатом заключается в том, что возвращаются значения 21474836, которых нет в исходном растре. Я понятия не имею, почему это происходит. Я подозреваю, что это где-то связано с небольшими числами (деление почти на 0), но возвращает неправильный результат.
ST_SummaryStats()
для № 1, но не уверен, как это сделать для № 2.Ответы:
Я нашел учебник по пересекающимся векторным полигонам с большим растровым покрытием, используя PostGIS WKT Raster . Может иметь ответы, которые вы ищете.
В учебном пособии использовались два набора данных: файл формы точки, который был буферизован для создания полигонов, и серия из 13 растров рельефа SRTM. Между ними было много шагов, но запрос, используемый для пересечения растра и вектора, выглядел так:
Значения были затем обобщены с использованием следующего:
На самом деле я не знаю достаточно PostGIS, чтобы объяснить это, но это звучит так, как будто вы пытались достичь. Учебник должен пролить свет на промежуточные этапы. Удачи :)
источник
Что касается пункта 2 в первоначальном вопросе - в нескольких выпусках разработки Postgis 2.0 использовалась версия библиотеки GDAL, которая преобразует числа с плавающей точкой в целые. Если в вашем растре есть значения с плавающей точкой, и вы использовали версию GDAL ниже 1.9.0 или предварительную версию PostGIS 2.0, которая неправильно вызывала GDALFPolygonize (), то вы можете столкнуться с этой ошибкой. Билеты в трекерах ошибок PostGIS и GDAL были зарегистрированы и закрыты. Эта ошибка была активна во время первоначального вопроса.
С точки зрения производительности вы обнаружите, что использование
ST_Intersects(raster, geom)
намного быстрее, чем использованиеST_Intersects(geom, raster)
. Первая версия растеризует геометрию и выполняет пересечение растрового пространства. Вторая версия векторизует геометрию и делает пересечение векторного пространства, что может быть гораздо более дорогим процессом.источник
Я также имеющие странные вопросы , используя
ST_SummaryStats
сST_Clip
. Запросы к данным по-разному говорят мне, что минимальное значение моего растра было 32, а затем максимум 300, ноST_SummaryStats
возвращалось значение -32700 для значений пикселей в моем целевом многоугольнике.Я закончил взломом вокруг вопроса, как это:
источник