我需要计算以下积分：

哪里 H(.) 是heaviside函数。
现在，我知道如何手动迭代数组并计算 I 这样做，但我显然想从xarray或者至少是numpy那里得到加速的好处。
这里一个巨大的复杂性是，我的数据集非常庞大，我必须及时进行数据块处理。但我在到达那里之前就失败了。
例如，我尝试从一个仅依赖于 x , y 及 z 在下面的mwe中，但我在呼叫时仍然收到一个错误 get_z_star() 因为显然 xr.apply_ufunc 正在传递numpy数组而不是dataarray，这正是我想要的：

import numpy as np
import xarray as xr

x = np.linspace(0, 1, 10)
b_tot = xr.DataArray(np.random.randn(10,10,10,10), dims=('x', 'y', 'z', 't'), coords=dict(x=x,
                                                                                          y=x,
                                                                                          z=x,
                                                                                          t=x))
def get_H_func_int(b, b_prime):
    """ b needs to be a float
        b_prime needs to be a DataArray but can't depend on time
    """
    H_func = np.heaviside(b_prime-b, 1/2)
    return H_func.integrate(('x', 'y', 'z'))

def get_z_star(b, b_prime,**kwargs):
    """ b here can be an dataArray but can't depend on time
        b_prime needs to be a DataArray but can't depend on time
    """
    I = xr.apply_ufunc(get_H_func_int, b, b_prime,**kwargs)
    return I

b = b_tot.min()
b_prime = b_tot.isel(t=0)
print(get_H_func_int(b, b_prime))
print(get_z_star(b_prime, b_prime))

有什么想法吗？