Python如何實現批量讀取HDF多波段柵格數據并繪制像元直方圖

本篇內容介紹了“Python如何實現批量讀取HDF多波段柵格數據并繪制像元直方圖”的有關知識，在實際案例的操作過程中，不少人都會遇到這樣的困境，接下來就讓小編帶領大家學習一下如何處理這些情況吧！希望大家仔細閱讀，能夠學有所成！

首先將本文所需代碼展示如下：

# -*- coding: utf-8 -*-
"""
Created on Tue Jul 20 11:05:31 2021

@author: fkxxgis
"""

import os
import numpy as np
from osgeo import gdal
import matplotlib.pyplot as plt

lai_file_path="G:/Postgraduate/LAI_Glass_RTlab/Test_DRT/h30v09.tif"
mcd_file_path="G:/Postgraduate/LAI_Glass_RTlab/Test_DRT/MCD12Q1.A2018001.h30v09.006.2019199233851.hdf"
pic_save_path="G:/Postgraduate/LAI_Glass_RTlab/Test_DRT/"

for veg_type in range(9):

    mcd_raster=gdal.Open(mcd_file_path)
    mcd_sub_dataset=mcd_raster.GetSubDatasets()
    hdf_band_num=len(mcd_sub_dataset)
    # for sub_dataset in mcd_sub_dataset:
    #     print(sub_dataset[1])
    # print(mcd_sub_dataset[2][1])
    mcd_sub_type=gdal.Open(mcd_sub_dataset[2][0])
    mcd_raster_array=mcd_sub_type.ReadAsArray()
    
    lai_raster=gdal.Open(lai_file_path)
    lai_raster_array=lai_raster.ReadAsArray()
    non_veg_type_lai_array=np.where(mcd_raster_array==veg_type+1,lai_raster_array,np.nan)
    plt.hist(non_veg_type_lai_array)
    plt.savefig(pic_save_path+"DRT_"+str(veg_type+1)+".png", dpi=300)
    plt.clf()
    plt.cla()

我們直接講解多波段HDF柵格圖像文件讀取部分的代碼：首先，多波段.hdf格式文件的讀取在一開始與單波段.hdf格式文件或.tif格式文件的讀取一致，即通過gdal.Open()函數實現；但隨后，需要額外借助len()函數獲取HDF文件對應的波段數量。

因為我們讀取的HDF文件是多波段，因此hdf_band_num肯定是大于1的，那么剛剛讀取進來的mcd_sub_dataset其實就是一個列表（List）；其中，這個列表的元素個數就是對應的多波段HDF文件波段數，列表的每一個元素則都是一個元組（tuple）；同時，每一個元組都有兩個元素，其每一個元素都是一個字符串；其中第一個元素為當前HDF文件的當前波段對應的文件路徑與部分提示信息，第二個元素作為當前HDF文件的當前波段對應的文件像素行列數、名稱與數據類型。

這么說可能不太明白，我們用一個實例來講解。例如，通過上述代碼讀取一景具有六個波段的MODIS LAI產品——MCD15A3H產品，其第一個波段為FPAR數據，第二個波段為LAI數據。那么讀取其后，得到的mcd_sub_dataset長這個樣子：

可以看到，是一個具有6個元素的列表。

點開列表，可以看到6個元素每一個都是一個具有2個元素的元組：

再點開第一個元組，可以看到其具有2個字符串格式的元素：

其第二個元素包含了該波段對應的數據行數與列數（即[2400×2400]）、數據名稱（即Fpar）、數據空間分辨率（即500m）、數據產品簡稱（即MOD_Grid_MCD15A3H），以及數據格式（即8-bit unsigned integer）；而第一個字符串沒有顯示完畢，我們可以點擊打開看看：

可以看到第一個元素則包含了該波段對應的數據路徑、文件全稱，以及部分與第二個元素重復的幾個數據信息參數。

有了上面的分析就比較清楚了，接下來再一次利用gdal.Open()函數讀取我們需要的波段，mcd_sub_dataset[2][0]表示第三個波段；其中，第三個波段卻用[2]來表示，是因為波段數量（也就是mcd_sub_dataset的Index）是從0開始計算的；而后面的[0]則表示元組中的第一個參數，也就是上面一幅圖中顯示的該波段對應的數據路徑。

隨后，再利用.ReadAsArray()函數將其讀取為Array即可。接下來的操作與本文開頭提及的那兩篇博客就一致，這里不再贅述~

“Python如何實現批量讀取HDF多波段柵格數據并繪制像元直方圖”的內容就介紹到這里了，感謝大家的閱讀。如果想了解更多行業相關的知識可以關注億速云網站，小編將為大家輸出更多高質量的實用文章！