datashaded image not zooming/panning properly

[jbednar]

Using https://daac.ornl.gov/daacdata/cms/CMS_DARTE_V2/data/onroad_2017.tif and the following code, I can show a geotiff:

import xarray as xr, holoviews as hv, colorcet as cc
from holoviews.operation.datashader import datashade

hv.extension('bokeh', logo=False)
emissions = xr.open_rasterio('CMS_DARTE_V2_1735/data/onroad_2017.tif').load()[0]
datashade(hv.Image(emissions, ['x','y']), cmap=cc.fire).opts(width=500, invert_yaxis=True)
#hv.Image(emissions, ['x','y']).opts(cmap=cc.fire, width=500, invert_yaxis=True)

However, zooming and panning don't work properly; the image jumps around and gets cut off, and in some scenarios shows only garbage:

Panning in the x direction seems to work fine, so I'm suspecting that it's to do with the inverted y axis; panning in y leaves the data in the same spot in the plot axes but the bounding box of the data going up. If I uncomment the non-datashaded hv.Image line, zooming and panning work fine, but the image is only barely visible and of course it's quite slow.

I'm on git master for the various HoloViz projects:

holoviews=1.13.0a10.post6+gc2f0c7bf6 datashader=0.8.0a1.post4+ge9a2889 geoviews=1.6.4a3 hvplot=0.5.0a1.post6+g8784838 xarray=0.13.0

[xavArtley]

Hello simple reproducer if it can help

import numpy as np
import holoviews as hv
from holoviews.operation.datashader import rasterize
hv.extension('bokeh')

ls = np.linspace(0, 10, 2000)
xx, yy = np.meshgrid(ls, ls)

bounds=(-1,-1,1,1)   # Coordinate system: (left, bottom, right, top)
img = hv.Image(np.sin(xx)*np.cos(yy), bounds=bounds)
rasterize(img).opts(hv.opts.Image(invert_axes=True))

[philippjfr]

At least the simple reproducer example appears to be working now.

[philippjfr]

Don't have access to the geotiff either, so please confirm @jbednar.

[rsignell-usgs]

With the latest holoviews 1.13.5 release, it seems we are back to the wacky zoom problems when rasterize=True. 😢

It's broken with the simple reproducer example above.

# Name                    Version                   Build  Channel
holoviews                 1.13.5             pyh9f0ad1d_0    conda-forge
datashader                0.11.1             pyh9f0ad1d_0    conda-forge
geoviews-core             1.8.2                      py_0    conda-forge
hvplot                    0.6.0              pyh9f0ad1d_0    conda-forge
bokeh                     2.2.3            py38h32f6830_0    conda-forge

It's also broken with my usual demos like https://nbviewer.jupyter.org/gist/rsignell-usgs/df96011574d2bc86b04a8e9f7c4393a2

[philippjfr]

That's super, bizarre, the latest release didn't touch any of that code. It really just fixed a very small number of things.

[philippjfr]

Looking into it now.

[philippjfr]

Can't reproduce the issue for some reason, but the example works if you "pre-project" it with project=True, it doesn't seem to like projecting the image after the fact.

ds.Significant_height_of_combined_wind_waves_and_swell_surface[-42:,:,:].hvplot(x='lon', y='lat', 
                            cmap='rainbow', rasterize=True, projection=crs, project=True, coastline=True,
                            widget_type='scrubber', widget_location='bottom')

[philippjfr]

Also I don't believe the simple reproducer ever worked correctly, as this issue was never closed, but your example does not use invert_axes or invert_yaxis.

[rsignell-usgs]

@philippjfr , okay, I didn't realize these were separate issues. It does look like the simple reproducer worked on June 7? https://github.com/holoviz/holoviews/issues/4052#issuecomment-640336446

[philippjfr]

@rsignell-usgs I don't believe what you are seeing is a regression from previous version, the example you posted never seems to have worked correctly, almost certainly due to projection issues in GeoViews. The only way of handling that correctly for the time being is to use project=True which projects the mesh before it is rasterized with datashader. I've gone back to master on the day I reported the simple reproducer working and it didn't in fact work, so I was lying :/

Handling inverted axes correctly does not seem like it would be hard so I'll assign to the next milestone.

[philippjfr]

Actually there's one mystery left, which is that I was not able to reproduce the axes jumping around like in your gif.

[rsignell-usgs]

@philippjfr, do you have access to https://staging.aws-uswest2.pangeo.io? This is the environment I'm testing in.

And one more thing: the hover in my global wave example never displayed lon/lat either?

[philippjfr]

I don't seem to have access, no:

OAuth error: User is not part of REQUIRED_GITHUB_TEAMS.

[jbednar]

The original geotiff file is from https://daac.ornl.gov/cgi-bin/dsviewer.pl?ds_id=1735 , which is freely accessible but requires registering and logging in, so I sent it privately to @philippjfr.

[rsignell-usgs]

@philippjfr sorry about that --@scottyhq can you allow Philipp access to the AWS uswest2 pangeo endpoint?

You can also duplicate the issue here with same environment:

[jbednar]

I've just tested the original issue again with holoviews=1.13.5, hvplot=0.6.0, datashader=0.10.0, and found no change; still pans properly in x but pans "invertedly" in y.

[philippjfr]

I'm suspecting this issue is actually a bug in datashader but ran out of time investigating it.

[oarcher]

using dask.array raise an error in datashader with invert_axes. When not using dask, no error is raised, but the output is incorrect.

This example is perhaps not strictly minimal, but it raise the error:

import xarray as xr
import dask
import holoviews as hv
hv.extension('bokeh')
from holoviews.operation.datashader import rasterize

da = xr.DataArray(dask.array.ones((838,1259)), dims=['atrack','xtrack'])
da = da.assign_coords(atrack=da.atrack.astype(float), xtrack = da.xtrack.astype(float))

