滑坡是一種地質過程，山體滑坡時重力會導致岩石、土壤、碎屑或其混合物沿下坡移動。山體滑坡會嚴重破壞建築物、道路及其他基礎設施，並破壞野生野生動物棲息地。一場野火之後，由於植被的突然損失，山體滑坡的風險通常更大。
2020年3月28日19時30分，四川省涼山州木裡縣喬瓦鎮鋤頭灣村與項腳蒙古族鄉項腳村交界處發生森林火災。

2020年7月1日四川省木裡縣項腳鄉因連日持續強降雨引發大股泥石流，造成項腳鄉阿牛窩子組境內山體滑坡、道路中斷、農作物大面積受災，村民房屋、集中安置點被淹沒。

本文探索根據ArcGIS Living Atlas of the World 的地形、影像和降雨柵格圖層生成滑坡敏感性表面，並根據此表面找到火災後山體滑坡風險最高的區域。

查詢資料

許多因素可導致山體滑坡風險增加，包括土壤組成、降雨、植被、坡度和坡向。本文將重點關注三個因素：植被密度、地形陡度和降雨量。
ArcGIS Living Atlas of the World是不斷髮展的精選和權威地理資訊的集合。該集合中可以找到我們需要的影像資料和地形等資料。
由於筆者並未找到合適的降雨量資料，選擇通過NOAA下載該區域部分測站的水文資料，通過差值形成降雨量資料。資料下載方法在ArcGIS pro通過聚合點建立時空立方體一文中有詳細說明。

1.通過美國NOAA網站下載中國測站的水文資料，點選https://gis.ncdc.noaa.gov/maps/，選擇日測資料

2. 定位到四川，框選木裡縣一定範圍內的監測站，Add to Cart

3. 選擇csv格式，資料日期本次選擇了6-8月的降雨量資料，Continue，後續選擇資料屬性帶有Precipitation降水量資料。

稍後會在郵箱中收到郵件，資料結果如下圖（需要注意將DATE轉換為yyyy-mm-dd格式）

降水量資料差值

1. 首先通過新增XY點資料，將降雨量資料新增到地圖中
   
   

2. 在工具中搜尋反距離權重法（IDW），將點插值成柵格表面
   

3. Z值欄位選擇降雨量欄位PRCP，點選執行。
   
   需要注意NOAA的PRCP單位為英寸，需要乘25.4轉為毫米
   

4. 得到插值結果
   
   拉伸後，降雨較多的區域以白色顯示
   

搜尋 Living Atlas of the World

獲取Landsat

在ArcGIS pro中依次點選檢視，目錄窗格，門戶。

選擇Living Atlas，搜尋landsat，選擇Multispectral Landsat，並單擊確定

在Mutispectral Landsat圖層屬性中，新建定義查詢，選擇 AcquisitionDate，選擇災後某一日的資料。（具體影像覆蓋資訊在USGS中檢視會比較快，這裡選擇5月6日-8日間的資料）

處理模板可以根據需要選擇，為了不影響分析，此次將處理模板更改為無，因此影象中的所有光譜波段均可用。

以下是幾種常用的處理模板對比

Bands shortwave IR-2, shortwave IR-1, red (7, 6, 4)

bands red, green, blue (4, 3, 2)

Bands near-IR, red, green (5, 4, 3)

*Bands shortwave IR-1, near-IR, blue (6, 5, 2) *

區劃過火範圍

在工程資料庫中新建要素類，命名為「木裡328火災」，以Landsat為底圖，勾繪出過火範圍（本次未精細區劃）。

獲取Terrain

在目錄窗格中，單擊門戶索引標籤，然後單擊 Living Atlas 索引標籤。

在搜尋欄中，輸入 terrain，然後按 Enter 鍵。
右鍵單擊名為 Terrain 的專案並選擇新增至當前地圖。

Terrain 影像圖層以預設灰度的視覺化效果顯示在地圖上。較高的高程以淺灰色顯示。

處理資料

