Re-projection problem when using DataArrays

[jlstevens]

Given a DataArray called data that looks like this:

<xarray.DataArray (y: 7941, x: 7821)>
dask.array<shape=(7941, 7821), dtype=float64, chunksize=(256, 256)>
Coordinates:
    band     int64 1
  * y        (y) float64 4.426e+06 4.426e+06 4.426e+06 4.426e+06 4.426e+06 ...
  * x        (x) float64 2.433e+05 2.433e+05 2.434e+05 2.434e+05 2.434e+05 ...

I get the following traceback if I try gv.Image(data, crs=ccrs.epsg(32611)):

--------------------------------------------------------------------------- ValueError Traceback (most recent call last) ~/miniconda3/envs/analytics/lib/python3.6/site-packages/geoviews/util.py in project_extents(extents, src_proj, dest_proj, tol) 48 try: ---> 49 geom_in_crs = dest_proj.project_geometry(geom_in_src_proj, src_proj) 50 except ValueError: ~/miniconda3/envs/analytics/lib/python3.6/site-packages/cartopy/crs.py in project_geometry(self, geometry, src_crs) 179 raise ValueError('Unsupported geometry ' --> 180 'type {!r}'.format(geom_type)) 181 return getattr(self, method_name)(geometry, src_crs) ValueError: Unsupported geometry type 'GeometryCollection' During handling of the above exception, another exception occurred: ValueError Traceback (most recent call last) ~/miniconda3/envs/analytics/lib/python3.6/site-packages/IPython/core/formatters.py in __call__(self, obj, include, exclude) 968 969 if method is not None: --> 970 return method(include=include, exclude=exclude) 971 return None 972 else: ~/Desktop/development/holoviews/holoviews/core/dimension.py in _repr_mimebundle_(self, include, exclude) 1263 combined and returned. 1264 """ -> 1265 return Store.render(self) 1266 1267 ~/Desktop/development/holoviews/holoviews/core/options.py in render(cls, obj) 1287 data, metadata = {}, {} 1288 for hook in hooks: -> 1289 ret = hook(obj) 1290 if ret is None: 1291 continue ~/Desktop/development/holoviews/holoviews/ipython/display_hooks.py in pprint_display(obj) 270 if not ip.display_formatter.formatters['text/plain'].pprint: 271 return None --> 272 return display(obj, raw_output=True) 273 274 ~/Desktop/development/holoviews/holoviews/ipython/display_hooks.py in display(obj, raw_output, **kwargs) 240 elif isinstance(obj, (CompositeOverlay, ViewableElement)): 241 with option_state(obj): --> 242 output = element_display(obj) 243 elif isinstance(obj, (Layout, NdLayout, AdjointLayout)): 244 with option_state(obj): ~/Desktop/development/holoviews/holoviews/ipython/display_hooks.py in wrapped(element) 140 try: 141 max_frames = OutputSettings.options['max_frames'] --> 142 mimebundle = fn(element, max_frames=max_frames) 143 if mimebundle is None: 144 return {}, {} ~/Desktop/development/holoviews/holoviews/ipython/display_hooks.py in element_display(element, max_frames) 186 return None 187 --> 188 return render(element) 189 190 ~/Desktop/development/holoviews/holoviews/ipython/display_hooks.py in render(obj, **kwargs) 63 renderer = renderer.instance(fig='png') 64 ---> 65 return renderer.components(obj, **kwargs) 66 67 ~/Desktop/development/holoviews/holoviews/plotting/bokeh/renderer.py in components(self, obj, fmt, comm, **kwargs) 268 # Bokeh has to handle comms directly in <0.12.15 269 comm = False if bokeh_version < '0.12.15' else comm --> 270 return super(BokehRenderer, self).components(obj,fmt, comm, **kwargs) 271 272 ~/Desktop/development/holoviews/holoviews/plotting/renderer.py in components(self, obj, fmt, comm, **kwargs) 327 plot = obj 328 else: --> 329 plot, fmt = self._validate(obj, fmt) 330 331 widget_id = None ~/Desktop/development/holoviews/holoviews/plotting/renderer.py in _validate(self, obj, fmt, **kwargs) 226 if isinstance(obj, tuple(self.widgets.values())): 227 return obj, 'html' --> 228 plot = self.get_plot(obj, renderer=self, **kwargs) 229 230 fig_formats = self.mode_formats['fig'][self.mode] ~/Desktop/development/holoviews/holoviews/plotting/bokeh/renderer.py in get_plot(self_or_cls, obj, doc, renderer) 153 curdoc().theme = self_or_cls.theme 154 doc.theme = self_or_cls.theme --> 155 plot = super(BokehRenderer, self_or_cls).get_plot(obj, renderer) 156 plot.document = doc 157 return plot ~/Desktop/development/holoviews/holoviews/plotting/renderer.py in get_plot(self_or_cls, obj, renderer) 213 init_key = tuple(v if d is None else d for v, d in 214 zip(plot.keys[0], defaults)) --> 215 plot.update(init_key) 216 else: 217 plot = obj ~/Desktop/development/holoviews/holoviews/plotting/plot.py in update(self, key) 511 def update(self, key): 512 if len(self) == 1 and ((key == 0) or (key == self.keys[0])) and not self.drawn: --> 513 return self.initialize_plot() 514 item = self.__getitem__(key) 515 self.traverse(lambda x: setattr(x, '_updated', True)) ~/Desktop/development/holoviews/holoviews/plotting/bokeh/element.py in initialize_plot(self, ranges, plot, plots, source) 729 self.handles['plot'] = plot 730 --> 731 self._init_glyphs(plot, element, ranges, source) 732 if not self.overlaid: 733 self._update_plot(key, plot, style_element) ~/Desktop/development/holoviews/holoviews/plotting/bokeh/element.py in _init_glyphs(self, plot, element, ranges, source) 680 else: 681 style = self.style[self.cyclic_index] --> 682 data, mapping, style = self.get_data(element, ranges, style) 683 current_id = element._plot_id 684 if source is None: ~/miniconda3/envs/analytics/lib/python3.6/site-packages/geoviews/plotting/bokeh/plot.py in get_data(self, element, ranges, style) 68 def get_data(self, element, ranges, style): 69 if self._project_operation and self.geographic and element.crs != DEFAULT_PROJ: ---> 70 element = self._project_operation(element) 71 return super(GeoPlot, self).get_data(element, ranges, style) 72 ~/Desktop/development/holoviews/holoviews/core/operation.py in __call__(self, element, **params) 161 operation=self, kwargs=params) 162 elif isinstance(element, ViewableElement): --> 163 processed = self._apply(element) 164 elif isinstance(element, DynamicMap): 165 if any((not d.values) for d in element.kdims): ~/Desktop/development/holoviews/holoviews/core/operation.py in _apply(self, element, key) 119 for hook in self._preprocess_hooks: 120 kwargs.update(hook(self, element)) --> 121 ret = self._process(element, key) 122 for hook in self._postprocess_hooks: 123 ret = hook(self, ret, **kwargs) ~/miniconda3/envs/analytics/lib/python3.6/site-packages/geoviews/operation.py in _process(self, img, key) 200 xn, yn = img.interface.shape(img, gridded=True)[:2] 201 px0, py0, px1, py1 = project_extents((x0, y0, x1, y1), --> 202 img.crs, proj) 203 src_ext, trgt_ext = (x0, x1, y0, y1), (px0, px1, py0, py1) 204 arrays = [] ~/miniconda3/envs/analytics/lib/python3.6/site-packages/geoviews/util.py in project_extents(extents, src_proj, dest_proj, tol) 54 'to %s projection. Ensure the coordinate ' 55 'reference system (crs) matches your data.' % ---> 56 (src_name, dest_name)) 57 else: 58 geom_in_crs = boundary_poly.intersection(domain_in_src_proj) ValueError: Could not project data from _EPSGProjection projection to Mercator projection. Ensure the coordinate reference system (crs) matches your data.

This happens whether or not I drop the irrelevant ``band`` dimension (as you might expect) but it does work correctly once expressed as a ``Dataset`` instead. It seems to me that we ought to be able to promote ``DataArray``s to ``DataSet``s that have a single data variable as all the necessary information is available from the ``DataArray`` (except for the name of the data variables which should match our default of 'z'). That way you can have a data array promoted to a dataset which you could redim from the 'z' default to give it a suitable name for the data variable. Of course if you specify a vdim that isn't 'z' and supply a ``DataArray``, you would get the usual error about data missing for the requested vdim. One motivation to support ``DataArray`` is that this is the format returned by ``xr.open_rasterio(f)`` for file ``f``.

[jlstevens]

The workaround I have found is as follows:

data.to_dataset(name='data')[['x','y', 'data']]

Using to_dataset to promote the array to a datset is fine, but having to set the ordering is rather annoying!

[jbednar]

I believe the ordering/naming issue has been discussed extensively somewhere, but I can't find it anywhere to link it in here. Maybe @philippjfr remembers.

[jlstevens]

One last comment: it is probably best to solve this at the holoviews level even though this particular problem manifests itself most obviously at the geoviews leve.