Двоичные / унарные функции предикатов, сравнивающие ВСЕ объекты со ВСЕМИ другими объектами в python - PullRequest
Новая
       39

Двоичные / унарные функции предикатов, сравнивающие ВСЕ объекты со ВСЕМИ другими объектами в python

0 голосов
/ 08 июля 2019

Я задавал очень похожий вопрос раньше. Поскольку решение arcpy очень громоздкое, сейчас я ищу в основном ту же функцию в geopandas. Вопрос в том, каков самый быстрый / лучший способ применения двоичной функции предикатов ( например, touches), где каждая функция x сравнивается с любой другой особенность либо x, либо другого набора данных y. Я ожидал бы вывод, аналогичный поведению по умолчанию в R:

Если y отсутствует, эффективно вызывается st_predicate(x, x) и возвращается квадратная матрица с диагональными элементами st_predicate(x[i], x[i]).

В качестве примера приведем некоторые фиктивные данные и функцию st_overlaps(): 

library(sf)

b0 = st_polygon(list(rbind(c(-1,-1), c(1,-1), c(1,1), c(-1,1), c(-1,-1))))
a0 = b0 * 0.8
a1 = a0 * 0.5 + c(2, 0.7)
a2 = a0 + 1
a3 = b0 * 0.5 + c(2, -0.5)
x = st_sfc(a0,a1,a2,a3)

plot(x)


st_overlaps(x)
#> Sparse geometry binary predicate list of length 4, where the predicate was `overlaps'
#>  1: 3
#>  2: 3
#>  3: 1, 2
#>  4: (empty)

Как мне добиться подобного поведения в python / geopandas? Очевидно, geopandas автоматически выравнивает x и x / y, и по элементам выполняется сравнение (см. этот вопрос SO и эту проблему на github ). В python выполнение x.overlaps(x) просто возвращает серию панд с четырьмя значениями True. 

import geopandas as gpd

x.overlaps(x)
0      True
1      True
2      True
3      True

Ответы [ 2 ]

1 голос
/ 08 июля 2019

Это определенно не самый быстрый способ, так как это простой итератор, но если ваши данные не велики, это может сработать. 

import geopandas as gpd
from shapely.geometry import Polygon

b0 = Polygon([(-1,-1), (1,-1), (1,1), (-1,1)])
a1 = Polygon([(1.5,0.2), (2.5,0.2), (2.5,1.2), (1.5,1.2)])
a2 = Polygon([(0,0), (2,0), (2,2), (0,2)])
a3 = Polygon([(1.5,-1), (2.5,-1), (2.5,-0.2), (1.5,-0.2)])

series = gpd.GeoSeries([b0, a1, a2, a3])

results = {}
for poly in series.iteritems():
    results[poly[0]] = []
    for poly2 in series.drop(poly[0]).iteritems():
        if poly[1].overlaps(poly2[1]):
            results[poly[0]].append(poly2[0])

Это даст вам контроль над вашими ценностями. 

{0: [2], 1: [2], 2: [0, 1], 3: []}

Однако имейте в виду, что он проверяет A-> B, а затем B-> A и что он также проверяет полигоны, даже если они явно находятся далеко. Чтобы ускорить его, вы можете использовать пространственный индекс rtree, чтобы проверять только тех, кто может перекрываться, вместо проверки каждого многоугольника друг на друга (дважды).

0 голосов
/ 09 июля 2019

Идиоматический способ Python выразить это заключается в понимании списка, например, чтобы создать список, состоящий из кортежа (index: (overlapping-index)), вы должны написать 

[ ( ind, 
    [ind2 for ind2, g2 in enumerate(series) if g.overlaps(g2)] 
  ) 
  for ind, g in enumerate(series) ]

Результат: 

[(0, [2]), (1, [2]), (2, [0, 1]), (3, [])]

Но, как отметил Мартинфлейс, это не очень эффективный способ сделать это, так как он не использует никакой пространственной индексации.

Вы можете повысить производительность, используя операцию наложения, см. http://geopandas.org/set_operations.html

...