ArcPy批次對大量遙感影像相減做差

  本文介紹基於Python中ArcPy模組，對大量柵格遙感影像檔案批次進行相減做差的方法。

  首先，我們來明確一下本文的具體需求。現有一個儲存有多張.tif格式遙感影像的資料夾，其中每一個遙感影像的檔名中都包含有該影象的成像年份，且每一個遙感影像的空間範圍、像元大小等都是一致的，可以直接進行柵格相減；且資料夾內除了.tif格式的遙感影像檔案外，還具有其它格式的檔案；如下圖所示。

  我們希望，對於同一年成像的兩景遙感影像分別進行做差處理。例如，將上圖中的2001.tif檔案減去2001_N.tif檔案，將2005.tif檔案減去2005_N.tif檔案，以此類推。

  明確了需求後，我們就可以開始具體的操作。首先，本文所需用到的程式碼如下。

# -*- coding: utf-8 -*-
"""
Created on Sun Apr 24 11:12:37 2022

@author: fkxxgis
"""

import arcpy

tif_file_path="E:/LST/Data/MODIS/16_True/"
dif_file_path="E:/LST/Data/MODIS/17_Difference/"
arcpy.env.workspace=tif_file_path

tif_file_name=arcpy.ListRasters("*","tif")
tif_file_year=tif_file_name[0][0:4]
one_year_tif_list=[]

for tif_file in tif_file_name:
    if tif_file[0:4]==tif_file_year:
        one_year_tif_list.append(tif_file)
        if tif_file==tif_file_name[len(tif_file_name)-1]:
            arcpy.gp.Minus_sa(one_year_tif_list[0],
                              one_year_tif_list[1],
                              dif_file_path+tif_file_year+"_Dif.tif")
    else:
        arcpy.gp.Minus_sa(one_year_tif_list[0],
                          one_year_tif_list[1],
                          dif_file_path+tif_file_year+"_Dif.tif")
        one_year_tif_list=[]
        one_year_tif_list.append(tif_file)
        tif_file_year=tif_file[0:4]

  其中，tif_file_path是原有計算平均值前遙感影象的儲存路徑，dif_file_path是我們新生成的求取平均值後遙感影像的儲存路徑，也就是結果儲存路徑。

  在這裡，和我們前期的部落格Python ArcPy批次拼接長時間序列柵格影象類似，需要首先在資源管理器中，將tif_file_path路徑下的各檔案以「名稱」排序的方式進行排序；隨後，利用arcpy.ListRasters()函數，獲取路徑下原有的全部.tif格式的影象檔案，並擷取第一個檔案的部分檔名，從而獲取其成像時間的具體年份。

  接下來，遍歷tif_file_path路徑下全部.tif格式影象檔案。其中，我們通過一個簡單的判斷語句if tif_file[0:4]==tif_file_year:，來確定某一年的遙感影像是否已經讀取完畢——如果已經讀取完畢，例如假如2001年成像的2幅遙感影像都已經遍歷過了，那麼就對這2景遙感影像做差，並開始對下一個年份（即2005年）成像的2景遙感影像繼續加以計算；如果還沒有讀取完畢，例如假如2001年成像的2幅遙感影像目前僅遍歷了第1幅，那麼就不做差，繼續往下遍歷，直到遍歷完2001年成像的2幅遙感影像。

  這裡相信大家也看到了為什麼我們要在前期先將資料夾中的檔案按照「名稱」排序——首先，是為了保證同一年成像的2景遙感影像都排列在一起，遍歷時只要遇到一個新的年份，程式就知道上一個年份的2張影象都已經遍歷完畢了，就可以將上一個年份的2張柵格影象加以做差；其次，是為了保證我們的被減數（例如2005.tif檔案）排在減數（例如2005_N.tif檔案）的前面，從而方便我們進行做差運算。

  在這裡，我們實現兩張柵格遙感影像相減操作的函數是arcpy.gp.Minus_sa()函數，其第一個引數是被減數，第二個引數是減數，第三個引數是結果儲存路徑與名稱。

  最後，通過if tif_file==tif_file_name[len(tif_file_name)-1]:這個判斷，來確認是否目前已經遍歷到資料夾中的最後一個影象檔案。如果是的話，就需要將當前成像年份的2景影象進行差值的求取，並宣告程式碼完成執行。

  在 IDLE (Python GUI) 中執行程式碼。程式碼執行完畢後，我們可以看到求取差值之後的遙感影像已經存在於我們的結果儲存路徑中了。

  至此，大功告成。