diff --git a/doc/source/api.rst b/doc/source/api.rst index 64bad58..0d27ae4 100644 --- a/doc/source/api.rst +++ b/doc/source/api.rst @@ -12,6 +12,6 @@ Indexing .. autosummary:: :toctree: generated/ - ShapelySTRTreeIndex.__init__ - ShapelySTRTreeIndex.from_variables - ShapelySTRTreeIndex.sel + GeometryIndex + GeometryIndex.crs + GeometryIndex.sindex diff --git a/doc/source/conf.py b/doc/source/conf.py index da50322..a49b0da 100644 --- a/doc/source/conf.py +++ b/doc/source/conf.py @@ -17,6 +17,7 @@ "sphinx.ext.autodoc", "numpydoc", "sphinx.ext.autosummary", + "sphinx.ext.intersphinx", "myst_nb", "sphinx_copybutton", ] @@ -24,6 +25,12 @@ templates_path = ["_templates"] exclude_patterns = [] +intersphinx_mapping = { + "python": ("https://docs.python.org/3", None), + "shapely": ("https://shapely.readthedocs.io/en/latest/", None), + "pyproj": ("https://pyproj4.github.io/pyproj/latest/", None), + "xarray": ("https://docs.xarray.dev/en/latest/", None), +} # -- Options for HTML output ------------------------------------------------- # https://www.sphinx-doc.org/en/master/usage/configuration.html#options-for-html-output diff --git a/doc/source/quickstart.ipynb b/doc/source/quickstart.ipynb index d1b2068..9c9ec83 100644 --- a/doc/source/quickstart.ipynb +++ b/doc/source/quickstart.ipynb @@ -30,7 +30,7 @@ "import xarray\n", "import numpy\n", "\n", - "from xvec import ShapelySTRTreeIndex" + "from xvec import GeometryIndex" ] }, { @@ -82,8 +82,19 @@ "source": [ "mode = [\"car\", \"bike\", \"foot\"]\n", "day = pandas.date_range(\"2015-01-01\", periods=100)\n", - "hours = range(24)\n", - "data = numpy.random.randint(1, 100, size=(3, 100, 24, 100, 100))" + "hours = range(24)" + ] + }, + { + "cell_type": "code", + "execution_count": 5, + "metadata": {}, + "outputs": [], + "source": [ + "from numpy.random import default_rng\n", + "\n", + "rng = default_rng(1)\n", + "data = rng.integers(1, 100, size=(3, 100, 24, 100, 100))" ] }, { @@ -95,7 +106,7 @@ }, { "cell_type": "code", - "execution_count": 5, + "execution_count": 6, "metadata": {}, "outputs": [ { @@ -465,93 +476,93 @@ " fill: currentColor;\n", "}\n", "
<xarray.DataArray (mode: 3, day: 100, time: 24, origin: 100, destination: 100)>\n",
-       "array([[[[[96, 29, 64, ..., 47, 80, 71],\n",
-       "          [79, 79, 96, ..., 92, 12, 90],\n",
-       "          [19, 85, 60, ..., 43, 28, 53],\n",
+       "array([[[[[47, 51, 75, ..., 15, 54, 82],\n",
+       "          [ 7, 68, 76, ..., 38,  3, 72],\n",
+       "          [50, 65, 16, ..., 53, 91, 15],\n",
        "          ...,\n",
-       "          [96,  8, 35, ..., 39,  7, 62],\n",
-       "          [82, 87, 18, ...,  1, 77, 38],\n",
-       "          [99, 94, 20, ..., 92, 16, 42]],\n",
+       "          [30,  1, 15, ...,  6, 36, 46],\n",
+       "          [19, 47,  3, ..., 75, 42, 69],\n",
+       "          [48, 94, 34, ..., 47, 54, 86]],\n",
        "\n",
-       "         [[79, 53, 28, ..., 30, 56, 14],\n",
-       "          [21, 88, 91, ..., 92, 73, 14],\n",
-       "          [74, 49, 88, ..., 30, 23, 28],\n",
+       "         [[67,  5, 73, ..., 95, 34, 44],\n",
+       "          [80, 22, 89, ..., 93,  5, 41],\n",
+       "          [47, 70, 65, ..., 44, 82, 54],\n",
        "          ...,\n",
-       "          [56, 76, 91, ..., 36, 94, 16],\n",
-       "          [80, 83,  4, ..., 91, 38,  3],\n",
-       "          [69, 61, 10, ..., 79, 15, 94]],\n",
+       "          [91, 33, 38, ..., 24, 40,  5],\n",
+       "          [78, 49,  9, ..., 45,  3,  9],\n",
+       "          [87, 73, 52, ..., 32, 58, 48]],\n",
        "\n",
-       "         [[91, 22, 13, ..., 45, 31, 17],\n",
-       "          [86, 48, 23, ..., 88, 55, 94],\n",
-       "          [56,  8, 55, ..., 68, 15, 48],\n",
+       "         [[22, 57, 42, ..., 24, 60,  1],\n",
+       "          [ 4,  3, 94, ..., 66, 40, 17],\n",
+       "          [87, 55, 87, ..., 41, 26, 55],\n",
        "          ...,\n",
        "...\n",
        "          ...,\n",
-       "          [27,  8, 29, ..., 96, 80, 60],\n",
-       "          [ 1, 98, 63, ...,  1, 73, 98],\n",
-       "          [77, 22, 88, ...,  6, 60,  7]],\n",
+       "          [67, 50, 96, ..., 28, 89, 55],\n",
+       "          [83, 98,  3, ..., 91, 25, 14],\n",
+       "          [81,  7, 45, ..., 18, 78,  9]],\n",
        "\n",
-       "         [[22, 40, 31, ..., 84, 62, 61],\n",
-       "          [27, 83, 95, ..., 81, 15, 54],\n",
-       "          [43, 45, 67, ..., 52,  2, 28],\n",
+       "         [[29, 94, 53, ...,  3, 52, 47],\n",
+       "          [93, 43, 11, ..., 94, 71, 37],\n",
+       "          [80, 34, 29, ..., 12, 67, 48],\n",
        "          ...,\n",
-       "          [91, 75, 49, ...,  1, 76, 65],\n",
-       "          [ 4, 78, 43, ..., 82, 20, 95],\n",
-       "          [16, 65, 36, ..., 27, 92, 43]],\n",
+       "          [31, 96, 58, ..., 71, 85, 90],\n",
+       "          [72, 95, 55, ..., 83, 17, 97],\n",
+       "          [72, 77, 15, ..., 42, 62,  9]],\n",
        "\n",
-       "         [[11, 54, 89, ..., 64, 56, 78],\n",
-       "          [ 1, 29, 49, ..., 46, 22, 95],\n",
-       "          [99,  4, 33, ..., 89, 48, 60],\n",
+       "         [[26, 76, 23, ..., 21, 92, 42],\n",
+       "          [60, 80, 37, ..., 49, 91, 88],\n",
+       "          [55, 64, 94, ..., 99, 91, 43],\n",
        "          ...,\n",
-       "          [18, 62, 84, ..., 12, 87, 52],\n",
-       "          [42, 79, 90, ...,  5, 24, 57],\n",
-       "          [51, 35, 71, ..., 55, 29, 11]]]]])\n",
+       "          [80, 85, 15, ..., 23, 42, 52],\n",
+       "          [66, 91, 73, ..., 10, 44, 22],\n",
+       "          [75, 27, 89, ..., 37, 57, 20]]]]])\n",
        "Coordinates:\n",
        "  * mode         (mode) <U4 'car' 'bike' 'foot'\n",
        "  * day          (day) datetime64[ns] 2015-01-01 2015-01-02 ... 2015-04-10\n",
        "  * time         (time) int64 0 1 2 3 4 5 6 7 8 9 ... 15 16 17 18 19 20 21 22 23\n",
        "  * origin       (origin) object MULTIPOLYGON (((-81.4727554321289 36.2343559...\n",
-       "  * destination  (destination) object MULTIPOLYGON (((-81.4727554321289 36.23...
    • " ], "text/plain": [ "\n", - "array([[[[[96, 29, 64, ..., 47, 80, 71],\n", - " [79, 79, 96, ..., 92, 12, 90],\n", - " [19, 85, 60, ..., 43, 28, 53],\n", + "array([[[[[47, 51, 75, ..., 15, 54, 82],\n", + " [ 7, 68, 76, ..., 38, 3, 72],\n", + " [50, 65, 16, ..., 53, 91, 15],\n", " ...,\n", - " [96, 8, 35, ..., 39, 7, 62],\n", - " [82, 87, 18, ..., 1, 77, 38],\n", - " [99, 94, 20, ..., 92, 16, 42]],\n", + " [30, 1, 15, ..., 6, 36, 46],\n", + " [19, 47, 3, ..., 75, 42, 69],\n", + " [48, 94, 34, ..., 47, 54, 86]],\n", "\n", - " [[79, 53, 28, ..., 30, 56, 14],\n", - " [21, 88, 91, ..., 92, 73, 14],\n", - " [74, 49, 88, ..., 30, 23, 28],\n", + " [[67, 5, 73, ..., 95, 34, 44],\n", + " [80, 22, 89, ..., 93, 5, 41],\n", + " [47, 70, 65, ..., 44, 82, 54],\n", " ...,\n", - " [56, 76, 91, ..., 36, 94, 16],\n", - " [80, 83, 4, ..., 91, 38, 3],\n", - " [69, 61, 10, ..., 79, 15, 94]],\n", + " [91, 33, 38, ..., 24, 40, 5],\n", + " [78, 49, 9, ..., 45, 3, 9],\n", + " [87, 73, 52, ..., 32, 58, 48]],\n", "\n", - " [[91, 22, 13, ..., 45, 31, 17],\n", - " [86, 48, 23, ..., 88, 55, 94],\n", - " [56, 8, 55, ..., 68, 15, 48],\n", + " [[22, 57, 42, ..., 24, 60, 1],\n", + " [ 4, 3, 94, ..., 66, 40, 17],\n", + " [87, 55, 87, ..., 41, 26, 55],\n", " ...,\n", "...\n", " ...,\n", - " [27, 8, 29, ..., 96, 80, 60],\n", - " [ 1, 98, 63, ..., 1, 73, 98],\n", - " [77, 22, 88, ..., 6, 60, 7]],\n", + " [67, 50, 96, ..., 28, 89, 55],\n", + " [83, 98, 3, ..., 91, 25, 14],\n", + " [81, 7, 45, ..., 18, 78, 9]],\n", "\n", - " [[22, 40, 31, ..., 84, 62, 61],\n", - " [27, 83, 95, ..., 81, 15, 54],\n", - " [43, 45, 67, ..., 52, 2, 28],\n", + " [[29, 94, 53, ..., 3, 52, 47],\n", + " [93, 43, 11, ..., 94, 71, 37],\n", + " [80, 34, 29, ..., 12, 67, 48],\n", " ...,\n", - " [91, 75, 49, ..., 1, 76, 65],\n", - " [ 4, 78, 43, ..., 82, 20, 95],\n", - " [16, 65, 36, ..., 27, 92, 43]],\n", + " [31, 96, 58, ..., 71, 85, 90],\n", + " [72, 95, 55, ..., 83, 17, 97],\n", + " [72, 77, 15, ..., 42, 62, 9]],\n", "\n", - " [[11, 54, 89, ..., 64, 56, 78],\n", - " [ 1, 29, 49, ..., 46, 22, 95],\n", - " [99, 4, 33, ..., 89, 48, 60],\n", + " [[26, 76, 23, ..., 21, 92, 42],\n", + " [60, 80, 37, ..., 49, 91, 88],\n", + " [55, 64, 94, ..., 99, 91, 43],\n", " ...,\n", - " [18, 62, 84, ..., 12, 87, 52],\n", - " [42, 79, 90, ..., 5, 24, 57],\n", - " [51, 35, 71, ..., 55, 29, 11]]]]])\n", + " [80, 85, 15, ..., 23, 42, 52],\n", + " [66, 91, 73, ..., 10, 44, 22],\n", + " [75, 27, 89, ..., 37, 57, 20]]]]])\n", "Coordinates:\n", " * mode (mode) 
    <xarray.DataArray (mode: 3, day: 100, time: 24, origin: 100, destination: 100)>\n",