現在已經收集了建模所需的三個柵格圖層，其中landsat和terrain來自線上的ArcGIS Living Atlas of the World，因此未將任何資料下載到本地計算機。

將影像隔離到一個感興趣的場景

Landsat 圖層包含許多覆蓋不同年份的世界場景，先前已經通過定義查詢將圖層更新為僅顯示2020年5月6日至2020年5月8日之間捕獲的場景，這是在木裡森林火災被完全遏制後的1個月後。你可能僅看到地圖上顯示的一個場景，但是那幾天世界上的其他地點也捕獲了數百個其他場景。接下來，需要使用空間查詢來進一步過濾圖層，使其僅顯示覆蓋木裡森林火災區域的一個場景。

• 單擊功能區上的地圖索引標籤。在選擇組中，單擊按位元置選擇。

• 在地理處理窗格中，對於的輸入要素，選擇 Multispectral Landsat。

• 對於選擇要素，選擇之前區劃的過火範圍**「木裡328火災「」**
  

• 單擊確定。所選場景的數量顯示在地圖檢視的右下方。
  

• 在內容窗格中，右鍵 Multispectral Landsat，指向選擇，然後單擊鎖定至選擇內容。
  

• 此操作將確保僅顯示覆蓋木裡森林火災區域的場景。

通過將定義查詢與空間查詢相結合，Multispectral Landsat 影像圖層已被過濾到一個滿足日期和位置標準的場景。

更改 Terrain 處理模板和符號系統

與 Multispectral Landsat 圖層相似，釋出的 Terrain 圖層也具有不同的處理模板。其中一些模板（例如 Slope_Degrees）可以將資料轉換為可用於分析的圖層；其他模板（例如 Slope_Degrees_Map）是處理後資料的可視表示，並且僅用於製圖。

• 在內容窗格中，關閉 Multispectral Landsat 圖層並開啟 Terrain 圖層。
• 選擇 Terrain 圖層後，單擊功能區上的資料索引標籤。
• 在處理組中，單擊處理模板，然後選擇 Slope_Degrees。
  
  Terrain 圖層在地圖上的外觀隨即發生變化。在內容窗格中，圖例也會更新以顯示 0 到 90 之間的值範圍。較亮的區域具有更陡的坡度，最高可達 90 度。最暗的區域是平坦的，坡度為 0 度。
• 右鍵單擊 Terrain 圖層，然後選擇符號系統。
• 在符號系統窗格中，單擊配色方案選單並選中顯示名稱。
• 向下捲動並選擇坡度配色方案
  

現在，已經以綠色-黃色-紅色配色方案描繪了坡度。下一步將統計資料更改為僅使用當前檢視中的資料，以適應感興趣區域。

• 在符號系統窗格中的統計資料中，對於統計資料，選擇 DRA。
  
  

更改降雨量圖層的符號系統

為了更好地檢視降雨量圖層值的變化，我們更改該圖層的符號系統。

• 對於降雨量圖層，在符號系統窗格中，將配色方案更改為降水量。
• 將統計資訊更改為 DRA。
  
  地圖的外觀會發生變化。藍色區域降雨最多，而紅色區域最少。
  

以上通過應用定義查詢、更改處理模板和修改符號系統，將所需的三個影像圖層進行了可用於分析的處理。

構建柵格函數模板

建立一個柵格函數模板，以分析山體滑坡敏感性。該模板中將包含柵格函數，以對輸入資料進行重分類併合並結果。該模板將生成一個輸出圖層，可用於標識更容易發生山體滑坡的區域。

柵格函數是一種更快、更輕量化的地理處理柵格工具替代方案。可將其動態應用於資料集的顯示畫素。柵格函數的輸出是記憶體柵格圖層。這樣可以縮短處理時間，並免去冗餘資料的建立和儲存。

計算植被指數NDVI

將柵格函數合併到柵格函數模板中，就像 ModelBuilder 模型一樣。

• 在影像索引標籤的分析組中，單擊函數編輯器。
  
  隨即顯示空白的函數編輯器檢視。

