基于開源GIS的活動斷層矢量數據融合方案初探

中國地震局地質研究所 中國地震局活動構造與火山重點實驗室 吳熙彥

基于開源GIS的活動斷層矢量數據融合方案初探

中國地震局地質研究所 中國地震局活動構造與火山重點實驗室 吳熙彥

活動斷層是確定區域內最大發震震級和破壞強度的重要資料。隨著計算機技術的發展，在國家支持下，地震活動斷層調查相關項目也獲得了大量活動斷層調查空間數據。將這些數據進行綜合，提取關鍵信息，形成一套格式、結構統一的數據，有利于提高效率和進一步的數據分析。

矢量數據融合；活動斷層；數據；pyshp

一、全國活動斷層數據融合意義

中國位于世界兩大地震帶——環太平洋地震帶與歐亞地震帶之間，地震頻發，造成的災害損失巨大。然而，地震短臨預報仍是世界難題，短期內無法攻克。目前，國際和國內防震減災的一項主要措施是通過確定某區域可能發生地震時造成的破壞強度，設定抗震設防烈度，規定建構筑物要達到的抗震要求，減少地震災害造成的損失?；顒訑鄬觿t是確定區域內最大發震震級和破壞強度的重要資料。隨著計算機技術的發展，在國家支持下，地震活動斷層調查相關項目也獲得了大量活動斷層調查空間數據。這些數據獲取的時間跨度大，來源不一。因此，數據的格式、結構和制作標準等并不一致。使用時需要調用和查看不同的數據集，調用不同的圖層，效率較低。同類數據分層分包存儲，也不利于空間查詢、統計等分析研究。將這些數據進行綜合，提取關鍵信息，形成一套格式、結構統一的數據，有利于提高效率和進一步的數據分析。

二、全國活動斷層數據融合流程設計

不同時期、不同來源的活動斷層GIS數據融合的內容包括空間圖形和屬性值[1]。將數據融合的思路分解為圖形的融合和屬性值的融合，以先圖形后屬性的方式進行流程設計。數據融合的步驟如下。

（1）轉換文件格式：原有數據采用了不同的文件格式。數據處理的第一步，是將所有數據統一轉換為shape file文件。Shapefile文件是被廣泛接受的一種GIS數據文件。常用的GIS軟件都能將GIS數據轉換為這種格式。

（2）轉換空間坐標系：不同來源的GIS數據可能在不同的空間坐標系下，應在處理之前轉換到統一的坐標系下。目前，我國最新的坐標系統為2000國家大地坐標系，它是國家基礎地理信息數據（基礎地理底圖）所采用的坐標系，由國務院批準自2008年7月1日啟用。為推廣2000國家大地坐標系的使用，國家測繪地理信息局公布了現有測繪成果轉換到該坐標系的技術指南。依據技術指南，可以將現存數據都轉換到2000國家大地坐標系。數據轉換到統一坐標系后，可實現多來源空間數據疊加，進行編輯和修改。

（3）分區處理數據：空間數據疊加后，可以看到部分區域內只有一層圖層，每個對象只有一種圖形表達，而部分區域有一層以上的重疊數據。對于只有唯一空間數據的區域，不需進行空間操作，直接保留。對于重疊區域，須進行圖形合并。因此，在圖形合并前，需將將數據集切塊，切出重疊區的部分待下一步處理。

（4）重疊區圖形數據合并：地圖合并技術一般包括匹配實體和將匹配的同名實體合并兩個過程[2][3]。實體自動匹配主要依靠幾何相似度、拓撲、數據模型和屬性等。同名實體自動合并主要根據預設權重和待合并圖形幾何形態生成新的圖形。然而，以線要素表達的活動斷層線的幾何形態并非地質對象的真實形態，即地殼斷裂帶的真實幾何形態，而是地質學家對地質對象的一種解釋和抽象。同一實體的不同解釋和抽象可能使幾何形態相似程度較低。因此，活動斷層的幾何形態的，并不具備很好地自動匹配和自動合并的條件。因此，采用以下三步合并數據。1）選取重疊區某一圖層為最終保留層（稱基準圖層），用人機交互的方式修改該圖層；2）為每一圖層增加“同名實體編號”字段。所有同名的實體賦予唯一的一個編號填入該字段；3）利用“同名實體編號”連接其他數據屬性表到基準圖層，將同名實體的屬性全部合并入基準圖層。最終重疊區僅保留基準圖層一層數據。

