pgrouting路径分析（站点与GRID中心点最短距离）

drop function if exists pro_calculateGridPathAll();
create function pro_calculateGridPathAll()
returns text as
$$
declare
v_i integer;
begin
for v_i in 8..29 loop
select pro_calculateGridPath(v_i);
end loop;
end;
$$ language plpgsql;




drop function if exists pro_calculateGridPath(argStationId integer);
create function pro_calculateGridPath(argStationId integer)
returns text as
$$
    declare
v_cursStation cursor(argStationId integer)
for select gid as stationId, startpoint, srcpulloutable, dstpulloutable, source, target, isnode, nodeid, lineid, srcsegline, dstsegline from station where gid=argStationId order by gid asc;
v_cursGrid cursor 
for select fid as gridFid, geom as endpoint, srcaccessable, dstaccessable, source, target, isnode, nodeid, lineid, srcsegline, dstsegline from central_point order by fid asc;


v_res geometry;   --dijkstra最短路径
v_resLen float;
v_shPath geometry;--最终结果
v_shPathLen float; --最短路径长度

v_startNodeIndex integer;
v_endNodeIndex integer;


src_stationId integer; 
startpoint geometry;
src_srcpulloutable boolean;
src_dstpulloutable boolean;
src_source integer;
src_target integer;
src_isnode boolean;
src_nodeid integer;
src_inefraction float;
src_lineid integer;
src_srcsegline geometry;
src_dstsegline geometry;


dst_gridFid integer; 
endpoint geometry;
dst_srcaccessable boolean;
dst_dstaccessable boolean;
dst_source integer;
dst_target integer;
dst_isnode boolean;
dst_nodeid integer;
dst_inefraction float;
dst_lineid integer;
dst_srcsegline geometry;
dst_dstsegline geometry;


v_sbuffer geometry;
v_ebuffer geometry; 
v_perStart float;
v_perEnd float;


v_hasPath boolean;


v_calculateIndex integer;
v_spheroid  spheroid;


v_line geometry;
v_count integer;
v_cursPgrPath cursor (isource integer, itarget integer) 
for select a.seq,a.id1,a.id2,a.cost,ST_LineMerge(b.geom) from 
pgr_dijkstra('select gid as id, source::integer, target::integer, 
cost::double precision AS cost, reverse_cost::double precision AS reverse_cost 
from changsharoad_connected', 
isource, itarget, true, true)a left join changsharoad_connected b 
ON (a.id2=b.gid) order by a.seq asc;

begin
v_calculateIndex := 0;
v_spheroid := 'SPHEROID["WGS84", 6378137, 298.257223563]' ;  --WGS84椭球体参数定义
open v_cursGrid;


<<loopGrid>>
loop
fetch v_cursGrid into dst_gridFid, endpoint, dst_srcaccessable, dst_dstaccessable, dst_source, dst_target, dst_isnode, dst_nodeid, dst_lineid, dst_srcsegline, dst_dstsegline;
if not found then exit;end if;


declare
v_endPointSeg geometry[];  --终点所在线两侧的折线
v_endPointSegLen float[];
v_endPointIn boolean[];
v_endNodeId integer[];
begin
endpoint := st_setsrid(endpoint, 4326);

v_endPointIn := array_append(v_endPointIn, dst_srcaccessable);
v_endPointIn := array_append(v_endPointIn, dst_dstaccessable);


v_endNodeId := array_append(v_endNodeId, dst_source);
v_endNodeId := array_append(v_endNodeId, dst_target);


dst_srcsegline := ST_SetSRID(dst_srcsegline, 4326);
dst_dstsegline := ST_SetSRID(dst_dstsegline, 4326);

v_endPointSeg := array_append(v_endPointSeg, dst_srcsegline);
v_endPointSeg := array_append(v_endPointSeg, dst_dstsegline);

open v_cursStation(argStationId);
<<loopStation>>
loop
fetch v_cursStation into src_stationId, startpoint, src_srcpulloutable, src_dstpulloutable, src_source, src_target,src_isnode, src_nodeid, src_lineid, src_srcsegline, src_dstsegline;
if not found then exit;end if;


startpoint := st_setsrid(startpoint, 4326);
v_shPathLen := 1000000;
v_hasPath := false;