• 在函數編輯器工具列上，單擊新增柵格變數按鈕。
  
  一個綠色塊隨即顯示在編輯器中。這是輸入柵格的預留位置。

• 右鍵單擊柵格變數塊，然後單擊重新命名。將現有文字替換為多光譜影像，然後按 Enter 鍵。

• 單擊功能區上的影像索引標籤。在分析組中，單擊柵格函數按鈕。
  隨即顯示柵格函數窗格。

• 在系統索引標籤上的分析組中，找到 NDVI 函數並將其拖動到函數編輯器上。
  

• 指向多光譜影像預留位置。單擊OUT標註，然後將其拖動到 NDVI 函數上的柵格標註。
  這兩個塊隨即會通過灰色箭頭連結在一起。現在，將多光譜影像設定為 NDVI 函數的輸入柵格。
  

• 雙擊 NDVI 函數以開啟其屬性視窗。
  NDVI 函數使用可見光的紅色波長和近紅外波長來計算綠色植被的密度。在多光譜 Landsat 影像圖層中，這些波長是通過波段 4 和 5 進行測量的。
  關於Landsat 8的波段組合定義可以點選這裡

• 在引數索引標籤的可見波段 ID 中，輸入 4。對於紅外波段 ID，輸入 5。

• 單擊變數索引標籤。在柵格行中，選中 IsPublic 框。
  公開柵格輸入引數，以便在執行模板時可以選擇不同的多光譜影像。用於具有不同資料的不同區域，而不僅限於針對此分析找到的特定 Landsat 場景。

• 單擊確定關閉 NDVI 屬性視窗。
  NDVI 函數將評估植被健康狀況，但是要在山體滑坡敏感性分析中使用這些值，則需要將其重對映。
• 在柵格函數窗格中，搜尋重對映，然後將重對映函數拖動到函數編輯器上。
• 從 NDVI 函數上的OUT標註拖動到重對映函數上的柵格標註以將兩者進行連線。
  
• 雙擊重對映函數。
  在引數索引標籤上，確保將重對映定義型別設定為列表，並根據下表更新列表：

最小值	最大值	輸出
0	50	5
50	100	4
100	133	3
133	166	2
166	202	1

該操作會將 NDVI 值重分類或重對映為五個新的指數值。由於植被稀疏會導致更高的山體滑坡風險，因此最低的 NDVI 值（幾乎沒有綠色植被）將接收最高重對映值 5。

分類的值可以根據區域特點進行修改，推薦使用自然間斷法把中斷值求算出來，然後通過表匯入到定義型別中，精度可能會更好一些。

• 單擊常規索引標籤。對於名稱，輸入 NDVI Remap。
• 單擊確定關閉重對映屬性視窗。
• 在函數編輯器工具列上，單擊儲存按鈕。
• 在儲存視窗的名稱中，鍵入山體滑坡敏感性分析。
• 對於描述，可以對該模型進行說明。
• 單擊確定。
  模板的名稱將在函數編輯器索引標籤上更新。
  

計算坡度指數值

與植被指數類似，將坡度圖層重對映為五個指數值。

• 將柵格變數新增到編輯器畫布。

• 將柵格變數重新命名為坡度。

• 從柵格函數窗格中，新增重對映函數並將坡度連線到該函數。
  

• 雙擊重對映函數以開啟重對映屬性視窗。

• 在引數索引標籤上，針對重對映定義型別，選擇列表，並根據下表更新列表：

最小值	最大值	輸出
0	10	1
10	25	2
25	40	3
40	55	4
55	91	5

由於較陡的坡度更容易發生山體滑坡，因此將較高的坡度值對映到較高的指數值。

• 單擊變數索引標籤。對於柵格引數，選中 IsPublic 核取方塊。
• 單擊常規索引標籤。修改名稱為Slope Remap。
• 單擊確定。

計算降雨量指數值

與坡度指數值計算一樣，將降雨量圖層重對映為指數值。

