Teile und herrsche – das gilt oft auch für Geodaten

Kürzlich bekam ich einen Fachartikel zur Begutachtung, der die kartographische Darstellung eines globalen Datensatzes auf einem digitalen Globus behandelte. Da es sich beim besagten Datensatz um eine große Menge von Punkten handelte, die jeweils das Auftreten eines Ereignisses markierten, hatten die Autoren die sinnvolle Idee, die Punkte in Zellen gleicher Größe zusammenzufassen.  In Folge sollte jede Zelle (genaugenommen nur deren Mittelpunkt) basierend auf der Anzahl der darin enthaltenen Ereignisse visualisiert werden. So weit so gut, wir alle kennen bzw. nutzen den Zugang der Tesselation, also der Unterteilung des Raumes in regelmäßige Kacheln, die dann als Aggregationseinheit für Menschen, Bäumen, Unfälle etc. dienen können. Während solche Kacheln früher nahezu immer Quadratrasterzellen waren, sind Hexagone mittlerweile ähnlich häufig anzutreffen.  Letztere haben den Vorteil, dass die Mittelpunkte aller benachbarter Zellen denselben Abstand zur Ausgangszelle haben. Dem steht jedoch der Nachteil gegenüber, dass in hierarchischen Systemen, die mehrere räumliche Auflösungsstufen umfassen, Hexagone höherer Ordnung nicht randscharf über Hexagone niedrigerer Ordnung passen. Bei Quadraten besteht dieses Problem nicht. Wird jede Quadratzelle in vier oder neun kleinere Zellen unterteilt, ist die eindeutige Zuordnung einer Zelle zur „Elternzelle“ auf höherem Aggregationsniveau problemlos möglich.

Zelleinteilungsmodelle Quadratraster Aperture 4 bzw. 9 und Hexagone Aperture 7

Hierarchische Tesselation auf Quadratraster- und Hexagonbasis. Wie beim Vergleich der beiden Quadratrastervarianten zu sehen, kann die Auflösung mit jeder Hierarchiestufe unterschiedlich stark zunehmen. Man spricht dabei auch von der „Aperture“ (dt. „Blendenöffnung“) eines Systems, die das Flächenverhältnis einer Kachel zur Kachel der nächsthöheren Auflösungsebene beschreibt. Dementsprechend hat das linke Quadratraster Aperture 4, das mittlere Quadratraster Aperture 9 und das Hexagonraster Aperture 7 (eigener Entwurf).

Kommen wir aber zurück zum erwähnten Manuskript. Dort wurde beschrieben, wie in der gängigen Web-Mercator Projektion ein hexagonales Raster über die Punktdaten gelegt wurde, um jedem Sechseck die Zahl der darin aufgetretenen Ereignisse zuzuweisen. Als ich anschließend las, dass jede Zelle eine Fläche von 30.000 Quadratkilometern abdeckt, stutzte ich: Ein regelmäßiges Raster auf Basis einer Web-Mercator Projektion mit ihrer enormen Flächenverzerrung in höheren Breiten geht einfach nicht mit einheitlichen Flächeninhalten je Zelle zusammen! Beim Blick auf den frei rotierbaren digitalen Globus wurde das auch deutlich – dort schrumpfte die Distanz der Zellmittelpunkte auf der geographischen Breite von Stockholm oder St. Petersburg auf etwa die Hälfte der entsprechenden Distanz am Äquator. Die real abgedeckte Fläche einer Zelle bei St. Petersburg hat sich also auf etwa ein Viertel der Fläche einer äquatornahen Zelle reduziert. Keine gute Voraussetzung für ein „Zählraster“, welches das Vorkommen von Objekten oder Ereignissen vergleichbar machen soll.

Zählraster auf Web-Merkator Basis und das Problem der kleineren Zellgröße in höheren Breiten

Ein Zählraster basierend auf einer Web-Mercator Projektion (links). Rechts die Darstellung dieses Rasters am digitalen Globus, der bei konstantem Maßstab (bzw. Betrachtungsabstand) so gedreht wurde, dass der Betrachter jeweils zentral über dem jeweiligen Ausschnitt (Skandinavien oben, Kongo unten) positioniert ist (eigener Entwurf).

Die gerechte Aufteilung der Erde?

