Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Improved handling of non-default central longitudes #281

Merged
merged 7 commits into from
Jan 18, 2019
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension


Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
6 changes: 4 additions & 2 deletions examples/gallery/bokeh/filled_contours.ipynb
Original file line number Diff line number Diff line change
Expand Up @@ -46,7 +46,9 @@
" return lons, lats, data\n",
"\n",
"lons, lats, data = sample_data()\n",
"contours = hv.operation.contours(gv.Image((lons, lats, data)), filled=True, levels=8).redim.range(z=(-1.1, 1.1))"
"\n",
"# Make sure we declare the central longitude\n",
"contours = gv.FilledContours((lons, lats, data), crs=crs.PlateCarree(central_longitude=180))"
]
},
{
Expand All @@ -63,7 +65,7 @@
"outputs": [],
"source": [
"(contours * gf.coastline).opts(\n",
" opts.Polygons(colorbar=True, width=600, height=300, color_levels=10, projection=crs.Mollweide(), cmap='nipy_spectral'))"
" opts.FilledContours(levels=8, color_levels=10, cmap='nipy_spectral', colorbar=True, width=600, height=300, projection=crs.Mollweide()))"
]
}
],
Expand Down
6 changes: 4 additions & 2 deletions examples/gallery/matplotlib/filled_contours.ipynb
Original file line number Diff line number Diff line change
Expand Up @@ -45,7 +45,9 @@
" return lons, lats, data\n",
"\n",
"lons, lats, data = sample_data()\n",
"contours = gv.FilledContours((lons, lats, data))"
"\n",
"# Make sure we declare the central longitude\n",
"contours = gv.FilledContours((lons, lats, data), crs=ccrs.PlateCarree(central_longitude=180))"
]
},
{
Expand All @@ -61,7 +63,7 @@
"metadata": {},
"outputs": [],
"source": [
"contours.opts(projection=ccrs.Mollweide(), cmap='nipy_spectral') * gf.coastline"
"contours.opts(colorbar=True, levels=8, color_levels=10, cmap='nipy_spectral', projection=ccrs.Mollweide()) * gf.coastline"
]
}
],
Expand Down
26 changes: 25 additions & 1 deletion examples/user_guide/Projections.ipynb
Original file line number Diff line number Diff line change
Expand Up @@ -174,7 +174,31 @@
"metadata": {},
"outputs": [],
"source": [
"features.opts(projection=crs.Orthographic(central_longitude=0, central_latitude=30), global_extent=True)"
"features.opts(projection=crs.Orthographic(central_longitude=-90, central_latitude=30), global_extent=True)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"This is particularly important when the data is specified in the longitudes of the data are defined in the range 0 to 360 rather than the default -180 to 180. In this case ensure that the ``crs`` of the element declares the actual ``central_longitude``, e.g. 180 degrees rather than the default 0 degrees. Let us draw the same polygon with the default ``crs`` (centered on 0) and a custom ``crs`` centered on 180 degrees:"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"poly_data = [(180, 0), (270, 45), (360, 0)]\n",
"features * gv.Polygons(poly_data) * gv.Polygons(poly_data, crs=crs.PlateCarree(central_longitude=180)).opts(global_extent=True)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"As we can see the default projection wraps the values above 180, while the custom projection automatically translates the polygon by the ``central_longitude``."
]
},
{
Expand Down
2 changes: 0 additions & 2 deletions examples/user_guide/Resampling_Grids.ipynb
Original file line number Diff line number Diff line change
Expand Up @@ -253,8 +253,6 @@
"outputs": [],
"source": [
"ds = xr.tutorial.open_dataset('rasm').load()\n",
"# Adjust the longitudes so they are in the range -180 to 180\n",
"ds.xc.values[ds.xc.values>180] -= 360 \n",
"ds"
]
},
Expand Down
26 changes: 10 additions & 16 deletions geoviews/operation/projection.py
Original file line number Diff line number Diff line change
@@ -1,13 +1,12 @@
import logging

import param
import shapely
import numpy as np

