From 56b2cad2cf43fa3171a308599dab88acefaefe9d Mon Sep 17 00:00:00 2001 From: Benoit Bovy Date: Thu, 12 Dec 2024 14:29:53 +0100 Subject: [PATCH 1/2] rename set_crs -> assign_crs --- xproj/accessor.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/xproj/accessor.py b/xproj/accessor.py index eac36b9..9ba220d 100644 --- a/xproj/accessor.py +++ b/xproj/accessor.py @@ -234,7 +234,7 @@ def crs(self) -> pyproj.CRS | None: else: return crs_index.crs - def set_crs( + def assign_crs( self, coord_name_crs: Mapping[Hashable, Any] | None = None, allow_override: bool = False, @@ -253,7 +253,7 @@ def set_crs( :meth:`pyproj.CRS.from_user_input() ` such as an authority string (e.g. ``"EPSG:4326"``), EPSG code (e.g. ``4326``) or a WKT string. - If the coordinate(s) doesn't exist they will be created. + If the coordinate(s) doesn't exist it will be created. Only one item is currently allowed. allow_override : bool, default False Allow to replace the index if the coordinates already have an index. From 5d293dab67adec83f3e1dc262b2160cb5a8a8195 Mon Sep 17 00:00:00 2001 From: Benoit Bovy Date: Thu, 12 Dec 2024 14:30:09 +0100 Subject: [PATCH 2/2] update demo notebook --- examples/demo.ipynb | 1200 ++++++++++++++++++++++++++++++++++++------- 1 file changed, 1005 insertions(+), 195 deletions(-) diff --git a/examples/demo.ipynb b/examples/demo.ipynb index 66f1ad4..bed7bdd 100644 --- a/examples/demo.ipynb +++ b/examples/demo.ipynb @@ -420,55 +420,35 @@ " stroke: currentColor;\n", " fill: currentColor;\n", "}\n", - "
<xarray.Dataset> Size: 400B\n",
-       "Dimensions:  (space: 50)\n",
+       "
<xarray.Dataset> Size: 80B\n",
+       "Dimensions:  (space: 10)\n",
        "Coordinates:\n",
-       "  * points   (space) object 400B POINT (72.14265652432513 -68.90724731796287)...\n",
+       "  * points   (space) object 80B POINT (171.61191601133118 -10.618898030367944...\n",
        "Dimensions without coordinates: space\n",
        "Data variables:\n",
        "    *empty*\n",
        "Indexes:\n",
-       "    points   GeometryIndex (crs=None)
    • " ], "text/plain": [ - " Size: 400B\n", - "Dimensions: (space: 50)\n", + " Size: 80B\n", + "Dimensions: (space: 10)\n", "Coordinates:\n", - " * points (space) object 400B POINT (72.14265652432513 -68.90724731796287)...\n", + " * points (space) object 80B POINT (171.61191601133118 -10.618898030367944...\n", "Dimensions without coordinates: space\n", "Data variables:\n", " *empty*\n", @@ -482,7 +462,8 @@ } ], "source": [ - "random_points = shapely.points(np.random.uniform(-180, 180, 50), np.random.uniform(-90, 90, 50))\n", + "rng = np.random.default_rng(1234)\n", + "random_points = shapely.points(rng.uniform(-180, 180, 10), rng.uniform(-90, 90, 10))\n", "\n", "ds_vec = xr.Dataset(coords={\"points\": (\"space\", random_points)}).set_xindex(\n", " \"points\", xvec.GeometryIndex\n", @@ -498,9 +479,9 @@ "source": [ "## Set the CRS\n", "\n", - "`.proj.set_crs()` can be used to set (or update) the CRS: it creates a new scalar coordinate (if it doesn't exist yet) with a `xproj.CRSIndex`.\n", + "`.proj.assign_crs()` can be used to set the CRS: it creates a new scalar coordinate (if it doesn't exist yet) with a `xproj.CRSIndex`.\n", "\n", - "Note: the name of the CRS coordinate is abritrary." + "Note: the name of the spatial reference coordinate is abritrary (\"spatial_ref\" is a common name)." ] }, { @@ -882,51 +863,31 @@ " stroke: currentColor;\n", " fill: currentColor;\n", "}\n", - "
    <xarray.Dataset> Size: 408B\n",
-       "Dimensions:      (space: 50)\n",
+       "
    <xarray.Dataset> Size: 88B\n",
+       "Dimensions:      (space: 10)\n",
        "Coordinates:\n",
-       "  * points       (space) object 400B POINT (72.14265652432513 -68.90724731796...\n",
+       "  * points       (space) object 80B POINT (171.61191601133118 -10.61889803036...\n",
        "  * spatial_ref  int64 8B 0\n",
        "Dimensions without coordinates: space\n",
        "Data variables:\n",
        "    *empty*\n",
        "Indexes:\n",
        "    points       GeometryIndex (crs=None)\n",
-       "    spatial_ref  CRSIndex (crs=EPSG:4326)
    • " ], "text/plain": [ - " Size: 408B\n", - "Dimensions: (space: 50)\n", + " Size: 88B\n", + "Dimensions: (space: 10)\n", "Coordinates:\n", - " * points (space) object 400B POINT (72.14265652432513 -68.90724731796...\n", + " * points (space) object 80B POINT (171.61191601133118 -10.61889803036...\n", " * spatial_ref int64 8B 0\n", "Dimensions without coordinates: space\n", "Data variables:\n", @@ -960,21 +921,11 @@ } ], "source": [ - "ds_wgs84 = ds_vec.proj.set_crs(spatial_ref=\"epsg:4326\")\n", + "ds_wgs84 = ds_vec.proj.assign_crs(spatial_ref=\"epsg:4326\")\n", "\n", "ds_wgs84" ] }, - { - "cell_type": "markdown", - "id": "63bf6ab4-e01b-4429-a196-bca36ce2d429", - "metadata": {}, - "source": [ - "## \"Link\" the CRS coordinate with other, indexed geospatial coordinates\n", - "\n", - "TODO." - ] - }, { "cell_type": "markdown", "id": "585608fd-eb1a-42a2-ad47-75e011dd3f8d", @@ -1078,7 +1029,7 @@ "metadata": {}, "outputs": [], "source": [ - "ds_pseudo_mercator = ds_wgs84.proj.set_crs(spatial_ref=\"epsg:3857\", allow_override=True)" + "ds_pseudo_mercator = ds_wgs84.proj.assign_crs(spatial_ref=\"epsg:3857\", allow_override=True)" ] }, { @@ -1119,18 +1070,16 @@ }, { "cell_type": "markdown", - "id": "5f4fc990-ab39-443f-8f9b-a0b5130706c2", + "id": "be12fc95-e5d3-422a-8e94-256e8f4d704c", "metadata": {}, "source": [ - "## Unset the CRS\n", - "\n", - "Just drop the CRS index and/or coordinate" + "This also works seamlessly with `xarray.concat()`:" ] }, { "cell_type": "code", "execution_count": 8, - "id": "83f828d1-bdac-42d3-949f-95c21830a21a", + "id": "175d913a-4000-4965-8a76-63f63a64019e", "metadata": {}, "outputs": [ { @@ -1506,60 +1455,47 @@ " stroke: currentColor;\n", " fill: currentColor;\n", "}\n", - "
    <xarray.Dataset> Size: 400B\n",
-       "Dimensions:  (space: 50)\n",
+       "
    <xarray.Dataset> Size: 32B\n",
+       "Dimensions:      (space: 3)\n",
        "Coordinates:\n",
-       "  * points   (space) object 400B POINT (72.14265652432513 -68.90724731796287)...\n",
+       "  * points       (space) object 24B POINT (171.61191601133118 -10.61889803036...\n",
+       "  * spatial_ref  int64 8B 0\n",
        "Dimensions without coordinates: space\n",
        "Data variables:\n",
        "    *empty*\n",
        "Indexes:\n",
-       "    points   GeometryIndex (crs=None)
    "
+       "    points       GeometryIndex (crs=None)\n",
+       "    spatial_ref  CRSIndex (crs=EPSG:4326)
    " ], "text/plain": [ - " Size: 400B\n", - "Dimensions: (space: 50)\n", + " Size: 32B\n", + "Dimensions: (space: 3)\n", "Coordinates:\n", - " * points (space) object 400B POINT (72.14265652432513 -68.90724731796287)...\n", + " * points (space) object 24B POINT (171.61191601133118 -10.61889803036...\n", + " * spatial_ref int64 8B 0\n", "Dimensions without coordinates: space\n", "Data variables:\n", " *empty*\n", "Indexes:\n", - " points GeometryIndex (crs=None)" + " points GeometryIndex (crs=None)\n", + " spatial_ref CRSIndex (crs=EPSG:4326)" ] }, "execution_count": 8, @@ -1568,23 +1504,21 @@ } ], "source": [ - "ds_no_crs = ds_wgs84.drop_vars(\"spatial_ref\")\n", - "\n", - "ds_no_crs" + "xr.concat([ds_wgs84.isel(space=[0, 1]), ds_wgs84.isel(space=[-1])], \"space\")" ] }, { "cell_type": "markdown", - "id": "fff6acec-9d29-443e-a29c-802833b4388d", + "id": "cb9a1bc2-2052-4b62-a486-551a2dc4ef37", "metadata": {}, "source": [ - "It is possible to combine datasets with and without a defined CRS. The resulting dataset will have the common CRS found among all datasets. " + "It is possible to combine heterogeneous geospatial Datasets (e.g., raster, vector, grid, mesh) as long as they all have a spatial reference coordinate with the same name and with a `CRSIndex`." ] }, { "cell_type": "code", "execution_count": 9, - "id": "d142186b-9f5a-4809-b775-66e35e7af175", + "id": "aaf7ed29-66c5-42f6-872f-12ea0792c63d", "metadata": {}, "outputs": [ { @@ -1960,51 +1894,36 @@ " stroke: currentColor;\n", " fill: currentColor;\n", "}\n", - "
    <xarray.Dataset> Size: 408B\n",
-       "Dimensions:      (space: 50)\n",
+       "
    <xarray.Dataset> Size: 328B\n",
+       "Dimensions:      (space: 10, lat: 10, lon: 20)\n",
        "Coordinates:\n",
-       "  * points       (space) object 400B POINT (72.14265652432513 -68.90724731796...\n",
+       "  * points       (space) object 80B POINT (171.61191601133118 -10.61889803036...\n",
        "  * spatial_ref  int64 8B 0\n",
+       "  * lat          (lat) float64 80B -90.0 -70.0 -50.0 -30.0 ... 50.0 70.0 90.0\n",
+       "  * lon          (lon) float64 160B -180.0 -161.1 -142.1 ... 142.1 161.1 180.0\n",
        "Dimensions without coordinates: space\n",
        "Data variables:\n",
        "    *empty*\n",
        "Indexes:\n",
        "    points       GeometryIndex (crs=None)\n",
-       "    spatial_ref  CRSIndex (crs=EPSG:4326)
  • lat
    PandasIndex
    PandasIndex(Index([-90.0, -70.0, -50.0, -30.0, -10.0, 10.0, 30.0, 50.0, 70.0, 90.0], dtype='float64', name='lat'))
  • lon
    PandasIndex
    PandasIndex(Index([             -180.0, -161.05263157894737, -142.10526315789474,\n",
+       "       -123.15789473684211, -104.21052631578948,  -85.26315789473685,\n",
+       "        -66.31578947368422,  -47.36842105263159,  -28.42105263157896,\n",
+       "         -9.47368421052633,     9.4736842105263,   28.42105263157893,\n",
+       "         47.36842105263156,   66.31578947368419,   85.26315789473682,\n",
+       "        104.21052631578948,  123.15789473684208,  142.10526315789468,\n",
+       "        161.05263157894734,               180.0],\n",
+       "      dtype='float64', name='lon'))
    • " ], "text/plain": [ - " Size: 408B\n", - "Dimensions: (space: 50)\n", + " Size: 328B\n", + "Dimensions: (space: 10, lat: 10, lon: 20)\n", "Coordinates:\n", - " * points (space) object 400B POINT (72.14265652432513 -68.90724731796...\n", + " * points (space) object 80B POINT (171.61191601133118 -10.61889803036...\n", " * spatial_ref int64 8B 0\n", + " * lat (lat) float64 80B -90.0 -70.0 -50.0 -30.0 ... 50.0 70.0 90.0\n", + " * lon (lon) float64 160B -180.0 -161.1 -142.1 ... 142.1 161.1 180.0\n", "Dimensions without coordinates: space\n", "Data variables:\n", " *empty*\n", @@ -2038,16 +1966,898 @@ } ], "source": [ - "ds_wgs84 + ds_no_crs" + "# lat-lon rectilinear grid\n", + "\n", + "ds_grid = xr.Dataset(coords={\"lat\": np.linspace(-90, 90, 10), \"lon\": np.linspace(-180, 180, 20)})\n", + "ds_grid_wgs84 = ds_grid.proj.assign_crs(spatial_ref=\"epsg:4326\")\n", + "\n", + "# merge the lat-lon grid with the vector data cube\n", + "xr.merge([ds_wgs84, ds_grid_wgs84])" + ] + }, + { + "cell_type": "markdown", + "id": "5f4fc990-ab39-443f-8f9b-a0b5130706c2", + "metadata": {}, + "source": [ + "## Unset the CRS\n", + "\n", + "Just drop the CRS index and/or coordinate" ] }, { "cell_type": "code", - "execution_count": null, - "id": "f51552bb-654d-4187-bf30-67c3fe96a118", + "execution_count": 10, + "id": "83f828d1-bdac-42d3-949f-95c21830a21a", "metadata": {}, - "outputs": [], - "source": [] + "outputs": [ + { + "data": { + "text/html": [ + "
    \n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "
    <xarray.Dataset> Size: 80B\n",
+       "Dimensions:  (space: 10)\n",
+       "Coordinates:\n",
+       "  * points   (space) object 80B POINT (171.61191601133118 -10.618898030367944...\n",
+       "Dimensions without coordinates: space\n",
+       "Data variables:\n",
+       "    *empty*\n",
+       "Indexes:\n",
+       "    points   GeometryIndex (crs=None)
    " + ], + "text/plain": [ + " Size: 80B\n", + "Dimensions: (space: 10)\n", + "Coordinates:\n", + " * points (space) object 80B POINT (171.61191601133118 -10.618898030367944...\n", + "Dimensions without coordinates: space\n", + "Data variables:\n", + " *empty*\n", + "Indexes:\n", + " points GeometryIndex (crs=None)" + ] + }, + "execution_count": 10, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "ds_no_crs = ds_wgs84.drop_vars(\"spatial_ref\")\n", + "\n", + "ds_no_crs" + ] + }, + { + "cell_type": "markdown", + "id": "fff6acec-9d29-443e-a29c-802833b4388d", + "metadata": {}, + "source": [ + "It is possible to combine datasets with and without a defined CRS. The resulting dataset will have the common CRS found among all datasets. " + ] + }, + { + "cell_type": "code", + "execution_count": 11, + "id": "d142186b-9f5a-4809-b775-66e35e7af175", + "metadata": {}, + "outputs": [ + { + "data": { + "text/html": [ + "
    \n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "
    <xarray.Dataset> Size: 88B\n",
+       "Dimensions:      (space: 10)\n",
+       "Coordinates:\n",
+       "  * points       (space) object 80B POINT (171.61191601133118 -10.61889803036...\n",
+       "  * spatial_ref  int64 8B 0\n",
+       "Dimensions without coordinates: space\n",
+       "Data variables:\n",
+       "    *empty*\n",
+       "Indexes:\n",
+       "    points       GeometryIndex (crs=None)\n",
+       "    spatial_ref  CRSIndex (crs=EPSG:4326)
    " + ], + "text/plain": [ + " Size: 88B\n", + "Dimensions: (space: 10)\n", + "Coordinates:\n", + " * points (space) object 80B POINT (171.61191601133118 -10.61889803036...\n", + " * spatial_ref int64 8B 0\n", + "Dimensions without coordinates: space\n", + "Data variables:\n", + " *empty*\n", + "Indexes:\n", + " points GeometryIndex (crs=None)\n", + " spatial_ref CRSIndex (crs=EPSG:4326)" + ] + }, + "execution_count": 11, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "ds_wgs84 + ds_no_crs" + ] } ], "metadata": {