issues

What is your issue?

I need to get from a DEM all flat areas and make a binary mask. I'm using xarray and xarray-spatial, below the code: ``` from xrspatial.aspect import aspect import xarray as xr

dem = xr.open_mfdataset( paths=dem_path, chunks={'x': 1000, 'y': 1000}, parallel=True, ) dem = dem.rename_vars({'band_data': 'dem'}) dem = dem.rename_dims({'band': 'elevation'})

dem_aspect = aspect(dem['dem'])

flat_dem = dem_aspect.where(dem_aspect > 0, other=0) flat_dem = flat_dem.where(flat_dem == 0, other=1) ```

At the end I need to save the mask as GeoTIFF: flat_dem.astype('float32').rio.to_raster(main_path.joinpath('flat_dem.tif'), driver="GTiff")

But when I try to save the mask I see the error below:

TypingError: Failed in nopython mode pipeline (step: nopython frontend) ('wrong tuple length for $46load_attr.20: ', 'expected 2, got 3') During: typing of exhaust iter at /home/max/.cache/pypoetry/virtualenvs/pygis_blog-iPs2SW5R-py3.10/lib/python3.10/site-packages/xrspatial/aspect.py (27)

File "../../../../../../../home/max/.cache/pypoetry/virtualenvs/pygis_blog-iPs2SW5R-py3.10/lib/python3.10/site-packages/xrspatial/aspect.py", line 27: def _run_numpy(data: np.ndarray): <source elided> out[:] = np.nan rows, cols = data.shape

``` TypingError Traceback (most recent call last) Cell In[4], line 17 14 flaat_dem = dem_aspect.where(dem_aspect > 0, other=0) 15 flaat_dem = flaat_dem.where(flaat_dem == 0, other=1) ---> 17 flaat_dem.astype('float32').rio.to_raster(main_path.joinpath('flat_dem.tif'), driver="GTiff")

File ~/.cache/pypoetry/virtualenvs/pygis_blog-iPs2SW5R-py3.10/lib/python3.10/site-packages/rioxarray/raster_array.py:1125, in RasterArray.to_raster(self, raster_path, driver, dtype, tags, windowed, recalc_transform, lock, compute, profile_kwargs) 1106 out_profile = { 1107 key: value 1108 for key, value in out_profile.items() (...) 1119 ) 1120 } 1121 rio_nodata = ( 1122 self.encoded_nodata if self.encoded_nodata is not None else self.nodata 1123 ) -> 1125 return RasterioWriter(raster_path=raster_path).to_raster( 1126 xarray_dataarray=self._obj, 1127 tags=tags, 1128 driver=driver, 1129 height=int(self.height), 1130 width=int(self.width), 1131 count=int(self.count), 1132 dtype=dtype, 1133 crs=self.crs, 1134 transform=self.transform(recalc=recalc_transform), 1135 gcps=self.get_gcps(), 1136 nodata=rio_nodata, 1137 windowed=windowed, 1138 lock=lock, 1139 compute=compute, 1140 out_profile, 1141 )

File ~/.cache/pypoetry/virtualenvs/pygis_blog-iPs2SW5R-py3.10/lib/python3.10/site-packages/rioxarray/raster_writer.py:286, in RasterioWriter.to_raster(self, xarray_dataarray, tags, windowed, lock, compute, **kwargs) 284 else: 285 out_data = xarray_dataarray --> 286 data = encode_cf_variable(out_data.variable).values.astype( 287 numpy_dtype 288 ) 289 if data.ndim == 2: 290 rds.write(data, 1, window=window)

File ~/.cache/pypoetry/virtualenvs/pygis_blog-iPs2SW5R-py3.10/lib/python3.10/site-packages/xarray/core/variable.py:614, in Variable.values(self) 611 @property 612 def values(self): 613 """The variable's data as a numpy.ndarray""" --> 614 return _as_array_or_item(self._data)

File ~/.cache/pypoetry/virtualenvs/pygis_blog-iPs2SW5R-py3.10/lib/python3.10/site-packages/xarray/core/variable.py:314, in _as_array_or_item(data) 300 def _as_array_or_item(data): 301 """Return the given values as a numpy array, or as an individual item if 302 it's a 0d datetime64 or timedelta64 array. 303 (...) 312 TODO: remove this (replace with np.asarray) once these issues are fixed 313 """ --> 314 data = np.asarray(data) 315 if data.ndim == 0: 316 if data.dtype.kind == "M":