Wie könnte nun aber ein globales Zählraster mit gleich großen Zellen definiert werden? Geographische Koordinaten scheiden hier aus, da die von Längen- und Breitenkreisen gebildeten Flächen (also etwa 1° Winkelabschnitte in Länge und Breite) zu den Polen hin immer kleiner werden. Datensätze wie das SRTM-Höhenmodell, die auf dieser Basis strukturiert sind, haben darüber hinaus den Nachteil, dass wir es in Folge mit trapezförmigen Rasterzellen variabler Breite zu tun haben, was für GIS-Analysen, die quadratische Raster voraussetzen, problematisch sein kann. Eine Möglichkeit wäre die Verwendung einer flächentreuen Projektion (z.B. Equal Earth) als Basis, über die dann ein regelmäßiges, gegebenenfalls hierarchisches Raster gelegt wird. Zwar hätte dann jede Rasterzelle denselben Flächeninhalt, die tatsächliche Form des von einer Zelle in der Natur abgedeckten Gebietes, würde in vielen Bereichen der Welt jedoch stark verzerrt sein. Für eine rein vergleichende, globale Vorkommensbilanzierung wäre das vielfach kein großes Problem, weiterführende Rasteranalysen, die in irgendeiner Weise auf einheitlichen Distanzannahmen im Raster basieren (z.B. „die Distanzen zu den diagonal angrenzenden Nachbarzellen rechts oben und links oben sind gleich“), wären jedoch unmöglich, bzw. würden zu falschen Ergebnissen führen.

Equal Earth Projektion mit Zählraster und Tissot Indikatrix

Ein Zählraster über eine flächentreue Projektion (Equal Earth) gelegt. Zwar ist die von einer Zelle abgedeckte Erdoberfläche jeweils gleich groß, die Form des abgedeckten Bereichs würde jedoch stark variieren, wie die orange dargestellte Tissotsche Indikatrix zeigt. Diese bildet real (d.h. auf der Kugeloberfläche liegende) gleich große Kreise im Kartenbild ab. So wie die aus der Projektion resultierenden Ellipsen in der Karte stark von der Kreisform abweichen, würden die Rasterzellen der Karte in der Natur von der Quadratform abweichen. Während sich die Abweichung im Zentrum der Projektion noch in Grenzen hält, nimmt sie zum Rand hin deutlich zu. (Quelle: modifiziert nach Justin Kunimune, CC0, via Wikimedia Commons)

Was wäre aber nun, wenn wir die Erdoberfläche nicht auf eine einzige Ebene projizieren, sondern stattdessen auf einen kugelähnlicheren, dreidimensionalen Körper mit ebenen Seiten (Polyeder), wie zum Beispiel einem über die Erde gelegten Würfel? Wenn die dann kleineren Abbildungsflächen der gekrümmten Erdoberfläche näher sind, müssten die durch die Verebnung entstehenden Verzerrungen deutlich geringer sein!

Polyeder aus Papier

Je kugelähnlicher der Polyeder (von oben: Oktaeder, Tetraeder, leicht beschädigter Ikosaeder) desto geringer die Verzerrungen bei der Projektion (…und desto einfacher die Jonglage ;-)).

Die Lösung heißt DGGS – Die Erde in Billionen von Teilen

Flächentreue Discrete Global Grid Systems (DGGS) verfolgen genau diesen Ansatz: Die Erdoberfläche wird flächentreu auf einen platonischen Körper wie etwa einen Würfel (gebildet aus 6 Quadraten) oder einen Ikosaeder (gebildet aus 20 gleichseitigen Dreiecken) projiziert. Basierend auf dieser, nicht nur flächentreuen, sondern auch halbwegs formerhaltenden Abbildung  wird dann ein hierarchisches Raster definiert. Dort kann auf jeder Hierarchie- bzw. Auflösungsebene jede Position auf der Erde einer eindeutigen Zelle zugewiesen werden. Im Gegensatz zur Visualisierung der gesamten Erde auf einen Blick, welche das „Aufschneiden und flache Ausbreiten des Polyedermantels“ erfordern würde, funktionieren Berechnungen in DGGS schnittfrei über die Polyederkanten hinweg und damit um die gesamte Erde herum.

Fuller Projektion

Die Idee, die Erde auf einen platonischen Körper zu projizieren, um Verzerrung gering zu halten ist nicht neu: 1954 präsentierten Buckminster Fuller und Shoji Sadao diese (allerdings nicht flächentreue) Projektion der Erde auf einen Ikosaeder, der anschließend so „geschnitten und aufgefaltet“ wurde, dass alle Schnitte in den Ozeanen liegen. (Quelle: Justin Kunimune, CC0, via Wikimedia Commons)

Kurz gesagt bieten DGGS also blattschnittfreie, erdumspannende Raster in hierarchisch organisierten Maßstabsebenen. Im Gegensatz zu punktorientierten, d.h. kontinuierlichen Längen- und Breitenangaben traditioneller, globaler Bezugssysteme diskretisieren DGGS die Erde in maßstäblich hierarchisch organisierte, flächige Kacheln.  Anwendungsvorteile ergeben sich dabei vor allem bei der Integration und Abfrage sehr großer, globaler Datenbestände.  Anstatt komplexe, sphärische Berechnungen in entsprechenden Vektordatenbeständen durchzuführen, werden diese zunächst in ein DGGS überführt, beziehungsweise dort auf passenden (oder sämtlichen) Auflösungsebenen flächig aggregiert. Abfragen und Analysen der nunmehr rasterbasiert vorliegenden Daten sind dank einer räumlichen Indizierung der Zellen sehr effizient möglich. Nicht ohne Grund wurde das derzeit prominenteste DGGS „H3“ vom Taxivermittler Uber entwickelt, um seine weltweiten Daten zu lokalem Angebot und Nachfrage zu analysieren. Ziel war dabei die Implementierung einer lokal differenzierten, dynamischen Preisgestaltung.

