告別ArcGis, 擁抱Python(1)：DEM輸入和顯示
前言：
樓主為在讀博士，先前無意間在網上搜到一個帖子，名為「Goodbye ArcGIS, hello Python」（http://pangea.stanford.edu/~samuelj/musings/dems-in-python-pt-1-intro.html），主要是介紹如何使用Python做簡單的DEM處理。作者為斯坦福大學的一名在讀博士生，可惜只更新了三章，就沒有了……Python因為其開源、簡單和豐富的庫在國外學術圈比較流行，但是在國內因為破解matlab成本較低，因此好像普及情況一直不是很高。樓主心血來潮，準備順著這位博士的思路，嘗試做個python和空間分析的小系列，希望能堅持下去。因為樓主沒有系統的學過python，如有問題或者改進的地方，請指出，謝謝！
題圖來自這個網站
————————————————————————————————————
開始之前：
1）這裡使用的是win10系統下的python2.7，沒使用3是因為2.7相對比較成熟，各種資源都是現成的，省去了不少麻煩。詳細的關於2.7和3.X的利弊可以參考這裡：
https://www.zhihu.com/question/20604044?rf=20605744
沒有使用macos是因為...沒錢...
（2）我現在使用的IDE是Anaconda中的Spyder，可以從這裡下載：下載鏈接
（3）學習過程中需要用的資料會通過OneDrive 分享，點擊鏈接即可下載
（4）樓主沒有系統的學習過GIS的知識（大學學費白交了），大部分都是自學，很多知識都很散，說得不好的請批評指出，讓樓主也學習進步！
之前使用過Enthought Canopy（因為導師推薦，用了大概一年） : Python Distribution and Integrated Analysis Environment
另外也有深受眾多大神追捧的pycharm，也是很好的選擇。
我個人比較喜歡anaconda和canopy的原因是，這兩個界面 IDE界面簡單，適合博士狗的科學計算…………
PS: 樓主學習壓力比較大……更新可能比較慢
————————————————————————————————————
今天是第一天，簡單的介紹下digital elevation mode(DEM)和如何在python中讀取和顯示DEM。
DEM是地學中（包括地理學、地形學、地貌學、大氣科學、水文學....以及時下最熱的氣候變化）最最最最最基本的數據。 DEM在wiki百科上是這麼解釋的：
A digital elevation model (DEM) is a digital model or 3D representation of a terrains surface — commonly for a planet (including Earth), moon, or asteroid — created from terrain elevation data.
翻譯過來就是用於表徵地表起伏的三維數字模型。樓主個人理解就是，含有高程信息（通常為海拔）的二階矩陣（圖1右，橫縱坐標通常為經緯度）。DEM的保存格式有很多，例如tiff，geoTiff，甚至可以用ascii表示。在這裡我個人比較喜歡喜歡用netcdf格式，原因是（1）在處理較大的範圍DEM時（例如，全球的dem）可以通過每次讀取一部分來避免內存不足的情況（當然如果有超算中心並且可以無限量使用的話，就無視吧）；（2）這個格式是很多氣候資料中非常常見的格式，很通用；（3）樓主就只懂這個格式。
netcdf的詳細介紹可以在這裡找到（不客氣）：NetCDF: Introduction and Overview
dem資料下載：https://1drv.ms/f/s!Avm_B4scAY3abHqNhJQCE5ZiBvU
dem是ASTER的數據，精度為1arc(~30 m)
好了，下面我們正式開始通過python讀取和顯示DEM。
————————————————————————————————————
第一步：讀取DEM文件
在這裡我們會使用到python中一個非常強大的模塊netcdf4，可以通過這個網站自學（有機會再講吧，各種挖坑）：netCDF4 API documentation
import netCDF4 as nc#python中強大的netcdf包import matplotlib.pylab as plt#畫圖模塊dem = nc.Dataset(example.nc)#讀取dem文件
dem為樓主自己轉的數據，的輸出結果為：
<type netCDF4._netCDF4.Dataset>
root group (NETCDF4_CLASSIC data model, file format HDF5):
description: Digital Earth Model with a resolution of 1 arc
source: Created by FreezeMe and for Python and GIS
dimensions(sizes): lon(3599), lat(3598)
variables(dimensions): float32 lon(lon), float32 lat(lat), float32 elevation(lat,lon)
groups:
翻譯：
1弧度（約30m精度）
數據大小為：3599（經度）*3598（緯度）
變數有：經度（lon），緯度（lat）和海拔（elevation），group信息為空
第二步：提取DEM中的地理信息
可以看到導入DEM文件後我們並沒有看到DEM和經緯度具體的信息，下面我們將提取這些信息。這裡我們就要用到netcdf裡面的變數命令，具體如下：
lons = dem.variables[lon][:]#經度值lats = dem.variables[lat][:]#緯度值ele = dem.variables[elevation][:]#海拔值
如果嘗試輸入dem.variables[lon]將會給出經度的基本信息：
<type netCDF4._netCDF4.Variable>
float32 lon(lon)
units: degree_east (demical)
unlimited dimensions:
current shape = (3599,)
filling on, default _FillValue of 9.96920996839e+36 used
第三步：繪製DEM
plt.figure(elevation)#創建名為elevation的圖窗plt.imshow(ele, cmap = rainbow)#繪製海拔圖，配色為rainbowplt.colorbar()#繪製海拔圖例dem.close()#關閉dem文件
下圖既是我們的DEM顯示效果
圖2 DEM圖
有興趣的朋友可以考慮下怎麼更改圖中的橫縱坐標，使其顯示經緯度信息。
下次我們將講解一下如何使用python得到坡度，坡向以及陰影圖
推薦閱讀：