python - 如何使用底图在靠近两极时绘制地球观测卫星的视野?
问题描述
我正在尝试绘制卫星沿其轨道的最大(理论)视野。我正在使用 Basemap,我想在其上绘制沿轨道的不同位置(带散点图),并且我想使用 tissot 方法(或等效方法)添加整个视野。下面的代码可以正常工作,直到纬度达到北纬 75 度,在 300 公里的高度轨道上。超出此代码输出 ValueError:“ValueError:未定义的逆测地线(可能是对映点)”
import matplotlib.pyplot as plt
from mpl_toolkits.basemap import Basemap
import math
earth_radius = 6371000. # m
fig = plt.figure(figsize=(8, 6), edgecolor='w')
m = Basemap(projection='cyl', resolution='l',
llcrnrlat=-90, urcrnrlat=90,
llcrnrlon=-180, urcrnrlon=180)
# draw the coastlines on the empty map
m.drawcoastlines(color='k')
# define the position of the satellite
position = [300000., 75., 0.] # altitude, latitude, longitude
# radius needed by the tissot method
radius = math.degrees(math.acos(earth_radius / (earth_radius + position[0])))
m.tissot(position[2], position[1], radius, 100, facecolor='tab:blue', alpha=0.3)
m.scatter(position[2], position[1], marker='*', c='tab:red')
plt.show()
需要注意的是,该代码在南极(纬度低于 -75)处工作正常。我知道这是一个已知的错误,是否有针对此问题的已知解决方法?谢谢你的帮助!
解决方案
您发现的是 Basemap 的一些限制。现在让我们切换到 Cartopy。工作代码会有所不同,但差别不大。
import matplotlib.pyplot as plt
import cartopy.crs as ccrs
import math
earth_radius = 6371000.
position = [300000., 75., 0.] # altitude (m), lat, long
radius = math.degrees(math.acos(earth_radius / (earth_radius + position[0])))
print(radius) # in subtended degrees??
fig = plt.figure(figsize=(12,8))
img_extent = [-180, 180, -90, 90]
# here, cartopy's' `PlateCarree` is equivalent with Basemap's `cyl` you use
ax = fig.add_subplot(1, 1, 1, projection = ccrs.PlateCarree(), extent = img_extent)
# for demo purposes, ...
# let's take 1 subtended degree = 112 km on earth surface (*** you set the value as needed ***)
ax.tissot(rad_km=radius*112, lons=position[2], lats=position[1], n_samples=64, \
facecolor='red', edgecolor='black', linewidth=0.15, alpha = 0.3)
ax.coastlines(linewidth=0.15)
ax.gridlines(draw_labels=False, linewidth=1, color='blue', alpha=0.3, linestyle='--')
plt.show()
使用上面的代码,结果图是:
现在,如果我们使用正交投影,(用这个替换相关的代码行)
ax = fig.add_subplot(1, 1, 1, projection = ccrs.Orthographic(central_longitude=0.0, central_latitude=60.0))
你得到这个情节:
推荐阅读
- javascript - 如何在jquery中选择当前文本框和其他动态加载的div元素?
- python - Pandas 为 Fiscar QTR/Year 创建日期列
- wordpress - woocommerce - 更改产品的海关属性,并面临 update_post_meta() 递增无唯一“_product_attributes”
- javascript - 从 vue.js 中的嵌套数组中过滤数据
- design-patterns - Task.Factory 和工厂模式
- javascript - Spotify web api 不适用于 ios 设备
- vba - 为什么 querySelectorAll 返回一个 JScriptTypeInfo 对象,然后可以将其视为字符串,而不是 nodeList?
- javascript - vscode 调试控制台自动完成
- r - gganimate:如何使堆积条形图从 x 轴向上平滑增长
- angular - 多个订阅完成/拆卸后如何执行操作