расстояние от точки до полигона (когда внутри)

Question

Я пытаюсь найти эффективный способ вычисления расстояния от точки до ближайшего края многоугольника в python. Я думал, что стройное было бы идеально для этого, но оно вычисляет только расстояние, когда точка находится за пределами многоугольника. Я думал о opencv, который имеет встроенную функцию, но для этого требуются целые координаты (я не могу дискретировать свои данные). Полигон не гарантированно выпуклый.

Answer 1

Используйте свойство [Multi] Polygon для возврата [Multi] LineString LinearRings, а затем используйте свойство .distance. Вы также можете использовать свойство .exterior, если хотите только наружное кольцо, а не любые внутренние кольца (если они существуют). Например.

from shapely import wkt 
poly = wkt.loads('POLYGON((30 10, 40 40, 20 40, 10 20, 30 10))') 
pt = wkt.loads('POINT(20 20)') 

# The point is in the polygon, so the distance will always be 0.0 
poly.distance(pt) # 0.0 

# The distance is calculated to the nearest boundary 
poly.boundary.distance(pt) # 4.47213595499958 

# or exterior ring 
poly.exterior.distance(pt) # 4.47213595499958 

В результате то же самое для .boundary и .exterior с этим примером, но может изменяться для других примеров.

Answer 2

У меня есть код для получения пересечения линий и полигонов (заданных ребрами) в cython. Если вы не хотите ждать здесь, это так. Я собираюсь спать сейчас, поэтому его единственный сырой код границы корабля судна и пересечения ватерлинии. Должно реализовать Что необходимо на Вашем собственном на сегодня: P ...

import numpy as np 
import pandas as pd 
import matplotlib.pyplot as plt 
import scipy.interpolate as sci_itp 
import scipy.integrate as sci_itg 
cimport numpy as np 
cimport cython 

from CROSS_SECTIONS_universal import * 

DTYPE = np.float64 
ctypedef np.float64_t DTYPE_t 

ITYPE = np.int 
ctypedef np.int_t ITYPE_t 

ITYPEP = np.intp 
ctypedef np.intp_t ITYPEP_t 

ITYPE1 = np.int64 
ctypedef np.int64_t ITYPE1_t 

BTYPE = np.uint8 
ctypedef np.uint8_t BTYPE_t 

ctypedef Py_ssize_t shape_t 

BOTYPE = np.bool_ 

cdef double[:,:] waterLine = np.array([[0,0],[0,0]],dtype=DTYPE) 
cdef double[:,:] crossLine = np.array([[0,0],[0,0]],dtype=DTYPE) 
cdef double[:,:] line 
cdef double[:] sectionLinePoint = np.full((2),np.nan,dtype=DTYPE) 
cdef double[:,:] sectionLengthPOINTS = np.full((2,2),np.nan,dtype=DTYPE) 
cdef double[:,:] sectionBreadthPOINTS = np.full((2,2),np.nan,dtype=DTYPE) 
cdef double sectionLength,sectionBreadth 


cdef point_side_of_line(double ax,double ay,double bx,double by,double cx,double cy): 
    """ Returns a position of the point c relative to the line going through a and b 
     Points a, b are expected to be different 
    """ 
    cdef double s10x,s10y,s20x,s20y,denom 
    cdef bint denomIsPositive,denomIsNegative 
    cdef int side 

    s10x,s10y=bx-ax,by-ay 
    s20x,s20y=cx-ax,cy-ay 

    denom = s20y*s10x - s10y*s20x 
    denomIsPositive,denomIsNegative=denom>0,denom<0 

    side=1 if denomIsPositive else (-1 if denomIsNegative else 0) 
    return side 

cdef is_point_in_closed_segment(double ax,double ay,double bx,double by,double cx,double cy): 
    """ Returns True if c is inside closed segment, False otherwise. 
     a, b, c are expected to be collinear 
    """ 
    cdef bint pointIn 

    if ax < bx: 
     pointIn=ax <= cx and cx <= bx 
     return pointIn 
    if bx < ax: 
     pointIn=bx <= cx and cx <= ax 
     return pointIn 
    if ay < by: 
     pointIn=ay <= cy and cy <= by 
     return pointIn 
    if by < ay: 
     pointIn=by <= cy and cy <= ay 
     return pointIn 

    pointIn=ax == cx and ay == cy 
    return pointIn 

