Notebooks from applied-cs/data-science@6dcc52e7

71b153b1 · Christof Kaufmann · b0999122 · 71b153b1 · 71b153b1 · 71b153b1
Commit 71b153b1 authored Jul 5, 2024 by Christof Kaufmann
--- a/10-geodaten/01-projections.ipynb
+++ b/10-geodaten/01-projections.ipynb
@@ -25,21 +25,28 @@
   },
   "outputs": [],
   "source": [
-    "from sklearn.datasets import fetch_openml\n",
+    "import numpy as np\n",
+    "import os\n",
    "import geopandas as gpd\n",
+    "from shapely import affinity\n",
+    "from sklearn.datasets import fetch_openml\n",
    "\n",
+    "path = 'geo_countries.geojson'\n",
+    "if not os.path.exists(path):\n",
    "    url = 'https://d2ad6b4ur7yvpq.cloudfront.net/naturalearth-3.3.0/ne_50m_admin_0_countries.geojson'\n",
    "    geo_countries = gpd.read_file(url)\n",
-    "path = 'geo_countries.geojson'\n",
-    "geo_countries.to_file(path)\n",
+    "    geo_countries.to_file(path)      # save local copy\n",
+    "\n",
+    "geo_countries = gpd.read_file(path)  # load from local copy\n",
    "geo_countries = geo_countries[['name', 'iso_a3', 'geometry']].set_index('iso_a3')\n",
    "\n",
    "countries_df, _ = fetch_openml('Countries-of-the-World', as_frame=True, return_X_y=True)\n",
    "countries_df = countries_df.rename(columns={'Area_(sq._mi.)': 'ref_area', 'Country': 'name'})  # unit was km² already\n",
    "countries_df = countries_df[['name', 'ref_area']]\n",
-    "countries_df['name'] = countries_df.name.str.strip()"
+    "countries_df['name'] = countries_df.name.str.strip()\n",
+    "countries_df.loc[len(countries_df)] = ['Antarctica', 14_200_000]  # from: https://en.wikipedia.org/wiki/Antarctica"
   ],
-   "id": "0003-4c6fa2e8089a613e6a3300cf4bd6fe01d88ed16280b83e1e378065b75a9"
+   "id": "0003-1224f13613ae889be0984402ff1578d533b9774aae88ca6305a774f6860"
  },
  {
   "cell_type": "markdown",
@@ -60,20 +67,63 @@
   "source": [
    "countries = geo_countries.merge(countries_df, on='name')\n",
    "\n",
+    "scaledown = False\n",
    "# epsg = 25832  # UTM Zone 32N, good for germany, does not work on global scale\n",
    "# epsg = 4087\n",
    "# epsg = 3395\n",
-    "epsg = 3857\n",
    "# epsg = 4326\n",
    "# ...\n",
    "\n",
+    "# EPSG:3857 is the original projection CRS of the background map and thus the map corresponds to the computed projection error. So in this case we can scale geometries down to show true sizes of countries.\n",
+    "epsg, scaledown = 3857, True\n",
+    "\n",
    "countries = countries.to_crs(f'EPSG:{epsg}')\n",
+    "\n",
+    "# calculate relative difference of projected and reference areas\n",
    "countries['proj_area'] = countries.area / 1e6  # in km²\n",
-    "countries['rel_diff_area'] = ((countries['proj_area'] / countries['ref_area'] - 1).round(2) * 100).astype(int)  # in percent\n",
+    "ratio = countries['proj_area'] / countries['ref_area']\n",
+    "countries['rel_diff_area'] = ((ratio - 1).round(2) * 100).astype(int)  # in percent\n",
+    "\n",
+    "if scaledown:\n",
+    "    # split MultiPolygons into separate Polygons, e. g. 1 MultiPolygon of Germany → 6 Polygons (6 rows with the same data)\n",
+    "    countries = countries.explode()\n",
+    "\n",
+    "    # scale down each poly\n",
+    "    scale = 1 / np.sqrt(ratio)\n",
+    "    countries.geometry = [affinity.scale(g, s, s) for g, s in zip(countries.geometry, scale[countries.index])]\n",
+    "\n",
+    "    # special handling for Antarctica, since it would be below the visible area\n",
+    "    max_antarctica_area = countries.loc[countries.name == 'Antarctica'].area.max()\n",
+    "    largest_antarctica_poly = countries.area == max_antarctica_area\n",
+    "    shifted_antarctica = affinity.translate(countries.loc[largest_antarctica_poly, 'geometry'].squeeze(), yoff=np.sqrt(max_antarctica_area)*7)\n",
+    "    countries.loc[largest_antarctica_poly, 'geometry'] = shifted_antarctica\n",
+    "\n",
+    "    # reverse splitting, i. e. combine Polygons with the same name into a single MultiPolygon\n",
+    "    countries = countries.dissolve('name')\n",
+    "\n",
+    "# plot\n",
    "vmax = countries['rel_diff_area'].abs().quantile(0.97)\n",
-    "countries.explore(column='rel_diff_area', tooltip=['name', 'rel_diff_area', 'proj_area', 'ref_area'], popup=True, tiles=\"CartoDB voyager\", vmin=-vmax, vmax=vmax, cmap='coolwarm', attribution=countries.crs.name)"
+    "countries.explore(\n",
+    "    column='rel_diff_area',\n",
+    "    vmin=-vmax, vmax=vmax, cmap='coolwarm',\n",
+    "    tooltip=['name', 'rel_diff_area', 'proj_area', 'ref_area'],\n",
+    "    popup=True,\n",
+    "    tiles=\"CartoDB voyager\",\n",
+    "    attribution=countries.crs.name\n",
+    ")"
+   ],
+   "id": "0005-87f4bcf5159c2d37428cabbcd55f47bf4cd95bbb86ed15eced56f1b465b"
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "Siehe auch:\n",
+    "\n",
+    "-   https://truesizeofcountries.com/\n",
+    "-   https://www.thetruesize.com"
   ],