Die Spezies und vier Exemplare

2017 veröffentlichte das Open Geospatial Consortium (OGC) die erste abstract specification zu DGGS, also eine technologieunabhängige Beschreibung der dahinterstehenden Konzepte. Während in dieser ersten Version neben dem Ansatz der hierarchischen Partition der gesamten Erde auch gleich große Kacheln innerhalb einer Hierarchieebene konstituierend für ein DGGS waren, ist das in der aktuellen Version 2.0 von 2021 nicht mehr der Fall. Flächentreue DGGS werden dort als eine Unterkategorie von DGGS beschrieben. Das liegt wohl nicht zuletzt daran, dass populäre DGGS wie das oben erwähnte H3 (https://h3geo.org/) oder Googles „S2 Geometry“ (https://s2geometry.io/) leichte Variationen im Flächeninhalt der Sechseck- (H3) bzw. Quadratkacheln (S2) aufweisen.

H3 (wie auch S2) ist open source, verfügt über eine eigene API, die in mehreren Sprachen adressiert werden kann und ist mittlerweile so verbreitet, dass klassische GI-Systeme wie QGIS oder ArcGIS die spezifische Kachelstruktur erzeugen können. Es basiert vorwiegend auf Hexagonen (Aperture 7) in 15 Auflösungsstufen, wobei eine Zelle in der höchsten Auflösung eine Fläche von unter einem Quadratmeter aufweist. „Vorwiegend“ deshalb, weil die Tesselation einer Ikosaederoberfläche ausschließlich mit Sechsecken nicht möglich ist, sondern immer auch 12 Fünfecke erfordert. Im Fall von H3 liegen diese meist im Meer, stören also kaum.

Ein klassischer Fußball besteht aus 12 Fünfecken und 20 Sechsecken. In H3 sind in jeder Auflösung ebenfalls 12 Fünfecke erforderlich, egal ob 110 (höchstes Aggregationsniveau) oder aber 569.707.381.193.150 Sechsecke (höchste Auflösung) hinzukommen. Die Pentagone sind deutlich kleiner, nach OGC-Spezifikation wird Flächengleichheit (für flächentreue DGGS, was für H3 ohnehin nicht zutrifft) aber nur innerhalb derselben Kachelgeometrie gefordert.

Beispiele für in der Wissenschaft etablierte, flächentreue DGGS sind ISEA3H (Sahr et al. 2003) oder das aus der Astronomie stammende HEALPix (Górski et al. 1998). ISEA3H verwendet flächentreu auf einen Ikosaeder projizierte Hexagone (Aperture 3), HEALPix quadratisch bis rautenförmige Kacheln, die auf einen Rhombendodekaeder projiziert sind. Durch ihre Flächentreue eignen sich beide Systeme zur Modellierung und statistischen Auswertung von globalen bio- und geowissenschaftlichen Daten, wie sie beispielsweise in der Habitat- oder Klimamodellierung vorkommen. Für die Aggregation der eingangs erwähnten Punktdaten, wären diese DGGS also bestens geeignet.

Górski, K. M., Hivon, E., & Wandelt, B. D. (1998) Analysis issues for large CMB data sets. https://doi.org/10.48550/arXiv.astro-ph/9812350.

Li, M., Stefanakis, E. (2020) Geospatial Operations of Discrete Global Grid Systems—a Comparison with Traditional GIS. Journal of Geovisualization and Spatial Analysis 4, 26. https://doi.org/10.1007/s41651-020-00066-3

Mechenich, M.F., Žliobaitė, I. (2023) Eco-ISEA3H, a machine learning ready spatial database for ecometric and species distribution modeling. Sci Data 10, 77. https://doi.org/10.1038/s41597-023-01966-x

Open Geospatial Consortium (2021) Topic 21 – Discrete Global Grid Systems – Part 1 Core Reference system and Operations and Equal Area Earth Reference System. Vers. 2.0 https://docs.ogc.org/as/20-040r3/20-040r3.html

Sahr, K., D. White & A. J. Kimerling (2003) Geodesic Discrete Global Grid Systems, Cartography and Geographic Information Science, 30:2, 121-134, DOI: 10.1559/152304003100011090