使用dask.array.map_overlap时dask输出的问题



我想使用dask.array.map_overlap来处理scipy插值函数。然而,我不断遇到我无法理解的错误,希望有人能回答我

如果我想运行.comute((.,这是我收到的错误消息

ValueError:无法将输入数组从形状(1070,0(广播到形状(1045,0(

为了解决这个问题,我开始使用.To_delayed((来检查每个分区的输出,这就是我发现的。


下面是我的python代码。

步骤1。通过Xarray加载netCDF文件,然后输出到具有块大小(400400(的dask.array

df = xr.open_dataset('./Brazil Sentinal2 Tile/' + data_file +'.nc')
lon, lat = df['lon'].data, df['lat'].data
slon = da.from_array(df['lon'], chunks=(400,400))
slat = da.from_array(df['lat'], chunks=(400,400))
data = da.from_array(df.isel(band=0).__xarray_dataarray_variable__.data, chunks=(400,400))

步骤2。为da.map_overlap声明一个函数,使用

def sumsum2(lon,lat,data,  hex_res=10):
hex_col = 'hex' + str(hex_res)
lon_max, lon_min = lon.max(), lon.min()
lat_max, lat_min = lat.max(), lat.min()

b = box(lon_min, lat_min, lon_max, lat_max, ccw=True)
b = transform(lambda x, y: (y, x), b)
b = mapping(b)

target_df = pd.DataFrame(h3.polyfill( b, hex_res), columns=[hex_col])    
target_df['lat'] = target_df[hex_col].apply(lambda x: h3.h3_to_geo(x)[0])
target_df['lon'] = target_df[hex_col].apply(lambda x: h3.h3_to_geo(x)[1])
tlon, tlat = target_df[['lon','lat']].values.T    
abc = lNDI(points=(lon.ravel(), lat.ravel()), 
values= data.ravel())(tlon,tlat)
target_df['out'] = abc
print(np.stack([tlon, tlat, abc],axis=1).shape)
return np.stack([tlon, tlat, abc],axis=1)

步骤3。应用da.map_overlap

b = da.map_overlap(sumsum2, slon[:1200,:1200], slat[:1200,:1200], data[:1200,:1200], depth=10, trim=True, boundary=None, align_arrays=False, dtype='float64', 
)

步骤4。使用to_delayerd((测试输出形状

print(b.to_delayed().flatten()[0].compute().shape, )
print(b.to_delayed().flatten()[1].compute().shape)

(1065,3(
(1045,0(
(1090,3(
(1070,0(

这意味着da.map_overlap的输出仅输出1-D维度(即(1045,0(和(1070,0((,而在da.map_overlap中,我准备的输出是2-D维度(为(1065,3(和(1090,3((。

此外,如果我关闭trim参数,即

c = da.map_overlap(sumsum2, 
slon[:1200,:1200], 
slat[:1200,:1200], 
data[:1200,:1200], 
depth=10,
trim=False,
boundary=None,
align_arrays=False,
dtype='float64', 
)
print(c.to_delayed().flatten()[0].compute().shape, )
print(c.to_delayed().flatten()[1].compute().shape)

输出变为

(1065,3(

这是说当trim=True时,我删掉了所有内容?

因为。。。

#-- print out the values 
b.to_delayed().flatten()[0].compute()[:10,:]

(1065,3(
数组([],shape=(1045,0(,dtype=float64(

while。。。

#-- print out the values
c.to_delayed().flatten()[0].compute()[:10,:]

数组([[-47.83683837、-18.98359832、1395.01848583]、
[47.8482856、-18.99038681、2663.68391094]、
[47.82800624、-18.99207069、1465.56517187]、
[47.81897323、-18.97919009、2769.91556363]、
[47.82066663、-19.00712956、1607.85927095]、
[47.82696896、-18.97167714、2110.7516765],
[47.81562653,-18.98302933,2662.72112163],
【-47.82176881、-18.98594465、2201.83205114】,
【-47.84567、-18.97512514、1283.20631652】,
【-47.8 4343568、-18.97270783、1282.92117225】(

对此有什么想法吗?

谢谢。

我想我得到了答案。如果我错了,请告诉我。

  1. 我不允许使用trim=True是因为我更改了输出数组的形状(上网后,我注意到输出数组的颜色应该与输入数组的颜色相同(。由于我改变了形状,dask不知道如何处理它,所以它将空数组返回给我(奇怪(。

  2. 由于我没有要求剪切缓冲区,所以现在可以输出返回值,而不是使用trim=False。(虽然我仍然不知道为什么dask不能连接分块的数组,但相信这也与形状有关(

  3. 解决方案是在da上使用延迟函数。concatenate,即

delayed(da.concatenate)([e.to_delayed().flatten()[idx] for idx in range(len(e.to_delayed().flatten()))])

在这种情况下,我们不依赖map_overlap中的concat函数,而是使用我们自己的concat来组合我们想要的输出。

最新更新