{ "cells": [ { "cell_type": "markdown", "id": "4f4d7519-6942-4ae4-a2f8-27aa167b6050", "metadata": {}, "source": [ "## Xarray: propagate bounds coordinate with an IntervalIndex in DataArray\n", "\n", "It should also work with any other coordinate associated with an Xarray index that shares at least one dimension with dataarray's dimensions." ] }, { "cell_type": "code", "execution_count": 1, "id": "779c11f2-2ec2-4c16-98e6-4b601bdf6c61", "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "import pandas as pd\n", "import xarray as xr\n", "\n", "xr.set_options(display_expand_indexes=True);" ] }, { "cell_type": "markdown", "id": "6e5dafab-e774-4380-9f2c-c37a9f3dd262", "metadata": {}, "source": [ "Example of an Xarray IntervalIndex that can be associated to a dimension coordinate and its CF bounds coordinate companion." ] }, { "cell_type": "code", "execution_count": 2, "id": "d5c647f4-b12f-4b22-9351-c201191a452c", "metadata": {}, "outputs": [], "source": [ "from collections.abc import Hashable\n", "\n", "from xarray.core.indexes import Index, PandasIndex\n", "from xarray import Variable\n", "\n", "\n", "class IntervalIndex(Index):\n", " # adapted from https://github.com/dcherian/xindexes/blob/main/interval-array.ipynb\n", "\n", " def __init__(self, index: PandasIndex, bounds_name: Hashable, bounds_dim: str):\n", " assert isinstance(index.index, pd.IntervalIndex)\n", " self._index = index\n", " self._bounds_name = bounds_name\n", " self._bounds_dim = bounds_dim\n", "\n", " @classmethod\n", " def from_variables(cls, variables, options):\n", " assert len(variables) == 2\n", "\n", " for k, v in variables.items():\n", " if v.ndim == 2:\n", " bounds_name, bounds = k, v\n", " elif v.ndim == 1:\n", " dim, _ = k, v\n", "\n", " bounds = bounds.transpose(..., dim)\n", " left, right = bounds.data.tolist()\n", " index = PandasIndex(pd.IntervalIndex.from_arrays(left, right), dim)\n", " bounds_dim = (set(bounds.dims) - set(dim)).pop()\n", " \n", " return cls(index, bounds_name, bounds_dim)\n", "\n", " @classmethod\n", " def concat(self, indexes, dim, positions=None):\n", " new_index = self._index.concat([idx._index for idx in indexes], dim, positions=positions)\n", "\n", " if indexes:\n", " bounds_name0 = indexes[0]._bounds_name\n", " bounds_dim0 = indexes[0]._bounds_dim\n", " if any(idx._bounds_name != bounds_name0 or idx._bounds_dim != bounds_dim0 for idx in indexes):\n", " raise ValueError(\n", " f\"Cannot concatenate along dimension {dim!r} indexes with different \"\n", " \"boundary coordinate or dimension names\"\n", " )\n", "\n", " return type(self)(new_index, self._bounds_name, self._bounds_dim)\n", " \n", " def create_variables(self, variables=None):\n", " empty_var = Variable((), 0)\n", " bounds_attrs = variables.get(self._bounds_name, empty_var).attrs\n", " mid_attrs = variables.get(self._index.dim, empty_var).attrs\n", "\n", " bounds_var = Variable(\n", " dims=(self._bounds_dim, self._index.dim),\n", " data=np.stack([self._index.index.left, self._index.index.right], axis=0),\n", " attrs=bounds_attrs,\n", " )\n", " mid_var = Variable(\n", " dims=(self._index.dim,),\n", " data=self._index.index.mid,\n", " attrs=mid_attrs,\n", " )\n", "\n", " return {self._index.dim: mid_var, self._bounds_name: bounds_var}\n", "\n", " def sel(self, labels, **kwargs):\n", " return self._index.sel(labels, **kwargs)\n", "\n", " def isel(self, indexers):\n", " new_index = self._index.isel(indexers)\n", " if new_index is not None:\n", " return type(self)(new_index, self._bounds_name, self._bounds_dim)\n", " else:\n", " return None\n", "\n", " def roll(self, shifts):\n", " new_index = self._index.roll(shifts)\n", " return type(self)(new_index, self._bounds_name, self._bounds_dim)\n", "\n", " def rename(self, name_dict, dims_dict):\n", " new_index = self._index.rename(name_dict, dims_dict)\n", "\n", " bounds_name = name_dict.get(self._bounds_name, self._bounds_name)\n", " bounds_dim = dims_dict.get(self._bounds_dim, self._bounds_dim)\n", " \n", " return type(self)(new_index, bounds_name, bounds_dim)\n", "\n", " def __repr__(self):\n", " string = f\"{self._index!r}\"\n", " return string\n" ] }, { "cell_type": "markdown", "id": "a4ff39b9-b36d-444f-add9-cfd77cea8b78", "metadata": {}, "source": [ "Create an example dataset with an IntervalIndex" ] }, { "cell_type": "code", "execution_count": 3, "id": "1c0ec9e9-234d-4a19-9a9c-23932f14d8a8", "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "
<xarray.Dataset> Size: 128B\n",
       "Dimensions:   (x: 4, bnds: 2)\n",
       "Coordinates:\n",
       "  * x         (x) float64 32B 1.0 2.0 3.0 4.0\n",
       "  * x_bounds  (bnds, x) float64 64B 0.5 1.5 2.5 3.5 1.5 2.5 3.5 4.5\n",
       "Dimensions without coordinates: bnds\n",
       "Data variables:\n",
       "    foo       (x) int64 32B 1 2 3 4\n",
       "Indexes:\n",
       "  ┌ x         IntervalIndex\n",
       "  └ x_bounds
" ], "text/plain": [ " Size: 128B\n", "Dimensions: (x: 4, bnds: 2)\n", "Coordinates:\n", " * x (x) float64 32B 1.0 2.0 3.0 4.0\n", " * x_bounds (bnds, x) float64 64B 0.5 1.5 2.5 3.5 1.5 2.5 3.5 4.5\n", "Dimensions without coordinates: bnds\n", "Data variables:\n", " foo (x) int64 32B 1 2 3 4\n", "Indexes:\n", " ┌ x IntervalIndex\n", " └ x_bounds" ] }, "execution_count": 3, "metadata": {}, "output_type": "execute_result" } ], "source": [ "left = np.arange(0.5, 3.6, 1)\n", "right = np.arange(1.5, 4.6, 1)\n", "bounds = np.stack([left, right])\n", "\n", "ds = xr.Dataset(\n", " {\"foo\": (\"x\", [1, 2, 3, 4])},\n", " coords={\"x\": [1, 2, 3, 4], \"x_bounds\": ((\"bnds\", \"x\"), bounds)},\n", ")\n", "\n", "newds = ds.drop_indexes(\"x\").set_xindex((\"x\", \"x_bounds\"), IntervalIndex)\n", "newds" ] }, { "cell_type": "markdown", "id": "60f3c29c-7b63-4666-9372-6acf62015de9", "metadata": {}, "source": [ "## Extract a DataArray from a Dataset with a bounds coordinate " ] }, { "cell_type": "code", "execution_count": 4, "id": "dcd26e55-17d9-417f-a3e7-e2ec4d31cb7b", "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "
<xarray.DataArray 'foo' (x: 4)> Size: 32B\n",
       "array([1, 2, 3, 4])\n",
       "Coordinates:\n",
       "  * x         (x) float64 32B 1.0 2.0 3.0 4.0\n",
       "  * x_bounds  (bnds, x) float64 64B 0.5 1.5 2.5 3.5 1.5 2.5 3.5 4.5\n",
       "Indexes:\n",
       "  ┌ x         IntervalIndex\n",
       "  └ x_bounds
" ], "text/plain": [ " Size: 32B\n", "array([1, 2, 3, 4])\n", "Coordinates:\n", " * x (x) float64 32B 1.0 2.0 3.0 4.0\n", " * x_bounds (bnds, x) float64 64B 0.5 1.5 2.5 3.5 1.5 2.5 3.5 4.5\n", "Indexes:\n", " ┌ x IntervalIndex\n", " └ x_bounds" ] }, "execution_count": 4, "metadata": {}, "output_type": "execute_result" } ], "source": [ "newds.foo" ] }, { "cell_type": "code", "execution_count": 5, "id": "c4b5cd0e-f053-4172-a9b4-2eaf480429fd", "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "
<xarray.DataArray 'x' (x: 4)> Size: 32B\n",
       "array([1., 2., 3., 4.])\n",
       "Coordinates:\n",
       "  * x         (x) float64 32B 1.0 2.0 3.0 4.0\n",
       "  * x_bounds  (bnds, x) float64 64B 0.5 1.5 2.5 3.5 1.5 2.5 3.5 4.5\n",
       "Indexes:\n",
       "  ┌ x         IntervalIndex\n",
       "  └ x_bounds
" ], "text/plain": [ " Size: 32B\n", "array([1., 2., 3., 4.])\n", "Coordinates:\n", " * x (x) float64 32B 1.0 2.0 3.0 4.0\n", " * x_bounds (bnds, x) float64 64B 0.5 1.5 2.5 3.5 1.5 2.5 3.5 4.5\n", "Indexes:\n", " ┌ x IntervalIndex\n", " └ x_bounds" ] }, "execution_count": 5, "metadata": {}, "output_type": "execute_result" } ], "source": [ "newds.x" ] }, { "cell_type": "markdown", "id": "4f84d4f7-517b-4757-b413-898ebd212a06", "metadata": {}, "source": [ "The bounds coordinate is not propagated if it is not associated with an index! The IntervalIndex is also dropped." ] }, { "cell_type": "code", "execution_count": 6, "id": "3607c2b9-72c3-49eb-8f29-f1563af0c57d", "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "
<xarray.DataArray 'foo' (x: 4)> Size: 32B\n",
       "array([1, 2, 3, 4])\n",
       "Coordinates:\n",
       "    x        (x) float64 32B 1.0 2.0 3.0 4.0
" ], "text/plain": [ " Size: 32B\n", "array([1, 2, 3, 4])\n", "Coordinates:\n", " x (x) float64 32B 1.0 2.0 3.0 4.0" ] }, "execution_count": 6, "metadata": {}, "output_type": "execute_result" } ], "source": [ "newds.drop_indexes([\"x\", \"x_bounds\"]).foo" ] }, { "cell_type": "markdown", "id": "3a8ecfd0-8db4-4f51-9a31-dbfb9b33f4e2", "metadata": {}, "source": [ "## Create a new DataArray with a bounds coordinate" ] }, { "cell_type": "code", "execution_count": 7, "id": "2a74463d-84a5-4320-9f8d-3f1214f8e04a", "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "
<xarray.DataArray 'foo' (x: 4)> Size: 32B\n",
       "array([1, 2, 3, 4])\n",
       "Coordinates:\n",
       "  * x         (x) float64 32B 1.0 2.0 3.0 4.0\n",
       "  * x_bounds  (bnds, x) float64 64B 0.5 1.5 2.5 3.5 1.5 2.5 3.5 4.5\n",
       "Indexes:\n",
       "  ┌ x         IntervalIndex\n",
       "  └ x_bounds
" ], "text/plain": [ " Size: 32B\n", "array([1, 2, 3, 4])\n", "Coordinates:\n", " * x (x) float64 32B 1.0 2.0 3.0 4.0\n", " * x_bounds (bnds, x) float64 64B 0.5 1.5 2.5 3.5 1.5 2.5 3.5 4.5\n", "Indexes:\n", " ┌ x IntervalIndex\n", " └ x_bounds" ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# note: the IntervalIndex associated with the \"x\" and \"x_bounds\" coordinates\n", "# is passed to the DataArray constructor below via the `coords` argument\n", "\n", "da = xr.DataArray(newds.foo, newds.coords, dims=\"x\")\n", "da" ] }, { "cell_type": "markdown", "id": "6cfb4343-47c8-43d8-ab31-5311514b718b", "metadata": {}, "source": [ "The bounds coordinate cannot be included if it is not associated with an index!" ] }, { "cell_type": "code", "execution_count": 8, "id": "e921a311-79df-4535-a462-df8e12d94c15", "metadata": {}, "outputs": [ { "ename": "ValueError", "evalue": "coordinate x_bounds has dimensions ('bnds', 'x'), but these are not a subset of the DataArray dimensions ('x',)", "output_type": "error", "traceback": [ "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m", "\u001b[0;31mValueError\u001b[0m Traceback (most recent call last)", "Cell \u001b[0;32mIn[8], line 2\u001b[0m\n\u001b[1;32m 1\u001b[0m \u001b[38;5;66;03m# note: convert `newds.coords` to a dict removes all indexes\u001b[39;00m\n\u001b[0;32m----> 2\u001b[0m \u001b[43mxr\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mDataArray\u001b[49m\u001b[43m(\u001b[49m\u001b[43mnewds\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mfoo\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[38;5;28;43mdict\u001b[39;49m\u001b[43m(\u001b[49m\u001b[43mnewds\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mcoords\u001b[49m\u001b[43m)\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mdims\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[38;5;124;43m\"\u001b[39;49m\u001b[38;5;124;43mx\u001b[39;49m\u001b[38;5;124;43m\"\u001b[39;49m\u001b[43m)\u001b[49m\n", "File \u001b[0;32m~/Git/github/benbovy/xarray/xarray/core/dataarray.py:522\u001b[0m, in \u001b[0;36mDataArray.__init__\u001b[0;34m(self, data, coords, dims, name, attrs, indexes, fastpath)\u001b[0m\n\u001b[1;32m 520\u001b[0m \u001b[38;5;28;01mif\u001b[39;00m \u001b[38;5;129;01mnot\u001b[39;00m \u001b[38;5;28misinstance\u001b[39m(coords, Coordinates):\n\u001b[1;32m 521\u001b[0m coords \u001b[38;5;241m=\u001b[39m create_coords_with_default_indexes(coords)\n\u001b[0;32m--> 522\u001b[0m \u001b[43mcheck_dataarray_coords\u001b[49m\u001b[43m(\u001b[49m\u001b[43mdata\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mshape\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mcoords\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mdims\u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m 523\u001b[0m indexes \u001b[38;5;241m=\u001b[39m \u001b[38;5;28mdict\u001b[39m(coords\u001b[38;5;241m.\u001b[39mxindexes)\n\u001b[1;32m 524\u001b[0m coords \u001b[38;5;241m=\u001b[39m {k: v\u001b[38;5;241m.\u001b[39mcopy() \u001b[38;5;28;01mfor\u001b[39;00m k, v \u001b[38;5;129;01min\u001b[39;00m coords\u001b[38;5;241m.\u001b[39mvariables\u001b[38;5;241m.\u001b[39mitems()}\n", "File \u001b[0;32m~/Git/github/benbovy/xarray/xarray/core/dataarray.py:167\u001b[0m, in \u001b[0;36mcheck_dataarray_coords\u001b[0;34m(shape, coords, dims)\u001b[0m\n\u001b[1;32m 165\u001b[0m raise_error \u001b[38;5;241m=\u001b[39m \u001b[38;5;28many\u001b[39m(d \u001b[38;5;129;01mnot\u001b[39;00m \u001b[38;5;129;01min\u001b[39;00m dims \u001b[38;5;28;01mfor\u001b[39;00m d \u001b[38;5;129;01min\u001b[39;00m var\u001b[38;5;241m.\u001b[39mdims)\n\u001b[1;32m 166\u001b[0m \u001b[38;5;28;01mif\u001b[39;00m raise_error:\n\u001b[0;32m--> 167\u001b[0m \u001b[38;5;28;01mraise\u001b[39;00m \u001b[38;5;167;01mValueError\u001b[39;00m(\n\u001b[1;32m 168\u001b[0m \u001b[38;5;124mf\u001b[39m\u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mcoordinate \u001b[39m\u001b[38;5;132;01m{\u001b[39;00mname\u001b[38;5;132;01m}\u001b[39;00m\u001b[38;5;124m has dimensions \u001b[39m\u001b[38;5;132;01m{\u001b[39;00mvar\u001b[38;5;241m.\u001b[39mdims\u001b[38;5;132;01m}\u001b[39;00m\u001b[38;5;124m, but these \u001b[39m\u001b[38;5;124m\"\u001b[39m\n\u001b[1;32m 169\u001b[0m \u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mare not a subset of the DataArray \u001b[39m\u001b[38;5;124m\"\u001b[39m\n\u001b[1;32m 170\u001b[0m \u001b[38;5;124mf\u001b[39m\u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mdimensions \u001b[39m\u001b[38;5;132;01m{\u001b[39;00mdims\u001b[38;5;132;01m}\u001b[39;00m\u001b[38;5;124m\"\u001b[39m\n\u001b[1;32m 171\u001b[0m )\n\u001b[1;32m 173\u001b[0m \u001b[38;5;66;03m# check dimension sizes\u001b[39;00m\n\u001b[1;32m 174\u001b[0m sizes \u001b[38;5;241m=\u001b[39m \u001b[38;5;28mdict\u001b[39m(\u001b[38;5;28mzip\u001b[39m(dims, shape, strict\u001b[38;5;241m=\u001b[39m\u001b[38;5;28;01mTrue\u001b[39;00m))\n", "\u001b[0;31mValueError\u001b[0m: coordinate x_bounds has dimensions ('bnds', 'x'), but these are not a subset of the DataArray dimensions ('x',)" ] } ], "source": [ "# note: convert `newds.coords` to a dict removes all indexes\n", "xr.DataArray(newds.foo, dict(newds.coords), dims=\"x\")" ] }, { "cell_type": "markdown", "id": "4da74fe4-5f17-42ac-8b38-5080e459ffe2", "metadata": {}, "source": [ "## Assign bounds coordinate to an existing DataArray" ] }, { "cell_type": "code", "execution_count": 9, "id": "d58ede39-6ddb-41fc-9669-e4830433a94e", "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "
<xarray.DataArray (x: 4)> Size: 32B\n",
       "array([0, 1, 2, 3])\n",
       "Coordinates:\n",
       "  * x         (x) float64 32B 1.0 2.0 3.0 4.0\n",
       "  * x_bounds  (bnds, x) float64 64B 0.5 1.5 2.5 3.5 1.5 2.5 3.5 4.5\n",
       "Indexes:\n",
       "  ┌ x         IntervalIndex\n",
       "  └ x_bounds
" ], "text/plain": [ " Size: 32B\n", "array([0, 1, 2, 3])\n", "Coordinates:\n", " * x (x) float64 32B 1.0 2.0 3.0 4.0\n", " * x_bounds (bnds, x) float64 64B 0.5 1.5 2.5 3.5 1.5 2.5 3.5 4.5\n", "Indexes:\n", " ┌ x IntervalIndex\n", " └ x_bounds" ] }, "execution_count": 9, "metadata": {}, "output_type": "execute_result" } ], "source": [ "xr.DataArray(range(4), dims=\"x\").assign_coords(newds.coords)" ] }, { "cell_type": "markdown", "id": "ef0b7e16-4919-45c9-b03c-e194cb1864b3", "metadata": {}, "source": [ "The bounds coordinate cannot be assigned if it is not associated with an index!" ] }, { "cell_type": "code", "execution_count": 10, "id": "e8b7f2bb-f955-4a8c-b62a-89d848db3697", "metadata": {}, "outputs": [ { "ename": "ValueError", "evalue": "coordinate x_bounds has dimensions ('bnds', 'x'), but these are not a subset of the DataArray dimensions ('x',)", "output_type": "error", "traceback": [ "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m", "\u001b[0;31mValueError\u001b[0m Traceback (most recent call last)", "Cell \u001b[0;32mIn[10], line 1\u001b[0m\n\u001b[0;32m----> 1\u001b[0m \u001b[43mxr\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mDataArray\u001b[49m\u001b[43m(\u001b[49m\u001b[38;5;28;43mrange\u001b[39;49m\u001b[43m(\u001b[49m\u001b[38;5;241;43m4\u001b[39;49m\u001b[43m)\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mdims\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[38;5;124;43m\"\u001b[39;49m\u001b[38;5;124;43mx\u001b[39;49m\u001b[38;5;124;43m\"\u001b[39;49m\u001b[43m)\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43massign_coords\u001b[49m\u001b[43m(\u001b[49m\u001b[38;5;28;43mdict\u001b[39;49m\u001b[43m(\u001b[49m\u001b[43mnewds\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mcoords\u001b[49m\u001b[43m)\u001b[49m\u001b[43m)\u001b[49m\n", "File \u001b[0;32m~/Git/github/benbovy/xarray/xarray/core/common.py:664\u001b[0m, in \u001b[0;36mDataWithCoords.assign_coords\u001b[0;34m(self, coords, **coords_kwargs)\u001b[0m\n\u001b[1;32m 661\u001b[0m \u001b[38;5;28;01melse\u001b[39;00m:\n\u001b[1;32m 662\u001b[0m results \u001b[38;5;241m=\u001b[39m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39m_calc_assign_results(coords_combined)\n\u001b[0;32m--> 664\u001b[0m \u001b[43mdata\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mcoords\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mupdate\u001b[49m\u001b[43m(\u001b[49m\u001b[43mresults\u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m 665\u001b[0m \u001b[38;5;28;01mreturn\u001b[39;00m data\n", "File \u001b[0;32m~/Git/github/benbovy/xarray/xarray/core/coordinates.py:603\u001b[0m, in \u001b[0;36mCoordinates.update\u001b[0;34m(self, other)\u001b[0m\n\u001b[1;32m 597\u001b[0m \u001b[38;5;66;03m# special case for PandasMultiIndex: updating only its dimension coordinate\u001b[39;00m\n\u001b[1;32m 598\u001b[0m \u001b[38;5;66;03m# is still allowed but depreciated.\u001b[39;00m\n\u001b[1;32m 599\u001b[0m \u001b[38;5;66;03m# It is the only case where we need to actually drop coordinates here (multi-index levels)\u001b[39;00m\n\u001b[1;32m 600\u001b[0m \u001b[38;5;66;03m# TODO: remove when removing PandasMultiIndex's dimension coordinate.\u001b[39;00m\n\u001b[1;32m 601\u001b[0m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39m_drop_coords(\u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39m_names \u001b[38;5;241m-\u001b[39m coords_to_align\u001b[38;5;241m.\u001b[39m_names)\n\u001b[0;32m--> 603\u001b[0m \u001b[38;5;28;43mself\u001b[39;49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43m_update_coords\u001b[49m\u001b[43m(\u001b[49m\u001b[43mcoords\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mindexes\u001b[49m\u001b[43m)\u001b[49m\n", "File \u001b[0;32m~/Git/github/benbovy/xarray/xarray/core/coordinates.py:978\u001b[0m, in \u001b[0;36mDataArrayCoordinates._update_coords\u001b[0;34m(self, coords, indexes)\u001b[0m\n\u001b[1;32m 974\u001b[0m obj_dims \u001b[38;5;241m=\u001b[39m \u001b[38;5;28mset\u001b[39m(\u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mdims)\n\u001b[1;32m 976\u001b[0m \u001b[38;5;28;01mif\u001b[39;00m coords_dims \u001b[38;5;241m>\u001b[39m obj_dims:\n\u001b[1;32m 977\u001b[0m \u001b[38;5;66;03m# need more checks\u001b[39;00m\n\u001b[0;32m--> 978\u001b[0m \u001b[43mcheck_dataarray_coords\u001b[49m\u001b[43m(\u001b[49m\u001b[38;5;28;43mself\u001b[39;49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43m_data\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mshape\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mCoordinates\u001b[49m\u001b[43m(\u001b[49m\u001b[43mcoords\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mindexes\u001b[49m\u001b[43m)\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[38;5;28;43mself\u001b[39;49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mdims\u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m 980\u001b[0m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39m_data\u001b[38;5;241m.\u001b[39m_coords \u001b[38;5;241m=\u001b[39m coords\n\u001b[1;32m 981\u001b[0m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39m_data\u001b[38;5;241m.\u001b[39m_indexes \u001b[38;5;241m=\u001b[39m indexes\n", "File \u001b[0;32m~/Git/github/benbovy/xarray/xarray/core/dataarray.py:167\u001b[0m, in \u001b[0;36mcheck_dataarray_coords\u001b[0;34m(shape, coords, dims)\u001b[0m\n\u001b[1;32m 165\u001b[0m raise_error \u001b[38;5;241m=\u001b[39m \u001b[38;5;28many\u001b[39m(d \u001b[38;5;129;01mnot\u001b[39;00m \u001b[38;5;129;01min\u001b[39;00m dims \u001b[38;5;28;01mfor\u001b[39;00m d \u001b[38;5;129;01min\u001b[39;00m var\u001b[38;5;241m.\u001b[39mdims)\n\u001b[1;32m 166\u001b[0m \u001b[38;5;28;01mif\u001b[39;00m raise_error:\n\u001b[0;32m--> 167\u001b[0m \u001b[38;5;28;01mraise\u001b[39;00m \u001b[38;5;167;01mValueError\u001b[39;00m(\n\u001b[1;32m 168\u001b[0m \u001b[38;5;124mf\u001b[39m\u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mcoordinate \u001b[39m\u001b[38;5;132;01m{\u001b[39;00mname\u001b[38;5;132;01m}\u001b[39;00m\u001b[38;5;124m has dimensions \u001b[39m\u001b[38;5;132;01m{\u001b[39;00mvar\u001b[38;5;241m.\u001b[39mdims\u001b[38;5;132;01m}\u001b[39;00m\u001b[38;5;124m, but these \u001b[39m\u001b[38;5;124m\"\u001b[39m\n\u001b[1;32m 169\u001b[0m \u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mare not a subset of the DataArray \u001b[39m\u001b[38;5;124m\"\u001b[39m\n\u001b[1;32m 170\u001b[0m \u001b[38;5;124mf\u001b[39m\u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mdimensions \u001b[39m\u001b[38;5;132;01m{\u001b[39;00mdims\u001b[38;5;132;01m}\u001b[39;00m\u001b[38;5;124m\"\u001b[39m\n\u001b[1;32m 171\u001b[0m )\n\u001b[1;32m 173\u001b[0m \u001b[38;5;66;03m# check dimension sizes\u001b[39;00m\n\u001b[1;32m 174\u001b[0m sizes \u001b[38;5;241m=\u001b[39m \u001b[38;5;28mdict\u001b[39m(\u001b[38;5;28mzip\u001b[39m(dims, shape, strict\u001b[38;5;241m=\u001b[39m\u001b[38;5;28;01mTrue\u001b[39;00m))\n", "\u001b[0;31mValueError\u001b[0m: coordinate x_bounds has dimensions ('bnds', 'x'), but these are not a subset of the DataArray dimensions ('x',)" ] } ], "source": [ "xr.DataArray(range(4), dims=\"x\").assign_coords(dict(newds.coords))" ] }, { "cell_type": "markdown", "id": "1638a449-3bc8-4f4e-84bf-6f3443285e7c", "metadata": {}, "source": [ "## Ensure common operations on DataArray propagate the bounds coordinate" ] }, { "cell_type": "code", "execution_count": 11, "id": "225aa06e-2641-463f-a46d-4f4243e6c5b9", "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "
<xarray.DataArray 'x_bounds' (bnds: 2, x: 4)> Size: 64B\n",
       "array([[0.5, 1.5, 2.5, 3.5],\n",
       "       [1.5, 2.5, 3.5, 4.5]])\n",
       "Coordinates:\n",
       "  * x         (x) float64 32B 1.0 2.0 3.0 4.0\n",
       "  * x_bounds  (bnds, x) float64 64B 0.5 1.5 2.5 3.5 1.5 2.5 3.5 4.5\n",
       "Dimensions without coordinates: bnds\n",
       "Indexes:\n",
       "  ┌ x         IntervalIndex\n",
       "  └ x_bounds
" ], "text/plain": [ " Size: 64B\n", "array([[0.5, 1.5, 2.5, 3.5],\n", " [1.5, 2.5, 3.5, 4.5]])\n", "Coordinates:\n", " * x (x) float64 32B 1.0 2.0 3.0 4.0\n", " * x_bounds (bnds, x) float64 64B 0.5 1.5 2.5 3.5 1.5 2.5 3.5 4.5\n", "Dimensions without coordinates: bnds\n", "Indexes:\n", " ┌ x IntervalIndex\n", " └ x_bounds" ] }, "execution_count": 11, "metadata": {}, "output_type": "execute_result" } ], "source": [ "da.x_bounds" ] }, { "cell_type": "code", "execution_count": 12, "id": "5c77bf21-cea8-4ae3-bec1-46e7edf63831", "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "
<xarray.Dataset> Size: 128B\n",
       "Dimensions:   (x: 4, bnds: 2)\n",
       "Coordinates:\n",
       "  * x         (x) float64 32B 1.0 2.0 3.0 4.0\n",
       "  * x_bounds  (bnds, x) float64 64B 0.5 1.5 2.5 3.5 1.5 2.5 3.5 4.5\n",
       "Dimensions without coordinates: bnds\n",
       "Data variables:\n",
       "    foo       (x) int64 32B 1 2 3 4\n",
       "Indexes:\n",
       "  ┌ x         IntervalIndex\n",
       "  └ x_bounds
" ], "text/plain": [ " Size: 128B\n", "Dimensions: (x: 4, bnds: 2)\n", "Coordinates:\n", " * x (x) float64 32B 1.0 2.0 3.0 4.0\n", " * x_bounds (bnds, x) float64 64B 0.5 1.5 2.5 3.5 1.5 2.5 3.5 4.5\n", "Dimensions without coordinates: bnds\n", "Data variables:\n", " foo (x) int64 32B 1 2 3 4\n", "Indexes:\n", " ┌ x IntervalIndex\n", " └ x_bounds" ] }, "execution_count": 12, "metadata": {}, "output_type": "execute_result" } ], "source": [ "da.to_dataset(name=\"foo\")" ] }, { "cell_type": "code", "execution_count": 13, "id": "bef7258d-bf2b-4eaa-bcb3-de19c74db305", "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "
<xarray.DataArray 'foo' (x: 2)> Size: 16B\n",
       "array([1, 3])\n",
       "Coordinates:\n",
       "  * x         (x) float64 16B 1.0 3.0\n",
       "  * x_bounds  (bnds, x) float64 32B 0.5 2.5 1.5 3.5\n",
       "Indexes:\n",
       "  ┌ x         IntervalIndex\n",
       "  └ x_bounds
" ], "text/plain": [ " Size: 16B\n", "array([1, 3])\n", "Coordinates:\n", " * x (x) float64 16B 1.0 3.0\n", " * x_bounds (bnds, x) float64 32B 0.5 2.5 1.5 3.5\n", "Indexes:\n", " ┌ x IntervalIndex\n", " └ x_bounds" ] }, "execution_count": 13, "metadata": {}, "output_type": "execute_result" } ], "source": [ "da.sel(x=[0.8, 3.25])" ] }, { "cell_type": "code", "execution_count": 14, "id": "df1b21b4-6e8f-4a8d-aca1-b6b5a1fb78e7", "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "
<xarray.DataArray 'foo' (x: 4)> Size: 32B\n",
       "array([2, 4, 6, 8])\n",
       "Coordinates:\n",
       "  * x         (x) float64 32B 1.0 2.0 3.0 4.0\n",
       "  * x_bounds  (bnds, x) float64 64B 0.5 1.5 2.5 3.5 1.5 2.5 3.5 4.5\n",
       "Indexes:\n",
       "  ┌ x         IntervalIndex\n",
       "  └ x_bounds
" ], "text/plain": [ " Size: 32B\n", "array([2, 4, 6, 8])\n", "Coordinates:\n", " * x (x) float64 32B 1.0 2.0 3.0 4.0\n", " * x_bounds (bnds, x) float64 64B 0.5 1.5 2.5 3.5 1.5 2.5 3.5 4.5\n", "Indexes:\n", " ┌ x IntervalIndex\n", " └ x_bounds" ] }, "execution_count": 14, "metadata": {}, "output_type": "execute_result" } ], "source": [ "da + da" ] }, { "cell_type": "code", "execution_count": 15, "id": "2942182d-551b-462f-8cd1-30ca1f0bad10", "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "
<xarray.DataArray 'foo' (x: 4)> Size: 32B\n",
       "array([4, 1, 2, 3])\n",
       "Coordinates:\n",
       "  * x         (x) float64 32B 4.0 1.0 2.0 3.0\n",
       "  * x_bounds  (bnds, x) float64 64B 3.5 0.5 1.5 2.5 4.5 1.5 2.5 3.5\n",
       "Indexes:\n",
       "  ┌ x         IntervalIndex\n",
       "  └ x_bounds
" ], "text/plain": [ " Size: 32B\n", "array([4, 1, 2, 3])\n", "Coordinates:\n", " * x (x) float64 32B 4.0 1.0 2.0 3.0\n", " * x_bounds (bnds, x) float64 64B 3.5 0.5 1.5 2.5 4.5 1.5 2.5 3.5\n", "Indexes:\n", " ┌ x IntervalIndex\n", " └ x_bounds" ] }, "execution_count": 15, "metadata": {}, "output_type": "execute_result" } ], "source": [ "da.roll(x=1, roll_coords=True)" ] }, { "cell_type": "code", "execution_count": 16, "id": "95c47a8d-1d53-4ff2-8672-c877d2b081b7", "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "
<xarray.DataArray 'foo' (y: 4)> Size: 32B\n",
       "array([1, 2, 3, 4])\n",
       "Coordinates:\n",
       "  * y         (y) float64 32B 1.0 2.0 3.0 4.0\n",
       "  * y_bounds  (bnds, y) float64 64B 0.5 1.5 2.5 3.5 1.5 2.5 3.5 4.5\n",
       "Indexes:\n",
       "  ┌ y         IntervalIndex\n",
       "  └ y_bounds
" ], "text/plain": [ " Size: 32B\n", "array([1, 2, 3, 4])\n", "Coordinates:\n", " * y (y) float64 32B 1.0 2.0 3.0 4.0\n", " * y_bounds (bnds, y) float64 64B 0.5 1.5 2.5 3.5 1.5 2.5 3.5 4.5\n", "Indexes:\n", " ┌ y IntervalIndex\n", " └ y_bounds" ] }, "execution_count": 16, "metadata": {}, "output_type": "execute_result" } ], "source": [ "da.rename(x=\"y\", x_bounds=\"y_bounds\")" ] }, { "cell_type": "code", "execution_count": 17, "id": "261187cf-c521-45d1-b032-0ad74ebeb646", "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "
<xarray.Dataset> Size: 128B\n",
       "Dimensions:   (x: 4, bnds: 2)\n",
       "Coordinates:\n",
       "  * x         (x) float64 32B 1.0 2.0 3.0 4.0\n",
       "  * x_bounds  (bnds, x) float64 64B 0.5 1.5 2.5 3.5 1.5 2.5 3.5 4.5\n",
       "Dimensions without coordinates: bnds\n",
       "Data variables:\n",
       "    foo       (x) int64 32B 1 2 3 4\n",
       "Indexes:\n",
       "  ┌ x         IntervalIndex\n",
       "  └ x_bounds
" ], "text/plain": [ " Size: 128B\n", "Dimensions: (x: 4, bnds: 2)\n", "Coordinates:\n", " * x (x) float64 32B 1.0 2.0 3.0 4.0\n", " * x_bounds (bnds, x) float64 64B 0.5 1.5 2.5 3.5 1.5 2.5 3.5 4.5\n", "Dimensions without coordinates: bnds\n", "Data variables:\n", " foo (x) int64 32B 1 2 3 4\n", "Indexes:\n", " ┌ x IntervalIndex\n", " └ x_bounds" ] }, "execution_count": 17, "metadata": {}, "output_type": "execute_result" } ], "source": [ "xr.merge([da, da])" ] }, { "cell_type": "markdown", "id": "7c65887b-43ec-4623-9180-5a4bc9e20bbb", "metadata": {}, "source": [ "`IntervalIndex.concat()` is implemented and should normally work but Xarray current alignment rules prevents it to be applied." ] }, { "cell_type": "code", "execution_count": 18, "id": "4a2da111-3bab-47f2-8851-0e2242dc68bb", "metadata": {}, "outputs": [ { "ename": "ValueError", "evalue": "cannot exclude dimension(s) x from alignment because these are used by an index together with non-excluded dimensions bnds", "output_type": "error", "traceback": [ "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m", "\u001b[0;31mValueError\u001b[0m Traceback (most recent call last)", "Cell \u001b[0;32mIn[18], line 1\u001b[0m\n\u001b[0;32m----> 1\u001b[0m \u001b[43mxr\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mconcat\u001b[49m\u001b[43m(\u001b[49m\u001b[43m[\u001b[49m\u001b[43mda\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mda\u001b[49m\u001b[43m]\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[38;5;124;43m\"\u001b[39;49m\u001b[38;5;124;43mx\u001b[39;49m\u001b[38;5;124;43m\"\u001b[39;49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mcoords\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[43m[\u001b[49m\u001b[38;5;124;43m\"\u001b[39;49m\u001b[38;5;124;43mx\u001b[39;49m\u001b[38;5;124;43m\"\u001b[39;49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[38;5;124;43m\"\u001b[39;49m\u001b[38;5;124;43mx_bounds\u001b[39;49m\u001b[38;5;124;43m\"\u001b[39;49m\u001b[43m]\u001b[49m\u001b[43m)\u001b[49m\n", "File \u001b[0;32m~/Git/github/benbovy/xarray/xarray/core/concat.py:264\u001b[0m, in \u001b[0;36mconcat\u001b[0;34m(objs, dim, data_vars, coords, compat, positions, fill_value, join, combine_attrs, create_index_for_new_dim)\u001b[0m\n\u001b[1;32m 259\u001b[0m \u001b[38;5;28;01mraise\u001b[39;00m \u001b[38;5;167;01mValueError\u001b[39;00m(\n\u001b[1;32m 260\u001b[0m \u001b[38;5;124mf\u001b[39m\u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mcompat=\u001b[39m\u001b[38;5;132;01m{\u001b[39;00mcompat\u001b[38;5;132;01m!r}\u001b[39;00m\u001b[38;5;124m invalid: must be \u001b[39m\u001b[38;5;124m'\u001b[39m\u001b[38;5;124mbroadcast_equals\u001b[39m\u001b[38;5;124m'\u001b[39m\u001b[38;5;124m, \u001b[39m\u001b[38;5;124m'\u001b[39m\u001b[38;5;124mequals\u001b[39m\u001b[38;5;124m'\u001b[39m\u001b[38;5;124m, \u001b[39m\u001b[38;5;124m'\u001b[39m\u001b[38;5;124midentical\u001b[39m\u001b[38;5;124m'\u001b[39m\u001b[38;5;124m, \u001b[39m\u001b[38;5;124m'\u001b[39m\u001b[38;5;124mno_conflicts\u001b[39m\u001b[38;5;124m'\u001b[39m\u001b[38;5;124m or \u001b[39m\u001b[38;5;124m'\u001b[39m\u001b[38;5;124moverride\u001b[39m\u001b[38;5;124m'\u001b[39m\u001b[38;5;124m\"\u001b[39m\n\u001b[1;32m 261\u001b[0m )\n\u001b[1;32m 263\u001b[0m \u001b[38;5;28;01mif\u001b[39;00m \u001b[38;5;28misinstance\u001b[39m(first_obj, DataArray):\n\u001b[0;32m--> 264\u001b[0m \u001b[38;5;28;01mreturn\u001b[39;00m \u001b[43m_dataarray_concat\u001b[49m\u001b[43m(\u001b[49m\n\u001b[1;32m 265\u001b[0m \u001b[43m \u001b[49m\u001b[43mobjs\u001b[49m\u001b[43m,\u001b[49m\n\u001b[1;32m 266\u001b[0m \u001b[43m \u001b[49m\u001b[43mdim\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[43mdim\u001b[49m\u001b[43m,\u001b[49m\n\u001b[1;32m 267\u001b[0m \u001b[43m \u001b[49m\u001b[43mdata_vars\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[43mdata_vars\u001b[49m\u001b[43m,\u001b[49m\n\u001b[1;32m 268\u001b[0m \u001b[43m \u001b[49m\u001b[43mcoords\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[43mcoords\u001b[49m\u001b[43m,\u001b[49m\n\u001b[1;32m 269\u001b[0m \u001b[43m \u001b[49m\u001b[43mcompat\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[43mcompat\u001b[49m\u001b[43m,\u001b[49m\n\u001b[1;32m 270\u001b[0m \u001b[43m \u001b[49m\u001b[43mpositions\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[43mpositions\u001b[49m\u001b[43m,\u001b[49m\n\u001b[1;32m 271\u001b[0m \u001b[43m \u001b[49m\u001b[43mfill_value\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[43mfill_value\u001b[49m\u001b[43m,\u001b[49m\n\u001b[1;32m 272\u001b[0m \u001b[43m \u001b[49m\u001b[43mjoin\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[43mjoin\u001b[49m\u001b[43m,\u001b[49m\n\u001b[1;32m 273\u001b[0m \u001b[43m \u001b[49m\u001b[43mcombine_attrs\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[43mcombine_attrs\u001b[49m\u001b[43m,\u001b[49m\n\u001b[1;32m 274\u001b[0m \u001b[43m \u001b[49m\u001b[43mcreate_index_for_new_dim\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[43mcreate_index_for_new_dim\u001b[49m\u001b[43m,\u001b[49m\n\u001b[1;32m 275\u001b[0m \u001b[43m \u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m 276\u001b[0m \u001b[38;5;28;01melif\u001b[39;00m \u001b[38;5;28misinstance\u001b[39m(first_obj, Dataset):\n\u001b[1;32m 277\u001b[0m \u001b[38;5;28;01mreturn\u001b[39;00m _dataset_concat(\n\u001b[1;32m 278\u001b[0m objs,\n\u001b[1;32m 279\u001b[0m dim\u001b[38;5;241m=\u001b[39mdim,\n\u001b[0;32m (...)\u001b[0m\n\u001b[1;32m 287\u001b[0m create_index_for_new_dim\u001b[38;5;241m=\u001b[39mcreate_index_for_new_dim,\n\u001b[1;32m 288\u001b[0m )\n", "File \u001b[0;32m~/Git/github/benbovy/xarray/xarray/core/concat.py:755\u001b[0m, in \u001b[0;36m_dataarray_concat\u001b[0;34m(arrays, dim, data_vars, coords, compat, positions, fill_value, join, combine_attrs, create_index_for_new_dim)\u001b[0m\n\u001b[1;32m 752\u001b[0m arr \u001b[38;5;241m=\u001b[39m arr\u001b[38;5;241m.\u001b[39mrename(name)\n\u001b[1;32m 753\u001b[0m datasets\u001b[38;5;241m.\u001b[39mappend(arr\u001b[38;5;241m.\u001b[39m_to_temp_dataset())\n\u001b[0;32m--> 755\u001b[0m ds \u001b[38;5;241m=\u001b[39m \u001b[43m_dataset_concat\u001b[49m\u001b[43m(\u001b[49m\n\u001b[1;32m 756\u001b[0m \u001b[43m \u001b[49m\u001b[43mdatasets\u001b[49m\u001b[43m,\u001b[49m\n\u001b[1;32m 757\u001b[0m \u001b[43m \u001b[49m\u001b[43mdim\u001b[49m\u001b[43m,\u001b[49m\n\u001b[1;32m 758\u001b[0m \u001b[43m \u001b[49m\u001b[43mdata_vars\u001b[49m\u001b[43m,\u001b[49m\n\u001b[1;32m 759\u001b[0m \u001b[43m \u001b[49m\u001b[43mcoords\u001b[49m\u001b[43m,\u001b[49m\n\u001b[1;32m 760\u001b[0m \u001b[43m \u001b[49m\u001b[43mcompat\u001b[49m\u001b[43m,\u001b[49m\n\u001b[1;32m 761\u001b[0m \u001b[43m \u001b[49m\u001b[43mpositions\u001b[49m\u001b[43m,\u001b[49m\n\u001b[1;32m 762\u001b[0m \u001b[43m \u001b[49m\u001b[43mfill_value\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[43mfill_value\u001b[49m\u001b[43m,\u001b[49m\n\u001b[1;32m 763\u001b[0m \u001b[43m \u001b[49m\u001b[43mjoin\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[43mjoin\u001b[49m\u001b[43m,\u001b[49m\n\u001b[1;32m 764\u001b[0m \u001b[43m \u001b[49m\u001b[43mcombine_attrs\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[43mcombine_attrs\u001b[49m\u001b[43m,\u001b[49m\n\u001b[1;32m 765\u001b[0m \u001b[43m \u001b[49m\u001b[43mcreate_index_for_new_dim\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[43mcreate_index_for_new_dim\u001b[49m\u001b[43m,\u001b[49m\n\u001b[1;32m 766\u001b[0m \u001b[43m\u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m 768\u001b[0m merged_attrs \u001b[38;5;241m=\u001b[39m merge_attrs([da\u001b[38;5;241m.\u001b[39mattrs \u001b[38;5;28;01mfor\u001b[39;00m da \u001b[38;5;129;01min\u001b[39;00m arrays], combine_attrs)\n\u001b[1;32m 770\u001b[0m result \u001b[38;5;241m=\u001b[39m arrays[\u001b[38;5;241m0\u001b[39m]\u001b[38;5;241m.\u001b[39m_from_temp_dataset(ds, name)\n", "File \u001b[0;32m~/Git/github/benbovy/xarray/xarray/core/concat.py:516\u001b[0m, in \u001b[0;36m_dataset_concat\u001b[0;34m(datasets, dim, data_vars, coords, compat, positions, fill_value, join, combine_attrs, create_index_for_new_dim)\u001b[0m\n\u001b[1;32m 513\u001b[0m \u001b[38;5;66;03m# Make sure we're working on a copy (we'll be loading variables)\u001b[39;00m\n\u001b[1;32m 514\u001b[0m datasets \u001b[38;5;241m=\u001b[39m [ds\u001b[38;5;241m.\u001b[39mcopy() \u001b[38;5;28;01mfor\u001b[39;00m ds \u001b[38;5;129;01min\u001b[39;00m datasets]\n\u001b[1;32m 515\u001b[0m datasets \u001b[38;5;241m=\u001b[39m \u001b[38;5;28mlist\u001b[39m(\n\u001b[0;32m--> 516\u001b[0m \u001b[43malign\u001b[49m\u001b[43m(\u001b[49m\n\u001b[1;32m 517\u001b[0m \u001b[43m \u001b[49m\u001b[38;5;241;43m*\u001b[39;49m\u001b[43mdatasets\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mjoin\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[43mjoin\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mcopy\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[38;5;28;43;01mFalse\u001b[39;49;00m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mexclude\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[43m[\u001b[49m\u001b[43mdim_name\u001b[49m\u001b[43m]\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mfill_value\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[43mfill_value\u001b[49m\n\u001b[1;32m 518\u001b[0m \u001b[43m \u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m 519\u001b[0m )\n\u001b[1;32m 521\u001b[0m dim_coords, dims_sizes, coord_names, data_names, vars_order \u001b[38;5;241m=\u001b[39m _parse_datasets(\n\u001b[1;32m 522\u001b[0m datasets\n\u001b[1;32m 523\u001b[0m )\n\u001b[1;32m 524\u001b[0m dim_names \u001b[38;5;241m=\u001b[39m \u001b[38;5;28mset\u001b[39m(dim_coords)\n", "File \u001b[0;32m~/Git/github/benbovy/xarray/xarray/core/alignment.py:883\u001b[0m, in \u001b[0;36malign\u001b[0;34m(join, copy, indexes, exclude, fill_value, *objects)\u001b[0m\n\u001b[1;32m 687\u001b[0m \u001b[38;5;250m\u001b[39m\u001b[38;5;124;03m\"\"\"\u001b[39;00m\n\u001b[1;32m 688\u001b[0m \u001b[38;5;124;03mGiven any number of Dataset and/or DataArray objects, returns new\u001b[39;00m\n\u001b[1;32m 689\u001b[0m \u001b[38;5;124;03mobjects with aligned indexes and dimension sizes.\u001b[39;00m\n\u001b[0;32m (...)\u001b[0m\n\u001b[1;32m 873\u001b[0m \n\u001b[1;32m 874\u001b[0m \u001b[38;5;124;03m\"\"\"\u001b[39;00m\n\u001b[1;32m 875\u001b[0m aligner \u001b[38;5;241m=\u001b[39m Aligner(\n\u001b[1;32m 876\u001b[0m objects,\n\u001b[1;32m 877\u001b[0m join\u001b[38;5;241m=\u001b[39mjoin,\n\u001b[0;32m (...)\u001b[0m\n\u001b[1;32m 881\u001b[0m fill_value\u001b[38;5;241m=\u001b[39mfill_value,\n\u001b[1;32m 882\u001b[0m )\n\u001b[0;32m--> 883\u001b[0m \u001b[43maligner\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43malign\u001b[49m\u001b[43m(\u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m 884\u001b[0m \u001b[38;5;28;01mreturn\u001b[39;00m aligner\u001b[38;5;241m.\u001b[39mresults\n", "File \u001b[0;32m~/Git/github/benbovy/xarray/xarray/core/alignment.py:572\u001b[0m, in \u001b[0;36mAligner.align\u001b[0;34m(self)\u001b[0m\n\u001b[1;32m 569\u001b[0m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mresults \u001b[38;5;241m=\u001b[39m (obj\u001b[38;5;241m.\u001b[39mcopy(deep\u001b[38;5;241m=\u001b[39m\u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mcopy),)\n\u001b[1;32m 570\u001b[0m \u001b[38;5;28;01mreturn\u001b[39;00m\n\u001b[0;32m--> 572\u001b[0m \u001b[38;5;28;43mself\u001b[39;49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mfind_matching_indexes\u001b[49m\u001b[43m(\u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m 573\u001b[0m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mfind_matching_unindexed_dims()\n\u001b[1;32m 574\u001b[0m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39massert_no_index_conflict()\n", "File \u001b[0;32m~/Git/github/benbovy/xarray/xarray/core/alignment.py:254\u001b[0m, in \u001b[0;36mAligner.find_matching_indexes\u001b[0;34m(self)\u001b[0m\n\u001b[1;32m 251\u001b[0m objects_matching_indexes \u001b[38;5;241m=\u001b[39m []\n\u001b[1;32m 253\u001b[0m \u001b[38;5;28;01mfor\u001b[39;00m obj \u001b[38;5;129;01min\u001b[39;00m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mobjects:\n\u001b[0;32m--> 254\u001b[0m obj_indexes, obj_index_vars \u001b[38;5;241m=\u001b[39m \u001b[38;5;28;43mself\u001b[39;49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43m_normalize_indexes\u001b[49m\u001b[43m(\u001b[49m\u001b[43mobj\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mxindexes\u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m 255\u001b[0m objects_matching_indexes\u001b[38;5;241m.\u001b[39mappend(obj_indexes)\n\u001b[1;32m 256\u001b[0m \u001b[38;5;28;01mfor\u001b[39;00m key, idx \u001b[38;5;129;01min\u001b[39;00m obj_indexes\u001b[38;5;241m.\u001b[39mitems():\n", "File \u001b[0;32m~/Git/github/benbovy/xarray/xarray/core/alignment.py:230\u001b[0m, in \u001b[0;36mAligner._normalize_indexes\u001b[0;34m(self, indexes)\u001b[0m\n\u001b[1;32m 228\u001b[0m excl_dims_str \u001b[38;5;241m=\u001b[39m \u001b[38;5;124m\"\u001b[39m\u001b[38;5;124m, \u001b[39m\u001b[38;5;124m\"\u001b[39m\u001b[38;5;241m.\u001b[39mjoin(\u001b[38;5;28mstr\u001b[39m(d) \u001b[38;5;28;01mfor\u001b[39;00m d \u001b[38;5;129;01min\u001b[39;00m exclude_dims)\n\u001b[1;32m 229\u001b[0m incl_dims_str \u001b[38;5;241m=\u001b[39m \u001b[38;5;124m\"\u001b[39m\u001b[38;5;124m, \u001b[39m\u001b[38;5;124m\"\u001b[39m\u001b[38;5;241m.\u001b[39mjoin(\u001b[38;5;28mstr\u001b[39m(d) \u001b[38;5;28;01mfor\u001b[39;00m d \u001b[38;5;129;01min\u001b[39;00m all_dims \u001b[38;5;241m-\u001b[39m exclude_dims)\n\u001b[0;32m--> 230\u001b[0m \u001b[38;5;28;01mraise\u001b[39;00m \u001b[38;5;167;01mValueError\u001b[39;00m(\n\u001b[1;32m 231\u001b[0m \u001b[38;5;124mf\u001b[39m\u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mcannot exclude dimension(s) \u001b[39m\u001b[38;5;132;01m{\u001b[39;00mexcl_dims_str\u001b[38;5;132;01m}\u001b[39;00m\u001b[38;5;124m from alignment because \u001b[39m\u001b[38;5;124m\"\u001b[39m\n\u001b[1;32m 232\u001b[0m \u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mthese are used by an index together with non-excluded dimensions \u001b[39m\u001b[38;5;124m\"\u001b[39m\n\u001b[1;32m 233\u001b[0m \u001b[38;5;124mf\u001b[39m\u001b[38;5;124m\"\u001b[39m\u001b[38;5;132;01m{\u001b[39;00mincl_dims_str\u001b[38;5;132;01m}\u001b[39;00m\u001b[38;5;124m\"\u001b[39m\n\u001b[1;32m 234\u001b[0m )\n\u001b[1;32m 236\u001b[0m key \u001b[38;5;241m=\u001b[39m (\u001b[38;5;28mtuple\u001b[39m(coord_names_and_dims), \u001b[38;5;28mtype\u001b[39m(idx))\n\u001b[1;32m 237\u001b[0m normalized_indexes[key] \u001b[38;5;241m=\u001b[39m idx\n", "\u001b[0;31mValueError\u001b[0m: cannot exclude dimension(s) x from alignment because these are used by an index together with non-excluded dimensions bnds" ] } ], "source": [ "xr.concat([da, da], \"x\", coords=[\"x\", \"x_bounds\"])" ] }, { "cell_type": "code", "execution_count": null, "id": "15538f65-4e78-4b28-9668-5d264101421b", "metadata": {}, "outputs": [], "source": [] } ], "metadata": { "kernelspec": { "display_name": "Python (Pixi)", "language": "python", "name": "pixi-kernel-python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.12.8" } }, "nbformat": 4, "nbformat_minor": 5 }