{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Alpine monatliche Niederschlagsdaten" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Alpine monatliche Niederschlagsdaten seit 1871, abgeleitet aus stationären Beobachtungen: Der LAPrec-Datensatz (Long-term Alpine Precipitation Reconstruction) bietet monatliche Niederschlagsdaten für die Alpenregion und basiert auf stationären Beobachtungen. Es gibt zwei Versionen:\n", "\n", "* LAPrec1871 (ab 1871, 85 Eingabereihen)\n", "* LAPrec1901 (ab 1901, 165 Eingabereihen)\n", " \n", "Der Datensatz erfüllt hohe klimatische Standards und ist eine wertvolle Grundlage für historische Klimaanalysen in den Alpen, einer Region, die stark vom Klimawandel betroffen ist. Es kombiniert die Datenquellen HISTALP (homogenisierte stationäre Niederschlagsdaten) und APGD (tägliche Rasterdaten von 1971–2008) und nutzt die Methode der Reduced Space Optimal Interpolation (RSOI). LAPrec wird alle zwei Jahre aktualisiert und wurde im Rahmen des Copernicus Climate Change Service in Zusammenarbeit mit den Wetterdiensten von Schweiz (MeteoSwiss) und Österreich (ZAMG) entwickelt.\n", "\n", "**Informationen zum Datensatz:**\n", "* Quelle: [Alpine Monthly Precipitation](https://cds.climate.copernicus.eu/datasets/insitu-gridded-observations-alpine-precipitation?tab=overview)\n", "* Author: T. Tewes (City of Konstanz)\n", "* Notebook-Version: 1.1 (Aktualisiert: Dezember 05, 2024)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## 1. Festlegen der Pfade und Arbeitsverzeichnisse" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "import os\n", "\n", "''' ---- Verzeichnisse hier angeben ---- '''\n", "download_folder = r\".\\data\\alpine-monthly-precipitation\\download\"\n", "working_folder = r\".\\data\\alpine-monthly-precipitation\\working\"\n", "geotiff_folder = r\".\\data\\alpine-monthly-precipitation\\geotiff\"\n", "csv_folder = r\".\\data\\alpine-monthly-precipitation\\csv\"\n", "output_folder = r\".\\data\\alpine-monthly-precipitation\\output\"\n", "''' ----- Ende der Angaben ---- '''\n", "\n", "os.makedirs(download_folder, exist_ok=True)\n", "os.makedirs(working_folder, exist_ok=True)\n", "os.makedirs(geotiff_folder, exist_ok=True)\n", "os.makedirs(csv_folder, exist_ok=True)\n", "os.makedirs(output_folder, exist_ok=True)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## 2. Herunterladen und Entpacken des Datensatzes" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### 2.1 Authentifizierung" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "import cdsapi\n", "\n", "def main():\n", " # API-Key für die Authentifizierung\n", " api_key = \"fdae60fd-35d4-436f-825c-c63fedab94a4\"\n", " api_url = \"https://cds.climate.copernicus.eu/api\"\n", "\n", " # Erstellung des CDS-API-Clients\n", " client = cdsapi.Client(url=api_url, key=api_key)\n", " return client" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### 2.2 Definieren die „request“ und laden Sie den Datensatz herunter" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "# Definition des Datensatzes und der Request-Parameter\n", "dataset = \"insitu-gridded-observations-alpine-precipitation\"\n", "request = {\n", " \"variable\": \"precipitation\",\n", " \"dataset_issue\": [\n", " \"laprec1871\",\n", " \"laprec1901\"\n", " ],\n", " \"version\": [\"1_2\"],\n", "}" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "ename": "NameError", "evalue": "name 'os' is not defined", "output_type": "error", "traceback": [ "\u001b[1;31m---------------------------------------------------------------------------\u001b[0m", "\u001b[1;31mNameError\u001b[0m Traceback (most recent call last)", "Cell \u001b[1;32mIn[3], line 16\u001b[0m\n\u001b[0;32m 13\u001b[0m \u001b[38;5;28mprint\u001b[39m(\u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mDatensatz bereits heruntergeladen.\u001b[39m\u001b[38;5;124m\"\u001b[39m)\n\u001b[0;32m 15\u001b[0m \u001b[38;5;28;01mif\u001b[39;00m \u001b[38;5;18m__name__\u001b[39m \u001b[38;5;241m==\u001b[39m \u001b[38;5;124m\"\u001b[39m\u001b[38;5;124m__main__\u001b[39m\u001b[38;5;124m\"\u001b[39m:\n\u001b[1;32m---> 16\u001b[0m \u001b[43mmain_retrieve\u001b[49m\u001b[43m(\u001b[49m\u001b[43m)\u001b[49m\n", "Cell \u001b[1;32mIn[3], line 4\u001b[0m, in \u001b[0;36mmain_retrieve\u001b[1;34m()\u001b[0m\n\u001b[0;32m 2\u001b[0m \u001b[38;5;28;01mdef\u001b[39;00m\u001b[38;5;250m \u001b[39m\u001b[38;5;21mmain_retrieve\u001b[39m():\n\u001b[0;32m 3\u001b[0m dataset_filename \u001b[38;5;241m=\u001b[39m \u001b[38;5;124mf\u001b[39m\u001b[38;5;124m\"\u001b[39m\u001b[38;5;132;01m{\u001b[39;00mdataset\u001b[38;5;132;01m}\u001b[39;00m\u001b[38;5;124m.zip\u001b[39m\u001b[38;5;124m\"\u001b[39m\n\u001b[1;32m----> 4\u001b[0m dataset_filepath \u001b[38;5;241m=\u001b[39m \u001b[43mos\u001b[49m\u001b[38;5;241m.\u001b[39mpath\u001b[38;5;241m.\u001b[39mjoin(download_folder, dataset_filename)\n\u001b[0;32m 6\u001b[0m \u001b[38;5;66;03m# Den Datensatz nur herunterladen, wenn er noch nicht heruntergeladen wurde\u001b[39;00m\n\u001b[0;32m 7\u001b[0m \u001b[38;5;28;01mif\u001b[39;00m \u001b[38;5;129;01mnot\u001b[39;00m os\u001b[38;5;241m.\u001b[39mpath\u001b[38;5;241m.\u001b[39misfile(dataset_filepath):\n\u001b[0;32m 8\u001b[0m \u001b[38;5;66;03m# Rufen Sie den CDS-Client nur auf, wenn der Datensatz noch nicht heruntergeladen wurde.\u001b[39;00m\n", "\u001b[1;31mNameError\u001b[0m: name 'os' is not defined" ] } ], "source": [ "# Führen Sie es aus, um den Datensatz herunterzuladen:\n", "def main_retrieve():\n", " dataset_filename = f\"{dataset}.zip\"\n", " dataset_filepath = os.path.join(download_folder, dataset_filename)\n", "\n", " # Den Datensatz nur herunterladen, wenn er noch nicht heruntergeladen wurde\n", " if not os.path.isfile(dataset_filepath):\n", " # Rufen Sie den CDS-Client nur auf, wenn der Datensatz noch nicht heruntergeladen wurde.\n", " client = main()\n", " # Den Datensatz mit den definierten Anforderungsparametern herunterladen\n", " client.retrieve(dataset, request, dataset_filepath)\n", " else:\n", " print(\"Datensatz bereits heruntergeladen.\")\n", "\n", "if __name__ == \"__main__\":\n", " main_retrieve()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### 2.3 Extrahieren die ZIP-Datei in Ordner" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "import zipfile\n", "\n", "# Erstellen des Dateinamens und des Dateipfads für die ZIP-Datei des Datensatzes\n", "dataset_filename = f\"{dataset}.zip\"\n", "dataset_filepath = os.path.join(download_folder, dataset_filename)\n", "\n", "# Definieren Sie einen Extraktionsordner für die ZIP-Datei, der dem Arbeitsordner entspricht\n", "extract_folder = working_folder\n", "\n", "# Entpacken der ZIP-Datei\n", "try:\n", " os.makedirs(extract_folder, exist_ok=True)\n", " \n", " if not os.listdir(extract_folder):\n", " # Versuchen Sie, die ZIP-Datei zu öffnen und zu extrahieren\n", " with zipfile.ZipFile(dataset_filepath, 'r') as zip_ref:\n", " zip_ref.extractall(extract_folder)\n", " print(f\"Dateien erfolgreich extrahiert nach: {extract_folder}\")\n", " else:\n", " print(\"Ordner ist nicht leer. Entpacken überspringen.\")\n", "except FileNotFoundError:\n", " print(f\"Fehler: Die Datei {dataset_filepath} wurde nicht gefunden.\")\n", "except zipfile.BadZipFile:\n", " print(f\"Fehler: Die Datei {dataset_filepath} ist keine gültige ZIP-Datei.\")\n", "except Exception as e:\n", " print(f\"Ein unerwarteter Fehler ist aufgetreten: {e}\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## 3. Untersuchen der Metadaten der NetCDF4-Datei" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# Liste der netCDF4-Dateien im Arbeits-/extrahierten Ordner drucken\n", "filename_list = os.listdir(extract_folder)\n", "print(filename_list)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Der Datensatz Monatliche Niederschläge in den Alpen enthält Aufzeichnungen ab 1871 oder 1901. In diesem Notizbuch verwenden wir den Datensatz von 1901 (LAPrec1901)." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "import netCDF4 as nc\n", "\n", "# Definieren Sie den Dateipfad für den ausgewählten NetCDF-Datensatz\n", "nc_filename = \"LAPrec1901.v1.2.nc\"\n", "nc_filepath = os.path.join(extract_folder, nc_filename)\n", "\n", "# Öffnen der NetCDF-Datei im Lesemodus\n", "nc_dataset = nc.Dataset(nc_filepath, mode='r')\n", "\n", "# Auflisten aller Variablen im Datensatz\n", "variables_list = list(nc_dataset.variables.keys())\n", "print(f\"Verfügbare Variablen: {list(variables_list)}\")" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "import pandas as pd\n", "\n", "# Variablennamen aus vorhandenen Variablen definieren und Variablendaten lesen\n", "variable_name = 'LAPrec1901'\n", "variable_data = nc_dataset[variable_name]\n", "\n", "# Erstellen einer Zusammenfassung der Hauptvariablen\n", "summary = {\n", " \"Variablename\": variable_name,\n", " \"Datentyp\": variable_data.dtype,\n", " \"Form\": variable_data.shape,\n", " \"Variableinfo\": f\"{variable_name}({', '.join(variable_data.dimensions)})\",\n", " \"Einheiten\": getattr(variable_data, \"units\", \"N/A\"),\n", " \"Langer Name\": getattr(variable_data, \"long_name\", \"N/A\"),\n", "}\n", "\n", "# Anzeigen der Zusammenfassung des Datensatzes als DataFrame zur besseren Visualisierung\n", "nc_summary = pd.DataFrame(list(summary.items()), columns=['Beschreibung', 'Bemerkungen'])\n", "\n", "# Anzeigen des Zusammenfassungs-DataFrames\n", "nc_summary" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# Drucken Sie eine Zusammenfassung aller Variablen des Datensatzes\n", "rows = []\n", "for variable in variables_list:\n", " try:\n", " var_obj = nc_dataset.variables[variable]\n", " unit = getattr(var_obj, 'units', 'N/A')\n", " shape = var_obj.shape\n", " rows.append({\n", " \"nc_variablen\": variable,\n", " \"einheit\": unit,\n", " \"form\": shape\n", " })\n", " except Exception as e:\n", " print(f\"Fehler bei der Verarbeitung der Variable {variable}: {e}\")\n", "\n", "# Erstelle ein DataFrame\n", "df = pd.DataFrame(rows)\n", "df" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## 5. Exportieren der Datensatz im CSV-Format" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### 4.1 Begrenzungsrahmen zum Filtern von Daten für Konstanz definieren" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# Zusätzliche Anforderungsfelder definieren, um sicherzustellen, dass die Anfrage innerhalb der Dateigrößenbegrenzung bleibt.\n", "# Diese Koordinaten wurden mit dem BBox Extractor Tool erhalten:\n", "# https://str-ucture.github.io/bbox-extractor/\n", "\n", "# Begrenzungsbox für die Region Konstanz (WGS84-Projektion):\n", "bbox_wgs84_konstanz = [47.9, 8.9, 47.6, 9.3] # Format: [Norden, Westen, Süden, Osten]\n", "bbox_wgs84_konstanz_standard = [9.0, 47.6, 9.3, 47.8] # Standardformat: [Westen, Süden, Osten, Norden]" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# Alternativ können Sie ein Shapefile für eine präzise geografische Filterung verwenden\n", "import geopandas as gpd\n", "import matplotlib.pyplot as plt\n", "from pyproj import Transformer\n", "\n", "# Beispiel: Shapefile von Konstanz laden (WGS84-Projektion)\n", "shapefile = r\"./shapefiles/kn_boundary.shp\"\n", "gdf = gpd.read_file(shapefile)\n", "\n", "# Extrahieren Sie den Begrenzungsrahmen des Shapefiles\n", "shapefile_epsg = gdf.crs.to_epsg()\n", "shapefile_bounds = gdf.total_bounds\n", "\n", "# Definieren Sie die Quell- und Zielprojektionen\n", "proj_shapefile = f\"epsg:{shapefile_epsg}\" # Quellprojektion (CRS der Shape-Datei)\n", "proj_laea = \"epsg:3035\" # Zielprojektion (LAEA für Europa)\n", "\n", "# Definieren Sie einen Transformator, um die Koordinaten von der Quell- in die Zielprojektion umzuwandeln\n", "transformer = Transformer.from_crs(proj_shapefile, proj_laea, always_xy=True)\n", "\n", "# Transformieren Sie die Koordinaten des Begrenzungsrahmens in die Zielprojektion (LAEA)\n", "x_min, y_min = transformer.transform(shapefile_bounds[0], shapefile_bounds[1]) # Untere linke Ecke\n", "x_max, y_max = transformer.transform(shapefile_bounds[2], shapefile_bounds[3]) # Obere rechte Ecke\n", "shapefile_bounds_laea = [x_min, y_min, x_max, y_max] # Begrenzungsrahmen in LAEA-Projektion\n", "\n", "# Anzeige des transformierten Begrenzungsrahmens\n", "print(\"Begrenzungsrahmen in LAEA:\", shapefile_bounds_laea)" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "print(f\"Koordinatensystem: EPSG({shapefile_epsg})\")\n", "print(f\"Bounding Box: {shapefile_bounds}\")\n", "\n", "# Schnelles Plotten Ihres Shapefiles\n", "gdf.plot()\n", "plt.title(\"Shapefile Darstellung\")\n", "plt.xlabel(\"Längengrad / X\")\n", "plt.ylabel(\"Breitengrad / Y\")\n", "plt.show()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "> Wichtig: Im Allgemeinen ist **kn_bounds_laea** ausreichend für die Extraktion von CSV- und GeoTIFF-Dateien. Allerdings können Anpassungen erforderlich sein, da:\n", "\n", "* **x** und **y** die Mittelpunkte der Rasterzellen darstellen, wodurch die Begrenzungsbox-Kanten möglicherweise nicht genau mit diesen Mittelpunkten übereinstimmen.\n", "* Ohne Anpassung könnten Rasterzellen nahe den Rändern ausgeschlossen werden, wenn die Begrenzungsbox-Kanten zwischen den Zellzentren liegen.\n", "* Eine Anpassung stellt sicher, dass alle relevanten Rasterzellen einbezogen werden.\n", "\n", "
\n", "
\n", " \n", "
Ausdehnung ohne BBox-Anpassung (kn_bounds_laea)
\n", "
\n", "
\n", " \n", "
Ausdehnung mit BBox-Anpassung (kn_bounds_laea_adjusted)
\n", "
\n", "
\n", "\n", "Es wird empfohlen, **kn_bounds_laea_adjusted** für die Extraktion von GeoTIFF-Dateien zu verwenden." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# Definiere den Variablennamen aus den verfügbaren Variablen und lese die Variablendaten ein\n", "variable_name = 'LAPrec1901'\n", "variable_data = nc_dataset[variable_name]\n", "\n", "x = nc_dataset['X'][:]\n", "y = nc_dataset['Y'][:]\n", "\n", "# Validierung der Begrenzungsbox-Dimensionen\n", "dx = x[1] - x[0] # X-Achsen-Auflösung\n", "dy = y[1] - y[0] # Y-Achsen-Auflösung\n", "\n", "# Anpassung der Begrenzungsbox an die Rasterkanten\n", "x_min_adjusted = max(x_min - dx / 2, x.min())\n", "x_max_adjusted = min(x_max + dx / 2, x.max())\n", "y_min_adjusted = max(y_min - dy / 2, y.min())\n", "y_max_adjusted = min(y_max + dy / 2, y.max())\n", "\n", "shp_bounds_laea_adjusted = [x_min_adjusted, y_min_adjusted, x_max_adjusted, y_max_adjusted]\n", "# print(shp_bounds_laea_adjusted)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### 4.2 Daten nach Begrenzungsbox filtern und als CSV exportieren" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "import xarray as xr\n", "\n", "# Funktion zur Konvertierung von NetCDF-Daten in ein Pandas DataFrame\n", "def netcdf_to_dataframe(nc_filepath, bounding_box=None):\n", "\n", " with xr.open_dataset(nc_filepath) as nc_dataset:\n", " # Zugriff auf die Variablendaten aus dem Dataset\n", " variable_data = nc_dataset[variable_name]\n", "\n", " # Falls eine Begrenzungsbox angegeben ist, Daten filtern\n", " if bounding_box:\n", " filtered_data = variable_data.where(\n", " (nc_dataset['X'] >= bounding_box[0]) & (nc_dataset['X'] <= bounding_box[2]) &\n", " (nc_dataset['Y'] >= bounding_box[1]) & (nc_dataset['Y'] <= bounding_box[3]),\n", " drop=True\n", " )\n", " else:\n", " filtered_data = variable_data\n", "\n", " # Konvertierung des xarray-Datasets in ein Pandas DataFrame\n", " df = filtered_data.to_dataframe().reset_index().set_index(['time', 'X', 'Y'])\n", "\n", " return df" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# Definiere den CSV-Dateinamen und den Dateipfad für die Ausgabe\n", "csv_filename = f\"{variable_name}_filtered_data.csv\"\n", "csv_filepath = os.path.join(csv_folder, csv_filename)\n", "\n", "# Exportiere das DataFrame als CSV, falls es noch nicht existiert\n", "if not os.path.isfile(csv_filepath):\n", " dataframe = netcdf_to_dataframe(nc_filepath=nc_filepath, bounding_box=shp_bounds_laea_adjusted)\n", " dataframe.to_csv(csv_filepath, sep=\",\", encoding='utf8')\n", "else:\n", " print(f\"Datei existiert bereits unter {csv_filepath}.\\nExport wird übersprungen.\")\n", " print(\"Lesen bestehende CSV-Datei ein...\")\n", " dataframe = pd.read_csv(csv_filepath).set_index(['time', 'Y', 'X'])\n", "\n", "# Zeige das DataFrame an\n", "dataframe" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### 4.3 Monatliche Mittelwerte berechnen und als CSV exportieren" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# Definiere den CSV-Dateinamen und den Dateipfad für die Ausgabe\n", "csv_filename_monthly_means = f\"{variable_name}_monthly_means.csv\"\n", "csv_filepath_monthly_means = os.path.join(csv_folder, csv_filename_monthly_means)\n", "\n", "# Exportiere die Pivot-Tabelle als CSV, falls sie noch nicht existiert\n", "if not os.path.isfile(csv_filepath_monthly_means):\n", " # Konvertiere die 'time'-Spalte in das Datetime-Format\n", " filtered_df_copy = dataframe.copy().reset_index()\n", " filtered_df_copy['time'] = pd.to_datetime(filtered_df_copy['time'])\n", "\n", " # Extrahiere Jahr und Monat aus der 'time'-Spalte\n", " filtered_df_copy['Year'] = filtered_df_copy['time'].dt.year\n", " filtered_df_copy['Month'] = filtered_df_copy['time'].dt.month\n", "\n", " # Gruppiere nach Jahr und Monat, berechne den Mittelwert für die angegebene Variable\n", " monthly_means = (\n", " filtered_df_copy.groupby(['Year', 'Month'])[variable_name]\n", " .mean()\n", " .reset_index()\n", " .round({variable_name: 2})\n", " )\n", "\n", " # Erstelle eine Pivot-Tabelle mit Jahren als Zeilen und Monaten als Spalten\n", " df_monthly_means = monthly_means.pivot(index='Year', columns='Month', values=variable_name)\n", "\n", " # Stelle sicher, dass alle Monate (1–12) enthalten sind, und füge jährliche Summen hinzu\n", " df_monthly_means = df_monthly_means.reindex(columns=range(1, 13))\n", " df_monthly_means[f\"Yearly_Sum\"] = df_monthly_means.sum(axis=1)\n", " df_monthly_means.to_csv(csv_filepath_monthly_means, index=True)\n", " print(f\"Monatliche Mittelwerte erfolgreich exportiert nach {csv_filepath}\")\n", "else:\n", " print(f\"Datei existiert bereits unter {csv_filepath}. Export wird übersprungen.\")\n", " print(\"Bestehende CSV-Datei wird eingelesen...\")\n", " df_monthly_means = pd.read_csv(csv_filepath_monthly_means).set_index(['Year'])\n", "\n", "# Zeige die Pivot-Tabelle an\n", "df_monthly_means" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### 4.4 Monatliche Durchschnittswerte für einen definierten Zeitraum darstellen" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "import matplotlib.pyplot as plt\n", "import matplotlib.ticker as ticker\n", "\n", "def plot_monthly_averages(\n", " start_year=None,\n", " end_year=None):\n", " \n", " \"\"\"\n", " Stellt die monatlichen Durchschnittswerte einer Variablen mit\n", " Fehlerbalken für einen bestimmten Zeitraum dar.\n", "\n", " Parameter:\n", " start_year (int): Startjahr für den Datensatz (optional).\n", " end_year (int): Endjahr für den Datensatz (optional).\n", " \"\"\"\n", " # Gefiltertes DataFrame aus der CSV-Datei einlesen\n", " filtered_df = pd.read_csv(csv_filepath)\n", " filtered_df['time'] = pd.to_datetime(filtered_df['time'])\n", "\n", " # Bestimmen Sie den Zeitraum\n", " min_year = min(filtered_df['time'].dt.year)\n", " max_year = max(filtered_df['time'].dt.year)\n", " \n", " # Start- und Endjahr an den verfügbaren Datenbereich anpassen\n", " if start_year and end_year:\n", " if start_year < min(filtered_df[\"time\"].dt.year):\n", " print(f\"Das angegebene Startjahr {start_year} liegt vor dem Datenbereich.\")\n", " print(f\"Startjahr wird auf {min_year} angepasst.\")\n", "\n", " if end_year > max(filtered_df[\"time\"].dt.year):\n", " print(f\"Das angegebene Endjahr {end_year} liegt nach dem Datenbereich.\")\n", " print(f\"Endjahr wird auf {max_year} angepasst.\")\n", "\n", " start_year = max(start_year, min_year)\n", " end_year = min(end_year, max_year)\n", " else:\n", " start_year = min_year\n", " end_year = max_year\n", " \n", " # Daten für den definierten Zeitraum filtern\n", " df_period = filtered_df[\n", " (filtered_df[\"time\"].dt.year >= start_year)\n", " & (filtered_df[\"time\"].dt.year <= end_year)\n", " ]\n", "\n", " # Monatliche Statistik berechnen: Mittelwert und Standardabweichung\n", " monthly_stats = df_period.groupby(df_period['time'].dt.month)[variable_name].agg(['mean', 'std'])\n", "\n", " # Diagramm erstellen\n", " fig, ax = plt.subplots(figsize=(12, 6), facecolor='#f1f1f1', edgecolor='k')\n", "\n", " # Balkendiagramm mit Fehlerbalken zeichnen\n", " ax.bar(\n", " monthly_stats.index,\n", " monthly_stats['mean'],\n", " yerr=monthly_stats['std'],\n", " capsize=5,\n", " color='skyblue',\n", " alpha=0.7,\n", " error_kw=dict(ecolor='black', lw=0.75),\n", " )\n", "\n", " # Dynamischen Y-Achsenbereich setzen\n", " y_min = max(0, monthly_stats['mean'].min() - monthly_stats['std'].max() - 0.5)\n", " y_min = y_min // 20 * 20\n", " y_max = monthly_stats['mean'].max() + monthly_stats['std'].max() + 0.5\n", " y_max = (y_max + 20) // 20 * 20\n", " ax.set_ylim(y_min, y_max)\n", "\n", " # Hilfslinien hinzufügen\n", " ax.grid(visible=True, color='#b0b0b0', linestyle='--', linewidth=0.8, alpha=0.6)\n", " ax.yaxis.set_major_formatter(ticker.FormatStrFormatter('%0.2f'))\n", "\n", " # Achsenbeschriftungen und Titel anpassen\n", " ax.set_xlabel('Monat', fontsize=12)\n", " ax.set_ylabel('Durchschnittlicher Niederschlag (mm)', fontsize=12)\n", " ax.set_title(\n", " f'Durchschnittlicher Niederschlag pro Monat (Zeitraum Jan {start_year} - Dez {end_year})',\n", " fontsize=14,\n", " fontweight='bold'\n", " )\n", "\n", " # X-Achse mit Monatsnamen beschriften\n", " month_labels = ['Jan', 'Feb', 'Mär', 'Apr', 'Mai', 'Jun', 'Jul', 'Aug', 'Sep', 'Okt', 'Nov', 'Dez']\n", " ax.set_xticks(range(1, 13))\n", " ax.set_xticklabels(month_labels, rotation=0)\n", "\n", " # Layout anpassen und Diagramm anzeigen\n", " plt.tight_layout()\n", " plt.show()" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# Geben Sie den Jahresbereich an, um den monatlichen Durchschnitt zu filtern und darzustellen.\n", "# Wenn kein Bereich angegeben wird, wird der komplette Datensatz verwendet.\n", "if __name__ == \"__main__\":\n", " plot_monthly_averages() # Plot mit vollständigem Datensatz\n", " plot_monthly_averages(2000,2020) # Plot mit angegebenem Jahresbereich" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## 5. Exportieren der NetCDF4-Dateien nach GeoTIFF" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### 5.1 Definieren eine Funktion zum Exportieren der NetCDF4-Datei als GeoTIFF-Datei(en)" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "from rasterio.transform import from_origin\n", "import rasterio\n", "\n", "from tqdm.notebook import tqdm\n", "\n", "def main_export_geotiff(\n", " nc_filepath,\n", " nc_file_variable,\n", " bounding_box=None,\n", " start_year=None,\n", " end_year=None,\n", " merged=None):\n", " \"\"\"\n", " Exportieren Sie Niederschlagsdaten als GeoTIFF-Dateien,\n", " entweder zusammengeführt in einer einzigen Multiband-Datei \n", " oder als Einzeldateien für jede Zeitscheibe, mit optionaler räumlicher Filterung.\n", " \n", " Parameters:\n", " nc_filepath (str): Dateipfad der NetCDF4-Datei.\n", " nc_file_variable (str): Variablenname der NetCDF4-Datei.\n", " bounding_box (list): [lon_min, lat_min, lon_max, lat_max] (optional).\n", " start_year (int): Startjahr für das Dataset (optional).\n", " end_year (int): Endjahr für das Dataset (optional).\n", " merged (bool): Gibt an, ob ein zusammengeführtes GeoTIFF oder einzelne GeoTIFFs erstellt werden sollen (optional).\n", " \"\"\"\n", "\n", " # Öffnet die NetCDF-Datei\n", " with nc.Dataset(nc_filepath, mode='r') as nc_dataset:\n", " x = nc_dataset['X'][:]\n", " y = nc_dataset['Y'][:]\n", "\n", " # Falls eine Begrenzungsbox angegeben wurde, filtere die Daten entsprechend\n", " if bounding_box:\n", " indices_x = np.where((x >= bounding_box[0]) & (x <= bounding_box[2]))[0]\n", " indices_y = np.where((y >= bounding_box[1]) & (y <= bounding_box[3]))[0]\n", " start_x, end_x = indices_x[0], indices_x[-1] + 1\n", " start_y, end_y = indices_y[0], indices_y[-1] + 1\n", " else:\n", " start_x, end_x = 0, len(x)\n", " start_y, end_y = 0, len(y)\n", " \n", " x = x[start_x:end_x]\n", " y = y[start_y:end_y]\n", "\n", " # Extrahiere die Zeitvariable und konvertiere sie in lesbare Datumsangaben\n", " time_var = nc_dataset.variables['time']\n", " time_units = nc_dataset.variables['time'].units\n", " time_calendar = getattr(time_var, \"calendar\", \"standard\")\n", " cftime = nc.num2date(time_var[:], units=time_units, calendar=time_calendar)\n", " \n", " # Berechnet die räumliche Auflösung und die Rastertransformation\n", " dx = abs(x[1] - x[0])\n", " dy = abs(y[1] - y[0])\n", " transform = from_origin(x.min() - dx / 2, y.min() - dy / 2, dx, -dy)\n", "\n", " # Bestimmt den verfügbaren Zeitraum\n", " min_year = cftime[0].year\n", " max_year = cftime[-1].year\n", "\n", " if start_year and end_year:\n", " # Passen Sie start_year und end_year basierend auf dem verfügbaren Zeitraum an\n", " if start_year < min_year:\n", " print(f\"Startjahr {start_year} liegt vor dem Datensatzbereich. Anpassen an {min_year}.\")\n", " if end_year > max_year:\n", " print(f\"Endjahr {end_year} überschreitet den Datensatzbereich. Anpassung an {max_year}.\")\n", " \n", " start_year = max(start_year, min_year)\n", " end_year = min(end_year, max_year)\n", " else:\n", " # Standardmäßig voller Datensatzbereich\n", " start_year, end_year = min_year, max_year\n", "\n", " # Findet die Indizes, die dem angegebenen Jahresbereich entsprechen\n", " start_index = next(i for i, dt in enumerate(cftime) if dt.year == start_year)\n", " end_index = next(i for i, dt in enumerate(cftime) if dt.year == end_year) + 12 # Ganzes Jahr (monatliche Daten)\n", "\n", " # Extrahiert die Niederschlagsdaten für den ausgewählten Bereich.\n", " precipitation_data = nc_dataset[variable_name]\n", " precipitation_data_subset = precipitation_data[start_index:end_index,start_y:end_y, start_x:end_x]\n", " \n", " if merged:\n", " # Erstellt ein zusammengeführtes GeoTIFF mit allen Zeitscheiben als separate Bänder\n", " output_folder = os.path.join(geotiff_folder, \"merged_geotiff\")\n", " os.makedirs(output_folder, exist_ok=True)\n", "\n", " output_filename = f\"{variable_name}_merged_{start_year}-{end_year}.tif\"\n", " output_filepath = os.path.join(output_folder, output_filename)\n", " \n", " # Erstellt eine GeoTIFF-Datei mit mehreren Bändern für jede Zeitscheibe\n", " with rasterio.open(\n", " output_filepath,\n", " \"w\",\n", " driver = \"GTiff\",\n", " dtype = str(precipitation_data_subset.dtype),\n", " width = precipitation_data_subset.shape[2],\n", " height = precipitation_data_subset.shape[1],\n", " count = precipitation_data_subset.shape[0],\n", " crs = \"EPSG:3035\",\n", " nodata = -9999,\n", " transform = transform, \n", " ) as dst:\n", " for year_index in tqdm(range(precipitation_data_subset.shape[0]),\n", " desc=f\"Exportiere zusammengeführte GeoTIFF-Datei von {start_year} bis {end_year}\"):\n", " band_data = precipitation_data_subset[year_index, :, :]\n", " dt = cftime[start_index + year_index]\n", " band_desc = f\"{dt.year:04d}-{dt.month:02d}-{dt.day:02d}\"\n", " \n", " # Schreibe jede Jahresscheibe als Band\n", " dst.write(band_data, year_index + 1)\n", " dst.set_band_description(year_index + 1, band_desc)\n", "\n", " else:\n", " # Export als einzelne GeoTIFF-Dateien\n", " output_folder = os.path.join(geotiff_folder, \"individual_geotiff\")\n", " os.makedirs(output_folder, exist_ok=True)\n", "\n", " for year_index in tqdm(range(precipitation_data_subset.shape[0]),\n", " desc=f\"Exportieren einzelner GeoTIFF-Dateien von {start_year} bis {end_year}\"):\n", " # Bestimmt das Datum für die aktuelle Zeitscheibe\n", " dt = cftime[start_index + year_index]\n", " dt_formatted = f\"{dt.year:04d}-{dt.month:02d}-{dt.day:02d}\"\n", "\n", " # Definiert den Speicherort der Ausgabe-GeoTIFF-Datei \n", " output_filename = f\"{variable_name}_{dt_formatted}.tif\"\n", " output_filepath = os.path.join(output_folder, output_filename)\n", "\n", " # Exportiert die aktuelle Zeitscheibe als GeoTIFF\n", " with rasterio.open(\n", " output_filepath,\n", " \"w\",\n", " driver=\"GTiff\",\n", " dtype=str(precipitation_data_subset.dtype),\n", " width=precipitation_data_subset.shape[2],\n", " height=precipitation_data_subset.shape[1],\n", " count=1,\n", " crs=\"EPSG:3035\",\n", " nodata=-9999,\n", " transform=transform,\n", " ) as dst:\n", " year_precipitation_data = precipitation_data_subset[year_index, :, :]\n", " dst.write(year_precipitation_data, 1)\n", " dst.set_band_description(1, f\"{dt.year:04d}-{dt.month:02d}-{dt.day:02d}\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### 5.2 Export selected NetCDF4 file(s) to GeoTIFF file(s)" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# Dateipfade und Variablennamen definieren\n", "nc_filepath = os.path.join(extract_folder, 'LAPrec1901.v1.2.nc')\n", "variable_name = 'LAPrec1901'\n", "\n", "if __name__ == \"__main__\":\n", " # Exportieren Sie alle NetCDF-Dateien als zusammengeführte GeoTIFF-Dateien.\n", " main_export_geotiff(\n", " nc_filepath = nc_filepath,\n", " nc_file_variable = variable_name,\n", " bounding_box = None,\n", " start_year = None,\n", " end_year = None,\n", " merged = True\n", " )\n", " \n", " # Alle NetCDF-Dateien als einzelne GeoTIFF-Dateien exportieren\n", " main_export_geotiff(\n", " nc_filepath = nc_filepath,\n", " nc_file_variable = variable_name,\n", " bounding_box = None,\n", " start_year = 1995,\n", " end_year = 2000,\n", " merged = False,\n", " )" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## 6. Analyse und Visualisierung Optionen" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### 6.1 Definieren Sie eine Funktion zur Visualisierung eines NetCDF4-Datensatzes als Heatmap" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "import matplotlib.pyplot as plt\n", "import cartopy.feature as cfeature\n", "import cartopy.crs as ccrs\n", "import numpy as np\n", "import geopandas as gpd\n", "\n", "def main_plt_plot(selected_year, selected_month, save):\n", " \n", " \"\"\"\n", " Plottet monatliche Niederschlagsdaten für ein bestimmtes Jahr und Monat.\n", "\n", " Parameter:\n", " selected_year (str): Das Jahr, das geplottet werden soll (z.B. '2020').\n", " selected_month (str): Der Monat, der geplottet werden soll (z.B. '08').\n", " save (bool): Gibt an, ob der Plot als PNG-Datei gespeichert werden soll.\n", " \"\"\"\n", " \n", " # Dateipfade und Variablennamen definieren\n", " nc_filepath = os.path.join(extract_folder, 'LAPrec1901.v1.2.nc')\n", " variable_name = 'LAPrec1901'\n", " \n", " # Öffnen des NetCDF-Datensatzes\n", " nc_dataset = nc.Dataset(nc_filepath, mode='r')\n", " lon = nc_dataset['lon'][:]\n", " lat = nc_dataset['lat'][:]\n", "\n", " # Zeitvariable extrahieren und in lesbare Daten umwandeln\n", " time_var = nc_dataset.variables['time']\n", " time_units = nc_dataset.variables['time'].units\n", " time_calendar = getattr(time_var, \"calendar\", \"standard\")\n", " cftime = nc.num2date(time_var[:], units=time_units, calendar=time_calendar)\n", "\n", " # Daten für das ausgewählte Datum extrahieren\n", " selected_date = f\"{selected_year}-{selected_month}-01\"\n", " target_index = np.where([t.strftime('%Y-%m-%d') == selected_date for t in cftime])[0][0]\n", " index_data = nc_dataset[variable_name][target_index, :, :]\n", " \n", " # Berechnung von vmin, vmax und bins\n", " vmin = np.floor(np.nanmin(index_data) / 25) * 25\n", " vmax = np.ceil(np.nanmax(index_data) / 25) * 25\n", " interval = 25\n", " bins = int((vmax - vmin) / interval)\n", "\n", " # Erstellen des Plots\n", " fig, ax = plt.subplots(figsize=(12, 8),\n", " facecolor='#f1f1f1',\n", " edgecolor='k',\n", " subplot_kw={'projection': ccrs.PlateCarree()})\n", "\n", " cmap = plt.get_cmap(\"rainbow_r\", bins)\n", " pcm = ax.pcolormesh(lon,\n", " lat,\n", " index_data,\n", " cmap=cmap,\n", " shading=\"auto\",\n", " vmin=vmin,\n", " vmax=vmax)\n", "\n", " # Hinzufügen der Shapefile\n", " konstanz_shp = r\".\\shapefiles\\kn_boundary.shp\"\n", " konstanz_boundary = gpd.read_file(konstanz_shp)\n", " konstanz_boundary = konstanz_boundary.to_crs(4326)\n", " konstanz_boundary.boundary.plot(ax=ax, edgecolor='red', linewidth=1.2)\n", "\n", " # Hinzufügen einer Farbskala\n", " ticks = np.linspace(vmin, vmax, num=bins + 1)\n", " cbar = plt.colorbar(pcm, ax=ax, orientation='horizontal', pad=0.06, shrink=0.8, ticks=ticks)\n", " cbar.set_label(\"Monthly Precipitation [mm]\", fontsize=12)\n", " cbar.ax.tick_params(labelsize=12)\n", "\n", " # Hinzufügen von Kartenmerkmalen\n", " ax.add_feature(cfeature.BORDERS, linestyle='-', linewidth=0.75, edgecolor='black')\n", " ax.add_feature(cfeature.COASTLINE, linewidth=0.75)\n", " ax.add_feature(cfeature.LAND, edgecolor='black', facecolor='none')\n", " # ax.add_feature(cfeature.RIVERS, linewidth=0.5, edgecolor='blue')\n", "\n", " # Hinzufügen von Gitterlinien\n", " gl = ax.gridlines(draw_labels=True,\n", " crs=ccrs.PlateCarree(),\n", " linewidth=0.5,\n", " color='gray',\n", " alpha=0.7,\n", " linestyle='--')\n", " gl.top_labels = False \n", " gl.right_labels = False\n", " gl.xlabel_style = {'size': 10, 'color': 'black'}\n", " gl.ylabel_style = {'size': 10, 'color': 'black'}\n", "\n", " # Hinzufügen von Titel und Achsenbeschriftungen\n", " fig.text(0.5, 0.175, 'Längengrad', ha='center', fontsize=14)\n", " fig.text(-0.01, 0.5, 'Breitengrad', va='center', rotation='vertical', fontsize=14)\n", " ax.set_title(f\"Monatlicher Niederschlag, {selected_date}\", fontsize=14, fontweight='bold')\n", " ax.set_aspect(\"equal\")\n", "\n", " # Plot speichern oder anzeigen\n", " plt.tight_layout()\n", " if save:\n", " temp_images_folder = os.path.join(output_folder, \"temp_images_folder\")\n", " os.makedirs(temp_images_folder, exist_ok=True)\n", " image_filepath = os.path.join(temp_images_folder, f\"{selected_year}-{selected_month}.png\")\n", " plt.savefig(image_filepath, format='png', bbox_inches='tight')\n", " print(f\"Plot gespeichert unter: {image_filepath}\")\n", " plt.show()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### 6.2 Erstellen einer Heatmap für den Monatsmittelwert" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# Beispielhafte Verwendung\n", "if __name__ == \"__main__\":\n", " main_plt_plot(\n", " selected_year='2020',\n", " selected_month='08',\n", " save=False,\n", " )" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### 6.3 Definieren Sie eine Funktion zur Visualisierung eines NetCDF4-Datensatzes als monatliche Heatmap" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "import seaborn as sns\n", "\n", "def plot_precipitation_heatmap(start_year=2010,\n", " end_year=None):\n", " \"\"\"\n", " Erstellt ein Heatmap der monatlichen durchschnittlichen Niederschläge aus einer CSV-Datei.\n", "\n", " Parameter:\n", " start_year (int, optional): Das Startjahr zur Filterung der Daten. Standardmäßig 2010.\n", " end_year (int, optional): Das Endjahr zur Filterung der Daten. Wenn nicht angegeben, wird bis zum maximalen Jahr gefiltert.\n", " \"\"\"\n", " \n", " # Daten laden und vorverarbeiten\n", " csv_path = os.path.join(csv_folder, 'LAPrec1901_monthly_means.csv')\n", " df = pd.read_csv(csv_path)\n", " df = df.drop(columns=['Yearly_Sum']) # Jahrensumme entfernen\n", " df = df.set_index('Year') # Jahr als Index setzen\n", " \n", " # Bestimmen des effektiven Endjahres\n", " max_year_in_data = df.index.max()\n", " if end_year is not None:\n", " end_year = min(max_year_in_data, end_year)\n", " else:\n", " end_year = max_year_in_data\n", " \n", " # Daten basierend auf dem Start- und Endjahr filtern\n", " df_filtered = df.loc[start_year:end_year]\n", " \n", " # Erstellen des Heatmaps\n", " fig, ax = plt.subplots(figsize=(12, 8))\n", " cmap = plt.get_cmap('Blues', 11)\n", " \n", " sns.heatmap(\n", " df_filtered,\n", " cmap=cmap,\n", " annot=True,\n", " annot_kws={\"fontsize\": 10},\n", " fmt=\".0f\",\n", " cbar_kws={\n", " \"label\": f\"Niederschlag ({summary.get('Units', 'mm')})\",\n", " \"ticks\": range(0, 221, 20),\n", " \"pad\": 0.015,\n", " },\n", " vmin=0,\n", " vmax=220\n", " )\n", " \n", " # X-Achse mit Monatsnamen anpassen\n", " month_labels = ['Jan', 'Feb', 'Mär', 'Apr', 'Mai', 'Jun', 'Jul', 'Aug', 'Sep', 'Okt', 'Nov', 'Dez'] \n", " ticks = range(1, 13)\n", " shifted_ticks = [tick - 0.5 for tick in ticks]\n", " # Setzen der Ticks und Labels\n", " ax.set_xticks(shifted_ticks)\n", " ax.set_xticklabels(month_labels, rotation=0, fontsize=10)\n", " \n", " # Hinzufügen von Beschriftungen und Titel\n", " ax.set_title(\"Monatlicher Durchschnittsniederschlag\", fontsize=16, fontweight='bold')\n", " ax.set_xlabel(\"Monat\", fontsize=12)\n", " ax.set_ylabel(\"Jahr\", fontsize=12)\n", "\n", " # Layout anpassen und Plot anzeigen\n", " plt.tight_layout()\n", " plt.show()" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# Beispielhafte Verwendung\n", "if __name__ == \"__main__\":\n", " plot_precipitation_heatmap(\n", " start_year=2010,\n", " end_year=2020\n", " )\n", " \n", " plot_precipitation_heatmap(\n", " start_year=2005,\n", " )" ] } ], "metadata": { "kernelspec": { "display_name": "cds_env", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.10.16" } }, "nbformat": 4, "nbformat_minor": 2 }