寫在前面（重要）
這學期有一門課有關于lisp對CAD二次開發的，老師上課時候放話誰編出了等高線自動繪制程式期末成績可以提檔，當時我下課就去問老師“什麼語言編都行嗎？”老師很肯定的說：随便！。我一聽樂了，這python調個matplotlib庫中的函數畫個等高線圖這不是十分到手了嗎？當老師把測試用的離散的點檔案發過來的時候我就迫不及待去試了。
後來發現，還真不是簡單用個庫就完事了。
contour函數所需要的參數X,Y,Z，Z是關于（X,Y）的函數。網上的文章的執行個體操作都挺千篇一律的。就是先用numpy的linspace生成一個清單x,和一個清單y，然後把x,y網格化成X,Y，接下來很關鍵，網上的文章就會說為了友善，這裡的Z就用f(X,Y）代替，把Z作為以（X,Y）為變量的函數值。這一聽很有道理！但是……假設你是要五個點的（x,y）坐标，你網格化之後就會有25個點，你采用Z=f(X,Y)的方式，你網格化之後有多少個點，他都能根據函數弄一個點的Z出來。但是老師給的文本檔案有兩百多個點，點号，點的橫縱坐标以及高程都給你安排得明明白白，是以當我按照網上的方法嘗試的時候瘋狂報錯，最後放棄，開始老老實實想怎麼用Lisp編出來。
但是，這學期上的地理資訊系統這門課讓我又有點想法了。老師在說到幾何校正這個内容的時候說起重采樣，重采樣是可以根據已知的點得到未知點的資訊，那我也可以對網格化後的x,y重采樣。接下來是過程。
具體過程
首先是讀取資料，讓我們先打開老師發過來的有關于點的文本檔案。
可以看到每一行資料，第一個是點号，第二第三個分别為X,Y坐标，第四個資料為高程。
f=open(r'C:\Users\58381\Desktop\等高線繪制資料.txt','r')
x=[]#建立三個清單用于存放點的資料
y=[]
h=[]
for line in f:
p=line[:-1]#去掉每行最後的換行符
a=p.split(',')#因為每行每個資料之間用","隔開，去掉逗号。
x.append(float(a[1]))
y.append(float(a[2]))
h.append(float(a[3]))
f.close()
然後就把根據x,y來做點表
import numpy as np
points=[]
for i in range(len(x)):
point=[]
point.append(x[i])
point.append(y[i])
points.append(point)
points=np.array(points)
再把這些點網格化
xi=np.linspace(min(x),max(x),300)
yi=np.linspace(min(y),max(y),300)
xi,yi=np.meshgrid(xi,yi)#網格化
接下來就是用這些點來重采樣了，也就是函數内插。python的scipy.interpolate庫裡就有很多關于内插的函數，我們可以采用griddata這個函數。有關griddata，官方的文檔是這樣說的：
scipy.interpolate.griddata(points, values, xi, method=‘linear’, fill_value=nan, rescale=False)
可以看到參數第一個是points，可以了解為已知的資料點的坐标，第二個參數為值，這裡是已知的點的高程，第三個為你要内插的點的坐标，也就是有着未知高程的點的坐标，第四個參數為方法。對于重采樣，方法有最近鄰元法，雙線性内插法，和三次卷積的方法，對應的名稱分别為"nearest",“linear”,“cubic”，對于不同的要求選擇不同的方法，三次卷積運作速度最慢，但是對于平滑的曲線内插的效果最好。
from scipy.interpolate import griddata
zi=griddata(points,h,(xi,yi),method='cubic')
得到重采樣後的值之後再畫等高線即可。
import matplotlib.pyplot as plt
plt.contour(xi,yi,zi)
plt.savefig(r'C:\Users\58381\Desktop\1.pdf') #把得到的圖儲存為pdf的格式
效果如下
檢驗
大家都知道由TIN生成等高線的方法，先是建構delaunay三角形，然後再在三角形的邊上内插等值點，然後追蹤等值點畫出等高線。下面我們就由已知點建構delaunay三角形。如圖所示
（藍色的點為已知點）可以看到TIN網建構出來了，然後再根據等高差得到内插的等值點（紅色）
然後我們可以發現，之前生成的等高線幾乎可以穿過紅色的點。
由于是通過格網資料生成的等高線，資料精度會下降。才疏學淺，歡迎批評指正！
寫在最後
其實這個方法生成格網的方式太過于粗暴，是以最後拟合出來的等高線并不完美。希望能有人告訴更好的方法，由TIN生成GRID的，感激不盡！