Python 柵格數(shù)據(jù)處理教程(二)
本文將介紹通過 ArcGIS Pro 的 Python 模塊(arcpy)對柵格數(shù)據(jù)進行柵格計算及數(shù)據(jù)統(tǒng)計的方法。
1 數(shù)據(jù)來源及介紹
本文使用的數(shù)據(jù)為國家青藏高原科學(xué)數(shù)據(jù)中心的 中國1km分辨率逐月降水量數(shù)據(jù)集基礎(chǔ)上通過《Python 柵格數(shù)據(jù)處理教程(一)》中的方法提取出的吉林省范圍降水量數(shù)據(jù)。該數(shù)據(jù)降水量單位為 0.1 mm,本文通過柵格計算器將數(shù)據(jù)單位批量轉(zhuǎn)換為 mm。根據(jù)各月份降水量疊加生成年降水量柵格及多年平均降水量柵格,同時統(tǒng)計歷年降水量數(shù)據(jù)信息,并整理為 Excel 表。
2 示例代碼
2.1 降水量單位轉(zhuǎn)換
import os
import arcpy
from arcpy import sa
from tqdm import tqdm # 進度條工具,需自行安裝
# 設(shè)置工作空間
arcpy.env.workspace = "Pre_Jilin.gdb"
# 如果在代碼同級文件夾中沒有結(jié)果數(shù)據(jù)庫,則創(chuàng)建該數(shù)據(jù)庫
if "Pre_Jilin1.gdb" not in os.listdir("."):
print("即將為您創(chuàng)建用于保存結(jié)果柵格的文件地理數(shù)據(jù)庫(Pre_Jilin1.gdb)……")
arcpy.CreateFileGDB_management(".", "Pre_Jilin1")
print("數(shù)據(jù)庫創(chuàng)建完成!")
# 遍歷工作空間內(nèi)所有柵格數(shù)據(jù)
for raster in tqdm(arcpy.ListRasters()):
# 通過柵格計算器將柵格值單位轉(zhuǎn)換為 mm,并保存到結(jié)果數(shù)據(jù)庫
sa.RasterCalculator([raster], [raster], f"{raster} * 0.1").save(f"Pre_Jilin1.gdb/{raster}")
print("計算完成!")
柵格數(shù)據(jù)轉(zhuǎn)換前如下所示:

柵格單位轉(zhuǎn)換后如下所示:

2.2 根據(jù)各月份降水量生成年降水量柵格
import os
import arcpy
from arcpy import sa
from tqdm import tqdm
# 設(shè)置工作空間
arcpy.env.workspace = "Pre_Jilin1.gdb"
# 創(chuàng)建柵格文件字典,以年份為鍵,值為該年份對應(yīng)的各月份降水量柵格列表
file_dict = dict()
for raster in arcpy.ListRasters():
# [4:8] 代表年份位于柵格數(shù)據(jù)名中的位置,需根據(jù)情況自行調(diào)整
if raster[4:8] not in file_dict:
file_dict[raster[4:8]] = [raster]
else:
file_dict[raster[4:8]].append(raster)
# 如果在代碼同級文件夾中沒有結(jié)果數(shù)據(jù)庫,則創(chuàng)建該數(shù)據(jù)庫
if "Pre_Jilin_Year.gdb" not in os.listdir("."):
print("即將為您創(chuàng)建用于保存結(jié)果柵格的文件地理數(shù)據(jù)庫(Pre_Jilin_Year.gdb)……")
arcpy.CreateFileGDB_management(".", "Pre_Jilin_Year")
print("數(shù)據(jù)庫創(chuàng)建完成!")
# 定義變量名,并構(gòu)建代數(shù)表達式,以對輸入柵格執(zhí)行空間分析
# 以 a-l 代表12個月份的降水量柵格,求全年總降水量
input_names = ["a", "b", "c", "d", "e", "f", "g", "h", "i", "j", "k", "l"]
expression = "a + b + c + d + e + f + g + h + i + j + k + l"
# 遍歷柵格文件字典
for year in tqdm(file_dict):
# 提取歷年降水量柵格,通過柵格計算器生成年降水量柵格,并保存到結(jié)果數(shù)據(jù)庫中
sa.RasterCalculator(file_dict[year], input_names, expression).save(f"Pre_Jilin_Year.gdb/Pre_{year}")
print("計算完成!")
2023年吉林省年降水量柵格如下圖所示:

