Localisation et caractérisation des têtes de bassin versant

Introduction[modifier]

Stage - A. Coudart 2018

Travail en cours (maj du 24 juillet 2018)

Bibliographie des TBV[modifier]

Lien vers la page publique

Fonctionnalités et services rendus par les têtes de bassin versant[modifier]

TODO

Contexte réglementaire[modifier]

La source de cette étude trouve son origine dans les textes réglementaires qui encadrent les actions territoriale en matière de gestion des milieux aquatiques. En premier lieu, rappelons que la France, comme tous les état membres de l'union européenne, a pour devoir de se conformer aux objectifs de la DCE (Directive Cadre sur l'Eau). Ces objectifs concernent une atteinte du "bon état écologique" des cours d'eau à différents pas de temps. Ainsi pour accompagner les acteurs de l'eau dans cette démarche des documents officiels font office de feuille de route. A l'échelle macro (Loire/Bretagne, territoire de l'Agence de l'Eau Loire Bretagne), nous trouverons le SDAGE (Schéma Directeur d'Aménagement et de Gestion des Eaux) et à l'échelle micro nous trouverons les SAGEs (Schéma d'Aménagement et de Gestion des Eaux) qui transcrivent le SDAGE dans un cadre local en l’occurrence la Sèvre Nantaise dans cette étude.

Ce que dit le SDAGE[modifier]

Le SDAGE Loire Bretagne 2016-2021, consacre un chapitre à part entière à la question des têtes de bassins versants. Mais cette thématique est présente dans d'autres chapitres. Ainsi dès la présentation générale du SDAGE il y fait référence dans les SRCE (Schéma Régionaux de Cohérence Ecologique), où elles sont décrites comme jouant un rôle de réservoirs de biodiversité.
Dans le chapitre 1, qui introduit la nécessaire refonte des aménagements de cours d'eau, présente les têtes de bassin versant comme des milieux sensibles, dont le bon fonctionnement est important pour l'ensemble du bassin et sur lesquels il est nécessaire de porter une attention particulière.
Dans le chapitre 8 concernant la préservation des zones humides et plus précisément dans la partie 8A-2, les têtes de bassins versant font parties des zones ciblées comme présentant un intérêt particulier en matière d'acquisition foncière.
Dans le chapitre 9 concernant la préservation de la biodiversité aquatique, il apparaît que les têtes de bassins versants peuvent faire l'objet de mesure particulière dans les plans départementaux pour la protection du milieu aquatique et la gestion des ressources piscicoles (PDPG).

Le chapitre 11 qui lui est spécialement consacré en 2 volets, 11A – Restaurer et préserver les têtes de bassin versant et 11B – Favoriser la prise de conscience et la valorisation des têtes de bassin versant. Elles y sont décrites comme "notre capital hydrologique". La disposition 11A-1 impose aux SAGEs un inventaire systématique des zones de têtes de bassin versant ainsi qu'une analyse des caractéristiques notamment écologiques et hydrologiques. C'est dans ce chapitre qu'une définition est donnée, elles sont décrites comme les bassins versants des cours d'eau dont le rang de Strahler est inférieur ou égal à 2 et dont la pente est supérieur à 1%. Ce critère de pente peut-être adapté localement pour les cours d'eau à faible puissance spécifique et présentant un risque de non atteinte des objectifs environnementaux. Ce dernier critère ne trouve pas d'appui dans la bibliographie scientifique. La disposition 11A-2 préconise que les SAGEs hiérarchisent les têtes de bassin versant en fonction des pressions et de l'état des masses d'eau pour définir des objectifs et principes de gestion adaptés.

Ce que dit le SAGE[modifier]

