0

それぞれが円上の 2 点で定義される 2 つの大円の交点 (緯度と経度) を計算しようとしています。ここで概説されている方法に従おうとしています。しかし、私が得た答えは間違っています。私のコードは以下の通りです。

import numpy as np
from numpy import cross
from math import cos, sin, atan2, asin, asinh

################################################
#### Intersection of two great circles.
# Points on great circle 1.
glat1 = 54.8639587
glon1 = -8.177818

glat2 = 52.65297082
glon2 = -10.78064876

# Points on great circle 2.
cglat1 = 51.5641564
cglon1 = -9.2754284

cglat2 = 53.35422063
cglon2 = -12.5767799

# 1. Put in polar coords.

x1 = cos(glat1) * sin(glon1)
y1 = cos(glat1) * cos(glon1)
z1 = sin(glat1)

x2 = cos(glat2) * sin(glon2)
y2 = cos(glat2) * cos(glon2)
z2 = sin(glat2)


cx1 = cos(cglat1) * sin(cglon1)
cy1 = cos(cglat1) * cos(cglon1)
cz1 = sin(cglat1)

cx2 = cos(cglat2) * sin(cglon2)
cy2 = cos(cglat2) * cos(cglon2)
cz2 = sin(cglat2)


# 2. Get normal to planes containing great circles.
#    It's the cross product of vector to each point from the origin.

N1 = cross([x1, y1, z1], [x2, y2, z2])
N2 = cross([cx1, cy1, cz1], [cx2, cy2, cz2])


# 3. Find line of intersection between two planes.
#    It is normal to the poles of each plane.

L = cross(N1, N2)


# 4. Find intersection points.

X1 = L / abs(L)
X2 = -X1


ilat = asin(X1[2]) * 180./np.pi
ilon = atan2(X1[1], X1[0]) * 180./np.pi

また、これは地球の表面上にあることにも言及する必要があります（球体を想定）。

4

2 に答える 2

1

上記のコメントの DSM からの解決策です。角度は度単位ですが、sin と cos はラジアンを期待しています。ラインも

X1 = L / abs(L)

あるべきです、

X1 = L / np.sqrt(L[0]**2 + L[1]**2 + L[2]**2)
于 2013-09-20T16:02:04.707 に答える
0

実行する必要があるもう 1 つの修正は、x および y 次元で "lon" の前の cos/sin を変更することです。

x = cos(lat) * cos(lon)
y = cos(lat) * sin(lon)
z = sin(lat)

これは、角度から球面システムへの元の変換が極/方位球角度を使用して行われ、緯度/経度角度と同じではないためです (wiki it https://en.wikipedia.org/wiki/Spherical_coordinate_system )。

于 2016-07-04T14:20:51.737 に答える