Path planning processing method based on pgrouting

For GIS business, path planning is a very basic business. Generally, if a company handles it, it will directly choose to call mature third-party interfaces, such as Gaode, Baidu, etc. Of course, there are actually many algorithms for path planning, such as the more famous Dijkstra, A* algorithm, etc. Of course, this article does not introduce the algorithm. The author of this article will use the Dijkstra algorithm that has been integrated with the postgresql database to process the shortest path.

1. Data processing

The core of path planning is data. The data is general road network data. However, after we get the road network data, we need to process the data. Since the idea of the algorithm is based on the principle of directed graphs, we first need to top-process the data. Through top-based we actually establish the vertex relationship between each road in the road network. The following are the main commands:

--Open the execution networktopoPlugins
create extension postgis;
create extension postgis_topology;
--Data creation topology
ALTER TABLE test_road ADD COLUMN source integer;
ALTER TABLE test_road ADD COLUMN target integer;
SELECT pgr_createTopology('test_road',0.00001, 'geom', 'gid');

where test_road is the table name that imports road network data into postgresql.

After processing the topo, it is basically enough. We can use the functions that come with pgrouting. In fact, there are many. We choose pgr_dijkstra

CREATE OR REPLACE FUNCTION public.pgr_dijkstra(
    IN edges_sql text,
    IN start_vid bigint,
    IN end_vid bigint,
    IN directed boolean,
    OUT seq integer,
    OUT path_seq integer,
    OUT node bigint,
    OUT edge bigint,
    OUT cost double precision,
    OUT agg_cost double precision)
  RETURNS SETOF record AS
$BODY$
DECLARE
BEGIN
    RETURN query SELECT *
    FROM _pgr_dijkstra(_pgr_get_statement($1), start_vid, end_vid, directed, false);
  END
$BODY$
  LANGUAGE plpgsql VOLATILE
  COST 100
  ROWS 1000;
ALTER FUNCTION public.pgr_dijkstra(text, bigint, bigint, boolean)
  OWNER TO postgres;

From the function input parameters, we can see that we need to query SQL, a starting point, an end point, and whether the direction is considered. After understanding the input parameters of the calling function, we will write this function.

2. Principle analysis

Generally, path planning basically enters a starting point position and an end point position and then directly plan. For us, if we want to apply the above function, we must find the starting point position target and the source of the end point position, and then plan based on the two top points found, call the above function to return the required results.

How to find the corresponding target based on the starting point? In fact, it is to find the target closest line to the starting point. Similarly, the source of the end point is to find the source closest line to the end point. Of course, after planning these two points, it is basically enough. However, in the end, you need to connect the target from the starting point to the nearest line of the starting point, and connect the source from the end point to the nearest line of the end point, so that the entire path planning is completed.

Let's take a look at the specific implementation stored procedures:

CREATE OR REPLACE FUNCTION public.pgr_shortest_road(
IN startx double precision,
IN starty double precision,
IN endx double precision,
IN endy double precision,
OUT road_name character varying,
OUT v_shpath character varying,
OUT cost double precision)
RETURNS SETOF record AS
$BODY$ 
declare 
v_startLine geometry;--The closest line to the starting point 
v_endLine geometry;--The closest line to the end 
v_startTarget integer;--The end point of the line closest to the starting point 
v_endSource integer;--The starting point of the closest line to the end 
v_statpoint geometry;--existv_startLineUp the closest point to the starting point 
v_endpoint geometry;--existv_endLineUp the closest point to the end 
v_res geometry;--Shortest path analysis results 
v_perStart float;--v_statpointexistv_resPercentage on 
v_perEnd float;--v_endpointexistv_resPercentage on 
v_rec record; 
first_name varchar;
end_name varchar;
first_cost double precision;
end_cost double precision;
begin 
--查询The closest line to the starting point 
execute 'select geom,target,name from china_road where 
ST_DWithin(geom,ST_Geometryfromtext(''point('|| startx ||' ' || starty||')''),0.01) 
order by ST_Distance(geom,ST_GeometryFromText(''point('|| startx ||' '|| starty ||')'')) limit 1' 
into v_startLine ,v_startTarget,first_name; 
--查询The closest line to the end 
execute 'select geom,source,name from china_road
where ST_DWithin(geom,ST_Geometryfromtext(''point('|| endx || ' ' || endy ||')''),0.01) 
order by ST_Distance(geom,ST_GeometryFromText(''point('|| endx ||' ' || endy ||')'')) limit 1' 
into v_endLine,v_endSource,end_name; 
--If the nearest line is not found，Just returnnull 
if (v_startLine is null) or (v_endLine is null) then 
return; 
end if ; 
select ST_ClosestPoint(v_startLine, ST_Geometryfromtext('point('|| startx ||' ' || starty ||')')) into v_statpoint; 
select ST_ClosestPoint(v_endLine, ST_GeometryFromText('point('|| endx ||' ' || endy ||')')) into v_endpoint; 