-       "array([[[[[96, 29, 64, ..., 47, 80, 71],\n",
-       "          [79, 79, 96, ..., 92, 12, 90],\n",
-       "          [19, 85, 60, ..., 43, 28, 53],\n",
+       "array([[[[[47, 51, 75, ..., 15, 54, 82],\n",
+       "          [ 7, 68, 76, ..., 38,  3, 72],\n",
+       "          [50, 65, 16, ..., 53, 91, 15],\n",
        "          ...,\n",
-       "          [96,  8, 35, ..., 39,  7, 62],\n",
-       "          [82, 87, 18, ...,  1, 77, 38],\n",
-       "          [99, 94, 20, ..., 92, 16, 42]],\n",
+       "          [30,  1, 15, ...,  6, 36, 46],\n",
+       "          [19, 47,  3, ..., 75, 42, 69],\n",
+       "          [48, 94, 34, ..., 47, 54, 86]],\n",
        "\n",
-       "         [[79, 53, 28, ..., 30, 56, 14],\n",
-       "          [21, 88, 91, ..., 92, 73, 14],\n",
-       "          [74, 49, 88, ..., 30, 23, 28],\n",
+       "         [[67,  5, 73, ..., 95, 34, 44],\n",
+       "          [80, 22, 89, ..., 93,  5, 41],\n",
+       "          [47, 70, 65, ..., 44, 82, 54],\n",
        "          ...,\n",
-       "          [56, 76, 91, ..., 36, 94, 16],\n",
-       "          [80, 83,  4, ..., 91, 38,  3],\n",
-       "          [69, 61, 10, ..., 79, 15, 94]],\n",
+       "          [91, 33, 38, ..., 24, 40,  5],\n",
+       "          [78, 49,  9, ..., 45,  3,  9],\n",
+       "          [87, 73, 52, ..., 32, 58, 48]],\n",
        "\n",
-       "         [[91, 22, 13, ..., 45, 31, 17],\n",
-       "          [86, 48, 23, ..., 88, 55, 94],\n",
-       "          [56,  8, 55, ..., 68, 15, 48],\n",
+       "         [[22, 57, 42, ..., 24, 60,  1],\n",
+       "          [ 4,  3, 94, ..., 66, 40, 17],\n",
+       "          [87, 55, 87, ..., 41, 26, 55],\n",
        "          ...,\n",
        "...\n",
        "          ...,\n",
-       "          [27,  8, 29, ..., 96, 80, 60],\n",
-       "          [ 1, 98, 63, ...,  1, 73, 98],\n",
-       "          [77, 22, 88, ...,  6, 60,  7]],\n",
+       "          [67, 50, 96, ..., 28, 89, 55],\n",
+       "          [83, 98,  3, ..., 91, 25, 14],\n",
+       "          [81,  7, 45, ..., 18, 78,  9]],\n",
        "\n",
-       "         [[22, 40, 31, ..., 84, 62, 61],\n",
-       "          [27, 83, 95, ..., 81, 15, 54],\n",
-       "          [43, 45, 67, ..., 52,  2, 28],\n",
+       "         [[29, 94, 53, ...,  3, 52, 47],\n",
+       "          [93, 43, 11, ..., 94, 71, 37],\n",
+       "          [80, 34, 29, ..., 12, 67, 48],\n",
        "          ...,\n",
-       "          [91, 75, 49, ...,  1, 76, 65],\n",
-       "          [ 4, 78, 43, ..., 82, 20, 95],\n",
-       "          [16, 65, 36, ..., 27, 92, 43]],\n",
+       "          [31, 96, 58, ..., 71, 85, 90],\n",
+       "          [72, 95, 55, ..., 83, 17, 97],\n",
+       "          [72, 77, 15, ..., 42, 62,  9]],\n",
        "\n",
-       "         [[11, 54, 89, ..., 64, 56, 78],\n",
-       "          [ 1, 29, 49, ..., 46, 22, 95],\n",
-       "          [99,  4, 33, ..., 89, 48, 60],\n",
+       "         [[26, 76, 23, ..., 21, 92, 42],\n",
+       "          [60, 80, 37, ..., 49, 91, 88],\n",
+       "          [55, 64, 94, ..., 99, 91, 43],\n",
        "          ...,\n",
-       "          [18, 62, 84, ..., 12, 87, 52],\n",
-       "          [42, 79, 90, ...,  5, 24, 57],\n",
-       "          [51, 35, 71, ..., 55, 29, 11]]]]])\n",
+       "          [80, 85, 15, ..., 23, 42, 52],\n",
+       "          [66, 91, 73, ..., 10, 44, 22],\n",
+       "          [75, 27, 89, ..., 37, 57, 20]]]]])\n",
        "Coordinates:\n",
        "  * mode         (mode) <U4 'car' 'bike' 'foot'\n",
        "  * day          (day) datetime64[ns] 2015-01-01 2015-01-02 ... 2015-04-10\n",
@@ -1519,48 +1530,48 @@
        "  * origin       (origin) object MULTIPOLYGON (((-81.4727554321289 36.2343559...\n",
        "  * destination  (destination) object MULTIPOLYGON (((-81.4727554321289 36.23...\n",
        "Indexes:\n",
-       "    origin       ShapelySTRTreeIndex\n",
-       "    destination  ShapelySTRTreeIndex
    • " ], "text/plain": [ "\n", - "array([[[[[96, 29, 64, ..., 47, 80, 71],\n", - " [79, 79, 96, ..., 92, 12, 90],\n", - " [19, 85, 60, ..., 43, 28, 53],\n", + "array([[[[[47, 51, 75, ..., 15, 54, 82],\n", + " [ 7, 68, 76, ..., 38, 3, 72],\n", + " [50, 65, 16, ..., 53, 91, 15],\n", " ...,\n", - " [96, 8, 35, ..., 39, 7, 62],\n", - " [82, 87, 18, ..., 1, 77, 38],\n", - " [99, 94, 20, ..., 92, 16, 42]],\n", + " [30, 1, 15, ..., 6, 36, 46],\n", + " [19, 47, 3, ..., 75, 42, 69],\n", + " [48, 94, 34, ..., 47, 54, 86]],\n", "\n", - " [[79, 53, 28, ..., 30, 56, 14],\n", - " [21, 88, 91, ..., 92, 73, 14],\n", - " [74, 49, 88, ..., 30, 23, 28],\n", + " [[67, 5, 73, ..., 95, 34, 44],\n", + " [80, 22, 89, ..., 93, 5, 41],\n", + " [47, 70, 65, ..., 44, 82, 54],\n", " ...,\n", - " [56, 76, 91, ..., 36, 94, 16],\n", - " [80, 83, 4, ..., 91, 38, 3],\n", - " [69, 61, 10, ..., 79, 15, 94]],\n", + " [91, 33, 38, ..., 24, 40, 5],\n", + " [78, 49, 9, ..., 45, 3, 9],\n", + " [87, 73, 52, ..., 32, 58, 48]],\n", "\n", - " [[91, 22, 13, ..., 45, 31, 17],\n", - " [86, 48, 23, ..., 88, 55, 94],\n", - " [56, 8, 55, ..., 68, 15, 48],\n", + " [[22, 57, 42, ..., 24, 60, 1],\n", + " [ 4, 3, 94, ..., 66, 40, 17],\n", + " [87, 55, 87, ..., 41, 26, 55],\n", " ...,\n", "...\n", " ...,\n", - " [27, 8, 29, ..., 96, 80, 60],\n", - " [ 1, 98, 63, ..., 1, 73, 98],\n", - " [77, 22, 88, ..., 6, 60, 7]],\n", + " [67, 50, 96, ..., 28, 89, 55],\n", + " [83, 98, 3, ..., 91, 25, 14],\n", + " [81, 7, 45, ..., 18, 78, 9]],\n", "\n", - " [[22, 40, 31, ..., 84, 62, 61],\n", - " [27, 83, 95, ..., 81, 15, 54],\n", - " [43, 45, 67, ..., 52, 2, 28],\n", + " [[29, 94, 53, ..., 3, 52, 47],\n", + " [93, 43, 11, ..., 94, 71, 37],\n", + " [80, 34, 29, ..., 12, 67, 48],\n", " ...,\n", - " [91, 75, 49, ..., 1, 76, 65],\n", - " [ 4, 78, 43, ..., 82, 20, 95],\n", - " [16, 65, 36, ..., 27, 92, 43]],\n", + " [31, 96, 58, ..., 71, 85, 90],\n", + " [72, 95, 55, ..., 83, 17, 97],\n", + " [72, 77, 15, ..., 42, 62, 9]],\n", "\n", - " [[11, 54, 89, ..., 64, 56, 78],\n", - " [ 1, 29, 49, ..., 46, 22, 95],\n", - " [99, 4, 33, ..., 89, 48, 60],\n", + " [[26, 76, 23, ..., 21, 92, 42],\n", + " [60, 80, 37, ..., 49, 91, 88],\n", + " [55, 64, 94, ..., 99, 91, 43],\n", " ...,\n", - " [18, 62, 84, ..., 12, 87, 52],\n", - " [42, 79, 90, ..., 5, 24, 57],\n", - " [51, 35, 71, ..., 55, 29, 11]]]]])\n", + " [80, 85, 15, ..., 23, 42, 52],\n", + " [66, 91, 73, ..., 10, 44, 22],\n", + " [75, 27, 89, ..., 37, 57, 20]]]]])\n", "Coordinates:\n", " * mode (mode) 
    <xarray.DataArray (mode: 3, origin: 2, destination: 2)>\n",
-       "array([[[14, 88],\n",
-       "        [ 9, 22]],\n",
-       "\n",
-       "       [[99, 46],\n",
-       "        [25, 62]],\n",
-       "\n",
-       "       [[20, 54],\n",
-       "        [44, 17]]])\n",
-       "Coordinates:\n",
-       "  * mode         (mode) <U4 'car' 'bike' 'foot'\n",
-       "    day          datetime64[ns] 2015-01-01\n",
-       "    time         int64 0\n",
-       "    origin       (origin) object MULTIPOLYGON (((-79.91995239257812 34.807918...\n",
-       "    destination  (destination) object MULTIPOLYGON (((-80.72651672363281 35.5...
    " - ], - "text/plain": [ - "\n", - "array([[[14, 88],\n", - " [ 9, 22]],\n", - "\n", - " [[99, 46],\n", - " [25, 62]],\n", - "\n", - " [[20, 54],\n", - " [44, 17]]])\n", - "Coordinates:\n", - " * mode (mode) \n", - "\n", - "\n", - "\n", - "\n", - "\n", - "\n", - "\n", - "\n", - "\n", - "\n", - "\n", - "\n", - "\n", - "\n", - "
    <xarray.DataArray (mode: 3, day: 100, time: 24, origin: 2, destination: 100)>\n",
+       "array([[[[[47, 51, 75, ..., 15, 54, 82],\n",
+       "          [ 7, 61, 46, ...,  7, 65, 17]],\n",
        "\n",
-       ".xr-header {\n",
-       "  padding-top: 6px;\n",
-       "  padding-bottom: 6px;\n",
-       "  margin-bottom: 4px;\n",
-       "  border-bottom: solid 1px var(--xr-border-color);\n",
-       "}\n",
+       "         [[67,  5, 73, ..., 95, 34, 44],\n",
+       "          [62, 80,  2, ..., 56, 38, 33]],\n",
        "\n",
-       ".xr-header > div,\n",
-       ".xr-header > ul {\n",
-       "  display: inline;\n",
-       "  margin-top: 0;\n",
-       "  margin-bottom: 0;\n",
-       "}\n",
+       "         [[22, 57, 42, ..., 24, 60,  1],\n",
+       "          [47,  3, 82, ..., 52, 56, 23]],\n",
        "\n",
