How do I slice an X-array netcdf dataset across multiple x, y and time values?

Question

I have created a function to select netcdf data for a single x, y and time using X-array but how can I input multiple x,y,time values (i.e. a route) and return the values from the netcdf? i.e. windspeed for a series of x,y,time arrays.

The only way I can think of is to loop through my x,y,time arrays and run the function each time but this seems very inefficient. Any suggestions on the best method?

My input dataset is:

| x y timestamp | | 0 77.49316 6.320333 2018-08-01 00:10:00+00:00 | 1 77.64450 6.287833 2018-08-01 00:44:00+00:00 | 2 77.98483 6.202500 2018-08-01 02:02:00+00:00

and I would like to return a netcdf variable as a column. Below is my code to select for a single location.

# Scan the folder for netcdf datasets and return variables
def listmodelvars(model):
    varlist = xr.open_mfdataset(os.path.join(model, '*.nc'))
    print(varlist.data_vars.keys())
    return varlist    

ds = listmodelvars(ncep_model)
KeysView(Data variables:
    tmp2m     (time, latitude, longitude) float64 dask.array<shape=(270, 311, 720), chunksize=(9, 311, 720)>
    rh2m      (time, latitude, longitude) float64 dask.array<shape=(270, 311, 720), chunksize=(9, 311, 720)>
    tmpsfc    (time, latitude, longitude) float64 dask.array<shape=(270, 311, 720), chunksize=(9, 311, 720)>
    pratesfc  (time, latitude, longitude) float64 dask.array<shape=(270, 311, 720), chunksize=(9, 311, 720)>
    prmslmsl  (time, latitude, longitude) float64 dask.array<shape=(270, 311, 720), chunksize=(9, 311, 720)>
    ugrd10m   (time, latitude, longitude) float64 dask.array<shape=(270, 311, 720), chunksize=(9, 311, 720)>
    vgrd10m   (time, latitude, longitude) float64 dask.array<shape=(270, 311, 720), chunksize=(9, 311, 720)>)
In [4]:

# Variables    
parameter = 'tmp2m'
start_date = '2018-08-01'
end_date = '2018-08-31'
lat = 0
lon = 0    

# Define a function to return time-series data for a single location    
def array_2d(model, parameter, start_date, end_date, lat, lon):
    data = xr.open_mfdataset(os.path.join(str(model), '*.nc'))[str(parameter)].sel(
        time=slice(str(start_date) + 'T00:00:00', str(end_date) + 'T00:00:00'), 
        latitude=lat, longitude=lon)
    return data    

data_2d = array_2d(ncep_model, parameter, start_date, end_date, lat, lon)