Vectorized functions

src/h3/_cy/CMakeLists.txt

@@ -19,3 +19,4 @@ add_cython_file(geo)
 add_cython_file(cells)
 add_cython_file(edges)
 add_cython_file(to_multipoly)
+add_cython_file(unstable_vect)

src/h3/_cy/h3lib.pxd

@@ -52,9 +52,9 @@ cdef extern from "h3api.h":
         int i
         int j
 
-    H3Index geoToH3(const GeoCoord *g, int res)
+    H3Index geoToH3(const GeoCoord *g, int res) nogil
 
-    void h3ToGeo(H3Index h3, GeoCoord *g)
+    void h3ToGeo(H3Index h3, GeoCoord *g) nogil
 
     void h3ToGeoBoundary(H3Index h3, GeoBoundary *gp)
 
@@ -82,9 +82,9 @@ cdef extern from "h3api.h":
 
     void destroyLinkedPolygon(LinkedGeoPolygon *polygon)
 
-    double degsToRads(double degrees)
+    double degsToRads(double degrees) nogil
 
-    double radsToDegs(double radians)
+    double radsToDegs(double radians) nogil
 
     double hexAreaKm2(int res)
 

src/h3/_cy/unstable_vect.pyx

@@ -0,0 +1,59 @@
+cimport h3lib
+from h3lib cimport H3int
+from .util cimport deg2coord
+
+from cython cimport boundscheck, wraparound
+from libc.math cimport sqrt, sin, cos, asin
+
+cdef double haversineDistance(double th1, double ph1, double th2, double ph2) nogil:
+    cdef:
+        double dx, dy, dz
+        double R = 6371.0088
+
+    ph1 -= ph2
+
+    dz = sin(th1) - sin(th2)
+    dx = cos(ph1) * cos(th1) - cos(th2)
+    dy = sin(ph1) * cos(th1)
+
+    return asin(sqrt(dx*dx + dy*dy + dz*dz) / 2)*2*R
+
+
+@boundscheck(False)
+@wraparound(False)
+cpdef void haversine_vect(
+    const H3int[:] a,
+    const H3int[:] b,
+         double[:] out
+) nogil:
+
+    cdef h3lib.GeoCoord p1, p2
+
+    with nogil:
+        # todo: add these back in when cython 3.0 comes out
+        #assert len(a) == len(b)
+        #assert len(a) <= len(out)
+
+        for i in range(len(a)):
+            h3lib.h3ToGeo(a[i], &p1)
+            h3lib.h3ToGeo(b[i], &p2)
+            out[i] = haversineDistance(
+                p1.lat, p1.lng,
+                p2.lat, p2.lng
+            )
+
+@boundscheck(False)
+@wraparound(False)
+cpdef void geo_to_h3_vect(
+    const double[:] lat,
+    const double[:] lng,
+    int res,
+    H3int[:] out
+) nogil:
+
+    cdef h3lib.GeoCoord c
+
+    with nogil:
+        for i in range(len(lat)):
+            c = deg2coord(lat[i], lng[i])
+            out[i] = h3lib.geoToH3(&c, res)

src/h3/_cy/util.pxd

@@ -1,7 +1,7 @@
 from .h3lib cimport H3int, H3str, GeoCoord
 
-cdef GeoCoord deg2coord(double lat, double lng)
-cdef (double, double) coord2deg(GeoCoord c)
+cdef GeoCoord deg2coord(double lat, double lng) nogil
+cdef (double, double) coord2deg(GeoCoord c) nogil
 
 cpdef H3int hex2int(H3str h) except? 0
 cpdef H3str int2hex(H3int x)

src/h3/_cy/util.pyx

@@ -5,7 +5,7 @@ from .h3lib cimport H3int, H3str, h3IsValid, h3UnidirectionalEdgeIsValid
 cimport h3lib
 
 
-cdef h3lib.GeoCoord deg2coord(double lat, double lng):
+cdef h3lib.GeoCoord deg2coord(double lat, double lng) nogil:
     cdef:
         h3lib.GeoCoord c
 
@@ -15,7 +15,7 @@ cdef h3lib.GeoCoord deg2coord(double lat, double lng):
     return c
 
 
-cdef (double, double) coord2deg(h3lib.GeoCoord c):
+cdef (double, double) coord2deg(h3lib.GeoCoord c) nogil:
     return (
         h3lib.radsToDegs(c.lat),
         h3lib.radsToDegs(c.lng)

src/h3/unstable/__init__.py

@@ -6,3 +6,4 @@
 )
 
 from . import v4
+from . import vect

src/h3/unstable/vect.py

@@ -0,0 +1,52 @@
+from .._cy import unstable_vect as _vect
+
+import numpy as np
+
+
+def geo_to_h3(lats, lngs, res):
+    """
+    Convert arrays describing lat/lng paris to cells.
+
+    Parameters
+    ----------
+    lats, lngs : arrays of floats
+
+    res: int
+        Resolution for output cells.
+
+    Returns
+    -------
+    array of H3Cells
+    """
+    assert len(lats) == len(lngs)
+
+    out = np.zeros(len(lats), dtype='uint64')
+
+    _vect.geo_to_h3_vect(lats, lngs, res, out)
+
+    return out
+
+
+def cell_haversine(a, b):
+    """
+    Compute haversine distance between the centers of cells given in
+    arrays a and b.
+
+
+    Parameters
+    ----------
+    a, b : arrays of H3Cell
+
+    Returns
+    -------
+    float
+        Haversine distance in kilometers
+    """
+
+    assert len(a) == len(b)
+
+    out = np.zeros(len(a), dtype='double')
+
+    _vect.haversine_vect(a, b, out)
+
+    return out