（5）確定屬性表結構：合并重疊區和非重疊區的空間數據，形成一層新的圖層。新數據的屬性表是原有數據所有屬性的集合。新屬性表存在重復信息，需要重新設計與規劃。根據屬性內容，重新設計屬性表結構，在不損失信息條件下精簡表結構。

（6）屬性融合：根據重新設計的表結構，為屬性表新增字段，準備填入信息融合后的屬性值。屬性處理有3種情況。1）屬性值直接賦給新表中對應屬性字段；2）多個信息重復的屬性值，合并到新表中的一個屬性字段；3）多個信息不同的屬性值，合并到新表中的一個屬性字段。根據新形成的數據和新設計的表結構，設計屬性遷移和合并規則。

（7）刪除冗余字段：將上一步形成的新數據中的冗余字段刪除，獲得最終成果數據。

圖1 流程圖

三、方案可行性分析

1.pyshp介紹

方案的第1-5步都較容易實現。方案實現的關鍵是第6步屬性融合。筆者采用了pyshp作為實現第6步的試驗工具。pyshp全稱Python Shapefile Library，是基于python語言的開源函數庫，提供了讀寫ESRI shapefile格式文件的方法。pyshp函數庫具有編譯工具簡單、免安裝、代碼開源、使用簡單的特點，并且提供了簡單易懂的使用說明和示例代碼，在ESRI shapefile編輯方面便利，非常適合地理數據融合工具。pyshp提供三個基礎類:Reader、Writer、Editor。Reader 提供了ESRI shapefile的文件頭信息、幾何圖形和屬性表的讀取方法；Writer提供了ESRI shapefile的文件頭信息、幾何圖形和屬性表的寫入方法；Editor是Writer子類，繼承了Writer的所有方法，主要用于已有文件的編輯。pyshp可以在任何pyshon平臺直接使用import命令導入模塊文件調用函數。本文使用的是pyshp 1.2.3。

2.關鍵步驟的實現

這里僅討論關鍵的兩種屬性處理。

（1）多個信息重復的屬性值，合并到新表中的一個屬性字段。例如，斷層活動年代的編碼方式不同，一部分代碼值為{a1,a2,a3,...ai}，另一部分代碼值為{b1,b2,b3,...ai}。實際上，a1和b1代表相同屬性。該屬性的新代碼值為{c1,c2,c3,...ci}。此時，設計一個代碼映射表，將代碼一一映射到新代碼集。核心python代碼如下：

iSetField=self.getEsfFieldIndex(esf, SetField)#新字段

iReferField=self.getEsfFieldIndex(esf, ReferField)#舊代碼值所在字段nRecord=0

for xRecord in esf.records:

#如果新代碼字段是空的，就將舊代碼值的對應新代碼填入

if xRecord[iSetField] is None or xRecord[iSetField]==0:

reSValue=str(u(xRecord[iReferField]))

xRecord[iSetField]=int(AgeDic.get(reSValue,0))#從代碼映射表獲得新代碼

nRecord+=1

return nRecord

（2）多個信息不同的屬性值，合并到新表中的一個屬性字段。例如，原數據有一個表示斷層是“出露”還是“隱伏”的字段，還有另一個表示“正斷”還是“逆斷”的字段。兩種信息可以合并，成為新的斷層性質字段。新編代碼的前幾位代表前一個字段值，后幾位代表后一個字段值。核心python代碼如下：

for xRecord in esf.records:

Sfirvalue=str(u(xRecord[iFirField])).strip() #獲得第一個字段值

xFir=NewCodeFirDic.get(Sfirvalue,0) #從代碼映射表獲得第一個字段值的新代碼

SSecvalue=str(u(xRecord[iSecField])).strip() #獲得第二個字段值

xSec=NewCodeFirDic.get(SSecvalue,0) #從代碼映射表獲得第二個字段值的新代碼

xRecord[iSetField]= xFir.join(xSec)

nRecord+=1

return nRecord

[1]唐文靜．多源地理空間矢量數據融合[M]．清華大學出版社,2014．

[2]張橋平,李德仁,龔健雅．地圖合并技術[J]．測繪通報,2001(1)：6-8．

[3]COBB M,CHUNG M,FOLEY H．A Rule-based Approach for the Conflation of Attributed Vector Data．GeoInformatica,1998,2(1)：7-35．

中國大陸主要地震構造帶活動斷層探察專項（IGCEA1621）子課題“郯廬斷裂帶三維地震構造分段模型與地震預測研究——數據庫建設”；中國地震局地震行業科研專項（20138001、201408023、200908001、201008001）。