маскировать и извлекать значения ячеек из файла vrt?

У меня есть растровые данные для застроенных территорий по всему миру с разрешением 40 м в виде файла vrt, скачать данные отсюда , и я пытаюсь обрезать данные по маске, а затем извлечь значение индекса цвета для каждой ячейки.

Примечание: существуют еще 2 файла с данными: vrt.clr и vrt.ovr

Вот пример данных: представление данных vrt в arcmap.

Мой вопрос: почему при обрезке по маске я получаю значения пустых ячеек?

Я пробовал следующее:

  • извлечь по маске с помощью панели инструментов arcmap

  • использование gdal в python 2.7

    import gdal ds = gdal.Open('input.vrt') ds = gdal.Translate('output.vrt', ds, projWin = [80.439,5.341,81.048,4.686]) ds = None

  • Я также попытался сохранить данные в формате tif

Кроме того, есть ли способ прочитать значение индекса цвета в заданных координатах (x, y) после маскирования данных?


python-2.7 raster arcgis gdal
person Jamal Dabbeek    schedule 07.03.2018    source источник

Ответы (1)


arrow_upward
2
arrow_downward

Данные представлены в псевдомеркаторской проекции (EPSG 3857). Поэтому вам следует либо указать экстент для projWin в этой системе координат, либо добавить projWinSRS, если вы хотите предоставить их в другой системе координат.

Кроме того, если вы хотите, чтобы gdal.Translate выводилось в файл VRT, вы должны добавить format='VRT. Потому что в вашем фрагменте кода выводится формат файла по умолчанию, то есть GeoTIFF.

Когда я предполагаю, что ваши координаты WGS84 (EPSG 4326), они определяют небольшой регион над океаном к югу от Шри-Ланки. Это не имеет особого смысла, учитывая характер данных.

Если вы хотите прочитать массив, заданный вашими координатами, вы можете использовать:

invrt = 'GHS_BUILT_LDSMT_GLOBE_R2015B_3857_38_v1_0.vrt'
outfile = '/vsimem/tmpfile'

ds = gdal.Translate(outfile, invrt, projWin=[80.439, 5.341, 81.048, 4.686], projWinSRS='EPSG:4326')
data = ds.ReadAsArray()
ds = None
gdal.Unlink(outfile)

Нанесенный массив выглядит следующим образом: введите здесь описание изображения

person Rutger Kassies    schedule 07.03.2018
comment
это сработало отлично !, я изменил координаты, чтобы они соответствовали (EPSG 3857). другой способ возврата: Не удалось обработать определение SRS: EPSG:4326. это мне очень помогло - спасибо - person Jamal Dabbeek; 08.03.2018
comment
Отлично, возможно, ваша переменная окружения GDAL_DATA установлена ​​неправильно? Вот где GDAL находит определения проекций, необходимые для преобразования. - person Rutger Kassies; 08.03.2018
comment
Готово!, установите путь к: (.../Anaconda2/Library/share/gdal) в системных переменных среды Windows с указанием на GDAL_DATA. - person Jamal Dabbeek; 08.03.2018