Растры с расстоянием до объекта

[Tifoso]

Добрый день. Помогите, пожалуйста, понять, в чём проблема.
Нужны растровые слои, где в значениях каналов будут расстояния в метрах до тех или иных объектов. Сначала до дороги. Основываюсь на алгоритме, описанном здесь:
http://www.mikemeredith.net/blog/1212_G ... n_QGIS.htm
Проект с СК WGS84/UTM54N.
ИЗ OSM беру слой дорог, обрезаю лишнее, пересохраняю слой в СК проекта - road.
Добавляю слой plots, чтобы ориентироваться на исследуемую площадь.
Растеризую слой road, вручную указывая границы искомого слоя - получаю road_raster.tif о значениями 1 по месту нахождения дороги, остальные пиксели со значением "0".
Далее иду Растр -> Анализ -> Карта близости. Исходный файл road_raster, целевой road_dist. Единицы измерения - GEO. Но прога выдаёт ошибку:

ERROR 6: GDALDriver::Create() ... no create method implemented for this format.
Traceback (most recent call last):
File "C:\OSGeo4W\bin\gdal_proximity.py", line 166, in
dst_ds.SetGeoTransform( src_ds.GetGeoTransform() )
AttributeError: 'NoneType' object has no attribute 'SetGeoTransform'

Если убрать параметр "-of GMT", то процесс проходит, но снова ошибка:

Warning 1: Pixels not square, distances will be inaccurate.

Но растр road_dist добавляется. Как и положено, на слое дорог значения "0", увеличиваются по мере удаления.
На первый взгляд, всё адекватно, но на проверку оказывается, что значения сильно разнятся с фактическими.
Точки, которые линейкой в том же QGISе или SAS.Планете имеют одни значения удалённости от дороги, имеют другие значения, если ориентироваться на значения пикселя слоя road_dist. Наиболее близкие точки имеют разницу в 10-15 м, наиболее отдалённые - в 2000 м (фактические и "линеечные" значения меньше, чем на полученном растре).
Понимаю, что проблема, скорее всего, кроется в выборе проекций/систем координат. Но пробовал всё то же самое с World Equidistant Cylindrical, WGS84/Plate Carree, WGS84/World Mercator, Pseudo Mercator - бесполезно.
Поэтому очень прошу подсказать, в чём моя ошибка, или что дополнительно нужно указать для решения проблемы?
Спасибо!

[Ariki]

Использовать нужно проекцию UTM, остальные из перечисленных сильно искажают расстояния.
Приведите вывод gdalinfo для исходного растра.

[Tifoso]

Значит, WGS 84 / UTM 54N должна подходить.

Вот вывод gdalinfo по итоговому растру с расстояниями (+ txt в приложении):

W>gdalinfo c:\road_dist.tif
Driver: GTiff/GeoTIFF
Files: c:\road_dist.tif
c:\road_dist.tif.aux.xml
Size is 3000, 3000
Coordinate System is:
PROJCS["WGS 84 / UTM zone 54N",
GEOGCS["WGS 84",
DATUM["WGS_1984",
SPHEROID["WGS 84",6378137,298.257223563,
AUTHORITY["EPSG","7030"]],
AUTHORITY["EPSG","6326"]],
PRIMEM["Greenwich",0],
UNIT["degree",0.0174532925199433],
AUTHORITY["EPSG","4326"]],
PROJECTION["Transverse_Mercator"],
PARAMETER["latitude_of_origin",0],
PARAMETER["central_meridian",141],
PARAMETER["scale_factor",0.9996],
PARAMETER["false_easting",500000],
PARAMETER["false_northing",0],
UNIT["metre",1,
AUTHORITY["EPSG","9001"]],
AUTHORITY["EPSG","32654"]]
Origin = (642675.000000000000000,5173418.000000000000000)
Pixel Size = (4.375000000000000,-3.301333333333334)
Metadata:
AREA_OR_POINT=Area
Image Structure Metadata:
INTERLEAVE=BAND
Corner Coordinates:
Upper Left ( 642675.000, 5173418.000) (142d51'58.46"E, 46d41'56.75"N)
Lower Left ( 642675.000, 5163514.000) (142d51'47.43"E, 46d36'36.05"N)
Upper Right ( 655800.000, 5173418.000) (143d 2'16.12"E, 46d41'46.21"N)
Lower Right ( 655800.000, 5163514.000) (143d 2' 4.09"E, 46d36'25.54"N)
Center ( 649237.500, 5168466.000) (142d57' 1.53"E, 46d39'11.25"N)
Band 1 Block=3000x1 Type=Float32, ColorInterp=Gray
Min=0.000 Max=12542.172
Minimum=0.000, Maximum=12542.172, Mean=4546.033, StdDev=3178.091
Metadata:
STATISTICS_MAXIMUM=12542.171875
STATISTICS_MEAN=4546.0333195243
STATISTICS_MINIMUM=0
STATISTICS_STDDEV=3178.0908246734

[Ariki]

Pixel Size = (4.375000000000000,-3.301333333333334)

Тут надо разбираться, почему на этапе растеризации у вас получились не квадратные пиксели. Судя по предупреждению GDAL, алгоритм на такой случай не рассчитан. Можно также, если вы уверены, что с привязкой растра всё в порядке, сделать ресемплинг растра, чтобы у него было одинаковое разрешение по горизонтали и вертикали. Вроде бы в QGIS при сохранении растра можно было задать разрешение. Ну или в консоли gdalwarp с параметром -tr

[Tifoso]

на этапе растеризации у вас получились не квадратные пиксели

Проблема прояснилась. Об этом предупреждал модуль Карта близости сообщением "Warning 1: Pixels not square, distances will be inaccurate".
Моя ошибка состояла в том, что при растеризации вектора, указав координаты крайних точек, я также задал параметр -ts вместо -tr, тем самым установил размер получаемого файла в пикселях вместо того, чтобы задать разрешение пикселя на единицу измерения СК (если я правильно уловил суть).
Большое спасибо, что натолкнули на решение!
Теперь все данные вполне адекватны.

Но всё равно, пока не удалю параметр -of GMT перед запуском Карты близости, так и выскакивает Error 6 (текст в первом посте).

[Ariki]

А зачем вам GMT? Параметр -of задаёт выходной формат, в данном случае - какой-то экзотический вариант NetCDF, который, судя по всему, не поддерживается GDAL на запись (ну или в вашу сборку QGIS не включена поддержка). По умолчанию же создаётся растр в формате GeoTIFF.