Dark xarray, созданный из операции 'isel', не может загрузить значения (слишком медленно) - PullRequest
0 голосов
/ 11 апреля 2020

У меня есть массив 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?

Спасибо!

1 Ответ

0 голосов
/ 13 апреля 2020

У меня нет набора данных для прохождения теста, и трудно сказать по ограниченной информации. Вот то, что моя первая догадка.

Когда вы используете xr.open_mfdataset для загрузки информации, вы фактически назначаете новое измерение «ocean_time» для топографической переменной h. Измерение "h" должно быть перенесено из [eta_rho, xi_rho] в [ocean_time, eta_rho, xi_rho]. Поэтому при работе с

depth_max=1e3  
deep = np.where( DS.h >= depth_max)

размеры глубины не являются [eta_rho, xi_rho]; на самом деле они также [ocean_time, eta_rho, xi_rho]. Я думаю, таким образом, проблема возникает здесь. Вы перемещаете систему координат из [ocean_time, s_rho, eta_rho, xi_rho] в [ocean_time, s_rho, gthmax], используя

eta_rho=xr.DataArray(deep[0],dims='gthmax')  
xi_rho=xr.DataArray(deep[1],dims='gthmax')  

Обратите внимание, что длина измерения «ocean_time» равна 1456, где она намного больше вашего горизонтального размера (eta_rho: 489, xi_rho: 655). Я считаю, что это запутает xarray / dask и замедлит процессы, когда вы захотите вызвать значение.


Поскольку я получил ваш набор данных, позвольте мне обновить мой ответ.

df = xr.open_mfdataset('./*.nc',  
                       combine='by_coords',  
                       chunks={'ocean_time':10}, 
                       parallel=True, 
)                                                                                                                 

temp = df[['temp']]                                                                                            
temp1 = temp.sel(xi_rho=xi_rho,eta_rho=eta_rho)                                                                   
%time temp1.isel(ocean_time=0,s_rho=46,z=100).compute()                                                          
CPU times: user 2.2 s, sys: 1.11 s, total: 3.31 s
Wall time: 3.05 s
...