【问题标题】:How does Python handle complicated calculations?Python如何处理复杂的计算?
【发布时间】:2015-08-25 08:16:09
【问题描述】:

这可能是一个非常简单的问题,但我似乎看不到它。 我有一个按顺时针顺序排列的点列表,并希望根据this 使用以下函数计算这些点(凸多边形)的质心:

def calculateCentroid(raLinks,raNodes, links, nodes):  

orderedPointsOfLinks = orderClockwise(raLinks,raNodes, links, nodes)

arg1 = 0
arg2 = 0
Xc = 0
Yc = 0
i = 0
for point in orderedPointsOfLinks:
    arg1 += point.Y*(orderedPointsOfLinks[i+1 if i+1<len(orderedPointsOfLinks) else 0].X)
    arg2 += (orderedPointsOfLinks[i+1 if i+1<len(orderedPointsOfLinks) else 0].Y)*point.X
    Xc += (point.X+(orderedPointsOfLinks[i+1 if i+1<len(orderedPointsOfLinks) else 0].X))*(((orderedPointsOfLinks[i+1 if i+1<len(orderedPointsOfLinks) else 0].Y)*point.X)-(point.Y*(orderedPointsOfLinks[i+1 if i+1<len(orderedPointsOfLinks) else 0].X)))
    Yc += (point.Y+(orderedPointsOfLinks[i+1 if i+1<len(orderedPointsOfLinks) else 0].Y))*(((orderedPointsOfLinks[i+1 if i+1<len(orderedPointsOfLinks) else 0].Y)*point.X)-(point.Y*(orderedPointsOfLinks[i+1 if i+1<len(orderedPointsOfLinks) else 0].X)))    
    i+=1

area = (arg1-arg2)*0.5
print area

X = -Xc/(6*area)
Y = -Yc/(6*area)
print X , "   ", Y

使用Arcpy计算面积和质心,上面的函数计算的面积是正确的,但是质心是错误的。

Xc 和 Yc 有什么问题我无法解决?

如果我按以下方式更改 for 循环,它会起作用:

    for point in orderedPointsOfLinks:
        y0 = point.Y
        x0 = point.X
        x1 = orderedPointsOfLinks[i+1 if i+1<len(orderedPointsOfLinks) else 0].X
        y1 = orderedPointsOfLinks[i+1 if i+1<len(orderedPointsOfLinks) else 0].Y
        a = x0*y1 - x1*y0
        area += a
        Xc += (x0+x1)*a
        Yc += (y0+y1)*a
         i+=1
    area *= 0.5
    print area

    X = Xc/(6*area)
    Y = Yc/(6*area)
    print X , "   ", Y

这里是检查代码的节点列表:

[(371623.876, 6159668.714),(371625.994, 6159661.094), (371624.319, 6159654.634), (371619.654, 6159649.86), (371614.194, 6159647.819), (371608.401, 6159648.449), (371601.544, 6159652.652), (371598.77, 6159658.058), (371599.318, 6159665.421), (371603.025, 6159671.805), (371611.372, 6159674.882 ), (371619.417, 6159673.065)]

【问题讨论】:

  • 就个人而言,如果我将第一个点的副本附加到点列表的末尾,我会发现我的代码更容易阅读,所以我没有每次我想参考第 (i+1) 个点时,都必须重新输入(正确)该环绕测试。这就是维基百科所做的。另外,为什么不只在列表的长度上迭代 i(如果将第一个点附加到末尾,则为 -1)并始终使用 i 进行索引?我试图从你的代码中解码你的表达,但放弃了。
  • @barny:你说得对,我可以用其他方式来做,将第一个点添加到列表末尾或在 for 循环结束后添加。但是使用列表推导可以更容易地理解代码(可能不会更有效),因为您不需要任何额外的步骤,只需要检查您是否在 lis 中,考虑第一个关闭循环。
  • 这是您的代码,但我的 0.01 美元说您的代码实际上难以辨认,而且难以辨认的代码难以理解、调试和维护。 “简化” for 循环以迭代点列表使您 a) 编写代码来初始化、使用和维护循环计数器 i,b) 每次需要 i+1 时都必须手动编码环绕到第一个节点- 八次你已经输入 'i+1 if i+1

标签: python list for-loop arcpy


【解决方案1】:

source

试试:

import numpy

tp = [(371623.876, 6159668.714),(371625.994, 6159661.094), (371624.319, 6159654.634), (371619.654, 6159649.86),\
      (371614.194, 6159647.819), (371608.401, 6159648.449), (371601.544, 6159652.652), (371598.77, 6159658.058), \
      (371599.318, 6159665.421), (371603.025, 6159671.805), (371611.372, 6159674.882 ), (371619.417, 6159673.065),(371623.876, 6159668.714)]



# cx = sigma (x[i]+x[i+1])*((x[i]*y[i+1]) - (x[i+1]*y[i] ))
# cy = sigma (y[i]+y[i+1])*((x[i]*y[i+1]) - (x[i+1]*y[i] ))
cx = 0
cy = 0

p = numpy.array(tp)

x = p[:, 0]
y = p[:, 1]

a = x[:-1] * y[1:]
b = y[:-1] * x[1:]

cx = x[:-1] + x[1:]
cy = y[:-1] + y[1:]



tp = tp[:-1] #dont need repeat


def area():
    tox=0
    toy=0
    for i in range(len(tp)):
        if i+1 == len(tp):
            tox += tp[-1][0]*tp[0][1]
        else:
            tox += tp[i][0]*tp[i+1][1]

    for i in range(len(tp)):
        if i+1 == len(tp):
            toy += tp[-1][1]*tp[0][0]
        else:
            toy += tp[i][1]*tp[i+1][0]
    return abs(tox-toy)*0.5

ar = area()
Cx = abs(numpy.sum(cx * (a - b)) / (6. * ar))
Cy = abs(numpy.sum(cy * (a - b)) / (6. * ar))
print Cx,Cy

警告!

tp[0] == tp[-1]

所以:第一个和最后一个坐标是相同的值...

【讨论】:

  • 我对该区域没有任何问题。问题在于计算坐标,因为它在区域正确时返回错误值。我怀疑python处理复杂方程的方式是问题所在,因为坐标值时不时地变化,如果点列表是固定的,这是不合理的!!!
  • @msc87 在立陶宛的一个地方?
  • 这是瑞典的一个地方......无论如何,感谢您的代码,但主要问题仍未得到解答......为什么我的第一个功能不起作用?而第二个是。这更多是为了理解您在任何其他复杂计算中可能遇到的问题,就像这样
  • 朋友我得到了很多计算错误,我忘记了首先和最后需要相同的值。并创建三角形形式的第一个坐标。 cx = x[:-1] + x[1:] 等于重复 x 但第一个和最后一个不能重复,因为已经包含!
猜你喜欢
  • 2019-08-21
  • 1970-01-01
  • 1970-01-01
  • 1970-01-01
  • 2019-08-11
  • 1970-01-01
  • 1970-01-01
  • 2023-03-27
  • 2010-09-13
相关资源
最近更新 更多