Erstellen Sie Großkreisbögen, die auf der Web Mercator-Karte gut aussehen?

10

Ich versuche, eine Karte mit Daten zu einigen Flügen zu erstellen und möchte Great Circle Arcs verwenden, um die Quelle und die Ziele zu verbinden.

Grundsätzlich möchte ich etwas Ähnliches wie die berühmte Facebook-Karte machen: Geben Sie hier die Bildbeschreibung ein

Ich habe die in diesem Beitrag angegebenen Funktionen verwendet: /gis//a/5205/442 (dh dieser Blog-Artikel: http://anitagraser.com/2011/08/20/visualizing-global-connections / ) und ich könnten Linien bekommen, aber sie überqueren die internationale Datumsgrenze und bündeln sich an den Polen:

Geben Sie hier die Bildbeschreibung ein

@underdark hat im verlinkten Blogpost erwähnt, dass diese Zeilen geteilt werden müssen, aber ich weiß nicht, wie ich sie in PostGIS automatisch teilen soll.

Zusätzlich muss auch das Bündeln der Linien in der Nähe der Pole gelöst werden.

Wie mache ich beides, wenn ich die Punktpositionen für Start und Ende dieser Flüge habe?

Devdatta Tengshe
quelle
Könnten Sie eine Polarprojektion verwenden?
Ian Turton
2
Ich habe hier einige Beispiele gegeben: gis.stackexchange.com/questions/133026/… . Beachten Sie, dass eqdc keine echten Großkreise liefert, aeqd jedoch.
AndreJ

Antworten:

5

Die Prinzipien in diesem Blog-Beitrag werden auf Allzweck-PostGIS übertragen.

http://blog.cartodb.com/jets-and-datelines/

Verwenden Sie ST_Segmentizeim Grunde genommen die Geografie und ein bisschen Magie, um Datumsgrenzen zu schneiden.

Paul Ramsey
quelle
8

Sie könnten die Geodäten berechnen. Wenn Sie die Geodät von A nach B anzeigen möchten, können Sie zuerst die Entfernung und den Azimut von A nach B berechnen (inverses geodätisches Problem) und dann Punkte von A bis zu mehreren Punkten zwischen A und B berechnen (direktes geodätisches Problem). Ich habe ein einfaches Skript in Python hinzugefügt, das GeographicLib verwendet und nur das Zeug in GeoJSON ausgibt:

from geographiclib.geodesic import Geodesic
from geojson import MultiLineString

def geodesic(lat1, lon1, lat2, lon2, steps):
    inverse = Geodesic.WGS84.Inverse(lat1, lon1, lat2, lon2)
    linestrings = []
    coordinates = []

    for i in range(0, steps + 1):
        direct = Geodesic.WGS84.Direct(inverse['lat1'], inverse['lon1'], inverse['azi1'], (i / float(steps)) * inverse['s12'])
        if len(coordinates) > 0:
            if (coordinates[-1][0] < -90 and direct['lon2'] > 90) or (coordinates[-1][0] > 90 and direct['lon2'] < -90):
                linestrings.append(coordinates)
                coordinates = []
        coordinates.append((direct['lon2'], direct['lat2']))

    linestrings.append(coordinates)
    geojson = MultiLineString(linestrings)
    return geojson

linestrings = []

# San Francisco: 37.7793, -122.4192
# Bangalore: 12.9, 77.616667
for linestring in geodesic(37.7793, -122.4192, 12.95, 77.616667, 100)['coordinates']:
    linestrings.append(linestring)

# Boston: 42.357778, -71.059444
# Bangalore: 12.9, 77.616667
for linestring in geodesic(42.357778, -71.059444, 12.95, 77.616667, 100)['coordinates']:
    linestrings.append(linestring)

print(MultiLineString(linestrings))

Das Ergebnis ist die wahre Geodät zwischen den Punkten in WGS-84. Natürlich können Sie dann die Koordinaten in jede gewünschte Projektion umwandeln. Das auf geojson.io visualisierte Ergebnis sieht folgendermaßen aus:

Geben Sie hier die Bildbeschreibung ein

sema
quelle
2
Unter geographiclib.sourceforge.net/html/python/… finden Sie eine alternative (und schnellere!) Möglichkeit, dies mit GeographicLib zu tun. Beachten Sie auch die Verwendung des LONG_UNROLL-Flags, um die Kontinuität der Längengrade sicherzustellen.
cffk