⚙️ How-To Gridded & climate data Raster & rasterio Gridded & Climate Data · step 4 of 20 Updated

How to Select a Time Range and Location from an xarray Dataset

Problem statement

.sel() is the first thing anyone does with a gridded dataset, and it fails in ways that do not look like failures. On the NCEP/NCAR Reanalysis 1 monthly surface air temperature file — latitude stored from 90 down to −90, longitude from 0 to 357.5 — four ordinary selections measured:

air.sel(lat=slice(30, 60))                        -> 0 latitudes, no error
air.sel(lat=40.71, lon=-74.01, method="nearest")  -> a cell off the coast of Spain, not New York
air.sel(lat=[40.71, 51.51], lon=[286.0, 0.0])     -> a 2 x 2 grid, not two points
air.sel(time="2024-01-20", method="nearest")      -> February 2024

None of them raises. The empty slice produces a dataset whose mean is NaN; the New York selection returns a temperature 7.61 °C too warm for January 2024; the list selection returns combinations of the coordinates that nobody asked for. Each is correct behaviour for a label-based index — and each is avoidable once you know which assumption it makes.

Quick answer 

import xarray as xr

air = xr.open_dataset("air.mon.mean.nc")["air"]

# a region: slice latitude in the order the file stores it (here 90 -> -90)
india_summer = air.sel(time=slice("2023-06", "2023-08"),
                       lat=slice(35, 5), lon=slice(68, 98))

# a point: match the grid's longitude convention, and refuse a distant match
new_york = air.sel(lat=40.71, lon=-74.01 % 360, method="nearest", tolerance=2.5)

print(dict(india_summer.sizes), float(new_york.lon))
{'time': 3, 'lat': 13, 'lon': 12} 285.0

Three rules cover most selections: read the coordinate order before slicing, convert point longitudes to the grid's convention, and give method="nearest" a tolerance so a wrong match raises instead of returning a number.

Four silent xarray selection failures, each paired with the change that prevents it.
All four return data. Only a check on the result, or a tolerance, turns them into errors.

Step-by-step solution

1. Read the coordinate order and range first 

for name in ("time", "lat", "lon"):
    index = air.indexes[name]
    order = "increasing" if index.is_monotonic_increasing else "decreasing"
    print(f"{name:5} {index[0]} -> {index[-1]}  ({order}, {len(index)} values)")
time  1948-01-01 00:00:00 -> 2026-02-01 00:00:00  (increasing, 938 values)
lat   90.0 -> -90.0  (decreasing, 73 values)
lon   0.0 -> 357.5  (increasing, 144 values)

This one print decides how every later slice has to be written. Descending latitude is common in reanalyses and in files converted from image formats; 0–360 longitude is common in global models.

2. Select time with partial date strings

xarray understands partial dates, and both ends of a label slice are inclusive:

year = air.sel(time="2023")                              # 12 months
summer = air.sel(time=slice("2020-06", "2020-08"))       # June, July, August 2020
print(year.sizes["time"], summer.time.dt.strftime("%Y-%m-%d").values)
12 ['2020-06-01' '2020-07-01' '2020-08-01']

A slice that runs past the end of the data is not an error either. slice("2025-01", "2030-12") returned 14 months — the 12 of 2025 and the two of 2026 that the file contains — so check the first and last timestamps whenever the period matters.

3. Slice latitude in the stored order

A label slice follows the order of the index. On a descending axis, slice(30, 60) asks for labels from 30 down to 60, and there are none:

print(air.sel(lat=slice(30, 60)).sizes["lat"])    # 0
print(air.sel(lat=slice(60, 30)).sizes["lat"])    # 13

To write code that works on either orientation, sort once, or build the slice from the data:

lat = air.indexes["lat"]
band = slice(60, 30) if lat.is_monotonic_decreasing else slice(30, 60)
mid_latitudes = air.sel(lat=band)

4. Put point longitudes in the grid's convention

A negative longitude on a 0–360 grid is simply a number smaller than every coordinate, so method="nearest" picks the first column, 0°. Measured for January 2024:

City Requested lon Cell picked Value Correct cell Value Error
New York −74.01 0.0° 9.89 °C 285.0° 2.28 °C +7.61 °C
Mexico City −99.13 0.0° 20.10 °C 260.0° 14.63 °C +5.47 °C
Los Angeles −118.24 0.0° 12.55 °C 242.5° 8.37 °C +4.18 °C
Rio de Janeiro −43.17 0.0° 22.54 °C 317.5° 23.98 °C −1.45 °C

The fix is lon % 360 for a 0–360 grid, or converting the grid itself — the trade-offs are in longitude conventions explained.

5. Give method="nearest" a tolerance

tolerance is the maximum distance, in coordinate units, that a nearest match may be from the request. With it, the New York mistake raises instead of returning a number:

try:
    air.sel(lat=40.71, lon=-74.01, method="nearest", tolerance=2.5)
except KeyError as error:
    print(error)
"not all values found in index 'lon'"

