I’m trying to compute the displacement of particules between two subsequent images and display the resulting displacement map as an overlay in Napari. I’m working with relatively large datasets and so preprocessing every frames in advance isn’t a very practical solution. The data is loaded as a 3d Dask array (2d + time) and I would like to lazily compute the displacement using the `extended_search_area_piv`

function from the `openpiv`

package in a similar fashion as demonstrated here. The openpiv function expects a pair of 2D frames as an argument. I attempted to use the `map_blocks()`

from `dask.array`

by rechunking the data along the time axis to single frames

`stack = stack.rechunk({0:1})`

and providing the piv function with two offset copies of the 2d+time array

```
area_search = partial(pyprocess.extended_search_area_piv, window_size=32, overlap=16, search_area_size=32)
piv_disp = da.map_blocks(area_search,
stack[1:],
stack[:-1],
chunks=(1, *x.shape), dtype='float64')
```

Once I try computing the result for a slice I get an error related to shape of the provided data.

```
~/.pyenv/versions/neuro_39/lib/python3.9/site-packages/OpenPIV-0.23.4-py3.9-macosx-11.2-x86_64.egg/openpiv/pyprocess.py in extended_search_area_piv()
753
754 if (window_size > frame_a.shape[0]) or (window_size > frame_a.shape[1]):
--> 755 raise ValueError("window size cannot be larger than the image")
756
757 # get field shape
ValueError: window size cannot be larger than the image
```

I’m assuming that this is a result of the 3rd dimension being added to the expected 2d frames. Is there a way to perform the computation on demand on 2d slices of the array?

Note that I already tried doing a stack of delayed function calls:

```
x, y = openpiv.pyprocess.get_coordinates(image_size=stack.shape[-2:], search_area_size=32, overlap=16)
@delayed
def displacement_map(frame_a, frame_b):
u, v, s2n = pyprocess.extended_search_area_piv(frame_a, frame_b, window_size=32, overlap=16, search_area_size=32)
displ = np.abs(np.dot(u, v))
return displ
UV = da.stack([da.from_delayed(displacement_map(stack[tt].astype(np.int32),
stack[tt+1].astype(np.int32)),
x.shape, dtype='float64')
for tt in range(stack.shape[0]-1)], axis=0)
```

however the ‘precompile’ time in this approach turned out to be extremely long.

If you have any suggestions let me know Thanks!