Digitale Ableitung von Relieffazetten

Eine automatisierte Ableitung von Relieffazetten bietet sich gegenüber einer manuellen Ableitung auf Grundlage des Höhenlinienbildes topografischer Karten an. Wesentliche Gründe hierfür sind:

• Eine weltweite Verfügbarkeit von digitaler Reliefinformation mit einer zumindest für Analysen des Meso- und Makroreliefs hinreichenden Genauigkeit
• Die Möglichkeit der mathematischen Formalisierung der Fazettierung bei gleichzeitiger Möglichkeit der Steuerung der Feinheit des Resultats über quantitative Parameter
• Beliebige Wiederholbarkeit und Übertragbarkeit auf andere Gebiete bei exakt reproduzierbaren Resultaten
• Einbindung der Fazettenberechnung in ein GIS, damit die Möglichkeit zur kartographischen Ausgabe der Resultate und zur Weiterbearbeitung z.B. zur Synthese komplexer Formen aus Relieffazetten.

Diese Versuche wurden bereits unternommen und liefern gute Ergebnisse. Ein bisher ungelöstes Problem besteht allerdings in der sehr langen Rechenzeit. Diese ist in der Einbindung des bestehenden Programmes in die GIS-Software Arc/Info begründet; eine eigenständige Implementierung des Programmes könnte jedoch die Rechenzeit erheblich reduzieren.

Im Folgenden werden die wesentlichen Schritte der digitalen Berechnung von Relieffazetten erläutert.

Digitale Reliefbeschreibung

Wie allgemein bekannt ist, existieren im Wesentlichen zwei digitale Formen der Reliefbeschreibung [Prechtel 2002]: die Dreiecksvermaschung aus irregulär über eine Fläche verteilten Punkten (im Englischen als Triangular Irregular Network, kurz TIN bezeichnet) sowie die regelmäßige Rasterstruktur, bei der jedem Mittelpunkt einer Zellmatrix ein diskreter Höhenwert zugeordnet ist. Obwohl die Dreiecksvermaschung einen wesentlichen Vorteil besitzt, nämlich eine flexible Stützpunktlage, die in schroffem Relief kleine Dreiecksmaschen liefert und in glattem Relief entsprechend großflächigere, sind Rasterhöhenmodelle weiter verbreitet. Sie lassen sich nämlich weitgehend automatisiert aus Messdaten der Fernerkundung erzeugen. Dabei kommen einmal Verfahren der digitalen Bildkorrelation in Betracht: Für jeden Bildpunkt einer Stereoaufnahme werden im Partnerbild über Korrelation homologe Punkte gesucht. Aus der relativen Lage von gemeinsamen Punkten innerhalb beider Bilder ergeben sich die so genannten Parallaxen. Sofern nun Lage der Kamera, ihre Brennweite und die Richtung der optischen Achse (Blickrichtung der Kamera) bekannt sind, lassen sich Bildpunkt für Bildpunkt die entsprechenden Geländehöhen rekonstruieren. Eine zweite Möglichkeit ist die der Höhenbestimmung aus Interferenzmustern kohärenter Strahlung. Dieses Aufnahmeverfahren wurde bei der SRTM-Mission genutzt, das auf zwei versetzt angebrachten Radarsensoren im Weltraum basiert. Wenn die Ergebnisse auch nicht völlig ohne Fehler sind, so liegt aus den Ergebnissen der Satellitenmission immerhin ein weltweites Rasterhöhenmodell für alle Landflächen in einer horizontalen Zellauflösung von 90m x 90m vor, welches für zahlreiche geomorphologische Aufgaben eine Datengrundlage liefern kann. Da diese Datenbasis aber noch sehr jung ist, fehlen noch repräsentative Auswertungen ihrer Höhengenauigkeit, vor allem im  Hochgebirgsrelief.

Grundlage unserer Auswertung war ein am Institut für Kartographie berechnetes Rasterhöhenmodell eines Gebietes im Russischen Altai, welches die Katunkette, die Uimon-Steppe und den nördlichen Rand der Terekta-Kette abdeckt. Die Fläche beträgt ca. 8500 km², die Auflösung des Originalmodells beträgt 100m x 100m. Für die Auswertung wurde allerdings ein feiner interpoliertes Modell mit Zellen der Größe 17m x 17m verwendet. Das Modell wurde aus einer Kombination von Karteninformation und photogrammetrischer Auswertung von Stereoaufnahmen der russischen Kamera MK4 erzeugt [Prechtel 2003] und ist in seiner Qualität heterogen. Einen Ausschnitt dieser Eingangsdaten für die Fazettenberechnung mit simulierter Beleuchtung zeigt Abbildung 1. Bei allen nachfolgenden Betrachtungen ist zu berücksichtigen, dass beliebige Ableitungen aus Reliefmodellen immer die Fehler und Ungenauigkeiten der Datengrundlage mitführen.

Berechnung initialer Fazetten

