Skip to content

Commit

Permalink
Fixes for projecting path data
Browse files Browse the repository at this point in the history
  • Loading branch information
philippjfr committed May 19, 2018
1 parent 16b7bec commit 8f4e644
Show file tree
Hide file tree
Showing 4 changed files with 64 additions and 11 deletions.
1 change: 1 addition & 0 deletions geoviews/element/geo.py
Original file line number Diff line number Diff line change
Expand Up @@ -285,6 +285,7 @@ class QuadMesh(_Element, HvQuadMesh):
def trimesh(self):
trimesh = super(QuadMesh, self).trimesh()
node_params = get_param_values(trimesh.nodes)
node_params['crs'] = self.crs
nodes = TriMesh.node_type(trimesh.nodes.data, **node_params)
return TriMesh((trimesh.data, nodes), crs=self.crs,
**get_param_values(trimesh))
Expand Down
28 changes: 19 additions & 9 deletions geoviews/operation/projection.py
Original file line number Diff line number Diff line change
Expand Up @@ -10,7 +10,7 @@
from ..element import (Image, Shape, Polygons, Path, Points, Contours,
RGB, Graph, Nodes, EdgePaths, QuadMesh, VectorField,
HexTiles, Labels)
from ..util import project_extents, geom_to_array
from ..util import project_extents, geom_to_array, wrap_path_data


class _project_operation(Operation):
Expand Down Expand Up @@ -65,19 +65,30 @@ def _process_element(self, element):
xs, ys = geom[xdim.name], geom[ydim.name]
if not len(xs):
continue
path = geom_type(np.column_stack([xs, ys]))
path = path.intersection(boundary_poly)
proj = self.p.projection.project_geometry(path, element.crs)
proj_arr = geom_to_array(proj)

proj_arr = self.p.projection.quick_vertices_transform(
np.column_stack([xs, ys]), element.crs)

if proj_arr is None:
vertices = np.column_stack([xs, ys])
if isinstance(element.crs, ccrs.PlateCarree):
vertices = wrap_path_data(vertices, element.crs, element.crs)
path = geom_type(vertices)
path = path.intersection(boundary_poly)
proj = self.p.projection.project_geometry(path, element.crs)
proj_arr = geom_to_array(proj)
geom[xdim.name] = proj_arr[:, 0]
geom[ydim.name] = proj_arr[:, 1]
projected.append(geom)
else:
# Handle iso-contour case
data = {k: vals[0] for k, vals in data.items()}
geom = path.geom()
if boundary_poly:
geom = geom.intersection(boundary_poly)
if isinstance(element.crs, ccrs.PlateCarree):
vertices = wrap_path_data(path.array([0, 1]), element.crs, element.crs)
geom = type(element)([vertices]).geom()
else:
geom = path.geom()

if not geom:
continue
proj = self.p.projection.project_geometry(geom, element.crs)
Expand All @@ -87,7 +98,6 @@ def _process_element(self, element):
return element.clone(projected, crs=self.p.projection)



class project_shape(_project_operation):
"""
Projects Shape Element from the source coordinate reference system
Expand Down
4 changes: 2 additions & 2 deletions geoviews/plotting/bokeh/plot.py
Original file line number Diff line number Diff line change
Expand Up @@ -98,7 +98,7 @@ def _postprocess_hover(self, renderer, source):
except:
CustomJSHover = None
if (not self.geographic or None in (hover, CustomJSHover) or
isinstance(hover.tooltips, basestring)):
isinstance(hover.tooltips, basestring) or self.projection is not GOOGLE_MERCATOR):
return
element = self.current_frame
xdim, ydim = [dimension_sanitizer(kd.name) for kd in element.kdims]
Expand Down Expand Up @@ -147,7 +147,7 @@ def get_extents(self, element, ranges):

def get_data(self, element, ranges, style):
proj = self.projection
if self._project_operation and self.geographic and element.crs != proj:
if self._project_operation and self.geographic:
element = self._project_operation(element, projection=proj)
return super(GeoPlot, self).get_data(element, ranges, style)

Expand Down
42 changes: 42 additions & 0 deletions geoviews/util.py
Original file line number Diff line number Diff line change
@@ -1,4 +1,6 @@
import numpy as np
import shapely.geometry as sgeom

from cartopy import crs as ccrs
from shapely.geometry import (MultiLineString, LineString, MultiPolygon, Polygon)

Expand Down Expand Up @@ -114,9 +116,19 @@ def polygon_to_geom(poly, multi=True):
else:
path = path
lines.append(path)
lines = [to_ccw(line) for line in lines]
return MultiPolygon(lines) if multi else lines


def to_ccw(geom):
"""
Reorients polygon to be wound counter-clockwise.
"""
if isinstance(geom, sgeom.Polygon) and not geom.exterior.is_ccw:
geom = sgeom.polygon.orient(geom)
return geom


def geom_to_arr(geom):
arr = geom.array_interface_base['data']
if (len(arr) % 2) != 0:
Expand Down Expand Up @@ -164,3 +176,33 @@ def geo_mesh(element):
zs = np.ma.concatenate([zs, zs[:, 0:1]], axis=1)
return xs, ys, zs


def wrap_data(vertices, src_crs, tgt_crs):
"""
Wraps path coordinates along the longitudinal axis.
"""
self_params = tgt_crs.proj4_params.copy()
src_params = src_crs.proj4_params.copy()
self_params.pop('lon_0'), src_params.pop('lon_0')

xs, ys = vertices[:, 0], vertices[:, 1]
potential = (self_params == src_params and
tgt_crs.y_limits[0] <= ys.min() and
tgt_crs.y_limits[1] >= ys.max())
if not potential:
return vertices

bboxes, proj_offset = tgt_crs._bbox_and_offset(src_crs)
mod = np.diff(src_crs.x_limits)[0]
x_lim = xs.min(), xs.max()
for poly in bboxes:
# Arbitrarily choose the number of moduli to look
# above and below the -180->180 range. If data is beyond
# this range, we're not going to transform it quickly.
for i in [-1, 0, 1, 2]:
offset = mod * i - proj_offset
if ((poly[0] + offset) <= x_lim[0] and
(poly[1] + offset) >= x_lim[1]):
vertices = vertices + [[-offset, 0]]
break
return vertices

0 comments on commit 8f4e644

Please sign in to comment.