Skip to content

Commit

Permalink
Fix handling wkt strings and add pyproj.CRS support (#1139)
Browse files Browse the repository at this point in the history 
* Fix handling wkt strings

* Remove incorrect assumption

* Fix tests

* Add not implemented error

* Remove check

* Update check

* Update hvplot/util.py

* Fix test

* Fix test

* Update hvplot/util.py

* Add tests, refactor ifs

* Add support for pyproj CRS and tests

* Update hvplot/util.py

* Update hvplot/tests/testutil.py

Co-authored-by: Simon Høxbro Hansen <[email protected]>

* Update hvplot/tests/testutil.py

* Apply suggestions from code review

---------

Co-authored-by: Simon Høxbro Hansen <[email protected]>
  • Loading branch information
@ahuang11 @hoxbro
ahuang11 and hoxbro authored Sep 26, 2023
1 parent 6aba047 commit e4a65f0
Show file tree
Hide file tree
Showing 2 changed files with 85 additions and 30 deletions.
54 changes: 51 additions & 3 deletions hvplot/tests/testutil.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -217,6 +217,14 @@ def test_proj_to_cartopy(self):

assert isinstance(crs, self.ccrs.CRS)

def test_proj_to_cartopy_wkt_string(self):
from ..util import proj_to_cartopy
crs = proj_to_cartopy('GEOGCRS["unnamed",BASEGEOGCRS["unknown",DATUM["unknown",ELLIPSOID["WGS 84",6378137,298.257223563,LENGTHUNIT["metre",1,ID["EPSG",9001]]]],PRIMEM["Greenwich",0,ANGLEUNIT["degree",0.0174532925199433],ID["EPSG",8901]]],DERIVINGCONVERSION["unknown",METHOD["PROJ ob_tran o_proj=latlon"],PARAMETER["o_lon_p",0,ANGLEUNIT["degree",0.0174532925199433,ID["EPSG",9122]]],PARAMETER["o_lat_p",37.5,ANGLEUNIT["degree",0.0174532925199433,ID["EPSG",9122]]],PARAMETER["lon_0",357.5,ANGLEUNIT["degree",0.0174532925199433,ID["EPSG",9122]]]],CS[ellipsoidal,2],AXIS["longitude",east,ORDER[1],ANGLEUNIT["degree",0.0174532925199433,ID["EPSG",9122]]],AXIS["latitude",north,ORDER[2],ANGLEUNIT["degree",0.0174532925199433,ID["EPSG",9122]]]]') # noqa: E501

assert isinstance(crs, self.ccrs.RotatedPole)
assert crs.proj4_params["lon_0"] == 357.5
assert crs.proj4_params["o_lat_p"] == 37.5


class TestDynamicArgs(TestCase):

Expand Down Expand Up @@ -283,20 +291,60 @@ def test_check_crs():

@pytest.mark.parametrize("input", [
"+init=epsg:26911",
])
def test_process_crs(input):
pytest.importorskip("pyproj")
ccrs = pytest.importorskip("cartopy.crs")
crs = process_crs(input)
assert isinstance(crs, ccrs.CRS)


def test_process_crs_pyproj_crs():
pyproj = pytest.importorskip("pyproj")
ccrs = pytest.importorskip("cartopy.crs")
crs = process_crs(pyproj.CRS.from_epsg(4326))
assert isinstance(crs, ccrs.PlateCarree)


def test_process_crs_pyproj_proj():
pyproj = pytest.importorskip("pyproj")
ccrs = pytest.importorskip("cartopy.crs")
crs = process_crs(pyproj.Proj(init='epsg:4326'))
assert isinstance(crs, ccrs.PlateCarree)


@pytest.mark.parametrize("input", [
"4326",
4326,
"epsg:4326",
"EPSG: 4326",
"+init=epsg:4326",
# Created with pyproj.CRS("EPSG:4326").to_wkt()
'GEOGCRS["WGS 84",ENSEMBLE["World Geodetic System 1984 ensemble",MEMBER["World Geodetic System 1984 (Transit)"],MEMBER["World Geodetic System 1984 (G730)"],MEMBER["World Geodetic System 1984 (G873)"],MEMBER["World Geodetic System 1984 (G1150)"],MEMBER["World Geodetic System 1984 (G1674)"],MEMBER["World Geodetic System 1984 (G1762)"],MEMBER["World Geodetic System 1984 (G2139)"],ELLIPSOID["WGS 84",6378137,298.257223563,LENGTHUNIT["metre",1]],ENSEMBLEACCURACY[2.0]],PRIMEM["Greenwich",0,ANGLEUNIT["degree",0.0174532925199433]],CS[ellipsoidal,2],AXIS["geodetic latitude (Lat)",north,ORDER[1],ANGLEUNIT["degree",0.0174532925199433]],AXIS["geodetic longitude (Lon)",east,ORDER[2],ANGLEUNIT["degree",0.0174532925199433]],USAGE[SCOPE["Horizontal component of 3D system."],AREA["World."],BBOX[-90,-180,90,180]],ID["EPSG",4326]]',
], ids=lambda x: str(x)[:20])
def test_process_crs_platecarree(input):
pytest.importorskip("pyproj")
ccrs = pytest.importorskip("cartopy.crs")
crs = process_crs(input)
assert isinstance(crs, ccrs.PlateCarree)


@pytest.mark.parametrize("input", [
"3857",
3857,
"epsg:3857",
"EPSG: 3857",
"+init=epsg:3857",
'PROJCS["WGS 84 / Pseudo-Mercator",GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AUTHORITY["EPSG","4326"]],PROJECTION["Mercator_1SP"],PARAMETER["central_meridian",0],PARAMETER["scale_factor",1],PARAMETER["false_easting",0],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AXIS["Easting",EAST],AXIS["Northing",NORTH],EXTENSION["PROJ4","+proj=merc +a=6378137 +b=6378137 +lat_ts=0 +lon_0=0 +x_0=0 +y_0=0 +k=1 +units=m +nadgrids=@null +wktext +no_defs"],AUTHORITY["EPSG","3857"]]',
# Created with pyproj.CRS("EPSG:3857").to_wkt()
'PROJCRS["WGS 84 / Pseudo-Mercator",BASEGEOGCRS["WGS 84",ENSEMBLE["World Geodetic System 1984 ensemble",MEMBER["World Geodetic System 1984 (Transit)"],MEMBER["World Geodetic System 1984 (G730)"],MEMBER["World Geodetic System 1984 (G873)"],MEMBER["World Geodetic System 1984 (G1150)"],MEMBER["World Geodetic System 1984 (G1674)"],MEMBER["World Geodetic System 1984 (G1762)"],MEMBER["World Geodetic System 1984 (G2139)"],ELLIPSOID["WGS 84",6378137,298.257223563,LENGTHUNIT["metre",1]],ENSEMBLEACCURACY[2.0]],PRIMEM["Greenwich",0,ANGLEUNIT["degree",0.0174532925199433]],ID["EPSG",4326]],CONVERSION["Popular Visualisation Pseudo-Mercator",METHOD["Popular Visualisation Pseudo Mercator",ID["EPSG",1024]],PARAMETER["Latitude of natural origin",0,ANGLEUNIT["degree",0.0174532925199433],ID["EPSG",8801]],PARAMETER["Longitude of natural origin",0,ANGLEUNIT["degree",0.0174532925199433],ID["EPSG",8802]],PARAMETER["False easting",0,LENGTHUNIT["metre",1],ID["EPSG",8806]],PARAMETER["False northing",0,LENGTHUNIT["metre",1],ID["EPSG",8807]]],CS[Cartesian,2],AXIS["easting (X)",east,ORDER[1],LENGTHUNIT["metre",1]],AXIS["northing (Y)",north,ORDER[2],LENGTHUNIT["metre",1]],USAGE[SCOPE["Web mapping and visualisation."],AREA["World between 85.06°S and 85.06°N."],BBOX[-85.06,-180,85.06,180]],ID["EPSG",3857]]'
'PROJCRS["WGS 84 / Pseudo-Mercator",BASEGEOGCRS["WGS 84",ENSEMBLE["World Geodetic System 1984 ensemble",MEMBER["World Geodetic System 1984 (Transit)"],MEMBER["World Geodetic System 1984 (G730)"],MEMBER["World Geodetic System 1984 (G873)"],MEMBER["World Geodetic System 1984 (G1150)"],MEMBER["World Geodetic System 1984 (G1674)"],MEMBER["World Geodetic System 1984 (G1762)"],MEMBER["World Geodetic System 1984 (G2139)"],ELLIPSOID["WGS 84",6378137,298.257223563,LENGTHUNIT["metre",1]],ENSEMBLEACCURACY[2.0]],PRIMEM["Greenwich",0,ANGLEUNIT["degree",0.0174532925199433]],ID["EPSG",4326]],CONVERSION["Popular Visualisation Pseudo-Mercator",METHOD["Popular Visualisation Pseudo Mercator",ID["EPSG",1024]],PARAMETER["Latitude of natural origin",0,ANGLEUNIT["degree",0.0174532925199433],ID["EPSG",8801]],PARAMETER["Longitude of natural origin",0,ANGLEUNIT["degree",0.0174532925199433],ID["EPSG",8802]],PARAMETER["False easting",0,LENGTHUNIT["metre",1],ID["EPSG",8806]],PARAMETER["False northing",0,LENGTHUNIT["metre",1],ID["EPSG",8807]]],CS[Cartesian,2],AXIS["easting (X)",east,ORDER[1],LENGTHUNIT["metre",1]],AXIS["northing (Y)",north,ORDER[2],LENGTHUNIT["metre",1]],USAGE[SCOPE["Web mapping and visualisation."],AREA["World between 85.06°S and 85.06°N."],BBOX[-85.06,-180,85.06,180]],ID["EPSG",3857]]',
], ids=lambda x: str(x)[:20])
def test_process_crs(input):
def test_process_crs_mercator(input):
pytest.importorskip("pyproj")
ccrs = pytest.importorskip("cartopy.crs")
crs = process_crs(input)
assert isinstance(crs, ccrs.Mercator)

assert isinstance(crs, ccrs.CRS)

def test_process_crs_rasterio():
pytest.importorskip("pyproj")
Expand Down
61 changes: 34 additions & 27 deletions hvplot/util.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -62,7 +62,7 @@ def check_crs(crs):
Returns
-------
A valid crs if possible, otherwise None
A valid crs if possible, otherwise None.
"""
import pyproj

Expand All @@ -79,8 +79,11 @@ class Dummy:
out = pyproj.Proj(crs.to_wkt(), preserve_units=True)
elif isinstance(crs, dict) or isinstance(crs, str):
if isinstance(crs, str):
# quick fix for https://github.com/pyproj4/pyproj/issues/345
crs = crs.replace(' ', '').replace('+', ' +')
try:
crs = pyproj.CRS.from_wkt(crs)
except RuntimeError:
# quick fix for https://github.com/pyproj4/pyproj/issues/345
crs = crs.replace(' ', '').replace('+', ' +')
try:
out = pyproj.Proj(crs, preserve_units=True)
except RuntimeError:
Expand Down Expand Up @@ -124,10 +127,10 @@ def proj_to_cartopy(proj):
except ImportError:
has_gdal = False

proj = check_crs(proj)

if proj_is_latlong(proj):
return ccrs.PlateCarree()
input_proj = proj
proj = check_crs(input_proj)
if proj is None:
raise ValueError(f"Invalid proj projection {input_proj!r}")

srs = proj.srs
if has_gdal:
Expand Down Expand Up @@ -167,19 +170,23 @@ def proj_to_cartopy(proj):
except:
pass
if k == 'proj':
if v == 'tmerc':
if v == "longlat":
cl = ccrs.PlateCarree
elif v == 'tmerc':
cl = ccrs.TransverseMercator
kw_proj['approx'] = True
if v == 'lcc':
elif v == 'lcc':
cl = ccrs.LambertConformal
if v == 'merc':
elif v == 'merc':
cl = ccrs.Mercator
if v == 'utm':
elif v == 'utm':
cl = ccrs.UTM
if v == 'stere':
elif v == 'stere':
cl = ccrs.Stereographic
if v == 'ob_tran':
elif v == 'ob_tran':
cl = ccrs.RotatedPole
else:
raise NotImplementedError(f'Unknown projection {v}')
if k in km_proj:
if k == 'zone':
v = int(v)
Expand Down Expand Up @@ -228,7 +235,9 @@ def process_crs(crs):
1. EPSG codes: Defined as string of the form "EPSG: {code}" or an integer
2. proj.4 string: Defined as string of the form "{proj.4 string}"
3. cartopy.crs.CRS instance
4. None defaults to crs.PlateCaree
3. pyproj.Proj or pyproj.CRS instance
4. WKT string: Defined as string of the form "{WKT string}"
5. None defaults to crs.PlateCaree
"""
missing = []
try:
Expand All @@ -248,29 +257,27 @@ def process_crs(crs):

if crs is None:
return ccrs.PlateCarree()
elif isinstance(crs, ccrs.CRS):
return crs
elif isinstance(crs, pyproj.CRS):
crs = crs.to_wkt()

errors = []
if isinstance(crs, str) and crs.lower().startswith('epsg'):
if isinstance(crs, (str, int)): # epsg codes
try:
crs = crs[5:].lstrip().rstrip()
return ccrs.epsg(crs)
crs = pyproj.CRS.from_epsg(crs).to_wkt()
except Exception as e:
errors.append(e)
if isinstance(crs, int):
try:
return ccrs.epsg(crs)
except Exception as e:
crs = str(crs)
errors.append(e)
if isinstance(crs, (str, pyproj.Proj)):
if isinstance(crs, (str, pyproj.Proj)): # proj4/wkt strings
try:
return proj_to_cartopy(crs)
except Exception as e:
errors.append(e)
if isinstance(crs, ccrs.CRS):
return crs

raise ValueError("Projection must be defined as a EPSG code, proj4 string, cartopy CRS or pyproj.Proj.") from Exception(*errors)
raise ValueError(
"Projection must be defined as a EPSG code, proj4 string, "
"WKT string, cartopy CRS, pyproj.Proj, or pyproj.CRS."
) from Exception(*errors)


def is_list_like(obj):
Expand Down

0 comments on commit e4a65f0

Please sign in to comment.