• 將柵格變數新增到函數編輯器中，並將其重新命名為降雨量。
• 新增重對映柵格函數，並將降雨量連線到該函數。
• 雙擊重對映以開啟其屬性視窗。
• 在引數索引標籤上，針對重對映定義型別，選擇列表，並根據下表更新列表：

最小值	最大值	輸出
0	3	1
3	6	2
6	9	3
9	15	4
15	22	5

高日均降水量值將重對映為高指數值。

• 單擊變數索引標籤。對於柵格引數，選中 IsPublic 核取方塊。
• 單擊常規索引標籤。修改名稱為Rainfall Remap。
• 單擊確定。

合併指數值以生成敏感性表面

使用加權總和函數將其合併。此函數將對每個輸入柵格進行加權，然後疊加所有輸入並計算其值的總和。

• 在柵格函數窗格中，搜尋加權總和函數並將其拖動到函數編輯器上。
  

• 將 NDVI Remap、Slope Remap 和 Rainfall Remap 連線到加權總和函數。

• 在函數編輯器工具列上，單擊自動佈局按鈕。
  
  接下來，需要為每個輸入分配權重，以指示哪些變數應被視為分析中最重要的變數。

• 雙擊加權總和柵格函數。
  在燒燬區域，植被缺乏程度將替代降雨或坡度成為評估山體滑坡風險的更為重要的因素。針對此分析，賦予 NDVI 輸入較高的權重，以指示植被指數對結果的影響是坡度或降雨輸入影響的兩倍。

• 在引數索引標籤上的NDVI Remap.OutputRaster錶行中，將權重值更改為 2。

• 對於像元大小型別，選擇最小值。
  這將確保輸出柵格採用最精細的輸入圖層的解析度而非最粗糙的解析度。
  

• 單擊確定，然後儲存柵格函數模板。

柵格函數模板的結果將是一個記憶體中的圖層，其值的範圍可能在 4 到 20 之間，較高的值表示山體滑坡敏感性更強。

指定柵格輸出的位大小

由於輸出值均為 4 到 20 之間的整數，結果可以用 8 位無符號資料型別表示。

-在柵格函數窗格中，搜尋轉為整型柵格函數並將其拖動到函數編輯器上。

• 加權總和函數拖動到轉為整型函數以將兩者進行連線。

• 雙擊轉為整型柵格函數以開啟其屬性視窗。

• 單擊常規索引標籤。對於輸出畫素型別，選擇 8 位無符號。
  

• 單擊確定。

將結果裁剪至感興趣區域

將結果裁剪至木裡森林火災邊界範圍。

• 在柵格函數窗格中，搜尋裁剪柵格函數並將其拖動到函數編輯器上。
  裁剪函數具有兩種可能的輸入連線：一種用於輸入柵格，一種用於裁剪幾何。本次需將轉為整型的輸出連線到裁剪的柵格輸入。

• 從轉為整型函數拖動到裁剪函數上的柵格標註以將兩者進行連線。
  

• 雙擊裁剪柵格函數以開啟其屬性視窗。

• 單擊變數索引標籤。對於 ClippingGeometry 和Extent引數，選中 IsPublic 核取方塊。
  

• 單擊確定。

• 儲存柵格函數模板並關閉函數編輯器。

• 儲存工程。

執行山體滑坡敏感性分析

執行柵格函數模板

在執行函數模板前，先將處理範圍設定到僅處理輸入資料範圍。避免無謂的計算。

• 在功能區上，單擊分析索引標籤。在地理處理組中，單擊環境。
  

• 對於處理範圍，選擇木裡328火災，範圍座標將自動更新以匹配圖層範圍。之後在此工程中執行任何地理處理工具或柵格函數時，將應用此環境設定。
  並行處理因子根據自身機器效能選擇，可以提高處理效率。
  

• 在柵格函數窗格中，單擊之前儲存的山體滑坡敏感性分析模型。

• 四個選項將顯示在引數索引標籤上。這些就是在構建柵格函數模板時設定為 IsPublic 的引數。依次選擇對應的資料。
  

