Clipping an xarray dataset with rio clip function produces empty xarray

[Simply-Adi]

I have an xarray dataset ds containing a variable 'pop'.
xarray.Dataset Dimensions: time: 1latitude: 5751longitude: 13197 Coordinates: time (time) datetime64[ns] 2000-12-31 latitude (latitude) float64 18.91 18.92 18.93 ... 71.38 71.39 longitude (longitude) float64 -179.1 -179.1 ... 179.8 179.8 spatial_ref () int32 0 Data variables: pop (time, latitude, longitude) float64 nan nan nan nan ... nan nan nan nan Indexes: (3) Attributes: (0)

I want to clip ds to the extent of California using geopandas dataframe 'aoi' with projection WGS84/EPSG4326. For this I tried:

`
from shapely.geometry import mapping

import rioxarray

import geopandas as gpd
shp_data = gpd.read_file('USA.shp') #USA.shp
aoi = shp_data[shp_data['state_code']=='CA'].reset_index()
ds = ds.rio.set_spatial_dims(x_dim="longitude", y_dim="latitude", inplace=True);
ds = ds.rio.write_crs("epsg:4326", inplace=True)
ds_clipped = ds.rio.clip( geometries = aoi.geometry.values, crs= aoi.crs, drop=True,invert=False);
#have tried using aoi.geometry.apply(mapping) as well
`
However, ds_clipped is empty.

To check if I am logically wrong, i did the same using regionmask as follows:

`
import regionmask

mask = regionmask.mask_3D_geopandas(aoi,ds.longitude,ds.latitude)

ds_clipped = ds.where(mask)
`
Surprisingly, this works. I am utterly confused about why this is happening. I would be thankful if somebody could correct me.

The links for ds and shp_data (whole USA) is here (https://drive.google.com/drive/folders/1NnvlKuL_pIZ4bPiLdcn1mvz5SbbKsINs?usp=share_link)

[snowman2]

I loaded it into QGIS and the California shapefile does not overlap the raster data where California should be:

[Simply-Adi]

Hi, thank you for your diagnosis. I cannot figure out what is causing the problem. I had matched the lat-lon ranges and they come under WGS84 system. I kind of suspect that I screwed up in creating the raster/xarray dataset (ds). I created it from a csv data in XYZ format. I have csv files for each year. So, to create ds, I did df.to_xarray() after setting year, lat, lob as indices. I also did sorting by indices. Do you have any suggestions for systrmatically approaching pandas df to xarray conversion to avoid coordinate system mismatch?

…

On Sat, 25 Feb, 2023, 8:03 am Alan D. Snow, ***@***.***> wrote: I loaded it into QGIS and the California shapefile does not overlap the raster data where California should be: [image: image] <https://user-images.githubusercontent.com/8699967/221331736-a2559f0b-ce15-461c-b4cc-a69c3b25922a.png> — Reply to this email directly, view it on GitHub <#638 (comment)>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/ASQFP5ZBMHPP7CRAXXB527DWZFVODANCNFSM6AAAAAAVEGRYOI> . You are receiving this because you authored the thread.Message ID: ***@***.***>

[snowman2]

You may be interested in geocube: https://github.com/corteva/geocube/