文章/答案/技术大牛

发布

Basemap工具函数

共 4 篇文章

Basemap工具函数(2)

Basemap工具函数(1)

Basemap工具函数(4)

Basemap工具函数(3)

清单首页Basemap工具函数文章详情

清单「Basemap工具函数」 04/04

Basemap工具函数(3)

bugsuse

rotate_vector

给定向量场的东西和南北方向分量以及经纬度点，然后对向量进行旋转，使向量场在地图投影上以适当的方向显示。

一些函数（比如 barbs，quiver，streamplot）使用的是向量数据，要求向量分量是地图坐标系（比如 u 是左右方向，v 是上下方向）。如果可用数据是地理学坐标系的（比如东西方向，南北方向），这些坐标必须进行转换，否则所绘制的向量方向会很怪异。这就是 rotate_vertor 方法的目的。

transform_scalar 方法也能完成同样的工作，但同时也会改变网格大小（插值）。

rotate_vector(uin, vin, lons, lats, returnxy=False)

uin 和 vin 是输入数据矩阵。方向是地理学坐标中的方向，u 是东西方向分量， v 是南北方向分量
lons 和 lats 是和 uin 及 vin 矩阵相对应的位置数组（二维）。也是地理学坐标
returnxy 会使方法返回重新投影到地图坐标的 lons 和 lats 矩阵

from mpl_toolkits.basemap import Basemap
import matplotlib.pyplot as plt
import numpy as np


fig=plt.figure(figsize=(9, 3))


map = Basemap(projection='sinu', 
              lat_0=0, lon_0=0)

lons = np.linspace(-180, 180, 10)
lats = np.linspace(-90, 90, 10)

lons, lats = np.meshgrid(lons, lats)

v10 = np.ones((lons.shape)) * 15
u10 = np.zeros((lons.shape))

u10_rot, v10_rot, x, y = map.rotate_vector(u10, v10, lons, lats, returnxy=True)

ax = fig.add_subplot(121)
ax.set_title('Without rotation')

map.drawmapboundary(fill_color='aqua')
map.fillcontinents(color='#cc9955', lake_color='aqua', zorder = 0)
map.drawcoastlines(color = '0.15')


map.barbs(x, y, u10, v10, 
    pivot='middle', barbcolor='#333333')

ax = fig.add_subplot(122)
ax.set_title('Rotated vectors')

map.drawmapboundary(fill_color='9999FF')
map.fillcontinents(color='#ddaa66', lake_color='9999FF', zorder = 0)
map.drawcoastlines(color = '0.15')

map.barbs(x, y, u10_rot, v10_rot, 
    pivot='middle', barbcolor='#ff7777')

plt.show()

lons 和 lats 被创建为覆盖全球的等间距网格
v10 和 u10 为风分量（v10 = 10， u10 = 0）
对创建的矩阵进行旋转并且计算在地图投影中的位置（returnxy = True）
使用 add_subplot 方法创建子图并绘制两个图 1）第一张图没有进行旋转。尽管假设风场为南北方向，但 barbs 显示出为由低到顶，而不是由南到北 2）第二张图展示的是旋转后的结果

set_axes_limits

此方法通常被用于改变 matplotlib axes以改变投影的形状。通常在内部进行调用。

set_axes_limits(ax=None)

ax 是需要改变范围的 axis 实例

大部分方法，比如 drawcountries, drawrivers, readshapefile等，都会在最后调用此方法。因此很难使用例子展示其是如何起作用的。

from mpl_toolkits.basemap import Basemap
import matplotlib.pyplot as plt
fig = plt.figure()

map = Basemap(projection='ortho', 
              lat_0=0, lon_0=0)

map.drawlsmask(land_color = "#ddaa66", 
               ocean_color="#7777ff",
               resolution = 'c')

x1, y1 = map(-60, -30)
x2, y2 = map(0, 0)
x3, y3 = map(45, 45)

plt.plot([x1, x2, x3], [y1, y2, y3], color='k', linestyle='-', linewidth=2)

ax1 = fig.add_axes([0.1, 0., 0.15, 0.15])
ax1.set_xticks([])
ax1.set_yticks([])

ax1.plot([x1, x2, x3], [y1, y2, y3], color='k', linestyle='-', linewidth=2)

map.set_axes_limits(ax=ax1)