--计算距离起点最近线上的点exist该线中的位置
select ST_Line_Locate_Point(st_linemerge(v_startLine), v_statpoint) into v_perStart;

select ST_Line_Locate_Point(st_linemerge(v_endLine), v_endpoint) into v_perEnd;

select ST_Distance_Sphere(v_statpoint,ST_PointN(ST_GeometryN(v_startLine,1), ST_NumPoints(ST_GeometryN(v_startLine,1)))) into first_cost;

select ST_Distance_Sphere(ST_PointN(ST_GeometryN(v_endLine,1),1),v_endpoint) into end_cost; 

if (ST_Intersects(st_geomfromtext('point('|| startx ||' '|| starty ||') '), v_startLine) and ST_Intersects(st_geomfromtext('point('|| endx ||' '|| endy ||') '), v_startLine)) then 
select ST_Distance_Sphere(v_statpoint, v_endpoint) into first_cost;

select ST_Line_Locate_Point(st_linemerge(v_startLine), v_endpoint) into v_perEnd;
for v_rec in 
select ST_Line_SubString(st_linemerge(v_startLine), v_perStart,v_perEnd) as point,COALESCE(end_name,'No Named Road') as name,end_cost as cost loop
v_shPath:= ST_AsGeoJSON(v_rec.point);
cost:= v_rec.cost;
road_name:= v_rec.name;
return next;
end loop;
return;
end if;
--Shortest path 
for v_rec in 
(select ST_Line_SubString(st_linemerge(v_startLine),v_perStart,1) as point,COALESCE(first_name,'No Named Road') as name,first_cost as cost
union all
SELECT st_linemerge() as point,COALESCE(,'No Named Road') as name, as cost
FROM pgr_dijkstra(
'SELECT gid as id, source, target, length as cost FROM china_road
where st_intersects(geom,st_buffer(st_linefromtext(''linestring('||startx||' ' || starty ||','|| endx ||' ' || endy ||')''),0.05))', 
v_startTarget, v_endSource , false 
) a, china_road b 
WHERE  = 
union all
select ST_Line_SubString(st_linemerge(v_endLine),0,v_perEnd) as point,COALESCE(end_name,'No Named Road') as name,end_cost as cost)
loop
v_shPath:= ST_AsGeoJSON(v_rec.point);
cost:= v_rec.cost;
road_name:= v_rec.name;
return next;
end loop; 
end; 
$BODY$
LANGUAGE plpgsql VOLATILE STRICT;

The above implementation is to return all query roads to a collection, and then the client merges each line. This method has a great impact on the final efficiency, so generally the roads will be merged into one road in the function. We can use the postgis st_union function to handle it. After a long time of experimentation, the editor has made many optimizations to the above stored procedures while ensuring efficiency and accuracy, and finally came to the following:

CREATE OR REPLACE FUNCTION public.pgr_shortest_road(
    startx double precision,
    starty double precision,
    endx double precision,
    endy double precision)
  RETURNS geometry AS
$BODY$   
declare  
    v_startLine geometry;--The closest line to the starting point  
    v_endLine geometry;--The closest line to the end 
    v_perStart float;--v_statpointexistv_resPercentage on  
    v_perEnd float;--v_endpointexistv_resPercentage on  
    v_shpath geometry;
    distance double precision;
    bufferInstance double precision;
    bufferArray double precision[];  