from cartopy import crs as ccrs
from cartopy.img_transform import warp_array, _determine_bounds
from holoviews.core.data import MultiInterface
from holoviews.core.util import cartesian_product, get_param_values, max_extents, pd
from holoviews.core.util import cartesian_product, get_param_values, pd
from holoviews.operation import Operation
from shapely.geometry import Polygon, MultiPolygon
from shapely.geometry.collection import GeometryCollection
Expand Down Expand Up @@ -49,37 +48,28 @@ class project_path(_project_operation):
supported_types = [Polygons, Path, Contours, EdgePaths]

def _process_element(self, element):
if not len(element):
if not bool(element):
return element.clone(crs=self.p.projection)

crs = element.crs
cylindrical = isinstance(crs, ccrs._CylindricalProjection)
proj = self.p.projection
if (isinstance(crs, ccrs.PlateCarree) and not isinstance(proj, ccrs.PlateCarree)
and crs.proj4_params['lon_0'] != 0):
element = self.instance(projection=ccrs.PlateCarree())(element)

if isinstance(proj, ccrs.CRS) and not isinstance(proj, ccrs.Projection):
raise ValueError('invalid transform:'
' Spherical contouring is not supported - '
' consider using PlateCarree/RotatedPole.')

boundary = Polygon(crs.boundary)
bounds = [round(b, 10) for b in boundary.bounds]
xoffset = round((boundary.bounds[2]-boundary.bounds[0])/2.)
if isinstance(element, Polygons):
geoms = polygons_to_geom_dicts(element, skip_invalid=False)
else:
geoms = path_to_geom_dicts(element, skip_invalid=False)

data_bounds = max_extents([g['geometry'].bounds for g in geoms])
total_bounds = tuple(round(b, 10) for b in data_bounds)

projected = []
for path in geoms:
geom = path['geometry']
if (cylindrical and total_bounds[0] >= (bounds[0]+xoffset) and
total_bounds[2] > (bounds[2]+xoffset//2) and
total_bounds[2] <= bounds[2]+xoffset):
# Offset if lon and not centered on 0 longitude
# i.e. lon_min > 0 and lon_max > 270
geom = shapely.affinity.translate(geom, xoff=-xoffset)

# Ensure minimum area for polygons (precision issues cause errors)
if isinstance(geom, Polygon) and geom.area < 1e-15:
Expand Down Expand Up @@ -313,6 +303,10 @@ def _process(self, img, key=None):
else:
projected, extents = arr, trgt_ext
arrays.append(projected)

if xn == 0 or yn == 0:
return img.clone([], bounds=extents, crs=proj)

xunit = ((extents[1]-extents[0])/float(xn))/2.
yunit = ((extents[3]-extents[2])/float(yn))/2.
xs = np.linspace(extents[0]+xunit, extents[1]-xunit, xn)
Expand Down
8 changes: 8 additions & 0 deletions geoviews/util.py
Original file line number Diff line number Diff line change
Expand Up @@ -28,6 +28,14 @@ def wrap_lons(lons, base, period):
def project_extents(extents, src_proj, dest_proj, tol=1e-6):
x1, y1, x2, y2 = extents

if (isinstance(src_proj, ccrs.PlateCarree) and
not isinstance(dest_proj, ccrs.PlateCarree) and
src_proj.proj4_params['lon_0'] != 0):
xoffset = src_proj.proj4_params['lon_0']
x1 = x1 - xoffset
x2 = x2 - xoffset
src_proj = ccrs.PlateCarree()

# Limit latitudes
cy1, cy2 = src_proj.y_limits
if y1 < cy1: y1 = cy1
Expand Down
2 changes: 1 addition & 1 deletion setup.py
Original file line number Diff line number Diff line change
Expand Up @@ -129,7 +129,7 @@ def package_assets(example_path):
_required = [
'bokeh >=1.0.0',
'cartopy >=0.16.0', # prevents pip alone (requires external package manager)
'holoviews >=1.11.0',
'holoviews >=1.11.1',
'numpy >=1.0',
'param >=1.6.1'
]
Expand Down