Context Navigation

back to first page..
back to CookBook Eda
back to RHT Go to CookBook RHT_UGA
See #80

Joining RHT with BDMAP

Joining the riversegment from RHT and stationsp2 (??)

-- creation de la table bdmap_rht, je rajoute un the_geom a la requete pour pouvoir voir la table dans Qgis  
DROP TABLE IF EXISTS rht.bdmap_rht;
CREATE TABLE rht.bdmap_rht as (
        SELECT distinct on (st_codecsp) st_codecsp, id_drain, min(distance) as distance, the_geom FROM (
               SELECT st_codecsp, id_drain ,CAST(distance(r.the_geom, s.the_geom) as  decimal(15,1)) as distance,s.the_geom 
               FROM bd_map.stationsp2 As s
               INNER JOIN  rht.rht r ON ST_DWithin(r.the_geom, s.the_geom,300)
               WHERE s.the_geom IS NOT NULL
               ORDER BY st_codecsp) AS sub 
        GROUP BY st_codecsp, distance,id_drain, the_geom
);
alter table rht.bdmap_rht add column id serial;
-- mise à jour de la table geometry_columns
INSERT INTO geometry_columns(f_table_catalog, f_table_schema, f_table_name, f_geometry_column, coord_dimension, srid, "type")
SELECT '', 'rht', 'bdmap_rht', 'the_geom', ST_CoordDim(the_geom), ST_SRID(the_geom), GeometryType(the_geom)
FROM rht.bdmap_rht LIMIT 1;

-- creation d'index, clé primaire, et constraintes qui vont bien
alter table rht.bdmap_rht add CONSTRAINT enforce_dims_the_geom CHECK (ndims(the_geom) = 2);
alter table rht.bdmap_rht add  CONSTRAINT enforce_geotype_the_geom CHECK (geometrytype(the_geom) = 'POINT'::text OR the_geom IS NULL);
alter table rht.bdmap_rht add  CONSTRAINT enforce_srid_the_geom CHECK (srid(the_geom) = 3035);
alter table rht.bdmap_rht ADD CONSTRAINT pk_id PRIMARY KEY(id);
CREATE INDEX indexbdmap_rht ON rht.bdmap_rht
  USING GIST ( the_geom GIST_GEOMETRY_OPS );

select count(*) from bd_map.stationsp2  --10193 lines
select count(*) from rht.rht  --114601 lines
select count(*) from rht.bdmap_rht  -- 9650 lines

Soit 94% des stations projetées sur le RHT

select count(*) from bd_map.bdmap_ccm_gid  --8188 lines

Soit 80% des stations projetées sur la CCM

Some station must be deleted or reproject and the id_drain must be changed:

DROP table if exists rht.bdmap_rht_id;
CREATE TABLE rht.bdmap_rht_id as (
select distinct on (st_codecsp) * from rht.bdmap_rht
);

---Tous les cours d'eau situés à plus de 200m ont été vérifiés.
delete from rht.bdmap_rht_id where st_codecsp in ('03760141','03270131','01620009','01620005','01620137','01590014','06890285','06210238','03890095','02570237','02670291','02670264','02670097','02670092','02670098','02670230','02670231','02670232','02680172','02680166','02680039','02680033','02680037','02680144','02680143');  -- 25 stations supprimées
delete from rht.bdmap_rht_id where st_codecsp in ('05400178','05400044','0547C020','06040033','06040206','06070217','06840018','03600120','03910037','03910021','04450005','04181015','04181000','03510116','06390254','06250304','06730233','02880208','06700020','02550268','02080064','02080065','03510111','03520110','02670261','02670202','02670203');  -- 27 stations supprimées
delete from rht.bdmap_rht_id where st_codecsp in ('04420483','04420401','04420181','04420190','04430154','03890102','03890135','04630032','04420437'); -- 9 stations supprimées
delete from rht.bdmap_rht_id where st_codecsp in ('06260192','06380273','02880051','02670109','02670055','02670176','02670177','02670143','02670178','02670141','02670082','02670039','02570237','02540023','02570203');  --14 lines
delete from rht.bdmap_rht_id where st_codecsp in ('06880149','02880048','04420458','02680085','02680039','02670004','02670053','06520051','03080085','02080027','06840052');  -- 11 stations