begin  
    execute 'select geom,
        case china_road.direction
        when ''3'' then
        source
        else
        target
        end 
        from china_road where  
        ST_DWithin(geom,ST_Geometryfromtext(''point('|| startx ||' ' || starty||')'',4326),0.05) 
        AND width::double precision &gt;= '||roadWidth||'
        order by ST_Distance(geom,ST_GeometryFromText(''point('|| startx ||' '|| starty ||')'',4326))  limit 1'  
    into v_startLine; 
     
    execute 'select geom,
        case china_road.direction
        when ''3'' then
        target
        else
        source
        end 
        from china_road
        where ST_DWithin(geom,ST_Geometryfromtext(''point('|| endx || ' ' || endy ||')'',4326),0.05)  
        AND width::double precision &gt;= '||roadWidth||'
        order by ST_Distance(geom,ST_GeometryFromText(''point('|| endx ||' ' || endy ||')'',4326))  limit 1'  
    into v_endLine;  
    
    if (v_startLine is null) or (v_endLine is null) then  
        return null;
    end if; 
         
    if (ST_equals(v_startLine,v_endLine)) then 
        select ST_LineLocatePoint(st_linemerge(v_startLine), ST_Geometryfromtext('point('|| startx ||' ' || starty ||')',4326)) into v_perStart;
        select ST_LineLocatePoint(st_linemerge(v_endLine), ST_Geometryfromtext('point('|| endx ||' ' || endy ||')',4326)) into v_perEnd;
        select ST_LineSubstring(st_linemerge(v_startLine),v_perStart,v_perEnd) into v_shPath;
        return v_shPath;
    end if;
    
    select ST_DistanceSphere(st_geomfromtext('point('|| startx ||' ' || starty ||')',4326),st_geomfromtext('point('|| endx ||' ' || endy ||')',4326)) into distance;
    if ((distance / 1000) &gt; 50) then
        bufferArray := ARRAY[0.1,0.2,0.3,0.5,0.8];    
    else
        bufferArray := ARRAY[0.02,0.05,0.08,0.1];
    end if;
    
    forEACH bufferInstance IN ARRAY bufferArray
    LOOP
    select _pgr_shortest_road(startx,starty,endx,endy,bufferInstance) into v_shPath;
    if (v_shPath is not null) then    
        return v_shPath;
    end if;
    end loop;
    
    end;  
    $BODY$
  LANGUAGE plpgsql VOLATILE STRICT
  COST 100;
ALTER FUNCTION public.pgr_shortest_road(double precision, double precision, double precision, double precision )
  OWNER TO postgres;

DROP FUNCTION public._pgr_shortest_road(double precision, double precision, double precision, double precision, double precision);

The above function can actually be basically satisfied for operations in most cases.

3. Efficiency optimization

In fact, in terms of data query, we use linear buffering between the starting point and the end point to improve efficiency, as follows:

SELECT gid as id, source, target, cost,rev_cost as reverse_cost FROM china_road
      where geom && st_buffer(st_linefromtext(''linestring('||startx||' ' || starty ||','|| endx ||' ' || endy ||')'',4326),'||bufferDistance||')

Of course, this is still good in most cases, and in some cases, it does not play a good role, because if the road is offset between the starting point and the end point (for example, when there are many mountains on the straight line, the road will be more circumferential), at this time, the buffer distance may be increased, and increasing the buffer distance will lead to the increase in the number of queries in some areas, which will affect efficiency. Therefore, we can actually consider using the mapid parameter. Where does this parameter come from? Generally, the road network data we get will have this field. We only need to generate a region table, and this region table has two fields, one is mapid and the other is the polygon range of this mapid. In this way, the query conditions above can be replaced by the following:

SELECT gid as id, source, target, cost,rev_cost as reverse_cost FROM china_road
      where mapid in (select mapid from maps where geom && st_buffer(st_linefromtext(''linestring('||startx||' ' || starty ||','|| endx ||' ' || endy ||')''),'||bufferDistance||'))

This will greatly improve efficiency.

4. Data bug handling

In fact, sometimes the road network data we get is not very accurate, or the input is flawed. What I encountered is the generated topo data. The target of the original road should coincide with the source of its adjacent road, but in fact it is different, which leads to problems in the final planning. Therefore, I simply wrote a function to deal with this problem.

CREATE OR REPLACE FUNCTION public.modity_road_data()
RETURNS void AS
$BODY$   
declare  
n integer;
begin  
    for n IN (select distinct(source) from china_road )  loop
        update china_road
        set geom = st_multi(st_addpoint(ST_geometryN(geom,1),
        (select st_pointn(ST_geometryN(geom,1),1) from china_road where source = n
        limit 1),
        st_numpoints(ST_geometryN(geom,1))))
        where target = n;
    end loop;
    end;  
    $BODY$
  LANGUAGE plpgsql VOLATILE STRICT
  COST 100;
ALTER FUNCTION public.modity_road_data()
OWNER TO postgres;

5. Follow-up planning

The above function has been verified in millions of data, and will verify the road network data of tens of millions of levels in the future. Of course, this level must be adjusted in strategy. For example, in the national road network recently tested, first plan the high-speed entrance from the starting point to the starting point, first plan the high-speed exit from the planning end point to the end point, and then plan the path from the high-speed entrance to the high-speed exit on the high-speed road network. This will find that the efficiency has been greatly improved. Of course, there are still many logic and services in it. After everything is verified, another edition of the article on the path planning of tens of millions of levels will be published.

This is the end of this article about path planning based on pgrouting. For more related path planning content for pgrouting, please search for my previous articles or continue browsing the related articles below. I hope everyone will support me in the future!