dasel = da.isel(atrack=slice(323, 519), xtrack=slice(529, 729))

rasterize(hv.Image(dasel).opts(invert_axes=False)) # OK

rasterize(hv.Image(dasel.compute()).opts(invert_axes=True))  # no error, incorrect output

rasterize(hv.Image(dasel).opts(invert_axes=True)) # ValueError: range() arg 3 must not be zero

WARNING:param.dynamic_operation: Callable raised "ValueError('range() arg 3 must not be zero')".
Invoked as dynamic_operation(height=400, scale=1.0, width=400, x_range=(528.5, 728.5), y_range=(322.5, 518.5))
WARNING:param.dynamic_operation: Callable raised "ValueError('range() arg 3 must not be zero')".
Invoked as dynamic_operation(height=400, scale=1.0, width=400, x_range=(528.5, 728.5), y_range=(322.5, 518.5))
Traceback (most recent call last):
  File "/home/oarcher/anaconda3/envs/xsar/lib/python3.9/site-packages/holoviews/plotting/util.py", line 275, in get_plot_frame
    return map_obj[key]
  File "/home/oarcher/anaconda3/envs/xsar/lib/python3.9/site-packages/holoviews/core/spaces.py", line 1341, in __getitem__
    val = self._execute_callback(*tuple_key)
  File "/home/oarcher/anaconda3/envs/xsar/lib/python3.9/site-packages/holoviews/core/spaces.py", line 1110, in _execute_callback
    retval = self.callback(*args, **kwargs)
  File "/home/oarcher/anaconda3/envs/xsar/lib/python3.9/site-packages/holoviews/core/spaces.py", line 714, in __call__
    ret = self.callable(*args, **kwargs)
  File "/home/oarcher/anaconda3/envs/xsar/lib/python3.9/site-packages/holoviews/util/__init__.py", line 1042, in dynamic_operation
    key, obj = resolve(key, kwargs)
  File "/home/oarcher/anaconda3/envs/xsar/lib/python3.9/site-packages/holoviews/util/__init__.py", line 1031, in resolve
    return key, map_obj[key]
  File "/home/oarcher/anaconda3/envs/xsar/lib/python3.9/site-packages/holoviews/core/spaces.py", line 1341, in __getitem__
    val = self._execute_callback(*tuple_key)
  File "/home/oarcher/anaconda3/envs/xsar/lib/python3.9/site-packages/holoviews/core/spaces.py", line 1110, in _execute_callback
    retval = self.callback(*args, **kwargs)
  File "/home/oarcher/anaconda3/envs/xsar/lib/python3.9/site-packages/holoviews/core/spaces.py", line 714, in __call__
    ret = self.callable(*args, **kwargs)
  File "/home/oarcher/anaconda3/envs/xsar/lib/python3.9/site-packages/holoviews/util/__init__.py", line 1043, in dynamic_operation
    return apply(obj, *key, **kwargs)
  File "/home/oarcher/anaconda3/envs/xsar/lib/python3.9/site-packages/holoviews/util/__init__.py", line 1035, in apply
    processed = self._process(element, key, kwargs)
  File "/home/oarcher/anaconda3/envs/xsar/lib/python3.9/site-packages/holoviews/util/__init__.py", line 1017, in _process
    return self.p.operation.process_element(element, key, **kwargs)
  File "/home/oarcher/anaconda3/envs/xsar/lib/python3.9/site-packages/holoviews/core/operation.py", line 194, in process_element
    return self._apply(element, key)
  File "/home/oarcher/anaconda3/envs/xsar/lib/python3.9/site-packages/holoviews/core/operation.py", line 141, in _apply
    ret = self._process(element, key)
  File "/home/oarcher/anaconda3/envs/xsar/lib/python3.9/site-packages/holoviews/operation/datashader.py", line 1477, in _process
    element = element.map(op, predicate)
  File "/home/oarcher/anaconda3/envs/xsar/lib/python3.9/site-packages/holoviews/core/data/__init__.py", line 205, in pipelined_fn
    result = method_fn(*args, **kwargs)
  File "/home/oarcher/anaconda3/envs/xsar/lib/python3.9/site-packages/holoviews/core/data/__init__.py", line 1222, in map
    return super(Dataset, self).map(*args, **kwargs)
  File "/home/oarcher/anaconda3/envs/xsar/lib/python3.9/site-packages/holoviews/core/dimension.py", line 710, in map
    return map_fn(self) if applies else self
  File "/home/oarcher/anaconda3/envs/xsar/lib/python3.9/site-packages/holoviews/core/operation.py", line 214, in __call__
    return self._apply(element)
  File "/home/oarcher/anaconda3/envs/xsar/lib/python3.9/site-packages/holoviews/core/operation.py", line 141, in _apply
    ret = self._process(element, key)
  File "/home/oarcher/anaconda3/envs/xsar/lib/python3.9/site-packages/holoviews/operation/datashader.py", line 939, in _process
    rarray = cvs.raster(xarr, upsample_method=interp,
  File "/windows_shared/Develop/datashader/datashader/core.py", line 1058, in raster
    data = resample_2d_distributed(
  File "/windows_shared/Develop/datashader/datashader/resampling.py", line 243, in resample_2d_distributed
    chunk_map = map_chunks(src.shape, (h, w), temp_chunks)
  File "/windows_shared/Develop/datashader/datashader/resampling.py", line 112, in map_chunks
    xchunks = list(range(0, outx, cxs)) + [outx]
ValueError: range() arg 3 must not be zero

[philippjfr]

@oarcher Thanks for the example, I can reproduce. For now I'd recommend just switching to an hv.QuadMesh instead of the Image but I'll make the appropriate fixes as well.