Создаём свою карту ДТП

Загружаем исходные данные

Источник данных: http://stat.gibdd.ru

Находим код региона.

Указываем диапазон дат.

Запускаем скрипт извлечения данных, по годам.

Результат: пачка XML с карточками ДТП.

Преобразовываем в табличный вид

С XML работать не удобно, поэтому используем следующий скрипт для преобразования XML в CSV.

Результат: CSV файл со всеми данными по ДТП.

Находим пересечения

Источник данных: http://data.nextgis.com Совсем свежие данные стоят 200-300р, старые можно скачать бесплатно. Так или иначе получаем свой регион. Нам нужен только один слой: highway-lines.shp (дорожная сеть).

Находим все пересечения и считаем в каждом количество примыканий, определяем таким образом Х-образные пересечения (у них это значение >= 4). Но можете оставить и другие, если нужно.

Скрипт детектирования перекрёстков в данных Openstreetmap

В OSM нет тега для перекрёстка, но их можно получить, взяв те точки, в которых пересекается несколько веев (более двух). Для пробы этот скрипт был написан на питоновских биндингах GDAL, и Москву он считает за ночь. Вероятно в PostGIS этот же алгоритм можно сделать более быстрым.

Если вы работаете на Windows или Mac OS, то установите NextGIS QGIS, найдите и запустите в меню программ "NextGIS command prompt". Если вы работаете на Linux, то установите GDAL с питоновскими биндингами.

Проверяем, что биндинги работают. Выполняем эту команду - она не должна выводить ничего.

python -c "import gdal"

Клонируем скрипт

git clone https://github.com/nextgis/data_processing_scripts.git
cd data_processing_scripts/nextgis_extracts_detect_crossings

При запуске скрипта нужно только указать путь к highway-line.shp, который распаковали на предыдущем шаге и задать фильтр по типам дорог - это то же самое, если бы вы отфильтровали highway-lines.shp в ogr2ogr. Этот фильтр вы задаёте в зависимости от конкретно вашего алгоритма анализа, в зависимости от его настройки можно получить аварии только на перекрёстках крупных улиц, или только на перекрёстках мелких улиц.

python detect_crossings.py -s highway-line.shp -d crossings.shp -f "HIGHWAY IN ('primary')"
python detect_crossings.py -s highway-line.shp -d crossings.shp -f "HIGHWAY IN ('motorway','motorway_link','trunk','trunk_link','primary','primary_link','secondary','secondary_link','tertiary','tertiary_link','unclassified','residential')"

Фильтр нужен обязательно, потому что без него у вас получатся перекрёстки тротуаров с тропинками в лесу

Скрипт создаст Shapefile с точками, в атрибуте WAYS_CNT будет количество дорог, подходящих к точке.

Создаём перекрестки

Перекресток - зона определенного радиуса вокруг пересечения (буферная зона).

Перекрестки получаются путем построения буферных зон вокруг пересечений. Соединенные при этом зоны - объединяются в единый объект.

На входе - точки (crossings), точки событий (dtp2014).

Расчёт делается в PostGIS. Хорошая инструкция по установке есть на https://www.digitalocean.com/community/tutorials/how-to-install-and-configure-postgis-on-ubuntu-14-04

Загружаем слои из файлов в PostGIS. Подставьте в эти команды логин и пароль к БД.

ogr2ogr -progress -overwrite PG:"host=localhost dbname=gis user= password=" ../moscow_crossings.gpkg

ogr2ogr -progress -overwrite -f GPKG dtp2014.gpkg --config SHAPE_ENCODING cp1251  dtp_post_14_all_Mow.shp -skipfailures 
ogr2ogr -progress -overwrite PG:"host=localhost dbname=gis user= password=" dtp2014.gpkg -skipfailures  -nln dtp2014