-   "id": "0005-b64aa750773482384137dd8fa7e5b5b069b9c692b3395c5c211650bf153"
+   "id": "0007-c166693366c20bec45d0a68d7bc44664690db85407f20dd66ceb12f7251"
  }
 ],
 "nbformat": 4,

 %% Cell type:markdown id:0002-a8f53080980ed874d387cc3bfb8c5d5ed65d7fdb4139eb4ed15105d383c tags:

 # Flächenverzerrung bei Projektionen

 In dieser Demo kann man sich die Flächenverzerrungen bei verschiedenen
 Projektionen auf einer interaktiven Karte anschauen.

 Zunächst laden wir eine GeoJSON-Datei mit groben Ländergrenzen herunter
 und speichern sie lokal. Und wir laden auch eine mehr oder weniger gute
 Datei mit den Referenzgrößen. Mit etwas Vorverarbeitung können wir die
 Datensätze in `countries` zusammenführen, was wir aber erst im nächsten
 Schritt machen.

-%% Cell type:code id:0003-4c6fa2e8089a613e6a3300cf4bd6fe01d88ed16280b83e1e378065b75a9 tags:
+%% Cell type:code id:0003-1224f13613ae889be0984402ff1578d533b9774aae88ca6305a774f6860 tags:

 ``` 
-from sklearn.datasets import fetch_openml
+import numpy as np
+import os
 import geopandas as gpd
+from shapely import affinity
+from sklearn.datasets import fetch_openml

-url = 'https://d2ad6b4ur7yvpq.cloudfront.net/naturalearth-3.3.0/ne_50m_admin_0_countries.geojson'
-geo_countries = gpd.read_file(url)
 path = 'geo_countries.geojson'
-geo_countries.to_file(path)
+if not os.path.exists(path):
+    url = 'https://d2ad6b4ur7yvpq.cloudfront.net/naturalearth-3.3.0/ne_50m_admin_0_countries.geojson'
+    geo_countries = gpd.read_file(url)
+    geo_countries.to_file(path)      # save local copy
+
+geo_countries = gpd.read_file(path)  # load from local copy
 geo_countries = geo_countries[['name', 'iso_a3', 'geometry']].set_index('iso_a3')

 countries_df, _ = fetch_openml('Countries-of-the-World', as_frame=True, return_X_y=True)
 countries_df = countries_df.rename(columns={'Area_(sq._mi.)': 'ref_area', 'Country': 'name'})  # unit was km² already
 countries_df = countries_df[['name', 'ref_area']]
 countries_df['name'] = countries_df.name.str.strip()
+countries_df.loc[len(countries_df)] = ['Antarctica', 14_200_000]  # from: https://en.wikipedia.org/wiki/Antarctica
 ```

 %% Cell type:markdown id:0004-5231db90b71a4c73fd2da36a3ca092775410c9eb4e21d28e7de691474ef tags:

 Hier kann ein EPSG-Code gewählt werden, der die Projektion definiert.
 Damit bestimmt sich wie groß die projizierte Fläche ist.