File ~/.cache/pypoetry/virtualenvs/pygis_blog-iPs2SW5R-py3.10/lib/python3.10/site-packages/dask/threaded.py:89, in get(dsk, keys, cache, num_workers, pool, kwargs) 86 elif isinstance(pool, multiprocessing.pool.Pool): 87 pool = MultiprocessingPoolExecutor(pool) ---> 89 results = get_async( 90 pool.submit, 91 pool._max_workers, 92 dsk, 93 keys, 94 cache=cache, 95 get_id=_thread_get_id, 96 pack_exception=pack_exception, 97 kwargs, 98 ) 100 # Cleanup pools associated to dead threads 101 with pools_lock:

File ~/.cache/pypoetry/virtualenvs/pygis_blog-iPs2SW5R-py3.10/lib/python3.10/site-packages/dask/local.py:511, in get_async(submit, num_workers, dsk, result, cache, get_id, rerun_exceptions_locally, pack_exception, raise_exception, callbacks, dumps, loads, chunksize, **kwargs) 509 _execute_task(task, data) # Re-execute locally 510 else: --> 511 raise_exception(exc, tb) 512 res, worker_id = loads(res_info) 513 state["cache"][key] = res

File ~/.cache/pypoetry/virtualenvs/pygis_blog-iPs2SW5R-py3.10/lib/python3.10/site-packages/dask/local.py:319, in reraise(exc, tb) 317 if exc.traceback is not tb: 318 raise exc.with_traceback(tb) --> 319 raise exc

File ~/.cache/pypoetry/virtualenvs/pygis_blog-iPs2SW5R-py3.10/lib/python3.10/site-packages/dask/local.py:224, in execute_task(key, task_info, dumps, loads, get_id, pack_exception) 222 try: 223 task, data = loads(task_info) --> 224 result = _execute_task(task, data) 225 id = get_id() 226 result = dumps((result, id))

File ~/.cache/pypoetry/virtualenvs/pygis_blog-iPs2SW5R-py3.10/lib/python3.10/site-packages/numba/core/dispatcher.py:468, in _DispatcherBase._compile_for_args(self, args, *kws) 464 msg = (f"{str(e).rstrip()} \n\nThis error may have been caused " 465 f"by the following argument(s):\n{args_str}\n") 466 e.patch_message(msg) --> 468 error_rewrite(e, 'typing') 469 except errors.UnsupportedError as e: 470 # Something unsupported is present in the user code, add help info 471 error_rewrite(e, 'unsupported_error')

File ~/.cache/pypoetry/virtualenvs/pygis_blog-iPs2SW5R-py3.10/lib/python3.10/site-packages/numba/core/dispatcher.py:409, in _DispatcherBase._compile_for_args.<locals>.error_rewrite(e, issue_type) 407 raise e 408 else: --> 409 raise e.with_traceback(None) ``` Here there is the DEM.

I'm using xarray 2023.8.0 with xarray-spatial 0.3.7

NB: I need to use aspect because in the next step I need to classify the exposition of my data avoiding the flat areas.

As I mentioned here, when I add a raster inside a geocube are produced a lot of nodata.

After a test @snowman2 has indicates that this is a possible bug of xarray libraries.

I was able to reproduce the issue you demonstrated. I think it has to do with decimal precision of the coordinates when adding the DataArray to the Dataset. That being said, I checked that the x coordinate and y coordinates were the exact same for raster_in and out_grid. This is likely an issue with xarray and not geocube.

For a workaround, this got it working without the gaps for me:

```python import geopandas import rioxarray from geocube.api.core import make_geocube

vector_in = geopandas.read_file("vector_data.shp") raster_in = rioxarray.open_rasterio("raster_data.tif", masked=True).sel(band=1).drop("band")

out_grid = make_geocube( vector_data=vector_in, measurements=["id"], like=raster_in, )

out_grid["process_value"] = (raster_in.dims, raster_in.values, raster_in.attrs, raster_in.encoding) ```