declare
v_startPointSeg geometry[];--起点所在线两侧的折线
v_startPointSegLen float[];
v_startPointOut boolean[];
v_startNodeId integer[]; 
begin


v_startPointOut := array_append(v_startPointOut, src_srcpulloutable);
v_startPointOut := array_append(v_startPointOut, src_dstpulloutable);


v_startNodeId := array_append(v_startNodeId, src_source);
v_startNodeId := array_append(v_startNodeId, src_target);


v_startPointSeg := array_append(v_startPointSeg, st_setsrid(src_srcsegline, 4326));
v_startPointSeg := array_append(v_startPointSeg, st_setsrid(src_dstsegline, 4326));


v_hasPath := false; --初始化为没有路径
v_line := null;
v_shPath := null;

raise notice  '【准备开始循环 -------------  gridFid: %, stationId: %  -------------】', dst_gridFid, src_stationId;

v_startNodeIndex := 1;


<<label1>>
while (v_startNodeIndex <= array_length(v_startPointOut, 1)) loop
v_endNodeIndex := 1;
--raise notice  'v_startPointOut[%]: %', v_startNodeIndex, v_startPointOut[v_startNodeIndex];
if (v_startPointOut[v_startNodeIndex] AND v_startNodeId[v_startNodeIndex] != v_endNodeId[v_endNodeIndex]) then
<<label2>>
while (v_endNodeIndex <= array_length(v_endPointIn, 1)) loop
declare 
v_path pgr_costResult;
v_shorts pgr_costResult[];
v_geoms geometry[];


begin
if (dst_isnode and v_endNodeId[v_endNodeIndex] != dst_nodeid) then
v_endNodeIndex := v_endNodeIndex + 1;
continue;
else
--raise notice  'v_endPointIn[%]: %', v_endNodeIndex, v_endPointIn[v_endNodeIndex];
--raise notice  '   组 -------------  %-->%  -------------', v_startNodeId[v_startNodeIndex], v_endNodeId[v_endNodeIndex];
if (v_endPointIn[v_endNodeIndex]) then
open v_cursPgrPath(v_startNodeId[v_startNodeIndex], v_endNodeId[v_endNodeIndex]);
v_line := null;
<<label3>>
loop
fetch v_cursPgrPath into v_path.seq,v_path.id1,v_path.id2,v_path.cost,v_line;
--exit when not found;-- 假如没有检索到(主表)数据，结束循环处理
----raise notice  'v_path.seq,v_path.id1,v_path.id2,v_path.cost: %,%,%,%', v_path.seq,v_path.id1,v_path.id2,v_path.cost;
if not found then exit;end if;


v_shorts := array_append(v_shorts,v_path);
v_geoms := array_append(v_geoms, ST_SetSRID(v_line, 4326));
end loop label3;
close v_cursPgrPath;
/*
SELECT ST_LineMerge(st_union(b.geom)) into v_line FROM (select * from pgr_dijkstra('     
SELECT gid AS id,
source::integer,
target::integer,
cost::double precision AS cost,
reverse_cost::double precision AS reverse_cost
FROM changsharoad_connected',
v_startNodeId[v_startNodeIndex], v_endNodeId[v_endNodeIndex], true, true) order by seq asc)
a LEFT JOIN changsharoad_connected b ON (a.id2 = b.gid);
*/
if( array_length(v_shorts, 1) > 1) then
v_hasPath := true;--找到一条最短路径
--raise notice  'v_hasPath := true;';
if (dst_lineid = src_lineid) then
v_sbuffer := ST_Buffer(startpoint::geography, 3);
if( false = ST_Intersects(v_geoms[1], v_sbuffer) ) then  --如果起点和起点线段没有相连，增加起始路段
----raise notice  '起点线段处理开始2';
--v_lstart := pro_truncation(v_lstart,v_geoms[1], v_statpoint,v_spheroid);
--v_geoms[1] :=  ST_Union(v_geoms[1],v_lstart);
v_geoms[1] := st_linemerge(ST_Union(array[v_startPointSeg[v_startNodeIndex], v_geoms[1], 
v_endPointSeg[v_endNodeIndex]]));
else    --如果起点和起点线段有相连，去除起始路段
select  ST_LineLocatePoint(v_geoms[1], startpoint) into v_perStart;  
select  ST_LineLocatePoint(v_geoms[1], endpoint) into v_perEnd;