update rht.bdmap_rht_id SET id_drain='7866' where st_codecsp='06010369';

Verification des stations projetées sur ccm et dont gid a été modifié --> changement si besoin sur le rht de l'id_drain Stations vérifiées et bien projetées :

st_codecsp IN ('06210010','05405071','05645242','05471015','05245023','03580002','03760112','05332023','05330018','06250227','05630038','05470066','03500175','06250317,'03580002','06690101','03760112')

BDMAP with data from 2009

---Create database BDMAP2009
d:
cd D:\CelineJouanin\BDMAP
C:\"Program Files"\PostgreSQL\8.4\bin\psql -d BDMAP2009 -h localhost -U postgres -p 5432 -f bdmap20110803.backup

d:
cd D:\CelineJouanin\export_schema
C:\"Program Files"\PostgreSQL\8.4\bin\pg_dump  -U postgres -p 5432 -t station_bdmapv2_final BDMAP2009> station_bdmapv2_final.sql
C:\"Program Files"\PostgreSQL\8.4\bin\psql -d eda2.0_RHT -h localhost -U postgres -p 5432 -f station_bdmapv2_final.sql 
C:\"Program Files"\PostgreSQL\8.4\bin\pg_dump  -U postgres -p 5432 -t station BDMAP2009> station.sql
C:\"Program Files"\PostgreSQL\8.4\bin\psql -d eda2.0_RHT -h localhost -U postgres -p 5432 -f station.sql 

---Create schema bdmap2009 into eda2.0_RHT
create schema bdmap2009;
Alter table station_bdmapv2_final set schema bdmap2009;
Alter table station set schema bdmap2009;

select * from bdmap2009.station_bdmapv2_final where codehydro is null   ---440 lines
select * from bdmap2009.station_bdmapv2_final where th_code is null   ---1 lines
select * from bdmap2009.station_bdmapv2_final where cgenelin is null  ---440 lines
select * from bdmap2009.station_bdmapv2_final where eh_codegen is null  ---0 lines
select count(*) from bdmap2009.station_bdmapv2_final ---8759 lines
select * from bdmap2009.station_bdmapv2_final stat inner join bd_carthage2011.troncon_hydrographique as bdc on bdc.code_hydro=stat.eh_codegen;
select * from bdmap2009.station_bdmapv2_final stat left join bd_carthage2011.troncon_hydrographique as bdc on bdc.code_hydro=stat.eh_codegen; 
select count(*) from bdmap2009.station_bdmapv2_final stat inner join bd_carthage2011.troncon_hydrographique as bdc on bdc.code_hydro=stat.eh_codegen; --58875 ca ne marche pas du tout

travail de projection des données

/* 
Travail brouillon intégration bd_map
cedric Briand 13-octobre 2011
*/

/* 
questions
Laurent

Est que les points ont été reprojettés. Est ce que comme lannee dernière le fichier correspond à la reprojection de 2008 seulement,
si non avez vous des shapes plus récents : => NON
Est ce que les données correspondent à station st_abcisse_l93 st_ordonnee_l93 ? => NON CES COLONNES SONT VIDES
Si c'est le cas à quoi sert le fichier bdmap2009.station_bdmapv2_final ? => OK VOUS N'AVEZ PAS RETRAVAILLE DESSUS MAIS IL FAUT RECUPERER LES FICHIERS DU RCS POUR AVOIR CERTAINES COORDONNEES MANQUANTES ... AU FINAL IL EN MANQUE PAS BEAUCOUP
Lors de cette reprojection avez vous une table avec les identifiants id_bdcarthage du cours d'eau.... c'est sur cette hypothèse que nous travaillions...
=> CELINE NE PEUT PAS AVANCER RAPIDEMENT AVEC LES DONNEES FOURNIES  !

Est ce que station_bdmapv2_final correspond à stationsp2 ? => NON IL FAUT LA REFAIRE
Pourquoi est ce que dans station_bdmapv2_final ? -- 8759 alors que céline select count(*) from bd_map.stationsp2  --10193 lines
=> PARCE QUE STATIONSP2 EST UNE JOINTURE RIGHT JOIN ENTRE DEUX TABLES
*/

-- quelques requêtes pour voir les tables
/*
select * from bdmap2009.station_bdmapv2_final limit 10;
select * from bdmap2009.station limit 100;
select count(*) from bdmap2009.station_bdmapv2_final -- 8759
select count(*) from bdmap2009.station --11379
select * from bd_carthage2011.troncon_hydrographique limit 10;
select distinct on (code_hydro) count(*) from bd_carthage2011.troncon_hydrographique group by code_hydro
-- pas de jointure avec la méthode proposée par Céline
select count(*) from bdmap2009.station_bdmapv2_final stat 
inner join bd_carthage2011.troncon_hydrographique as bdc on bdc.code_hydro=stat.th_code --12325
*/

-- reprojection de bdmap2009.station_bdmapv2_final

alter table bdmap2009.station_bdmapv2_final drop constraint enforce_srid_the_geom;
update bdmap2009.station_bdmapv2_final set the_geom=st_transform(the_geom,3035);
alter table bdmap2009.station_bdmapv2_final add constraint enforce_srid_the_geom CHECK (st_srid(the_geom) = 3035);

--select st_transform(PointFromText('POINT(' || st_abcisse || ' ' || st_ordonnee || ')',27572),3035) from bdmap2009.station;
--select * from bdmap2009.station where st_abcisse_l93>0 -- 0 lines on en peut pas utiliser

SELECT AddGeometryColumn('bdmap2009', 'station','the_geom', 3035,'POINT',2); 
UPDATE bdmap2009.station SET the_geom=st_transform(PointFromText('POINT(' || st_abcisse || ' ' || st_ordonnee || ')',27572),3035);
ALTER TABLE bdmap2009.station SET WITH OIDS;
CREATE INDEX indexStations ON bdmap2009.station
  USING GIST ( the_geom GIST_GEOMETRY_OPS ); 

-- il faut aller chercher les stations qui manquent dans bdmap2009.station_bdmapv2_final

DROP TABLE IF EXISTS bdmap2009.stationsp2;
CREATE TABLE bdmap2009.stationsp2 WITH OIDS AS(
  select 
  bdmap2009.station.st_altitude,
  bdmap2009.station.st_abcisse,
  bdmap2009.station.st_codecsp,
  bdmap2009.station.st_datearret,
  bdmap2009.station.st_datecreation,
  bdmap2009.station.st_distancesource,
  bdmap2009.station.st_distancemer,
  bdmap2009.station.st_finalite,
  bdmap2009.station.st_imageign,
  bdmap2009.station.st_imagedept,
  bdmap2009.station.st_lieudit,
  bdmap2009.station.st_limites,
  bdmap2009.station.st_localisation,
  bdmap2009.station.st_longueur,
  bdmap2009.station.st_moduleia,
  bdmap2009.station.st_cd_naturecourseau,
  bdmap2009.station.st_ordonnee,
  bdmap2009.station.st_penteign,
  bdmap2009.station.st_pkaval,
  bdmap2009.station.st_raisremp,
  bdmap2009.station.st_sbv,
  bdmap2009.station.st_t_janvier,
  bdmap2009.station.st_t_juillet,
  bdmap2009.station.st_cd_typecourseau,
  bdmap2009.station.st_cd_tet,
  bdmap2009.station.st_st_id,
  bdmap2009.station.st_cm_id,
  bdmap2009.station.st_cx_id,
  bdmap2009.station.st_th_id,
  bdmap2009.station.st_eh_id,
  bdmap2009.station.st_uh_id,
  bdmap2009.station.st_dt_cre,
  bdmap2009.station.st_dt_maj,
  bdmap2009.station.st_qi_maj,
  bdmap2009.station.st_masseeau,
  CASE WHEN bdmap2009.station_bdmapv2_final.the_geom IS NULL THEN station.the_geom
       ELSE bdmap2009.station_bdmapv2_final.the_geom
       END AS the_geom
  FROM bdmap2009.station_bdmapv2_final right join bdmap2009.station on station_bdmapv2_final.st_codecsp=bdmap2009.station.st_codecsp);--11379
COMMENT ON TABLE bdmap2009.stationsp2 is 'Données BDMAP avec l année 2009 transmises par Laurent Beaulaton bdmap20110803.backup et mise à jour avec les coordonnées reprojetées par Hélène'

-- je vire les stations hors de France, par rapport à la dernière fois deux sont rentrées en France avec des coordonnées correctes.
delete from bdmap2009.stationsp2 where st_codecsp IN ('03890162','03270062','03270061','03270063','03270064','06260126','06420055','06070241')




CREATE INDEX indexStationsp2 ON bdmap2009.stationsp2
  USING GIST ( the_geom GIST_GEOMETRY_OPS );


-- creation de la table bdmap2009.bdmap_rht, je rajoute un the_geom a la requete pour pouvoir voir la table dans Qgis  
DROP TABLE IF EXISTS bdmap2009.bdmap_rht;
CREATE TABLE bdmap2009.bdmap_rht as (
        SELECT distinct on (st_codecsp) st_codecsp, id_drain, min(distance) as distance, the_geom FROM (
               SELECT st_codecsp, id_drain ,CAST(distance(r.the_geom, s.the_geom) as  decimal(15,1)) as distance,s.the_geom 
               FROM bdmap2009.stationsp2 As s
               INNER JOIN  rht.rht r ON ST_DWithin(r.the_geom, s.the_geom,300)
               WHERE s.the_geom IS NOT NULL
               ORDER BY st_codecsp) AS sub 
        GROUP BY st_codecsp, distance,id_drain, the_geom
);
alter table bdmap2009.bdmap_rht add column id serial;
-- mise à jour de la table geometry_columns
INSERT INTO geometry_columns(f_table_catalog, f_table_schema, f_table_name, f_geometry_column, coord_dimension, srid, "type")
SELECT '', 'bdmap2009', 'bdmap_rht', 'the_geom', ST_CoordDim(the_geom), ST_SRID(the_geom), GeometryType(the_geom)
FROM bdmap2009.bdmap_rht LIMIT 1;

-- creation d'index, clé primaire, et constraintes qui vont bien
alter table bdmap2009.bdmap_rht add CONSTRAINT enforce_dims_the_geom CHECK (ndims(the_geom) = 2);
alter table bdmap2009.bdmap_rht add  CONSTRAINT enforce_geotype_the_geom CHECK (geometrytype(the_geom) = 'POINT'::text OR the_geom IS NULL);
alter table bdmap2009.bdmap_rht add  CONSTRAINT enforce_srid_the_geom CHECK (srid(the_geom) = 3035);
alter table bdmap2009.bdmap_rht ADD CONSTRAINT pk_id PRIMARY KEY(id);
CREATE INDEX indexbdmap_rht ON bdmap2009.bdmap_rht
  USING GIST ( the_geom GIST_GEOMETRY_OPS );

-- je vire les deux anciennes tables de Céline qui ne servent plus puisque je viens de refaire le travail de projection

DROP TABLE IF EXISTS rht.bdmap_rht_id;
DROP TABLE IF EXISTS rht.bdmap_rht;

-- creation de la table bdmap_bd_carthage, je rajoute un the_geom a la requete pour pouvoir voir la table dans Qgis  
-- normalement on projette sur la bd_carthage et on utilise la jointure Rht bd_carthage pour récupérer les tronçons.


DROP table if exists bdmap2009.bdmap_bd_carthage;
CREATE TABLE bdmap2009.bdmap_bd_carthage as (
        SELECT distinct on (st_codecsp) st_codecsp,  id_bdcarth, min(distance) as distance, the_geom FROM (
               SELECT st_codecsp, id_bdcarth ,CAST(distance(r.the_geom, s.the_geom) as  decimal(15,1)) as distance,s.the_geom 
               FROM bdmap2009.stationsp2 As s
               INNER JOIN  bd_carthage2011.troncon_hydrographique r ON ST_DWithin(r.the_geom, s.the_geom,300)
               WHERE s.the_geom IS NOT NULL
               AND nature !='Aqueduc, conduite forcée' and num_superp !=1 -- je vire les tronçons superposés.
               ORDER BY st_codecsp) AS sub 
        GROUP BY st_codecsp, distance, id_bdcarth, the_geom
);--11181
alter table bdmap2009.bdmap_bd_carthage add column id serial;
-- mise à jour de la table geometry_columns
INSERT INTO geometry_columns(f_table_catalog, f_table_schema, f_table_name, f_geometry_column, coord_dimension, srid, "type")
SELECT '', 'bdmap2009', 'bdmap_bd_carthage', 'the_geom', ST_CoordDim(the_geom), ST_SRID(the_geom), GeometryType(the_geom)
FROM bdmap2009.bdmap_bd_carthage LIMIT 1;

-- creation d'index, clé primaire, et constraintes qui vont bien
alter table bdmap2009.bdmap_bd_carthage add CONSTRAINT enforce_dims_the_geom CHECK (ndims(the_geom) = 2);
alter table bdmap2009.bdmap_bd_carthage add  CONSTRAINT enforce_geotype_the_geom CHECK (geometrytype(the_geom) = 'POINT'::text OR the_geom IS NULL);
alter table bdmap2009.bdmap_bd_carthage add  CONSTRAINT enforce_srid_the_geom CHECK (srid(the_geom) = 3035);
CREATE INDEX indexbdmap_ccm ON bdmap2009.bdmap_bd_carthage
  USING GIST ( the_geom GIST_GEOMETRY_OPS );



  -- vérifications il manque peu de stations
select count(*) from bdmap2009.station --11379
select count(*) from bdmap2009.stationsp2 --11371 0.07 %
select * from bdmap2009.stationsp2 where the_geom is null; -- 110 stations ou il manque la geométrie
select count(*) from bdmap2009.bdmap_bd_carthage; --11181 sur les tronçons hydro 1.74 %
select count(*) from bdmap2009.bdmap_bd_carthage t join  rht.rht_bdcarthage r on r.id_bdcarth=t.id_bdcarth; -- 10998 sur les tronçons du RHT 3.34 %

Travail Laurent -- récupération du tronçon à partir des coordonnées des stations

-- script qui fonctionne uniquement chez Laurent (lien à la table "SIG"
UPDATE station_geography
        SET id_bdcarthage_troncon = id_bdcarth, 
        id_bdcarthage_cours_d_eau = c_hyd_cdo, 
        nom_cours_d_eau = toponyme1
        FROM "SIG".troncon_hydrographique, 
        (SELECT sg_id, MIN(ST_Distance(troncon_hydrographique.the_geom, station_geography.the_geom)) AS min 
        FROM station_geography, "SIG".troncon_hydrographique
         WHERE ST_Expand(station_geography.the_geom, 300) && troncon_hydrographique.the_geom GROUP BY sg_id) min_dist
        WHERE station_geography.sg_id = min_dist.sg_id 
        AND ST_Expand(station_geography.the_geom, 300) && troncon_hydrographique.the_geom 
        AND ST_Distance(troncon_hydrographique.the_geom, station_geography.the_geom) = min_dist.min
-- test de cédric
/*
e:
cd E:\IAV\eda\bdmap
psql -h localhost -U postgres --verbose  --file "station_geography.sql" eda2.0_RHT
alter table station_geography set schema bdmap2009;
*/
SET search_path TO bdmap2009, public; -- vous la connaissez celle là ?
SELECT * FROM station_geography; --11379
select  * from station_geography sg join station st on st.st_codecsp=sg.st_codecsp; --11379 c'est effectivement la même chose....

BDMAP - Projection Laurent

d:
cd D:\CelineJouanin\EDA20RHT\BDMAP\geolocalisation
C:\"Program Files"\PostgreSQL\8.4\bin\psql -d eda2.0_RHT -h localhost -U postgres -p 5432 -f station_geography.sql

alter table station_geography set schema bdmap2009;
---Some station outside of France
---Changement du système de projection  RGF93_Lambert_93 (srid 2154) --> ETRS LAEA (srid 3035)
alter table bdmap2009.station_geography drop constraint enforce_srid_the_geom;
drop trigger tr_station_geography on  bdmap2009.station_geography ;
update bdmap2009.station_geography set the_geom =ST_Transform(ST_SetSRID(the_geom,2154),3035);

alter table bdmap2009.station_geography add CONSTRAINT enforce_dims_the_geom CHECK (ndims(the_geom) = 2);
alter table bdmap2009.station_geography add  CONSTRAINT enforce_geotype_the_geom CHECK (geometrytype(the_geom) = 'POINT'::text OR the_geom IS NULL);
alter table bdmap2009.station_geography add  CONSTRAINT enforce_srid_the_geom CHECK (srid(the_geom) = 3035);
alter table bdmap2009.station_geography ADD CONSTRAINT pk_id PRIMARY KEY(id);
CREATE INDEX indexstatgeo ON bdmap2009.station_geography
  USING GIST ( the_geom GIST_GEOMETRY_OPS );

-- SCRIPT POUR CEDRIC
/*
E:\IAV\eda\bdmap\geolocalisation
E:
CD E:\IAV\eda\bdmap\geolocalisation
psql -d eda2.0_RHT -h localhost -U postgres -p 5432 -f station_geography.sql
*/
-- psql erreur la fonction geom_update n'existe pas (Céline ca tu n'avais pas=> je n'ai pas le trigger)
alter table station_geography set schema bdmap2009;
---Changement du système de projection  RGF93_Lambert_93 (srid 2154) --> ETRS LAEA (srid 3035)
alter table bdmap2009.station_geography drop constraint enforce_srid_the_geom;
---update bdmap2009.station_geography set the_geom =ST_Transform(the_geom,3035);
update bdmap2009.station_geography set the_geom =ST_Transform(ST_SetSRID(the_geom,2154),3035);--11379 lines
-- pour moi geom_update n'a pas marché, il doit bien s'agir d'un problem de trigger, fait drop trigger
alter table bdmap2009.station_geography add  CONSTRAINT enforce_srid_the_geom CHECK (srid(the_geom) = 3035);
CREATE INDEX index_station_geography ON bdmap2009.station_geography
  USING btree ( id_bdcarthage_troncon);
select * from bdmap2009.station_geography where a_conserver = true

alter table bdmap2009.station_geography add column id_drain integer;
update bdmap2009.station_geography set id_drain=sub.id_drain from (select id_drain, id_bdcarth from rht.rht_bdcarthage) as sub where id_bdcarthage_troncon=id_bdcarth;  ---11070 lines

Last modified 7 years ago Last modified on Jun 1, 2018 5:57:54 PM

Download in other formats:

Plain Text