0

I have a bunch of daily netcdf files that cover the same spatial extent. I can open just one file with xr.open_dataset:

ds = xr.open_datset(filename, decode_cf=True)

This results in

<xarray.Dataset>
Dimensions:        (time: 1, bnds: 2, y: 292, x: 278, soil: 4)
Coordinates:
  * time           (time) datetime64[ns] 2021-01-01T12:00:00
    longitude      (y, x) float32 ...
    latitude       (y, x) float32 ...
Dimensions without coordinates: bnds, y, x, soil
Data variables:
    time_bnds      (time, bnds) datetime64[ns] ...
    arr1           (time, soil, y, x) float32 ...
    arr2           (time, y, x) float32 ...
    arr3           (time, y, x) float32 ...
    arr4           (time, y, x) float32 ...
    arr5           (time, y, x) float32 ...
    arr6           (time, y, x) float32 ...
Attributes:
    CDI:          Climate Data Interface version 1.9.8 (https://mpimet.mpg.de...
    Conventions:  CF-1.6

Now, I would like to select data by location (and time), but

ds.sel(longitude=0.4, latitude=9.4, method="nearest")

fails with a KeyError: 'no index found for coordinate longitude' error, which I guess is linked to the Dimensions without coordinates bit.

Ideally, I would like to open all my netcdf files (using xr.open_mfdataset) and be able to extract a time series for a given location in arr5 above.

The only way I found to do this is calculating distances to points:

lon = ds.coords["longitude"].values * 1
lat = ds.coords["latitude"].values * 1
dist = ((lon - sel_lon) ** 2 + (lat + sel_lat) ** 2)
locs = np.where(dist==np.min(dist))

I am sure that there has to be a better way of doing this!

Jose
  • 1,987
  • 2
  • 21
  • 27
  • I am not sure there is a better way using `xarray`. See discussion https://stackoverflow.com/questions/57145557/select-data-by-latitude-and-longitude – Robert Davy Feb 24 '22 at 00:33

0 Answers0