ogr2ogr -progress -overwrite -f GPKG dtp2015.gpkg --config SHAPE_ENCODING cp1251  dtp_post_15_all_Mow.shp -skipfailures 
ogr2ogr -progress -overwrite PG:"host=localhost dbname=gis user= password=" dtp2015.gpkg -skipfailures  -nln dtp2015

ogr2ogr -progress -overwrite -f GPKG dtp2016.gpkg --config SHAPE_ENCODING cp1251  dtp_post_16_all_Mow.shp -skipfailures 
ogr2ogr -progress -overwrite PG:"host=localhost dbname=gis user= password=" dtp2016.gpkg -skipfailures  -nln dtp2016

Для перекрёстков строится 3 буфера, они объединяются. Делается spatial self join, и подсчитывается количество и разные суммы атрибутов аварий у перекрёстков. На выходе - полигональный слой с кружками вокруг пеерекрёстков и разными атрибутами.

Подключитесь к БД в pgAdmin или Valentina Studio, выполните этот запрос. Если вы работаете не с Москвой, то вместо 32637 подставьте номер системы координат UTM на нужную местность.

-- генерируем три слоя буферов с разным расстоянием
-- Укажите тут нужную СК, если это не Москва.
CREATE TEMPORARY TABLE buffers100 ON COMMIT DROP AS 
SELECT st_buffer(ST_Transform(wkb_geometry,32637),100) AS wkb_geometry,100 AS buffer_size FROM moscow_crossings WHERE ways_cnt = 4;
CREATE TEMPORARY TABLE buffers50 ON COMMIT DROP AS 
SELECT st_buffer(ST_Transform(wkb_geometry,32637),50) AS wkb_geometry,50 AS buffer_size FROM moscow_crossings WHERE ways_cnt = 4;
CREATE TEMPORARY TABLE buffers20 ON COMMIT DROP AS 
SELECT st_buffer(ST_Transform(wkb_geometry,32637),20) AS wkb_geometry,20 AS buffer_size FROM moscow_crossings WHERE ways_cnt = 4;



DROP TABLE IF EXISTS crossingsarea_disolved;
CREATE TABLE crossingsarea_disolved ( 
    wkb_geometry geometry,
    buffer integer);
    
--для 100 метров
DROP TABLE IF EXISTS crossingsarea_disolved_onebuffer ;
CREATE 
TEMPORARY 
TABLE crossingsarea_disolved_onebuffer 
ON COMMIT DROP 
AS 
WITH
 clusters(wkb_geometry) AS 
      (SELECT ST_CollectionExtract(unnest(ST_ClusterIntersecting(wkb_geometry)), 3) 
         FROM buffers100),
 multis(id, wkb_geometry) AS 
      (SELECT row_number() over() as id, wkb_geometry FROM clusters)
 SELECT ST_UNION(wkb_geometry) AS wkb_geometry FROM 
      (SELECT id, (ST_DUMP(wkb_geometry)).geom AS wkb_geometry FROM multis) d GROUP BY id;
INSERT INTO  crossingsarea_disolved SELECT wkb_geometry, 100 as buffer FROM crossingsarea_disolved_onebuffer;  

--для 50 метров
DROP TABLE IF EXISTS crossingsarea_disolved_onebuffer ;
CREATE 
TEMPORARY 
TABLE crossingsarea_disolved_onebuffer 
ON COMMIT DROP 
AS 
WITH
 clusters(wkb_geometry) AS 
      (SELECT ST_CollectionExtract(unnest(ST_ClusterIntersecting(wkb_geometry)), 3) 
         FROM buffers50),
 multis(id, wkb_geometry) AS 
      (SELECT row_number() over() as id, wkb_geometry FROM clusters)
 SELECT ST_UNION(wkb_geometry) AS wkb_geometry FROM 
      (SELECT id, (ST_DUMP(wkb_geometry)).geom AS wkb_geometry FROM multis) d GROUP BY id;
INSERT INTO  crossingsarea_disolved SELECT wkb_geometry, 50 as buffer FROM crossingsarea_disolved_onebuffer;   

