Postgis实现最短路径查询（Dijkstra算法）

最新推荐文章于 2024-07-23 22:06:37 发布

就爱安逸啊

最新推荐文章于 2024-07-23 22:06:37 发布

阅读量1.1k

点赞数 2

分类专栏：空间数据库文章标签：算法数据库 postgresql 学习 javascript

本文链接：https://blog.csdn.net/Webgisera/article/details/131876866

版权

空间数据库专栏收录该内容

1 篇文章 0 订阅

订阅专栏

文章目录

前言

前言

这篇文章是我的webgis系统中的最后一个功能，也是我花了好长时间才明白一点的。内容借鉴博主：不睡觉的怪叔叔。

GIS数据处理

首先获取了osm的shp数据，在Arcgis中将其转化为web墨卡托投影（EPSG：3857），并使用gis工具中的要素转线，可以将相交的折线在交点处打断。
在这里插入图片描述

数据导入

现在pgAdmin中建立一个数据库，并安装postgis和pgrouting插件。
sql：create extension postgis 和 create extension pgrouting
测试是否安装成功：select postgis_full_version(),pg_version()
最后将数据导入对应的数据库，可以通过postgis自带的用来导入shp数据的工具。
注意设置SRID：3857

设置路径成本

首先为shenzhen_roads表（右键点击查询工具）添加四个字段（后面涉及到对应字段时再详细解释）：

source —— 用于保存路径起始顶点的id
target —— 用于保存路径终止顶点的id
cost —— 用于保存路径正向的成本（或者代价）
reverse_cost —— 用于保存路径反向的成本（或者代价）
SQL语句：

ALTER TABLE shenzhen_roads
ADD COLUMN source INTEGER,
ADD COLUMN target INTEGER,
ADD COLUMN cost DOUBLE PRECISION,
ADD COLUMN reverse_cost DOUBLE PRECISION;

OSM数据包含一个oneway列，它的值可能是以下三个值之一：

‘F’ —— 表示该路径是单向的，且路径方向是正向的。
‘T’ —— 表示该路径是单向的，且路径方向是反向的。
‘B’ —— 表示该路径是双向的。
在这里插入图片描述
后续是根据路径方向计算路径成本以及创建路网拓：将路径的长度作为成本，如果路径是正向的，则cost值为路径的长度，reverse_cost值为-1；如果路径是反向的，则reverse_cost值为路径的长度，cost值为-1；如果路径是双向的，则cost和reverse_cost的值都为路径的长度。
路网拓朴：

SELECT pgr_createTopology(
	'shenzhen_roads', 
	0.001,
	'geom',
	'gid',
	'source',
	'target'
);

上面的六个参数分别表示：

路网表
路径之间的容差，两条路径的距离大于这个容差值，就表示它们不相交，否则就是相交。
路网表中包含空间信息的列
路网表的主码列
保存路径起始顶点的id的列
保存路径终止顶点的id的列

之后多处一张表shenzhen_roads_vertices_pgr
shenzhen_roads中source和target的值即为这里保存的id

路径优先级

这是无优先级
现在添加一个penalty列用于保存每条路径的优先级信息（penalty值越小，优先级越高）：

ALTER TABLE shenzhen_roads ADD COLUMN penalty DOUBLE PRECISION;

把每条路径的penalty的值都设为1.0，表示所有路径都没有优先级之分（penalty值都为1.0）。

UPDATE shenzhen_roads SET penalty = 1.0;

有优先级
车辆路径可以让路径的优先级符合以下规则：
禁止在pedestrian、footway、steps、cycleway、track这些路径上行驶。
不鼓励在rsidential路径上行驶。
鼓励在primary、primary_link路径上行驶。

修改路径的优先级具体如下：

UPDATE shenzhen_roads SET penalty = -1.0
WHERE fclass IN ('steps','footway','pedestrian', 'cycleway', 'track', 'track_grade2');

UPDATE shenzhen_roads SET penalty = 5.0
WHERE fclass IN ('residential');

UPDATE shenzhen_roads SET penalty = 0.8
WHERE fclass IN ('tertiary', 'tertiary_link', 'motorway', 'motorway_link', 'living_street');
UPDATE shenzhen_roads SET penalty = 0.5
WHERE fclass IN ('secondary', 'secondary_link');

UPDATE shenzhen_roads SET penalty = 0.3
WHERE fclass IN ('primary', 'primary_link', 'trunk', 'trunk_link');

这时候再进行查询就会得到不同的路径成本

SELECT * FROM pgr_dijkstra(
	'SELECT gid AS id,
		source, target,
		cost * penalty AS cost, 
	    reverse_cost
	 FROM shenzhen_roads
	',
	50013, 19688, 
	directed := true
);

创建数据库视图

车辆不允许在行人道路上行驶，因此：

创建包含车辆可以行驶的道路的数据库视图。
道路成本设置为时间（单位分钟）。
验证以上过程减少了多少路段。

