Как преобразовать в картографическую проекцию из географического, как в базовой карте? - PullRequest
1 голос
/ 07 июня 2019

Я хочу преобразовать lon / lat (в градусах) в координаты проекции карты x / y (в метрах), но с использованием cartopy + pyplot вместо базовой карты.

скажем, что это код базовой карты:

>>> from mpl_toolkits.basemap import Basemap
>>> import numpy as np
>>> import matplotlib.pyplot as plt
>>> # read in topo data (on a regular lat/lon grid)
>>> etopo = np.loadtxt('etopo20data.gz')
>>> lons  = np.loadtxt('etopo20lons.gz')
>>> lats  = np.loadtxt('etopo20lats.gz')
>>> # create Basemap instance for Robinson projection.
>>> m = Basemap(projection='robin',lon_0=0.5*(lons[0]+lons[-1]))
>>> # compute map projection coordinates for lat/lon grid.
>>> x, y = m(*np.meshgrid(lons,lats))

Я хочу эмулировать аналогичную функциональность в cartopy, как я могу это сделать?

1 Ответ

1 голос
/ 07 июня 2019

Шаги для достижения точек сетки сетки, подходящих для построения с помощью Cartopy, отличаются и, насколько я знаю, более сложны.

Вот рабочий код, использующий Cartopy:

import matplotlib.pyplot as plt
import cartopy
import cartopy.crs as ccrs
import numpy as np

# create arrays of values for long and lat
lons = np.linspace(0,160,10)
lats = np.linspace(0,70,5)

# create meshgrid of points
x, y = np.meshgrid(lons, lats)

# select a CRS/projection to tranform/plot points for demo
use_proj = ccrs.Robinson();

# transform all the meshgrid points arrays ..
# .. from geodetic long/lat to Robinson x/y/z
out_xyz = use_proj.transform_points(ccrs.Geodetic(), x, y)
# out_xyz.shape -> (5, 10, 3)

# separate x_array, y_array from the result(x,y,z) above
x_array = out_xyz[:,:,0]
y_array = out_xyz[:,:,1]

# setup fig/axis and plot the meshgrid of points
fig = plt.figure()
ax = fig.add_axes([0, 0, 1, 1], projection=use_proj)
ax.add_feature(cartopy.feature.LAND, facecolor='lightgray')
ax.scatter(x_array, y_array, s=25, c="r", zorder=10)
ax.set_global()
plt.show()

Выходной график будет:

enter image description here

...