四、绘制上海到芝加哥大圆航线
不知道坐过国际航班的你是否也产生过这样的疑问：为啥国内出发的去美国的飞机不按照地图上两点之间直线最短拉一条线段从太平洋上飞而要越过西伯利亚穿过白令海峡在北极圈兜一圈再缓缓绕过加拿大最后抵达美国？
并不是害怕飞机掉海里...(虽然迫降死亡率最高确实在海面，在陆地上还有微乎其微的生还几率）
这得从地图的绘制来讲。
首先我们要明确... 地球是一个近似球体...
我们日常所见到的世界地图是基于墨卡托投影的：于是粘一段维基百科，麦卡托投影法(Mercator projection)，又称麦卡托投影法、正轴等角圆柱投影，是一种等角的圆柱形地图投影法。在以此投影法绘制的地图上，经线纬线于任何位置皆垂直相交，世界地图可以绘制在一个长方形上。由于可显示任两点间的正确方位，航海用途的海图、航路图大都以此方式绘制。在该投影中线型比例尺在图中任意一点周围都保持不变，从而可以保持大陆轮廓投影后的角度和形状不变（即等角投影)）；但麦卡托投影会使面积产生变形，极点的比例甚至达到了无穷大。
所以我们日常看见的世界地图南极被拉成了一长条，战斗民族国家显得幅员辽阔（虽然国土面积仍然是世界第一但没有地图上看起来的那么大）
但其实遵从我们高中所学的（还是大学所学的空间解析几何？）地球上两点之间最短的距离并不是在地图上画一条连线，而是通过连接球上两点作出过圆心的大弧~~~
贴一个计算两点中点经纬度的公式：
不高兴了，居然不支持Latex,还得我自己在tex弄完保存成图片传，简书真low，Latex中文支持的包也不好使了，屋漏偏逢连夜雨
哼，求得中点经纬线后喜闻乐见的调用basemap画出上海到芝加哥的大圆航线：
#把地球看做一个球体，通过地面上任意两点和地心做一平面，平面与地球表面橡胶看到的圆周就是大圆。
#两点之间的大圆劣弧线是两点在地面上的最短距离。沿着这一段大圆弧线航行时的航线称为大圆航线。
#由于大圆航线是两点之间的最短航线，故有时称为最经济航线。
chilat = 41.52; chilon = -87.37
# 芝加哥经纬度
shlat = 31.22; shlon = 121.48
#上海经纬度
m = Basemap(llcrnrlon=-140.,llcrnrlat=10.,urcrnrlon=140.,urcrnrlat=70.,\
rsphere=(6378137.00,6356752.3142),\
resolution='l',projection='merc',\
lat_0=(chilat+shlat)/2.,lon_0=(chilon+shlon)/2.,lat_ts=0.)
#设置图像边界四点经纬度坐标，中心点经纬度参数坐标， 地球半径，分辨率设置为低
# 麦卡托式地图投影法设置 projection = 'merc'
# draw great circle route between NY and London
m.drawgreatcircle(shlon,shlat,chilon,chilat,linewidth=2,color='b')
#画出上海芝加哥两地间的大圆航线，设置线宽及颜色
m.drawcoastlines()
#海岸线
m.fillcontinents()
#填充大陆
m.drawparallels(np.arange(10,90,20),labels=[1,1,0,1])
# 画出纬线， 在北纬10度到90度区间内以20度为单位， 纬度标记在图形左右和下测
m.drawmeridians(np.arange(-180,180,30),labels=[1,1,0,1])
# # 画出经线， 从西经180度到东经180度区间内以30度为单位， 经度标记在图形左右和下测
ax.set_title('Great Circle from Chicago to Shanghai')
#图像命名为芝加哥到上海的大圆航线
plt.show()
芝加哥-上海大圆航线.png
这里其实有一个疑惑，我设置的投影中心点没有按照上文提到的球面中心点计算（因为算了以后生成的地图视角特别奇怪...主要显现的是欧洲大陆和南亚次大陆，而且是我们常见世界地图反转的...可是我的上下左右经纬度设置没错啊...）
奇怪的视角.png
五、其他看起来很炫酷的地图
偷下懒如果没有新导入的库，import部分就不写了。
5.1 微信开始界面的蓝色弹珠（Blue Marble）
看起来很高大上其实就加了一行代码...
plt.figure(figsize=(8, 8))
#设新图
m = Basemap(projection='ortho', resolution=None, lat_0=30, lon_0=120)
#正射投影，投影原点设在了上海周边
m.bluemarble(scale=0.5);
#图像原始分辨率是5400*2700，设置scale = 0.5以后分辨率为2700*1350,如此作图
#迅速不少也不那么占用内存了
蓝色弹珠-上海投影视角.png
5.2 据说原模型是用来测量储水量的EPOTO1图
查了下好像很高大上的样子...这个模型很能反映地表信息，主要用来看冰层和基岩的。输入芝加哥经纬度标记下芝加哥的位置
fig = plt.figure(figsize=(8, 8))
m = Basemap(projection='lcc', resolution=None,
width=8E6, height=8E6,
lat_0=45, lon_0=-100,)
#这里投影模式llc是圆锥投影
m.etopo(scale=0.5, alpha=0.5)
#其实也就加了这一句...scale原理同上文bluemarble中的scale，alpha是透明度图层
# Map (long, lat) to (x, y) for plotting
x, y = m(-87.37, 41.52)
plt.plot(x, y, 'ok', markersize=5)
#这里不是打印ok啊...是设置标记点为黑色圆点
plt.text(x, y, ' Chicago', fontsize=12);
#其实我在想瞄成红点加一个launched字符是不是更有冲击力（不带这么诅咒自己的）
OPOTO1_Chicago.png
五、绘制即时晨昏线
晨昏圈，又称晨昏线，或是曙暮光区是一条虚拟的线，它在行星的表面画出了白天和黑夜的交界线（也称为灰线）。晨昏圈由晨线和昏线组成，晨线和昏线各是一个半圆弧，晨线的东边是昼半球，昏线的西边是夜半球。
其实这个图没什么好说的，稍微值得补充的部分都放在代码的注释部分，直接上代码：
import numpy as np
from mpl_toolkits.basemap import Basemap
import matplotlib.pyplot as plt
from datetime import datetime
#python标准库 日期与时间
from dateutil import tz
#这是一个第三方库，拓展了标准时间库，tz = time zone
fig=plt.figure()
ax=fig.add_axes([0.1,0.1,0.8,0.8])
# 新图对象
map = Basemap(projection='mill',lon_0=180)
#米勒投影法，此投影与墨卡托投影类似，只是极点区域的面积变形不如后者大。向极点靠近时，两
#条纬线的间距比墨卡托投影的小。这样就降低了面积变形程度，但这会导致局部形状和方向发生变形。
map.drawcoastlines()
#画海岸线
map.drawparallels(np.arange(-90,90,30),labels=[1,0,0,0])
map.drawmeridians(np.arange(map.lonmin,map.lonmax+30,60),labels=[0,0,0,1])
#numpy的作用只是提供了数组
map.drawmapboundary(fill_color='aqua')
#边界aqua(水）色
map.fillcontinents(color='coral',lake_color='aqua')
#大陆珊瑚色，湖泊水色
date = datetime.now()
#这个是local时间，我所在的芝加哥地区是(CST)central standard time
#如果想取格林尼治时间（Coordinated Universal Time），以上这句改为
#date = datetime.utcnow()即可
CS=map.nightshade(date)
#basemap中的nightshade函数用作将夜晚区域覆盖作阴影，具体参数不扩展
#有兴趣了解搜索basemap.nightshade即可
plt.title('Day/Night Map for %s (CST)' % date.strftime("%d %b %Y %H:%M:%S"))
#设置图名， datetime格式日/月/年/时/分/秒
plt.show()
即时晨昏线.png
把这张图放在最后也是为了秀一下夜已经深了。
如此python关于basemap的地图可视化整理告一段落，不排除以后会添加新应用实例或者补充诸如如何释放内存的解决方案（画了几个高级图遇到了MemoryError,gc.collect()也不好使...容我再研究下）。晚安！
网友评论