-       ".xr-obj-type,\n",
-       ".xr-array-name {\n",
-       "  margin-left: 2px;\n",
-       "  margin-right: 10px;\n",
-       "}\n",
+       "         ...,\n",
        "\n",
-       ".xr-obj-type {\n",
-       "  color: var(--xr-font-color2);\n",
-       "}\n",
+       "         [[ 6, 83, 84, ..., 50, 77, 13],\n",
+       "          [66, 67, 24, ..., 60, 14, 91]],\n",
        "\n",
-       ".xr-sections {\n",
-       "  padding-left: 0 !important;\n",
-       "  display: grid;\n",
-       "  grid-template-columns: 150px auto auto 1fr 20px 20px;\n",
-       "}\n",
+       "         [[55, 85, 99, ..., 28, 13, 39],\n",
+       "          [76, 92, 33, ..., 18, 31, 75]],\n",
        "\n",
-       ".xr-section-item {\n",
-       "  display: contents;\n",
-       "}\n",
+       "         [[86, 47, 79, ..., 70, 99, 93],\n",
+       "          [86, 87, 78, ..., 66, 26, 80]]],\n",
        "\n",
-       ".xr-section-item input {\n",
-       "  display: none;\n",
-       "}\n",
+       "...\n",
        "\n",
-       ".xr-section-item input + label {\n",
-       "  color: var(--xr-disabled-color);\n",
-       "}\n",
+       "        [[[52, 34, 75, ...,  7, 84, 13],\n",
+       "          [34, 11, 38, ..., 70, 89, 56]],\n",
        "\n",
-       ".xr-section-item input:enabled + label {\n",
-       "  cursor: pointer;\n",
-       "  color: var(--xr-font-color2);\n",
-       "}\n",
+       "         [[ 9, 52, 95, ..., 82, 89, 44],\n",
+       "          [84, 94, 76, ..., 54, 34, 49]],\n",
        "\n",
-       ".xr-section-item input:enabled + label:hover {\n",
-       "  color: var(--xr-font-color0);\n",
-       "}\n",
+       "         [[ 9, 64, 85, ..., 28, 44, 22],\n",
+       "          [15, 18, 14, ..., 47, 81, 18]],\n",
        "\n",
-       ".xr-section-summary {\n",
-       "  grid-column: 1;\n",
-       "  color: var(--xr-font-color2);\n",
-       "  font-weight: 500;\n",
-       "}\n",
+       "         ...,\n",
        "\n",
-       ".xr-section-summary > span {\n",
-       "  display: inline-block;\n",
-       "  padding-left: 0.5em;\n",
-       "}\n",
+       "         [[ 9, 82, 60, ..., 64, 56, 98],\n",
+       "          [30, 70, 34, ..., 17, 34, 86]],\n",
        "\n",
-       ".xr-section-summary-in:disabled + label {\n",
-       "  color: var(--xr-font-color2);\n",
-       "}\n",
+       "         [[29, 94, 53, ...,  3, 52, 47],\n",
+       "          [56, 99, 33, ..., 81, 88, 85]],\n",
        "\n",
-       ".xr-section-summary-in + label:before {\n",
-       "  display: inline-block;\n",
-       "  content: 'â–º';\n",
-       "  font-size: 11px;\n",
-       "  width: 15px;\n",
-       "  text-align: center;\n",
-       "}\n",
+       "         [[26, 76, 23, ..., 21, 92, 42],\n",
+       "          [44, 47, 95, ..., 84,  5, 86]]]]])\n",
+       "Coordinates:\n",
+       "  * mode         (mode) <U4 'car' 'bike' 'foot'\n",
+       "  * day          (day) datetime64[ns] 2015-01-01 2015-01-02 ... 2015-04-10\n",
+       "  * time         (time) int64 0 1 2 3 4 5 6 7 8 9 ... 15 16 17 18 19 20 21 22 23\n",
+       "  * origin       (origin) object MULTIPOLYGON (((-81.4727554321289 36.2343559...\n",
+       "  * destination  (destination) object MULTIPOLYGON (((-81.4727554321289 36.23...\n",
+       "Indexes:\n",
+       "    origin       GeometryIndex (crs=EPSG:4267)\n",
+       "    destination  GeometryIndex (crs=EPSG:4267)
    " + ], + "text/plain": [ + "\n", + "array([[[[[47, 51, 75, ..., 15, 54, 82],\n", + " [ 7, 61, 46, ..., 7, 65, 17]],\n", + "\n", + " [[67, 5, 73, ..., 95, 34, 44],\n", + " [62, 80, 2, ..., 56, 38, 33]],\n", + "\n", + " [[22, 57, 42, ..., 24, 60, 1],\n", + " [47, 3, 82, ..., 52, 56, 23]],\n", + "\n", + " ...,\n", + "\n", + " [[ 6, 83, 84, ..., 50, 77, 13],\n", + " [66, 67, 24, ..., 60, 14, 91]],\n", + "\n", + " [[55, 85, 99, ..., 28, 13, 39],\n", + " [76, 92, 33, ..., 18, 31, 75]],\n", + "\n", + " [[86, 47, 79, ..., 70, 99, 93],\n", + " [86, 87, 78, ..., 66, 26, 80]]],\n", + "\n", + "...\n", + "\n", + " [[[52, 34, 75, ..., 7, 84, 13],\n", + " [34, 11, 38, ..., 70, 89, 56]],\n", + "\n", + " [[ 9, 52, 95, ..., 82, 89, 44],\n", + " [84, 94, 76, ..., 54, 34, 49]],\n", + "\n", + " [[ 9, 64, 85, ..., 28, 44, 22],\n", + " [15, 18, 14, ..., 47, 81, 18]],\n", + "\n", + " ...,\n", + "\n", + " [[ 9, 82, 60, ..., 64, 56, 98],\n", + " [30, 70, 34, ..., 17, 34, 86]],\n", + "\n", + " [[29, 94, 53, ..., 3, 52, 47],\n", + " [56, 99, 33, ..., 81, 88, 85]],\n", + "\n", + " [[26, 76, 23, ..., 21, 92, 42],\n", + " [44, 47, 95, ..., 84, 5, 86]]]]])\n", + "Coordinates:\n", + " * mode (mode) \n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "
    <xarray.DataArray (mode: 3, origin: 2, destination: 2)>\n",
+       "array([[[47, 23],\n",
+       "        [93, 84]],\n",
+       "\n",
+       "       [[35, 46],\n",
+       "        [18, 68]],\n",
+       "\n",
+       "       [[62, 69],\n",
+       "        [49, 32]]])\n",
+       "Coordinates:\n",
+       "  * mode         (mode) <U4 'car' 'bike' 'foot'\n",
+       "    day          datetime64[ns] 2015-01-01\n",
+       "    time         int64 0\n",
+       "  * origin       (origin) object MULTIPOLYGON (((-79.91995239257812 34.807918...\n",
+       "  * destination  (destination) object MULTIPOLYGON (((-80.72651672363281 35.5...\n",
+       "Indexes:\n",
+       "    origin       GeometryIndex (crs=EPSG:4267)\n",
+       "    destination  GeometryIndex (crs=EPSG:4267)
    " + ], + "text/plain": [ + "\n", + "array([[[47, 23],\n", + " [93, 84]],\n", + "\n", + " [[35, 46],\n", + " [18, 68]],\n", + "\n", + " [[62, 69],\n", + " [49, 32]]])\n", + "Coordinates:\n", + " * mode (mode) \n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "
    <xarray.DataArray (mode: 3, day: 100, time: 24, origin: 41, destination: 100)>\n",
