2016-12-08 155 views
1

我有一个函数,它使用一个生成器来循环大型二维python浮点坐标列表,以便创建表示坐标之间距离的整数平面列表。Cython - 计算二维坐标之间的距离数组

point_input = {"x": -8081441.0, "y": 5685214.0} 
output = [-8081441, 5685214] 

polyline_input = {"paths" : [[-8081441.0, 5685214.0], [-8081446.0, 5685216.0], [-8081442.0, 5685219.0], [-8081440.0, 5685211.0], [-8081441.0, 5685214.0]]} 
output = [[-8081441, 5685214, 5, -2, -4, -3, -2, 8, 1, -3]] 

polygon_input = {"rings" : [[-8081441.0, 5685214.0], [-8081446.0, 5685216.0], [-8081442.0, 5685219.0], [-8081440.0, 5685211.0], [-8081441.0, 5685214.0]]} 
output = [[-8081441, 5685214, 5, -2, -4, -3, -2, 8, 1, -3]] 

纯Python:

def geometry_to_distance(geometry, geometry_type): 
    def calculate_distance(coords): 
     iterator = iter(coords) 
     previous_x, previous_y = iterator.next() 
     yield int(previous_x) 
     yield int(previous_y) 
     for current_x, current_y in iterator: 
      yield int(previous_x - current_x) 
      yield int(previous_y - current_y) 
      previous_x, previous_y = current_x, current_y 

    if geometry_type == "POINT": 
     distance_array = [int(geometry["x"]), int(geometry["y"])] 
    elif geometry_type == "POLYLINE": 
     distance_array = [list(calculate_distance(path)) for path in geometry["paths"]] 
    elif geometry_type == "POLYGON": 
     distance_array = [list(calculate_distance(ring)) for ring in geometry["rings"]] 
    else: 
     raise Exception("{} geometry type not supported".format(geometry_type)) 

    return distance_array 

对于速度的表现,我想用同样的功能的实现用Cython。我在calculate_distance函数中使用整型变量的类型声明。

用Cython实现:

def geometry_to_distance(geometry, geometry_type): 
    def calculate_distance(coords): 
     cdef int previous_x, previous_y, current_x, current_y 
     iterator = iter(coords) 
     previous_x, previous_y = iterator.next() 
     yield previous_x 
     yield previous_y 
     for current_x, current_y in iterator: 
      yield previous_x - current_x 
      yield previous_y - current_y 
      previous_x, previous_y = current_x, current_y 

    if geometry_type == "POINT": 
     distance_array = [geometry["x"], geometry["y"]] 
    elif geometry_type == "POLYLINE": 
     distance_array = [list(calculate_distance(path)) for path in geometry["paths"]] 
    elif geometry_type == "POLYGON": 
     distance_array = [list(calculate_distance(ring)) for ring in geometry["rings"]] 
    else: 
     raise Exception("{} geometry type not supported".format(geometry_type)) 

    return distance_array 

这里可以用来基准功能的脚本:

import time 
from functools import wraps 
import numpy as np 
import geometry_converter as gc 

def timethis(func): 
    '''Decorator that reports the execution time.''' 
    @wraps(func) 
    def wrapper(*args, **kwargs): 
     start = time.time() 
     result = func(*args, **kwargs) 
     end = time.time() 
     print(func.__name__, end-start) 
     return result 
    return wrapper 


def prepare_data(featCount, size): 
    ''' Create arrays of polygon geometry (see polygon_input above)''' 
    input = [] 
    for i in xrange(0, featCount): 
     polygon = {"rings" : []} 
     #random x,y coordinates inside a quadrant of the world bounding box in a spherical mercator (epsg:3857) projection 
     ys = np.random.uniform(-20037507.0,0,size).tolist() 
     xs = np.random.uniform(0,20037507.0,size).tolist() 
     polygon["rings"].append(zip(xs,ys)) 
     input.append(polygon) 
    return input 

@timethis 
def process_data(data): 
    output = [gc.esriJson_to_CV(x, "POLYGON") for x in data] 
    return output 

data = prepare_data(100, 100000) 
process_data(data) 

是否有改进,可在用Cython实现提高性能?也许通过使用2D cython数组或carrays?

+0

为什么不使用'numpy.diff'来获取X和Y坐标的第一个差值? – pbreach

+0

因为看起来从巨大的2D Python列表创建numpy.array太慢了。 –

+0

你也会遇到和cython或c数组一样的问题。列表不存储在连续的内存中,而(同类)numpy,cython和c数组。因此,不管这些方法如何,转换都需要一些时间。我很惊讶'numpy.diff'并不比使用生成器和列表的cython实现快。 – pbreach

回答

1

Python的,没有发电机改写,是

In [362]: polyline_input = {"paths" : [[-8081441.0, 5685214.0], [-8081446.0, 568 
    ...: 5216.0], [-8081442.0, 5685219.0], [-8081440.0, 5685211.0], [-8081441.0 
    ...: , 5685214.0]]} 
In [363]: output=polyline_input['paths'][0][:] # copy 
In [364]: i0,j0 = output 
    ...: for i,j in polyline_input['paths'][1:]: 
    ...:  output.extend([i0-i, j0-j][:]) 
    ...:  i0,j0 = i,j 
    ...:  
In [365]: output 
Out[365]: [-8081441.0, 5685214.0, 5.0, -2.0, -4.0, -3.0, -2.0, 8.0, 1.0, -3.0] 

我只是想表达,虽然计算的替代方法。我本可以使用append来取代平面列表的配对列表。

阵列当量:

In [375]: arr=np.array(polyline_input['paths']) 
In [376]: arr[1:,:]=arr[:-1,:]-arr[1:,:] 
In [377]: arr.ravel().tolist() 
Out[377]: [-8081441.0, 5685214.0, 5.0, -2.0, -4.0, -3.0, -2.0, 8.0, 1.0, -3.0] 

忽略列表转换为阵列的成本,看起来像一个高效numpy的操作。为了在cython中改进它,我希望你不得不将数组转换为memoryview,并且在值对上迭代c样式。

我忘了你为什么要切换到这种距离格式。你想保存一些文件空间吗?或者加快一些下游计算?

+0

谢谢您的回答!是的,目标是获得几何图形的轻量级表示以提高网络共享的速度 –