```CREATE VIEW vehicle_net AS
	SELECT gid, source, target,
		   cost / 1.3 / 60 AS cost,
		   reverse_cost / 1.3 / 60 AS reverse_cost, -- 假设行人每秒走1.3m
		   geom
	FROM shenzhen_roads 
	WHERE fclass NOT IN ('steps', 'cycleway', 'footway', 'track_grade2', 'track');
	
-- 验证
SELECT count(*) FROM shenzhen_roads;
SELECT count(*) FROM vehicle_net;

在这里插入图片描述

将车辆限制在只能再一定矩形区域中行驶，该区域的范围是：（12673569, 2573130, 12680410, 2578570）坐标系ESPG:3857。

车辆只能在区域（12673569, 2573130, 12680410, 2578570）中行驶。
创建一个保存该区域路段的数据库视图。
验证缺少的路段数

CREATE VIEW little_net AS
	SELECT *
	FROM vehicle_net
	WHERE vehicle_net.geom && ST_MakeEnvelope(12673569, 2573130, 12680410, 2578570, 3857);

-- 验证
SELECT COUNT(*) FROM little_net;

在这里插入图片描述
ST_MakeEnvelope创建一个矩形多边形用于筛选数据

创建一个函数

我们创建的函数具有如下特点：

可以在车辆道路的全部区域进行搜索。
将数据库视图作为参数。
结果包含道路名、路径的几何信息（二进制格式）、路径各个路段的大地方位角（单位十进制度）等等。

CREATE OR REPLACE FUNCTION wrk_dijkstra(
    IN edges_subset regclass,	-- 视图作为参数
    IN source BIGINT,
    IN target BIGINT,
    OUT seq INTEGER,
    OUT gid BIGINT,
    OUT name TEXT,
    OUT cost FLOAT,
    OUT azimuth FLOAT,
    OUT route_readable TEXT,
    OUT route_geom geometry
) RETURNS SETOF record AS 
$BODY$
    WITH
    dijkstra AS (
        SELECT * FROM pgr_dijkstra(
            -- 使用参数化的视图
            'SELECT gid AS id, * FROM ' || $1,
            $2, $3)
    ),
    get_geom AS (
        SELECT dijkstra.*, ways.name,
            CASE
                WHEN dijkstra.node = ways.source THEN geom
                ELSE ST_Reverse(geom)
            END AS route_geom
        FROM dijkstra JOIN shenzhen_roads AS ways ON (edge = gid)
        ORDER BY seq)
    SELECT
        seq,
        edge,
        name,
        cost,
        degrees(ST_azimuth(ST_StartPoint(route_geom), ST_EndPoint(route_geom))) AS azimuth,
        ST_AsText(route_geom),
        route_geom
    FROM get_geom
    ORDER BY seq;
$BODY$
LANGUAGE 'sql';

创建函数wrk_fromAtoB：

CREATE OR REPLACE FUNCTION wrk_fromAtoB(
	IN edges_subset regclass,
	IN x1 NUMERIC, IN y1 numeric,
	IN x2 NUMERIC, IN y2 NUMERIC,
	OUT seq INTEGER,
	OUT gid BIGINT,
	OUT name TEXT,
	OUT costs FLOAT,
	OUT azimuth FLOAT,
	OUT geom GEOMETRY
)
RETURNS SETOF record AS
$BODY$
DECLARE 
	final_query TEXT;
BEGIN
	final_query := 
		FORMAT($$
			WITH
			vertices AS (
				SELECT * FROM shenzhen_roads_vertices_pgr
				WHERE id IN (
					SELECT source FROM %1$I
					UNION
					SELECT target FROM %1$I
				)
			),
			dijkstra AS (
				SELECT * 
				FROM wrk_dijkstra(
					'%1$I',
					-- source
					(SELECT id FROM vertices
                        ORDER BY the_geom <-> ST_SetSRID(ST_Point(%2$s, %3$s), 3857) LIMIT 1),
					-- target
					(SELECT id FROM vertices
                        ORDER BY the_geom <-> ST_SetSRID(ST_Point(%4$s, %5$s), 3857) LIMIT 1)
				)
			)
			SELECT 
			   seq,
			   dijkstra.gid,
			   dijkstra.name,
			   dijkstra.cost / 1000.0 AS costs,
			   azimuth,
			   route_geom AS geom
			FROM dijkstra 
			JOIN shenzhen_roads ways USING(gid);$$,
		edges_subset, x1,y1,x2,y2
		);
	RAISE notice '%', final_query; -- 执行该函数时显示函数的逻辑代码信息
	RETURN QUERY EXECUTE final_query;
END;
$BODY$
LANGUAGE 'plpgsql';

此函数作用即为计算从一个坐标点到另一个坐标点的路径信息

Geoserver发布接口

首先连接postgis数据库，点击postgis数据源将数据导入，之后点击图层。
在创建新的SQL视图页面，填写视图名称为"shenzhen"，并填写如下SQL语句：

SELECT ST_MakeLine(route.geom) FROM (
    SELECT geom FROM wrk_fromAtoB('vehicle_net', %x1%, %y1%, %x2%, %y2%) ORDER BY seq
) AS route;

以上SQL语句中，使用两个百分号表示中间的变量是参数。

然后再点击"从SQL猜想的参数"，页面中就会显示两个百分号包含的参数，设置四个参数的默认值为0，并且配置"验证的正则表达式"都为"^-?[\d.]+$"，这样参数只能接受数字值（为了防止SQL注入）。

接着再点击属性下面的"刷新"按钮，将会出现一个名称为st_makeline的属性，设置该属性的类型为LineString，设置该属性的SRID为3857。
最后点击"保存"按钮，保存我们新配置的SQL视图。
这样就会返回到编辑图层页面，我们需要在该页面设置"边框"参数，点击从数据计算和computed…
打开浏览器，发送一个WFS的GetFeature请求：

http://localhost:8080/geoserver/shenzhen/ows?service=WFS&version=1.0.0&request=GetFeature&typeName=shenzhen:shenzhen&maxFeatures=50&outputFormat=application/json&viewparams=x1:12677354.9;y1:2578172.3;
	x2:12677441.2;y2:2577908.3

在这里插入图片描述
返回一个geojson文件表示成功！！！

使用openlayers实现路径规划可视化

  import Feature from 'ol/Feature'
  import VectorLayer from 'ol/layer/Vector';
  import VectorSource from 'ol/source/Vector';
  import {
  	map
  } from './baseMap.js';
  import ImageWMS from 'ol/source/ImageWMS';
  import {
  	Image as ImageLayer,
  	Tile as TileLayer
  } from 'ol/layer.js';
  import {
  	Point
  } from 'ol/geom';
  // 初始化起始点要素和目标点要素
  var startPoint = new Feature();
  var destPoint = new Feature();

  // 用于包含起始点要素和目标点要素的图层
  const vectorLayer = new VectorLayer({
  	source: new VectorSource({
  		features: [startPoint, destPoint]
  	})
  });
  map.addLayer(vectorLayer);

  // 注册一个鼠标点击事件，用户点击地图时就会触发
  var result = null;
  let startButton = document.getElementById('one_start')
  startButton.addEventListener('click',()=>{
		  map.on('click', function(event) {
		  	if (startPoint.getGeometry() == null) {
		  		// 设置起始点要素的坐标信息
		  		startPoint.setGeometry(new Point(event.coordinate));
		  	} else if (destPoint.getGeometry() == null) {
		  		// 设置目标点要素的坐标信息
		  		destPoint.setGeometry(new Point(event.coordinate));
		  		var startCoord = startPoint.getGeometry().getCoordinates();
		  		var destCoord = destPoint.getGeometry().getCoordinates();
		  		// 设置GeoServer的SQL视图的请求参数
		  		var viewparams = [
		  			'x1:' + startCoord[0], 'y1:' + startCoord[1],
		  			'x2:' + destCoord[0], 'y2:' + destCoord[1]
		  		];
		  		viewparams = viewparams.join(';');
		  		// 向GeoServer发送WMS请求，并将结果路径渲染出来
		  		result = new ImageLayer({
		  			source: new ImageWMS({
		  				url: 'http://localhost:8080/geoserver/shenzhen/wms',
		  				params: {
		  					LAYERS: 'shenzhen:shenzhen',
		  					FORMAT: 'image/png',
		  					viewparams: viewparams
		  				}
		  			})
		  		});
		  		map.addLayer(result);
		  	}
		  });
	 
  })
  // 用户点击clear按钮，就会触发事件
  const clearButton = document.getElementById('one_clear');
  clearButton.addEventListener('click', function(event) {
  	// 将起始点要素和目标点要素的坐标信息清空
  	startPoint.setGeometry(null);
  	destPoint.setGeometry(null);
  	// 移除结果路径图层
  	map.removeLayer(result);
  });
  export default {
  	startButton,clearButton
  }

总结

具体的实现细节还需要仔细看文章开头那位博主写的，还有一些练习很不错。最后还想说暑假回到家后和家人商量了一下要准备考研，后续就是等秋招投简历找工作以及同时备考，暂时就不会再更新了，能考上研是最好，考不上的话找到工作就去了，反正也是备考的计算机，学的专业知识也有用。好了就这样吧。🤭

就爱安逸啊

关注

2
点赞
踩
10

收藏

觉得还不错? 一键收藏
1
评论
Postgis实现最短路径查询（Dijkstra算法）

这篇文章是我的webgis系统中的最后一个功能，也是我花了好长时间才明白一点的。不睡觉的怪叔叔。具体的实现细节还需要仔细看文章开头那位博主写的，还有一些练习很不错。最后还想说暑假回到家后和家人商量了一下要准备考研，后续就是等秋招投简历找工作以及同时备考，暂时就不会再更新了，能考上研是最好，考不上的话找到工作就去了，反正也是备考的计算机，学的专业知识也有用。好了就这样吧。🤭。
复制链接

扫一扫

专栏目录