--для 20 метров
DROP TABLE IF EXISTS crossingsarea_disolved_onebuffer ;
CREATE 
TEMPORARY 
TABLE crossingsarea_disolved_onebuffer 
ON COMMIT DROP 
AS 
WITH
 clusters(wkb_geometry) AS 
      (SELECT ST_CollectionExtract(unnest(ST_ClusterIntersecting(wkb_geometry)), 3) 
         FROM buffers20),
 multis(id, wkb_geometry) AS 
      (SELECT row_number() over() as id, wkb_geometry FROM clusters)
 SELECT ST_UNION(wkb_geometry) AS wkb_geometry FROM 
      (SELECT id, (ST_DUMP(wkb_geometry)).geom AS wkb_geometry FROM multis) d GROUP BY id;
INSERT INTO  crossingsarea_disolved SELECT wkb_geometry, 20 as buffer FROM crossingsarea_disolved_onebuffer;   


--добавление столбцов и рассчёт числовых значений
ALTER TABLE crossingsarea_disolved ADD COLUMN count_2014 integer;
ALTER TABLE crossingsarea_disolved ADD COLUMN count_2015 integer;
ALTER TABLE crossingsarea_disolved ADD COLUMN count_2016 integer;

ALTER TABLE crossingsarea_disolved ADD COLUMN cp_2014 integer;
ALTER TABLE crossingsarea_disolved ADD COLUMN cp_2015 integer;
ALTER TABLE crossingsarea_disolved ADD COLUMN cp_2016 integer;

ALTER TABLE crossingsarea_disolved ADD COLUMN cps_2014 integer;
ALTER TABLE crossingsarea_disolved ADD COLUMN cps_2015 integer;
ALTER TABLE crossingsarea_disolved ADD COLUMN cps_2016 integer;

ALTER TABLE crossingsarea_disolved ADD COLUMN cp_sum integer;
COMMENT ON COLUMN crossingsarea_disolved.cp_sum IS 'КОЛИЧЕСТВО АВАРИЙ С ТРУПАМИ';
ALTER TABLE crossingsarea_disolved ADD COLUMN cps_sum integer;
COMMENT ON COLUMN crossingsarea_disolved.cp_sum IS 'КОЛИЧЕСТВО ТРУПОВ';

UPDATE crossingsarea_disolved SET wkb_geometry=ST_SetSRID(wkb_geometry,32637);
UPDATE crossingsarea_disolved SET count_2014=0;
UPDATE crossingsarea_disolved SET count_2015=0;
UPDATE crossingsarea_disolved SET count_2016=0;
UPDATE crossingsarea_disolved SET cp_2014=0;
UPDATE crossingsarea_disolved SET cp_2015=0;
UPDATE crossingsarea_disolved SET cp_2016=0;
UPDATE crossingsarea_disolved SET cps_2014=0;
UPDATE crossingsarea_disolved SET cps_2015=0;
UPDATE crossingsarea_disolved SET cps_2016=0;

CREATE TEMPORARY TABLE medium2014 ON COMMIT DROP AS 
SELECT
dtp2014.wkb_geometry, crossingsarea_disolved.wkb_geometry AS crossing_geometry,corpse
FROM 
dtp2014 , crossingsarea_disolved
WHERE 
ST_Within(dtp2014.wkb_geometry, crossingsarea_disolved.wkb_geometry);


CREATE TEMPORARY TABLE medium2014_cnt ON COMMIT DROP AS 
SELECT medium2014.crossing_geometry, COUNT(*) AS cnt  FROM medium2014 GROUP BY crossing_geometry;

UPDATE crossingsarea_disolved SET count_2014 = medium2014_cnt.cnt
FROM medium2014_cnt
WHERE medium2014_cnt.crossing_geometry = crossingsarea_disolved.wkb_geometry;
DROP TABLE medium2014_cnt;

