python gdal教學之：幾何形狀geometry與投影projection

建立空的geometry物件：ogr.Geometry

定義各種不同的geometry使用的方法是不一樣的(point, line, polygon, etc)

新建點point，使用方法AddPoint( , , [])。其中的z座標一般是省略的，預設值是0

例如：

point = ogr.Geometry(ogr.wkbPoint)

point.AddPoint(10,20)

line

使用AddPoint(, , [])加入點

使用SetPoint(, , , [])更改點的座標

例如下面這段程式碼，更改了0號點的座標：

line = ogr.Geometry(ogr.wkbLineString)

line.AddPoint(10,10)

line.AddPoint(20,20)

line.SetPoint(,30 )   #(10,10) -> (30,30)

統計所有點的數目

print line.GetPointCount()

讀取0號點的x座標和y坐標

print

讀取0號點的x座標和y坐標

print line.GetX(000line. )

print line.GetY(0)

新建多邊形，首先要新建環(ring)，然後把環加到多邊形物件中。

如何創建一個ring？先新建一個ring對象，然後向裡面逐一增加點。

ring = ogr.Geometry(ogr.wkbLinearRing)

ring.AddPoint(0,0)

ring.AddPoint(100,0)

ring.AddPoint(100,100)

結束的時候，用CloseRings關閉ring，或是將最後一個點的座標設定為與第一個點相同。

ring.CloseRings()

ring.AddPoint(0,0)

下面舉一個例子，建立一個方框。這是個polygon對象，又例外兩層ring構成。

outring = ogr.Geometry(ogr.wkbLinearRing)

outring.AddPoint(0,0)

outring.AddPoint(100,0)

outring.AddPoint(0,0)

inring = ogr.Geometry(ogr.wkbLinearRing)inring = ogr.Geometry(ogr.wkbLinearRing)

inring.AddPoint(25,25)

與

.AddGeometry(inring)

最後三句話比較重要，就是先建立一個polygon對象，然後添加外層ring和內層ring

下面這句話可以幫你數數你的polygon能有幾個ring

print polygon.GetGeometryCount()

從polygon讀取ring，index的順序和建立polygon時加入ring的順序相同

outring = polygon.GetGeometryRef(0)

inring = polygon.GetGeometryRef(0)

inring = polygon.Getome建立複合幾何形狀multi geometry

例如MultiPoint, MultiLineString, MultiPolygon

用AddGeometry把普通的幾何形狀加到複合幾何形狀中，例如：

multipoint = ogr.Geometry(ogr.wkbMultiPoint). Geometry(ogr.wkbPoint)point = ogr.Geometry(ogr.wkbPoint)

point.AddPoint(10,10)

multipoint.AddGeometry(point)

Addpoint.AddPoint(20,20)

. point)

讀取MultiGeometry中的Geometry，方法和從Polygon讀取ring是一樣的，可以說Polygon是一種內建的MultiGeometry。

不要刪除一個已存在的Feature的Geometry，會把python搞崩潰的

只能刪除一個已存在的Feature的Geometry，比方說手工創建出來的，或者調用其他函數自動創建的。就算這個Geometry已經用來創建別的Feature，你還是可以刪除它。

例如：Polygon.Destroy()

關於投影Projections，使用SpatialReference物件

多樣的Projections，GDAL支援WKT, PROJ.4, ESPG, USGS, ESRI.prjections，GDAL支援WKT, PROJ.4, ESPG, USGS, ESRI.prjections，GDAL支援WKT, PROJ.4, ESPG, USGS, ESRI.prjections，GDMrome已讀可從內讀可從寫取Projections，例如：

spatialRef = layer.GetSpatialRef()

spatialRef = geom.GetSpatialReference()

Projection：

首先導入osr庫，之後使用osr.SpatialReference()建立SpatialReference物件

之後用下列語句向SpatialReference物件匯入投影資訊

)  

    •ImportFromProj4()

    •ImportFromESRI()

       •ImportFromUSGS(, )

•ImportFromXML()

導出Projection，使用下面的語句可以導出為字串

    •ExportToWkt()

  •ExportToPCI()

    •ExportToUSGS ()

•ExportToXML()

對一個幾何形狀Geometry進行投影變換，要先初始化兩個Projection，然後建立一個CoordinateTransformation對象，用它來做變換

sourceSR = osr.SpatialReference()

sourceSR.ImportFromEPSG(32612)

targetSR = osr.SpatialReference()

targetSR.ImportFromEPSG(4326) #Geo WGS84

coordTrans = osr.CoordinateTransformation(sourceSR, targetSRgets

（ ！在windows裡面跑不通。老外的論壇裡面有討論，說在linux裡面沒問題，windows死活不行，哎。 。 。

另外還有幾個要注意的地方：

要在適當的時候編輯Geometry，投影變換之後最好就不要再動了吧。

對一個資料來源DataSource裡面的所有Geometry做投影變換，你得一個一個來。用個循環吧。

將你的投影寫入.prj文件，其實很簡單。首先MorphToESRI()，轉成字串，然後開個文字檔往裡面寫就行了。例如：

targetSR.MorphToESRI()

file = open('test.prj', 'w')

file.write(targetSR.ExportToWkt())

ffile.close()教學課程：幾何形狀geometry與投影projection的內容，更多相關內容請關注PHP中文網（www.php.cn）！