• 單擊建立新圖層按鈕。
  柵格函數模板完成後，新的記憶體柵格圖層將新增到地圖。該圖層包含感興趣區域的山體滑坡敏感性值。該結果是通過動態處理線上資料建立的記憶體圖層。
  

匯出柵格

• 在內容窗格中，右鍵單擊 Landslide Susceptibility Analysis，指向資料，然後選擇匯出柵格。
  匯出柵格窗格隨即出現。
• 像元大小預設設定為 0.25 米。對於如此大的區域，此解析度很小。這裡設定10，將其增大，以使柵格可以更快地匯出。

修改柵格符合系統

匯出的新圖層隨即會自動新增到地圖。
圖層具有黑色背景。可以通過符號系統，掩膜，勾選顯示背景值去除黑色背景
圖層以灰度繪製。白色地區發生山體滑坡的風險最高。

為了使其在視覺上更清晰，重新修改柵格符號系統。

• 在內容窗格中，右鍵單擊 LandslideSusceptibility.tif，然後選擇符號系統。
• 對於配色方案，選擇 地獄紅。
  當前，將以最淺的顏色顯示山體滑坡敏感性最高的區域。
• 重新開啟配色方案選單，然後單擊格式化配色方案。
• 在反向配色方案中，單擊反向配色方案按鈕，然後單擊確定。

當前，將以最深的顏色顯示山體滑坡敏感性最高的區域。

使用分析結果查詢存在風險的位置

通過地理處理工具，可以將道路要素按等距離間隔分配敏感性值。由於沒有該區域道路資訊，因此選擇建立隨機點的方式來評估整個區域範圍內的敏感性位置。

將山體滑坡敏感性值分配至隨機點

• 在資料管理工具中，選擇取樣工具列，選擇建立隨機點。本次建立2000個隨機樣點。
  

• 在地理處理窗格中，搜尋並開啟多值提取至點工具。
  該工具將為點要素新增新屬性。此屬性的值將等於每個點位置的山體滑坡敏感性分析輸出柵格中的像元值。

• 對於輸入點要素，選擇Q2000。

• 對於輸入柵格，選擇Landslide Analysis.tif

• 對於輸出欄位名稱，輸入 Susceptibility_Score。
  

• 單擊執行。每個點將根據基礎柵格分配敏感性得分。
  

具有最高山體滑坡風險的位置

使用定義查詢和分級符號來隔離山體滑坡風險最高的道路點。

• 在內容窗格中右鍵單擊Q2000並選擇屬性。
• 單擊定義查詢索引標籤，然後單擊新建定義查詢。
• 使用下拉式選單構建子句 Where Susceptibility_Score 大於 12。
  
  現在，該地圖僅顯示山體滑坡風險最高的道路點。
• 右鍵單擊Q2000 圖層，然後選擇符號系統。
  分級符號將為具有較大敏感性得分的點繪製較大的符號。
  
• 切換為熱力圖顯示效果如下
  
  縮放到區域性。
  

森林火災的負面影響之一是破壞了健康的綠色植被。當大雨落在已被砍掉植被的陡坡上時，即存在很高的山體滑坡風險。基於植被密度、坡度和降水量，對木裡森林火災地區的山體滑坡敏感性分析迅速確定了具有較高山體滑坡風險的潛在區域。後續可根據地質災害風險相關模型豐富引數，研究結果可用於主動計劃最關鍵地區的山體滑坡減災工作。

鳴謝

本文參考了Mark Gilbert和Aileen Buckley的相關文章，使用的所有資料均來自 ArcGIS Living Atlas of the World 、ArcGIS Online和NOAA。

• 多光譜 Landsat 衛星：Esri、USGS、AWS、NASA
• 地形：Airbus、USGS、NGA、NASA、CGIAR、NLS、OS、NMA、Geodatastyrelsen、GSA、GSI 以及 GIS 使用者社群
• 降雨量：NOAA