Dans le SAGE du bassin de la Sèvre Nantaise validé par arrêté préfectoral en 2015, les préconisation du SDAGE sont reprises et font l'objet de la disposition 46-1 : Améliorer la connaissance des têtes de bassin versant. Cette disposition reprend la travail de pré localisation validée par la CLE (Commission Locale de l'Eau) le 11 juillet 2013 et prévoit une étude précisant la localisation, la caractérisation des fonctionnalités aux regard des enjeux associés (qualité d'eau, quantité, biodiversité, morphologie), ainsi que le niveau de vulnérabilité et les modes de gestion. Il est également précisé que la méthode de travail doit-être établie en concertation avec les bassins versants limitrophe dans un soucis d'homogénéité.

Les leviers réglementaires[modifier]

La délimitation des têtes de bassin versant et leur caractérisation n'a pas pour finalité de rajouter un nouveau cadre réglementaire (MAMAN L., 2017 -- AUGIER A., 2017 -- Entretiens avec la DDTM 35 et l'Agence de l'Eau Loire Bretagne, délégations Centre et Bretagne -- Retours d'expérience TBV v1.8 de l'AFB). Comme énoncé plus tôt, leur objectif est de proposer un découpage fin du territoire pour aider les différents acteurs de l'aménagement du territoire à rationaliser leur actions visant la reconquête de la qualité de l'eau. Les têtes de bassins versant, que ce soit les cours d'eau ou leurs zones humides associées sont déjà concernées par la législation en vigueur (Loi sur l'eau, LEMA, 2006). Toutefois, leur protection est associée à leur cartographie, hors ce sont précisément sur ces cours d'eau de faibles dimensions et intermittents que la cartographie fait défaut (NGUYEN Van R., 2012). L'état mène actuellement par le biais de ces services décentralisés (DDT), une campagne de cartographie des cours d'eau. Cette dernière à une portée juridique forte puisque c'est sur cette base que la protection lié à la LEMA s'appliquera ou non. Dans le cas ou un rang 1 ne serait pas cartographié, il perdrait ainsi sa protection juridique. De plus les rangs 0 ou zones de sources (BENDA et AL, 2005), hors périmètre de zone humide inscrit au PLU(i) ou autre protection communautaire (Natura 2000...), ne sont aujourd'hui pas concernés par les mesures de protection alors qu'ils sont des éléments indissociables des cours d'eau en tête de bassin (LE BIHAN, 2009).

En conséquence et au regard des enjeux concernant la protection de ces milieux et des pressions dont ils font l'objet, il apparaît intéressant de pouvoir proposer au pouvoirs publics des pistes politiques. En matière d'urbanisme, les communes et intercommunalités possèdent un levier majeur avec les PLU(i), la délimitation et la caractérisation des têtes de bassins versant peuvent apporter des éléments de connaissance pouvant faciliter la prise de décisions, comme par exemple des sites sensibles à protéger de l’imperméabilisation, ou des sites de réservoirs de biodiversité. Ainsi sans rentrer dans une politique de protection globale de tête de bassin versant, il est tout à fait possible de jouer sur des leviers conditionnant l'aménagement du territoire à la prise en compte des services rendus par ces espaces.

Quelques exemples :

• Urbanisation possible avec traitement des eaux pluviales en infiltration à la parcelle.
• Compensation environnementale si destruction des services rendus par une tête de bassin versant.
• Identification de zone de source
• Protection d'éléments bocager en secteur à fort risque érosif
• ...

Solutions et incitation réglementaires expérimentées sur le bassin Loire Bretagne[modifier]

• La charte de l'agriculture et de l'urbanisme du Morbihan -- par le département du Morbihan (depuis 2008). Cette charte incite les communes et intercommunalités à traduire dans leur documents d'urbanisme un règlement "interdisant toutes constructions à 35 m des cours d’eau et tout comblement, affouillement, exhaussement de terrains (sauf création de retenues d’irrigation autorisées dans le respect de la loi sur l’eau et mise aux normes environnementales autorisée)."

• Le SCOT 2011 -- Communauté d'agglomération de la presqu'île de Guérande (CAP Atlantique). L'objectif 1.2.2 Protéger les milieux humides et les abords des cours d'eau, contiens une prescription concernant les zones de sources (objet d'une étude interne consistant à élaborer une méthodologie d'identification de l'origine des rangs 1) : "Les collectivités préciseront les zones de sources effectives afin de :
  • les intégrer dans la trame verte et bleue des documents d'urbanisme locaux
  • s'assurer que les travaux d'aménagements ne suppriment pas ces zones ni ne les déconnectent des cours d'eau ou, en cas de reconfiguration des zones de sources (notamment en espace urbain existant ou programmé), que ces travaux intègrent les moyens de préserver l'alimentation des cours d'eau ou d'améliorer leur fonctionnement hydraulique. "

• Le SAGE Couesnon -- (depuis 2013). La délimitation des têtes de bassin versant à été inscrite dans le règlement et le PAGD du SAGE Couesnon. Il y est indiqué que les travaux sur les cours de tête de bassin versant sont interdits sur le bassin versant du Couesnon (rangs 1 et 2 avec 1% de pente). Plus précisément :
  • La mise en place d’ouvrage faisant obstacle à l’écoulement de l’eau et à la migration des poissons (chute supérieure à 20 cm).
  • Le reprofilage du cours d’eau.
  • Le busage d’un cours d’eau.
  • L’artificialisation des berges (à partir de 20 m).

C'est à ce jour, le seul exemple connu d'intégration effective des têtes de bassin versant dans un règlement de SAGE.

• Les zones humides prioritaires -- Étude menée par le Forum des Marais (FMA). Étude et éléments de compréhension sur la portée réglementaire des dispositifs ZHIEP (Zone Humide d'Intérêt Environnemental Particulier) et ZSGE (Zones humides Stratégiques pour la Gestion de l'Eau). A consulter ici.
  

Localisation des têtes de bassin versant[modifier]

La première étape de l'étude des têtes de bassin versant est logiquement de les localiser, donc définir leur périmètre et leur exutoire. Selon la directive cadre sur l'eau, un bassin versant est l'espace drainé par un cours d'eau et ses affluents. L'ensemble des eaux qui tombent dans cet espace convergent vers un même point de sortie appelé exutoire : cours d'eau, lac, mer, océan, etc..

Référentiels hydrographiques[modifier]

Le choix du référentiel hydrographique est une donnée majeur dans le travail d'étude des têtes de bassin versant car de fait, il conditionne leur localisation et le travail de recoupement de données qui en découle. Son choix est donc, à juste titre, la première pierre à poser pour édifier une cartographie cohérente des têtes de bassin versant.

Les référentiels disponibles :

• Bd Carthage
• Bd Topo V2.1
• Bd Topo 151
• Inventaire DDTM
• Tracé(s) théorique(s)

Tableau comparatif du linéaire hydrographique en fonction des référentiel[modifier]

Chaque référentiel possède ces points forts et points faibles. L'échelle de numérisation et les objectifs n'étant pas les mêmes le linéaire visible est très disparate. Comparatif sur le bassin versant de la Sèvre Nantaise avec un zoom sur le département 49 :

Décryptage du tableau :[modifier]

La Bd Carthage est une donnée réalisé à petite échelle sur tout le territoire français (Validité au 1:100 000), à l'échelle locale sa fiabilité est donc toute relative, c'est pourquoi c'est le référentiel qui nous donne le plus petit niveau de linéaire de CE (Cours d'Eau).
La différence entre les 2 tracés théoriques s'explique par une différence notable entre les paramètres utilisés. Pour le travail de Moreau, 2016, aucune valeur d'accumulation dans les cellules n'avait été retenue. Pour le travail de Coudart, 2018, une valeur minimale de cellule d'accumulation de 4 000 à été retenue. Cette valeur équivaut à 10 hectares sur le terrain, au dessous, il est considéré que la puissance spécifique n'est pas à même de structuré un cours d'eau (Cirou, 2017).
Le flux de la bd topo, tronçon cours d'eau, est la donnée la plus récente du référentiel IGN à grande échelle. Le référentiel qui a été fournit aux services de l'état dans le cadre de la démarche d'inventaire des cours d'eau (la version bd topo 151) se rapprochait plus de la donnée actuel (flux) que de la version 2.1. A titre de comparaison le calcul sur la version 151 nous donnait 3387 km de cours d'eau.

Script de calcul SQL[modifier]

Le calcul à été réalisé par requêtes SQL sur une base PostreSQL / Postgis.

select 'Carthage' as couche 
	, 'BVSN' as zone
	, (sum(ST_Length(xx_99_utils.st_intersection_btr(a.geom, b.geom, false)))/1000)::int as km 
from r020_territoire_hydrographie.t_troncon_hydrographique_actu_carthage a
	, r020_territoire_hydrographie.t_bvsn_millesime_2014_carthage b 
where st_intersects(a.geom, b.geom)
union all 
select 'Flux IGN troncon cours deau bd_topo' as couche
	, 'BVSN' as zone
	, (sum(ST_Length(xx_99_utils.st_intersection_btr(a.geom, b.geom, false)))/1000)::int as km 
from xx_99_50_10_import_wfs.t_ign_bdtopo_troncon_ce a
	, r020_territoire_hydrographie.t_bvsn_millesime_2014_carthage b 
where st_intersects(a.geom, b.geom)
union all
select 'IGN v2.1' as couche
	, 'BVSN' as zone
	, (sum(ST_Length(xx_99_utils.st_intersection_btr(a.geom, b.geom, false)))/1000)::int as km 
from r020_territoire_hydrographie.t_hydro_troncon_cours_eau_ign_bdtopo_v2_1 a
	, r020_territoire_hydrographie.t_bvsn_millesime_2014_carthage b 
where st_intersects(a.geom, b.geom)
union all
select 'Tracé théorique Colin' as couche
	, 'BVSN' as zone
	, (sum(ST_Length(xx_99_utils.st_intersection_btr(a.geom, b.geom, false)))/1000)::int as km 
from r020_territoire_hydrographie.t_troncon_cours_eau_theorique_mnt_5x5 a
	, r020_territoire_hydrographie.t_bvsn_millesime_2014_carthage b 
where st_intersects(a.geom, b.geom)
union all
select 'Tracé théorique Anthony' as couche
	, 'BVSN' as zone
	, (sum(ST_Length(xx_99_utils.st_intersection_btr(a.wkb_geometry, b.geom, false)))/1000)::int as km 
from _agent_acoudart.t_mnt_sc_sk_stream_4000 a
	, r020_territoire_hydrographie.t_bvsn_millesime_2014_carthage b 
where  a.strahler::int >= 1 and st_intersects(a.wkb_geometry, b.geom)

union all

select 'Carthage' as couche
	, 'BVSN 49' as zone 
	, (sum(ST_Length(xx_99_utils.st_intersection_btr(a.geom, d.geom, false)))/1000)::int as km 
from r020_territoire_hydrographie.t_troncon_hydrographique_actu_carthage a
	, r020_territoire_hydrographie.t_bvsn_millesime_2014_carthage b 
	, r010_administratif_reglementaire.t_departements_2xxx_bd_topo d
where si_admin_metad_version_date = '2016-10-20'::date and d.depart = 'MAINE-ET-LOIRE' and st_intersects(a.geom, d.geom) and st_intersects(a.geom, b.geom)
union all
select 'Flux IGN troncon cours deau bd_topo' as couche
	, 'BVSN 49' as zone
	, (sum(ST_Length(xx_99_utils.st_intersection_btr(a.geom, d.geom, false)))/1000)::int as km 
from xx_99_50_10_import_wfs.t_ign_bdtopo_troncon_ce a
	, r020_territoire_hydrographie.t_bvsn_millesime_2014_carthage b
	, r010_administratif_reglementaire.t_departements_2xxx_bd_topo d 
where si_admin_metad_version_date = '2016-10-20'::date and d.depart = 'MAINE-ET-LOIRE' and st_intersects(a.geom, d.geom)and st_intersects(a.geom, b.geom)
union all
select 'IGN v2.1' as couche 
	, 'BVSN 49' as zone 
	, (sum(ST_Length(xx_99_utils.st_intersection_btr(a.geom, d.geom, false)))/1000)::int as km 
from r020_territoire_hydrographie.t_hydro_troncon_cours_eau_ign_bdtopo_v2_1 a
	, r020_territoire_hydrographie.t_bvsn_millesime_2014_carthage b 
	, r010_administratif_reglementaire.t_departements_2xxx_bd_topo d
where si_admin_metad_version_date = '2016-10-20'::date and d.depart = 'MAINE-ET-LOIRE' and st_intersects(a.geom, d.geom) and st_intersects(a.geom, b.geom)
union all
select 'DDT 49' as couche 
	, 'BVSN 49' as zone 
	, (sum(ST_Length(xx_99_utils.st_intersection_btr(a.geom, d.geom, false)))/1000)::int as km 
from r020_territoire_hydrographie.t_cours_eau_police_de_l_eau_ddtm a
	, r020_territoire_hydrographie.t_bvsn_millesime_2014_carthage b 
	, r010_administratif_reglementaire.t_departements_2xxx_bd_topo d
where si_admin_metad_version_date = '2016-10-20'::date and d.depart = 'MAINE-ET-LOIRE' and st_intersects(a.geom, d.geom) and st_intersects(a.geom, b.geom)
union all
select 'Tracé théorique Colin' as couche 
	, 'BVSN 49' as zone 
	, (sum(ST_Length(xx_99_utils.st_intersection_btr(a.geom, b.geom, false)))/1000)::int as km 
from r020_territoire_hydrographie.t_troncon_cours_eau_theorique_mnt_5x5 a
	, r020_territoire_hydrographie.t_bvsn_millesime_2014_carthage b 
	, r010_administratif_reglementaire.t_departements_2xxx_bd_topo d
where si_admin_metad_version_date = '2016-10-20'::date and d.depart = 'MAINE-ET-LOIRE' and st_intersects(a.geom, d.geom) and st_intersects(a.geom, b.geom)
union all
select 'Tracé théorique Anthony' as couche
	, 'BVSN 49' as zone
	, (sum(ST_Length(xx_99_utils.st_intersection_btr(a.wkb_geometry, d.geom, false)))/1000)::int as km 
from _agent_acoudart.t_mnt_sc_sk_stream_4000 a
	, r020_territoire_hydrographie.t_bvsn_millesime_2014_carthage b
	, r010_administratif_reglementaire.t_departements_2xxx_bd_topo d 
where si_admin_metad_version_date = '2016-10-20'::date and d.depart = 'MAINE-ET-LOIRE' and st_intersects(a.wkb_geometry, d.geom)  and st_intersects(a.wkb_geometry, b.geom) and a.strahler::int >= 1

Choix méthodologiques[modifier]

Comme le montre la bibliographie existante sur le sujet. De nombreux outils permettent de réaliser les traitements nécessaires à la correction ou à la création de référentiels hydrographique. Le choix a été fait de recourir à la technologie GRASS. D'une part GRASS est un système d'information géographique historique, robuste, fiable, il a largement fait ces preuves. D'autres parts l'ensemble des traitements dont nous avons besoin sont disponible sous forme de fonctions GRASS, ce qui nous permet d'effectuer l'ensemble de la chaîne de traitements avec la même technologie. Enfin, GRASS pouvait être monter sur le serveur distant de la structure afin d'interagir avec la base de donnée Postgres (elle aussi présente sur le serveur distant). L'alliance GRASS / Postgres est interfacée sur des pages web en php ayant pour résultat une interaction ergonomique et surtout une capacité de traitement basé sur celle du serveur (64 go RAM) et non sur les machines locales de la structure.

GRASS GIS 7.4 (ou antérieures) peu être téléchargé sur cette page.

L'étude des têtes de basin versant suit la logique méthodologique suivante :

• Localisation des têtes de bassin versant
  • Préparation du MNT
  • Préparation du référentiel hydrographique
  • Ordination des cours d'eau selon le rang de Strahler
  • Création des bassins versants

• Caractérisation des têtes de bassin versant

Les extensions Grass[modifier]

Les extensions sont des modules complémentaires qui sont téléchargeables très simplement et permettent d'ajouter des fonctions au logiciel de base. r.stream.order et r.stream.basins en font parties.
L'ensemble des modules complémentaires sont disponibles ici

Installation d'une extension[modifier]

Installer une extension est très simple, g.extension est la fonction qui permet cette opération. Il est possible d'installer une extension dans l'interface graphique de Grass ou via la console. La fonction télécharge et installe l'extension appelée en paramètre, à la manière d'un 'APT GET' sur linux.

Exemple d'appel de la fonction g.extension dans la console GRASS[modifier]

g.extension extension=r.stream.order

Les extensions "Stream"[modifier]

Il existe une série d'extensions "stream" qui regroupe plusieurs fonctions d'ordre hydrographique.

Introduction à GRASS GIS[modifier]

Avant de rentrer dans l'utilisation d'un outil il est toujours intéressant d'en connaître l'identité. GRASS GIS est un logiciel de SIG relativement atypique, et mérite que l'on s'y attarde un peu. L'origine du projet remonte à 1982 et a été institué par différentes agences gouvernementales et universités américaines, le cœur de l'application étant développé par l'armée américaine elle-même. L'objectif initial était d'aider à l'aménagement du territoire. Sans rentrer dans les détails techniques, on peu préciser que le développement est continu depuis 1982, la dernière version stable en date étant la 7.4 (1 fev. 2018), il est open source, libre et construit sur une base modulaire, ce qui signifie que chaque fonction est un module à part entière. Si GRASS offre la possibilité de travailler avec des vecteurs, sa spécialité principale reste les rasters.
Si on cherchait à comparer ces fonctionnalités, il se rapprocherait de Spatial Analyst ou 3D Analyst d'Arcgis, son gros point positif étant qu'il est gratuit et tout aussi efficace. En revanche sa prise en main peu demander une phase d'adaptation. Car si une partie de sa boîte à outils est disponible via Qgis, pour profiter pleinement de toutes les fonctions il est nécessaire de passer directement par sa stand-alone ou sa console.

GRASS diffère de la plupart des autres SIG dans le fait qu'il est nécessaire de constituer une base de données propre pour travailler différents fichiers. Il n'est pas possible de glisser déposer un .tif ou un .shp à la volée par exemple. Toutefois une fois ce pré-requis intégrer, il devient d'une grande efficacité. Pour aider à la compréhension, un tutoriel en français est disponible sur le site portail sig, bien que datant de 2012, l'article reste pertinent pour appréhender l'outil. Les modules complémentaires ajouter au fur et à mesure son parfaitement décrit dans la documentation de GRASS et permettent de rester à jour.

Les 2 éléments à retenir pour démarrer avec GRASS sont la 'LOCATION' et le 'MAPSET'. La LOCATION correspond à la zone géographique de travail et à un système de projection (SRID) particulier, on pourra donc changer de LOCATION au choix en changeant de projet ou en changeant de zone géographique d'étude. Le MAPSET correspond lui au jeux de données. Par défaut un MAPSET PERMANENT est créer à la création d'une LOCATION.

Il fois la LOCATION crée et le MAPSET choisit, il est possible de travailler via l'interface graphique ou par la console. Ci-dessous vous trouver les lignes de commandes permettant d'appeler les fonctions via la console de GRASS. (ces lignes de commandes commençant par le nom du module il sera facile, le cas échéant, de les retrouver dans l'interface graphique).

Attention ! les exemples devront être adaptés à votre LOCATION et vos MAPSET.

Démarrage de Grass[modifier]

Pour les modalités de démarrage de l'application se reporter au tutoriel disponible sur le portail sig précédemment cité.

Préparation du MNT[modifier]

La deuxième étape du travail de localisation des têtes de bassin versant, après le choix et le nettoyage topologique du référentiel est la préparation du Modèle Numérique de Terrain. Cette modélisation 3D du territoire nous sera utile pour les calculs de pente, d'accumulation, d'écoulement... Celui-ci étant issu d'une construction algorithmique il est indispensable de le retravailler pour qu'il soit valide. Les différentes étapes sont décrites ci-dessous.

Creusement du MNT[modifier]

Dans l'optique de construire un réseau d'écoulement théorique basé sur les tracés IGN des cours d'eau, il est nécessaire de passé par une série d'étapes. Parmi lesquelles le surcreusement du MNT. Le module r.mapcalc sera utilisé pour la partie calcul de creusement et de surélévation, v.to.rast sera utilisé pour transformer le vecteur en raster, r.null sera utilisé pour transformer les valeurs NULL du raster en 0 afin de ne pas interférer dans le creusement ultérieur.

Attention ! Si l'étape de "lissage" des dépressions est indispensable, l'étape de creusement est tout à fait optionnel. Dans notre étude présente, son utilité réside dans son utilisation ultérieure. Ce réseau théorique forcé en IGN nous permettra, comme nous le verrons plus loin, de combler les zones de déconnexion du référentiel IGN. De même l'étape de surélévation avant l'étape de creusement permet dans des zones proches du niveau de la mer de ne pas générer des valeurs négatives lors du surcreusement. Sur le territoire de la Sèvre Nantaise, cette étape est indispensable au surcreusement.

• Transformation du vecteur IGN en raster (inscription de la valeur 50m en dur dans les données du raster grâce à l'option "value=50"):

v.to.rast input="ign_bd_topo_151_bvsn@PERMANENT" layer=ign_bd_topo_151_bvsn type="line" output="ign_bd_topo_151_bvsn" use="val" value=50

• Transformation des valeur NULL du raster cours d'eau :

r.null map="ign_bd_topo_151_bvsn@PERMANENT" null=0

• Surélévation du MNT brut (pour éviter les valeurs négatives) :

r.mapcalc expression="mnt_sangueze_sureleve = mnt_sangueze_brut@PERMANENT + 50" --overwrite

• Creusement du MNT surélevé par rapport au référentiel hydrographique (-50m) :

r.mapcalc expression="mnt_sangueze_creuse = mnt_sangueze_sureleve@PERMANENT - ign_bd_topo_151_bvsn@PERMANENT" --overwrite

Comblement des dépressions[modifier]

Pour rappel le MNT est un raster construit par interpolation des points composant le nuage de points originel, ce sont donc des algorithmes qui lient les points entre eux, créant ainsi de potentiel zones de surélévation ou de dépressions issu de l'interprétation. Afin de corriger cet état de fait et d'obtenir un MNT topologiquement correct, ou du moins exploitable en vue de calcul des écoulements, il est nécessaire de combler les dépressions. Le module r.fill.dir sera utilisé. Dans l'exemple ci-dessous nous nous servirons du MNT creusé précédemment. Le module boucle sur le MNT jusqu'à épuisement des dépressions.

Attention ! Il est à noter que le module ne peu faire la différence entre zone de dépressions réelle ou issu de l'interpolation. Il est possible dans les options de rajouter une couche de dépressions réelles (lac, étang, zones humides...)

• Comblement des zones de dépressions :

r.fill.dir input=mnt_sangueze_creuse@PERMANENT output=mnt_sangueze_creuse_elevation direction=mnt_sangueze_creuse_direction --overwrite

Construction d'un réseau théorique d'écoulement[modifier]

Un réseau théorique d'écoulement est une donnée générée à partir des informations d'élévation d'un MNT (Modèle Numérique de Terrain). Il simule le chemin de l'eau en fonction des valeurs d'accumulation de chaque cellule.
Grass Gis offre de nombreux modules hydrographiques performants permettant de réaliser les opérations qui suivent. De plus, il est possible dans appeler les fonctions au travers de traitements Postgrès.
L'outil GRASS r.stream.extract permet de générer un raster de direction de flux, un raster de cours d'eau orienté de l'amont vers l'aval et segmenter entre chaque confluence, un vecteur de cours d'eau. Pour générer le vecteur de cours d'eau, en plus d'une donnée d'élévation en entrée, une donnée d'accumulation est nécessaire. Cette dernière peut-être obtenue à partir de r.watershed.

MNT[modifier]

• RGE Alti, IGN, 5m x 5m
• RGE Alti sinkless12 (MNT corrigé en supprimant les dépressions avec l'outil GRASS r.sink)

Fonctions GRASS[modifier]

• r.stream.extract
• r.watershed

Construction du flux d'accumulation[modifier]

• Exemple d'appel de la fonction r.watershed dans la console GRASS

r.watershed -b elevation=RGE_Alti_5mx5m_60cm_2015_mnt_sinkless_test_ouin@PERMANENT accumulation=_watershed_accu_test_ouin_1

Construction du réseau d'écoulement théorique[modifier]

• Exemple d'appel de la fonction r.stream.extract dans la console GRASS

r.stream.extract elevation=RGE_Alti_5mx5m_60cm_2015_mnt_sinkless_test_ouin accumulation=_watershed_accu_test_ouin_1 threshold=4000 stream_raster=rast_stream_ouin_2 stream_vector=vect_stream_ouin_2 direction=rast_dir_ouin_2

• Le paramètre "threshold" dont la valeur à été initialisé à 4000 dans l'exemple ci-dessous renvoi à la valeur d'accumulation minimale pour qu'un cours d'eau soit créé (en l’occurrence, 4000 cellules doivent s'écouler vers un points avant de former un cours d'eau, étant donné qu'une cellule mesure sur le terrain 25 m2, cela équivaut à un bassin versant minimal de 10 ha (cette valeur à été identifié dans le travail de Julien Cirou et validé par l'AFB)

Export d'un raster de Grass vers GDAL[modifier]

• Pour la suite des traitements, il est indispensable de passer les données dans une base Postgres.

r.out.gdal input=test_bvsn_MNT_sinkless output=/var/www/observ/donnees-brutes/IGN_MNT_RGE/RGE_Alti_5mx5m_60cm_2015_mnt_sinkless_grass.tif

Ordination des cours d'eau selon le rang de Strahler[modifier]

Les rangs de Strahler sont obtenus grâce à une méthode d’arborescence. Cette classification est très populaire pour les réseaux hydrographiques et permet d'identifier la complexité des réseaux d'affluents et de sous affluents.
Un segment de cours d'eau est représenté par une ligne et des nœuds qui le lie au segment suivant et/ou précédent. Ainsi le cours d'eau de rang 1 est celui qui trouve son origine dans une source jusqu'à la première confluence à l'aval. Par la suite :

• Si il rencontre un cours d'eau de rang 1, ils se transformerons en cours d'eau de rang 2.
• Si il rencontre un cours d'eau de rang 2 ou supérieur, les segments suivants ne seront pas impactés
  

Voir schéma ci-dessous (TODO)

MNT[modifier]

• RGE Alti, IGN, 5m x 5m
• RGE Alti sinkless12 (MNT corrigé en supprimant les dépressions avec l'outil GRASS r.sink)

Fonction GRASS[modifier]

• r.stream.order

Exemple d'appel de la fonction r.stream.order dans la console GRASS[modifier]

• "stream_vect=" correspondant à la sortie vecteur du calcul, "strahler=" correspondant à la sortie Strahlerisée (en raster) du calcul

r.stream.order stream_rast=unique_stream_from_stream_extract_3 direction=flow_direction_from_stream_extract_3 stream_vect=stream_vect_from_stream_order strahler=riverorder_strahler_test_1

Correction d'un référentiel hydrographique[modifier]

A Revoir

Si le réseau théorique nous donne une information intéressante en matière de cours d'eau supposé et malgré sa construction basé sur la physique, il reste peu fiable car ne prend pas en compte les logiques pédologiques, géologiques, les divers aménagements anthropiques comme les fossés, enterrement de cours d'eau, dérivation, etc, mais seulement les logiques gravitaires (et selon des modèles d'accumulation). Il doit donc être forcément mis en relation avec les autres référentiels hydrologique, comme celui de la BD topo ou les inventaires de cours d'eau locaux.

Si l'on souhaite utiliser la BD topo ou tout autres référentiels local, il sera indispensable de corriger ces derniers afin de les rendre exploitables par un traitement SIG. En effet, les erreurs topologiques, les sens d'écoulement aberrants, ou les "boucles" (ex : cours d'eau anastomosés (ou tressés), deltas ou bief de moulin...) seront autant de biais qui rendront le traitement impossible ou pire, invalide.

Ces questions sont au cœur de toutes les démarches analogues, jusqu'à aujourd'hui l'automatisation de la correction des référentiesl n'était pas connu. Chacun réalisant des corrections mi-manuelles, mi-automatisées, le travail réalisé pendant cette étude permet à partir de script SQL de corriger automatiquement 99 % des erreurs.

La logique des corrections est la suivante :

• Comblement des "trous" dans le référentiel
• Orientation amont -> aval du référentiel
• Suppression des "boucles" du référentiel
• Ordination du référentiel

Correction topologique avec postgis[modifier]

Avant de lancer le script de correction automatique, il est indispensable d'effectuer quelques vérifications. Le référentiel "tronçon cours d'eau" de la BD topo de l'IGN comporte de nombreuses erreurs. Idem pour les segments routiers, Makina corpus propose un tutoriel pour corriger la topologie d'un graph avec postgis.

• Le script ci-dessous est un exemple d'adaptation au cours d'eau sur la base de la Sèvre nantaise :

drop table  a000_temp.ign_ce ;

CREATE TABLE a000_temp.ign_ce as

 

SELECT gid, ogc_fid, gml_id, id, prec_plani, prec_alti, artif, fictif,

       franchisst, nom, pos_sol, regime, z_ini, z_fin, geom

  FROM rx000_data_referentiel_mutualise.t_ign_bdtopo_troncon_ce_bbox;

 

 

grant select on a000_temp.ign_ce TO grp_db_ext_mut_02_layon_structure;

 

--SELECT topology.DropTopoGeometryColumn('a000_temp', 't_ocs_epsn_v01', 'topo_geom');

select topology.DropTopology('ign_ce_topo');

SELECT topology.CreateTopology('ign_ce_topo', 2154);

SELECT topology.AddTopoGeometryColumn('ign_ce_topo', 'a000_temp', 'ign_ce', 'topo_geom', 'LINESTRING');

 

 

GRANT USAGE ON SCHEMA ign_ce_topo TO grp_db_ext_mut_02_layon_structure;

 

ALTER DEFAULT PRIVILEGES IN SCHEMA ign_ce_topo

    GRANT INSERT, SELECT, UPDATE, DELETE, TRUNCATE, REFERENCES, TRIGGER ON TABLES

    TO grp_db_ext_mut_02_layon_sysma_admin;

 

 

DO $$DECLARE r record;

BEGIN

  FOR r IN SELECT * FROM a000_temp.ign_ce LOOP

    BEGIN

      UPDATE a000_temp.ign_ce SET topo_geom = topology.toTopoGeom(st_force2d(geom), 'ign_ce_topo', 1, 1.0)

      WHERE gid = r.gid;

    EXCEPTION

      WHEN OTHERS THEN

        RAISE WARNING 'Loading of record % failed: %', r.gid, SQLERRM;

    END;

  END LOOP;

END$$;

 

--alter table a000_temp.ign_ce drop column geom_from_topogeom;

alter table a000_temp.ign_ce add column geom_from_topogeom geometry(multilinestring, 2154);

UPDATE  a000_temp.ign_ce set geom_from_topogeom = st_multi(topo_geom)::geometry(multilinestring, 2154)

 

);

Création de bassins versants[modifier]

A revoir

Les fonctions du module r.stream.basins de Grass permettent de passer à la commande de multiples paramètres qui nous permettent d'atteindre différents objectifs. Par exemple si on utilise tous les paramètres par défaut, le module nous renverra une cartographie des bassins-versants de chaque segment de cours d'eau. Si on souhaite préciser son analyse il sera nécessaire de lui passer de nouvelles informations qui seront obtenus à partir d'autres modules.

MNT[modifier]

• RGE Alti, IGN, 5m x 5m
• RGE Alti sinkless12 (MNT corrigé en supprimant les dépressions avec l'outil GRASS r.sink)

Fonctions GRASS[modifier]

• r.stream.basins
• r.reclass
• r.mapcalc
• r.to.vect

Exemple d'appel de la fonction r.stream.basins dans la console GRASS - [DEFAUT] micro-bv[modifier]

r.stream.basins direction=flow_direction_from_stream_extract_3  stream_rast=unique_stream_from_stream_extract_3  basins=bas_basins_elem_test_ouin_1

• l'appel ci-dessus prend tous les paramètres par défaut, il renvoi donc l'ensemble des micros bassins versants affiliés à chaque segment de cours d'eau.

Exemple d'appel de la fonction r.stream.basins dans la console GRASS - [Rang 2] Têtes de BV[modifier]

echo 2 = 2 > test2.txt | echo * = NULL >> test2.txt
r.reclass input=riverorder_strahler_test_1 output=riverorder_strahler_2_ouin_test_2_shell rules=test2.txt
r.stream.basins direction=flow_direction_from_stream_extract_3  stream_rast=unique_stream_from_stream_extract_3  basins=bas_basins_elem_test_ouin_1

• L'appel ci-dessus prend comme paramètre une variable supplémentaire qui permet de reclasser grâce au module r.reclass le raster en fonction de sa valeur, ici nous n'avons sélectionné que les cours d'eau de rangs 2 correspondant à la définition du SDAGE Loire Bretagne.
• Le résultat est un raster ne distinguant pas les id des différents cours d'eau.
• Si l'on souhaite un résultat qui différencie les différents cours d'eau les uns des autres alors il convient d'employé la méthode ci-dessous, qui ajoute l'option -c à la commande précédente et utilise le module r.mapcalc:
  

r.mapcalc "sel_strahler_cat_2 = if(riverorder_strahler_2_ouin_test_2_shell==2,riverorder_strahler_test_1,null())"
r.stream.basins -c dir=flow_direction_from_stream_extract_3 stream=sel_strahler_cat_2 basins=bas_basin_strahler_2_test_ouin_2

Exemple d'appel de la fonction r.stream.basins dans la console GRASS - [Rang 1] Têtes de BV[modifier]

• Calcul des bassins versants des cours d'eau de rang1

r.mapcalc "riverorder_strahler_ouin_s1_from_rmapcalc = if(riverorder_strahler_ouin_all_from_rstreamorder==1,riverorder_strahler_ouin_all_from_rstreamorder,null())"
r.stream.basins -c dir=direction_ouin_from_rstreamextract stream=riverorder_strahler_ouin_s1_from_rmapcalc basins=basin_stra1_ouin_from_rbasin

Transformation des bassins raster en vecteurs[modifier]

• Le module Grass r.to.vect prend en paramètre la couche raster en entée, le nom du fichier en sortie ainsi que le type (point, line, area).

r.to.vect input=basin_stra1_ouin_from_rbasin output=basin_stra1_ouin_vect type=area
r.to.vect input=basin_stra2_ouin_from_rbasin output=basin_stra2_ouin_vect type=area

Intégration des bassins vecteurs dans postgres[modifier]

L'étape d'assemblage des bassins versants étant réalisée grâce à une requête SQL, il est nécessaire d'importer les bassins versants dans la base Postgres. Pour le moment cette étape d'import nécessite une transformation des bassins raster en vecteurs shape, mais cette dernière devra à terme être remplacé par une fonction d'import direct dans la base sans passer par la forme shape.

Assemblage des bassins versants de rangs 1 et 2[modifier]

L'assemblage des bassins versants suit la logique suivante, prise en compte de tous les rangs 2 et union des rangs 1 si ils ne sont pas reliés directement à des rangs 2.

• Exemple d'appel de la fonction d'unification des bassins versant de rang 2 et 1 :

DROP table if exists _agent_acoudart.t_mnt_sc_sk_strahler_4000_bas_rk_1_2; 
CREATE table _agent_acoudart.t_mnt_sc_sk_strahler_4000_bas_rk_1_2 AS
 
WITH w_prepare as 
( 
SELECT ta.* 
FROM _agent_acoudart.t_mnt_sc_sk_strahler_4000_rank2_bas ta 
WHERE ta.value IS NOT NULL 

UNION ALL 

SELECT DISTINCT tb.* 
FROM _agent_acoudart.t_mnt_sc_sk_strahler_4000_rank1_bas tb 
LEFT JOIN _agent_acoudart.t_mnt_sc_sk_strahler_4000_rank2_bas ta ON ( ST_Within(tb.wkb_geometry,ta.wkb_geometry) AND ta.value IS NOT NULL) 
WHERE ta.ogc_fid IS NULL AND tb.value IS NOT NULL
) 
SELECT Row_number() over() ::bigint as gid_all, * from w_prepare ; 

Interface Web et fonction php[modifier]

Afin de rendre les traitements le plus intuitif possible, une surcouche d'appel via une interface web a été réalisée, elle a été implantée sur le serveur web distant de la structure processing.observatoire (accès limité aux agents de l'EPTB). Les commandes GRASS vues au-dessus ont été transformées en fonctions php pour permettre leur appel par une page Web. Ceci permet également de réaliser tous les traitements directement sur le serveur (avec sa puissance de calcul) plutôt que de passer par les machines locales.

• Exemple d'une fonction GRASS passée en php.
• Cette fonction, grass_r_stream_basins, correspond au module GRASS r.stream.basins.
  • Dans un premier temps on déclare toutes les variables qui seront appelées et qui sont paramétrables par l'utilisateur.
  • Puis une commande est lancé à la console. Cette commande se connecte à GRASS et une LOCATION spécifique
  • Ensuite la commande exécute la commande du module que nous avions préalablement testée et validée.
• Pour visualiser le code cliquer sur "Afficher".

<?php

/**
 * grass_r_stream_basins
 * 
 * @param string $grassUser TextType DEFAULT=GRASS_USER
 * @param string $grassMemoryMB TextType DEFAULT=GRASS_MEMORY
 * @param string $grassDataDir TextType DEFAULT='/var/www/observ/out_processing/grass_data'
 * @param string $grassLocation TextType DEFAULT='test_13'
 * @param string $grassMapset TextType DEFAULT='PERMANENT'
 * @param string $grassInputRasterDirectionName TextType DEFAULT='MNT_sinkless_direction'
 * @param string $grassInputRasterStreamCategoriseName TextType DEFAULT='sel_strahler_cat_2'
 * @param string $grassOutputRasterBasinsName TextType DEFAULT='MNT_sinkless_basins' 
 * 
 */
function grass_r_stream_basins($grassUser,$grassMemoryMB, $grassDataDir, $grassLocation, $grassMapset, $grassInputRasterDirectionName, $grassInputRasterStreamCategoriseName, $grassOutputRasterBasinsName) {
//$process_uuid=printf("uniqid(): %s\r\n", uniqid());


    $command = GRASS_EXEC . ' ' . $grassDataDir . '/' . $grassLocation . '/' . $grassMapset . ' --exec' .
            ' r.stream.basins dir=' . $grassInputRasterDirectionName . '@' . $grassMapset . ' stream=' . $grassInputRasterStreamCategoriseName . '@' . $grassMapset .
            ' basins=' . $grassOutputRasterBasinsName .
            ' memory='.$grassMemoryMB.
            ' --overwrite';

    $command = '"' . $command . '"';
    $command = 'sudo /bin/su - ' . $grassUser . ' -c ' . $command;

  //  echo '<BR>' . $command . '<BR>';

    return array('command' => $command, 'type' => 'exec');
}

Une fois les fonctions unitaire crées, des process des traitements permettent de réalisé un ensemble de traitement, ci-dessous la fonction hydro_strahler_rank_from_rast_mnt permet de créer un réseau d'écoulement théorique à partir d'un simple MNT non traité.

• La logique est la suivante :
  • On commence par déclarer l'ensemble des variable nécessaire au traitement qui seront autant de variables paramétrables par l'utilisateur.
  • Puis on appelle les fonctions une par une.
• Pour visualiser le code cliquer sur "Afficher".

<?php

/**
 * hydro_strahler_rank_from_rast_mnt
 * 
 * @param SilexApp $app
 * @param string $pghost TextType DEFAULT=PGHOST
 * @param string $pgbase TextType DEFAULT=PGBASE
 * @param int $pgPort NumberType DEFAULT=PGPORT
 * @param string $pgUser TextType DEFAULT=PGUSER
 * @param string $pgPass TextType DEFAULT=PGPASS
 * @param string $destSchemaName TextType DEFAULT='_agent_acoudart'
 * @param integer $epsg NumberType DEFAULT='2154'
 * @param string $grassUser TextType DEFAULT=GRASS_USER
 * @param string $grassMemoryMB TextType DEFAULT=GRASS_MEMORY
 * @param string $grassDataDir TextType DEFAULT='/var/www/observ/out_processing/grass_data'
 * @param string $grassLocation TextType DEFAULT='bvsn_surcreuse'
 * @param string $grassMapset TextType DEFAULT='PERMANENT'
 * @param string $grassInputRastMNT TextType DEFAULT='/var/www/observ/donnees-brutes/IGN_MNT_RGE/RGE_Alti_5mx5m_60cm_2015.tif' 
 * @param string $grassRastMNTAlreadyInGrassLocation TextType DEFAULT='mnt_sc' 
 * @param string $grassInputStreamThreshold NumberType DEFAULT='500'
 */

function hydro_strahler_rank_from_rast_mnt($app, $pghost, $pgbase, $pgPort, $pgUser, $pgPass,$destSchemaName, $epsg, 
        $grassUser, $grassMemoryMB,$grassDataDir,$grassLocation,$grassMapset,$grassInputRastMNT,$grassRastMNTAlreadyInGrassLocation,$grassInputStreamThreshold) {

        $result = null ;
        if ($grassRastMNTAlreadyInGrassLocation=='') { // il n'y a pas encore de raster MNT dans GRASS , il faut donc créer une location et ajouter le raster
		$LOC_rastMNTName 		=$grassLocation.'_MNT';  // LOC_ pour variable locale
                
                // Création de la location et import du MNT
                $result .= $app['twig']->render('/app/process-res3.twig',[ 'res' => execute($app, 'grass_location_init_from_rast',
                            [$grassUser, $grassDataDir,$grassLocation,$grassMapset,$grassInputRastMNT])]);
                
                $result .= $app['twig']->render('/app/process-res3.twig',[ 'res' => execute($app, 'grass_rast_import',
                            [$grassUser, $grassDataDir,$grassLocation,$grassMapset,$grassInputRastMNT, $LOC_rastMNTName])]);
             
        } else   {$LOC_rastMNTName=$grassRastMNTAlreadyInGrassLocation; }
           
		
        
            $LOC_rastMNTSinkLessName 	=$LOC_rastMNTName.'_sk';
            $LOC_rastAccumulationName 	=$LOC_rastMNTSinkLessName.'_accumulation';
            $LOC_rastStreamName 		=$LOC_rastMNTSinkLessName.'_stream_'.$grassInputStreamThreshold;
            $LOC_vectStreamName 		=$LOC_rastMNTSinkLessName.'_stream_'.$grassInputStreamThreshold;
            $LOC_rastDirectionName 		=$LOC_rastMNTSinkLessName.'_direction';
            $LOC_rastStrahlerName 		=$LOC_rastMNTSinkLessName.'_strahler_'.$grassInputStreamThreshold;
		

                // Comblement des dépressions
            $result .= $app['twig']->render('/app/process-res3.twig',[ 'res' => execute($app, 'grass_r_fill_dir',    
                            [$grassUser, $grassDataDir,$grassLocation,$grassMapset,$LOC_rastMNTName, $LOC_rastMNTSinkLessName, $LOC_rastDirectionName])]);
            
                // Création du raster d'accumulation            
            $result .= $app['twig']->render('/app/process-res3.twig',[ 'res' => execute($app, 'grass_r_watershed',
                            [$grassUser, $grassMemoryMB,$grassDataDir,$grassLocation,$grassMapset, $LOC_rastMNTSinkLessName, $LOC_rastAccumulationName])]);
            
                // Construction du reseau théorique d'écoulement            
            $result .= $app['twig']->render('/app/process-res3.twig',[ 'res' => execute($app, 'grass_r_stream_extract',
                            [$grassUser, $grassMemoryMB,$grassDataDir,$grassLocation,$grassMapset,$LOC_rastMNTSinkLessName, $LOC_rastAccumulationName,$grassInputStreamThreshold, $LOC_rastStreamName, $LOC_rastStreamName, $LOC_rastDirectionName])]);
                
                // Ordination du réseau hydrologique selon Strahler            
            $result .= $app['twig']->render('/app/process-res3.twig',[ 'res' => execute($app, 'grass_r_stream_order',
                            [$grassUser, $grassMemoryMB,$grassDataDir,$grassLocation,$grassMapset,$LOC_rastStreamName,$LOC_rastDirectionName, $LOC_rastMNTSinkLessName, $LOC_rastAccumulationName, $LOC_rastStrahlerName, $LOC_vectStreamName])]);
                
            /*
             * export du $LOC_vectStreamName en vecteur pg
             */

            $result .= $app['twig']->render('/app/process-res3.twig',[ 'res' => execute($app, 'grass_v_out_ogr',
                        [$grassUser, $grassDataDir,$grassLocation,$grassMapset,$LOC_vectStreamName, 'line', 'ESRI_Shapefile',$LOC_vectStreamName])]);


            $result .= $app['twig']->render('/app/process-res3.twig',[ 'res' => execute($app, 'gdal_import_data_from_file_to_pg',
                        [$pghost, $pgbase, $pgPort, $pgUser, $pgPass, 'SHP', $grassDataDir.'/'.$grassLocation.'/'.$LOC_vectStreamName.'/'.$LOC_vectStreamName.'.shp', 'UTF-8', $destSchemaName, pg_object_define_valid_alias('t_'.$LOC_vectStreamName), 'LINESTRING', $epsg])]);

           echo $result ;
           die ;
}

Localisation des TBV pas à pas[modifier]

A. Coudart - MAJ 2020-01-20

Le référentiel hydrographique[modifier]

Le référentiel hydrographique doit être préparé en amont selon les règles suivantes :

• Table PG
  • gid serial
  • geom (linestring, 2154)

• Pour gommer les tronçons de cours d'eau que l'on souhaite supprimer du référentiel de travail (cours d'eau isolé, fossé connu...), on privilégie des masques polygones qui déduiront les tronçons s'y trouvant.

select  ce.geom
from r020_territoire_hydrographie.t_ce_bdtopo ce
where ce.gid not in (
select ce.gid
from r020_territoire_hydrographie.t_ce_bdtopo ce
, m060_milieux_continuite_tetes_de_bv.t_mask_ce_out mask
where st_intersects(mask.geom,ce.geom)

• En cas de réseau hydrographique à multiple exutoires, on privilégie la construction d'un exutoire fictif sur lequel chaque exutoire réel sera rattaché. Cet exutoire sera construit à l'amont avec une ramification suffisante pour que le premier exutoire réel se rattache à un rang de strahler supérieur ou égal à 3.

  Insert into m060_milieux_continuite_tetes_de_bv.ce_without_mask_v2
  (
  geom
  )
  select geom
  from m060_milieux_continuite_tetes_de_bv.t_ce_for_work

• La correction topologique est une étape indispensable se référer à la partie la traitant

• Pour créer les tables de masques et de cours d'eau fictif, lancer le script suivant dans l'éditeur SQL du SGBD

------------Création table des masques

create table a000_temp.t_ce_mask(
gid serial 
, geom geometry (polygon,2154)
);

alter table a000_temp.t_ce_mask
	add constraint t_ce_mask_pk primary key (gid);

CREATE INDEX t_ce_mask_geom_idx
  ON a000_temp.t_ce_mask
  USING gist
  (geom);

------------Création table des cours d'eau fictifs
  

create table a000_temp.t_ce_fictif(
gid serial 
, geom geometry (linestring,2154)
);

alter table a000_temp.t_ce_fictif
	add constraint t_ce_fictif_pk primary key (gid);

CREATE INDEX t_ce_fictif_geom_idx
  ON a000_temp.t_ce_fictif
  USING gist
  (geom);

Le modèle numérique de terrain[modifier]

Le modèle numériques de terrain doit être préparé en suivant les étapes suivantes :

• Assemblage des dalles du territoires
• Définition d'une transparence sur les pixels 0 (noirs) (ne pas créer de bande de transparence. Le MNT ne doit avoir qu'une seule bande)
• Découpage de l'assemblage avec un buffer de 1000 à 2000 m autour du territoire
• Dépôt du MNT sur le serveur

Le surcreusement du MNT[modifier]

Cette étape de préparation du MNT permet de rendre le modèle numérique conforme au chemin de l'eau réel sur le terrain.

Fonction processing = Processes/Tr_Hydro > Hydro_mnt_surcreusement_from_pg_hydro_vect

• Hydro_mnt_surcreusement_from_pg_hydro_vect
  • Utiliser le référentiel de cours d'eau préparé avant
  • Utiliser le MNT préparé avant

Constitution du réseau théorique d'écoulement[modifier]

Cette étape créer un réseau théorique d'écoulement à partir du MNT adapté au réel. Ce référentiel sera utilisé plus tard comme donnée de comblement en cas de "trous" dans le réseau.

Fonction processing = Processes/Tr_Hydro > hydro_strahler_rank_from_rast_mnt

• hydro_strahler_rank_from_rast_mnt
  • Utiliser le MNT surcreusé
  • La valeur d'accumulation de 500 permet d'avoir un réseau théorique important qui sera fort utile pour combler les trous éventuels dans le réseau.

1ère étape Délimitation[modifier]

Cette étape effectue la correction topologique et à la suppression des tronçons isolés

Fonction Processing = Processes/Tr_TBV/Delimitation > tbv_delim_step_010_reseau_hydro_correction_topologie_et_seg_isole 

• tbv_delim_step_010_reseau_hydro_correction_topologie_et_seg_isole
  • Utiliser le référentiel de cours d'eau préparé avant
  • Pour le choix de l'exutoire, pointé sur Qgis un point à l'aval du dernier exutoire (tolérance de 100 mètres, modulable)
  • Pour le nom de la table de destination, le suffixe "_corr" à ajouter au nom de la table d'origine est un bon compromis pour s'y retrouver
  • Ne pas oublier le table de réseau théorique calculée précédemment.

2ème étape Délimitation[modifier]

TODO AR : voir à utiliser postgis topology sur réseau hydro car semble très bien détecter aussi les anastomoses cf. réflexion topology.

Cette étape crée un graph à partir du référentiel hydro et corrige les anastomoses

Fonction Processing = Processes/Tr_TBV/Delimitation > tbv_delim_step_020_reseau_hydro_creation_de_graph_avec_correction_anastomoses  

• tbv_delim_step_020_reseau_hydro_creation_de_graph_avec_correction_anastomoses
  • Utiliser la table des cours d'eau avec le suffixe "_corr"

3ème étape Délimitation[modifier]

Cette étape calcule les rang de strahler

Fonction Processing = Processes/Tr_TBV/Delimitation > tbv_delim_step_030_reseau_hydro_calcul_rang_de_strahler_et_calc_rang_avant_et_apres_exutoire

• tbv_delim_step_030_reseau_hydro_calcul_rang_de_strahler_et_calc_rang_avant_et_apres_exutoire
  • Utiliser la table avec les suffixe "_corr_graph"

Attention !!! Des erreurs peuvent faire échouer le traitement !

ERROR:  more than one row returned by a subquery used as an expression
********** Erreur **********

****** Exemple ****** 
Si doublons faire le test suivant :
with lst as (
select t_up.gid_graph, t_up.strahler_rank
,t_in.gid_graph as dest_gid_graph
-- (select gid_graph from m060_milieux_continuite_tetes_de_bv.t_ce_ref_bdtopo_corr_graph t_in where
-- st_intersects(t_up.geom_graph_ori,t_in.geom_graph_ori) AND st_intersects(st_endpoint(t_up.geom_graph_ori),st_startpoint(t_in.geom_graph_ori))) as t2
from m060_milieux_continuite_tetes_de_bv.t_ce_ref_bdtopo_corr_graph t_up
left join m060_milieux_continuite_tetes_de_bv.t_ce_ref_bdtopo_corr_graph t_in on (st_intersects(t_up.geom_graph_ori,t_in.geom_graph_ori) AND st_intersects(st_endpoint(t_up.geom_graph_ori),st_startpoint(t_in.geom_graph_ori)))
)

select gid_graph, count(dest_gid_graph), string_agg(dest_gid_graph::text,','::text),
'Attention le gid_grah : ' || gid_graph || ' à ' ||count(dest_gid_graph) ||' gid_graph de destination (' || string_agg(dest_gid_graph::text,','::text)
|| ') ! Il faut résoudre manuellement.'
from lst
group by gid_graph
having count(dest_gid_graph)>1

4ème étape Délimitation[modifier]

Cette étape calcule l'ensemble des micros BV pour chaque tronçon de cours d'eau.

Fonction Processing = Processes/Tr_TBV/Delimitation > tbv_delim_step_040_reseau_hydro_delimitation_des_micro_bv

• tbv_delim_step_040_reseau_hydro_delimitation_des_micro_bv
  • Utiliser la table avec les suffixe "_corr_graph_grass_geom"
  • La règle de nommage "t_micro_bv_v01" est un bon compromis

5ème et dernière étape Délimitation[modifier]

Cette étape sélectionne les micros BV inférieurs ou égal à 2 et les assemblent pour créer la couche TBV définitive.

Fonction Processing = Processes/Tr_TBV/Delimitation > tbv_delim_step_050_reseau_hydro_extraction_des_tetes_de_bv

• tbv_delim_step_050_reseau_hydro_extraction_des_tetes_de_bv
  • Utiliser la table des cours d'eau avec les suffixe "_corr_graph"
  • Utiliser la table des micros BV précédemment crée
  • La règle de nommage "t_tbv_without_carac_v01" est un bon compromis

AntoineR : Détermination des Têtes de bassins versants[modifier]

Introduction[modifier]

TODO AR : revoir intro : ajout étapes de détermination. cours d'eau DDTM + IGN ->assemblage -> correction -> prolongation des segments isolés -> graph -> strahler -> micro bv -> têtes de bv

Selon le SAGE les têtes de bassins versants correspondent aux surfaces de bassins versants des cours d'eau de rang de Strahler <=2. Le référentiel hydrographique retenu pour la localisation des têtes de bassins versants et un référentiel composite constitué des données DDMT (cartographie des cours d'eau) préférentiellement et des données BDTOPO IGN à défaut (cours d'eau permanents + intermittents). En Aout 2018 : 44 et 49 inventaires DDTM réalisés, 85 trop partiel et 79 non réalisé.

Les couches sur DDTM ne sont pas jointives entre département (ou BDTOPO IGN) et possèdent des chevauchements ou incohérences en limites de départements. De plus, à l'heure actuelle, il n'existe pas d'identifiant unique permettant d'identifier l'origine d'un tracé DDTM par rapport au référentiel BDTOPO. Il faudra donc réaliser différentes opérations pour obtenir un référentiel hydrographique topologique pouvant être strahlerisé afin de déterminer les surfaces de têtes de bassins versants.

Grandes étapes :

1. Assemblage des référentiels de cours d'eau dans une couche composite
2. Correction de la couche composite
3. Délimitation des micro-bassins versants
4. Extraction/fusion des micro-bassins versant pour obtenir les géométries des têtes de bassins versants

Automatisation pour EPTBSN (Processing) : TODO AC : réaliser la chaine de traitement qui fait tout.

Automatisation pour autre structure : Dérouler pas à pas les étapes suivantes.

Assemblage des couches DDTM et IGN[modifier]

MAJ AR (2018_09_12): Les données DDTM 85 sont intégrées maintenant.

TODO AR : Il faut refaire les calculs avec ces nouvelles données (V03).

L'assemblage s'effectue en 2 étapes:

• Constitution d'une couche DDTM (r020_territoire_hydrographie.t_cours_eau_police_de_l_eau_ddtm_v02) par insertion sélective(et découpage) par département des données DDTM.
  • En Aout 2018 :
    • Buffer - 10m = coverage49 (sur le département) pour l'import des données DDTM 49 (linéaire de cours d'eau retenus) avec découpage des géométries pour ne conserver que celles qui sont dans le buffer.
    • Buffer + 10m =coverage44 (sur le département) pour l'import des données DDTM 44 (linéaire de cours d'eau retenus) avec découpage des géométries pour ne conserver que celles qui sont dans le buffer.
    • 85 et 79 pas d'import, données respectivement trop partielles et inexistantes.

• Ajout des données BDTOPO IGN (cours d'eau intermittents et permanents) pour créer la couche DDTM + IGN :"r020_territoire_hydrographie"."t_cours_ref_ddtm_et_ign"
  • En Aout 2018 :
    • Import des données BDTOPO IGN hors du coverageDDTM (= coverage49+coverage44) avec découpage des géométries pour ne conserver que celles qui sont en dehors du coverageDDTM.

Automatisation pour EPTBSN (Processing) : TODO AC : Scripte de traitement qui fait l'assemblage (puis indiquer ici le lien html).

Automatisation pour autre structure : Emplacement du scripte SQL :

• O:\TRAVAUX\ariviere\2017_RISQUES_EROSIFS_INONDATIONS_TRANSFERTS\70_Work_EPTBSN\11_COURS_EAU_THEORIQUES\Anastomoses\Creation_couche_ddtm_et_ddtm_ign.sql
• O:\TRAVAUX\ariviere\2017_RISQUES_EROSIFS_INONDATIONS_TRANSFERTS\70_Work_EPTBSN\11_COURS_EAU_THEORIQUES\Anastomoses\Creation_couche_ddtm_et_ddtm_ign_03.sql (nouvelle version)

IMPORTANT : : Pour un calcul de Strahler, cette couche n'est pas correcte :

• sur le plan topologique (segments qui ne se touchent pas)
• au niveau de la cohérence hydrologique (cours d'eau isolés)
• pour le calcul du Strahler (orientation non rigoureuse des segments et anastomoses)

Il va falloir la corriger/modifier pour opérer le calcul de Strahler

Correction de la couche composite d'assemblage DDTM et IGN pour créer une table graph[modifier]

TODO AR : Il faut refaire les calculs avec ces nouvelles données (V03).

Traitement préparatoire : création d'une couche d’écoulements théoriques forcée sur le référentiel[modifier]

Calcul des écoulements théoriques sur la base d'un MNT surcreusé (50 m) par le référentiel de cours eau (couche composite DDTM + BDTO IGN, cad r020_territoire_hydrographie.t_cours_ref_ddtm_et_ign). Cette donnée va permettre d'obtenir un tracé d'écoulement pour les cours d'eau isolés. La table ainsi créée (t_mnt_sc_ddtm_sk_stream_500) sera utilisée pour la correction de la couche composite (prolongation des segments isolés).

Paramètres :

• MNT source : MNT RGE IGN.
• Couche vectorielle hydro qui servira à surcreuser : r020_territoire_hydrographie.t_cours_ref_ddtm_et_ign
• surcreusement : 50m
• Threshold : 500

Automatisation pour EPTBSN (Processing) :

• http://processing.observatoire.sevre-nantaise.com/process?file=/var/www/observ/processing/src/App/Processes/Tr_Hydro/hydro_mnt_surcreusement_from_pg_hydro_vect.php&from_last_run=2726
• http://processing.observatoire.sevre-nantaise.com/process?file=/var/www/observ/processing/src/App/Processes/Tr_Hydro/hydro_strahler_rank_from_rast_mnt.php&from_last_run=2734

Automatisation pour autre structure : TODO AC lien vers section du wiki correspondantes :

• grass surcreusement
• réseau théorique avec stream_order

Résultats :

• table geographique avec différentes informations de hiérarchisation (Strahler, shreve (qui sera utilisé pour la prolongation des cours d'eau isolés), ... )  : m060_milieux_continuite_tetes_de_bv.t_mnt_sc_ddtm_sk_stream_500

Correction de la couche composite DDTM IGN (création d'un réseau graphé) à l'aide d'une couche d’écoulements théoriques forcée sur ce référentiel[modifier]

ATTENTION Opération manuelle à réaliser en amont : suppression d'un cours d'eau entre 2 plans d'eau 422385.3,6624714.4 dans la table r020_territoire_hydrographie.t_cours_ref_ddtm_et_ign.

mémo AR : correction + prolongation seg isolés + graph + strahlerisation + création bv dans grass

Etapes du traitement :

1. Correction de la topologie et des ségments isolés
   • fonction pg :xx_99_utils.f_reseau_hydro_correction_topologie_et_seg_isoles_v01_04).
   • Cette fonction a besoin de connaitre l’extrémité aval du réseau.
   1. Dans un premier temps le réseau va être remonté à partir du point aval et tous les segments vont être connecté à ce réseau s'ils sont à une distance <= à la valeur _buffer_snap (défaut 15m).
   2. [optionnel] : Les segments isolés restants vont être rattachés au réseau principal grâce aux données contenues dans le couche d'écoulements théorique basée sur un MNT surcreusé (cf. partie précédente )
   • -> résultat :
     • t_cours_ref_ddtm_et_ign.m060_milieux_continuite_tetes_de_bv.t_cours_ref_ddtm_et_ign_corr (cours d'eau corrigés, ne conserve pas les données de la table corrigée (perte des id d'objets...). La colonne correction_result permet d'identifier les tracés qui restent isolés après le traitement (correction_result like 'ISOLE')
     • t_cours_ref_ddtm_et_ign.m060_milieux_continuite_tetes_de_bv.t_cours_ref_ddtm_et_ign_corr_detail (information sur les opérations réalisées)
2. Création d'un réseau de graph avec suppression des anastomoses :
   • fonction pg : xx_99_utils.f_reseau_correction_anastomoses_v04_01
   • -> résultat :
     • table identique à la table en entrée avec ajout colonne graph : m060_milieux_continuite_tetes_de_bv.t_cours_ref_ddtm_et_ign_corr_graph
     • table avec géométries d'anastomoses détectées (utile pour visualiser les anastomoses détectées et pour le traitement ultérieur ): m060_milieux_continuite_tetes_de_bv.t_cours_ref_ddtm_et_ign_corr_anast
3. Calcul des rang de strahler [premier passage] (fonction pg : xx_99_utils.f_reseau_hydro_strahler_rank_v01_04 ) . Nécessite de connaitre le point bas du réseau. Cette function va orienter tous les segments et calculer leur rang de strahler.
   • -> résultat :
     • m060_milieux_continuite_tetes_de_bv.t_cours_ref_ddtm_et_ign_corr_graph avec ajout de colonnes : geom_graph_ori (orienté) et strahler_rank
4. Ajout des rangs de strahler avant confluence (strahler_before_outlet) et après confluence (strahler_after_outlet) et gid du tronçon suivant (graph_gid_next) [premier passage] à chaque tronçon de Cours d'eau dans la Table des Cours d'eau graph et strahlerisé (t_cours_ref_ddtm_et_ign_corr_graph)
   • -> résultat :
     • m060_milieux_continuite_tetes_de_bv.t_cours_ref_ddtm_et_ign_corr_graph avec ajout de colonnes : strahler_before_outlet, strahler_after_outlet graph_gid_next.
5. suppression des mini troncons provenant de problèmes de digit ou lors des opérations de correction correction
   • rang de strahler =1 qui génère un strahler >= 2 et dont la longueur est inférieure une distance légèrement supérieure à la diagonale de 2 pixels de MNT (14.5 ). pour mémoire diagonale de 2 pixel de MNT (5mx5m) = 14.14 m
   • -> résultat :
     • table graph nettoyée des petites erreurs : m060_milieux_continuite_tetes_de_bv.t_cours_ref_ddtm_et_ign_corr_graph
6. Calcul des rang de strahler [deuxième passage]
   • toujours la même fonction : xx_99_utils.f_reseau_hydro_strahler_rank_v01_04
   • -> résultat :
     • table graph strahlerisée : m060_milieux_continuite_tetes_de_bv.t_cours_ref_ddtm_et_ign_corr_graph
7. Ajout des rangs de strahler avant confluence (strahler_before_outlet) et après confluence (strahler_after_outlet) et gid du tronçon suivant (graph_gid_next) [deuxième passage] à chaque tronçon de Cours d'eau dans la Table des Cours d'eau graph et strahlerisé (t_cours_ref_ddtm_et_ign_corr_graph)
   • -> résultat :
     • m060_milieux_continuite_tetes_de_bv.t_cours_ref_ddtm_et_ign_corr_graph avec ajout de colonnes : strahler_before_outlet, strahler_after_outlet graph_gid_next.
8. Orientation de la géometrie d'origine (geom) en geom_ori sur la base l'orientation de la geom_graph_ori. Cette donnée servira lors de l'extraction des points haut/bas du chemin le plus long dans chaque bv.
   • -> résultat :
     • m060_milieux_continuite_tetes_de_bv.t_cours_ref_ddtm_et_ign_corr_graph avec ajout de colonne : geom_ori

Paramètres :

• buffer_snap : 15
• ...

Automatisation pour EPTBSN (Processing) :

• TODO AC à intégrer dans un scripte php + lien vers le last run
  • http://processing.observatoire.sevre-nantaise.com/process?file=/var/www/observ/processing/src/App/Processes/Tr_TBV/Delimitation/tbv_delim_step_010_reseau_hydro_correction_topologie_et_seg_isole.php
  • http://processing.observatoire.sevre-nantaise.com/process?file=/var/www/observ/processing/src/App/Processes/Tr_TBV/Delimitation/tbv_delim_step_020_reseau_hydro_creation_de_graph_avec_correction_anastomoses.php
  • http://processing.observatoire.sevre-nantaise.com/process?file=/var/www/observ/processing/src/App/Processes/Tr_TBV/Delimitation/tbv_delim_step_030_reseau_hydro_calcul_rang_de_strahler_et_calc_rang_avant_et_apres_exutoire.php

Automatisation pour autre structure : Lien vers le fichier SQL

O:\TRAVAUX\ariviere\2017_RISQUES_EROSIFS_INONDATIONS_TRANSFERTS\70_Work_EPTBSN\11_COURS_EAU_THEORIQUES\Anastomoses\Correction_et_strahlerisation_CE_DDTM_TOPO_03.sql

Résultats :

• couche des cours d'eau corrigés (snap, prolongation des segments isolés, mais conserves des artefacts de digitaliséation et mauvaise orientation des cours d'eau) : t_cours_ref_ddtm_et_ign.m060_milieux_continuite_tetes_de_bv.t_cours_ref_ddtm_et_ign_corr
• couche informative des opérations réalisées en étapes de correction  : t_cours_ref_ddtm_et_ign.m060_milieux_continuite_tetes_de_bv.t_cours_ref_ddtm_et_ign_corr_detail
• couche des cours d'eau corrigée et graphés, strahlerisés, avec suppression des artefacts  :m060_milieux_continuite_tetes_de_bv.t_cours_ref_ddtm_et_ign_corr_graph
• couche de position des anastomoses :m060_milieux_continuite_tetes_de_bv.t_cours_ref_ddtm_et_ign_corr_anast

Délimitation des micro-bassins versants[modifier]

La couche de cours corrigée est utilisée pour calculer les micro-bassins versants.

Paramètres :

• correction_result not like 'ISOLE'
• Threshold : 4000 (m2)

Automatisation pour EPTBSN (Processing) : TODO AC scripte + lien vers last run:

mémo interne : fichier SQL + commandes GRASS à ne pas diffuser : O:\TRAVAUX\ariviere\2017_RISQUES_EROSIFS_INONDATIONS_TRANSFERTS\70_Work_EPTBSN\11_COURS_EAU_THEORIQUES\Anastomoses\Calcul_microBV_DDTM_IGN_v01.sql

1. Réalisé dans le traitement précédant !!! . SQL : Copie de la table de graph dans une table qui ne contiendra qu'une seule géométrie : geom (pour te traitement grass suivant) et qui ne contient pas les segments isolés ('ISOLE'), car ce sont des segments soit : hors du bassin versant de la Sèvre Nantaise, soit trop isolés pour être rattachés au réseau hydrographique principal.

-- Copie des enregistrements filtrés
DROP TABLE IF EXISTS   m060_milieux_continuite_tetes_de_bv.t_cours_ref_ddtm_et_ign_corr_graph_grass_geom ;
   CREATE TABLE m060_milieux_continuite_tetes_de_bv.t_cours_ref_ddtm_et_ign_corr_graph_grass_geom as select 
   * from m060_milieux_continuite_tetes_de_bv.t_cours_ref_ddtm_et_ign_corr_graph
   where correction_result not like 'ISOLE'
   ;
-- Suppression des colonnes géométriques inutiles (notre version de Grass ne peut pas charger une table avec plusieurs colonnes géométriques)
ALTER TABLE m060_milieux_continuite_tetes_de_bv.t_cours_ref_ddtm_et_ign_corr_graph_grass_geom DROP COLUMN geom_graph;
ALTER TABLE m060_milieux_continuite_tetes_de_bv.t_cours_ref_ddtm_et_ign_corr_graph_grass_geom DROP COLUMN geom_graph_ori;
ALTER TABLE m060_milieux_continuite_tetes_de_bv.t_cours_ref_ddtm_et_ign_corr_graph_grass_geom DROP COLUMN geom_ori;

2. Surcreusement du MNT (utilisation de la colonne geom de la table précédente ) : http://processing.observatoire.sevre-nantaise.com/process?file=/var/www/observ/processing/src/App/Processes/Tr_Hydro/hydro_mnt_surcreusement_from_pg_hydro_vect.php&from_last_run=3242

Launched the 2018-08-08 @ 12:26:39 by ariviere

3. Calcul du rast direction  : à partir du t_cours_ref_ddtm_et_ign_corr_graph_grass_geom réalise le sinklesss et créé le raster de direction ( TODO : méthode lourde pour direction : il faudrait s'arrêter à r.stream.extract) := mnt_sc_ddtm_cor_sk_direction : http://processing.observatoire.sevre-nantaise.com/process?file=/var/www/observ/processing/src/App/Processes/Tr_Hydro/hydro_strahler_rank_from_rast_mnt.php&from_last_run=3250

Launched the 2018-08-08 @ 12:37:14 by ariviere

4. import des données vect pg strahlerisées => rast grass (note l'import a déjà été réalisé en 2. donc voir à simplifier cette étape si possible, mais attention à générer le bon raster pour la création des bvs (pas de fusion sur le même id ..).

Directement dans GRASS (location : bvsn_surcreuse_02_cor_v02) :

v.in.ogr input='PG:dbname=observatoire host=sevre-nantaise.com port=5432 user=XXXX password=XXXXX layer=m060_milieux_continuite_tetes_de_bv.t_cours_ref_ddtm_et_ign_corr_graph_grass_geom output=stream_ddtm_ign type=line geometry=geom --overwrite

5. v to rast :

Directement dans GRASS (location : bvsn_surcreuse_02_cor_v02) :

v.to.rast --overwrite input=stream_ddtm_ign@PERMANENT type=line output=stream_ddtm_ign use=cat --overwrite

6. Stream basin :

Directement dans GRASS (location : bvsn_surcreuse_02_cor_v02) :

r.stream.basins dir=mnt_sc_ddtm_cor_sk_direction@PERMANENT stream=stream_ddtm_ign@PERMANENT basins=stream_ddtm_ign_bv memory=10000 --overwrite

7. r to vect :

Directement dans GRASS (location : bvsn_surcreuse_02_cor_v02) :

r.to.vect input=stream_ddtm_ign_bv@PERMANENT output=stream_ddtm_ign_bv type=area --overwrite

8. v out ogr (=>shp) :

Directement dans GRASS (location : bvsn_surcreuse_02_cor_v02) :

v.out.ogr input=stream_ddtm_ign_bv@PERMANENT output=/var/www/observ/out_processing/grass_data/bvsn_surcreuse_02_cor_v02/t_micro_bv_ign_ddtm.shp format=ESRI_Shapefile --overwrite

9. Import du Shp dans PG:

TODO AC : utiliser la fonction php import shp dans pg? TODO AR faire le trie entre gdal_import_data_from_file_to_pg.php et pg_import_data_from_file_to_pg.php

/*
Commentaire sur la table : 
*/
COMMENT ON TABLE r020_territoire_hydrographie.t_micro_bv_ddtm_ign
  IS '{"metadata"={
        "intitule": "Couche de délimitation des micros Bassins versants",
        "producteur": "EPTBSN",
        "descripteur": "AR",
        "description": "Couche de délimitation des micros BV (surface min 4000 (m2)) sur la base des données cours d''eau corrigée (topologique) des données DDTM et IGN  ",
        "tags": [
            "Bassin versant", "micro bassin versant"],
        "import_date": "2018-08-20 18:12:56.301708+02",
        "import_type": ["FONCTION_PHP"],
        "intitule_version": "2017_01",
        "procedure_import": "http://wiki.sevre-nantaise.com/index.php/Localisation_et_caract%C3%A9risation_des_t%C3%AAtes_de_bassin_versant#D.C3.A9limitation_des_micro-bassins_versantsl",
        "procedure_collecte": "http://wiki.sevre-nantaise.com/index.php/Localisation_et_caract%C3%A9risation_des_t%C3%AAtes_de_bassin_versant#D.C3.A9limitation_des_micro-bassins_versantsl",
        "intitule_couverture": "2017",
        "donnees_sources_type": "DIV",
        "donnees_sources_encoding": "UTF8",  
        "donnees_sources_localisation_loc": "",
        "donnees_sources_localisation_web": "",
        "donnees_sources_metadonnees_localisation_loc": "",
        "donnees_sources_metadonnees_localisation_web": "",
        "json_metedata_update_date": "2018-09-11 10:42:56.301708+02"
     "},

    "old_comment={  
        "value": "" 
    }
}';

Etapes de traitement pas à pas pour autre Structure :

TODO AC : voir les étapes pour processing à copier ici = compléments à ajouter

Résultats :

• table pg des micro bv  :m60_milieux_continuite_tetes_de_bv.t_micro_bv_ign_ddtm

Création de la couche des têtes de basin versants[modifier]

Les micro BV sont fusionnés jusqu'au rang de strahler =2.

TODO AR : Finir de commenter le SQL.

Paramètres :

• NA

Automatisation pour EPTBSN (Processing) :

• TODO AC : SQL à mettre en php + lien vers last run

Automatisation pour autre structure : Lien vers le fichier SQL:

• O:\PROJETS\2018_Tetes_BV\2018_tetes_bv_acoudart\01_donnees\01_sql_query\create_TBV_Final_without_carac_v03.sql

Résultats :

• Table des têtes de bassins versants : m060_milieux_continuite_tetes_de_bv.t_tbv_01_without_caract

Caractérisation des têtes de bassin versant[modifier]

La caractérisation des têtes de bassin versant consiste à évaluer l'état et les pressions qui leur sont associés. Cette étape comprend 3 grands volets qui seront étudier en 2 temps selon 3 enjeux majeurs.

Les 3 grands volets :

• La caractérisation physique et morphologique du bassin versant : les critères définisse la structure physique pure et via des formules mathématiques éprouvées permettent d'évaluer les capacités de réaction et de résilience du BV.
• La caractérisation du lit mineur et de la bande riveraine : les critères définissent l'état et les pressions qui s'exercent à l'échelle du cours d'eau. Les résultats permettront de cibler certains type de travaux de restauration ou mesure de protection spécifique au lit mineur ou au berges...
• La caractérisation du bassin versant : les critères sont peu ou prou les mêmes que pour le lit mineur et la bande riveraine, avec des éléments supplémentaires comme la densité de prélèvement et la dynamique d'urbanisation. Les résultats permettront également de cibler directement des typologies d'actions à l'échelle du bassin versant, inscription dans les PLU(i), sensibilisation des élus, accompagnement agricole...

L'analyse se fait en 2 temps :

• Les résultats bruts : Ils sont stockés et conservés dans les tables des TBV mais il ne seront pas analysés sous cette forme
• Les résultats discrétisés : Ils sont stockés, conservés et analysés et seront le point d'appui des pondérations appliqués dans la phase de priorisation des têtes de bassins versant à enjeux.

Chaque critères peu ou non intervenir dans la grille de lecture des enjeux, ces derniers sont le reflet des différents volets composant l'objectif final de reconquête de la qualité de l'eau portée par la directive cadre sur l'eau (DCE) du 23 octobre 2000 (directive 2000/60) :

• Enjeux qualité de l'eau (QE)
• Enjeux quantité de l'eau (QTE)
• Enjeux qualité des milieux (QM)

La liste des données indicateurs est disponible en interne dans ce tableau Calc : "O:\PROJETS\2018_Tetes_BV\2018_tetes_bv_acoudart\00_final_all_files\05_analyse_données\donnees_caracterisation.ods".

Caractérisation physique et morphologique[modifier]

La caractérisation physique et morphologique de la tête de bassin versant consiste à éprouvé sa capacité de réaction et donc de résilience. La dynamique hydraulique d'un cours d'eau pentu et réactif étant pour importante, il aura plus de faciliter à s'auto restaurer qu'un cours d'eau de plaine à faible débit.

Grille d'indicateurs[modifier]

Attention, la prise en compte des indicateurs de manière positive ou négative est à revoir !

Préparation des données MNT rasters nécessaires aux calculs[modifier]

• calcul du raster de slope (pente en pourcentage = percent) avec grass; Utiliser en entrée le MNT RGE de l'iIGN en brut : http://processing.observatoire.sevre-nantaise.com/process?file=/var/www/observ/processing/src/App/Processes/Generics/grass/grass_r_slope_aspect.php
• création d'un raster 3 bandes : http://processing.observatoire.sevre-nantaise.com/process?file=/var/www/observ/processes_master/Processes/Generics/grass/grass_i_group.php

i.group group=topo input=bvsn_brut_MNT,bvsn_brut_MNT_aspect,bvsn_brut_MNT_slope

• export tif : r.out.gdal input=topo output=/var/www/observ/donnees-brutes/IGN_MNT_RGE/RGE_Alti_5mx5m_60cm_2015_aspect_slope_percent_grass.tif
• export to pg rast (attention il faut encore saisir le mdp d'accès à pg après avoir copié la ligne de commande dans un shell) :

http://processing.observatoire.sevre-nantaise.com/process?file=/var/www/observ/processing/src/App/Processes/Generics/ImportExport/pg_import_raster_from_file_to_pg_rast.php

• conversion du rast en table polygon :
  

/*
 Création de la table à partir du raster pg
 MAJ 2020/10 v02:
- optimisation de la durée du traitement (<> 1 heure pour vie et jaunay)

 */
set work_mem ='32GB';
SET client_min_messages = ERROR;

/*
alti temp
*/
drop table if exists r020_territoire_physique.t_mnt_rge_r1_alti_temp  ;

Create table r020_territoire_physique.t_mnt_rge_r1_alti_temp  as  (
 with 
 gv as (
   SELECT r.rid, ST_PixelAsPolygons(r.rast,1) gv  
   FROM r020_territoire_physique.r_mnt_rge_alti_aspect_slope_percent_grass r
 )
 select 
    rid
  , (gv).x
  , (gv).y 
  , (gv).val::real as alti

  , (gv).geom::geometry(polygon,2154) as geom
  from gv 
 )
 ; 

/*
aspect temp
*/
  
set work_mem ='32GB';
drop table if exists r020_territoire_physique.t_mnt_rge_r2_aspect_temp;
Create table r020_territoire_physique.t_mnt_rge_r2_aspect_temp  as  (
 with gv as (
   SELECT r.rid, ST_PixelAsPolygons(r.rast,2) gv  
   FROM r020_territoire_physique.r_mnt_rge_alti_aspect_slope_percent_grass r
  
  )
 select rid,(gv).val::real as aspect, (gv).x, (gv).y 
  from gv 
 )
 ; 

set work_mem ='32GB';
  CREATE INDEX idx_t_mnt_rge_r2_aspect_temp
  ON r020_territoire_physique.t_mnt_rge_r2_aspect_temp  
  USING btree
  (rid, x, y);

/*
slope temp
*/
  
set work_mem ='32GB';
DROP  table IF EXISTS r020_territoire_physique.t_mnt_rge_r3_slope_temp ;
Create table r020_territoire_physique.t_mnt_rge_r3_slope_temp  as  (
 with gv as (
   SELECT r.rid, ST_PixelAsPolygons(r.rast,3) gv  
   FROM r020_territoire_physique.r_mnt_rge_alti_aspect_slope_percent_grass r
  
  )
 select rid,(gv).val::real as slope_percent, (gv).x, (gv).y  from gv 
 )
 ; 

set work_mem ='32GB';
  CREATE INDEX idx_t_mnt_rge_r3_slope_temp
  ON r020_territoire_physique.t_mnt_rge_r3_slope_temp  
  USING btree
  (rid, x, y);

/*
Assemblage :
*/
set work_mem ='32GB';
drop table if exists r020_territoire_physique.t_mnt_rge_alti_aspect_slope_percent_grass  ;

Create table r020_territoire_physique.t_mnt_rge_alti_aspect_slope_percent_grass  as  (
   select row_number() over() ::bigint as gid
   , r1.alti::real 
   , r2.aspect::real
   , r3.slope_percent::real
   , r1.geom::geometry(polygon,2154) as geom
   FROM r020_territoire_physique.t_mnt_rge_r1_alti_temp r1
   LEFT JOIN r020_territoire_physique.t_mnt_rge_r2_aspect_temp r2 on (r1.rid=r2.rid and r1.x=r2.x AND r1.y=r2.y)
   LEFT JOIN r020_territoire_physique.t_mnt_rge_r3_slope_temp r3 on (r1.rid=r3.rid and r1.x=r3.x AND r1.y=r3.y)
);

/*
AJOUT des commentaires
*/

COMMENT ON TABLE r020_territoire_physique.t_mnt_rge_alti_aspect_slope_percent_grass
  IS 'Données du RGE alti de l''IGN (5mx5m)
date de création : 10/2020

Méthode de calcul et de création de la couche :

http://wiki.sevre-nantaise.com/index.php/Localisation_et_caract%C3%A9risation_des_t%C3%AAtes_de_bassin_versant#Pr.C3.A9paration_des_donn.C3.A9es_rasters_n.C3.A9cessaires_aux_calculs

Référence grass :
https://grass.osgeo.org/grass70/manuals/r.slope.aspect.html

';
COMMENT ON COLUMN r020_territoire_physique.t_mnt_rge_alti_aspect_slope_percent_grass.alti IS 'altitude en m';
COMMENT ON COLUMN r020_territoire_physique.t_mnt_rge_alti_aspect_slope_percent_grass.slope_percent IS 'pente en poucentage (calcul grass https://grass.osgeo.org/grass70/manuals/r.slope.aspect.html)';
COMMENT ON COLUMN r020_territoire_physique.t_mnt_rge_alti_aspect_slope_percent_grass.aspect IS 'Direction de pente  = orientation (calcul grass : https://grass.osgeo.org/grass70/manuals/r.slope.aspect.html)

Description des valeurs :
The aspect output raster map indicates the direction that slopes are facing. The aspect categories represent the number degrees of east. Category and color table files are also generated for the aspect raster map. The aspect categories represent the number degrees of east and they increase counterclockwise: 90 degrees is North, 180 is West, 270 is South 360 is East.
  90 : Nord
180 : Ouest
270 : Sud
360 : Est';



/*
AJOUT pk et index geom
*/
set work_mem ='32GB';


ALTER TABLE r020_territoire_physique.t_mnt_rge_alti_aspect_slope_percent_grass  
  ADD CONSTRAINT t_mnt_alti_aspect_slope_percent_grass_pkey PRIMARY KEY(gid);

CREATE INDEX sidx_t_mnt_alti_aspect_slope_percent_grass_geom
  ON r020_territoire_physique.t_mnt_rge_alti_aspect_slope_percent_grass  
  USING gist
  (geom);


/*
Supression des tables temporaires
*/
drop table if exists r020_territoire_physique.t_mnt_rge_r1_alti_temp  ;
drop table if exists r020_territoire_physique.t_mnt_rge_r2_aspect_temp;
DROP table IF EXISTS r020_territoire_physique.t_mnt_rge_r3_slope_temp ;

Surface de BV[modifier]

Calcul visant à connaitre la taille du BV

select a.*, st_area(a.wkb_geometry)/10000 as surface_bv
from _agent_acoudart._t_v_test_sevre_amont_mnt_sinkless_basins_2et1 a

Périmètre de BV[modifier]

Calcul nécessaire pour les calculs suivant.

select a.*, st_perimeter(a.wkb_geometry) as perimetre_bv
from _agent_acoudart._t_v_test_sevre_amont_mnt_sinkless_basins_2et1 a

Chemin le plus long (écoulement et non cours d'eau)[modifier]

TODO AR : A implémenter en chaine. + adapter tbv_passini.php + tbv_pente_moy_longest_path.php

• o GRASS pour le calcul de stream sur la base du MNT surcreusé avec Threshold = 1 : r.fill.dir + r.watershed + r.stream.extract + r.stream.order , la sortie vectorielle de r.stream.order contient une colonne avec cum_length = longeur cumulée en m.

http://processing.observatoire.sevre-nantaise.com/process?file=/var/www/observ/processing/src/App/Processes/Tr_Hydro/hydro_strahler_rank_from_rast_mnt.php&from_last_run=3659

Launched the 2018-09-19 @ 14:18:55 by ariviere

table en sortie  : m060_milieux_continuite_tetes_de_bv.t_mnt_sc_ddtm_cor_sk_stream_1

• SQL parcours de Graph à récupérer ici + à adapter :

"O:\PROJETS\2018_Tetes_BV\2018_tetes_bv_acoudart\00_final_all_files\04_scripts_sql\Annexe_6 create_TBV_Final_without_carac_v03.sql"

Pente moyenne du chemin le plus long[modifier]

• cf : http://wiki.sevre-nantaise.com/index.php?title=Localisation_et_caract%C3%A9risation_des_t%C3%AAtes_de_bassin_versant#Chemin_le_plus_long_.28.C3.A9coulement_et_non_cours_d.27eau.29

Temps de concentration des pluies (Formule de Passini)[modifier]

• cf : http://wiki.sevre-nantaise.com/index.php?title=Localisation_et_caract%C3%A9risation_des_t%C3%AAtes_de_bassin_versant#Chemin_le_plus_long_.28.C3.A9coulement_et_non_cours_d.27eau.29

Indice de Gravélius (ou indice de compacité)[modifier]

L'indice de Gravélius est un indice de compacité utilisé dans les modèles hydrauliques, il permet de connaitre la forme générale d'un bassin versant. Un indice de 1 correspond à un bassin versant rond peu réactif au précipitation car très probablement composé d'un réseau de cours d'eau dense, à l'inverse plus l'indice est élevé et plus la forme du bassin versant s'allonge, traduisant une plus grande réactivité et assimilé à un réseau hydraulique principal majoritaire alimenté par un réseau peu dense de cours d'eau perpendiculaire. Plus le Bassin versant est allongé plus le débit de pointe des crue sera important.

Drop table if exists _agent_acoudart._test_indice_gravelius_sangueze;
create table _agent_acoudart._test_indice_gravelius_sangueze as

with w_prepare as(

select a.ogc_fid as id
, a.wkb_geometry as geom
,st_area(a.wkb_geometry) as surface_bv
,st_perimeter(a.wkb_geometry) as perimetre_bv
from _agent_acoudart._t_test_sangueze_mnt_sinkless_basins_2et1 a
)

select id
, geom
, (perimetre_bv/(2*(|/pi()*surface_bv))) as indice_gravelius
from _agent_acoudart._t_test_sangueze_mnt_sinkless_basins_2et1 a, w_prepare
group by id, geom, indice_gravelius

Indice de densité de chevelu[modifier]

La densité de chevelu nous informe sur les capacité d'auto épuration du bassin versant. En effet plus le réseau est dense et plus les zones de connexion entre l'eau et la terre (zones hyporéïques) seront nombreuses, augmentant ainsi les zones de dénitrification et les capacités d'auto épuration du cours d'eau.

drop table if exists _agent_acoudart._test_densite_drainage_sangueze;
create table _agent_acoudart._test_densite_drainage_sangueze as
with w_prepare as 
(
select b.ogc_fid   as gid
, b.wkb_geometry as geom
, st_area(b.wkb_geometry)/10000. as area_ha
, st_multi(ST_union(st_collectionextract(xx_99_utils.st_intersection_btr(b.wkb_geometry, a.geom),2)))::geometry(multilinestring,2154) as geom_ce
from _agent_acoudart._t_test_sangueze_mnt_sinkless_basins_2et1 b
left join _agent_acoudart."_r_theo_sangueze_5ha_from_IGN_force" a on ( st_intersects(b.wkb_geometry,a.geom))
GROUP  by b.ogc_fid
, b.wkb_geometry
)
select 
gid
, geom
, area_ha
, st_length (geom_ce)  as long_ce_ml
, st_length (geom_ce)/area_ha as densite_chevelu

from w_prepare

Densité du réseau de points bas[modifier]

Au contraire de la densité du réseau de chevelu, cet indice porte sur le réseau d'écoulement théorique créé grâce au MNT. Cet élément apporte une information en terme de capacité de ruissellement du bassin versant.

Pente moyenne du bassin versant[modifier]

La pente moyenne du bassin versant nous donne un indice en terme de relief, plus la valeur sera élevée et plus le bassin versant sera accidenté, augmentant ainsi sa réactivité.

with bv as (
select bv.*
, ST_buffer(bv.wkb_geometry,-1) geom_b_neg
from
_agent_acoudart._t_test_sangueze_mnt_sinkless_basins_2et1 bv
)
,
a as (
SELECT
bv.ogc_fid
, bv.wkb_geometry
,r.slope_percent  as rval
FROM
bv
LEFT JOIN r020_territoire_physique.t_mnt_rge_slope_grass_percent r
on (st_intersects(bv.geom_b_neg, r.geom))
--t_intersects(bv.wkb_geometry, r.rast, 1) AND
)

select

 a.ogc_fid
, count(rval)
, avg(rval)
, stddev(rval)
, min(rval)
, max(rval)

  from a
  group by a.ogc_fid

Position de la TBV dans le réseau hydrographique[modifier]

Cet indice nous apporte une information en terme de localisation de la TBV, si sa position est apicale (connectée à extrémité amont d'un rang 3) alors la bibliographie montre que son incidence sur le milieu sera plus importante qu'une position tributaire latéral (connectée à des rangs 3 et 4) ou tributaire aval (connectée à des rangs >5). Cela s'explique par les capacités de dilution plus importante de rangs supérieurs, ainsi une tête de bassin versant apicale dégradée aura un impacté négatif fort sur les milieux à l'aval contrairement à une tête de bassin versant tributaire aval qui verra son rôle minimisé.

Caractérisation du lit mineur et de la bande riveraine[modifier]

Grille d'indicateurs[modifier]

Caractérisation du bassin versant[modifier]

Grille d'indicateurs[modifier]

Indice de pression[modifier]

Les indices de pressions ci-dessous ont été validés par des comités d'experts sur des bases bibliographique et scientifique. Cette expertise comporte inévitablement des biais (non connaissance des parcelles en agriculture biologique, surfaces imperméabilisées comprenant les parcs et jardins des villes, rotation des cultures [Une année en fourrage, l'autre en culture], pas de différence entre les typologies de cultures...). Toutefois, les informations apportées sont jugées exploitables en l'état permettent de caractériser convenablement l'occupation du sol du territoire.

Tableau de correspondance indice de pression / occupation du sol[modifier]

Construction de l'occupation du sol de la tête de bassin versant[modifier]

L'occupation du sol est une donnée majeure en matière d'état écologique. Une surface artificialisée ne permettra plus aux milieux naturels d'assurer leur fonctionnalités originelles. Ainsi pour caractériser l'état des têtes de bassin versant il conviens d'en apprécier leur couverture. Pour cela Il nous a fallu croiser les données d'occupation urbaine, des tronçons routiers, des tronçons de voies ferrées, les zones de cultures, les zones de prairies permanentes, les surfaces boisées, les surfaces en eau...

Il existe plusieurs sources de données d'occupation du sol, la plus connues à l'échelon européen est Corine Land Cover. Cette donnée possède l'avantage d'être homogène et cohérente à cette échelle européenne mais son échelle d'utilisation est à 1/100 000e, ce qui lui confère peu d'intérêt en matière d'analyse locale.

La carte d'occupation du sol du pôle Theia propose une couverture à l'échelle nationale avec une résolution de 10 mètres pour les années 2016 et 2017. Ce pôle est une infrastructure de données satellitaire sous la gouvernance d'un comité directeur et d'un comité scientifique composés par les grandes institutions scientifique nationale, comme le CNES, l'IGN, l'INRA, le CNRS, l'IRSTEA, Météo-France, le CEREMA. Les données sont produites à partir d'images satellite du programme Copernicus Sentinel 2 (résolution à 10m) pour les dernières années (2016 / 2017) et Landsat 5 (résolution à 30m) pour les millésimes de 2009 - 2011 et 2014.
Les données seront plus précises qu'avec CLC, mais toujours relativement peu précise à grande échelle, par exemple pour l'étude d'impact sur une bande riveraine d'un cours cours sur 20 mètres de large, la précision n'est pas suffisante.

L'OCS GE, ou le référentiel à grande échelle concernant l'occupation du sol, est éditer par l'IGN et fournit une donnée précise, homogène sur le territoire national, en cohérence avec le RGE et millésimée. Toutefois, elle n'est pas encore disponible sur l'ensemble du territoire, en ce qui concerne le bassin versant de la Sèvre Nantaise, le département des Deux Sèvres manque à l'appel.
C'est la source de donnée la plus précise, bien que partielle.

Enfin il est également possible de construire son propre référentiel d'occupation du sol, c'est ce que nous avons entrepris.

Construction table : zones urbaines[modifier]

Pour construire la couche zone urbanisée, nous avons utilisées différentes données issues de la Bd Topo de l'IGN :

• Bâti industriel
• Bâti remarquable
• Bâti indifférencié
• Gare
• Surface route (parkings, places, les routes étaient décrites dans une autre source de donnée qui sera exploitée ultérieurement...)
  

La méthode utilisée est celle éditée par l'AGRESTE de Franche Comté et que l'on retrouve dans de nombreuses études sur les têtes de bassins versant (cf. SYLOA, Vilaine...). La logique suivit elle celle de la dilatation des éléments (buffer 50m), de leur agrégation puis de l'érosion de l'élément composite (buffer -50m).

Répertoire de travail: O:\TRAVAUX\ariviere\2018_ASSEMBLAGE_REFERENTIELS\2018_URBAIN

• Pour afficher le script SQL cliquer sur "Afficher"

-- 
-- Etape 1 : buffers de 50m sur les différentes couches de bati et union dans une table unique

drop table if exists _agent_acoudart.__bati_buff_bvsn_step01 ;
create table _agent_acoudart.__bati_buff_bvsn_step01 as 

with 
w_prepare as (
select bi.geom, st_multi(st_buffer(bi.geom,50,'join=mitre mitre_limit=2'))::geometry(multipolygon, 2154) as geom_b50
from xx_99_50_10_import_wfs.t_ign_bdtopo_bati_indifferencie bi

union all
select  bindus.geom, st_multi(st_buffer(bindus.geom,50,'join=mitre mitre_limit=2'))::geometry(multipolygon, 2154) as geom_b50
from xx_99_50_10_import_wfs.t_ign_bdtopo_bati_industriel bindus

union all
select br.geom, st_multi(st_buffer(br.geom,50,'join=mitre mitre_limit=2'))::geometry(multipolygon, 2154) as geom_b50
from xx_99_50_10_import_wfs.t_ign_bdtopo_bati_remarquable br

union all
select  g.geom, st_multi(st_buffer(g.geom,50,'join=mitre mitre_limit=2'))::geometry(multipolygon, 2154) as geom_b50
from xx_99_50_10_import_wfs.t_ign_bdtopo_gare g

union all
select  sur_rte.geom, st_multi(st_buffer(sur_rte.geom,50,'join=mitre mitre_limit=2'))::geometry(multipolygon, 2154) as geom_b50
from xx_99_50_10_import_wfs.t_ign_bdtopo_surface_route sur_rte

) 
select row_number()over()::bigint as gid,*
from w_prepare
;

DROP INDEX if exists _agent_acoudart.sidx___bati_buff_bvsn_geom_b50;

CREATE INDEX sidx___bati_buff_bvsn_geom_b50
  ON _agent_acoudart.__bati_buff_bvsn
  USING gist
  (geom_b50);


-- Etape 2 : fusion des buffers qui se touchent

drop table if exists  _agent_acoudart.__bati_buff_bvsn_step02 ;
create table _agent_acoudart.__bati_buff_bvsn_step02 as 
with
w_prepare as (
select (st_dump(st_collectionextract(st_union(geom_b50),3))).geom::geometry(polygon,2154) as geom
from _agent_acoudart.__bati_buff_bvsn_step01  bb
)
select row_number()over()::bigint as gid , * from w_prepare
;


-- Etape 3 : erosion du buffer de -50m

drop table if exists  _agent_acoudart.__bati_buff_bvsn ;
create table _agent_acoudart.__bati_buff_bvsn as 
with 
w_prepare as (
select st_multi(st_collectionextract(st_buffer(bb.geom,-50,'join=mitre mitre_limit=1'),3))::geometry(multipolygon,2154) as geom
from _agent_acoudart.__bati_buff_bvsn_step02  bb
)
select row_number()over()::bigint as gid , * from w_prepare
;

Construction table : Zone urbaine perméable[modifier]

Les zones zones urbaines perméables ont été constituées avec la même logique que pour la table précédente en utilisant les données piste aérodrome et terrain de sport de la bd topo de l'IGN.

Note :Attention utilisation de l'OCE ING (r020_territoire_occupation_du_sol.t_ocsge_ign_v1_1_2013) qui ne couvre pas tout le territoire !! Voir avec A.Coudart pourquoi ce choix.

Répertoire de travail: O:\TRAVAUX\ariviere\2018_ASSEMBLAGE_REFERENTIELS\2018_URBAIN

• Pour afficher le script SQL cliquer sur "Afficher"

/*
Note Attention utilisation de l'OCE ING (r020_territoire_occupation_du_sol.t_ocsge_ign_v1_1_2013)  qui ne couvre pas tout le territoire !!
Voir avec A.Coudart pourquoi ce choix.
*/

/*
Duration 50 min sur 4 deptartements
*/


-- Etape 1 : buffers de 50m sur les différentes couches de bati et union dans une table unique


set   work_mem='20000MB' ;
drop  table if exists __bati_buff_perm_bvsn_step01 ;
create TEMP table __bati_buff_perm_bvsn_step01 as 


with 
w_prepare as (
select pa.geom, st_multi(st_buffer(pa.geom,50,'join=mitre mitre_limit=2'))::geometry(multipolygon, 2154) as geom_b50
from xx_99_50_10_import_wfs.t_ign_bdtopo_piste_aerodrome pa

union all

select ts.geom, st_multi(st_buffer(ts.geom,50,'join=mitre mitre_limit=2'))::geometry(multipolygon, 2154) as geom_b50
from xx_99_50_10_import_wfs.t_ign_bdtopo_terrain_sport ts

)
 
select row_number()over()::bigint as gid,*
from w_prepare
;

CREATE INDEX sidx___bati_buff_perm_bvsn_step01
  ON __bati_buff_perm_bvsn_step01
  USING gist
  (geom);

-- Etape 2 : fusion des buffers qui se touchent


drop  table if exists __bati_buff_perm_bvsn_step02 ;
create TEMP table __bati_buff_perm_bvsn_step02 as 
 
with
w_prepare as (
select
(st_dump(st_collectionextract(st_union(geom_b50),3))).geom::geometry(polygon,2154) as geom
 from __bati_buff_perm_bvsn_step01  bb
)
select row_number()over()::bigint as gid , * from w_prepare
;

CREATE INDEX sidx___bati_buff_perm_bvsn_step02
  ON __bati_buff_perm_bvsn_step02
  USING gist
  (geom);


-- Etape 3 : erosion du buffer de -50m


drop table if exists  r020_territoire_occupation_du_sol.t_zone_batie_permeable_eptbsn_v01 ;
create table r020_territoire_occupation_du_sol.t_zone_batie_permeable_eptbsn_v01 as 

with 
w_prepare as (
select
st_multi(st_collectionextract(st_buffer(bb.geom,-50,'join=mitre mitre_limit=1'),3))::geometry(multipolygon,2154) as geom
from _agent_acoudart.__bati_buff_perm_bvsn_step02  bb

union all

select 
st_multi(st_intersection(ocs.geom, urb.geom))::geometry(multipolygon,2154) as geom
from 
r020_territoire_occupation_du_sol.t_zone_batie_eptbsn_v01 urb
, r020_territoire_occupation_du_sol.t_ocsge_ign_v1_1_2013 ocs
where
st_intersects(ocs.geom, urb.geom)
AND ocs.code_cs ilike 'CS2.2.1'

)
select row_number()over()::bigint as gid , * from w_prepare
;

CREATE INDEX sidx___bati_permeable_bvsn
  ON r020_territoire_occupation_du_sol.t_zone_batie_permeable_eptbsn_v01
  USING gist
  (geom);

Réseaux routiers et voies ferrées[modifier]

Le réseaux routier et les voies ferrées proviennent toutes les deux de la bd topo de l'IGN. Initialement ces données "linestring" ont été converties en polygones. Les largeurs appliquées sont celles contenus dans la table attributaire, si la largeur des routes n'est pas indiquée, une convention de 2 m à été appliqué pour les chemin et routes empierrées et de 5 m pour toutes les autres routes. Les routes de nature "bac" ont été retirée car ne elles n'ont pas de représentation physique sur le terrain
Les voies ferrées sur le territoire ont toutes une largeur dite normale, c'est à dire environ 1.5m qui correspondent à l'écartement des rails. sachant que les voies sont dans la majorité des cas au moins au nombre de deux et en comptant les ballastes et poteaux supportant les caténaires et après mesure à différents endroits via photo satellite il apparait qu'une largeur moyenne de 15 mètres semble pertinente.

Boisements[modifier]

Les boisements sont issus de la bd foret v.2, ils sont intégrés tel quel dans la table d'occupation du sol locale.

Surfaces agricoles[modifier]

Les éléments agricoles sont issus du RPG 2016, en particulier le registre parcellaire (et non par îlots) de ce dernier. Ces données sont issues des déclarations des exploitants agricoles à la PAC (Politique Agricole Commune), n'y seront donc présente que les parcelles qui font l'objet d'une déclaration, par exemple les usages agricoles "de loisir" seront absent (prairie équine de particulier, micro exploitation, exploitation agricole ne demandant pas d'aides européennes...).
La donnée d'origine comprend 246 typologies différentes regroupées en 23 groupes, le tableau ci-dessous propose un regroupement en grandes catégories pour en facilité l'exploitation ultérieure.

Codification occupation du sol OCS Théia[modifier]

Pour compléter les données issu de l'IGN, nous avons fait le choix d'utiliser le travail du pôle Théia. La codification de l'OCS Théia est consultable ici.

Assemblage des différents type d'occupation du sol[modifier]

Ces données d'occupation du sol pouvant être utile pour de nombreux besoins, il est convenu de les réunirent en une seule et même table. La table a été construite avec les références citées ci-dessus et en comblant les manques avec l'OCS du pôle Theia :

La construction suit une règle de priorité, chaque couche ajouter se découpant aux précédentes suivant l'ordre suivant :

1. Réseau routier et ferré
2. Surface en eau
3. Boisement
4. RPG
5. Urbain perméable
6. Urbain imperméable
7. OCS Théia

La table est composé des champs de base suivant :

• gid : identifiant unique, clé primaire
• data_theme : Thème principal (ex : Urbain)
• data_subtheme : Sous thème (ex : Perméable)
• data_utility : Utilité de la donnée (ex : Etude TBV)
• source_schema : schéma source donnée origine
• source_table : table source donnée origine
• source_id : identifiant source donnée origine
• source_data : Json regroupant l'intégralité des données de la source d'origine (champs, valeur)
• data_create_query : Requête SQL de création de donnée
• data_insert_date : Date d'insertion originel de la donnée dans l'OCS général
• data_update_date : Date de mise à jour de la donnée
• geom : Géométrie de la couche

Cette organisation nous permet de tracer la donnée d'origine et le requêtage en est facilité.

Assemblage final de la table OCS bvsn[modifier]

Ancien fichier : O:\PROJETS\2018_Tetes_BV\2018_tetes_bv_acoudart\00_final_all_files\04_scripts_sql\Annexe_9 construction_bd_ocs_local_union_all_v03.sql"

-- AJOUT DE LA COLONNE niveau_01 :
ALTER TABLE r020_territoire_occupation_du_sol.t_ocs_epsn_v01 ADD COLUMN niveau_01 character varying(50);
-- MISE A JOUR DE LA COLONNE niveau_01
update r020_territoire_occupation_du_sol.t_ocs_epsn_v01  set niveau_01
=case   data_theme 
     when 'Agriculture' then 'SURFACE AGRICOLE'
     when 'bâti dense' then 'SURFACE BATIE'
   when   'bâti diffus' then 'SURFACE BATIE'
   when   'Boisement' then 'SURFACE BOISEE'
   when   'culture été' then 'SURFACE AGRICOLE'
   when   'culture hiver' then 'SURFACE AGRICOLE'
   when   'eau' then 'SURFACE EN EAU'
   when   'forêt conifères' then 'SURFACE BOISEE'
   when   'forêt feuillus' then 'SURFACE BOISEE'
   
   when    'lande ligneuse' then 'SURFACE AGRICOLE'
   when   'pelouses' then 'SURFACE AGRICOLE'
   
   when   'plages et dunes' then 'SURFACE PLAGES ET DUNES'
   when   'prairie' then 'SURFACE AGRICOLE'
   when   'Surface en eau' then 'SURFACE EN EAU'
   when   'surface route' then 'SURFACE VOIERIE'
   when   'Transport' then 'SURFACE VOIERIE'
   when   'Urbain' then 'SURFACE BATIE'

   when   'verger' then 'SURFACE AGRICOLE'
   when    'vigne' then 'SURFACE AGRICOLE'
   when     'zones ind et com' then 'SURFACE INDUSTRIELLE ET COMMERCIALE'
     END


   ;select distinct niveau_01 from r020_territoire_occupation_du_sol.t_ocs_epsn_v01 

</interne>

Pour afficher le script d'assemblage de la table d'occupation du sol cliquez sur "afficher"

/*
attention pour eviter des erreurs de topology et la rapidité de création,  cette couche n'est calculée que sur un buffer de 2000 autour du bvsn
*/
set work_mem = '150000';

drop table if exists m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01;
create table m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01

(

gid serial primary key
, data_theme character varying(50)
, data_subtheme character varying(50)
, data_utility character varying(50)[]
, source_schema character varying(255)
, source_table character varying(255)
, source_id character varying(255)
, source_data json
, data_create_query text 
, data_insert_date date default now()
, data_update_date date default now()
, geom geometry(Multipolygon,2154)
, source_geom geometry(Multipolygon,2154)

)
;
CREATE INDEX sidx_t_ocs_epsn_v01_geom
  ON m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01
  USING gist
  (geom);

CREATE INDEX sidx_t_ocs_epsn_v01_source_geom
  ON m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01
  USING gist
  (source_geom);
/*
Import des données route topo
*/

INSERT INTO m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01(
             data_theme
            , data_subtheme
            , data_utility
            , source_schema
            , source_table
            , source_id
            , source_data
            , data_create_query
            , data_update_date
            , source_geom) 

with 
w_routes as(
select *
, (case when nature not in ('Chemin', 'Route empierrée') AND largeur = 0 then 5 
	when nature in ('Chemin', 'Route empierrée', 'Escalier') AND largeur = 0 then 2 
	else (largeur)::numeric end) as w_largeur
from xx_99_50_10_import_wfs.t_ign_bdtopo_routes
where nature not in ('Bac piéton')
AND st_intersects(geom, (select st_buffer(geom,2000)  from r020_territoire_hydrographie.t_bvsn_actu_carthage))
)

select 
  'Transport'
, null
, null
, 'xx_99_50_10_import_wfs'
, 't_ign_bdtopo_routes'
, id
, to_json(a) as source_data
, current_query()
, null
, ST_Multi(st_Collectionextract(st_makevalid(ST_RemoveRepeatedPoints(st_snaptogrid(st_buffer(geom, (w_largeur)::numeric,'quad_segs=2'),0.01))),3)) 

from w_routes a
;


-- update m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01 t_up set 
-- geom =NULL WHERE source_table = 't_ign_bdtopo_routes';

update m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01 t_up set 
geom = (SELECT 
   case when st_intersects(t_up.source_geom, st_union(t_in.source_geom)) THEN   st_multi(st_collectionextract(st_difference(t_up.source_geom, st_union(t_in.source_geom)),3))
   ELSE st_multi(t_up.source_geom)
   END :: geometry
from m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01  t_in
WHERE ST_intersects(t_up.source_geom, t_in.source_geom) and  t_up.gid > t_in.gid)
where geom is null;




/*
Import des données voie_ferrée topo
*/

INSERT INTO m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01(
             data_theme
            , data_subtheme
            , data_utility
            , source_schema
            , source_table
            , source_id
            , source_data
            , data_create_query
            , data_update_date
            , source_geom) 


select 
'Transport'
, null
, null 
, 'xx_99_50_10_import_wfs'
, 't_ign_bdtopo_troncon_voie_ferree'
, id
, to_json(a) as source_data
, current_query()
, null

, ST_Multi(st_Collectionextract(st_makevalid(ST_RemoveRepeatedPoints(st_snaptogrid(st_buffer(geom, 15, 'quad_segs=2'),0.01))),3)) 


from xx_99_50_10_import_wfs.t_ign_bdtopo_troncon_voie_ferree a
WHERE  st_intersects(geom, (select st_buffer(geom,2000)  from r020_territoire_hydrographie.t_bvsn_actu_carthage))
;


-- update m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01 t_up set 
-- geom =NULL WHERE source_table = 't_ign_bdtopo_troncon_voie_ferree';

update m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01 t_up set 
geom = (SELECT 
   case when st_intersects(t_up.source_geom, st_union(t_in.source_geom)) THEN   st_multi(st_collectionextract(st_difference(t_up.source_geom, st_union(t_in.source_geom)),3))
   ELSE st_multi(t_up.source_geom)
   END :: geometry
from m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01  t_in
WHERE ST_intersects(t_up.source_geom, t_in.source_geom) and  t_up.gid > t_in.gid)
where geom is null;






/*
Import des données surface en eau topo
*/

INSERT INTO m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01(
             data_theme
            , data_subtheme
            , data_utility
            , source_schema
            , source_table
            , source_id
            , source_data
            , data_create_query
            , data_update_date
            , source_geom) 


select 
'Surface en eau'
, null
, null 
, 'xx_99_50_10_import_wfs'
, 't_ign_bdtopo_surface_eau'
, id
, to_json(a) as source_data
, current_query()
, null
, ST_Multi(st_Collectionextract(st_makevalid(ST_RemoveRepeatedPoints(st_snaptogrid(geom,0.01))),3)) 

from xx_99_50_10_import_wfs.t_ign_bdtopo_surface_eau a
WHERE  st_intersects(geom, (select st_buffer(geom,2000)  from r020_territoire_hydrographie.t_bvsn_actu_carthage))
;


-- update m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01 t_up set 
-- geom =NULL WHERE source_table = 't_ign_bdtopo_surface_eau';

update m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01 t_up set 
geom = (SELECT 
   case when st_intersects(t_up.source_geom, st_union(t_in.source_geom)) THEN   st_multi(st_collectionextract(st_difference(t_up.source_geom, st_union(t_in.source_geom)),3))
   ELSE st_multi(t_up.source_geom)
   END :: geometry
from m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01  t_in
WHERE ST_intersects(t_up.source_geom, t_in.source_geom) and  t_up.gid > t_in.gid)
where geom is null;

/*
Import des données boisement topo
*/

INSERT INTO m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01(
             data_theme
            , data_subtheme
            , data_utility
            , source_schema
            , source_table
            , source_id
            , source_data
            , data_create_query
            , data_update_date
            , source_geom) 


select 
'Boisement'
, null
, null 
, 'xx_99_50_10_import_wfs'
, 't_ign_bdforet_v2'
, id
, to_json(a) as source_data
, current_query()
, null
, ST_Multi(st_Collectionextract(st_makevalid(ST_RemoveRepeatedPoints(st_snaptogrid(geom,0.01))),3)) 

from xx_99_50_10_import_wfs.t_ign_bdforet_v2 a
WHERE  st_intersects(geom, (select st_buffer(geom,2000)  from r020_territoire_hydrographie.t_bvsn_actu_carthage))
;


-- update m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01 t_up set 
-- geom =NULL WHERE source_table = 't_ign_bdforet_v2';

update m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01 t_up set 
geom = (SELECT 
   case when st_intersects(t_up.source_geom, st_union(t_in.source_geom)) THEN   st_multi(st_collectionextract(st_difference(t_up.source_geom, st_union(t_in.source_geom)),3))
   ELSE st_multi(t_up.source_geom)
   END :: geometry
from m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01  t_in
WHERE ST_intersects(t_up.source_geom, t_in.source_geom) and  t_up.gid > t_in.gid)
where geom is null;


/*
Import des données RPG 
--- Attention réalisation d'un ST_buffer negatif et positif lors de l'injection de la source_geom
*/
set work_mem = '150000';
INSERT INTO m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01(
             data_theme
            , data_subtheme
            , data_utility
            , source_schema
            , source_table
            , source_id
            , source_data
            , data_create_query
            , data_update_date
            , source_geom) 


select 
'Agriculture'
, null
, null 
, 'r080_agriculture'
, 't_rpg_2016_parcelles_graphiques'
, id_parcel
, to_json(a) as source_data
, current_query()
, null
, st_multi(st_collectionextract(st_makevalid(st_buffer(st_buffer(st_buffer(st_buffer(geom,-0.01,'join=mitre mitre_limit=5.0'),0.01,'join=mitre mitre_limit=5.0'),0.01,'join=mitre mitre_limit=5.0'),-0.01,'join=mitre mitre_limit=5.0')),3))

from r080_agriculture.t_rpg_2016_parcelles_graphiques a
WHERE  st_intersects(geom, (select st_buffer(geom,2000)  from r020_territoire_hydrographie.t_bvsn_actu_carthage))
;
update m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01 t_up set 
source_geom =ST_Multi(st_Collectionextract(st_makevalid(ST_RemoveRepeatedPoints(st_snaptogrid(source_geom,0.01))),3)) where  source_table = 't_rpg_2016_parcelles_graphiques';

---delete from m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01 WHERE source_table = 't_rpg_2016_parcelles_graphiques'
-- update m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01 t_up set 
-- geom =NULL WHERE source_table = 't_rpg_2016_parcelles_graphiques';

set work_mem = '150000';
update m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01 t_up set 
geom = (SELECT 
   case when st_intersects(t_up.source_geom, st_union(t_in.source_geom)) THEN   st_multi(st_collectionextract(st_difference(t_up.source_geom, st_union(t_in.source_geom)),3))
   ELSE st_multi(t_up.source_geom)
   END :: geometry
from m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01  t_in
WHERE ST_intersects(t_up.source_geom, t_in.source_geom) and  t_up.gid > t_in.gid)
where geom is null;

/*
Import des données Urbain permeable 
*/
set work_mem = '150000';
INSERT INTO m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01(
             data_theme
            , data_subtheme
            , data_utility
            , source_schema
            , source_table
            , source_id
            , source_data
            , data_create_query
            , data_update_date
            , source_geom) 


select 
'Urbain'
, 'Permeable'
, null 
, '_agent_acoudart'
, '__bati_buff_perm_bvsn_step03'
, gid
, to_json(a) as source_data
, current_query()
, null
,  ST_Multi(st_Collectionextract(st_makevalid(ST_RemoveRepeatedPoints(st_snaptogrid(geom,0.01))),3))

from _agent_acoudart.__bati_buff_perm_bvsn_step03 a
WHERE  st_intersects(geom, (select st_buffer(geom,2000)  from r020_territoire_hydrographie.t_bvsn_actu_carthage))
;

-- delete from m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01 WHERE source_table = '__bati_buff_perm_bvsn_step03'
-- update m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01 t_up set 
-- geom =NULL WHERE source_table = '__bati_buff_perm_bvsn_step03';

update m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01 t_up set 
geom = (SELECT 
   case when st_intersects(t_up.source_geom, st_union(t_in.source_geom)) THEN   st_multi(st_collectionextract(st_difference(t_up.source_geom, st_union(t_in.source_geom)),3))
   ELSE st_multi(t_up.source_geom)
   END :: geometry
from m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01  t_in
WHERE ST_intersects(t_up.source_geom, t_in.source_geom) and  t_up.gid > t_in.gid)
where geom is null;

/*
Import des données Urbain 
*/
set work_mem = '150000';
INSERT INTO m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01(
             data_theme
            , data_subtheme
            , data_utility
            , source_schema
            , source_table
            , source_id
            , source_data
            , data_create_query
            , data_update_date
            , source_geom) 


select 
'Urbain'
, null
, null 
, '_agent_acoudart'
, '__bd_urb_v02'
, gid
, to_json(a) as source_data
, current_query()
, null
,  ST_Multi(st_Collectionextract(st_makevalid(ST_RemoveRepeatedPoints(st_snaptogrid(geom,0.01))),3))

from _agent_acoudart.__bd_urb_v02 a
WHERE  st_intersects(geom, (select st_buffer(geom,2000)  from r020_territoire_hydrographie.t_bvsn_actu_carthage))
;



-- delete from m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01 WHERE source_table = '__bd_urb_v02'
-- update m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01 t_up set 
-- geom =NULL WHERE source_table = '__bd_urb_v02';

update m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01 t_up set 
geom = (SELECT 
   case when st_intersects(t_up.source_geom, st_union(t_in.source_geom)) THEN   st_multi(st_collectionextract(st_difference(t_up.source_geom, st_union(t_in.source_geom)),3))
   ELSE st_multi(t_up.source_geom)
   END :: geometry
from m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01  t_in
WHERE ST_intersects(t_up.source_geom, t_in.source_geom) and  t_up.gid > t_in.gid)
where geom is null;



/*
Import de l'OCS Théia
*/


set work_mem = '150000';
INSERT INTO m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01(
             data_theme
            , data_subtheme
            , data_utility
            , source_schema
            , source_table
            , source_id
            , source_data
            , data_create_query
            , data_update_date
            , source_geom) 

with 
-- grille sur les 4 deparements
   grille as (
      select * from  xx_99_utils.f_geo_makegrid2d_table(
          ---  4 depts
         -- (select st_union(geom) from r010_administratif_reglementaire.t_departements_2xxx_bd_topo 
         --             where "si_admin_metad_version_date" =(select max(si_admin_metad_version_date) from r010_administratif_reglementaire.t_departements_2xxx_bd_topo) 
         --          AND "si_admin_datageom_is_in_bvsn" IS TRUE
         --   )
       (select st_buffer(geom,2000)  from r020_territoire_hydrographie.t_bvsn_actu_carthage)
      ,1000,2154)
      
   )
-- decoupe theia avec grille
, decoupe_par_grille as (
    select 
         ocs.*
      ,  ST_intersection(ocs.geom, grille.geom) as geom_intersection_grille 
      from 
      r020_territoire_occupation_du_sol.t_ocs_ces_2017 ocs
      , grille
      WHERE ST_intersects (grille.geom, ocs.geom)
      
   )

select 
   CASE a.classe 
      WHEN 11 THEN 'culture été'
      WHEN 12 THEN 'culture hiver'
      WHEN 31 THEN 'forêt feuillus'
      WHEN 32 THEN 'forêt conifères'
      WHEN 34 THEN 'pelouses'
      WHEN 36 THEN 'lande ligneuse'
      WHEN 41 THEN 'bâti dense'
      WHEN 42 THEN 'bâti diffus'
      WHEN 43 THEN 'zones ind et com'
      WHEN 44 THEN 'surface route'
      WHEN 45 THEN 'surfaces minérales'
      WHEN 46 THEN 'plages et dunes'
      WHEN 51 THEN 'eau'
      WHEN 53 THEN 'glaciers ou neige'
      WHEN 211 THEN 'prairie'
      WHEN 221 THEN 'verger'
      WHEN 222 THEN 'vigne'
      ELSE '!!! OCS THEIA non classe !!!'
   END ::text
   , null
   , null 
   , 'r020_territoire_occupation_du_sol'
   , 't_ocs_ces_2017'
   , ogc_fid
   , to_json(a) as source_data
   , current_query()
   , null
   , st_multi(st_collectionextract(st_makevalid(a.geom_intersection_grille),3))::geometry(multipolygon,2154)
from decoupe_par_grille a
;


-- update m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01 
-- set data_theme = (SELECT CASE data_theme::int
--       WHEN 11 THEN 'culture été'
--       WHEN 12 THEN 'culture hiver'
--       WHEN 31 THEN 'forêt feuillus'
--       WHEN 32 THEN 'forêt conifères'
--       WHEN 34 THEN 'pelouses'
--       WHEN 36 THEN 'lande ligneuse'
--       WHEN 41 THEN 'bâti dense'
--       WHEN 42 THEN 'bâti diffus'
--       WHEN 43 THEN 'zones ind et com'
--       WHEN 44 THEN 'surface route'
--       WHEN 45 THEN 'surfaces minérales'
--       WHEN 46 THEN 'plages et dunes'
--       WHEN 51 THEN 'eau'
--       WHEN 53 THEN 'glaciers ou neige'
--       WHEN 211 THEN 'prairie'
--       WHEN 221 THEN 'verger'
--       WHEN 222 THEN 'vigne'
--       ELSE '!!! OCS THEIA non classe !!!'
--    END ::text
--    )
--    WHERE source_table = 't_ocs_ces_2017'
  
-- update m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01 t_up set geom =NULL WHERE source_table = 't_ocs_ces_2017';

/*
update m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01 t_up set 
geom = st_collectionextract(st_makevalid(geom),3)
WHERE geom is not null;

DELETE FROM m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01
WHERE source_table = 't_ocs_ces_2017';

*/

-- attention pour eviter des erreurs de topology cette couche n'est calculée que sur un buffer de 2000 autour du bvsn
/*set work_mem = '150000';
update m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01 t_up set 
geom = (
   SELECT 
      case when st_intersects(t_up.source_geom, st_collectionextract(st_makevalid(st_union(st_makevalid(t_in.source_geom))),3)) THEN   st_multi(st_collectionextract(xx_99_utils.st_difference_safe(t_up.source_geom, st_makevalid(st_collect(st_makevalid(t_in.source_geom)))),3))
            ELSE st_multi(t_up.source_geom)
            END :: geometry
         from m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01  t_in
         WHERE ST_intersects(t_up.source_geom, t_in.source_geom) and  t_up.gid > t_in.gid
   )
where geom is null ;
*/

-- update m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01 t_up set 
-- geom =ST_Multi(st_Collectionextract(st_makevalid(ST_RemoveRepeatedPoints(st_snaptogrid(source_geom,0.01))),3)) where  source_table = 't_ign_bdtopo_routes';
-- 
-- update m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01 t_up set 
-- geom =ST_Multi(st_Collectionextract(st_makevalid(ST_RemoveRepeatedPoints(st_snaptogrid(source_geom,0.01))),3)) where  source_table = 't_ign_bdtopo_routes';

-- SET postgis.backend = sfcgal;
-- set work_mem = '150000';
-- update m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01 t_up set 
-- geom = (
--    SELECT 
--       case when st_intersects(t_up.source_geom, st_collectionextract(st_makevalid(st_union(st_makevalid(t_in.geom))),3)) THEN   st_multi(st_collectionextract(xx_99_utils.st_difference_safe(t_up.source_geom, st_makevalid(st_collect(st_makevalid(t_in.geom)))),3))
--             ELSE st_multi(t_up.source_geom)
--             END :: geometry
--          from m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01  t_in
--          WHERE ST_intersects(t_up.source_geom, t_in.geom) and  t_up.gid > t_in.gid
--    )
-- where geom is null ;


-- update m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01 t_up set 
-- geom = (SELECT 
--    case when st_intersects(t_up.source_geom, st_union(t_in.source_geom)) THEN   st_multi(st_collectionextract(st_difference(t_up.source_geom, st_union(t_in.source_geom)),3))
--    ELSE st_multi(t_up.source_geom)
--    END :: geometry
-- from m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01  t_in
-- WHERE ST_intersects(t_up.source_geom, t_in.source_geom) and  t_up.gid > t_in.gid)
-- where geom is null;

-- SET postgis.backend = sfcgal;
--update m060_milieux_continuite_tetes_de_bv.t_ocs_epsn_v01 t_up set geom =NULL WHERE source_table = 't_ocs_ces_2017';


create OR  REPLACE function m060_milieux_continuite_tetes_de_bv.f_decoupe_geom_source_pour_combler_ocs_geom(
_schemaname name
, _tablename name
, _pkname NAME
, _source_geomname name
, _dest_geomname NAME
)
returnS text as 

$BODY$
   DECLARE
   /*
   Découpe geom_source_geom pour occuper l'esopace libre dans l'ensemble de la table. (utile pour couche de couverture de type occupation du sol)

   Exemple d'appel :

   SELECT m060_milieux_continuite_tetes_de_bv.f_decoupe_geom_source_pour_combler_ocs_geom (
        _schemaname := 'm060_milieux_continuite_tetes_de_bv'    -- nom du schema de la table de couverture
      , _tablename := 't_ocs_epsn_v01'                        -- nom  de la table de couverture
      , _pkname:='gid'                                         -- nom de la colonne primary key
      , _source_geomname := 'source_geom'                      -- nom de la colonne geom à découper                               
      , _dest_geomname := 'geom'                               -- nom de la colonne geom qui contiendra le résultat du découpage                                             

      )
      
   notes de version :

   - V01 : first version

   */
   mode_debug BOOLEAN := TRUE ; ---true;  -- TRUE  renvoi des notices

   lst_gids record;
   gid_to_up bigint;
   i int;
   str_sql_lst_gids text;
   str_sql_lst_gids_count TEXT;
   str_sql_up_geom text;
   return_txt text :='';
   nb_gid_to_up bigint;



   
   BEGIN

   -- 
   -- 1 . Liste des gids à découper :
   -- 
   str_sql_lst_gids  := format (
      $SQL$
      SELECT %I as gid FROM %I.%I WHERE %I IS NULL
      
      $SQL$
      ,_pkname, _schemaname , _tablename , _dest_geomname
      )
      ;
      
   str_sql_lst_gids_count  := format (
      $SQL$
      SELECT count(%I) as nb FROM %I.%I WHERE %I IS NULL
      
      $SQL$
      ,_pkname, _schemaname , _tablename , _dest_geomname
      )
      ;
   Execute str_sql_lst_gids_count into nb_gid_to_up;
   i:=0;
   For lst_gids in EXECUTE str_sql_lst_gids LOOP
   i:=i+1;
   gid_to_up := lst_gids.gid;

   str_sql_up_geom := format (
      $SQL$
      update %I.%I t_up set 
      %I= (
         SELECT 
            case when st_intersects(t_up.%I, st_collectionextract(st_makevalid(st_union(st_makevalid(t_in.%I))),3)) 
               THEN   st_multi(st_collectionextract(St_difference(t_up.%I, st_collectionextract(st_makevalid(st_union(st_makevalid(t_in.%I))),3)),3))
                  ELSE st_multi(t_up.%I)
                  END :: geometry
               from %I.%I  t_in
               WHERE ST_intersects(t_up.%I, t_in.%I) 
                  and  t_up.%I > t_in.%I
         )
      where %I = %s;
      $SQL$
       , _schemaname , _tablename 
       , _dest_geomname
       , _source_geomname, _source_geomname
       , _source_geomname, _source_geomname
       , _source_geomname
       , _schemaname , _tablename 
       , _source_geomname, _source_geomname
       , _pkname, _pkname
       , _pkname, gid_to_up
       
      );
      --RAISE NOTICE 'str_sql_up_geom : %', str_sql_up_geom;
      RAISE NOTICE '% / % -> gid_to_up : % ', i, nb_gid_to_up, gid_to_up;
      begin
         execute str_sql_up_geom;
         EXCEPTION WHEN OTHERS THEN 
            str_sql_up_geom := format (
               $SQL$
               update %I.%I t_up set 
               %I= (
                  SELECT 
                     case when st_intersects(t_up.%I, st_collectionextract(st_makevalid(st_union(st_makevalid(ST_Buffer(t_in.%I, 0.00001)))),3)) 
                        THEN   st_multi(st_collectionextract(St_difference(t_up.%I, st_collectionextract(st_makevalid(st_union(st_makevalid(ST_Buffer(t_in.%I,0.00001)))),3)),3))
                           ELSE st_multi(t_up.%I)
                           END :: geometry
                        from %I.%I  t_in
                        WHERE ST_intersects(t_up.%I, t_in.%I) 
                           and  t_up.%I > t_in.%I
                  )
               where %I = %s;
               $SQL$
                , _schemaname , _tablename 
                , _dest_geomname
                , _source_geomname, _source_geomname
                , _source_geomname, _source_geomname
                , _source_geomname
                , _schemaname , _tablename 
                , _source_geomname, _source_geomname
                , _pkname, _pkname
                , _pkname, gid_to_up
                
               );
            begin
               execute str_sql_up_geom;
               EXCEPTION WHEN OTHERS THEN 
               
               return_txt := return_txt || ',' || gid_to_up::text;
               END;
      end ;
      
   END LOOP;
   
   
   

   RETURN return_txt;
   END;
$BODY$
  LANGUAGE plpgsql VOLATILE
  COST 100;


  

   SELECT m060_milieux_continuite_tetes_de_bv.f_decoupe_geom_source_pour_combler_ocs_geom (
        _schemaname := 'm060_milieux_continuite_tetes_de_bv'    -- nom du schema de la table de couverture
      , _tablename := 't_ocs_epsn_v01'                        -- nom  de la table de couverture
      , _pkname:='gid'                                         -- nom de la colonne primary key
      , _source_geomname := 'source_geom'                      -- nom de la colonne geom à découper                               
      , _dest_geomname := 'geom'                               -- nom de la colonne geom qui contiendra le résultat du découpage                                             

      )

Construction référentiel ripisylve ( _bd_haie_vege)[modifier]

Dans un soucis de cohérence et d'homogénéité territoriale et pour analyser la densité de ripisylve nous avons fait le choix de ne pas utiliser les seules données DEC. Nous avons donc agréger les données végétation de la bd topo aux données DEC quand elles étaient disponibles.

1ère étape : Découper les DEC avec l'ensemble de la végétation bd topo et buffer de 4 mètres (8mètres de houppier)

2ème étape : assembler ce découpage avec la bd topo végétation

Il est possible de rechercher les haies seules sur cette couche au moyen d'une requête sur l'attribut 'nature'. where nature ilike 'haie' OR nature is NULL.

NULL = données DEC (couche de haie, il n'y a que des haies dedans). haie = donnée bd topo végétation avec filtre sur les haies.

/* Création de la table haie issu des DEC découpés suivant la bd topo zone végétation et bufferisé à 4 m */

drop table if exists _agent_acoudart._haies_buff_dec_difference_bdtopo;
create table _agent_acoudart._haies_buff_dec_difference_bdtopo as

with w_prepare as (
select *, st_multi(st_buffer(geom, 4))::geometry(multipolygon, 2154) as geom_buff
from 
m060_milieux_continuite_dec_agreg.t_haies_agreg
)

SELECT
xx_99_utils.st_difference_geom_from_table (
       _geom := geom_buff                                                               -- Geom de base                                                                                                                               
      , _schemaname := 'xx_99_50_10_import_wfs'                         -- Schema de la table contenant les infos à intersecter avec la base
      , _tablename:= 't_ign_bdtopo_zone_vegetation'                      -- Table contenant les infos à intersecter avec la base
      , _geomname := 'geom'                                                            -- Geom des infos à intersecter avec la base
      , _whereclause := NULL                                                            -- Type geom de base (POLYGON, LINESTRING, POINT)
      , _return_geom_typeid := 3
      , _use_st_safe := FALSE
      , _st_makevalid := TRUE                                                           -- [optionnel] : Clause Where (TODO AR : non implémenté)
      ) ::geometry(multipolygon,2154) as geom_diff
      ,*
   FROM 
w_prepare

/* Création de la table haie avec agrégat bd topo végétation (prioritaire) et la table haie DEC créée précédemment (secondaire) */

drop table if exists _agent_acoudart._bd_haie_vege;
create table _agent_acoudart._bd_haie_vege as


WITH w_prepare as(

SELECT
gid2 as id_source
, geom_diff as geom
, '' as nature 
FROM
_agent_acoudart._haies_buff_dec_difference_bdtopo

UNION ALL

SELECT 
gid as id_source
, geom
, nature
FROM
xx_99_50_10_import_wfs.t_ign_bdtopo_zone_vegetation
)

SELECT
row_number()over()::bigint as gid
, *
FROM
w_prepare 

Suite des traitements[modifier]

• Calcul des différents indices nécessaires pour le calcul de l'indice de pression

http://processing.observatoire.sevre-nantaise.com/processing-chain/5 (indice_pression_tbv)

• Calcul de tous les paramètres nécessaires pour l'évaluation des pressions sur les enjeux

http://processing.observatoire.sevre-nantaise.com/processing-chain/8 (indicateurs_caracterisation_tbv)

dans cette chaine le process /tbv_scenario_calcul.php permet de calculer l'indice de pression (qui est lui même issu d'une série d'indicateurs pondérés calculés à l'étape précédente.

• Centrage réduction des résultats des paramètres de pression

http://processing.observatoire.sevre-nantaise.com/process?file=/var/www/observ/processing/src/App/Processes/Tr_TBV/Caracterisation/tbv_discretisation.php

• Calcul des notes de pression qe_syloa, qte_sylo et qm_syloa et de la sensibilité : sensibilite_syloa

Dans le cas des ruisseaux salmonicoles on utilise des pondérations spécifiques (liste d'indicateur plus complète) pour les enjeux : qm_eptbsn_rsalmo_01, qe_eptbsn_rsalmo_01, qte_eptbsn_rsalmo_01,

http://processing.observatoire.sevre-nantaise.com/process?file=/var/www/observ/processing/src/App/Processes/Tr_TBV/Caracterisation/tbv_scenario_calcul.php pour le qe_syloa, qte_sylo et qm_syloa et pour sensibilite_syloa (qm_eptbsn_rsalmo_01, qe_eptbsn_rsalmo_01, qte_eptbsn_rsalmo_01, sensibilite_syloa pour rsalmo)

• Centrage réduction des notes de pression et sensibilité

http://processing.observatoire.sevre-nantaise.com/process?file=/var/www/observ/processing/src/App/Processes/Tr_TBV/Caracterisation/tbv_discretisation.php pour les 1360,1361,1362,1363

• Assignation des notes de vulnérabilité

http://processing.observatoire.sevre-nantaise.com/process?file=/var/www/observ/processing/src/App/Processes/Tr_TBV/Caracterisation/tbv_vulnerabilite.php par enjeu

• Rafraichissement de la VM des résultats

/*
Pour mémoire code VM + code création index geo pour tous les id_param
*/
create materialized view m060_milieux_continuite_tetes_de_bv.vm_tbv_01_with_result_geom
as select r.* , p.nom_param,
g.geom
from 
 m060_milieux_continuite_tetes_de_bv.t_tbv_01_with_result r
LEFT JOIN m060_milieux_continuite_tetes_de_bv.t_tbv_01_param_carac p on (p.id_param=r.id_param)
LEFT JOIN m060_milieux_continuite_tetes_de_bv.t_tbv_01_without_caract g on (g.gid = r.gid_tbv)


/*
a jouer et executer le retour :
*/

select 'CREATE index if not exists idx_vm_tbv_01_with_result_geom'|| id_param ||' ON m060_milieux_continuite_tetes_de_bv.vm_tbv_01_with_result_geom using gist (geom) WHERE id_param = ' || id_param || ';'


FROM m060_milieux_continuite_tetes_de_bv.t_tbv_01_param_carac
union all 

select 'CREATE index if not exists idx_vm_tbv_01_with_result_geom_2_'|| id_param ||' ON m060_milieux_continuite_tetes_de_bv.vm_tbv_01_with_result_geom using gist (geom) WHERE nom_param = ''' || nom_param || ''';'


FROM m060_milieux_continuite_tetes_de_bv.t_tbv_01_param_carac

TODO AR : à organiser : projet qgis : O:\PROJETS\2018_Tetes_BV\2018_tetes_bv_acoudart\04_travail\prj\tbv_result.qgs

zone urbaine: "O:\PROJETS\2018_Tetes_BV\2018_tetes_bv_acoudart\05_rapport\annexes\Scripts SQL\Annexe_8 construction_zones_urbaines.sql"

haie :