У меня есть массив 4d dask, размеры которого соответствуют (время, глубина, широта, долгота). К вашему сведению, это набор данных Oceani c.
# Create xarray dataset object for ocean temperature (time x depth x lat x lon)
DS=xr.open_mfdataset([outdir + s for s in flist],combine='by_coords',chunks={'xi_rho':25,'eta_rho':25})
DS.temp
Вывод:
xarray.DataArray
'temp' (ocean_time: 1456, s_rho: 50, eta_rho: 489, xi_rho: 655)
dask.array<chunksize=(1456, 50, 25, 25), meta=np.ndarray>
Когда я пытаюсь загрузить 1d-вектор из этого массива dask, проблем нет .
T=DS.temp
%time
T.isel(ocean_time=0,eta_rho=100,xi_rho=500).values
Вывод (я опускаю значения ниже):
CPU times: user 0 ns, sys: 0 ns, total: 0 ns
Wall time: 5.96 µs
Теперь я хочу выбрать только те (широта, долгота), где дно океана глубже, чем, скажем, 1000 м
depth_max=1e3;
deep=np.where(DS.h >=depth_max); # DS.h is the depth of the ocean bottom.
# Number of locations where the ocean is deeper than depth_max
xy_num=len(deep[0])
Это дает мне кортеж 'deep'
, первый элемент которого (deep[0])
является списком всех значений 'eta_rho'
(индекс широты), удовлетворяющих условию. Второй элемент кортежа (deep[1])
представляет собой список соответствующих значений 'xi_rho'
(индекс долготы). Итак, я могу построить индексные пары (lat,lon)
, используя (deep[0][0],deep[1][0]), (deep[0][1],deep[1][1]), (deep[0][2],deep[1][2]), (deep[0][3],deep[1][3])
, et c.
Это здорово, потому что я хочу создать одну координату, которая будет проходить через пары (lat,lon)
, удовлетворяющие вышеуказанному условию. Цель состоит в том, чтобы получить от (time,depth,lat,lon) -> (time,depth,gthmax)
, где gthmax
является новой координатой, описанной abvoe. Я делаю это так:
# Picking only those (lat,lon) where the condition is satisfied
T_deep=T.isel(eta_rho=xr.DataArray(deep[0],dims='gthmax'),xi_rho=xr.DataArray(deep[1],dims='gthmax'))
# Explicitly assign the new coordinate
T_deep=T_deep.assign_coords({"gthmax":range(0,xy_num)})
# Create chunks along this new coordinate
T_deep=T_deep.chunk({'gthmax':1000})
T_deep
Вывод (просто показывает размеры):
xarray.DataArray 'temp' (ocean_time: 1456, s_rho: 50, gthmax: 133446)
dask.array<chunksize=(1456, 50, 1000), meta=np.ndarray>
Вот где возникает проблема. Когда я пытаюсь получить доступ к значениям из этого нового массива 3d dask, даже в одной точке, приведенная ниже команда никогда не завершает завершение, и мне приходится заканчивать тем, что убил ядро. Я также пытался использовать load()
и compute()
, но безрезультатно.
T_deep.isel(ocean_time=0,s_rho=46,gthmax=100).values
Есть мысли о том, где я облажался при преобразовании исходного массива 4d dask в массив 3d dask?
Спасибо!