BasinDelineation.py — Einzugsgebiets-Definition¶

Quelldatei

src/BasinDelineation.py

Identifiziert alle Rasterzellen innerhalb eines Einzugsgebiets. Nutzt die Abflussrichtungsdatei (G_OUTFLC.UNF4) für das stromaufwärts-Routing oder alternativ einen Shapefile-Overlay mit dem WaterGAP-Referenzgitter.

import pandas as pd
from typing import List
from BinaryFileHandler import ReadBin
import Paths_and_params as PP

Parametrisierung aus Konfigurationsdatei (config.yaml)

continent_index: int = PP.continent_index
outflowpath: str = PP.outflow_file
MostDownstreamCell: List[int] = PP.downstream_cells

def DelineateBasin(
    MostDownstreamCell: List[int],
    outflowpath: str,
    ng: List[int],
    continent_index: int
) -> List[int]:

Rekursive Bestimmung des Einzugsgebiets mittels Routing-Verfolgung.

Ausgehend von der Auslass-Zelle werden iterativ alle oberstromigen Zellen identifiziert, indem die Routing-Tabelle rückwärts durchsucht wird. Der Algorithmus terminiert, wenn keine neuen oberstromigen Zellen mehr gefunden werden (Zustand der Konvergenz erreicht).

Parameter

Parameter	Typ	Beschreibung
`MostDownstreamCell`	`List[int]`	Liste der Auslass-Zellen-IDs, die als Startpunkt der Rückwärts-Verfolgung dienen.
`outflowpath`	`str`	Pfad zur Binärdatei (G_OUTFLC.UNF4) mit Routing-Informationen.
`ng`	`List[int]`	Liste mit der Anzahl von Rasterzellen pro Kontinent.
`continent_index`	`int`	Index des Kontinents in den ng-Array.

Rückgabe

List[int]

Liste aller Rasterzellen-IDs innerhalb des Einzugsgebiets.

Routing-Daten aus Binärdatei laden

    OutflowData = ReadBin(outflowpath, ng[continent_index])

Mapping-Tabelle erstellen: Rasterzellen-ID → Oberstrom-Zelle

    DataTable = {'GCRC': list(range(1,ng[continent_index]+1)), 'DownstreamCell': OutflowData}
    df2 = pd.DataFrame(DataTable)
    df2.set_index('GCRC')

Initialisierung: Basismenge mit Auslass-Zelle(n) (Kopie erstellen)

    Basincells = list(MostDownstreamCell)
    NewDownstreamcells = list(MostDownstreamCell)
    temp = []

Iterative Rückwärts-Routing bis Konvergenz

    while temp != NewDownstreamcells:
        upstreamcells = [df2[df2['DownstreamCell'] == i]['GCRC'].values.tolist() for i in NewDownstreamcells]
        flat_list_upstreamcells = [item for sublist in upstreamcells for item in sublist]
        Basincells.extend(flat_list_upstreamcells)
        temp = NewDownstreamcells
        NewDownstreamcells = flat_list_upstreamcells
    return Basincells

Routing-basierte Einzugsgebietsdelineation durchführen

List_of_cells_in_the_basin = DelineateBasin(MostDownstreamCell, outflowpath, PP.ng, continent_index)
df6 = pd.DataFrame({'CellsInBasin': List_of_cells_in_the_basin})

Ergebnisse in Excel-Datei speichern

df6.to_excel("BasinCells.xlsx", header=False)

import geopandas as gpd
import matplotlib.pyplot as plt
import numpy as np

Shapefile-Overlay-Methode zur Einzugsgebietsdelineation

WaterGAP-Referenzgitter laden

mother_grid = gpd.read_file(PP.mother_grid_shapefile)

ANMERKUNG: Die Pfade zu Einzugsgebiet- und Speicherbecken-Shapefiles müssen projektspezifisch manuell angepasst werden. Sie liegen nicht in der Standard-Konfigurationsdatei vor. basin_shp_file = gpd.read_file("pfad/zum/einzugsgebiet.shp") reservoir_shp_file = gpd.read_file("pfad/zum/reservoir.shp")

Implementierungsbeispiel für Shapefile-Overlay-Verfahren

Räumliches Overlay durchführen: Schnittmenge von Gitter und Einzugsgebiet berechnen intersection = gpd.overlay(mother_grid, basin_shp_file, how='intersection')

Koordinatenreferenzsystem in flächentreue Projektion (Equal Area) transformieren¶

eqArea_intersection = intersection.to_crs(epsg=6933)

Flächen in Quadratkilometern berechnen¶

areas = eqArea_intersection.area / 1000000 intersection['Actual Area (Km2)'] = areas intersection['Actual Area (Km2)'].round(decimals=6)

Flächenanteil pro Rasterzelle berechnen (in Prozent)¶

intersection['Portion of cell in the basin'] = intersection['Actual Area (Km2)'] / intersection['AREA_KM2'] * 100

Ergebnistabelle strukturieren und als CSV-Datei speichern¶

columns_names = ['Cell_ID', 'Portion of Cell in Basin (%)'] Basin_Cells = pd.DataFrame(columns=columns_names) Basin_Cells['Cell_ID'] = intersection['ARCID'] Basin_Cells['Portion of Cell in Basin (%)'] = intersection['Portion of cell in the basin'] Basin_Cells.to_csv("List_of_cells_in_basin.csv")