NetCDF 文件批量轉柵格並導出柵格各波段

兩年前，我曾發佈過一篇名為《導出 NetCDF 柵格圖層的各個波段》的公眾號推文，講述了通過網路中的ArcGIS工具將單個 NetCDF 文件的各個波段分別導出為 tif 文件的方法。該工具提供了 arcpy 源代碼，我們以該代碼為基礎，將其轉換為 ArcGIS Pro 環境下的 Python 3 代碼，並使程式可對多個文件進行批處理，本文將對該代碼及其在 ArcGIS Pro 中的配置過程進行介紹。

1 數據來源及介紹

本文所使用的數據為國家青藏高原科學數據中心的 中國1km解析度逐月降水量數據集。對下載得到的 nc 文件通過 Panoply 軟體進行讀取可知，其X、Y及波段維度名稱分別為lon、lat、time，nc文件中共12個波段，分別存儲該年度1-12月每個月份的降水量。

2 腳本配置及相關代碼

2.1 NetCDF 轉換為 tif

在 ArcGIS Pro 地圖工程的預設工具箱中新建腳本工具，命名為 NetCDF to TIF，並仿照 ArcGIS 的 創建 NetCDF 柵格圖層 工具設置其參數。該工具用於讀取文件夾中的所有 NetCDF 文件，根據其 X、Y、波段等維度將其轉換為 TIF 文件並逐一輸出至 TIF 文件夾中。由於單波段 NetCDF 文件不存在波段維度，故該參數設置為可選類型。

右鍵單擊腳本工具，在彈出菜單中單擊編輯，將編輯視窗中的預設代碼替換為下方 Python 代碼：

import os
import arcpy

# 腳本工具核心代碼
def script_tool(nc_folder, variable, x_dimension, y_dimension, band_dimension, tif_folder):
    # 基於工作空間為 tif_folder，後續輸出文件都將存放在該文件夾中
    arcpy.env.workspace = tif_folder
    # 遍歷文件夾中所有尾碼名為 .nc 的文件，將其逐一轉換為 tif 文件
    for file in os.listdir(nc_folder):
        if file[-3:] == ".nc":
            # 以原始文件名稱為基礎，得到不包含尾碼名部分的 tif 文件名
            tif_file_name = file[:-3]
            # 判斷 band_dimension 參數是否存在輸入值，有值則在後續函數中指定 band_dimension 參數
            if band_dimension:
                # 基於 NetCDF 文件創建柵格圖層
                arcpy.MakeNetCDFRasterLayer_md(f"{nc_folder}/{file}", variable, x_dimension, y_dimension, tif_file_name, band_dimension)
            else:
                arcpy.MakeNetCDFRasterLayer_md(f"{nc_folder}/{file}", variable, x_dimension, y_dimension, tif_file_name)
            # 基於柵格圖層創建 tif 文件
            arcpy.CopyRaster_management(tif_file_name, tif_file_name + ".tif")
            arcpy.AddMessage(tif_file_name + " " + "conversion successful!")
    return

# 讀取輸入參數，執行程式
if __name__ == "__main__":
    # 獲取 ArcGIS 腳本工具參數
    nc_folder = arcpy.GetParameterAsText(0)
    variable = arcpy.GetParameterAsText(1)
    x_dimension = arcpy.GetParameterAsText(2)
    y_dimension = arcpy.GetParameterAsText(3)
    band_dimension = arcpy.GetParameterAsText(4)
    tif_folder = arcpy.GetParameterAsText(5)
    # 執行腳本工具核心代碼
    script_tool(nc_folder, variable, x_dimension, y_dimension, band_dimension, tif_folder)

保存後可在右鍵菜單的工具屬性中看到下方界面：

配置完成後，將 Panoply 讀取到的參數輸入腳本，如下圖所示：

程式運行完成後輸出的柵格數據如下圖所示：

2.2 遍歷並導出柵格各波段

新建腳本工具，命名為 Band Split，並設置參數。

將腳本工具的預設代碼替換為下方 Python 代碼：

import os
import arcpy

def script_tool(input_folder, output_folder):
    # 遍歷文件夾中所有 tif 文件
    for file in os.listdir(input_folder):
        if file[-4:] == ".tif":
            # 根據文件名稱獲取對應年份，此處為倒數第5位至倒數第8位，可根據實際情況調整
            year = file[-8:-4]
            # 創建迴圈，遍歷各月份對應波段，將波段逐一導出為單獨的 tif 文件並保存到輸出文件夾
            for i in range(12):
                input_band = f"{input_folder}/{file}/Band_{i+1}"
                output_raster = f"{output_folder}/Pre_{year}_{i + 1}.tif"
                arcpy.CopyRaster_management(input_band, output_raster)
                arcpy.AddMessage(f"Pre_{year}_{i+1} conversion successful!")
    return

if __name__ == "__main__":
    # 獲取 ArcGIS 腳本工具參數
    input_folder = arcpy.GetParameterAsText(0)
    output_folder = arcpy.GetParameterAsText(1)
    # 執行腳本工具核心代碼
    script_tool(input_folder, output_folder)

配置完成後，將上一步得到的 tif 文件夾作為輸入文件夾參數輸入腳本，如下圖所示：

程式運行完成後輸出的柵格數據如下圖所示：