Set it to one grid spacing. The same parameter protects time selection: monthly timestamps sit on the first of the month, so sel(time="2024-01-20", method="nearest") returned 1 February, which is 12 days away, rather than 1 January, which is 19. An exact sel(time="2024-01-20") raises KeyError and suggests method='nearest' — the suggestion is the trap.

6. Select several points with DataArray indexers

Passing lists selects every combination of the values — orthogonal indexing:

lats = [40.71, 51.51, 35.68]
lons = [285.99, 359.87, 139.69]

grid = air.sel(lat=lats, lon=lons, method="nearest")
points = air.sel(lat=xr.DataArray(lats, dims="city"),
                 lon=xr.DataArray(lons, dims="city"), method="nearest")
print(dict(grid.sizes), dict(points.sizes))
{'time': 938, 'lat': 3, 'lon': 3} {'time': 938, 'city': 3}

Nine series for three cities, six of them for places like "New York's latitude at Tokyo's longitude". Wrapping the indexers in DataArrays that share a new dimension selects the three pairs instead.

7. Choose between the nearest cell and interpolation

sel(method="nearest") returns a cell value; interp() blends the surrounding cells. On a 2.5° grid the difference is not small. For July 2023:

            nearest   interp   difference
New York     24.63     23.80     -0.83
Tokyo        24.96     24.69     -0.27
Denver       19.03     19.29     +0.26
London       15.25       nan

London is NaN because its longitude, 359.87° on this grid, lies beyond the last column at 357.5° and interp does not wrap around the globe. A cell value is what the data contains; an interpolated value is an estimate. Use nearest when you will compare with other gridded values, and interpolation only on a grid whose longitude seam is not in the way.

Bar chart of January 2024 temperature errors when negative longitudes are matched on a 0 to 360 grid: +7.61 °C for New York, +5.47 for Mexico City, +4.18 for Los Angeles and −1.45 for Rio de Janeiro.
Every one of these picked longitude 0°, so the error is just the difference between each city and a point on the Greenwich meridian.

Code examples

Example 1 — a box selection that works on any orientation or convention 

import numpy as np
import pandas as pd
import xarray as xr


def select_box(da, south, north, west, east, lat="lat", lon="lon"):
    """Select a latitude/longitude box whatever the storage order or longitude convention.

    A box with west > east crosses the grid's seam (for example 150 to -150).
    """
    ascending = da.indexes[lat].is_monotonic_increasing
    da = da.sel({lat: slice(south, north) if ascending else slice(north, south)})

    if float(da[lon].max()) > 180:                      # 0-360 grid
        west, east = west % 360, east % 360
    else:                                               # -180-180 grid
        west, east = (west + 180) % 360 - 180, (east + 180) % 360 - 180

    if west <= east:
        return da.sel({lon: slice(west, east)})
    return da.where((da[lon] >= west) | (da[lon] <= east), drop=True)
pacific = select_box(air.sel(time="2023"), -30, 30, 150, -150)
print(dict(pacific.sizes), round(float(pacific.mean()), 2))
{'time': 12, 'lat': 25, 'lon': 25} 25.58

Example 2 — point extraction that wraps longitude and refuses distant matches 

def extract_points(da, sites, max_distance=None, lat="lat", lon="lon"):
    """Nearest grid cell for each named site, along a new 'site' dimension.

    Longitude distance is measured around the globe, so -0.13 matches 0.0 on a
    0-360 grid rather than 357.5. Raises if any site is further than
    max_distance degrees (default: one grid spacing) from its cell.
    """
    names = list(sites)
    plat = np.array([sites[n][0] for n in names], dtype=float)
    plon = np.array([sites[n][1] for n in names], dtype=float)
    glat = da[lat].values.astype(float)
    glon = da[lon].values.astype(float)

    ilat = np.abs(glat[None, :] - plat[:, None]).argmin(axis=1)
    dlon = np.abs((glon[None, :] - plon[:, None] + 180) % 360 - 180)
    ilon = dlon.argmin(axis=1)

    limit = max_distance or max(np.abs(np.diff(glat)).max(), np.abs(np.diff(glon)).max())
    far = (np.abs(glat[ilat] - plat) > limit) | (dlon[np.arange(len(names)), ilon] > limit)
    if far.any():
        raise ValueError(f"no cell within {limit} degrees of {np.array(names)[far].tolist()}")

    picked = da.isel({lat: xr.DataArray(ilat, dims="site"),
                      lon: xr.DataArray(ilon, dims="site")})
    return picked.assign_coords(site=names)
cities = {"New York": (40.71, -74.01), "London": (51.51, -0.13), "Tokyo": (35.68, 139.69)}
july = extract_points(air.sel(time="2023-07-01"), cities)
for site, cell_lon, value in zip(july.site.values, july.lon.values, july.values):
    print(f"{site:9} lon {cell_lon:6.1f}  {value:5.2f} °C")

It works unchanged on a grid already converted to −180–180, because the distance calculation never assumes a convention.

Example 3 — a period selection that fails loudly when the data do not cover it 