CREATE TEMPORARY TABLE medium2014_cnt ON COMMIT DROP AS 
SELECT medium2014.crossing_geometry, COUNT(*) AS cnt, SUM(corpse) AS corpses FROM medium2014 WHERE corpse > 0 GROUP BY crossing_geometry;

UPDATE crossingsarea_disolved SET cp_2014 = medium2014_cnt.cnt, cps_2014 = medium2014_cnt.corpses
FROM medium2014_cnt
WHERE medium2014_cnt.crossing_geometry = crossingsarea_disolved.wkb_geometry;
DROP TABLE medium2014_cnt;


CREATE TEMPORARY TABLE medium2015 ON COMMIT DROP AS 
SELECT
dtp2015.wkb_geometry, crossingsarea_disolved.wkb_geometry AS crossing_geometry,corpse
FROM 
dtp2015 , crossingsarea_disolved
WHERE 
ST_Within(dtp2015.wkb_geometry, crossingsarea_disolved.wkb_geometry);

CREATE TEMPORARY TABLE medium2015_cnt ON COMMIT DROP AS 
SELECT medium2015.crossing_geometry, COUNT(*) AS cnt FROM medium2015 GROUP BY crossing_geometry;

UPDATE crossingsarea_disolved SET count_2015 = medium2015_cnt.cnt
FROM medium2015_cnt
WHERE medium2015_cnt.crossing_geometry = crossingsarea_disolved.wkb_geometry;
DROP TABLE medium2015_cnt;

CREATE TEMPORARY TABLE medium2015_cnt ON COMMIT DROP AS 
SELECT medium2015.crossing_geometry, COUNT(*) AS cnt, SUM(corpse) AS corpses FROM medium2015 WHERE corpse > 0 GROUP BY crossing_geometry;

UPDATE crossingsarea_disolved SET cp_2015 = medium2015_cnt.cnt, cps_2015 = medium2015_cnt.corpses
FROM medium2015_cnt
WHERE medium2015_cnt.crossing_geometry = crossingsarea_disolved.wkb_geometry;
DROP TABLE medium2015_cnt;


CREATE TEMPORARY TABLE medium2016 ON COMMIT DROP AS 
SELECT
dtp2016.wkb_geometry, crossingsarea_disolved.wkb_geometry AS crossing_geometry,corpse
FROM 
dtp2016 , crossingsarea_disolved
WHERE 
ST_Within(dtp2016.wkb_geometry, crossingsarea_disolved.wkb_geometry);

CREATE TEMPORARY TABLE medium2016_cnt ON COMMIT DROP AS 
SELECT medium2016.crossing_geometry, COUNT(*) AS cnt FROM medium2016 GROUP BY crossing_geometry;

UPDATE crossingsarea_disolved SET count_2016 = medium2016_cnt.cnt
FROM medium2016_cnt
WHERE medium2016_cnt.crossing_geometry = crossingsarea_disolved.wkb_geometry;
DROP TABLE medium2016_cnt;

CREATE TEMPORARY TABLE medium2016_cnt ON COMMIT DROP AS 
SELECT medium2016.crossing_geometry, COUNT(*) AS cnt, SUM(corpse) AS corpses FROM medium2016 WHERE corpse > 0 GROUP BY crossing_geometry;

UPDATE crossingsarea_disolved SET cp_2016 = medium2016_cnt.cnt, cps_2016 = medium2016_cnt.corpses
FROM medium2016_cnt
WHERE medium2016_cnt.crossing_geometry = crossingsarea_disolved.wkb_geometry;
DROP TABLE medium2016_cnt;


UPDATE crossingsarea_disolved SET cp_sum = cp_2014+cp_2015+cp_2016;
UPDATE crossingsarea_disolved SET cps_sum = cps_2014+cps_2015+cps_2016;

Выгружаем слой в файл

ogr2ogr -progress -overwrite -f GPKG crossings_dtp.gpkg PG:"host=localhost dbname=gis " crossingsarea_disolved  -skipfailures