Welche Werkzeuge in Python stehen für die Erstellung von großen Kreisabständen und Linien zur Verfügung?

20

Ich muss Python verwenden, um eine große Kreisentfernung zu erstellen - sowohl eine Zahl als auch vorzugsweise eine Art 'Kurve', mit der ich eine clientseitige Karte zeichnen kann. Das Format der Kurve ist mir egal - sei es WKT oder ein Satz von Koordinatenpaaren - aber ich möchte nur die Daten herausholen.

Welche Tools gibt es da draußen? Was soll ich benutzen?

Christopher Schmidt
quelle

Antworten:

8

Die Antworten anderer sind etwas eleganter, aber hier ist ein ultraschnelles, etwas unpythonisches Stück Python, das die Grundlagen liefert. Die Funktion verwendet zwei Koordinatenpaare und eine benutzerdefinierte Anzahl von Segmenten. Es ergibt sich eine Reihe von Zwischenpunkten entlang eines großen Kreiswegs. Ausgabe: schreibfertiger Text als KML. Vorsichtsmaßnahmen: Der Code berücksichtigt keine Antipoden und geht von einer sphärischen Erde aus.

Code von Alan Glennon http://enj.com Juli 2010 (Der Autor stellt diesen Code öffentlich zur Verfügung. Die Verwendung erfolgt auf eigenes Risiko.)

-

def tweensegs (Längengrad1, Breitengrad1, Längengrad2, Breitengrad2, Anzahl_der_Segmente):

import math

ptlon1 = longitude1
ptlat1 = latitude1
ptlon2 = longitude2
ptlat2 = latitude2

numberofsegments = num_of_segments
onelessthansegments = numberofsegments - 1
fractionalincrement = (1.0/onelessthansegments)

ptlon1_radians = math.radians(ptlon1)
ptlat1_radians = math.radians(ptlat1)
ptlon2_radians = math.radians(ptlon2)
ptlat2_radians = math.radians(ptlat2)

distance_radians=2*math.asin(math.sqrt(math.pow((math.sin((ptlat1_radians-ptlat2_radians)/2)),2) + math.cos(ptlat1_radians)*math.cos(ptlat2_radians)*math.pow((math.sin((ptlon1_radians-ptlon2_radians)/2)),2)))
# 6371.009 represents the mean radius of the earth
# shortest path distance
distance_km = 6371.009 * distance_radians

mylats = []
mylons = []

# write the starting coordinates
mylats.append([])
mylons.append([])
mylats[0] = ptlat1
mylons[0] = ptlon1 

f = fractionalincrement
icounter = 1
while (icounter <  onelessthansegments):
        icountmin1 = icounter - 1
        mylats.append([])
        mylons.append([])
        # f is expressed as a fraction along the route from point 1 to point 2
        A=math.sin((1-f)*distance_radians)/math.sin(distance_radians)
        B=math.sin(f*distance_radians)/math.sin(distance_radians)
        x = A*math.cos(ptlat1_radians)*math.cos(ptlon1_radians) + B*math.cos(ptlat2_radians)*math.cos(ptlon2_radians)
        y = A*math.cos(ptlat1_radians)*math.sin(ptlon1_radians) +  B*math.cos(ptlat2_radians)*math.sin(ptlon2_radians)
        z = A*math.sin(ptlat1_radians) + B*math.sin(ptlat2_radians)
        newlat=math.atan2(z,math.sqrt(math.pow(x,2)+math.pow(y,2)))
        newlon=math.atan2(y,x)
        newlat_degrees = math.degrees(newlat)
        newlon_degrees = math.degrees(newlon)
        mylats[icounter] = newlat_degrees
        mylons[icounter] = newlon_degrees
        icounter += 1
        f = f + fractionalincrement

# write the ending coordinates
mylats.append([])
mylons.append([])
mylats[onelessthansegments] = ptlat2
mylons[onelessthansegments] = ptlon2

# Now, the array mylats[] and mylons[] have the coordinate pairs for intermediate points along the geodesic
# My mylat[0],mylat[0] and mylat[num_of_segments-1],mylat[num_of_segments-1] are the geodesic end points

# write a kml of the results
zipcounter = 0
kmlheader = "<?xml version=\"1.0\" encoding=\"UTF-8\"?><kml xmlns=\"http://www.opengis.net/kml/2.2\"><Document><name>LineString.kml</name><open>1</open><Placemark><name>unextruded</name><LineString><extrude>1</extrude><tessellate>1</tessellate><coordinates>"
print kmlheader
while (zipcounter < numberofsegments):
        outputstuff = repr(mylons[zipcounter]) + "," + repr(mylats[zipcounter]) + ",0 "
        print outputstuff
        zipcounter += 1
kmlfooter = "</coordinates></LineString></Placemark></Document></kml>"
print kmlfooter
Glennon
quelle
8

GeographicLib hat eine Python-Oberfläche :

Dies kann Computer-Geodäten auf einem Ellipsoid (Abflachen auf Null setzen, um große Kreise zu erhalten) und Zwischenpunkte auf einer Geodät erzeugen (siehe die "Linien" -Befehle im Beispiel).

So drucken Sie Punkte auf der geodätischen Linie von JFK zum Flughafen Changi (Singapur) aus:

from geographiclib.geodesic import Geodesic
geod = Geodesic.WGS84

g = geod.Inverse(40.6, -73.8, 1.4, 104)
l = geod.Line(g['lat1'], g['lon1'], g['azi1'])
num = 15  # 15 intermediate steps

for i in range(num+1):
    pos = l.Position(i * g['s12'] / num)
    print(pos['lat2'], pos['lon2'])

->
(40.60, -73.8)
(49.78, -72.99)
(58.95, -71.81)
(68.09, -69.76)
(77.15, -65.01)
(85.76, -40.31)
(83.77, 80.76)
(74.92, 94.85)
...
cffk
quelle
Der Python - Port von GeographicLib ist jetzt bei pypi.python.org/pypi/geographiclib
cffk
Siehe auch dieses Papier: CFF Karney, Algorithmen für Geodäsie, J. Geod, DOI: dx.doi.org/10.1007/s00190-012-0578-z
cffk
7

pyproj verfügt über die Funktion Geod.npts , die eine Reihe von Punkten entlang des Pfades zurückgibt . Beachten Sie, dass die Anschlusspunkte nicht im Array enthalten sind, sodass Sie sie berücksichtigen müssen:

import pyproj
# calculate distance between points
g = pyproj.Geod(ellps='WGS84')
(az12, az21, dist) = g.inv(startlong, startlat, endlong, endlat)

# calculate line string along path with segments <= 1 km
lonlats = g.npts(startlong, startlat, endlong, endlat,
                 1 + int(dist / 1000))

# npts doesn't include start/end points, so prepend/append them
lonlats.insert(0, (startlong, startlat))
lonlats.append((endlong, endlat))
scruss
quelle
Vielen Dank! Lösung von bekannten & massiv genutzten Bibliothek hier zur Verfügung gestellt :)
tdihp
-2

Ich habe dieses Paket nicht verwendet, aber es scheint interessant und eine mögliche Lösung: http://trac.gispython.org/lab/wiki/Shapely

Hugo Estrada
quelle
7
AFAIK, formschön rechnet nicht sphäroidal:Shapely is a Python package for set-theoretic analysis and manipulation of **planar** features
fmark 25.07.10