-%% Cell type:code id:0005-b64aa750773482384137dd8fa7e5b5b069b9c692b3395c5c211650bf153 tags:
+%% Cell type:code id:0005-87f4bcf5159c2d37428cabbcd55f47bf4cd95bbb86ed15eced56f1b465b tags:

 ``` 
 countries = geo_countries.merge(countries_df, on='name')

+scaledown = False
 # epsg = 25832  # UTM Zone 32N, good for germany, does not work on global scale
 # epsg = 4087
 # epsg = 3395
-epsg = 3857
 # epsg = 4326
 # ...

+# EPSG:3857 is the original projection CRS of the background map and thus the map corresponds to the computed projection error. So in this case we can scale geometries down to show true sizes of countries.
+epsg, scaledown = 3857, True
+
 countries = countries.to_crs(f'EPSG:{epsg}')
+
+# calculate relative difference of projected and reference areas
 countries['proj_area'] = countries.area / 1e6  # in km²
-countries['rel_diff_area'] = ((countries['proj_area'] / countries['ref_area'] - 1).round(2) * 100).astype(int)  # in percent
+ratio = countries['proj_area'] / countries['ref_area']
+countries['rel_diff_area'] = ((ratio - 1).round(2) * 100).astype(int)  # in percent
+
+if scaledown:
+    # split MultiPolygons into separate Polygons, e. g. 1 MultiPolygon of Germany → 6 Polygons (6 rows with the same data)
+    countries = countries.explode()
+
+    # scale down each poly
+    scale = 1 / np.sqrt(ratio)
+    countries.geometry = [affinity.scale(g, s, s) for g, s in zip(countries.geometry, scale[countries.index])]
+
+    # special handling for Antarctica, since it would be below the visible area
+    max_antarctica_area = countries.loc[countries.name == 'Antarctica'].area.max()
+    largest_antarctica_poly = countries.area == max_antarctica_area
+    shifted_antarctica = affinity.translate(countries.loc[largest_antarctica_poly, 'geometry'].squeeze(), yoff=np.sqrt(max_antarctica_area)*7)
+    countries.loc[largest_antarctica_poly, 'geometry'] = shifted_antarctica
+
+    # reverse splitting, i. e. combine Polygons with the same name into a single MultiPolygon
+    countries = countries.dissolve('name')
+
+# plot
 vmax = countries['rel_diff_area'].abs().quantile(0.97)
-countries.explore(column='rel_diff_area', tooltip=['name', 'rel_diff_area', 'proj_area', 'ref_area'], popup=True, tiles="CartoDB voyager", vmin=-vmax, vmax=vmax, cmap='coolwarm', attribution=countries.crs.name)
+countries.explore(
+    column='rel_diff_area',
+    vmin=-vmax, vmax=vmax, cmap='coolwarm',
+    tooltip=['name', 'rel_diff_area', 'proj_area', 'ref_area'],
+    popup=True,
+    tiles="CartoDB voyager",
+    attribution=countries.crs.name
+)
 ```
+
+%% Cell type:markdown id:0007-c166693366c20bec45d0a68d7bc44664690db85407f20dd66ceb12f7251 tags:
+
+Siehe auch:
+
+-   https://truesizeofcountries.com/
+-   https://www.thetruesize.com