def select_period(da, start, end, time="time"):
    """Slice start..end (partial dates allowed) and raise if the data fall short."""
    out = da.sel({time: slice(start, end)})
    if out.sizes[time] == 0:
        raise ValueError(f"no timesteps between {start} and {end}")
    first = pd.Timestamp(out[time].values[0])
    last = pd.Timestamp(out[time].values[-1])
    if first > pd.Period(start).end_time or last < pd.Period(end).start_time:
        raise ValueError(f"asked for {start} to {end}; data cover "
                         f"{first:%Y-%m-%d} to {last:%Y-%m-%d}")
    return out
print(select_period(air, "2020-06", "2020-08").sizes["time"])
select_period(air, "2025-01", "2030-12")
3
ValueError: asked for 2025-01 to 2030-12; data cover 2025-01-01 to 2026-02-01

Explanation

Why a slice follows the stored order

.sel() with a slice is label-based: it finds the position of the start label and the position of the stop label in the index and returns everything between them, in storage order. On a decreasing latitude index, the label 30 comes after the label 60, so a slice from 30 to 60 spans nothing.

Returning an empty result is deliberate — an empty selection is a legitimate answer for a range that holds no data. It becomes a problem only because the next step, usually a mean, silently turns the empty array into NaN.

Why nearest does not know that longitude wraps

To the index, longitude is just a sorted list of numbers. The distance from −74.01 to 0.0 is 74.01; to 285.0 it is 359.01. Nearest picks 0.0, correctly by arithmetic and absurdly by geography.

The seam causes the same problem in the other direction. London at −0.13° converted to 359.87° is 0.13° from the column at 0.0° around the globe, but the index measures 359.87 from 0.0 and 2.37 from 357.5, so sel(method="nearest") picked the column at 357.5° instead. A distance measured modulo 360, as in Example 2, is the only selection that is right on both sides of the seam.

Why a list selects a grid and a DataArray selects points

xarray keeps dimensions independent unless told otherwise. A list for lat and a list for lon are two separate selections on two separate dimensions, and the result is their outer product — three latitudes by three longitudes.

When both indexers are DataArrays with the same dimension name, xarray treats them as paired coordinates and walks along that shared dimension, the vectorised behaviour NumPy calls fancy indexing. The dimension name you choose becomes a real dimension of the result, which is why naming it city or site is worth the extra words.

Why nearest in time picked the following month

Monthly means are stamped with a single instant, the first of the month here, even though they describe the whole month. "Nearest" compares instants, not periods, so any day after the 16th of a 31-day month is nearer to the next month's stamp.

For monthly data, select with a partial string — sel(time="2024-01") — which matches on the period rather than the instant. Use method="nearest" in time only for data whose timestamps are the instants they describe, such as 6-hourly analyses.

Two panels contrasting list indexers, which select a three by three grid of cells, with DataArray indexers, which select three paired points.
The shape of the result is the check: a site dimension means pairs, separate lat and lon dimensions mean a grid.

Edge cases or notes

  • An unsorted coordinate breaks slicing silently. After reassigning longitudes without sorting, sel(lon=slice(-10, 30)) returned 0 columns; call sortby after any coordinate change.
  • Scalar selection drops the dimension. sel(lat=50) returns data without a lat dimension but keeps lat as a scalar coordinate; pass a one-element list to keep the dimension.
  • isel is positional and ignores order. air.isel(time=-1, lat=0, lon=0) and the matching sel gave the same value; isel avoids label questions entirely when you know the positions.
  • Chunk layout decides what a selection costs. This file stores one month per compressed chunk, so a single-cell time series decompresses all 938 of them.
  • Tolerance is in coordinate units. For time it is a timedelta such as "15D"; for latitude and longitude it is degrees.
  • Label slices are inclusive at both ends. Unlike Python's positional slices, slice("2020-06", "2020-08") includes August.
  • Model calendars change the rules for strings. On a noleap cftime axis, partial strings work, but a date such as 29 February does not exist.
  • Irregular regions need a mask, not a slice. Use .where() with a condition, or clip to a polygon.

FAQ

Why does xarray sel return an empty result for a latitude slice?

Because the latitude coordinate is stored in decreasing order and the slice was written increasing. On this file slice(30, 60) returned 0 latitudes and slice(60, 30) returned 13.

How do I select the nearest grid point to a longitude such as −74?

Convert it to the grid's convention first — −74.01 % 360 is 285.99 — and pass a tolerance. Without either, the NCEP grid matched longitude 0 and returned a value 7.61 °C too warm.

How do I select several points rather than a grid?

Wrap the latitude and longitude lists in DataArrays that share a dimension name, such as city. Plain lists select every combination, which turned three cities into nine series.

Why did selecting 20 January give me February?

Monthly data are stamped on the first of the month, and 1 February is nearer to 20 January than 1 January is. Select monthly data with a partial string such as 2024-01 instead of method nearest.

Should I use sel with nearest or interp for a point?

Nearest returns a real cell value; interp returns a blend of neighbouring cells. On a 2.5° grid they differed by up to 0.83 °C, and interp returned NaN near the longitude seam.

Does a time slice beyond the end of the data raise an error?

No. A slice from 2025-01 to 2030-12 returned the 14 months that exist, so compare the first and last timestamps of the result with the period you asked for.