等高线图与Cartopy - 静止卫星视图【图、代码、尝试】

1 投票
1 回答
71 浏览
提问于 2025-04-14 18:08

我有一些数据,这些数据是从球坐标转换过来的,现在是以X,Y坐标的形式。我想把这些数据的等高线绘制在我从Cartopy上获取的地球仪上。但是一直没成功,我已经尝试了几个小时在网上查资料和进行尝试!如果有任何建议,我会非常感激。

我已经修改了代码,加入了一些虚拟数据,这些数据是我想在地球仪上绘制等高线的。

下面的链接里也有一些示例图。

import numpy as np
import matplotlib.pyplot as plt
# import cartopy.feature as cf
import cartopy.crs as ccrs
import cartopy.feature as cf

plt.close('all')


#%% Co-ordindate System
resolution = 201
theta = np.linspace(-90, 90, num=resolution) * np.pi/180
phi = np.linspace(0, 360, num=resolution) * np.pi/180
theta_grid, phi_grid = np.meshgrid(theta,phi)

u_grid = np.sin(theta_grid) * np.cos(phi_grid)
v_grid = np.sin(theta_grid) * np.sin(phi_grid)


#%% Dummy Data - End result called 'AF' is what is important, please ignore how its done. I cant get 'AF' to plot ontop of a globe.
x_num = 10
y_num = 10
num_elements = x_num * y_num
dx = 2
dy = 2
[array_x, array_y] = np.mgrid[0:x_num, 0:y_num]
array_x = (array_x - (x_num- 1) / 2)
array_y = np.flipud((array_y - (y_num - 1) / 2))
array_x = np.transpose(array_x) * dx
array_y = np.flip(np.transpose(array_y)) * dy  
array_x = np.reshape(array_x, (x_num * y_num, -1)).flatten()
array_y = np.reshape(array_y, (x_num * y_num, -1)).flatten()

AF = np.zeros(np.size(u_grid),dtype = np.complex128) #Make sure u_grid and v_grid are the same size
for n in range(num_elements):

    kx = array_x[n] * u_grid.ravel()
    ky = array_y[n] * v_grid.ravel()
    kr = kx+ky
    AF = AF  + np.exp(1j*kr)
    
AF = np.resize(AF, np.shape(u_grid))
AF = AF / np.max(abs(AF))
AF = 20*np.log10(abs(AF))
AF[AF< -50] =-50  


#%% Plots

c = np.linspace(np.min(AF),np.max(AF),50)
)

fig = plt.figure('UV') 
ax = fig.add_subplot(111)
plt.contourf(u_grid,v_grid,AF,c,extend="both", cmap = 'Spectral_r')  #plasma #gnuplot2
plt.colorbar(label ='Directivity (dBi)')
ax.set_xlabel('u = sin(theta)cos(phi)')
ax.set_ylabel('v = sin(theta)sin(phi)')
ax.legend(loc=1)


projection=ccrs.Geostationary(central_longitude=0.0, satellite_height=35785831)

fig = plt.figure('No Overlay') 
ax = plt.axes(projection=projection)
ax.add_feature(cf.COASTLINE)
plt.show()


fig = plt.figure('Attempt to Overlay') 
ax = plt.axes(projection=projection)
plt.contourf(u_grid,v_grid,AF,c,extend="both", cmap = 'Spectral_r', alpha = 0.2)  #plasma #gnuplot2
ax.add_feature(cf.COASTLINE)
plt.show()

这是数据的示例图,以及等高线无法叠加到地球仪上的情况。

在XY坐标系中转换的球坐标数据等高线图:

在XY坐标系中转换的球坐标数据等高线图

我想在上面绘制等高线的地球仪:

我想在上面绘制等高线的地球仪

我尝试了很多次在网上查找这个问题,查看Cartopy网站和stackexchange上的帖子,但还是没有找到解决办法……任何建议我都会非常感激。

1 个回答

1

因为你想把等高线图和地理地图一起展示,所以你需要使用合适的坐标,也就是经度和纬度,单位是度。这和很多Cartopy的例子一样,比如在这个例子中,他们把原来的弧度值转换成了度。在你的代码中,可以这样做:theta_deg = np.rad2deg(theta)phi_deg = np.rad2deg(phi)。另外,你还需要在matplotlib.pyplot.contourf()中设置transform参数。这里有一个使用你数据的简化示例:

import numpy as np
import matplotlib.pyplot as plt
import cartopy.crs as ccrs

#co-ordindate system
resolution = 201
lats = np.linspace(-90, 90, num=resolution)
lons = np.linspace(0, 360, num=resolution)
lons, lats = np.meshgrid(lons, lats)

# load numpy array of data
AF = np.load('Plotted_data.npy')

projection=ccrs.Geostationary(central_longitude=0.0, satellite_height=35785831)

fig = plt.figure()
ax = fig.add_subplot(1, 1, 1, projection=projection)
# show coastlines
ax.coastlines()
# add color filled contours
filled_c = ax.contourf(lons, lats, AF, transform=ccrs.PlateCarree(), extend="both", cmap='Spectral_r', alpha=0.2)

plt.show()

结果如下:

enter image description here

可能结果看起来和你想要的并不完全一样,特别是那个在本初子午线附近的突兀过渡,看起来很难看。但这就是Plotted_data.npy的数据在使用地理坐标时,如何映射到地球上的样子。你可能需要检查一下数据,或者进行某种数据转换。

撰写回答