python如何處理衛(wèi)星云圖獲取亮溫值

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

這篇文章主要介紹了python如何處理衛(wèi)星云圖獲取亮溫值問題,具有很好的參考價(jià)值,希望對(duì)大家有所幫助,如有錯(cuò)誤或未考慮完全的地方,望不吝賜教

根據(jù)最佳路徑數(shù)據(jù)集，獲取經(jīng)緯度和時(shí)間。
根據(jù)時(shí)間找對(duì)應(yīng)的hdf文件，獲取其中的經(jīng)緯度，然后根據(jù)hdf文件中的經(jīng)度確定是哪顆衛(wèi)星的。
在全圓盤標(biāo)對(duì)稱文件經(jīng)緯度對(duì)照表dat文件獲取經(jīng)緯度時(shí)加上相應(yīng)的經(jīng)度。
最后根據(jù)最佳路徑的經(jīng)緯度獲取dat文件中最近格點(diǎn)的經(jīng)緯度索引，再通過hdf文件獲得亮溫值。

1.處理hdf文件

import h5py
hdfFile = h5py.File('F:/新建文件夾/Satellite_Imagery/IR_data/韋森特/FY2D_FDI_ALL_NOM_20120723_1945.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']
a = hdfFile['NomFileInfo'][()]
print(a[0][3])  # 獲取緯度
print(a[0][4])  # 獲取經(jīng)度

對(duì)于風(fēng)云衛(wèi)星的數(shù)據(jù)格式，可以查閱相應(yīng)的資料，獲取靜止氣象衛(wèi)星的數(shù)據(jù)格式。這里不再詳細(xì)描述。

亮溫值的獲取分為四個(gè)IR通道和一個(gè)可見光通道，不同通道的用處不同，根據(jù)自己的需求選取數(shù)據(jù)。

在這里插入圖片描述

對(duì)于獲取衛(wèi)星經(jīng)緯度等數(shù)據(jù)，需要，根據(jù)NomFileInfo的數(shù)據(jù)集分析。

a = hdfFile['NomFileInfo'][()]
print(a)
[(b'FY2D', b'NOM', b'Nomalized Projection Full Disc Image', 0., 86.5, 35785864., b'fit', 2012, 7, 10, 12, 34, 0, 2012, 7, 23, 19, 45, 1, 2012, 7, 23, 19, 58, 58, 2012, 7, 23, 20, 0, 27, 6378137., 3.5e-05, 0.00014, 298.25722356, b'NOM Fit HDF5')]

對(duì)應(yīng)信息見下圖數(shù)據(jù)集。

獲取亮溫值：

獲取通道1的亮溫值，代碼如下：

hw = hw1[()]  # 通道IR1，定標(biāo)表
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]
print(tb)
tb = tb.T  # 注意需要轉(zhuǎn)置

2. 處理dat文件

對(duì)于FY2氣象衛(wèi)星而言，經(jīng)緯度查找表文件（NOM_ITG_2288_2288(0E0N)_LE.zip）可以從網(wǎng)上下載。

數(shù)據(jù)解壓之后，里面有三個(gè)文件，分別是：

在數(shù)據(jù)說明文件中，詳細(xì)說明了經(jīng)緯度查找表數(shù)據(jù)的使用方法，同時(shí)針對(duì)FY2系列衛(wèi)星的經(jīng)度進(jìn)行了說明，由于上述經(jīng)緯度數(shù)據(jù)是以中心點(diǎn)為0度經(jīng)度生成的，所以對(duì)于FY2的數(shù)據(jù)，需要在經(jīng)度數(shù)據(jù)上加上衛(wèi)星所在的經(jīng)度。

5	FY-2G	104.5°E
6	FY-2H	79°E

對(duì)于不同的星可以根據(jù)hdf文件中的經(jīng)度，判斷對(duì)應(yīng)的查找表的經(jīng)度。

下面以H星為例，讀取dat文件，獲取網(wǎng)格經(jīng)緯度

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]

這里需要注意的是，需要先取出來經(jīng)度，再取出緯度，這是經(jīng)緯度在dat文件中的存儲(chǔ)順序，然后加上衛(wèi)星的經(jīng)度，這里是79，故經(jīng)度加79，獲取經(jīng)緯度。

# 經(jīng)度表
[[379. 379. 379. ... 379. 379. 379.]
 [379. 379. 379. ... 379. 379. 379.]
 [379. 379. 379. ... 379. 379. 379.]
 ...
 [379. 379. 379. ... 379. 379. 379.]
 [379. 379. 379. ... 379. 379. 379.]
 [379. 379. 379. ... 379. 379. 379.]]
# 緯度表 
 [[300. 300. 300. ... 300. 300. 300.]
 [300. 300. 300. ... 300. 300. 300.]
 [300. 300. 300. ... 300. 300. 300.]
 ...
 [300. 300. 300. ... 300. 300. 300.]
 [300. 300. 300. ... 300. 300. 300.]
 [300. 300. 300. ... 300. 300. 300.]]

這樣就獲取了dat文件中的網(wǎng)格經(jīng)緯度，即2288*2288的對(duì)照表，這樣根據(jù)最佳路徑數(shù)據(jù)集中的經(jīng)緯度，找到網(wǎng)格中最近的經(jīng)緯度。

注：因?yàn)榫W(wǎng)格存儲(chǔ)的經(jīng)緯度不是等間距的，所以不能直接近似找最近的點(diǎn)，需要根據(jù)角度計(jì)算距離，找到網(wǎng)格上離所求經(jīng)緯度最近的網(wǎng)格點(diǎn)，找到對(duì)應(yīng)的索引值，再根據(jù)索引值找到亮溫值。