wiki:import layers postgis

Context Navigation

Version 26 (modified by cedric, 7 years ago) (diff)
--

back to start page ..

Import layer postgis

first fixing a problem with the database

At first the following failed and I had to fix a problem in the database

~~spatial_ref_sys is empty in the eda.0 database => collecting it from another database~~

pg_dump -U postgres -f "spatial_ref_sys.sql" --data-only --table public.spatial_ref_sys --inserts --verbose eda2.2
psql -U postgres -f spatial_ref_sys.sql eda2.0

this does not work, it's a known bug https://trac.osgeo.org/postgis/ticket/1815

I'm using the sql file from the postgis installation. For this we will try on a subset in the oria

psql -U postgres -f "C:/Program Files/PostgreSQL/9.6/share/contrib/postgis-2.4/spatial_ref_sys.sql" eda2.0

First trial with two layers

Our objective is to test the best way to route the river network

-- first create a new schema for spain
psql -U postgres -c "create schema spain" eda2.0

-- shp2pgsql with create index option
-- Subcuencas de ríos completos clasificadas según Pfafstetter modificado
shp2pgsql -s 4258 -I A_cuencas_rios_Atl_Norte.shp spain.a_cuencas_rios_atl_norte > a_cuencas_rios_atl_norte.sql
psql -U postgres -f a_cuencas_rios_atl_norte.sql eda2.0

-- Tramos de ríos de España clasificados según Pfafstetter modificado
shp2pgsql -s 4258 -W LATIN1 -I A_Rios_v2.shp spain.a_rios_v2 > a_rios_v2.sql
psql -U postgres -f a_rios_v2.sql eda2.0

in fact the system of projection indicated on the inspire layers is wrong it is 25830 http://www.mapama.gob.es/es/cartografia-y-sig/ide/directorio_datos_servicios/caracteristicas_wms.aspx

select UpdateGeometrySRID('spain','a_cuencas_rios_atl_norte', 'geom', 25830 );
select UpdateGeometrySRID('spain','a_rios_v2', 'geom', 25830 );

Now we need to select only the Oria and neighbouring rivers for our trial

set search_path to public, spain ;-- using this I will not have to say that the tables are in the spain scheme
select  * from a_cuencas_rios_atl_norte limit 10
select  * from a_cuencas_rios_atl_norte where nom_rio_1 like '%ORIA%';

802F802

"A"

172988

"1003714"

"ORIA IBAIA"

2707

"A"

175359

"10037294"

...

10037 corresponds to the Biscay coast of Spain, 1003714 correspond à Oria

select  * from a_cuencas_rios_atl_norte where pfafrio like '10037%';
select  * from a_cuencas_rios_atl_norte where pfafrio like '1003714%';

Creating a table corresponding to ORIA

Table of basins

drop table if exists spain.oria_cuenca;
create table spain.oria_cuenca as select * from a_cuencas_rios_atl_norte where pfafrio like '1003714%';--422

CREATE INDEX oria_cuenca_geom_idx
  ON oria_cuenca
  USING gist
  (geom);

Table of rivers

drop table if exists spain.oria_a_rios_v2;
create table spain.oria_a_rios_v2 as select * from a_rios_v2 where pfafrio like '1003714%';--422

CREATE INDEX oria_a_rios_v2_geom_idx
  ON oria_a_rios_v2
  USING gist
  (geom);

Trials for the pfafsteter system

To be able to explore the topology we could use pfafstetter

