【发布时间】:2019-04-22 20:09:52
【问题描述】:
我的问题包括以下内容:给我两对角(在球坐标中),它们由两部分组成——方位角和纬度角。如果我们无限扩展两个角度(从而增加它们各自的半径)以形成一条指向这对角度给定方向的长线,那么我的目标是确定
- 如果它们相交或彼此非常接近并且
- 它们相交的确切位置。
目前我尝试了几种方法:
最明显的方法是迭代比较每个半径,直到两者之间存在匹配或足够小的距离。 (当我说比较每个半径时,我指的是将每个球坐标转换为笛卡尔坐标,然后找到两者之间的欧几里得距离)。但是,这个运行时是 $O(n^{2})$,如果我试图扩展这个程序,这会非常慢
第二种最明显的方法是使用优化包来找到这个距离。不幸的是,我不能迭代地优化包,并且在一个实例之后优化算法重复相同的答案,这是没有用的。
最不明显的方法是直接从角度计算(使用微积分)精确的半径。虽然这是一种快速的方法,但它并不是非常准确。
注意:虽然交叉点始终位于零原点 (0,0,0) 似乎很简单,但情况并非总是如此。有些点永远不会相交。
方法 (1) 的代码
def match1(azimuth_recon_1,colatitude_recon_1,azimuth_recon_2, colatitude_recon_2,centroid_1,centroid_2 ):
# Constants: tolerance factor and extremely large distance
tol = 3e-2
prevDist = 99999999
# Initialize a list of radii to loop through
# Checking iteravely for a solution
for r1 in list(np.arange(0,5,tol)):
for r2 in list(np.arange(0,5,tol)):
# Get the estimates
estimate_1 = np.array(spher2cart(r1,azimuth_recon_1,colatitude_recon_1)) + np.array(centroid_1)
estimate_2 = np.array(spher2cart(r2,azimuth_recon_2,colatitude_recon_2))+ np.array(centroid_2)
# Calculate the euclidean distance between them
dist = np.array(np.sqrt(np.einsum('i...,i...', (estimate_1 - estimate_2), (estimate_1 - estimate_2)))[:,np.newaxis])
# Compare the distance to this tolerance
if dist < tol:
if dist == 0:
return estimate_1, [], True
else:
return estimate_1, estimate_2, False
## If the distance is too big break out of the loop
if dist > prevDist:
prevDist = 9999999
break
prevDist = dist
return [], [], False
方法(3)的代码
def match2(azimuth_recon_1,colatitude_recon_1,azimuth_recon_2, colatitude_recon_2,centriod_1,centroid_2):
# Set a Tolerance factor
tol = 3e-2
def calculate_radius_2(azimuth_1,colatitude_1,azimuth_2,colatitude_2):
"""Return radius 2 using both pairs of angles (azimuth and colatitude). Equation is provided in the document"""
return 1/((1-(math.sin(azimuth_1)*math.sin(azimuth_2)*math.cos(colatitude_1-colatitude_2))
+math.cos(azimuth_1)*math.cos(azimuth_2))**2)
def calculate_radius_1(radius_2,azimuth_1,colatitude_1,azimuth_2,colatitude_2):
"""Returns radius 1 using both pairs of angles (azimuth and colatitude) and radius 2.
Equation provided in document"""
return (radius_2)*((math.sin(azimuth_1)*math.sin(azimuth_2)*math.cos(colatitude_1-colatitude_2))
+math.cos(azimuth_1)*math.cos(azimuth_2))
# Compute radius 2
radius_2 = calculate_radius_2(azimuth_recon_1,colatitude_recon_1,azimuth_recon_2,colatitude_recon_2)
#Compute radius 1
radius_1 = calculate_radius_1(radius_2,azimuth_recon_1,colatitude_recon_1,azimuth_recon_2,colatitude_recon_2)
# Get the estimates
estimate_1 = np.array(spher2cart(radius_1,azimuth_recon_1,colatitude_recon_1))+ np.array(centroid_1)
estimate_2 = np.array(spher2cart(radius_2,azimuth_recon_2,colatitude_recon_2))+ np.array(centroid_2)
# Calculate the euclidean distance between them
dist = np.array(np.sqrt(np.einsum('i...,i...', (estimate_1 - estimate_2), (estimate_1 - estimate_2)))[:,np.newaxis])
# Compare the distance to this tolerance
if dist < tol:
if dist == 0:
return estimate_1, [], True
else:
return estimate_1, estimate_2, False
else:
return [], [], False
我的问题有两个:
是否有更快、更准确的方法来查找两者的半径 积分?
如果是,我该怎么做?
编辑:我正在考虑创建两个半径的两个 numpy 数组,然后通过 numpy 布尔逻辑比较它们。但是,我仍然会迭代地比较它们。有没有更快的方法来进行这种比较?
【问题讨论】:
-
我不明白交点总是(0,0,0)。你能解释一下“如果我们无限扩展两个角度(从而增加它们各自的半径)”的意思吗
-
啊,对不起,我应该更清楚。每对角度都给我们一个方向。所以,当我说“无限延伸两个角度(从而增加它们各自的半径)”时,我基本上是说要画一条长线,指向这对角度给定的方向。
-
另外,交点并不总是在点 (0,0,0)。有些角度对永远不会相交或彼此接近
-
那个长线原点是 (0,0,0) ?
-
不,线的原点(角度 1)是某个预先确定的点(centroid_1: x1,y1,z1)——我在上面的代码中编辑过——其中 x1, y1 , 和 z1 > 0。同样,线的原点(角度 2)是某个预定点 (centroid_2:x2,y2,z2),其中 x2、y2 和 z2 > 0
标签: python-3.x distance numpy-ndarray numpy-einsum