Relieffazetten stellen eine Abstraktion des Reliefs dar. Einen mathematisch glatten, also weder in der Horizontalebene noch in Höhenlinienrichtung aufgespannten Vertikalebene völlig krümmungsfreien Geländeausschnitt, gibt es nahezu nicht. Gleichzeitig existiert keine exakte Datengrundlage für das Relief, weder wenn sie aus der klassischen topografischen Aufnahme noch wenn sie aus Fernerkundungsverfahren abgeleitet wurde. Daher ist die Berechnung von Fazetten nicht abgeschlossen, indem man für jede Zelle eines Höhenmodells lokale Neigung und Hangrichtung, oder gleichbedeutend die Flächennormale bestimmt und anschließend alle Zellen mit gleichen Resultaten zu Gruppen zusammenfasst. Die Zahl der Fazetten wäre so, da (siehe oben) in den Modellen keine mathematisch glatten Hänge existieren, nahezu gleich der Anzahl der Zellen des Höhenmodells.

Der entscheidende Punkt ist aber die sinnvolle Zusammenfassung von Initialfazetten zu größeren Einheiten, denen eine wesentliche Bedeutung zur geomorphologischen, aber auch darüber hinaus naturräumlichen Ordnung einer Landschaft zukommt. Einfach zu berechnen, aber im Resultat nicht befriedigend, wäre eine globale Gruppierung der Normalenvektoren in eine definierte Anzahl von schematischen Klassen. Sofort einsichtig ist, dass diese Vorgehensweise in vielen Fällen zu völlig willkürlichen Fazettengrenzen führen muss (Abbildung 2).

Abbildung 2: Höhenmodell mit 9 Zellen (links) mit Darstellung der Flächennormale in Polarkoordinaten (rechts). Schematische Gruppierung in Richtungssektoren führt zu keiner sinnvollen Fazettenbildung. Eine sinnvollere Fazettenbildung ist durch Gruppierung in beiden Diagrammen über verschiedene Grauwerte angedeutet.

Kern des neu entwickelten Algorithmus [Pollack 2003] ist eine Erzeugung von Initialfazetten durch eine iterative Berechnung von Trendflächen erster Ordnung (Gleichung 1) aus dem Originalhöhenmodell. Dieser Ansatz geht auf eine Anregung durch einem Vortrag zur Ermittlung von Kluftflächen zurück [Menz et al. 2002]. Zunächst wird eine Ausgleichsebene für das gesamte Modell berechnet. Die Ebene erfüllt die Bedingung, dass die quadrierten Höhendifferenzen zwischen Originalmodell und der Trendfläche minimiert sind. Aus dem zweiten Schritt erhält man für jede Zelle Differenzhöhen zu dieser Ausgleichsebene (Gleichung 2). Geländekanten bilden sich nun als Grate (positiv) und Talwege (negativ) im Differenzmodell ab. Nutzt man jedoch den Betrag der Differenz, dann erscheinen beide Typen von Strukturlinien (also Grate und Talwege) als positive Modellkanten oberhalb der Trendebene. Die Flächen zwischen diesen Kanten, im morphologischen Sinne Becken, erhält man durch Bestimmung von Einzugsgebieten (im hydrologischen Sinn) aus dem Differenzmodell. Nun setzt man die Berechnung der Trendflächen für jedes Einzugsgebiet einzeln fort und bestimmt wiederum die Differenzhöhen zum Originalmodell. Gebiete, bei denen die Differenzhöhen unterhalb eines vorab eingegebenen Schwellwertes verbleiben, werden nicht mehr bearbeitet, alle anderen Einzugsgebiete werden der nächsten Iteration zugeführt. Ebenso endet die Zerlegung der Fläche spätestens dann, wenn ein verbleibendes Einzugsgebiet weniger als drei Zellen umfasst, da sich dann die Ebenengleichung mit ihren drei Unbekannten nicht mehr lösen lässt.

(1) Gleichung der Trendfläche:      f(x,y) = a₀ + a₁x + a₂y

(2) Absolute Differenzhöhe einer Zelle zur Trendfläche:Dz=|z_o–z_t|*cos a

                                  z_o: Höhenwert einer Zelle im Originalmodell
                                  z_t: Höhenwert einer Zelle in der Trendfläche
                                   a: Zenitwinkel der Neigung

Die so erzeugten Initialfazetten sind jedoch zu fein und müssen im nächsten Schritt zu größeren Einheiten aggregiert werden. Dabei wird aus den genannten Gründen aber kein globaler Ansatz zur Gruppierung verwendet, sondern stets die Nachbarschaften der Initialfazetten betrachtet.

Neugruppierung zu endgültigen Fazetten