if (v_perStart < v_perEnd) then
v_geoms[1] := ST_LineSubstring(v_geoms[1], v_perStart, v_perEnd);
else 
v_geoms[1] := ST_LineSubstring(v_geoms[1], v_perEnd, v_perStart);
v_geoms[1] := ST_REVERSE(v_geoms[1]);
end if;
end if;
else
--起点线段处理
--raise notice  '起点线段处理开始1';
v_sbuffer := ST_Buffer(startpoint::geography, 3);
if( false = ST_Intersects(v_geoms[1], v_sbuffer) ) then  --如果起点和起点线段没有相连，增加起始路段
--raise notice  '起点线段处理开始2, %', st_astext(v_startPointSeg[v_startNodeIndex]);
v_geoms[1] := st_linemerge(ST_Union(array[v_startPointSeg[v_startNodeIndex],  st_setsrid(v_geoms[1], 4326)]));
else    --如果起点和起点线段有相连，去除起始路段
if (v_startNodeIndex = 1) then
----raise notice  '起点线段处理开始3';
v_geoms[1] := v_startPointSeg[2];
else
----raise notice  '起点线段处理开始4';
v_geoms[1] := v_startPointSeg[1];
end if;
end if;
--raise notice  '起点线段处理完成5';


--raise notice 'dst_isnode: %', dst_isnode;
--终点线段处理
if (dst_isnode <> true or dst_isnode is null) then
--raise notice  '终点线段处理开始1';
v_count := array_length(v_geoms, 1);
v_ebuffer := ST_Buffer(endpoint::geography, 3);
if (false = ST_Intersects(v_geoms[v_count - 1], v_ebuffer)) then  --如果终点和终点线段没有相连
--raise notice  '终点线段处理开始2';
v_geoms[v_count - 1] := st_linemerge(ST_Union(array[v_geoms[v_count - 1], v_endPointSeg[v_endNodeIndex]]));
else
if (v_endNodeIndex = 1) then
----raise notice  '终点线段处理开始3';
v_geoms[v_count - 1] := v_endPointSeg[2];
else
----raise notice  '终点线段处理开始4';
v_geoms[v_count - 1] := v_endPointSeg[1];
end if;
end if;
--raise notice  '终点线段处理完成5';
end if;
end if;

v_res := ST_SetSRID(st_linemerge(ST_Union(v_geoms)), 4326);
----raise notice  '合并所有路段，用来计算总长';
--长度判断
v_resLen := ST_LengthSpheroid(v_res, v_spheroid);
----raise notice  '合并总长：%', v_resLen;
if (v_resLen < v_shPathLen) then
v_shPath := v_res;
/*delete from dijkstra_res;
insert into dijkstra_res(gid, seq, geom) values(1, 1, ST_GeomFromText(st_astext(v_res), 4326));*/
v_shPathLen := v_resLen;
end if;
end if;
end if;


v_endNodeIndex := v_endNodeIndex + 1;


if (dst_isnode and v_endNodeId[v_endNodeIndex] = dst_nodeid) then
exit label2;
end if;
end if;
end;
end loop label2;
end if;
v_startNodeIndex := v_startNodeIndex + 1;
end loop label1;


end;

if (v_hasPath) then
v_calculateIndex := dst_gridFid * 29 + src_stationId - 1;
--raise notice '%, 找到最短路径', src_stationId;
INSERT INTO kxjc_gridstationresponse(id,stationid,gridfid,timecost,distance,schemeid, geom)
VALUES(v_calculateIndex,src_stationId,dst_gridFid,v_shPathLen / 666.7,v_shPathLen,1, v_shPath)
ON conflict(id) 
DO update set 
id = v_calculateIndex,
stationid = src_stationId,
gridfid = dst_gridFid,
timecost = v_shPathLen / 666.7,
distance = v_shPathLen,
schemeid = 1,
geom = v_shPath;
end if;
end loop loopStation;
close v_cursStation;


end;
end loop loopGrid;
close v_cursGrid;


return 'ok';
end;
$$ language plpgsql;