2.3 統(tǒng)計研究區(qū)年降水量的最值、均值、標(biāo)準(zhǔn)差等信息
import os
import arcpy
from arcpy import sa
from tqdm import tqdm
# 設(shè)置工作空間
arcpy.env.workspace = "Pre_Jilin_Year.gdb"
# 如果在代碼同級文件夾中沒有結(jié)果數(shù)據(jù)庫,則創(chuàng)建該數(shù)據(jù)庫
if "Pre_Jilin_Year_Info.gdb" not in os.listdir("."):
print("即將為您創(chuàng)建用于保存結(jié)果柵格的文件地理數(shù)據(jù)庫(Pre_Jilin_Year_Info.gdb)……")
arcpy.CreateFileGDB_management(".", "Pre_Jilin_Year_Info")
print("數(shù)據(jù)庫創(chuàng)建完成!")
# 遍歷工作空間內(nèi)所有柵格數(shù)據(jù)
for raster in tqdm(arcpy.ListRasters()):
# 統(tǒng)計吉林省范圍內(nèi)年降水量的最值、均值、標(biāo)準(zhǔn)差等信息,并以表格的形式保存到結(jié)果數(shù)據(jù)庫中
sa.ZonalStatisticsAsTable("Data/Jilin.shp", "FID", raster, f"Pre_Jilin_Year_Info.gdb/{raster}")
print("統(tǒng)計完成!")
2.4 將統(tǒng)計信息表格導(dǎo)出為 txt 文本文件
import os
import arcpy
from tqdm import tqdm
arcpy.env.workspace = "Pre_Jilin_Year_Info.gdb"
# 如果在代碼同級文件夾中沒有結(jié)果文件夾,則創(chuàng)建該文件夾
if "Pre_Jilin_Year_Table" not in os.listdir("."):
print("即將為您創(chuàng)建用于保存轉(zhuǎn)換結(jié)果的文件夾(Pre_Jilin_Year_Table)……")
os.mkdir("Pre_Jilin_Year_Table")
print("文件夾創(chuàng)建完成!")
# 遍歷工作空間內(nèi)所有表格,將表格逐一導(dǎo)出為 txt 文本文件
for table in tqdm(arcpy.ListTables()):
arcpy.ExportTable_conversion(table, f"Pre_Jilin_Year_Table/{table}.txt")
print("轉(zhuǎn)換完成!")
2023 年的統(tǒng)計信息文本如下所示:

2.5 將歷年的 txt 文件整理為 Excel 表
import os
import pandas as pd
pre_dict = dict()
# 遍歷文件列表
for fileName in os.listdir("Pre_Jilin_Year_Table"):
# 通過后綴名判斷是否為 txt 文件
if fileName[-4:] == ".txt":
# 基于 UTF-8 編碼讀取 txt 文件內(nèi)容
with open(f"Pre_Jilin_Year_Table/{fileName}", "r", encoding="utf-8") as file:
file_txt = file.readlines()
# 文本內(nèi)容的第一行是標(biāo)題,第二行是數(shù)據(jù)
# 讀取第二行數(shù)據(jù)第5個數(shù)值開始的所有數(shù)據(jù)值,即 MIN,MAX,RANGE,MEAN,STD,SUM,MEDIAN,PCT90
# 注:ArcGIS Pro 不同版本轉(zhuǎn)換得到的 txt 數(shù)據(jù)表字段可能不同,需要根據(jù)實際情況調(diào)整代碼
val_list = file_txt[1].replace("\n", "").split(",")[4:]
# 將所有數(shù)據(jù)值從字符串類型轉(zhuǎn)換為數(shù)值類型
for i in range(len(val_list)):
val_list[i] = eval(val_list[i])
# 以數(shù)值類型年份為鍵,數(shù)據(jù)值列表為值,將數(shù)據(jù)存入字典
pre_dict[eval(fileName[4:8])] = val_list
# 將字典轉(zhuǎn)換為 DataFrame 二維數(shù)據(jù)表
df = pd.DataFrame(pre_dict, index="MIN,MAX,RANGE,MEAN,STD,SUM,MEDIAN,PCT90".split(",")).T
print(df)
# 將 DataFrame 導(dǎo)出為 Execl 表格
df.to_excel("Pre_Jilin_Year_Table.xlsx")
print("導(dǎo)出表格完成!")
簡單調(diào)整 Excel 表的格式后,其展示效果如下圖所示:

2.6 基于歷年降水量生成多年平均降水量柵格
import os
import arcpy
from arcpy import sa
# 設(shè)置工作空間,并根據(jù)數(shù)據(jù)庫中的柵格數(shù)據(jù)數(shù)量確定年份數(shù)量
arcpy.env.workspace = "Pre_Jilin_Year.gdb"
year_num = len(arcpy.ListRasters())
# 定義變量名,并構(gòu)建代數(shù)表達式,以對輸入柵格執(zhí)行空間分析
# input_names = ["Pre_2001", "Pre_2002", "Pre_2003", ......, "Pre_2021", "Pre_2022", "Pre_2023"]
# expression = "(Pre_2001 + Pre_2002 + Pre_2003 + ...... + Pre_2021 + Pre_2022 + Pre_2023) / 23"
input_names, expression = [], ""
for raster in arcpy.ListRasters():
input_names.append(raster)
if expression == "":
expression += raster
else:
expression += f" + {raster}"
expression = f"({expression}) / {year_num}"
# 如果在代碼同級文件夾中沒有結(jié)果數(shù)據(jù)庫,則創(chuàng)建該數(shù)據(jù)庫
if "Pre_Jilin_Mean.gdb" not in os.listdir("."):
print("即將為您創(chuàng)建用于保存結(jié)果柵格的文件地理數(shù)據(jù)庫(Pre_Jilin_Mean.gdb)……")
arcpy.CreateFileGDB_management(".", "Pre_Jilin_Mean")
print("數(shù)據(jù)庫創(chuàng)建完成!")
# 計算多年平均降水量柵格并保存到結(jié)果數(shù)據(jù)庫中
sa.RasterCalculator(arcpy.ListRasters(), input_names, expression).save("Pre_Jilin_Mean.gdb/Pre_Mean")
print("計算完成!")
運行完成后生成的結(jié)果柵格如下圖所示:

3 降水?dāng)?shù)據(jù)參考文獻格式
3.1 數(shù)據(jù)的引用
彭守璋. (2020). 中國1km分辨率逐月降水量數(shù)據(jù)集(1901-2023). 國家青藏高原數(shù)據(jù)中心. https://doi.org/10.5281/zenodo.3114194.
Peng, S. (2020). 1-km monthly precipitation dataset for China (1901-2023). National Tibetan Plateau / Third Pole Environment Data Center. https://doi.org/10.5281/zenodo.3114194.
3.2 文章的引用
1、Peng, S.Z., Ding, Y.X., Wen, Z.M., Chen, Y.M., Cao, Y., & Ren, J.Y. (2017). Spatiotemporal change and trend analysis of potential evapotranspiration over the Loess Plateau of China during 2011-2100. Agricultural and Forest Meteorology, 233, 183-194. https://doi.org/10.1016/j.agrformet.2016.11.129
2、Ding, Y.X., & Peng, S.Z. (2020). Spatiotemporal trends and attribution of drought across China from 1901–2100. Sustainability, 12(2), 477.
3、Peng, S.Z., Ding, Y.X., Liu, W.Z., & Li, Z. (2019). 1 km monthly temperature and precipitation dataset for China from 1901 to 2017. Earth System Science Data, 11, 1931–1946. https://doi.org/10.5194/essd-11-1931-2019
4、Peng, S., Gang, C., Cao, Y., & Chen, Y. (2017). Assessment of climate change trends over the loess plateau in china from 1901 to 2100. International Journal of Climatology.
浙公網(wǎng)安備 33010602011771號