Базовая оценка транспортной доступности средствами GRASS GIS и QGIS

Материал из GIS-Lab
(Различия между версиями)
Перейти к: навигация, поиск
(Загрузка данных OSM)
м
 
(не показаны 26 промежуточных версий 2 участников)
Строка 1: Строка 1:
{{Статья|Черновик}}
+
{{Статья|Опубликована|isochrone-map-grass-qgis}}
  
{{Аннотация|Построение карт транспортной доступности на основе данных OpenStreetMap средствами открытых геоинформационных систем GRASS GIS 7 и QGIS}}
+
{{Аннотация|Построение карт транспортной доступности на основе данных OpenStreetMap средствами открытых геоинформационных систем GRASS GIS и QGIS}}
  
Картографирование транспортной доступности на основе данных о дорожной сети - одна из классических задач ГИС. Наиболее распространенным способом моделирования транспортной доступности является построение изохрон - линий равных затрат времени на преодоление пространства из относительно заданных точек. В представленной статье обсуждается алгоритм построения изохрон по данным [http://www.openstreetmap.org/ OpenStreetMap] с использоваем открытых ГИС [https://grass.osgeo.org/ GRASS GIS 7] и [http://qgis.org/ru/site/ QGIS]. В QGIS будет осуществляться подготовка данных и картографическое представление результатов, а в GRASS собственно моделирование. Всю работу можно выполнить целиком в GRASS, но, по мнению автора, общие манипуляции геоданными и представление картографических материалов удачнее и удобнее реализованы в QGIS.
+
Картографирование транспортной доступности на основе данных о дорожной сети - одна из классических задач ГИС. Наиболее распространенным способом моделирования транспортной доступности является построение изохрон - линий равных затрат времени на преодоление пространства относительно заданных точек. В представленной статье обсуждается алгоритм построения изохрон по данным [http://www.openstreetmap.org/ OpenStreetMap] с использованием открытых ГИС [https://grass.osgeo.org/ GRASS GIS] и [http://qgis.org/ru/site/ QGIS]. В QGIS будет осуществляться подготовка данных и картографическое представление результатов, а в GRASS собственно моделирование. Всю работу можно выполнить целиком в GRASS, но, по мнению автора, общие манипуляции геоданными и представление картографических материалов удачнее и удобнее реализованы в QGIS.
  
Описываемые в статье действия выполнялись в средах Ubuntu 16.04 и Windows 8.1 (x64), версия GRASS 7.2, версия QGIS 2.14.
+
Описываемые в статье действия выполнялись в средах Ubuntu 14.04 LTS и Windows 8.1 (x64), GRASS GIS 7.2, QGIS 2.14.
  
 
==Получение и подготовка данных==
 
==Получение и подготовка данных==
  
Для осуществления расчётов нам потребуется набор векторных линейных геоданных, содержащий информацию о дорожно-транспортной сети исследуемой территории. Заполучить подобную информацию можно различными способами: приобрести у специализированных поставщиков, найти в одном из источников открытых данных, оцифровать атлас автомобильных дорог, нарисовать по космическому снимку и так далее (не забывайте про лицензии данных!). В данном случае мы воспользуемся данными OSM как достаточно качественными и подробными, и, что важно, доступными бесплатно и легально. Начать знакомство с OSM вы можете с [http://wiki.openstreetmap.org/wiki/RU:%D0%9E_%D0%BF%D1%80%D0%BE%D0%B5%D0%BA%D1%82%D0%B5 этой страницы]. Этап подготовки данных являются ключевым: помимо технических этапов здесь нам предстоит определить транспортные характеристики, которые будут использованы для моделирования.
+
Для осуществления расчётов нам потребуется набор векторных линейных геоданных, содержащий информацию о дорожно-транспортной сети исследуемой территории. Заполучить подобную информацию можно различными способами: приобрести у специализированных поставщиков, найти в одном из источников открытых данных, оцифровать атлас автомобильных дорог, нарисовать по космическому снимку, и так далее (не забывайте про лицензии данных!). В данном случае мы воспользуемся данными OSM как достаточно качественными и подробными, и, что важно, доступными бесплатно и легально. Начать знакомство с OSM вы можете с [http://wiki.openstreetmap.org/wiki/RU:%D0%9E_%D0%BF%D1%80%D0%BE%D0%B5%D0%BA%D1%82%D0%B5 этой страницы]. Этап подготовки данных являются ключевым: помимо технических этапов, здесь нам предстоит определить транспортные характеристики, которые будут использованы для моделирования.
  
Для демонстрации принципов базового транспортного моделирования рассмотрим Тосненский район Ленинградской области.
+
Для демонстрации принципов базового транспортного моделирования рассмотрим Тосненский район Ленинградской области. Транспорт: автомобиль.
  
 
===Загрузка данных OSM===
 
===Загрузка данных OSM===
  
Для загрузки данных OSM существует множество возможностей. Одним из наиболее простых способов является загрузка уже подготовленных в формате ESRI Shapefile наборов данных по слоям. Популярные сервисы в сети: [http://download.geofabrik.de/ Geofabrik], где вы можете найти комплекты данных на весь мир, и [http://data.nextgis.com/osmshp/ NextGIS], где вы найдёте вместе с данными оформленные QGIS-проекты по субъектам РФ и странам СНГ. Также удобный способ быстро получить данные OSM - воспользоваться одним из соответствующих плагинов для QGIS, например [http://plugins.qgis.org/plugins/OSMDownloader/ OSMDownloader]. При его установке в интерфейс QGIS добавится кнопка выделения прямоугольной области, для которой и будут загружены данные в указанное вами расположение. Для демонстрации в этой статье используется набор данных OSM от NextGIS на территорию Ленинградской области от 13 ноября 2016 года, доступный бесплатно ([https://drive.google.com/file/d/0BztNh0nXzyjbTUpNWVl4bk15Sjg/view?usp=sharing загрузить]).
+
Для загрузки данных OSM существует множество возможностей. Одним из наиболее простых способов является загрузка уже подготовленных наборов данных по слоям в формате ESRI Shapefile. Популярные сервисы в сети: [http://download.geofabrik.de/ Geofabrik], где вы можете найти комплекты данных на весь мир, и [http://data.nextgis.com/osmshp/ NextGIS], где вы найдёте вместе с данными оформленные QGIS-проекты по субъектам РФ и странам СНГ. Также удобный способ быстро получить данные OSM - воспользоваться одним из соответствующих плагинов для QGIS, например, [http://plugins.qgis.org/plugins/OSMDownloader/ OSMDownloader]. При его установке в интерфейс QGIS добавится кнопка выделения прямоугольной области, для которой и будут загружены данные в указанное вами расположение. Для демонстрации в этой статье используется набор данных OSM от NextGIS на территорию Ленинградской области от 13 ноября 2016 года, доступный бесплатно ([https://drive.google.com/file/d/0BztNh0nXzyjbTUpNWVl4bk15Sjg/view?usp=sharing загрузить]).
  
===Первичная подготовка данных===
+
Для нашей задачи понадобятся только данные о дорожной сети (в случае работы с выгрузками NextGIS - слой ''highway-line'' из директории ''data''; в случае работы с OSMDownloader - все загруженные объекты линейного типа, у которых в атрибуте "Highway" содержится какое-либо значение). Для оформительских задач также может быть использован набор геоданных с границами административно-территориального деления (например, для обрезки изохрон по границам района).
  
===Определение транспортных характеристик данных о дорожной сети===
+
===Первичная подготовка данных в QGIS===
 +
 
 +
Для начала подготовим данные по дорожной сети: обрежем их по границам района интереса, спроецируем куда нужно, отфильтруем лишнее. Распаковываем архив ''RU-LEN.7z'' в удобное расположение, лучше без кириллицы и пробелов в путях. Затем запускаем QGIS и добавляем два векторных слоя из папки .../RU-LEN/data/:
 +
* highway-line
 +
* boundary-polygon
 +
 
 +
Нужно проявить осторожность с кодировками, и добавлять слои с кодировкой UTF-8 (кодировку можно выбрать в меню выбора векторного слоя при его добавлении). Выглядит это как-то так (цвета могут быть произвольными):
 +
 
 +
[[Файл:Grass_qgis_isochrones_map_raw.png|600px|thumb|center|]]
 +
 
 +
В первую очередь, избавимся от всего лишнего. Найдём полигон Тосненского района и сохраним его в отдельный слой. Для этого в таблице атрибутов слоя ''boundary-polygon'' (доступно через контекстное меню слоя) найдём объект с соответствующим содержанием поля NAME. Для этого можно воспользоваться поиском по запросу, либо просто отсортировать колонку и выделить нужную строку. Затем через контекстное меню слоя выбираем "Сохранить как", определяя новое расположение для файла с границами Тосненского района (в моём случае Tosnensky_boundary.shp), и не забываем активировать флаг "Сохранить только выделенные объекты".
 +
 
 +
 
 +
[[Файл:Grass_qgis_isochrones_table_selection.png|400px|thumb|center|]]
 +
[[Файл:Grass_qgis_isochrones_save_as.png|400px|thumb|center|]]
 +
 
 +
Затем необходимо обрезать набор геоданных по дорожной сети по этим границам. Для этого воспользуемся инструментом обрезки (Вектор - Инструменты геообработки - Обрезка). В качестве исходного слоя выбираем ''highway-line'', в качестве слоя обрезки ''Tosnensky_bondary'', результат сохраняем в удобное новое расположение (в моём случае Tosnensky_roads.shp).
 +
 
 +
[[Файл:Grass_qgis_isochrones_clip.png|400px|thumb|center|]]
 +
 
 +
Теперь, чтобы избежать возможных проблем с расчётом длин участков дорог (и некоторых других) в будущем, перепроецируем результирующий слой ''Tosnensky_roads'' в подходящую систему координат. Для расчётных задач рекомендуется использовать поперечную цилиндрическую проекцию (UTM или Гаусса-Крюгера) для вашей зоны. В данном случае используем WGS 84 / UTM zone 36N (EPSG:32636), для этого через меню "Сохранить как" создаем копию слоя ''Tosnensky_roads'', в выпадающем списке "Система координат", выбрав UTM zone 36N и задав удобное расположение для нового файла (Tosnensku_roads_utm).
 +
 
 +
[[Файл:Grass_qgis_isochrones_reproject_roads.png|400px|thumb|center|]]
 +
 
 +
Аналогичным образом перепроецируем слой с границами Тосненского района, на всякий случай (получаем дополнительно Tosnensky_boundary_utm). Данные почти готовы! Теперь нужно отфильтровать геоданные по дорогам с учётом требований к моделированию. Это первый тематический шаг подготовки: '''отбор только тех объектов дорожной сети, которые могут быть использованы для перемещения предполагаемыми транспортными средствами'''. Здесь вы должны решить, по объектам каких типов будет разрешено перемещаться, например, для моделирования пешеходного движения можно оставить все типы дорог; для легковых автомобилей исключить пешеходные дорожки, тропы, просеки; для грузовых автомобилей исключить дороги с соответствующими ограничениями на массу, и так далее. Шаг этот очень ответственный. При этом подход к фильтрации сильно зависит от качества исходных данных и их детализации. Вспоминаем, что имеем дело с OSM, где различия между объектами существуют на уровне тегов, которые записаны в таблице атрибутов в поле HIGHWAY. Посмотрим на список уникальных значений для нашего набора данных:
 +
<pre>
 +
path steps footway residential primary service unclassified secondary road track construction tertiary raceway tertiary_link trunk secondary_link proposed pedestrian bridleway primary_link living_street trunk_link
 +
</pre>
 +
 
 +
Все эти теги [http://wiki.openstreetmap.org/wiki/Key:highway описаны на ресурсах OSM], даже с картинками. Нам ещё предстоит вернуться к этим описаниям, а сейчас важно определить, по объектам каких типов может нормально передвигаться легковой автомобиль, а по каким нет. Очевидно, что нужно исключить объекты следующих типов (пригодны только для пешеходов, проектируются, строятся и т.д.):
 +
<pre>
 +
path, steps, footway, construction, proposed, pedestrian, bridleway
 +
</pre>
 +
 
 +
Также много объектов с типом ''unclassified''. При осуществлении серьезной работы, конечно, стоит прояснить природу каждого из таких объектов по вспомогательным данным (и внести информацию в OSM, конечно же). В данном случае оставим их как подходящие для работы. Объекты всех отобранных типов нужно удалить. Для этого переходим в таблицу атрибутов слоя Tosnensky_roads_utm, выбираем инструмент "Выбрать по выражению" и формируем запрос вида:
 +
<pre>
 +
"HIGHWAY" IN (  'path' ,  'steps' ,  'footway' ,  'construction' ,  'proposed' ,  'pedestrian' ,  'bridleway' )
 +
</pre>
 +
 
 +
отмечая в списке, соответственно, выбранные вами типы дорог, которые должны быть проигнорированы. Нажимаем кнопку "выбрать" и видим на карте лишние дороги.
 +
 
 +
[[Файл:Grass_qgis_isochrones_for_deleting.png|600px|thumb|center|]]
 +
 
 +
Переходим в режим редактирования и удаляем выбранные объекты, сохраняем изменения. Если что, все исходные данные остаются в других слоях. Общая подготовка закончена.
 +
 
 +
 
 +
===Определение транспортных характеристик данных о дорожной сети в QGIS===
 +
 
 +
Второй этап подготовки является самым важным и трудоёмким в реальных условиях: необходимо каждому участку дороги назначить некоторую среднюю скорость перемещения по нему. Для примера мы используем простую логику: каждому семейству объектов (по тегу) назначим общую ожидаемую скорость, не вникая в каждый отдельный объект. При серьезном исследовании, опять же, стоит более внимательно отнестись к этому этапу, по возможности используя дополнительные материалы. Итак, рассмотрим те категории, которые у нас остались после фильтрации:
 +
 
 +
* trunk - важнейшие и крупнейшие дороги, например, для нашей территории, Московское шоссе. Ожидаемая скорость 90 км/ч
 +
* primary - крупные шоссе, следующий уровень после trunk. Ожидаемая скорость 90 км/ч
 +
* secondary - относительно крупные дороги, следующий уровень после primary. Ожидаемая скорость 60 км/ч
 +
* tertiary - обычные автомобильные дороги между небольшими населенными пунктами. Ожидаемая скорость 60 км/ч
 +
* living_street - жилые зоны, где у пешеходов явное преимущество в праве передвижения. Ожидаемая скорость 15 км/ч
 +
* residential - автомобильные дороги в жилых кварталах. Ожидаемая скорость около 40 км/ч
 +
* service - сервисные подъезды, въезды и проч. Ожидаемая скорость 30 км/ч
 +
* road - автомобильная дорога неизвестного типа. Примем скорость 60 км/ч
 +
* track - грунтовые дороги, обычно для сельхоз-техники. Примем скорость 30 км/ч
 +
* raceway - дороги для автомобильных видов спорта. Примем скорость 90 км/ч
 +
* tertiary_link - участки, соединяющие tertiary с другими tertiary или дорогами других типов. Примем скорость 40 км/ч
 +
* secondary_link - участки, соединяющие secondary с другими secondary или дорогами других типов. Примем скорость 40 км/ч
 +
* primary_link - участки, соединяющие primary с другими primary или дорогами других типов. Примем скорость 40 км/ч
 +
* trunk_link - участки, соединяющие trunk с другими trunk или дорогами других типов. Примем скорость 40 км/ч
 +
* unclassified - дороги без тега. Примем скорость 40 км/ч
 +
 
 +
Принятые значения достаточно условны, но позволят продемонстрировать принципы дальнейшей работы. Назначаемые скорости довольно существенно зависят от конкретной территории, местного законодательства и так далее. Важно заметить, что у некоторых участков (обычно их немного) в свойствах заполнено значение MAXSPEED, оно может быть как численным (в км/ч), так и вида RU:urban. Расшифровку этих обозначений можно найти, к примеру, [http://wiki.openstreetmap.org/wiki/RU:Key:source:maxspeed здесь] и использовать при определении скоростей. '''В целом этот этап является ключевым, и в боевых условиях здесь важно внимательно и тщательно определить характеристики дорожной сети'''.
 +
 
 +
Для того, чтобы применить принятые характеристики, воспользуемся калькулятором полей, запустив его из таблицы атрибутов слоя ''Tosnensky_roads_utm''. Создадим новый атрибут SPEED как целочисленный и заполним его следующим выражением:
 +
<pre>
 +
CASE
 +
WHEN  "HIGHWAY" = 'trunk' THEN 90
 +
WHEN  "HIGHWAY" = 'primary' THEN 90
 +
WHEN  "HIGHWAY" = 'secondary' THEN 60
 +
WHEN  "HIGHWAY" = 'tertiary' THEN 60
 +
WHEN  "HIGHWAY" = 'living_street' THEN 15
 +
WHEN  "HIGHWAY" = 'residential' THEN 40
 +
WHEN  "HIGHWAY" = 'service' THEN 30
 +
WHEN  "HIGHWAY" = 'road' THEN 60
 +
WHEN  "HIGHWAY" = 'track' THEN 30
 +
WHEN  "HIGHWAY" = 'raceway' THEN 90
 +
WHEN  "HIGHWAY" = 'tertiary_link' THEN 40
 +
WHEN  "HIGHWAY" = 'secondary_link' THEN 40
 +
WHEN  "HIGHWAY" = 'primary_link' THEN 40
 +
WHEN  "HIGHWAY" = 'trunk_link' THEN 40
 +
WHEN  "HIGHWAY" = 'unclassified' THEN 40
 +
END
 +
</pre>
 +
 
 +
Здесь мы проверяем принадлежность каждого объекта слоя к одной из категорий и назначаем скорость в соответствии с определенным ранее значениями. Очень удобно использовать условный оператор CASE.
 +
 
 +
[[Файл:Grass_qgis_isochrones_speed_calc.png|600px|thumb|center|]]
 +
 
 +
После этой операции у каждого объекта слоя есть характеристика скорости. Теперь, при необходимости, можно поправить её для отдельных объектов в таблице атрибутов с учётом изучения реальных характеристик территории. Мы этот шаг пропустим и приступим к следующему - расчёту времени, - которое потребуется для преодоления каждого отдельного участка дорожной сети. Для этого нужно сначала рассчитать длины всех объектов слоя. В калькуляторе атрибутов создаем новое поле LENGTH как десятичное число по формуле
 +
<pre>
 +
$length/1000
 +
</pre>
 +
 
 +
[[Файл:Grass_qgis_isochrones_length_calc.png|400px|thumb|center|]]
 +
 
 +
По умолчанию расчёт длин и площадей производится в единицах проекции, т.е. в нашем случае в метрах (квадратных метрах для площади). Делением на 1000 мы приводим расстояния к километрам. Следующий шаг - расчёт времени (предпочтительно в минутах), которое потребуется для преодоления каждого отдельного участка дороги с учётом его длины и ожидаемой скорости. В калькуляторе атрибутов создаем новое поле TIME как десятичное число по формуле
 +
<pre>
 +
"LENGTH"/"SPEED" * 60
 +
</pre>
 +
 
 +
[[Файл:Grass_qgis_isochrones_time_calc.png|400px|thumb|center|]]
 +
 
 +
 
 +
=== Подготовка данных через командную строку и утилиты GDAL/OGR ===
 +
 
 +
Альтернативным подходом при подготовке данных является использование утилит командной строки. Такой способ удобно использовать, к примеру, для автоматизации пакетной обработки большого количества наборов данных. Для осуществления вышеописанных операций понадобятся две утилиты: [http://www.gdal.org/ogr2ogr.html ogr2ogr] и [http://www.gdal.org/ogrinfo.html ogrinfo]. Если вы работаете в среде Linux, при установленных QGIS и/или GRASS GIS, эти утилиты доступны вам из штатной системной командной строки (терминала). При работе в Windows наиболее простой способ доступа к ним - использование терминала OSGeo4W Shell, который автоматически устанавливается вместе с QGIS. Итак, после запуска терминала в Linux или OSGeo4W Shell в Windows выполняем следующие команды:
 +
 
 +
1. Переходим в директорию, содержащую данные OSM, с помощью команды [https://en.wikipedia.org/wiki/Cd_(command) cd]:
 +
<pre>cd <путь до папки с распакованным архивом>/RU-LEN/data</pre>
 +
 
 +
2. С помощью утилиты ''ogr2ogr'' обрезаем слой дорог по границам Тосненского района, одновременно отфильтровывая ненужные объекты дорожной сети и проецируем результат в UTM 36N:
 +
<pre>ogr2ogr -skipfailures -t_srs "EPSG:32636" Tosnensky_roads_utm.shp  highway-line.shp -clipsrc boundary-polygon.shp -clipsrcwhere "NAME = 'Тосненский район'" -nlt "LINESTRING" -sql "SELECT * FROM \"highway-line\" WHERE OGR_GEOMETRY='LINESTRING'" -sql "SELECT * FROM \"highway-line\" WHERE HIGHWAY NOT IN ('path','steps','footway','construction','proposed','pedestrian','bridleway')" -lco ENCODING=UTF-8</pre>
 +
 
 +
3. Конвертируем набор данных из ESRI Shapefile в SQLite для оптимизации дальнейшей работы:
 +
<pre>ogr2ogr -explodecollections -f "SQLite" -dsco SPATIALITE=YES Tosnensky_roads_utm.sqlite Tosnensky_roads_utm.shp</pre>
 +
 
 +
4. Добавляем в атрибуты колонки SPEED, LENGTH, TIME:
 +
<pre>
 +
ogrinfo  Tosnensky_roads_utm.sqlite -sql "ALTER TABLE tosnensky_roads_utm ADD COLUMN SPEED integer"
 +
ogrinfo  Tosnensky_roads_utm.sqlite -sql "ALTER TABLE tosnensky_roads_utm ADD COLUMN LENGTH real"
 +
ogrinfo  Tosnensky_roads_utm.sqlite -sql "ALTER TABLE tosnensky_roads_utm ADD COLUMN TIME real"
 +
</pre>
 +
 
 +
5. Рассчитываем значения скоростей для объектов исходя из атрибута "highway":
 +
<pre>
 +
ogrinfo  Tosnensky_roads_utm.sqlite -sql "UPDATE tosnensky_roads_utm SET SPEED = CASE WHEN  "highway" = 'trunk' THEN 90 WHEN  "highway" = 'primary' THEN 90 WHEN  "highway" = 'secondary' THEN 60 WHEN  "highway" = 'tertiary' THEN 60 WHEN  "highway" = 'living_street' THEN 15 WHEN  "highway" = 'residential' THEN 40 WHEN  "highway" = 'service' THEN 30 WHEN  "highway" = 'road' THEN 60 WHEN  "highway" = 'track' THEN 30 WHEN  "highway" = 'raceway' THEN 90 WHEN  "highway" = 'tertiary_link' THEN 40 WHEN  "highway" = 'secondary_link' THEN 40 WHEN  "highway" = 'primary_link' THEN 40 WHEN  "highway" = 'trunk_link' THEN 40 WHEN  "highway" = 'unclassified' THEN 40 END"
 +
</pre>
 +
 
 +
6. Рассчитываем длины объектов дорожной сети
 +
<pre>ogrinfo  Tosnensky_roads_utm.sqlite -sql "UPDATE tosnensky_roads_utm SET LENGTH = ST_length(GEOMETRY)/1000"</pre>
 +
 
 +
7. Рассчитываем, исходя из скоростей и длин, затраты времени на преодоление объектов дорожной сети:
 +
<pre>ogrinfo  Tosnensky_roads_utm.sqlite -sql "UPDATE tosnensky_roads_utm SET TIME = LENGTH/SPEED*60"</pre>
 +
 
 +
 
 +
Теперь у нас есть для каждого участка дорожной сети оценка времени, требуемая для его преодоления. На этом подготовка данных завершена, можно приступать к моделированию и переходить в GRASS. Если вы готовили данные в QGIS, не забудьте сохранить изменения в слое и выйти из режима редактирования.