使用PostGIS完成兩點(diǎn)間的河流軌跡及流經(jīng)長度的計(jì)算(推薦)
基礎(chǔ)準(zhǔn)備工作
1.PostGIS 的安裝
在安裝PostGIS前首先必須安裝PostgreSQL,然后再安裝好的Stack Builder中選擇安裝PostGIS組件。具體安裝步驟可參照PostGIS的安裝與初步使用
2.加載Post GIS擴(kuò)展
選中指定數(shù)據(jù)庫,執(zhí)行加載擴(kuò)展語句
–添加支持 CREATE EXTENSION postgis; --添加postgis擴(kuò)展 CREATE EXTENSION pgrouting; --添加pgrouting擴(kuò)展 CREATE EXTENSION postgis_topology; CREATE EXTENSION fuzzystrmatch; CREATE EXTENSION postgis_tiger_geocoder;
在做兩點(diǎn)間河流軌跡及流經(jīng)長度計(jì)算過程中,需要加載postgis和pgrouting兩個擴(kuò)展
可以通過查看加載擴(kuò)展的版本驗(yàn)證擴(kuò)展加載是否成功
–查看postgresql版本 show server_version; –查看postgis版本 SELECT PostGIS_full_version(); –查看pgrouting版本 select pgr_version();
3.河流矢量圖層轉(zhuǎn)成單線格式
河流包括各種匯入和匯出,為了實(shí)現(xiàn)流經(jīng)流域的計(jì)算,河流水系矢量數(shù)據(jù)需要一個河流一個ID的方式,可以在河流交匯點(diǎn)處將河流進(jìn)行打段處理。
4.河流矢量數(shù)據(jù)導(dǎo)入PostgreSQL數(shù)據(jù)庫
打開位于“開始>所有程序>PostGIS 2.3 bundle for PostgreSQL”之中的PostGIS Shapefile Import/Export Manager。
首先單擊"View connection details"按鈕,打開"PostGIS connection"對話框,輸入用戶名"postgres"及其對應(yīng)的密碼,設(shè)置連接的數(shù)據(jù)庫,如下圖所示:
連接數(shù)據(jù)庫之后,單擊"Add file"按鈕,加入***.shp文件,并將其SRID設(shè)置為"4326",如下圖所示。這一步絕對不能省略,否則不能正確導(dǎo)入數(shù)據(jù)。
5.河流數(shù)據(jù)拓?fù)涮幚?/h3>
在數(shù)據(jù)分析過程中,使用到了pgrouting擴(kuò)展中的 pgr_dijkstra 算法
Dijkstra算法(迪杰斯特拉算法),由荷蘭計(jì)算機(jī)科學(xué)家Edsger Dijkstra于1956年提出。它是一種圖搜索算法,它解決了非負(fù)代價(jià)邊路徑圖的最短路徑問題,即從起始頂點(diǎn)(start_vid)到結(jié)束頂點(diǎn)(end_vid)的最短路徑。此算法可以與有向圖或無向圖一起使用。
函數(shù)的簽名摘要:
在實(shí)際使用中,需要先明確所有的頂點(diǎn),并為所有頂點(diǎn)分配唯一的編號,函數(shù)的 start_vid 和 end_vid 都是整型數(shù)值,函數(shù)使用edges_sql參數(shù)(sql腳本)篩選出和頂點(diǎn)相鄰的所有邊信息(即河流信息)。
所以,在使用pgr_dijkstra方法前,需要
- 對找到河流的所有頂點(diǎn)信息,并做唯一整型值編號
- 在數(shù)據(jù)庫中為每條河流設(shè)置好起始頂點(diǎn)和結(jié)束頂點(diǎn)
--篩選出所有頂點(diǎn)信息,st_dump函數(shù)主要是將MultiLineString類型 調(diào)整成 LineString類型 select st_astext(st_startpoint((ST_Dump(geom)).geom)) from singleriver union select st_astext(st_endpoint((ST_Dump(geom)).geom)) from singleriver
將查詢結(jié)果在Excel中進(jìn)行整型值編號,再導(dǎo)入到postgresql中的新建表distinctpoint 中,然后關(guān)聯(lián)河流數(shù)據(jù)表,更新河流的開始頂點(diǎn)(source)和結(jié)束頂點(diǎn)編號(target)
--更新起始頂點(diǎn)編號 update singleriver q set source=tt.sourcepoint from singleriver s, (select gid,p.id as sourcepoint from (select gid,st_astext(st_startpoint((ST_Dump(geom)).geom)) as startpoint, st_astext(st_endpoint((ST_Dump(geom)).geom)) as endpoint from singleriver )s left join distinctpoint p on s.startpoint=p.point) tt where q.gid=tt.gid
--插入結(jié)束頂點(diǎn)編號 update singleriver q set target=tt.endpoint from singleriver s, (select gid,p.id as endpoint from (select gid,st_astext(st_startpoint((ST_Dump(geom)).geom)) as startpoint, st_astext(st_endpoint((ST_Dump(geom)).geom)) as endpoint from singleriver )s left join distinctpoint p on s.endpoint=p.point) tt where q.gid=tt.gid
至此,河流拓?fù)鋽?shù)據(jù)處理完成
PG分析處理函數(shù)
1.函數(shù)編寫
CREATE OR REPLACE FUNCTION "public"."pgr_shortest_river"(IN "startx" float8, IN "starty" float8, IN "endx" float8, IN "endy" float8, OUT "river_name" varchar, OUT "v_shpath" varchar, OUT "cost" float8) RETURNS SETOF "pg_catalog"."record" AS $BODY$ declare v_startLine geometry;--離起點(diǎn)最近的線 v_endLine geometry;--離終點(diǎn)最近的線 v_startTarget integer;--距離起點(diǎn)最近線的終點(diǎn) v_endSource integer;--距離終點(diǎn)最近線的起點(diǎn) v_statpoint geometry;--在v_startLine上距離起點(diǎn)最近的點(diǎn) v_endpoint geometry;--在v_endLine上距離終點(diǎn)最近的點(diǎn) v_res geometry;--最短路徑分析結(jié)果 v_perStart float;--v_statpoint在v_res上的百分比 v_perEnd float;--v_endpoint在v_res上的百分比 v_rec record; first_name varchar; end_name varchar; first_cost double precision; end_cost double precision; begin --查詢離起點(diǎn)最近的線 execute 'select (st_dump(geom)).geom as geom,target as target,name from singleriver 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; raise notice '起點(diǎn)線段%',v_startLine; raise notice '起點(diǎn)位置%',v_startTarget; raise notice '河流名稱%',first_name; --查詢離終點(diǎn)最近的線 execute 'select (st_dump(geom)).geom as geom,"source" as source,name from singleriver 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; --如果沒找到最近的線,就返回null 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; --計(jì)算距離起點(diǎn)最近線上的點(diǎn)在該線中的位置 select st_linelocatepoint(st_linemerge(v_startLine), v_statpoint) into v_perStart; select st_linelocatepoint(st_linemerge(v_endLine), v_endpoint) into v_perEnd; select st_distancesphere(v_statpoint,ST_PointN(ST_GeometryN(v_startLine,1), ST_NumPoints(ST_GeometryN(v_startLine,1)))) into first_cost; select st_distancesphere(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_distancesphere(v_statpoint, v_endpoint) into first_cost; select st_linelocatepoint(st_linemerge(v_startLine), v_endpoint) into v_perEnd; for v_rec in select st_linesubstring(st_linemerge(v_startLine), v_perStart,v_perEnd) as point,COALESCE(end_name,'無名河流') as name,end_cost as cost loop v_shPath:= ST_AsGeoJSON(v_rec.point); cost:= v_rec.cost; river_name:= v_rec.name; return next; end loop; return; end if; --最短路徑 for v_rec in (select st_linesubstring(st_linemerge(v_startLine),v_perStart,1) as point,COALESCE(first_name,'無名河流') as name,first_cost as cost union all SELECT st_linemerge(b.geom) as point,COALESCE(b.name,'無名河流') as name,st_length(geom, false) as cost FROM pgr_dijkstra( 'SELECT gid as id, source, target, st_length(geom, false) as cost FROM singleriver where st_intersects(geom,st_buffer(st_linefromtext(''linestring('||startx||' ' || starty ||','|| endx ||' ' || endy ||')''),0.05))', v_startTarget, v_endSource , false ) a, singleriver b WHERE a.edge = b.gid union all select st_linesubstring(st_linemerge(v_endLine),0,v_perEnd) as point,COALESCE(end_name,'無名河流') as name,end_cost as cost) loop v_shPath:= ST_AsGeoJSON(v_rec.point); cost:= v_rec.cost; river_name:= v_rec.name; return next; end loop; end; $BODY$ LANGUAGE plpgsql VOLATILE STRICT COST 100 ROWS 1000
2.參數(shù)說明
輸入?yún)?shù):開始點(diǎn)和結(jié)束點(diǎn)的經(jīng)緯度坐標(biāo)
輸出結(jié)果:river_name:河流名稱;v_shppath:流經(jīng)的河流路徑; cost:河流流經(jīng)長度
3.內(nèi)部調(diào)用函數(shù)說明
函數(shù)調(diào)用過程,根據(jù)postgis不同版本,函數(shù)名稱可能會有偏差,有版本展示形式為st_linesubstring ,有版本展示形式為st_line_substring
4.輸出結(jié)果驗(yàn)證
為了驗(yàn)證河流輸出結(jié)果是否正確,流經(jīng)河流路徑是否連通,可以通過在線geojson地圖(geojson.io)呈現(xiàn)出來驗(yàn)證。
在右側(cè)json-features-geometry 中填充函數(shù)輸出的v_shppath參數(shù)內(nèi)容(按照行單獨(dú)輸入,可以輸入多個,注意需要增加json屬性)
到此這篇關(guān)于使用PostGIS完成兩點(diǎn)間的河流軌跡及流經(jīng)長度的計(jì)算的文章就介紹到這了,更多相關(guān)PostGIS兩點(diǎn)間的河流軌跡計(jì)算內(nèi)容請搜索本站以前的文章或繼續(xù)瀏覽下面的相關(guān)文章希望大家以后多多支持本站!
版權(quán)聲明:本站文章來源標(biāo)注為YINGSOO的內(nèi)容版權(quán)均為本站所有,歡迎引用、轉(zhuǎn)載,請保持原文完整并注明來源及原文鏈接。禁止復(fù)制或仿造本網(wǎng)站,禁止在非www.sddonglingsh.com所屬的服務(wù)器上建立鏡像,否則將依法追究法律責(zé)任。本站部分內(nèi)容來源于網(wǎng)友推薦、互聯(lián)網(wǎng)收集整理而來,僅供學(xué)習(xí)參考,不代表本站立場,如有內(nèi)容涉嫌侵權(quán),請聯(lián)系alex-e#qq.com處理。