Der folgende Schritt folgt dem Verfahren des Regionenwachstums (englisch „region growing“). Für jede Fazette aus der initialen Berechnung – ein zusammenhängender Cluster von Zellen – werden zunächst Polygone gebildet. Die Parameter der Ebenengleichung werden in Normalvektoren umgerechnet und bilden Attribute der Flächen. Da ein Mosaik von Flächen vorliegt, sind zu jeder Fläche auch die jeweiligen Nachbarflächen bekannt, da sie eine gemeinsame Kante besitzen. Die Gruppierung zu endgültigen Fazetten kann dann so erfolgen, dass jeweils die Differenzwinkel der Normalvektoren benachbarter Flächen verglichen werden. Die erste Vereinigung erfolgt zwischen Flächen mit der geringsten Winkeldifferenz. Nach Vereinigung wird der nächst kleinste Differenzwinkel gesucht. Bei jeder Vereinigung wird zudem ein Öffnungswinkel errechnetet oder aktualisiert, der einen Kegelmantel beschreibt, welcher sämtliche Originalvektoren einschließt. Mit diesem Öffnungswinkel kann der Vereinigungsvorgang zwischen den Originalvektoren gesteuert werden. Liegen keine Flächen mehr vor, die ohne Überschreitung des maximal zulässigen Öffnungswinkels vereinigt werden können, ist die Berechnung beendet.

Darstellung als morphometrische Karte

Die Untersuchung ist ein Teil von vielfältigen physisch-geographischen Untersuchungen im zentralen Altai. Daher ist auch ein großer strukturierter GIS-Datenbestand mit zahlreichen topografischen und thematischen Informationsebenen vorhanden [Prechtel 2002]. In Anlehnung an den Maßstab der Eingangsdaten (1: 100 000 bis 1: 200 000) und das Detailniveau wurde das GIS als „Altai 100“ bezeichnet. Die Resultate der Fazettierung konnten zusammen mit Orientierungselementen daher automatisiert als „Morphometrische Karte“ ausgegeben werden. Aus den Orientierungsparametern der Fazetten konnten für die Kartendarstellung Hangneigungsklassen (Tab. 1) gebildet werden. Diese wurden mit Farben kodiert. Ebenso ließen sich konkave und konvexe Fazettengrenzen ausscheiden, die durch unterschiedlichen Strichtyp gekennzeichnet sind (konvex: durchgezogen, konkav: gerissen). Der Orientierung in der Karte dienen weiterhin Höhenlinien, Höhenpunkte sowie das Gewässernetz. Zudem wurde der Abbildung noch eine Reliefschummerung unterlegt. Einen Ausschnitt aus der Karte zeigt die Abbildung 3.

Weiterentwicklungsmöglichkeiten

Die erfolgreichen Versuche mit dem bestehenden Programm sollen in zwei Richtungen weiterverfolgt werden: Erstens könnte man die Rechenzeit durch eine Programmierung außerhalb des GIS wesentlich verbessern und die Resultate erst sekundär dem GIS übergeben. Zweitens könnte die Fazettenbildung noch geringfügig verbessert werden, indem zur Vereinigung von Initialfazetten neben der Winkeldifferenz noch die Kompaktheit der Vereinigungsfläche (Verhältnis Umfang zur Fläche) als ein zusätzlicher Steuerfaktor hinzugenommen würde.

Literatur

Prechtel, Nikolas (2002): „Digitales Höhenmodell“. Lexikonstichwort in: Bollmann, J. u. Koch, W.-G. (Hrsg.): Lexikon der Kartographie und Geomatik. Spektrum Akademischer Verlag, Heidelberg.

Pollack, Mario (2003): Methodischer Beitrag zur GIS-basierten mittelmaßstäbigen geomorphologischen Klassifizierung und Kartierung des Zentral-Altai auf Grundlage von Digitalem Höhenmodell und Geländeerhebungen. Diplomarbeit am Institut für Kartographie der TU Dresden, Betreuer: N. Prechtel und M. Buchroithner, 64 S.

Prechtel, Nikolas (2003): GIS-Aufbau für den Naturschutz im Russischen Altai. In: Geoinformationssysteme – Theorie, Anwendungen, Problemlösungen. Kartographische Bausteine, Bd. 21, Inst. für Kartographie, TU Dresden, S. 82–100.

Menz, Joachim, Martienßen, Thomas und Nguyen, Van Dien (2002): Photogrammetrishe Erfassung der Klüftigkeit des Gebirges und Zugänge zur weitgehend automatisierten Bereitstellung der Ergebnisse. Vortrag zum Ehrenkolloquium aus Anlass des 90. Geburtstages von Prof. Dr.-Ing. Werner Rüger, TU Dresden 18.1.2002, Pdf-Dokument im Internet: http://www1.tu-freiberg.de/~wwwmage/forschung/hase/rueger.pdf, 21 S.

Autoren:

Dr. habil. Wladimir Butwilowski, Institut für Polymerforschung, Dresden_ Deutschland, Wladimirbutwilowski@gmail.com

Dr. Nikolas Prechtel, Institut für Kartographie, TU Dresden, D-01062 Dresden, Deutschland, Nikolas.Prechtel@mailbox.tu-dresden.de.