将Basemap转换为Cartopy，是否有Basemap的shiftgrid()等等效函数？

Question

我正在将使用 matplotlib 的工具包 Basemap 的应用程序转换为使用 Cartopy 以准备从 Python 2 迁移到 Python 3。 我在 Cartopy 中为 Basemap 的“addcyclic()”和“maskoceans()”找到了类似的功能， 但是，对于 Basemap 的 shiftgrid() 函数，我在 numpy 或 Cartopy 中找不到类似的东西。

这是使用底图的代码： '''

    import matplotlib.pyplot as plt
    from mpl_toolkits.basemap import Basemap
    import cartopy
    import cartopy.crs as ccrs
    import cartopy.feature as cfeature
    import numpy as np
    from mpl_toolkits.basemap import shiftgrid

    bmap = Basemap(projection='ortho', lat_0=0, lon_0=0)
    lons = np.arange(30, 410, 30)
    lons[1] = 70
    lats = np.arange(0, 100, 10)

    data = np.indices((lats.shape[0], lons.shape[0]))
    data = data[0] + data[1]

    data, lons = shiftgrid(180., data, lons, start=False)

    llons, llats = np.meshgrid(lons, lats)
    x, y = bmap(llons, llats)
    bmap.contourf(x, y, data)
    bmap.drawcoastlines()

'''

初始数据： 数据 '''

    [[ 0  1  2  3  4  5  6  7  8  9 10 11 12]
     [ 1  2  3  4  5  6  7  8  9 10 11 12 13]
     [ 2  3  4  5  6  7  8  9 10 11 12 13 14]
     [ 3  4  5  6  7  8  9 10 11 12 13 14 15]
     [ 4  5  6  7  8  9 10 11 12 13 14 15 16]
     [ 5  6  7  8  9 10 11 12 13 14 15 16 17]
     [ 6  7  8  9 10 11 12 13 14 15 16 17 18]
     [ 7  8  9 10 11 12 13 14 15 16 17 18 19]
     [ 8  9 10 11 12 13 14 15 16 17 18 19 20]
     [ 9 10 11 12 13 14 15 16 17 18 19 20 21]]

     lons

     [ 30  70  90 120 150 180 210 240 270 300 330 360 390]

     After the 'data, lons = shiftgrid(180., data, lons, start=False)':
     data

     [[ 5  6  7  8  9 10 11 12  1  2  3  4  5]
      [ 6  7  8  9 10 11 12 13  2  3  4  5  6]
      [ 7  8  9 10 11 12 13 14  3  4  5  6  7]
      [ 8  9 10 11 12 13 14 15  4  5  6  7  8]
      [ 9 10 11 12 13 14 15 16  5  6  7  8  9]
      [10 11 12 13 14 15 16 17  6  7  8  9 10]
      [11 12 13 14 15 16 17 18  7  8  9 10 11]
      [12 13 14 15 16 17 18 19  8  9 10 11 12]
      [13 14 15 16 17 18 19 20  9 10 11 12 13]
      [14 15 16 17 18 19 20 21 10 11 12 13 14]]

      lons

      [-180 -150 -120  -90  -60  -30    0   30   70   90  120  150  180]

''' 我尝试了以下 cartopy 代码来重新创建 Basemap shiftgrid 所做的事情。 这是 Cartopy 代码，当我一次尝试时，有些东西被注释掉了： '''

    DATA_CRS = ccrs.PlateCarree()
    lons = np.arange(30, 410, 30)
    lons[1] = 70
    lats = np.arange(0, 100, 10)

    data = np.indices((lats.shape[0], lons.shape[0]))
    data = data[0] + data[1]
    # data2 = np.roll(data, -5)
    # lons2 = np.mod(lons2 - 180.0, 360.0) - 180.0
    cm_lon = 0
    #llons, llats = np.meshgrid(lons2, lats)
    llons, llats = np.meshgrid(lons, lats)
    PROJECTION = ccrs.Orthographic(central_longitude=cm_lon)
    fig1 = plt.figure(num=1, figsize=(11, 8.5), dpi=150)
    ax = plt.axes(projection=PROJECTION)
    ax.add_feature(cfeature.COASTLINE, linewidths=0.7)
    ax.add_feature(cfeature.BORDERS, edgecolor='black', linewidths=0.7)
    ax.contourf(llons, llats, data, transform=ccrs.PlateCarree())

'''

数据和经度为原始数据，我只是在投影中使用了“central_longitude”。 底图图像显示整个地球，但 Cartopy 图像仅显示赤道以上。 除了最右侧之外，数据的颜色似乎相似，所以我担心 Cartopy 中的数据与 Basemap 中的映射不同。

所以，问题是......有什么相当于 Basemap 的 shiftgrid() 或者我需要找出类似于 Basemap 的 shiftgrid() 的东西还是只在投影中使用“central_longitude”？ 我似乎无法粘贴 .png 文件。 非常感谢任何帮助。 我在网上搜索了相同的功能，但没有找到 shiftgrid() 的功能。 谢谢。

Answer 1

