Re-projection problem when using DataArrays
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)):
~/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 DataArrays to DataSets 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.
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!
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.
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.