ax2 = fig.add_axes([0.3, 0., 0.15, 0.15])
ax2.set_xticks([])
ax2.set_yticks([])

ax2.plot([x1, x2, x3], [y1, y2, y3], color='k', linestyle='-', linewidth=2)

ax3 = fig.add_axes([0.5, 0., 0.15, 0.15])
ax3.set_xticks([])
ax3.set_yticks([])

map.plot([x1, x2, x3], [y1, y2, y3], color='k', linestyle='-', linewidth=2, ax=ax3)

plt.show()

shiftdata

添加经度值从而正确适应地图。改变数据数组也能适应新的地图。有时候，精度数据给定的间隔范围是从 -180 到 180，或是从 0 到 360。要正确绘图，数据必须要进行转换。

shiftdata(lonsin, datain=None, lon_0=None)

lonsin 是原始经度数据。必须是一维或二维数组，但总是规则的 lat/lon 网格
datain 是数据数组。如果 longitudes 必须要分割（比如，-180 变为 360），那么数据的顺序必要要改变。如果给定了数据，则执行此操作。
lon_0 是地图的原点.。经度为[lon_0-180, lon_0 + 180]区间

from mpl_toolkits.basemap import Basemap
import matplotlib.pyplot as plt
import numpy as np

map = Basemap(projection='sinu', 
              lat_0=0, lon_0=0)

lons = np.arange(30, 390, 30)
lats = np.arange(0, 100, 10)

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

print data

print lons

lons, data = map.shiftdata(lons, datain = data, lon_0=0)

print lons

llons, llats = np.meshgrid(lons, lats)

x, y = map(llons, llats)

map.contourf(x, y, data)

map.drawcoastlines()
plt.show()

坐标和数据不是严格意义上的经纬度和有意义的数据。坐标是简单的范围，而数据是 x 和 y 坐标的和。因此，左下角值最小，右上角值最大。

注意：经度以 30 为起始点。

使用 shiftdata 转换经度和数据
使用 meshgrid 方法转换lons 和 lats 为二维数组，并且使用 basemap 实例对地图投影进行重新投影
绘制填充等值线图
最终会出现一块白色区域，可使用 addcyclic 方法解决

应用此方法到数据后的结果

未使用此方法的结果。注意：蓝色（值为0）区域起始点不是经度=0

未使用此方法的结果。注意结果和预期不同，因为经度在正确区间外

shiftgrid

此函数和 shiftdata 方法类似。向西或向东以东所有经度和数据。

basemap.shiftgrid(lon0, datain, lonsin, start=True, cyclic=360.0)

注意这不是 basemap实例方法，而是basemap模块中的函数
lon_0 是最终格点的起始或终止经度。必须在输入经度区间内，因此有时要设置为起始点，有时作为终止点（设置 start 参数）
datain 是输入数据数组，将根据经度进行重新排序
lonsin 是要进行转换的输入经度，是一维数组，但不一定是均匀间隔
start 参数，默认是True。如果 lon0是经度的初始点或是输出数组的最后一个点，要设置 start
cyclic 设置经度值，经度在[0 lon0]区间
函数返回被排序后的数据和转换后的经度

注意：

和 shiftdata 的主要差异是：投影不需要是 cylindrical，因为这不是 basemap 实例方法，而且经度也不需要是均匀间隔的。

from mpl_toolkits.basemap import Basemap
from mpl_toolkits.basemap import shiftgrid
import matplotlib.pyplot as plt
import numpy as np

map = Basemap(projection='sinu', 
              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]

print data
print lons

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

print data
print lons

llons, llats = np.meshgrid(lons, lats)

x, y = map(llons, llats)

map.contourf(x, y, data)
map.drawcoastlines()
plt.show()

坐标和数据不是严格意义上的经纬度和有意义的数据。坐标是简单的范围，而数据是 x 和 y 坐标的和。因此，左下角值最小，右上角值最大。

注意：经度以 30 为起始点。

第二个元素由60改为70，以说明经度可以不是均匀间隔

使用 shiftgrid 转换经度和数据 1）lon0 设置为 180，start = False，因此最终经度为 180 2）不可能设置 lon0 为 -180，因为这在经度区间之外，因此必须要使用上述方法。
使用 meshgrid 方法转换lons 和 lats 为二维数组，并且使用 basemap 实例对地图投影进行重新投影
绘制填充等值线图