--- a/10-geodaten/02-bike-track-sol.ipynb
+++ b/10-geodaten/02-bike-track-sol.ipynb
+{
+ "cells": [
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "# Bike Track\n",
+    "\n",
+    "Es ist ein GPS-Track von einer Fahrradtour gegeben. Dieser, sowie ein\n",
+    "GeoDataFrame von Kreisen und kreisfreien Städten wird schon im\n",
+    "Start-Code geladen. Finden Sie heraus durch welche Kreise und kreisfreie\n",
+    "Städte die Tour lief.\n",
+    "\n",
+    "*Bonus: Stellen Sie die Kreise und kreisfreie Städte, durch die die\n",
+    "Fahrradtour lief, in einer interaktiven Karte da. Ggf. stellen Sie\n",
+    "zusätzlich den GPS-Track dar. Siehe\n",
+    "[hier](https://geopandas.org/en/stable/docs/user_guide/interactive_mapping.html)\n",
+    "um mehrere GeoDataFrames in einer Karte darzustellen.*"
+   ],
+   "id": "0002-96cd7cef045240a7b1b0c57d7e374c0dd9716f4dc6f170b31077d61a462"
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {
+    "style": "python"
+   },
+   "outputs": [],
+   "source": [
+    "import geopandas as gpd\n",
+    "import zstandard\n",
+    "\n",
+    "cities_path = 'kreise-geo-mittel.geojson.zst'\n",
+    "with open(cities_path, 'rb') as fh:\n",
+    "    dctx = zstandard.ZstdDecompressor()\n",
+    "    with dctx.stream_reader(fh) as reader:\n",
+    "        grid = gpd.read_file(reader).set_index('ID_3')\n",
+    "\n",
+    "track_path = 'bike_track_2024-05-25.geojson.zst'\n",
+    "with open(track_path, 'rb') as fh:\n",
+    "    dctx = zstandard.ZstdDecompressor()\n",
+    "    with dctx.stream_reader(fh) as reader:\n",
+    "        track = gpd.read_file(reader)"
+   ],
+   "id": "0003-f5d67afb904d0f49edc541c7a33c18d0985cd0356b8747150b8687a8e18"
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## Lösung\n",
+    "\n",
+    "Mit `sjoin` können wir die Polygone auswählen, in denen die Punkte des\n",
+    "GPS-Track liegen. Dabei wird unter anderem eine Spalte `index_right`\n",
+    "erzeugt, die wir nutzen um die Kreise und kreisfreien Städte\n",
+    "entsprechend des Tracks zu sortieren. Aber eigentlich brauchen wir nur\n",
+    "den Namen `NAME_3` und die Geometrie zum Plotten. Da für jeden Punkt im\n",
+    "GPS-Track eine Zeile hinzugefügt wird entfernen wir die Duplikate."
+   ],
+   "id": "0005-cd644122221c1cbcbcf18caa67556650fe84a458610ba7d11da6109e7b1"
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {
+    "style": "python"
+   },
+   "outputs": [],
+   "source": [
+    "cities = grid.sjoin(track).sort_values('index_right')\n",
+    "cities = cities[['NAME_3', 'geometry']]\n",
+    "cities = cities.drop_duplicates()"
+   ],
+   "id": "0006-2c45647891a1292c7e47ee77fd9ed51784bffeba4a292c7a25555e9e529"
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [
+    {
+     "output_type": "stream",
+     "name": "stdout",
+     "text": [
+      "['Mettmann', 'Ennepe-Ruhr', 'Bochum Städte', 'Dortmund Städte', 'Unna', 'Soest', 'Hochsauerlandkreis', 'Waldeck-Frankenberg', 'Höxter', 'Kassel', 'Kassel Städte']"
+     ]
+    }
+   ],
+   "source": [
+    "cities['NAME_3'].to_list()"
+   ],
+   "id": "0007-4a30b87cb18cc8bb69afce6b37dc80d502773b474ca99b60625b42e2ccf"
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "Für die Bonusaufgabe ist nun nicht mehr viel zu tun. Aber wir geben uns\n",
+    "etwas Mühe mit dem GPS-Track und Plotten die Farbe entsprechend der\n",
+    "Geschwindigkeit und jeden Segmentanfang (beim Anhalten pausiert die\n",
+    "Aufnahme) etwas größer."
+   ],
+   "id": "0008-5ee78ebd0383ea7c29c4bade37d4ab2368b2ed24fb9719f4e7d070eb1d3"
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {
+    "style": "python"
+   },
+   "outputs": [],
+   "source": [
+    "m = cities.explore()\n",
+    "track.explore(\n",
+    "    column='speed',\n",
+    "    style_kwds={\"style_function\": lambda x: {\"radius\": (x[\"properties\"][\"track_seg_point_id\"] == 0) * 3 + 1}},\n",
+    "    m=m\n",
+    ")"
+   ],
+   "id": "0009-bcaada0a444065b3b854e5a440df01ed30acd0e553f55216d4f44cc7e51"
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "Die vorbereiteten Daten zu den Kreisen und kreisfreien Städten stammen\n",
+    "von\n",
+    "[hier](https://github.com/isellsoap/deutschlandGeoJSON/tree/main/4_kreise)\n",
+    "(3_mittel.geo.json).\n",
+    "\n",
+    "Wenn man nicht nur Kreise und kreisfreie Städte haben will, sondern eine\n",
+    "höhere Auflösung, kann man z. B.\n",
+    "[Gemeindegrenzen](https://hub.arcgis.com/datasets/esri-de-content::gemeindegrenzen-2022/about)\n",
+    "verwenden."
+   ],
+   "id": "0011-2a24f7ded33e3832d9f7e06b4143c0b22b4763564f9cc12dec28965e820"
+  }
+ ],
+ "nbformat": 4,
+ "nbformat_minor": 5,
+ "metadata": {}
+}
+%% Cell type:markdown id:0002-96cd7cef045240a7b1b0c57d7e374c0dd9716f4dc6f170b31077d61a462 tags:
+
+# Bike Track
+
+Es ist ein GPS-Track von einer Fahrradtour gegeben. Dieser, sowie ein
+GeoDataFrame von Kreisen und kreisfreien Städten wird schon im
+Start-Code geladen. Finden Sie heraus durch welche Kreise und kreisfreie
+Städte die Tour lief.
+
+*Bonus: Stellen Sie die Kreise und kreisfreie Städte, durch die die
+Fahrradtour lief, in einer interaktiven Karte da. Ggf. stellen Sie
+zusätzlich den GPS-Track dar. Siehe
+[hier](https://geopandas.org/en/stable/docs/user_guide/interactive_mapping.html)
+um mehrere GeoDataFrames in einer Karte darzustellen.*
+
+%% Cell type:code id:0003-f5d67afb904d0f49edc541c7a33c18d0985cd0356b8747150b8687a8e18 tags:
+
+``` 
+import geopandas as gpd
+import zstandard
+
+cities_path = 'kreise-geo-mittel.geojson.zst'
+with open(cities_path, 'rb') as fh:
+    dctx = zstandard.ZstdDecompressor()
+    with dctx.stream_reader(fh) as reader:
+        grid = gpd.read_file(reader).set_index('ID_3')
+
+track_path = 'bike_track_2024-05-25.geojson.zst'
+with open(track_path, 'rb') as fh:
+    dctx = zstandard.ZstdDecompressor()
+    with dctx.stream_reader(fh) as reader:
+        track = gpd.read_file(reader)
+```
+
+%% Cell type:markdown id:0005-cd644122221c1cbcbcf18caa67556650fe84a458610ba7d11da6109e7b1 tags:
+
+## Lösung
+
+Mit `sjoin` können wir die Polygone auswählen, in denen die Punkte des
+GPS-Track liegen. Dabei wird unter anderem eine Spalte `index_right`
+erzeugt, die wir nutzen um die Kreise und kreisfreien Städte
+entsprechend des Tracks zu sortieren. Aber eigentlich brauchen wir nur
+den Namen `NAME_3` und die Geometrie zum Plotten. Da für jeden Punkt im
+GPS-Track eine Zeile hinzugefügt wird entfernen wir die Duplikate.
+
+%% Cell type:code id:0006-2c45647891a1292c7e47ee77fd9ed51784bffeba4a292c7a25555e9e529 tags:
+
+``` 
+cities = grid.sjoin(track).sort_values('index_right')
+cities = cities[['NAME_3', 'geometry']]
+cities = cities.drop_duplicates()
+```
+
+%% Cell type:code id:0007-4a30b87cb18cc8bb69afce6b37dc80d502773b474ca99b60625b42e2ccf tags:
+
+``` 
+cities['NAME_3'].to_list()
+```
+
+%% Output
+
+    ['Mettmann', 'Ennepe-Ruhr', 'Bochum Städte', 'Dortmund Städte', 'Unna', 'Soest', 'Hochsauerlandkreis', 'Waldeck-Frankenberg', 'Höxter', 'Kassel', 'Kassel Städte']
+
+%% Cell type:markdown id:0008-5ee78ebd0383ea7c29c4bade37d4ab2368b2ed24fb9719f4e7d070eb1d3 tags:
+
+Für die Bonusaufgabe ist nun nicht mehr viel zu tun. Aber wir geben uns
+etwas Mühe mit dem GPS-Track und Plotten die Farbe entsprechend der
+Geschwindigkeit und jeden Segmentanfang (beim Anhalten pausiert die
+Aufnahme) etwas größer.
+
+%% Cell type:code id:0009-bcaada0a444065b3b854e5a440df01ed30acd0e553f55216d4f44cc7e51 tags:
+
+``` 
+m = cities.explore()
+track.explore(
+    column='speed',
+    style_kwds={"style_function": lambda x: {"radius": (x["properties"]["track_seg_point_id"] == 0) * 3 + 1}},
+    m=m
+)
+```
+
+%% Cell type:markdown id:0011-2a24f7ded33e3832d9f7e06b4143c0b22b4763564f9cc12dec28965e820 tags:
+
+Die vorbereiteten Daten zu den Kreisen und kreisfreien Städten stammen
+von
+[hier](https://github.com/isellsoap/deutschlandGeoJSON/tree/main/4_kreise)
+(3_mittel.geo.json).
+
+Wenn man nicht nur Kreise und kreisfreie Städte haben will, sondern eine
+höhere Auflösung, kann man z. B.
+[Gemeindegrenzen](https://hub.arcgis.com/datasets/esri-de-content::gemeindegrenzen-2022/about)
+verwenden.
--- a/10-geodaten/02-bike-track.ipynb
+++ b/10-geodaten/02-bike-track.ipynb
+{
+ "cells": [
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "# Bike Track\n",
+    "\n",
+    "Es ist ein GPS-Track von einer Fahrradtour gegeben. Dieser, sowie ein\n",
+    "GeoDataFrame von Kreisen und kreisfreien Städten wird schon im\n",
+    "Start-Code geladen. Finden Sie heraus durch welche Kreise und kreisfreie\n",
+    "Städte die Tour lief.\n",
+    "\n",
+    "*Bonus: Stellen Sie die Kreise und kreisfreie Städte, durch die die\n",
+    "Fahrradtour lief, in einer interaktiven Karte da. Ggf. stellen Sie\n",
+    "zusätzlich den GPS-Track dar. Siehe\n",
+    "[hier](https://geopandas.org/en/stable/docs/user_guide/interactive_mapping.html)\n",
+    "um mehrere GeoDataFrames in einer Karte darzustellen.*"
+   ],
+   "id": "0002-96cd7cef045240a7b1b0c57d7e374c0dd9716f4dc6f170b31077d61a462"
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {
+    "style": "python"
+   },
+   "outputs": [],
+   "source": [
+    "import geopandas as gpd\n",
+    "import zstandard\n",
+    "\n",
+    "cities_path = 'kreise-geo-mittel.geojson.zst'\n",
+    "with open(cities_path, 'rb') as fh:\n",
+    "    dctx = zstandard.ZstdDecompressor()\n",
+    "    with dctx.stream_reader(fh) as reader:\n",
+    "        grid = gpd.read_file(reader).set_index('ID_3')\n",
+    "\n",
+    "track_path = 'bike_track_2024-05-25.geojson.zst'\n",
+    "with open(track_path, 'rb') as fh:\n",
+    "    dctx = zstandard.ZstdDecompressor()\n",
+    "    with dctx.stream_reader(fh) as reader:\n",
+    "        track = gpd.read_file(reader)"
+   ],
+   "id": "0003-f5d67afb904d0f49edc541c7a33c18d0985cd0356b8747150b8687a8e18"
+  }
+ ],
+ "nbformat": 4,
+ "nbformat_minor": 5,
+ "metadata": {}
+}
+%% Cell type:markdown id:0002-96cd7cef045240a7b1b0c57d7e374c0dd9716f4dc6f170b31077d61a462 tags:
+
+# Bike Track
+
+Es ist ein GPS-Track von einer Fahrradtour gegeben. Dieser, sowie ein
+GeoDataFrame von Kreisen und kreisfreien Städten wird schon im
+Start-Code geladen. Finden Sie heraus durch welche Kreise und kreisfreie
+Städte die Tour lief.
+
+*Bonus: Stellen Sie die Kreise und kreisfreie Städte, durch die die
+Fahrradtour lief, in einer interaktiven Karte da. Ggf. stellen Sie
+zusätzlich den GPS-Track dar. Siehe
+[hier](https://geopandas.org/en/stable/docs/user_guide/interactive_mapping.html)
+um mehrere GeoDataFrames in einer Karte darzustellen.*
+
+%% Cell type:code id:0003-f5d67afb904d0f49edc541c7a33c18d0985cd0356b8747150b8687a8e18 tags:
+
+``` 
+import geopandas as gpd
+import zstandard
+
+cities_path = 'kreise-geo-mittel.geojson.zst'
+with open(cities_path, 'rb') as fh:
+    dctx = zstandard.ZstdDecompressor()
+    with dctx.stream_reader(fh) as reader:
+        grid = gpd.read_file(reader).set_index('ID_3')
+
+track_path = 'bike_track_2024-05-25.geojson.zst'
+with open(track_path, 'rb') as fh:
+    dctx = zstandard.ZstdDecompressor()
+    with dctx.stream_reader(fh) as reader:
+        track = gpd.read_file(reader)
+```
--- a/10-geodaten/03-grid-transform.ipynb
+++ b/10-geodaten/03-grid-transform.ipynb
+{
+ "cells": [
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "# Gitter-Transformation\n",
+    "\n",
+    "Installieren Sie ggf. fehlende Pakete:"
+   ],
+   "id": "0001-9d10715da27eea0563e559964a2d386bf45e78b85abc9cb22a230f2e35c"
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {
+    "style": "python"
+   },
+   "outputs": [],
+   "source": [
+    "!pip install wetterdienst"
+   ],
+   "id": "0002-2b71ddfdefd192f1b26f0e96bf6da9ba9144a57b5b0e0e4e81bbe24bb52"
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "1.  Laden Sie die Wetterdaten mit dem vorbereiteten Code und schauen Sie\n",
+    "    sich die Daten grob an.\n",
+    "2.  Laden Sie das TK 25 Gitter (schon vorbereitet) und die Kreise (und\n",
+    "    kreisfreie Städte) analog.\n",
+    "3.  Interpolieren Sie die Temperaturen mit $k$-NN in einem der Gitter.\n",
+    "    Als Koordinaten können Sie die `centroid`-Koordinaten verwenden.\n",
+    "4.  Transformieren Sie die Daten in das andere Gitter."
+   ],
+   "id": "0003-d9bedf792849d3240f20100ae2d5727c2843229862be9f3803725199e01"
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {
+    "style": "python"
+   },
+   "outputs": [],
+   "source": [
+    "import geopandas as gpd\n",
+    "import pandas as pd\n",
+    "\n",
+    "import zstandard\n",
+    "from wetterdienst.provider.dwd.observation import (DwdObservationRequest,\n",
+    "                                                   DwdObservationResolution)\n",
+    "\n",
+    "\n",
+    "cache_path = '.'\n",
+    "\n",
+    "tk25_path = 'tk25-grid-utm32s.geojson.zst'\n",
+    "with open(tk25_path, 'rb') as fh:\n",
+    "    dctx = zstandard.ZstdDecompressor()\n",
+    "    with dctx.stream_reader(fh) as reader:\n",
+    "        grid = gpd.read_file(reader).set_index('id')\n",
+    "\n",
+    "\n",
+    "def retrieve_weatherdata(start_date=\"2005-01-01\", end_date=\"2023-08-31\", drop_cache=False):\n",
+    "    try:\n",
+    "        if drop_cache:\n",
+    "            raise FileNotFoundError('Dropping cache')\n",
+    "\n",
+    "        dwd_stations = pd.read_parquet(cache_path + f'/dwd_stations-{end_date}.parquet.zst')\n",
+    "        dwd_data = pd.read_parquet(cache_path + f'/dwd_data-{end_date}.parquet.zst')\n",
+    "    except FileNotFoundError:\n",
+    "        all_station_request = DwdObservationRequest(\n",
+    "            parameter=[\"temperature_air_mean_200\",  # monthly mean temperature 2 m above ground\n",
+    "                       \"precipitation_height\"],     # monthly sum of precipitation height\n",
+    "            resolution=DwdObservationResolution.MONTHLY,\n",
+    "            start_date=start_date,\n",
+    "            end_date=end_date,\n",
+    "        )\n",
+    "\n",
+    "        dwd_stations = all_station_request.all().df\n",
+    "        if not isinstance(dwd_stations, pd.DataFrame):\n",
+    "            # convert from polars to pandas\n",
+    "            dwd_stations = dwd_stations.to_pandas()\n",
+    "\n",
+    "        have_data = (dwd_stations.end_date > start_date) & (dwd_stations.start_date < end_date)\n",
+    "        dwd_stations = dwd_stations[have_data]\n",
+    "\n",
+    "        request = all_station_request.filter_by_station_id(dwd_stations.station_id.astype(int))\n",
+    "        dwd_data = request.values.all().df\n",
+    "        if not isinstance(dwd_data, pd.DataFrame):\n",
+    "            # convert from polars to pandas\n",
+    "            dwd_data = dwd_data.to_pandas()\n",
+    "\n",
+    "        dwd_data = dwd_data.dropna().drop(columns=['dataset', 'quality'])\n",
+    "\n",
+    "        dwd_stations.to_parquet(cache_path + f'/dwd_stations-{end_date}.parquet.zst', engine='pyarrow', compression='zstd')\n",
+    "        dwd_data.to_parquet(cache_path + f'/dwd_data-{end_date}.parquet.zst', engine='pyarrow', compression='zstd')\n",
+    "\n",
+    "    return dwd_stations, dwd_data"
+   ],
+   "id": "0004-bf8c1d8b475519dc7862b9af3394e1c98d67f29981feff16de470d18939"
+  }
+ ],
+ "nbformat": 4,
+ "nbformat_minor": 5,
+ "metadata": {}
+}
+%% Cell type:markdown id:0001-9d10715da27eea0563e559964a2d386bf45e78b85abc9cb22a230f2e35c tags:
+
+# Gitter-Transformation
+
+Installieren Sie ggf. fehlende Pakete:
+
+%% Cell type:code id:0002-2b71ddfdefd192f1b26f0e96bf6da9ba9144a57b5b0e0e4e81bbe24bb52 tags:
+
+``` 
+!pip install wetterdienst
+```
+
+%% Cell type:markdown id:0003-d9bedf792849d3240f20100ae2d5727c2843229862be9f3803725199e01 tags:
+
+1.  Laden Sie die Wetterdaten mit dem vorbereiteten Code und schauen Sie
+    sich die Daten grob an.
+2.  Laden Sie das TK 25 Gitter (schon vorbereitet) und die Kreise (und
+    kreisfreie Städte) analog.
+3.  Interpolieren Sie die Temperaturen mit $k$-NN in einem der Gitter.
+    Als Koordinaten können Sie die `centroid`-Koordinaten verwenden.
+4.  Transformieren Sie die Daten in das andere Gitter.
+
+%% Cell type:code id:0004-bf8c1d8b475519dc7862b9af3394e1c98d67f29981feff16de470d18939 tags:
+
+``` 
+import geopandas as gpd
+import pandas as pd
+
+import zstandard
+from wetterdienst.provider.dwd.observation import (DwdObservationRequest,
+                                                   DwdObservationResolution)
+
+
+cache_path = '.'
+
+tk25_path = 'tk25-grid-utm32s.geojson.zst'
+with open(tk25_path, 'rb') as fh:
+    dctx = zstandard.ZstdDecompressor()
+    with dctx.stream_reader(fh) as reader:
+        grid = gpd.read_file(reader).set_index('id')
+
+
+def retrieve_weatherdata(start_date="2005-01-01", end_date="2023-08-31", drop_cache=False):
+    try:
+        if drop_cache:
+            raise FileNotFoundError('Dropping cache')
+
+        dwd_stations = pd.read_parquet(cache_path + f'/dwd_stations-{end_date}.parquet.zst')
+        dwd_data = pd.read_parquet(cache_path + f'/dwd_data-{end_date}.parquet.zst')
+    except FileNotFoundError:
+        all_station_request = DwdObservationRequest(
+            parameter=["temperature_air_mean_200",  # monthly mean temperature 2 m above ground
+                       "precipitation_height"],     # monthly sum of precipitation height
+            resolution=DwdObservationResolution.MONTHLY,
+            start_date=start_date,
+            end_date=end_date,
+        )
+
+        dwd_stations = all_station_request.all().df
+        if not isinstance(dwd_stations, pd.DataFrame):
+            # convert from polars to pandas
+            dwd_stations = dwd_stations.to_pandas()
+
+        have_data = (dwd_stations.end_date > start_date) & (dwd_stations.start_date < end_date)
+        dwd_stations = dwd_stations[have_data]
+
+        request = all_station_request.filter_by_station_id(dwd_stations.station_id.astype(int))
+        dwd_data = request.values.all().df
+        if not isinstance(dwd_data, pd.DataFrame):
+            # convert from polars to pandas
+            dwd_data = dwd_data.to_pandas()
+
+        dwd_data = dwd_data.dropna().drop(columns=['dataset', 'quality'])
+
+        dwd_stations.to_parquet(cache_path + f'/dwd_stations-{end_date}.parquet.zst', engine='pyarrow', compression='zstd')
+        dwd_data.to_parquet(cache_path + f'/dwd_data-{end_date}.parquet.zst', engine='pyarrow', compression='zstd')
+
+    return dwd_stations, dwd_data
+```
--- a/10-geodaten/bike_track_2024-05-25.geojson.zst
+++ b/10-geodaten/bike_track_2024-05-25.geojson.zst
--- a/10-geodaten/dwd_data-2023-08-31.parquet.zst
+++ b/10-geodaten/dwd_data-2023-08-31.parquet.zst
--- a/10-geodaten/dwd_stations-2023-08-31.parquet.zst
+++ b/10-geodaten/dwd_stations-2023-08-31.parquet.zst
--- a/10-geodaten/kreise-geo-mittel.geojson.zst
+++ b/10-geodaten/kreise-geo-mittel.geojson.zst
--- a/10-geodaten/tk25-grid-utm32s.geojson.zst
+++ b/10-geodaten/tk25-grid-utm32s.geojson.zst