快捷導(dǎo)航

python如何繪制登陸時(shí)的衛(wèi)星云圖(TBB)

更新時(shí)間：2023年08月16日 10:40:47 作者：小朱小朱絕不服輸

這篇文章主要介紹了python如何繪制登陸時(shí)的衛(wèi)星云圖(TBB),具有很好的參考價(jià)值,希望對(duì)大家有所幫助,如有錯(cuò)誤或未考慮完全的地方,望不吝賜教

1.讀取hdf文件

詳見上篇文章，這里不再贅述，再貼一下代碼。

主要是從全圓盤標(biāo)對(duì)稱文件經(jīng)緯度對(duì)照表獲取經(jīng)緯度和從hdf文件中獲取亮溫值。

# 從全圓盤標(biāo)對(duì)稱文件經(jīng)緯度對(duì)照表獲取經(jīng)緯度
def getLatLonFromdat():
    lonlatfile = 'F:/Satellite_Imagery/Code/NOM_ITG_2288_2288(0E0N)_LE.dat'
    with open(lonlatfile, 'rb') as f:
        #lon_fy = np.fromfile(f, count=2288 * 2288, dtype='float32') + 79  # 先存經(jīng)度，根據(jù)衛(wèi)星的不同加上對(duì)應(yīng)的經(jīng)度值
        #lat_fy = np.fromfile(f, count=2288 * 2288, dtype='float32')  # 再存緯度
        data = np.fromfile(f, dtype='float32')
        data = data.reshape([2288, 2288, 2], order='F')
    #lon = lon_fy.reshape([2288, 2288], order='F')
    #lat = lat_fy.reshape([2288, 2288], order='F')
    lon = data[:, :, 0] + 104.5
    lat = data[:, :, 1]
    return lon, lat

# 從hdf文件中獲取亮溫值
def getTBBFromhdf():
    hdfFile = h5py.File('F:/IR/Satellite_Imagery/IR_data/利奇馬/FY2G_FDI_ALL_NOM_20190811_1200.hdf', 'r')
    db1 = hdfFile['/CALChannelIR1']
    hw1 = hdfFile['/NOMChannelIR1']
    # db2 = hdfFile['/CALChannelIR2']
    # hw2 = hdfFile['/NOMChannelIR2']
    # db3 = hdfFile['/CALChannelIR3']
    # hw3 = hdfFile['/NOMChannelIR3']
    # db4 = hdfFile['/CALChannelIR4']
    # hw4 = hdfFile['/NOMChannelIR4']
    infoh = hdfFile['/NomFileInfo']
    # 查看衛(wèi)星的經(jīng)緯度
    lat_hdf = infoh[0][3]
    lon_hdf = infoh[0][4]
    # print(lat_hdf)
    # print(lon_hdf)
    hw = hw1[()]
    db = db1[()]  # 獲取定標(biāo)表的值
    tb = np.zeros(shape=(2288, 2288))  # 2288*2288的圖像每個(gè)具體的亮溫值
    for i in range(2288):
        for j in range(2288):
            if hw[i][j] == 65535 or hw[i][j] == 65534:
                tb[i][j] = 0
            else:
                a = hw[i][j]
                tb[i][j] = db[0][a]
    tb = tb.T
    return tb

2. 畫圖

fig = plt.figure(figsize=(8, 6))
m = Basemap(projection='cyl', llcrnrlat=10, llcrnrlon=110, urcrnrlat=40,  urcrnrlon=140)  # 使用Basemap繪制地圖，這里可以讀取對(duì)應(yīng)的地圖shp文件。
m.drawcoastlines(color='black')
m.drawstates(color='black')
m.drawcountries(color='black')
x, y = m(lon, lat)  # 將lats / lons轉(zhuǎn)換為地圖投影坐標(biāo)
# 繪制輪廓圖
# 這里data就是計(jì)算的亮溫值，x,y就是經(jīng)緯度投影的坐標(biāo)
cf = m.contourf(x, y, data, levels=np.linspace(180, 301, 400), cmap='jet')
cbar = m.colorbar(cf, location='right', size='5%', pad='2%')
font = {'family': 'serif',
        'color': 'darkred',
        'weight': 'normal',
        'size': 16,
        }
cbar.set_label('Brightness Temperature ( K )', fontdict=font)
m.drawmeridians(np.arange(110, 140, 5), labels=[0, 0, 0, 1])
m.drawparallels(np.arange(10, 40, 5), labels=[1, 0, 0, 0])
# 將最佳路徑集上的經(jīng)緯度映射到地圖上，再把該點(diǎn)繪制在地圖上
best_lon, best_lat = m(best_lon, best_lat)
m.plot(best_lon, best_lat, 'o', color='fuchsia', ms=5)
plt.xticks(fontsize=20)
plt.yticks(fontsize=20)
plt.title('Brightness Temperature(2019081112) ', fontdict=font)
plt.savefig('test.png')
plt.show()

注：讀取地圖文件時(shí)，可以使用readshapefile函數(shù)。

以2019081112時(shí)刻為例，對(duì)應(yīng)地圖10-40，110-140的區(qū)域，讀取對(duì)應(yīng)的hdf文件，繪制該點(diǎn)的衛(wèi)星云圖的亮溫值。