cdef closed_segment_intersect(double ax,double ay,double bx,double by,double cx,double cy,double dx,double dy): 
    """ Verifies if closed segments a, b, c, d do intersect. 
    """ 
    cdef int s1,s2 

    if (ax is bx) and (ay is by): 
     assert (ax is cx and ay is cy) or (ax is dx and ay is dy) 
    if (cx is dx) and (cy is dy): 
     assert (cx is ax and cy is ay) or (cx is bx and cy is by) 

    s1,s2 = point_side_of_line(ax,ay,bx,by,cx,cy),point_side_of_line(ax,ay,bx,by,dx,dy) 

    # All points are collinear 
    if s1 is 0 and s2 is 0: 
     assert \ 
      is_point_in_closed_segment(ax,ay,bx,by,cx,cy) or \ 
      is_point_in_closed_segment(ax,ay,bx,by,dx,dy) or \ 
      is_point_in_closed_segment(cx,cy,dx,dy,ax,ay) or \ 
      is_point_in_closed_segment(cx,cy,dx,dy,bx,by) 


    # No touching and on the same side 
    if s1 and s1 is s2: 
     assert False 

    s1,s2 = point_side_of_line(cx,cy,dx,dy,ax,ay),point_side_of_line(cx,cy,dx,dy,bx,by) 

    # No touching and on the same side 
    if s1 and s1 is s2: 
     assert False 

cdef find_intersection(double[:,:] L1,double[:,:] L2, double[:] intersectionPoint) : 
    cdef double ax,ay,bx,by,cx,cy,dx,dy 
    cdef double s10x,s10y,s32x,s32y,s02x,s02y 
    cdef double tNumer,t,denom 


    ax,ay=L1[0] 
    bx,by=L1[1] 
    cx,cy=L2[0] 
    dx,dy=L2[1] 

    closed_segment_intersect(ax,ay,bx,by,cx,cy,dx,dy) 

    s10x,s10y=bx-ax,by-ay 
    s02x,s02y=ax-cx,ay-cy 
    s32x,s32y=dx-cx,dy-cy 

    denom = s10x * s32y - s32x * s10y 
    tNumer = s32x * s02y - s32y * s02x 

    t = tNumer/denom 
    intersectionPoint[0]=ax+t*s10x 
    intersectionPoint[1]=ay+t*s10y 

cdef section_length_breadth(double[:,:] planeSectionCLUSTER,double x0,double x1, double y0, double y1): 

    cdef int counterL=0,counterB=0 

    x01=(x0+x1)/2. 

    crossLine[0,0]=x01 
    crossLine[1,0]=x01 
    crossLine[0,1]=y0-0.1 
    crossLine[1,1]=y1+0.1 
    waterLine[0,0]=x0-0.1 
    waterLine[1,0]=x1+0.1 


    for i in range(1,len(planeSectionCLUSTER)): 
     line=planeSectionCLUSTER[i-1:i+1] 
     plt.plot(line[:,0],line[:,1],'c-',lw=3) 
     try: 
      try: 
       find_intersection(line,waterLine,sectionLinePoint) 
       assert sectionLinePoint[0] not in sectionLengthPOINTS[:,0] 
       sectionLengthPOINTS[counterL]=sectionLinePoint 
       counterL+=1 
      except AssertionError:   
       find_intersection(line,crossLine,sectionLinePoint) 
       assert sectionLinePoint[1] not in sectionBreadthPOINTS[:,1] 
       sectionBreadthPOINTS[counterB]=sectionLinePoint 
       counterB+=1 
     except AssertionError: 
      pass  

    print '#' 
    print [[j for j in i] for i in sectionLengthPOINTS] 
    print [[j for j in i] for i in sectionBreadthPOINTS] 

    plt.plot(waterLine[:,0],waterLine[:,1],'b--',lw=3) 
    plt.plot(sectionLengthPOINTS[:,0],sectionLengthPOINTS[:,1],'ko-',lw=3) 
    plt.plot(crossLine[:,0],crossLine[:,1],'b--',lw=3) 
    plt.plot(sectionBreadthPOINTS[:,0],sectionBreadthPOINTS[:,1],'ro-',lw=3) 
    return 0,0