+       "array([[[[[77, 84, 80, ..., 31, 71, 18],\n",
+       "          [30, 62, 58, ..., 26,  9, 61],\n",
+       "          [96, 86, 84, ..., 31, 69, 57],\n",
+       "          ...,\n",
+       "          [98, 86, 56, ..., 94, 22, 75],\n",
+       "          [85, 42, 89, ..., 23, 69, 10],\n",
+       "          [98, 90, 56, ..., 35, 26, 39]],\n",
+       "\n",
+       "         [[ 7, 45, 63, ..., 82, 61, 99],\n",
+       "          [88, 96, 84, ...,  7, 87, 18],\n",
+       "          [67, 72, 77, ..., 27, 54,  5],\n",
+       "          ...,\n",
+       "          [86, 52, 62, ..., 50, 51, 28],\n",
+       "          [50,  8, 17, ..., 63, 93, 77],\n",
+       "          [95, 50, 26, ..., 56, 85, 46]],\n",
+       "\n",
+       "         [[27, 67, 87, ..., 88, 63, 37],\n",
+       "          [27, 46, 96, ..., 26, 29, 10],\n",
+       "          [60, 86, 15, ...,  9, 51, 29],\n",
+       "          ...,\n",
+       "...\n",
+       "          ...,\n",
+       "          [86,  7, 13, ..., 66, 74, 38],\n",
+       "          [59, 33, 37, ..., 85, 92, 18],\n",
+       "          [90, 87, 39, ..., 15, 99, 71]],\n",
+       "\n",
+       "         [[28, 62,  2, ..., 41,  2, 84],\n",
+       "          [70, 99,  7, ..., 55, 28, 78],\n",
+       "          [27, 78, 16, ..., 73, 89, 24],\n",
+       "          ...,\n",
+       "          [15, 80, 70, ..., 20, 12, 10],\n",
+       "          [77,  7, 58, ..., 90, 12, 78],\n",
+       "          [67, 89, 62, ..., 39,  8, 68]],\n",
+       "\n",
+       "         [[87, 17, 38, ..., 41, 73, 36],\n",
+       "          [58, 93, 31, ..., 13, 78, 28],\n",
+       "          [30, 47, 25, ..., 92, 47, 74],\n",
+       "          ...,\n",
+       "          [12,  8, 60, ..., 86, 25, 44],\n",
+       "          [88, 56, 41, ..., 39, 19, 93],\n",
+       "          [29, 88, 17, ..., 54, 93, 57]]]]])\n",
+       "Coordinates:\n",
+       "  * mode         (mode) <U4 'car' 'bike' 'foot'\n",
+       "  * day          (day) datetime64[ns] 2015-01-01 2015-01-02 ... 2015-04-10\n",
+       "  * time         (time) int64 0 1 2 3 4 5 6 7 8 9 ... 15 16 17 18 19 20 21 22 23\n",
+       "  * origin       (origin) object MULTIPOLYGON (((-78.11376953125 34.720985412...\n",
+       "  * destination  (destination) object MULTIPOLYGON (((-81.4727554321289 36.23...\n",
+       "Indexes:\n",
+       "    origin       GeometryIndex (crs=EPSG:4267)\n",
+       "    destination  GeometryIndex (crs=EPSG:4267)
    " + ], + "text/plain": [ + "\n", + "array([[[[[77, 84, 80, ..., 31, 71, 18],\n", + " [30, 62, 58, ..., 26, 9, 61],\n", + " [96, 86, 84, ..., 31, 69, 57],\n", + " ...,\n", + " [98, 86, 56, ..., 94, 22, 75],\n", + " [85, 42, 89, ..., 23, 69, 10],\n", + " [98, 90, 56, ..., 35, 26, 39]],\n", + "\n", + " [[ 7, 45, 63, ..., 82, 61, 99],\n", + " [88, 96, 84, ..., 7, 87, 18],\n", + " [67, 72, 77, ..., 27, 54, 5],\n", + " ...,\n", + " [86, 52, 62, ..., 50, 51, 28],\n", + " [50, 8, 17, ..., 63, 93, 77],\n", + " [95, 50, 26, ..., 56, 85, 46]],\n", + "\n", + " [[27, 67, 87, ..., 88, 63, 37],\n", + " [27, 46, 96, ..., 26, 29, 10],\n", + " [60, 86, 15, ..., 9, 51, 29],\n", + " ...,\n", + "...\n", + " ...,\n", + " [86, 7, 13, ..., 66, 74, 38],\n", + " [59, 33, 37, ..., 85, 92, 18],\n", + " [90, 87, 39, ..., 15, 99, 71]],\n", + "\n", + " [[28, 62, 2, ..., 41, 2, 84],\n", + " [70, 99, 7, ..., 55, 28, 78],\n", + " [27, 78, 16, ..., 73, 89, 24],\n", + " ...,\n", + " [15, 80, 70, ..., 20, 12, 10],\n", + " [77, 7, 58, ..., 90, 12, 78],\n", + " [67, 89, 62, ..., 39, 8, 68]],\n", + "\n", + " [[87, 17, 38, ..., 41, 73, 36],\n", + " [58, 93, 31, ..., 13, 78, 28],\n", + " [30, 47, 25, ..., 92, 47, 74],\n", + " ...,\n", + " [12, 8, 60, ..., 86, 25, 44],\n", + " [88, 56, 41, ..., 39, 19, 93],\n", + " [29, 88, 17, ..., 54, 93, 57]]]]])\n", + "Coordinates:\n", + " * mode (mode) 1\u001b[0m \u001b[43marr\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43msel\u001b[49m\u001b[43m(\u001b[49m\u001b[43morigin\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[43m[\u001b[49m\u001b[43mshapely\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mPoint\u001b[49m\u001b[43m(\u001b[49m\u001b[38;5;241;43m-\u001b[39;49m\u001b[38;5;241;43m80\u001b[39;49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[38;5;241;43m35\u001b[39;49m\u001b[43m)\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mshapely\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mPoint\u001b[49m\u001b[43m(\u001b[49m\u001b[38;5;241;43m-\u001b[39;49m\u001b[38;5;241;43m76\u001b[39;49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[38;5;241;43m35.6\u001b[39;49m\u001b[43m)\u001b[49m\u001b[43m]\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mmethod\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[38;5;124;43m\"\u001b[39;49m\u001b[38;5;124;43mintersects\u001b[39;49m\u001b[38;5;124;43m\"\u001b[39;49m\u001b[43m)\u001b[49m\n", + "File \u001b[0;32m~/miniconda3/envs/xvec_dev/lib/python3.11/site-packages/xarray/core/dataarray.py:1526\u001b[0m, in \u001b[0;36mDataArray.sel\u001b[0;34m(self, indexers, method, tolerance, drop, **indexers_kwargs)\u001b[0m\n\u001b[1;32m 1416\u001b[0m \u001b[38;5;28;01mdef\u001b[39;00m \u001b[38;5;21msel\u001b[39m(\n\u001b[1;32m 1417\u001b[0m \u001b[38;5;28mself\u001b[39m: T_DataArray,\n\u001b[1;32m 1418\u001b[0m indexers: Mapping[Any, Any] \u001b[38;5;241m=\u001b[39m \u001b[38;5;28;01mNone\u001b[39;00m,\n\u001b[0;32m (...)\u001b[0m\n\u001b[1;32m 1422\u001b[0m \u001b[38;5;241m*\u001b[39m\u001b[38;5;241m*\u001b[39mindexers_kwargs: Any,\n\u001b[1;32m 1423\u001b[0m ) \u001b[38;5;241m-\u001b[39m\u001b[38;5;241m>\u001b[39m T_DataArray:\n\u001b[1;32m 1424\u001b[0m \u001b[38;5;124;03m\"\"\"Return a new DataArray whose data is given by selecting index\u001b[39;00m\n\u001b[1;32m 1425\u001b[0m \u001b[38;5;124;03m labels along the specified dimension(s).\u001b[39;00m\n\u001b[1;32m 1426\u001b[0m \n\u001b[0;32m (...)\u001b[0m\n\u001b[1;32m 1524\u001b[0m \u001b[38;5;124;03m Dimensions without coordinates: points\u001b[39;00m\n\u001b[1;32m 1525\u001b[0m \u001b[38;5;124;03m \"\"\"\u001b[39;00m\n\u001b[0;32m-> 1526\u001b[0m ds \u001b[38;5;241m=\u001b[39m \u001b[38;5;28;43mself\u001b[39;49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43m_to_temp_dataset\u001b[49m\u001b[43m(\u001b[49m\u001b[43m)\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43msel\u001b[49m\u001b[43m(\u001b[49m\n\u001b[1;32m 1527\u001b[0m \u001b[43m \u001b[49m\u001b[43mindexers\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[43mindexers\u001b[49m\u001b[43m,\u001b[49m\n\u001b[1;32m 1528\u001b[0m \u001b[43m \u001b[49m\u001b[43mdrop\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[43mdrop\u001b[49m\u001b[43m,\u001b[49m\n\u001b[1;32m 1529\u001b[0m \u001b[43m \u001b[49m\u001b[43mmethod\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[43mmethod\u001b[49m\u001b[43m,\u001b[49m\n\u001b[1;32m 1530\u001b[0m \u001b[43m \u001b[49m\u001b[43mtolerance\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[43mtolerance\u001b[49m\u001b[43m,\u001b[49m\n\u001b[1;32m 1531\u001b[0m \u001b[43m \u001b[49m\u001b[38;5;241;43m*\u001b[39;49m\u001b[38;5;241;43m*\u001b[39;49m\u001b[43mindexers_kwargs\u001b[49m\u001b[43m,\u001b[49m\n\u001b[1;32m 1532\u001b[0m \u001b[43m \u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m 1533\u001b[0m \u001b[38;5;28;01mreturn\u001b[39;00m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39m_from_temp_dataset(ds)\n", + "File \u001b[0;32m~/miniconda3/envs/xvec_dev/lib/python3.11/site-packages/xarray/core/dataset.py:2554\u001b[0m, in \u001b[0;36mDataset.sel\u001b[0;34m(self, indexers, method, tolerance, drop, **indexers_kwargs)\u001b[0m\n\u001b[1;32m 2493\u001b[0m \u001b[38;5;124;03m\"\"\"Returns a new dataset with each array indexed by tick labels\u001b[39;00m\n\u001b[1;32m 2494\u001b[0m \u001b[38;5;124;03malong the specified dimension(s).\u001b[39;00m\n\u001b[1;32m 2495\u001b[0m \n\u001b[0;32m (...)\u001b[0m\n\u001b[1;32m 2551\u001b[0m \u001b[38;5;124;03mDataArray.sel\u001b[39;00m\n\u001b[1;32m 2552\u001b[0m \u001b[38;5;124;03m\"\"\"\u001b[39;00m\n\u001b[1;32m 2553\u001b[0m indexers \u001b[38;5;241m=\u001b[39m either_dict_or_kwargs(indexers, indexers_kwargs, \u001b[38;5;124m\"\u001b[39m\u001b[38;5;124msel\u001b[39m\u001b[38;5;124m\"\u001b[39m)\n\u001b[0;32m-> 2554\u001b[0m query_results \u001b[38;5;241m=\u001b[39m \u001b[43mmap_index_queries\u001b[49m\u001b[43m(\u001b[49m\n\u001b[1;32m 2555\u001b[0m \u001b[43m \u001b[49m\u001b[38;5;28;43mself\u001b[39;49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mindexers\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[43mindexers\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mmethod\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[43mmethod\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mtolerance\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[43mtolerance\u001b[49m\n\u001b[1;32m 2556\u001b[0m \u001b[43m\u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m 2558\u001b[0m \u001b[38;5;28;01mif\u001b[39;00m drop:\n\u001b[1;32m 2559\u001b[0m no_scalar_variables \u001b[38;5;241m=\u001b[39m {}\n", + "File \u001b[0;32m~/miniconda3/envs/xvec_dev/lib/python3.11/site-packages/xarray/core/indexing.py:183\u001b[0m, in \u001b[0;36mmap_index_queries\u001b[0;34m(obj, indexers, method, tolerance, **indexers_kwargs)\u001b[0m\n\u001b[1;32m 181\u001b[0m results\u001b[38;5;241m.\u001b[39mappend(IndexSelResult(labels))\n\u001b[1;32m 182\u001b[0m \u001b[38;5;28;01melse\u001b[39;00m:\n\u001b[0;32m--> 183\u001b[0m results\u001b[38;5;241m.\u001b[39mappend(\u001b[43mindex\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43msel\u001b[49m\u001b[43m(\u001b[49m\u001b[43mlabels\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[38;5;241;43m*\u001b[39;49m\u001b[38;5;241;43m*\u001b[39;49m\u001b[43moptions\u001b[49m\u001b[43m)\u001b[49m) \u001b[38;5;66;03m# type: ignore[call-arg]\u001b[39;00m\n\u001b[1;32m 185\u001b[0m merged \u001b[38;5;241m=\u001b[39m merge_sel_results(results)\n\u001b[1;32m 187\u001b[0m \u001b[38;5;66;03m# drop dimension coordinates found in dimension indexers\u001b[39;00m\n\u001b[1;32m 188\u001b[0m \u001b[38;5;66;03m# (also drop multi-index if any)\u001b[39;00m\n\u001b[1;32m 189\u001b[0m \u001b[38;5;66;03m# (.sel() already ensures alignment)\u001b[39;00m\n", + "File \u001b[0;32m~/Git/github/benbovy/xvec/xvec/index.py:229\u001b[0m, in \u001b[0;36mGeometryIndex.sel\u001b[0;34m(self, labels, method, tolerance)\u001b[0m\n\u001b[1;32m 222\u001b[0m \u001b[38;5;28;01mreturn\u001b[39;00m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39m_index\u001b[38;5;241m.\u001b[39msel(labels)\n\u001b[1;32m 223\u001b[0m \u001b[38;5;28;01melse\u001b[39;00m:\n\u001b[1;32m 224\u001b[0m \u001b[38;5;66;03m# We reuse here `method` and `tolerance` options of\u001b[39;00m\n\u001b[1;32m 225\u001b[0m \u001b[38;5;66;03m# `xarray.indexes.PandasIndex` as `predicate` and `distance`\u001b[39;00m\n\u001b[1;32m 226\u001b[0m \u001b[38;5;66;03m# options when `labels` is a single geometry.\u001b[39;00m\n\u001b[1;32m 227\u001b[0m \u001b[38;5;66;03m# Xarray currently doesn't support custom options\u001b[39;00m\n\u001b[1;32m 228\u001b[0m \u001b[38;5;66;03m# (see https://github.com/pydata/xarray/issues/7099)\u001b[39;00m\n\u001b[0;32m--> 229\u001b[0m \u001b[38;5;28;01mreturn\u001b[39;00m \u001b[38;5;28;43mself\u001b[39;49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43m_sel_sindex\u001b[49m\u001b[43m(\u001b[49m\u001b[43mlabels\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mmethod\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mtolerance\u001b[49m\u001b[43m)\u001b[49m\n", + "File \u001b[0;32m~/Git/github/benbovy/xvec/xvec/index.py:182\u001b[0m, in \u001b[0;36mGeometryIndex._sel_sindex\u001b[0;34m(self, labels, method, tolerance)\u001b[0m\n\u001b[1;32m 180\u001b[0m \u001b[38;5;28;01mif\u001b[39;00m method \u001b[38;5;241m!=\u001b[39m \u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mnearest\u001b[39m\u001b[38;5;124m\"\u001b[39m:\n\u001b[1;32m 181\u001b[0m \u001b[38;5;28;01mif\u001b[39;00m \u001b[38;5;129;01mnot\u001b[39;00m \u001b[38;5;28misinstance\u001b[39m(label, shapely\u001b[38;5;241m.\u001b[39mGeometry):\n\u001b[0;32m--> 182\u001b[0m \u001b[38;5;28;01mraise\u001b[39;00m \u001b[38;5;167;01mValueError\u001b[39;00m(\n\u001b[1;32m 183\u001b[0m \u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mselection with another method than nearest only supports \u001b[39m\u001b[38;5;124m\"\u001b[39m\n\u001b[1;32m 184\u001b[0m \u001b[38;5;124m\"\u001b[39m\u001b[38;5;124ma single geometry as input label.\u001b[39m\u001b[38;5;124m\"\u001b[39m\n\u001b[1;32m 185\u001b[0m )\n\u001b[1;32m 187\u001b[0m \u001b[38;5;28;01mif\u001b[39;00m \u001b[38;5;28misinstance\u001b[39m(label, DataArray):\n\u001b[1;32m 188\u001b[0m label_array \u001b[38;5;241m=\u001b[39m label\u001b[38;5;241m.\u001b[39m_variable\u001b[38;5;241m.\u001b[39m_data\n", + "\u001b[0;31mValueError\u001b[0m: selection with another method than nearest only supports a single geometry as input label." + ] + } + ], + "source": [ + "arr.sel(origin=[shapely.Point(-80, 35), shapely.Point(-76, 35.6)], method=\"intersects\")" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Other operations like `isel`, `concat` and alignment work just like with any default (pandas) index:" + ] + }, + { + "cell_type": "code", + "execution_count": 14, + "metadata": {}, + "outputs": [ + { + "data": { + "text/html": [ + "
    \n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "
    <xarray.DataArray (mode: 3, day: 100, time: 24, origin: 41, destination: 100)>\n",
-       "array([[[[[41, 61, 90, ..., 86, 69, 39],\n",
-       "          [38, 45, 93, ..., 95, 53, 80],\n",
-       "          [53, 31,  5, ..., 74, 48, 96],\n",
-       "          ...,\n",
-       "          [98, 48, 42, ..., 41, 81, 69],\n",
-       "          [87, 18, 10, ..., 97, 15, 33],\n",
-       "          [38, 68, 19, ..., 71,  9, 70]],\n",
+       "
    <xarray.DataArray (mode: 3, day: 100, time: 24, origin: 2, destination: 2)>\n",
+       "array([[[[[108, 164],\n",
+       "          [  6, 144]],\n",
+       "\n",
+       "         [[ 68,  88],\n",
+       "          [ 10,  82]],\n",
+       "\n",
+       "         [[120,   2],\n",
+       "          [ 80,  34]],\n",
+       "\n",
+       "         ...,\n",
        "\n",
-       "         [[30,  3, 51, ..., 60, 70,  6],\n",
-       "          [40,  3, 36, ..., 32, 22, 32],\n",
-       "          [74,  8, 25, ..., 91, 65, 71],\n",
-       "          ...,\n",
-       "          [72, 63, 10, ..., 54, 84, 63],\n",
-       "          [12, 77, 74, ..., 10, 30, 32],\n",
-       "          [27, 73, 36, ..., 63, 16, 57]],\n",
+       "         [[154,  26],\n",
+       "          [104,  94]],\n",
+       "\n",
+       "         [[ 26,  78],\n",
+       "          [  8, 152]],\n",
+       "\n",
+       "         [[198, 186],\n",
+       "          [ 40, 114]]],\n",
        "\n",
-       "         [[68, 46, 26, ..., 71, 71, 33],\n",
-       "          [21, 92, 65, ..., 94, 44, 88],\n",
-       "          [62, 13, 97, ..., 10, 42,  2],\n",
-       "          ...,\n",
        "...\n",
-       "          ...,\n",
-       "          [53, 77, 85, ...,  9, 56, 44],\n",
-       "          [17, 72, 69, ..., 34, 63, 43],\n",
-       "          [17,  3,  1, ..., 40, 55, 74]],\n",
        "\n",
-       "         [[ 1, 27, 55, ..., 74, 12, 45],\n",
-       "          [88, 31, 98, ..., 24,  4, 54],\n",
-       "          [25, 94, 82, ..., 44,  6, 42],\n",
-       "          ...,\n",
-       "          [97, 21, 97, ..., 25, 99, 15],\n",
-       "          [49, 44,  3, ..., 83, 13, 14],\n",
-       "          [57, 69, 97, ..., 27, 32, 77]],\n",
+       "        [[[168,  26],\n",
+       "          [ 78,  10]],\n",
        "\n",
-       "         [[76, 34, 34, ..., 62, 81, 11],\n",
-       "          [59, 57, 36, ..., 86, 47, 48],\n",
-       "          [71, 58, 18, ..., 26, 33, 70],\n",
-       "          ...,\n",
-       "          [55,  9, 44, ..., 63, 60,  5],\n",
-       "          [44, 28,  3, ..., 61, 81,  9],\n",
-       "          [51, 19,  7, ..., 39, 92, 50]]]]])\n",
+       "         [[178,  88],\n",
+       "          [114,  20]],\n",
+       "\n",
+       "         [[ 88,  44],\n",
+       "          [192, 184]],\n",
+       "\n",
+       "         ...,\n",
+       "\n",
+       "         [[112, 196],\n",
+       "          [168, 120]],\n",
+       "\n",
+       "         [[104,  94],\n",
+       "          [142,  74]],\n",
+       "\n",
+       "         [[184,  84],\n",
+       "          [182, 176]]]]])\n",
        "Coordinates:\n",
        "  * mode         (mode) <U4 'car' 'bike' 'foot'\n",
        "  * day          (day) datetime64[ns] 2015-01-01 2015-01-02 ... 2015-04-10\n",
        "  * time         (time) int64 0 1 2 3 4 5 6 7 8 9 ... 15 16 17 18 19 20 21 22 23\n",
-       "    origin       (origin) object MULTIPOLYGON (((-78.11376953125 34.720985412...\n",
-       "  * destination  (destination) object MULTIPOLYGON (((-81.4727554321289 36.23...\n",
+       "  * origin       (origin) object MULTIPOLYGON (((-81.4727554321289 36.2343559...\n",
+       "  * destination  (destination) object MULTIPOLYGON (((-77.96073150634766 34.1...\n",
        "Indexes:\n",
-       "    destination  ShapelySTRTreeIndex
  • origin
    GeometryIndex (crs=EPSG:4267)
    <xvec.index.GeometryIndex object at 0x16ab6ca50>
  • destination
    GeometryIndex (crs=EPSG:4267)
    <xvec.index.GeometryIndex object at 0x16b513210>
    • " ], "text/plain": [ - "\n", - "array([[[[[41, 61, 90, ..., 86, 69, 39],\n", - " [38, 45, 93, ..., 95, 53, 80],\n", - " [53, 31, 5, ..., 74, 48, 96],\n", - " ...,\n", - " [98, 48, 42, ..., 41, 81, 69],\n", - " [87, 18, 10, ..., 97, 15, 33],\n", - " [38, 68, 19, ..., 71, 9, 70]],\n", + "\n", + "array([[[[[108, 164],\n", + " [ 6, 144]],\n", "\n", - " [[30, 3, 51, ..., 60, 70, 6],\n", - " [40, 3, 36, ..., 32, 22, 32],\n", - " [74, 8, 25, ..., 91, 65, 71],\n", - " ...,\n", - " [72, 63, 10, ..., 54, 84, 63],\n", - " [12, 77, 74, ..., 10, 30, 32],\n", - " [27, 73, 36, ..., 63, 16, 57]],\n", + " [[ 68, 88],\n", + " [ 10, 82]],\n", + "\n", + " [[120, 2],\n", + " [ 80, 34]],\n", + "\n", + " ...,\n", + "\n", + " [[154, 26],\n", + " [104, 94]],\n", + "\n", + " [[ 26, 78],\n", + " [ 8, 152]],\n", + "\n", + " [[198, 186],\n", + " [ 40, 114]]],\n", "\n", - " [[68, 46, 26, ..., 71, 71, 33],\n", - " [21, 92, 65, ..., 94, 44, 88],\n", - " [62, 13, 97, ..., 10, 42, 2],\n", - " ...,\n", "...\n", - " ...,\n", - " [53, 77, 85, ..., 9, 56, 44],\n", - " [17, 72, 69, ..., 34, 63, 43],\n", - " [17, 3, 1, ..., 40, 55, 74]],\n", "\n", - " [[ 1, 27, 55, ..., 74, 12, 45],\n", - " [88, 31, 98, ..., 24, 4, 54],\n", - " [25, 94, 82, ..., 44, 6, 42],\n", - " ...,\n", - " [97, 21, 97, ..., 25, 99, 15],\n", - " [49, 44, 3, ..., 83, 13, 14],\n", - " [57, 69, 97, ..., 27, 32, 77]],\n", + " [[[168, 26],\n", + " [ 78, 10]],\n", "\n", - " [[76, 34, 34, ..., 62, 81, 11],\n", - " [59, 57, 36, ..., 86, 47, 48],\n", - " [71, 58, 18, ..., 26, 33, 70],\n", - " ...,\n", - " [55, 9, 44, ..., 63, 60, 5],\n", - " [44, 28, 3, ..., 61, 81, 9],\n", - " [51, 19, 7, ..., 39, 92, 50]]]]])\n", + " [[178, 88],\n", + " [114, 20]],\n", + "\n", + " [[ 88, 44],\n", + " [192, 184]],\n", + "\n", + " ...,\n", + "\n", + " [[112, 196],\n", + " [168, 120]],\n", + "\n", + " [[104, 94],\n", + " [142, 74]],\n", + "\n", + " [[184, 84],\n", + " [182, 176]]]]])\n", "Coordinates:\n", " * mode (mode) 
    <xarray.DataArray (mode: 3, day: 100, time: 24, points: 2, destination: 100)>\n",
-       "array([[[[[31, 69, 46, ..., 85, 87, 41],\n",
-       "          [63, 34, 15, ..., 62, 36, 23]],\n",
+       "
    <xarray.DataArray (mode: 3, day: 100, time: 24, origin: 4, destination: 100)>\n",
+       "array([[[[[47, 51, 75, ..., 15, 54, 82],\n",
+       "          [ 7, 68, 76, ..., 38,  3, 72],\n",
+       "          [19, 47,  3, ..., 75, 42, 69],\n",
+       "          [48, 94, 34, ..., 47, 54, 86]],\n",
        "\n",
-       "         [[79, 46, 26, ..., 75, 89, 12],\n",
-       "          [99, 67, 18, ..., 93, 63,  9]],\n",
+       "         [[67,  5, 73, ..., 95, 34, 44],\n",
+       "          [80, 22, 89, ..., 93,  5, 41],\n",
+       "          [78, 49,  9, ..., 45,  3,  9],\n",
+       "          [87, 73, 52, ..., 32, 58, 48]],\n",
        "\n",
-       "         [[ 7, 22, 63, ..., 78, 89, 64],\n",
-       "          [60, 33, 13, ..., 39, 54, 33]],\n",
+       "         [[22, 57, 42, ..., 24, 60,  1],\n",
+       "          [ 4,  3, 94, ..., 66, 40, 17],\n",
+       "          [81, 36, 30, ...,  5, 51,  6],\n",
+       "          [54, 70, 25, ..., 30, 18, 37]],\n",
        "\n",
        "         ...,\n",
        "\n",
-       "         [[22, 59, 96, ..., 54, 84, 18],\n",
-       "          [39, 52, 31, ..., 40, 45, 43]],\n",
-       "\n",
-       "         [[15, 37, 54, ..., 30,  2, 57],\n",
-       "          [29, 90,  1, ..., 60, 21, 29]],\n",
-       "\n",
-       "         [[97, 47, 98, ..., 29, 76, 77],\n",
-       "          [43, 42, 23, ..., 68, 56, 42]]],\n",
-       "\n",
+       "         [[ 6, 83, 84, ..., 50, 77, 13],\n",
+       "          [44, 84, 52, ..., 38, 52, 47],\n",
+       "          [65, 26, 50, ..., 12, 21, 25],\n",
        "...\n",
-       "\n",
-       "        [[[53, 70, 33, ...,  6, 78, 35],\n",
-       "          [62, 58, 58, ..., 35, 88, 12]],\n",
-       "\n",
-       "         [[27, 77,  3, ..., 62, 44, 34],\n",
-       "          [95, 88, 66, ..., 67, 63, 83]],\n",
-       "\n",
-       "         [[43, 78, 66, ..., 22, 16, 56],\n",
-       "          [15, 82, 12, ..., 20, 26, 43]],\n",
+       "          [44, 67, 26, ..., 29, 96, 92],\n",
+       "          [38, 63, 87, ..., 32, 89, 62],\n",
+       "          [ 1,  5, 28, ..., 88, 84, 58]],\n",
        "\n",
        "         ...,\n",
        "\n",
-       "         [[ 6, 65, 35, ..., 46, 17, 12],\n",
-       "          [21, 22, 69, ..., 45, 58, 11]],\n",
+       "         [[ 9, 82, 60, ..., 64, 56, 98],\n",
+       "          [79, 62, 12, ..., 47, 84, 60],\n",
+       "          [83, 98,  3, ..., 91, 25, 14],\n",
+       "          [81,  7, 45, ..., 18, 78,  9]],\n",
        "\n",
-       "         [[75, 78, 24, ..., 13, 37, 66],\n",
-       "          [12, 90, 93, ..., 17, 52, 47]],\n",
+       "         [[29, 94, 53, ...,  3, 52, 47],\n",
+       "          [93, 43, 11, ..., 94, 71, 37],\n",
+       "          [72, 95, 55, ..., 83, 17, 97],\n",
+       "          [72, 77, 15, ..., 42, 62,  9]],\n",
        "\n",
-       "         [[49, 62, 62, ..., 74, 70, 19],\n",
-       "          [32,  8, 19, ..., 79, 16, 87]]]]])\n",
+       "         [[26, 76, 23, ..., 21, 92, 42],\n",
+       "          [60, 80, 37, ..., 49, 91, 88],\n",
+       "          [66, 91, 73, ..., 10, 44, 22],\n",
+       "          [75, 27, 89, ..., 37, 57, 20]]]]])\n",
        "Coordinates:\n",
        "  * mode         (mode) <U4 'car' 'bike' 'foot'\n",
        "  * day          (day) datetime64[ns] 2015-01-01 2015-01-02 ... 2015-04-10\n",
        "  * time         (time) int64 0 1 2 3 4 5 6 7 8 9 ... 15 16 17 18 19 20 21 22 23\n",
-       "    origin       (points) object MULTIPOLYGON (((-79.91995239257812 34.807918...\n",
+       "  * origin       (origin) object MULTIPOLYGON (((-81.4727554321289 36.2343559...\n",
        "  * destination  (destination) object MULTIPOLYGON (((-81.4727554321289 36.23...\n",
-       "Dimensions without coordinates: points\n",
        "Indexes:\n",
-       "    destination  ShapelySTRTreeIndex
  • origin
    GeometryIndex (crs=EPSG:4267)
    <xvec.index.GeometryIndex object at 0x16ab95890>
  • destination
    GeometryIndex (crs=EPSG:4267)
    <xvec.index.GeometryIndex object at 0x16c4de650>
    • " ], "text/plain": [ - "\n", - "array([[[[[31, 69, 46, ..., 85, 87, 41],\n", - " [63, 34, 15, ..., 62, 36, 23]],\n", - "\n", - " [[79, 46, 26, ..., 75, 89, 12],\n", - " [99, 67, 18, ..., 93, 63, 9]],\n", - "\n", - " [[ 7, 22, 63, ..., 78, 89, 64],\n", - " [60, 33, 13, ..., 39, 54, 33]],\n", + "\n", + "array([[[[[47, 51, 75, ..., 15, 54, 82],\n", + " [ 7, 68, 76, ..., 38, 3, 72],\n", + " [19, 47, 3, ..., 75, 42, 69],\n", + " [48, 94, 34, ..., 47, 54, 86]],\n", + "\n", + " [[67, 5, 73, ..., 95, 34, 44],\n", + " [80, 22, 89, ..., 93, 5, 41],\n", + " [78, 49, 9, ..., 45, 3, 9],\n", + " [87, 73, 52, ..., 32, 58, 48]],\n", + "\n", + " [[22, 57, 42, ..., 24, 60, 1],\n", + " [ 4, 3, 94, ..., 66, 40, 17],\n", + " [81, 36, 30, ..., 5, 51, 6],\n", + " [54, 70, 25, ..., 30, 18, 37]],\n", "\n", " ...,\n", "\n", - " [[22, 59, 96, ..., 54, 84, 18],\n", - " [39, 52, 31, ..., 40, 45, 43]],\n", - "\n", - " [[15, 37, 54, ..., 30, 2, 57],\n", - " [29, 90, 1, ..., 60, 21, 29]],\n", - "\n", - " [[97, 47, 98, ..., 29, 76, 77],\n", - " [43, 42, 23, ..., 68, 56, 42]]],\n", - "\n", + " [[ 6, 83, 84, ..., 50, 77, 13],\n", + " [44, 84, 52, ..., 38, 52, 47],\n", + " [65, 26, 50, ..., 12, 21, 25],\n", "...\n", - "\n", - " [[[53, 70, 33, ..., 6, 78, 35],\n", - " [62, 58, 58, ..., 35, 88, 12]],\n", - "\n", - " [[27, 77, 3, ..., 62, 44, 34],\n", - " [95, 88, 66, ..., 67, 63, 83]],\n", - "\n", - " [[43, 78, 66, ..., 22, 16, 56],\n", - " [15, 82, 12, ..., 20, 26, 43]],\n", + " [44, 67, 26, ..., 29, 96, 92],\n", + " [38, 63, 87, ..., 32, 89, 62],\n", + " [ 1, 5, 28, ..., 88, 84, 58]],\n", "\n", " ...,\n", "\n", - " [[ 6, 65, 35, ..., 46, 17, 12],\n", - " [21, 22, 69, ..., 45, 58, 11]],\n", + " [[ 9, 82, 60, ..., 64, 56, 98],\n", + " [79, 62, 12, ..., 47, 84, 60],\n", + " [83, 98, 3, ..., 91, 25, 14],\n", + " [81, 7, 45, ..., 18, 78, 9]],\n", "\n", - " [[75, 78, 24, ..., 13, 37, 66],\n", - " [12, 90, 93, ..., 17, 52, 47]],\n", + " [[29, 94, 53, ..., 3, 52, 47],\n", + " [93, 43, 11, ..., 94, 71, 37],\n", + " [72, 95, 55, ..., 83, 17, 97],\n", + " [72, 77, 15, ..., 42, 62, 9]],\n", "\n", - " [[49, 62, 62, ..., 74, 70, 19],\n", - " [32, 8, 19, ..., 79, 16, 87]]]]])\n", + " [[26, 76, 23, ..., 21, 92, 42],\n", + " [60, 80, 37, ..., 49, 91, 88],\n", + " [66, 91, 73, ..., 10, 44, 22],\n", + " [75, 27, 89, ..., 37, 57, 20]]]]])\n", "Coordinates:\n", " * mode (mode) 1\u001b[0m arr\u001b[39m.\u001b[39;49msel(origin\u001b[39m=\u001b[39;49m[shapely\u001b[39m.\u001b[39;49mPoint(\u001b[39m-\u001b[39;49m\u001b[39m80\u001b[39;49m, \u001b[39m35\u001b[39;49m), shapely\u001b[39m.\u001b[39;49mPoint(\u001b[39m-\u001b[39;49m\u001b[39m76\u001b[39;49m, \u001b[39m35.6\u001b[39;49m)], method\u001b[39m=\u001b[39;49m\u001b[39m\"\u001b[39;49m\u001b[39mintersects\u001b[39;49m\u001b[39m\"\u001b[39;49m)\n", - "File \u001b[0;32m~/mambaforge/envs/xvec_dev/lib/python3.10/site-packages/xarray/core/dataarray.py:1526\u001b[0m, in \u001b[0;36mDataArray.sel\u001b[0;34m(self, indexers, method, tolerance, drop, **indexers_kwargs)\u001b[0m\n\u001b[1;32m 1416\u001b[0m \u001b[39mdef\u001b[39;00m \u001b[39msel\u001b[39m(\n\u001b[1;32m 1417\u001b[0m \u001b[39mself\u001b[39m: T_DataArray,\n\u001b[1;32m 1418\u001b[0m indexers: Mapping[Any, Any] \u001b[39m=\u001b[39m \u001b[39mNone\u001b[39;00m,\n\u001b[0;32m (...)\u001b[0m\n\u001b[1;32m 1422\u001b[0m \u001b[39m*\u001b[39m\u001b[39m*\u001b[39mindexers_kwargs: Any,\n\u001b[1;32m 1423\u001b[0m ) \u001b[39m-\u001b[39m\u001b[39m>\u001b[39m T_DataArray:\n\u001b[1;32m 1424\u001b[0m \u001b[39m\"\"\"Return a new DataArray whose data is given by selecting index\u001b[39;00m\n\u001b[1;32m 1425\u001b[0m \u001b[39m labels along the specified dimension(s).\u001b[39;00m\n\u001b[1;32m 1426\u001b[0m \n\u001b[0;32m (...)\u001b[0m\n\u001b[1;32m 1524\u001b[0m \u001b[39m Dimensions without coordinates: points\u001b[39;00m\n\u001b[1;32m 1525\u001b[0m \u001b[39m \"\"\"\u001b[39;00m\n\u001b[0;32m-> 1526\u001b[0m ds \u001b[39m=\u001b[39m \u001b[39mself\u001b[39;49m\u001b[39m.\u001b[39;49m_to_temp_dataset()\u001b[39m.\u001b[39;49msel(\n\u001b[1;32m 1527\u001b[0m indexers\u001b[39m=\u001b[39;49mindexers,\n\u001b[1;32m 1528\u001b[0m drop\u001b[39m=\u001b[39;49mdrop,\n\u001b[1;32m 1529\u001b[0m method\u001b[39m=\u001b[39;49mmethod,\n\u001b[1;32m 1530\u001b[0m tolerance\u001b[39m=\u001b[39;49mtolerance,\n\u001b[1;32m 1531\u001b[0m \u001b[39m*\u001b[39;49m\u001b[39m*\u001b[39;49mindexers_kwargs,\n\u001b[1;32m 1532\u001b[0m )\n\u001b[1;32m 1533\u001b[0m \u001b[39mreturn\u001b[39;00m \u001b[39mself\u001b[39m\u001b[39m.\u001b[39m_from_temp_dataset(ds)\n", - "File \u001b[0;32m~/mambaforge/envs/xvec_dev/lib/python3.10/site-packages/xarray/core/dataset.py:2554\u001b[0m, in \u001b[0;36mDataset.sel\u001b[0;34m(self, indexers, method, tolerance, drop, **indexers_kwargs)\u001b[0m\n\u001b[1;32m 2493\u001b[0m \u001b[39m\"\"\"Returns a new dataset with each array indexed by tick labels\u001b[39;00m\n\u001b[1;32m 2494\u001b[0m \u001b[39malong the specified dimension(s).\u001b[39;00m\n\u001b[1;32m 2495\u001b[0m \n\u001b[0;32m (...)\u001b[0m\n\u001b[1;32m 2551\u001b[0m \u001b[39mDataArray.sel\u001b[39;00m\n\u001b[1;32m 2552\u001b[0m \u001b[39m\"\"\"\u001b[39;00m\n\u001b[1;32m 2553\u001b[0m indexers \u001b[39m=\u001b[39m either_dict_or_kwargs(indexers, indexers_kwargs, \u001b[39m\"\u001b[39m\u001b[39msel\u001b[39m\u001b[39m\"\u001b[39m)\n\u001b[0;32m-> 2554\u001b[0m query_results \u001b[39m=\u001b[39m map_index_queries(\n\u001b[1;32m 2555\u001b[0m \u001b[39mself\u001b[39;49m, indexers\u001b[39m=\u001b[39;49mindexers, method\u001b[39m=\u001b[39;49mmethod, tolerance\u001b[39m=\u001b[39;49mtolerance\n\u001b[1;32m 2556\u001b[0m )\n\u001b[1;32m 2558\u001b[0m \u001b[39mif\u001b[39;00m drop:\n\u001b[1;32m 2559\u001b[0m no_scalar_variables \u001b[39m=\u001b[39m {}\n", - "File \u001b[0;32m~/mambaforge/envs/xvec_dev/lib/python3.10/site-packages/xarray/core/indexing.py:183\u001b[0m, in \u001b[0;36mmap_index_queries\u001b[0;34m(obj, indexers, method, tolerance, **indexers_kwargs)\u001b[0m\n\u001b[1;32m 181\u001b[0m results\u001b[39m.\u001b[39mappend(IndexSelResult(labels))\n\u001b[1;32m 182\u001b[0m \u001b[39melse\u001b[39;00m:\n\u001b[0;32m--> 183\u001b[0m results\u001b[39m.\u001b[39mappend(index\u001b[39m.\u001b[39;49msel(labels, \u001b[39m*\u001b[39;49m\u001b[39m*\u001b[39;49moptions)) \u001b[39m# type: ignore[call-arg]\u001b[39;00m\n\u001b[1;32m 185\u001b[0m merged \u001b[39m=\u001b[39m merge_sel_results(results)\n\u001b[1;32m 187\u001b[0m \u001b[39m# drop dimension coordinates found in dimension indexers\u001b[39;00m\n\u001b[1;32m 188\u001b[0m \u001b[39m# (also drop multi-index if any)\u001b[39;00m\n\u001b[1;32m 189\u001b[0m \u001b[39m# (.sel() already ensures alignment)\u001b[39;00m\n", - "File \u001b[0;32m~/Git/xvec/xvec/strtree.py:40\u001b[0m, in \u001b[0;36mShapelySTRTreeIndex.sel\u001b[0;34m(self, labels, method, tolerance)\u001b[0m\n\u001b[1;32m 38\u001b[0m \u001b[39mif\u001b[39;00m method \u001b[39mis\u001b[39;00m \u001b[39mnot\u001b[39;00m \u001b[39mNone\u001b[39;00m \u001b[39mand\u001b[39;00m method \u001b[39m!=\u001b[39m \u001b[39m\"\u001b[39m\u001b[39mnearest\u001b[39m\u001b[39m\"\u001b[39m:\n\u001b[1;32m 39\u001b[0m \u001b[39mif\u001b[39;00m \u001b[39mnot\u001b[39;00m \u001b[39misinstance\u001b[39m(label, shapely\u001b[39m.\u001b[39mGeometry):\n\u001b[0;32m---> 40\u001b[0m \u001b[39mraise\u001b[39;00m \u001b[39mValueError\u001b[39;00m(\n\u001b[1;32m 41\u001b[0m \u001b[39m\"\u001b[39m\u001b[39mselection with another method than nearest only supports \u001b[39m\u001b[39m\"\u001b[39m\n\u001b[1;32m 42\u001b[0m \u001b[39m\"\u001b[39m\u001b[39ma single geometry as input label.\u001b[39m\u001b[39m\"\u001b[39m\n\u001b[1;32m 43\u001b[0m )\n\u001b[1;32m 45\u001b[0m \u001b[39mif\u001b[39;00m \u001b[39misinstance\u001b[39m(label, xarray\u001b[39m.\u001b[39mDataArray):\n\u001b[1;32m 46\u001b[0m label_array \u001b[39m=\u001b[39m label\u001b[39m.\u001b[39m_variable\u001b[39m.\u001b[39m_data\n", - "\u001b[0;31mValueError\u001b[0m: selection with another method than nearest only supports a single geometry as input label." - ] - } - ], - "source": [ - "arr.sel(origin=[shapely.Point(-80, 35), shapely.Point(-76, 35.6)], method=\"intersects\")" - ] + "outputs": [], + "source": [] } ], "metadata": { "kernelspec": { - "display_name": "Python 3.10.6 ('xvec_dev')", + "display_name": "Python [conda env:xvec_dev]", "language": "python", - "name": "python3" + "name": "conda-env-xvec_dev-py" }, "language_info": { "codemirror_mode": { @@ -3925,9 +5272,8 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.10.6" + "version": "3.11.0" }, - "orig_nbformat": 4, "vscode": { "interpreter": { "hash": "0d73ffc659e32f9becefa85ab414fb538900c004037e66c85673d8fcbf126471" @@ -3935,5 +5281,5 @@ } }, "nbformat": 4, - "nbformat_minor": 2 + "nbformat_minor": 4 } diff --git a/pyproject.toml b/pyproject.toml index b14b7de..1333517 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -28,6 +28,7 @@ classifiers = [ requires-python = ">=3.8" dependencies = [ "xarray >= 2022.11.0", + "pyproj >= 3.0.0", "shapely >= 2.0b1", ] diff --git a/xvec/__init__.py b/xvec/__init__.py index 5aca488..15ad2cb 100644 --- a/xvec/__init__.py +++ b/xvec/__init__.py @@ -1,9 +1,9 @@ from importlib.metadata import PackageNotFoundError, version -from .strtree import ShapelySTRTreeIndex # noqa +from .index import GeometryIndex # noqa try: - __version__ = version("package-name") + __version__ = version("xvec") except PackageNotFoundError: # noqa # package is not installed pass diff --git a/xvec/index.py b/xvec/index.py new file mode 100644 index 0000000..1c4eee6 --- /dev/null +++ b/xvec/index.py @@ -0,0 +1,265 @@ +from __future__ import annotations + +import warnings +from typing import Any, Hashable, Iterable, Mapping, Sequence + +import numpy as np +import pandas as pd +import shapely +from pyproj import CRS +from xarray import DataArray, Variable, get_options +from xarray.core.indexing import IndexSelResult +from xarray.indexes import Index, PandasIndex + + +def _format_crs(crs: CRS | None, max_width: int = 50) -> str: + if crs is not None: + srs = crs.to_string() + else: + srs = "None" + + return srs if len(srs) <= max_width else " ".join([srs[:max_width], "..."]) + + +def _get_common_crs(crs_set: set[CRS | None]): + # code taken from geopandas (BSD-3 Licence) + + crs_not_none = [crs for crs in crs_set if crs is not None] + names = [crs.name for crs in crs_not_none] + + if len(crs_not_none) == 0: + return None + if len(crs_not_none) == 1: + if len(crs_set) != 1: + warnings.warn( + "CRS not set for some of the concatenation inputs. " + f"Setting output's CRS as {names[0]} " + "(the single non-null CRS provided)." + ) + return crs_not_none[0] + + raise ValueError( + f"cannot determine common CRS for concatenation inputs, got {names}. " + # "Use `to_crs()` to transform geometries to the same CRS before merging." + ) + + +class GeometryIndex(Index): + """An CRS-aware, Xarray-compatible index for vector geometries. + + This index can be set from any 1-dimensional coordinate of + (shapely 2.0) :class:`shapely.Geometry` elements. + + It provides all the basic functionality of an + :class:`xarray.indexes.PandasIndex`. In addition, it allows spatial + filtering based on geometries (powered by :class:`shapely.STRtree`). + + Parameters + ---------- + index : :class:`xarray.indexes.PandasIndex` + An Xarray (pandas) index built from an array-like of + :class:`shapely.Geometry` objects. + crs : :class:`pyproj.crs.CRS` or any, optional + The coordinate reference system. Any value accepted by + :meth:`pyproj.crs.CRS.from_user_input`. + + """ + + _index: PandasIndex + _sindex: shapely.STRtree | None + _crs: CRS | None + + def __init__(self, index: PandasIndex, crs: CRS | Any | None = None): + if not np.all(shapely.is_geometry(index.index)): + raise ValueError("array must contain shapely.Geometry objects") + + if crs is not None: + crs = CRS.from_user_input(crs) + + self._crs = crs + self._index = index + self._sindex = None + + @property + def crs(self) -> CRS | None: + """Returns the coordinate reference system of the index as a + :class:`pyproj.crs.CRS` object. + """ + return self._crs + + @property + def sindex(self) -> shapely.STRtree: + """Returns the spatial index, i.e., a :class:`shapely.STRtree` object. + + It may build the index before returning it if it hasn't been built before. + """ + if self._sindex is None: + self._sindex = shapely.STRtree(self._index.index) + return self._sindex + + def _check_crs(self, other_crs: CRS | None, allow_none: bool = False) -> bool: + """Check if the index's projection is the same than the given one. + If allow_none is True, empty CRS is treated as the same. + """ + if allow_none: + if self.crs is None or other_crs is None: + return True + if not self.crs == other_crs: + return False + return True + + def _crs_mismatch_raise( + self, other_crs: CRS | None, warn: bool = False, stacklevel: int = 3 + ): + """Raise a CRS mismatch error or warning with the information + on the assigned CRS. + """ + srs = _format_crs(self.crs, max_width=50) + other_srs = _format_crs(other_crs, max_width=50) + + # TODO: expand message with reproject suggestion + msg = ( + "CRS mismatch between the CRS of index geometries " + "and the CRS of input geometries.\n" + f"Index CRS: {srs}\n" + f"Input CRS: {other_srs}\n" + ) + + if warn: + warnings.warn(msg, UserWarning, stacklevel=stacklevel) + else: + raise ValueError(msg) + + @classmethod + def from_variables( + cls, + variables: Mapping[Any, Variable], + *, + options: Mapping[str, Any], + ): + # TODO: try getting CRS from coordinate attrs or GeometryArray or SRID + + index = PandasIndex.from_variables(variables, options={}) + return cls(index, crs=options.get("crs")) + + @classmethod + def concat( + cls, + indexes: Sequence[GeometryIndex], + dim: Hashable, + positions: Iterable[Iterable[int]] | None = None, + ) -> GeometryIndex: + crs_set = {idx.crs for idx in indexes} + crs = _get_common_crs(crs_set) + + indexes_ = [idx._index for idx in indexes] + index = PandasIndex.concat(indexes_, dim, positions) + return cls(index, crs) + + def create_variables( + self, variables: Mapping[Any, Variable] | None = None + ) -> dict[Hashable, Variable]: + return self._index.create_variables(variables) + + def to_pandas_index(self) -> pd.Index: + return self._index.index + + def isel(self, indexers: Mapping[Any, Any]): + index = self._index.isel(indexers) + + if index is not None: + return type(self)(index, self.crs) + else: + return None + + def _sel_sindex(self, labels, method, tolerance): + # only one entry expected + assert len(labels) == 1 + label = next(iter(labels.values())) + + if method != "nearest": + if not isinstance(label, shapely.Geometry): + raise ValueError( + "selection with another method than nearest only supports " + "a single geometry as input label." + ) + + if isinstance(label, DataArray): + label_array = label._variable._data + elif isinstance(label, Variable): + label_array = label._data + elif isinstance(label, shapely.Geometry): + label_array = np.array([label]) + else: + label_array = np.array(label) + + if not np.all(shapely.is_geometry(label_array)): + raise ValueError("labels must be shapely.Geometry objects") + + if method == "nearest": + indices = self.sindex.nearest(label_array) + else: + indices = self.sindex.query(label, predicate=method, distance=tolerance) + + # attach dimension names and/or coordinates to positional indexer + if isinstance(label, Variable): + indices = Variable(label.dims, indices) + elif isinstance(label, DataArray): + indices = DataArray(indices, coords=label._coords, dims=label.dims) + + return IndexSelResult({self._index.dim: indices}) + + def sel( + self, labels: dict[Any, Any], method=None, tolerance=None + ) -> IndexSelResult: + if method is None: + return self._index.sel(labels) + else: + # We reuse here `method` and `tolerance` options of + # `xarray.indexes.PandasIndex` as `predicate` and `distance` + # options when `labels` is a single geometry. + # Xarray currently doesn't support custom options + # (see https://github.com/pydata/xarray/issues/7099) + return self._sel_sindex(labels, method, tolerance) + + def equals(self, other: Index) -> bool: + if not isinstance(other, GeometryIndex): + return False + if not self._check_crs(other.crs, allow_none=True): + return False + return self._index.equals(other._index) + + def join( + self: GeometryIndex, other: GeometryIndex, how: str = "inner" + ) -> GeometryIndex: + if not self._check_crs(other.crs, allow_none=True): + self._crs_mismatch_raise(other.crs) + + index = self._index.join(other._index, how=how) + return type(self)(index, self.crs) + + def reindex_like( + self, other: GeometryIndex, method=None, tolerance=None + ) -> dict[Hashable, Any]: + if not self._check_crs(other.crs, allow_none=True): + self._crs_mismatch_raise(other.crs) + + return self._index.reindex_like( + other._index, method=method, tolerance=tolerance + ) + + def roll(self, shifts: Mapping[Any, int]) -> GeometryIndex: + index = self._index.roll(shifts) + return type(self)(index, self.crs) + + def rename(self, name_dict, dims_dict): + index = self._index.rename(name_dict, dims_dict) + return type(self)(index, self.crs) + + def _repr_inline_(self, max_width: int): + # TODO: remove when fixed in XArray + if max_width is None: + max_width = get_options()["display_width"] + + srs = _format_crs(self.crs, max_width=max_width) + return f"{self.__class__.__name__} (crs={srs})" diff --git a/xvec/strtree.py b/xvec/strtree.py deleted file mode 100644 index 592e8e1..0000000 --- a/xvec/strtree.py +++ /dev/null @@ -1,72 +0,0 @@ -import numpy -import shapely -import xarray -from xarray.core.indexes import IndexSelResult -from xarray.indexes import Index - - -class ShapelySTRTreeIndex(Index): - def __init__(self, array, dim, crs): - assert numpy.all(shapely.is_geometry(array)) - - # only support 1-d coordinate for now - assert len(array.shape) == 1 - - self._tree = shapely.STRtree(numpy.ravel(array)) - self.dim = dim - self.crs = crs - - @classmethod - def from_variables(cls, variables, *, options): - # only supports one coordinate of shapely geometries - assert len(variables) == 1 - var = next(iter(variables.values())) - - return cls(var._data, var.dims[0], options["crs"]) - - def sel(self, labels, method=None, tolerance=None): - # We reuse here `method` and `tolerance` options of - # `xarray.indexes.PandasIndex` as `predicate` and `distance` - # options when `labels` is a single geometry. - # Xarray currently doesn't support custom options - # (see https://github.com/pydata/xarray/issues/7099) - - # only one coordinate supported - assert len(labels) == 1 - label = next(iter(labels.values())) - - if method is not None and method != "nearest": - if not isinstance(label, shapely.Geometry): - raise ValueError( - "selection with another method than nearest only supports " - "a single geometry as input label." - ) - - if isinstance(label, xarray.DataArray): - label_array = label._variable._data - elif isinstance(label, xarray.Variable): - label_array = label._data - elif isinstance(label, shapely.Geometry): - label_array = numpy.array([label]) - else: - label_array = numpy.array(label) - - # check for possible CRS of geometry labels - # (by default assume same CRS than the index) - if hasattr(label_array, "crs") and label_array.crs != self.crs: - raise ValueError("conflicting CRS for input geometries") - - assert numpy.all(shapely.is_geometry(label_array)) - - if method is None or method == "nearest": - indices = self._tree.nearest(label_array) - else: - indices = self._tree.query(label, predicate=method, distance=tolerance) - - # attach dimension names and/or coordinates to positional indexer - if isinstance(label, xarray.Variable): - indices = xarray.Variable(label.dims, indices) - elif isinstance(label, xarray.DataArray): - indices = xarray.DataArray(indices, coords=label._coords, dims=label.dims) - - return IndexSelResult({self.dim: indices}) diff --git a/xvec/tests/test_index.py b/xvec/tests/test_index.py new file mode 100644 index 0000000..efd9fdd --- /dev/null +++ b/xvec/tests/test_index.py @@ -0,0 +1,209 @@ +import numpy as np +import pytest +import shapely +import xarray as xr +from pyproj import CRS + +from xvec import GeometryIndex + + +@pytest.fixture(scope="session") +def geom_array(): + return np.array([shapely.Point(1, 2), shapely.Point(3, 4)]) + + +@pytest.fixture(scope="session") +def geom_dataset_no_index(geom_array): + # a dataset with a geometry coordinate but no index + ds = xr.Dataset(coords={"geom": geom_array}) + return ds.drop_indexes("geom") + + +@pytest.fixture(scope="session") +def geom_dataset(geom_dataset_no_index): + # a dataset with a geometry coordinate baked by a GeometryIndex + crs = CRS.from_user_input(26915) + return geom_dataset_no_index.set_xindex("geom", GeometryIndex, crs=crs) + + +@pytest.fixture(scope="session") +def geom_dataset_no_crs(geom_dataset_no_index): + # a dataset with a geometry coordinate baked by a GeometryIndex (no CRS) + return geom_dataset_no_index.set_xindex("geom", GeometryIndex) + + +@pytest.fixture(scope="session") +def first_geom_dataset(geom_dataset, geom_array): + return ( + xr.Dataset(coords={"geom": [geom_array[0]]}) + .drop_indexes("geom") + .set_xindex("geom", GeometryIndex, crs=geom_dataset.xindexes["geom"].crs) + ) + + +def test_set_index(geom_dataset_no_index): + crs = CRS.from_user_input(26915) + ds = geom_dataset_no_index.set_xindex("geom", GeometryIndex, crs=crs) + + # test properties + assert isinstance(ds.xindexes["geom"], GeometryIndex) + assert ds.xindexes["geom"].crs == crs + np.testing.assert_array_equal(ds.xindexes["geom"].sindex.geometries, ds.geom.values) + + # test `GeometryIndex.create_variables` + assert ds.geom.variable._data.array is ds.xindexes["geom"]._index.index + + no_crs_ds = geom_dataset_no_index.set_xindex("geom", GeometryIndex) + assert no_crs_ds.xindexes["geom"].crs is None + + no_geom_ds = xr.Dataset(coords={"no_geom": ("x", [0, 1, 2])}) + with pytest.raises(ValueError, match="array must contain shapely.Geometry objects"): + no_geom_ds.set_xindex("no_geom", GeometryIndex, crs=crs) + + +def test_concat(geom_dataset, geom_array, geom_dataset_no_index, geom_dataset_no_crs): + expected = ( + xr.Dataset(coords={"geom": np.concatenate([geom_array, geom_array])}) + .drop_indexes("geom") + .set_xindex("geom", GeometryIndex, crs=geom_dataset.xindexes["geom"].crs) + ) + actual = xr.concat([geom_dataset, geom_dataset], "geom") + xr.testing.assert_identical(actual, expected) + + # different CRS + crs = CRS.from_user_input(4267) + geom_dataset_alt = geom_dataset_no_index.set_xindex("geom", GeometryIndex, crs=crs) + + with pytest.raises(ValueError, match="cannot determine common CRS"): + xr.concat([geom_dataset, geom_dataset_alt], "geom") + + # no CRS + expected = ( + xr.Dataset(coords={"geom": np.concatenate([geom_array, geom_array])}) + .drop_indexes("geom") + .set_xindex("geom", GeometryIndex) + ) + actual = xr.concat([geom_dataset_no_crs, geom_dataset_no_crs], "geom") + xr.testing.assert_identical(actual, expected) + + # mixed CRS / no CRS + with pytest.warns( + UserWarning, match="CRS not set for some of the concatenation inputs" + ): + xr.concat([geom_dataset, geom_dataset_no_crs], "geom") + + +def test_to_pandas_index(geom_dataset): + index = geom_dataset.xindexes["geom"] + assert index.to_pandas_index() is index._index.index + + +def test_isel(geom_dataset, first_geom_dataset): + actual = geom_dataset.isel(geom=[0]) + xr.testing.assert_identical(actual, first_geom_dataset) + + # scalar selection + actual = geom_dataset.isel(geom=0) + assert len(actual.geom.dims) == 0 + assert "geom" not in actual.xindexes + + +def test_sel_strict(geom_dataset, geom_array, first_geom_dataset): + actual = geom_dataset.sel(geom=[geom_array[0]]) + xr.testing.assert_identical(actual, first_geom_dataset) + + +@pytest.mark.parametrize( + "label", + [ + shapely.Point(1, 1), + [shapely.Point(1, 1)], + xr.Variable("geom", [shapely.Point(1, 1)]), + xr.DataArray([shapely.Point(1, 1)], dims="geom"), + ], +) +def test_sel_nearest(geom_dataset, first_geom_dataset, label): + actual = geom_dataset.sel(geom=label, method="nearest") + xr.testing.assert_identical(actual, first_geom_dataset) + + +def test_sel_nearest_error(geom_dataset): + with pytest.raises(ValueError, match="labels must be shapely.Geometry objects"): + geom_dataset.sel(geom=[0], method="nearest") + + +def test_sel_query(geom_dataset, first_geom_dataset): + actual = geom_dataset.sel(geom=shapely.box(0, 0, 2, 2), method="intersects") + xr.testing.assert_identical(actual, first_geom_dataset) + + +def test_equals( + geom_dataset, geom_dataset_no_index, first_geom_dataset, geom_dataset_no_crs +): + assert geom_dataset.xindexes["geom"].equals(geom_dataset.xindexes["geom"]) + + # different index types + other = xr.Dataset(coords={"geom": [0, 1]}) + assert not geom_dataset.xindexes["geom"].equals(other.xindexes["geom"]) + + # different CRS + crs = CRS.from_user_input(4267) + other = geom_dataset_no_index.set_xindex("geom", GeometryIndex, crs=crs) + assert not geom_dataset.xindexes["geom"].equals(other.xindexes["geom"]) + + # no CRS + assert geom_dataset_no_crs.xindexes["geom"].equals( + geom_dataset_no_crs.xindexes["geom"] + ) + + # different geometries + assert not geom_dataset.xindexes["geom"].equals(first_geom_dataset.xindexes["geom"]) + + +def test_align( + geom_dataset, first_geom_dataset, geom_dataset_no_index, geom_dataset_no_crs +): + # test both GeometryIndex's `join` and `reindex_like` + aligned = xr.align(geom_dataset, first_geom_dataset, join="inner") + assert all([ds.identical(first_geom_dataset) for ds in aligned]) + + # test conflicting CRS + crs = CRS.from_user_input(4267) + geom_dataset_alt = geom_dataset_no_index.set_xindex("geom", GeometryIndex, crs=crs) + + with pytest.raises(ValueError, match="CRS mismatch"): + xr.align(geom_dataset_alt, first_geom_dataset, join="inner") + + with pytest.raises(ValueError, match="CRS mismatch"): + first_geom_dataset.reindex_like(geom_dataset_alt) + + # no CRS + first_geom_dataset_no_crs = geom_dataset_no_crs.isel(geom=[0]) + aligned = xr.align(geom_dataset_no_crs, first_geom_dataset_no_crs, join="inner") + assert all([ds.identical(first_geom_dataset_no_crs) for ds in aligned]) + + +def test_roll(geom_dataset, geom_array): + expected = ( + xr.Dataset(coords={"geom": np.roll(geom_array, 1)}) + .drop_indexes("geom") + .set_xindex("geom", GeometryIndex, crs=geom_dataset.xindexes["geom"].crs) + ) + actual = geom_dataset.roll(geom=1, roll_coords=True) + xr.testing.assert_identical(actual, expected) + + +def test_rename(geom_dataset): + ds = geom_dataset.rename_vars(geom="renamed") + assert ds.xindexes["renamed"]._index.index.name == "renamed" + assert ds.xindexes["renamed"]._index.dim == "geom" + + +def test_repr_inline(geom_dataset, geom_dataset_no_crs): + actual = geom_dataset.xindexes["geom"]._repr_inline_(70) + expected = "GeometryIndex (crs=EPSG:26915)" + assert actual == expected + + actual = geom_dataset_no_crs.xindexes["geom"]._repr_inline_(70) + expected = "GeometryIndex (crs=None)" + assert actual == expected diff --git a/xvec/tests/test_strtree.py b/xvec/tests/test_strtree.py deleted file mode 100644 index 7e5d728..0000000 --- a/xvec/tests/test_strtree.py +++ /dev/null @@ -1,2 +0,0 @@ -def test_import_dummy(): - from xvec import ShapelySTRTreeIndex # noqa