python如何绘制极坐标轮廓图contourf
作者：小朱小朱绝不服输
这篇文章主要介绍了python如何绘制极坐标轮廓图contourf问题,具有很好的参考价值,希望对大家有所帮助,如有错误或未考虑完全的地方,望不吝赐教
python绘制极坐标轮廓图contourf
任务：将每个象限的风速，距离，角度绘制成极坐标轮廓图。
使用极坐标画等直线图，可以用极地图的 ax.contour 或 ax.contourf 。
1.变量计算
每个象限的风速，距离就不再说怎么画了，这里说下角度的计算。
两个经纬度之间的夹角（与正北方向的夹角）：
# 点1到点2方向沿逆时针方向转到正北方向的夹角
def LatLng2Degree(LatZero,LngZero,Lat,Lng):
"""
Args:
point p1(latA, lonA)
point p2(latB, lonB)
Returns:
bearing between the two GPS points,
default: the basis of heading direction is north
"""
radLatA = math.radians(LatZero)
radLonA = math.radians(LngZero)
radLatB = math.radians(Lat)
radLonB = math.radians(Lng)
dLon = radLonB - radLonA
y = math.sin(dLon) * math.cos(radLatB)
x = math.cos(radLatA) * math.sin(radLatB) - math.sin(radLatA) * math.cos(radLatB) * math.cos(dLon)
brng = math.degrees(math.atan2(y, x))
brng = (brng + 360) % 360
return brng四个象限的距离，风速，角度一一对应。
2.极坐标绘制ax.contourf
import numpy as np import matplotlib.pyplot as plt #-- Generate Data ----------------------------------------- # Using linspace so that the endpoint of 360 is included... azimuths = np.radians(np.linspace(0, 360, 20)) zeniths = np.arange(0, 70, 10) r, theta = np.meshgrid(zeniths, azimuths) values = np.random.random((azimuths.size, zeniths.size)) #-- Plot... ------------------------------------------------ fig, ax = plt.subplots(subplot_kw=dict(projection='polar')) ax.contourf(theta, r, values) plt.show()
首先，创建极坐标轴，绘制轮廓图。
fig, ax = subplots(subplot_kw=dict(projection='polar')) cax = ax.contourf(theta, r, values, nlevels)
参数分别为theta:角度，r:半径，values:轮廓的实际值，nlevels：要绘制的轮廓图的层数。
这里theta， r， values 都需要是二维数组，需要通过列表转化为列。
相当于需要先构造网格矩阵，对于坐标轴，一般使用
r, theta = np.meshgrid(zeniths, azimuths)
np.meshgrid函数去构造。
一般情况下，theta， r， values的对应关系是：
可以通过reshape将值转化为二维数组
这里可以使用：
theta = np.radians(azimuths) zeniths = np.array(zeniths) values = np.array(values) values = values.reshape(len(azimuths), len(zeniths)) r, theta = np.meshgrid(zeniths, np.radians(azimuths))
但当theta， r， values一一对应时，需要插值实现。
这里记得将角度转化为弧度。
r = np.array(r) v = np.array(v) angle = np.array(angle) angle = np.radians(angle) angle1, r1 = np.meshgrid(angle, r) # r和angle构造网格矩阵 v_new = griddata((angle, r), v, (angle1, r1), method='linear') cax = ax.contourf(angle1, r1, v_new, 20, cmap='jet')
这里python插值可以参考：python插值（scipy.interpolate模块的griddata和Rbf）
可以解释一下使用Rbf插值会报错的原因，是因为有的方位角相同，比如多行是0度，这样在使用Rbf插值的时候，要求矩阵可逆，则会出现报错的情况。
插值完以后，这里还有一个点，就是把0度调整为正北方向。
ax.set_theta_zero_location("N")
ax.set_theta_direction(-1)运行结果：
可以看到四个象限的风速就用极坐标轮廓图描述出来了。
总结
以上为个人经验，希望能给大家一个参考，也希望大家多多支持脚本之家。