请参阅上面的问题,只想从虹膜立方体获取Google地图的经纬度...
cube.coord_system()
LambertAzimuthalEqualArea(latitude_of_projection_origin = 54.9,经度_projection_origin = -2.5,false_easting = 0.0,false_northing = 0.0,椭球= GeogCS(semi_major_axis = 6378137.0,semi_minor_axis = 6356752.314140356))
import iris
import numpy as np
首先,获取本机投影中的坐标点
proj_x = cube.coord("projection_x_coordinate").points
proj_y = cube.coord("projection_y_coordinate").points
接下来,制作一对形状相同的二维数组
xx, yy = np.meshgrid(proj_x, proj_y)
然后提取原始投影并将其转换为cartopy
投影:
cs_nat = cube.coord_system()
cs_nat_cart = cs_nat.as_cartopy_projection()
[下一步,创建目标投影,例如标准椭球投影
cs_tgt = iris.coord_systems.GeogCS(iris.analysis.cartography.DEFAULT_SPHERICAL_EARTH_RADIUS)
# Again, convert it to a cartopy projection
cs_tgt_cart = cs_tgt.as_cartopy_projection()
最后,使用cartopy
的transform方法将本机投影中的2D坐标数组转换为目标投影中的坐标。
lons, lats, _ = cs_tgt_cart.transform_points(cs_nat_cart, xx, yy).T
# Note the transpose at the end.
还要注意,上面的函数总是返回z
坐标数组,但在这种情况下为0。
然后您可以使用lons
和lats
在Google Maps或其他应用程序中绘制数据。请记住,这些新坐标是曲线的,因此实际上必须是2D数组。但是,如果要在
iris
(和matplotlib
)中绘制数据,则可以这样进行:
import cartopy.crs as ccrs import matplotlib.pyplot as plt proj_x = cube.coord("projection_x_coordinate").points proj_y = cube.coord("projection_y_coordinate").points cs_nat_cart = cube.coord_system().as_cartopy_projection() fig = plt.figure() ax = fig.add_subplot(111, projection=ccrs.PlateCarree()) ax.pcolormesh(proj_x, proj_y, cube.data, transform=cs_nat_cart) ax.coastlines()
希望这会有所帮助。