Как рассчитать площадь самопересекающихся многоугольников на Земле

Question

Я пытаюсь найти область самопересекающихся (сложных) многоугольников на Земле. Есть ли библиотека, которая правильно выполняет вычисления площадей для таких геометрий?

До сих пор я пробовал geographiclib и Polygon2 (на основе gpc), но они дают неправильные результаты для сложных полигонов.

В альтернативном варианте есть простой способ преобразования самопересекающихся многоугольников в множество простых многоугольников, так что я могу использовать алгоритмы для простых геометрий, а затем суммировать все полигоны? Я знаю, что решением было бы реализовать Bentley-Ottman algorithm, чтобы найти точки пересечения, а затем разделить многоугольник на множество простых полигонов, но если есть библиотека, которую я мог бы использовать, я бы с удовольствием избегал изобретать колесо (и, возможно, вводить ошибки).

Обновление: Координаты, которые я использую, не содержат явные точки пересечения. Поэтому библиотеки отсечения, похоже, не обрабатывают их должным образом. Например, в следующем коде coords_1 и coords_2 определяют один и тот же полигон. Использование библиотеки обрезки Polygon2.

>>> coords_1 = [(0,0), (1,1), (1,0), (0,1)] 
>>> coords_2 = [(0,0), (0.5, 0.5), (1,0), (1,1), (0.5,0.5), (0,1)] 
>>> Polygon(coords_1).area() 
0 
>>> Polygon(coords_2).area() 
0.5 

Я бы хотел получить второй результат, используя coords_1. Я также пробовал другие библиотеки отсечения, но пока не повезло.

Answer 1

Приблизительное решение было бы отображать многоугольники в плоскость и использовать библиотеку геометрии плоскости для выполнения вычислительной геометрии. Это приближение должно работать очень хорошо для многоугольников, которые малы относительно всей сферы; но даже для больших полигонов вы можете быть произвольно точными, достаточно усовершенствуя свой полигон.

Я не думаю, что, к сожалению, существует большая функциональность библиотеки для операций с сферической сферой.

Edit: здесь есть стройное решение, которое я думаю, что делает то, что вы хотите:

from shapely.geometry.polygon import LinearRing, Polygon 

coords_1 = [(0,0), (1,1), (1,0), (0,1)] 
##coords_1 = [(0,0), (0,1), (1,1), (1,0)] 

lr = LinearRing(coords_1) 

if lr.is_valid: 
    print Polygon(lr).area 
else: 
    print Polygon(coords_1).buffer(0).area + Polygon(coords_1[::-1]).buffer(0).area 

Answer 2

Преобразования широты/долготы многоугольники, чтобы играть хорошо на евклидовой плоскости. Затем используйте что-то вроде: http://sourceforge.net/projects/polyclipping/.