/*-----------------------------------------
Function to get all segments downstream from one one in the same subbasin
Problem when the subbasin changes it is difficult to go downstream
-----------------------------------------*/
CREATE OR REPLACE FUNCTION downstream_segments(_id text) 
RETURNS TABLE(pfafcuen character varying(254), lngtramo_m numeric) AS
$$
BEGIN
        RAISE NOTICE 'downstream segment for %', _id;                                 -- prints the selected segment
        RETURN QUERY                                                                  -- RETURN QUERY allows to return a table without loop
             SELECT r.pfafcuen, r.lngtramo_m from  oria_a_rios_v2 r 
                WHERE
                 substring(_id::text,1,char_length(_id::text)-1)= 
                   substring(r.pfafcuen::text,1,char_length(r.pfafcuen::text)-1)     -- the selection is done on a subset eg. 10011 10013  10015  10017  will be selected using string 1001
                AND
                 _id>r.pfafcuen                                                      -- we select only segments downstream from the current one
                AND
                 (r.pfafcuen::numeric % 2) = 1                                       -- even number     
               order by r.pfafcuen::numeric DESC;       
END
$$
  LANGUAGE plpgsql VOLATILE;
select * from downstream_segments('10037146457');

this query returns the three odd nodes dowstream from '10037146457'

pfafcuen	lngtramo_m
10037146451	747
10037146453	304
10037146455	437

Problem from that step we need to go to another number

10037146451 => 1003714643

-- this query returns nothing
with d1 as (SELECT    pfafcuen, 
                  lngtramo_m, 
                  (substring(pfafcuen::text,1,char_length(pfafcuen::text)-1)::integer-2)::character varying(254) 
                  as pfafcuen_downstream 
  from downstream_segments('10037146457'))
select * from downstream_segments((select min(pfafcuen_downstream) FROM d1 limit 1));
--because this query returns nothing
select downstream_segments('1003714643');

Before trying too much the following problem arises
where 1037146(47) is immediately upstream from 1037146(457)

pg_routing

https://docs.pgrouting.org/2.5/en/pgRouting-concepts.html#description-of-the-edges-sql-query-for-dijkstra-like-functions

-------------------------------------
-- pgrouting installation and configuration on table
CREATE EXTENSION pgrouting;
ALTER TABLE spain.oria_a_rios_v2 DROP COLUMN source;
ALTER TABLE spain.oria_a_rios_v2 DROP COLUMN target;
ALTER TABLE spain.oria_a_rios_v2 ADD COLUMN source integer;
ALTER TABLE spain.oria_a_rios_v2 ADD COLUMN target integer;
CREATE INDEX oria_a_rios_v2_source_idx ON spain.oria_a_rios_v2 (source);
CREATE INDEX oria_a_rios_v2_target_idx ON spain.oria_a_rios_v2 (target);

---------------------------------------------------------------------
-- createTopology also creates spain.oria_a_rios_v2_vertices_pgr
------------------------------------------------------------------
set search_path to spain,public;
SELECT pgr_createTopology('oria_a_rios_v2', 0.0001, 'geom', 'gid');--1.4s
---------------------------------------------------------------------
-- routing this is done using nodes, FALSE is necessery otherwise fails (using a directed path=
------------------------------------------------------------------
SELECT seq, node, edge, agg_cost, gid, pfafcuen
FROM pgr_dijkstra(
   'SELECT gid as id, source, target, st_length(geom) as cost FROM oria_a_rios_v2',
   5,
   540,
  FALSE
)  pt
JOIN oria_a_rios_v2 r ON pt.edge = r.gid order by seq;

"seq"	"node"	"edge"	"agg_cost"	"gid"	"pfafcuen"
1	5	774	0	774	"1003714113"
2	28	783	340.239109235141	783	"1003714115"
3	33	916	434.71909815275	916	"1003714117"
4	61	920	1381.42299559645	920	"100371413"
5	62	956	1931.41349601216	956	"10037141511"
6	54	992	2751.07786782504	992	"10037141513"
7	53	1023	3006.69704002442	1023	"10037141515"
8	16	1028	3923.51604452287	1028	"10037141517"
9	15	962	4550.98432547035	962	"10037141531"
10	47	879	5241.61915518433	879	"10037141533"
11	39	880	5777.24646955956	880	"10037141535"

Download in other formats:

Plain Text