我不知道有任何 shiftgrid 等价物。可能值得在CartoPy issue tracker 请求此类功能上打开一个问题。提及一个可靠的用例来帮助驱动功能将有助于这样做。

Answer 2

这一定是最不优雅的解决方案，但我一直在做的一些 Basemap 的有用功能（还没有？）在 cartopy 中，只是从 Basemap 的源代码中复制功能定义。它工作正常。比如shiftgrid：

def shiftgrid(lon0,datain,lonsin,start=True,cyclic=360.0):
"""
Shift global lat/lon grid east or west.
.. tabularcolumns:: |l|L|
==============   ====================================================
Arguments        Description
==============   ====================================================
lon0             starting longitude for shifted grid
                 (ending longitude if start=False). lon0 must be on
                 input grid (within the range of lonsin).
datain           original data with longitude the right-most
                 dimension.
lonsin           original longitudes.
==============   ====================================================
.. tabularcolumns:: |l|L|
==============   ====================================================
Keywords         Description
==============   ====================================================
start            if True, lon0 represents the starting longitude
                 of the new grid. if False, lon0 is the ending
                 longitude. Default True.
cyclic           width of periodic domain (default 360)
==============   ====================================================
returns ``dataout,lonsout`` (data and longitudes on shifted grid).
"""
if np.fabs(lonsin[-1]-lonsin[0]-cyclic) > 1.e-4:
    # Use all data instead of raise ValueError, 'cyclic point not included'
    start_idx = 0
else:
    # If cyclic, remove the duplicate point
    start_idx = 1
if lon0 < lonsin[0] or lon0 > lonsin[-1]:
    raise ValueError('lon0 outside of range of lonsin')
i0 = np.argmin(np.fabs(lonsin-lon0))
i0_shift = len(lonsin)-i0
if ma.isMA(datain):
    dataout  = ma.zeros(datain.shape,datain.dtype)
else:
    dataout  = np.zeros(datain.shape,datain.dtype)
if ma.isMA(lonsin):
    lonsout = ma.zeros(lonsin.shape,lonsin.dtype)
else:
    lonsout = np.zeros(lonsin.shape,lonsin.dtype)
if start:
    lonsout[0:i0_shift] = lonsin[i0:]
else:
    lonsout[0:i0_shift] = lonsin[i0:]-cyclic
dataout[...,0:i0_shift] = datain[...,i0:]
if start:
    lonsout[i0_shift:] = lonsin[start_idx:i0+start_idx]+cyclic
else:
    lonsout[i0_shift:] = lonsin[start_idx:i0+start_idx]
dataout[...,i0_shift:] = datain[...,start_idx:i0+start_idx]
return dataout,lonsout

Answer 3

我找到了底图的 shiftgrid 功能 here 您可以将其与 cartopy 一起作为单独的函数调用。

import numpy as np
import numpy.ma as ma
def shiftgrid(lon0,datain,lonsin,start=True,cyclic=360.0):
    """
    Shift global lat/lon grid east or west.
    .. tabularcolumns:: |l|L|
    ==============   ====================================================
    Arguments        Description
    ==============   ====================================================
    lon0             starting longitude for shifted grid
                     (ending longitude if start=False). lon0 must be on
                     input grid (within the range of lonsin).
    datain           original data with longitude the right-most
                     dimension.
    lonsin           original longitudes.
    ==============   ====================================================
    .. tabularcolumns:: |l|L|
    ==============   ====================================================
    Keywords         Description
    ==============   ====================================================
    start            if True, lon0 represents the starting longitude
                     of the new grid. if False, lon0 is the ending
                     longitude. Default True.
    cyclic           width of periodic domain (default 360)
    ==============   ====================================================
    returns ``dataout,lonsout`` (data and longitudes on shifted grid).
    """
    if np.fabs(lonsin[-1]-lonsin[0]-cyclic) > 1.e-4:
        # Use all data instead of raise ValueError, 'cyclic point not included'
        start_idx = 0
    else:
        # If cyclic, remove the duplicate point
        start_idx = 1
    if lon0 < lonsin[0] or lon0 > lonsin[-1]:
        raise ValueError('lon0 outside of range of lonsin')
    i0 = np.argmin(np.fabs(lonsin-lon0))
    i0_shift = len(lonsin)-i0
    if ma.isMA(datain):
        dataout  = ma.zeros(datain.shape,datain.dtype)
    else:
        dataout  = np.zeros(datain.shape,datain.dtype)
    if ma.isMA(lonsin):
        lonsout = ma.zeros(lonsin.shape,lonsin.dtype)
    else:
        lonsout = np.zeros(lonsin.shape,lonsin.dtype)
    if start:
        lonsout[0:i0_shift] = lonsin[i0:]
    else:
        lonsout[0:i0_shift] = lonsin[i0:]-cyclic
    dataout[...,0:i0_shift] = datain[...,i0:]
    if start:
        lonsout[i0_shift:] = lonsin[start_idx:i0+start_idx]+cyclic
    else:
        lonsout[i0_shift:] = lonsin[start_idx:i0+start_idx]
    dataout[...,i0_shift:] = datain[...,start_idx:i0+start_idx]
    return dataout,lonsout