拂衣天氣（微天氣）— 行政區劃資料（二）

前言

此前提到微天氣應用程式需要使用到行政區劃資料，不過上一章所使用的資料來源於網路，或多或少都可以考慮一下是否還有其他獲取的方式，所以也就有了本文的內容。

在這裡，將基於全國1:100萬基礎地理資訊資料進行行政區劃資料的提取。本文用於記錄使用程式實現全國1:100萬基礎地理資訊資料合併的全過程。

當然，本文產生的最重要原因其實並不是受到微天氣的啟發，更多的是個人想試一試能不能用。

資料說明

全國1：100萬公眾版基礎地理資訊資料（2021）覆蓋全國陸地範圍和包括臺灣島、海南島、釣魚島、南海諸島在內的主要島嶼及其臨近海域，共77幅1:100萬圖幅，該資料集整體現勢性為2019年。資料採用2000國家大地座標系，1985國家高程基準，經緯度座標。

為滿足廣大社會群眾對地理資訊資料資源的需求，經自然資源部授權，全國地理資訊資源目錄服務系統提供全國1:100萬全圖層要素免費下載的服務。下載資料採用1：100萬標準圖幅分發，內容包括水系、居民地及設施、交通、管線、境界與政區、地貌與土質、植被、地名及註記9個資料集，且儲存要素間空間關係和相關屬性資訊。

💡 提供下載的是向量資料，不是最終地圖，與符號化後的地圖再視覺化表達上存在一定差異。使用者利用此資料編制地圖，應當嚴格執行《地圖管理條例》有關規定；編制的地圖如需向社會公開的，還應當依法履行地圖稽核程式。

成果規格

• 分幅編號及範圍

  1：100萬公眾版基礎地理資訊資料（2021）的圖幅總數為77幅，分幅資料按照GB/T 13989-2012《國家基本比例尺地形圖分幅和編號》執行。空間儲存單元為6°（經差）×4°（緯差）。

• 座標系統

  • 平面座標系： 2000國家大地座標系。
  • 高程基準：1985國家高程基準。
  • 地圖投影：分幅資料採用地理座標，座標單位為度。

• 幾何精度

  更新後地物點對於附近野外控制點的平面位置和高程中誤差符合下表的要求，以兩倍中誤差值為最大誤差。

  地物點誤差	最小	最大
  平面位置	100	500
  高程	50	200

• 現勢性

  1:100萬地形資料現勢性與更新使用的資料來源的現勢性一致，資料整體現勢性達到2019年。

成果資料組織

全國1：100萬公眾版地形資料（2021）內容包括水系、居民地及設施、交通、管線、境界與政區、地貌與土質、植被、地名及註記9個資料集。

資料分層的命名採用四個字元，第一個字元代表資料分類，第二三個字元是資料內容的縮寫，第四個字元代表幾何型別。

目標

• 實現分圖層合併（處理圖幅合併時接邊問題）
  • 水系（暫緩）
  • 交通（暫緩）
  • 境界與政區（國、省、市、縣）
  • 地名及註記
• 根據合併後的行政區劃與地名註記，製作行政區劃資料庫
• 資料庫使用PostgreSQL（PostGIS）

設計

圖層解析程式

Java程式， 使用GDAL 讀取GDB

邏輯

根據《國家基本比例尺地形圖分幅和編號》規定可知網格範圍，透過網格範圍動態生成對應網格的分幅編號，並以該編號進行資料檢索。如果命中則根據成果資料組織規格以及相關標準對資料進行解析，如果未命中則跳過，直至網格掃描完畢。

經度：72~138（E），43~53（11）

緯度：0~56（N），A~N（14）

即，網格範圍:[43,53] x [A,N]

資料的處理流程可使用責任鏈模式進行，後續也方便加入其他的處理流程。即，整體的執行框架為策略模式+責任鏈模式。為統一策略選擇模型，在此提出圖層定義（LayerDefinition）的概念。

LayerDefinition由如下幾個關鍵要素組成：

• 圖層資料來源（LayerSource）
  • driver：驅動，參考java.sql.UnWrapper實現
    • instance
    • name
  • catalog：
  • schema：
  • table：
  • commonDefinitionKey ：常規定義快取的key
  • fieldDefinitionKey ：欄位定義快取的key
  • featureCarrierKey ：要素載體快取的key
• origin：資料來源（分幅檔案路徑）
• scale：比例尺（如：1000000，表示1:1000000）
• sourceSpatialRef：源座標系（表現形式可為標準ID、PROJ Text、WKT Text）
• sinkSpatialRef：目標座標系（表現形式可為標準ID、PROJ Text、WKT Text）
• featureCode：要素分類碼，對應成果資料組織中的要素分類（如：C、B）
• name：圖層命名，也是圖層分類碼
• layerCode：圖層分類碼
• release：釋放格式（比如：WMTS、Shapefile、GDB…）

描述字元為：比例尺:源座標系:目標座標系:要素分類碼:圖層名稱:釋放格式，在描述字串中，座標系僅使用標準ID表示

💡 源座標系、要素分類碼、圖層名稱均來自於源資料。當且僅當原始資料中無法獲取到源座標系定義時，源座標系定義可來自於配置。

基本流程

1. 掃描網格定義（可配置）

2. GDB讀取（支援zip，並推薦使用zip）

   /home/fuyi/GeoDatabase/China_Basic_1-100/J50.gdb
   
   /home/fuyi/GeoDatabase/China_Basic_1-100/J50.gdb.zip
   

3. 策略管理器（StrategyManager），職責：根據配置選擇具體策略方案，並呼叫執行

   策略（Strategy），持有相應的處理器物件，責任鏈由此開始啟動

   總體上分為三個階段

   1、圖層歸一（資料、位置、座標系）
   2、圖層合併
   3、圖層釋放（service、dump、nothing）
   
   Layer normalization (data, location, coordinate system)
   Layer merging
   Layer release (service, dump, nothing)
   

   • 歸一化處理階段（應保證空間參考歸一、圖層資料儲存歸一，即應保持整體一致）
   • 合併處理階段（或者稱之為中間處理階段，在這裡不僅可以處理合併，還可以做其他的事情）
   • 釋放階段（資料的處理總歸是有目的存在的，那麼本階段就是確定資料最終存在或服務的形式，比如匯出為Shapefile檔案、寫入PostGIS資料庫或透過WMTS服務對外提供資料訪問等）

   GDAL策略：預設策略，即全程基本上使用GDAL處理

4. 匹配獲取指定圖層處理器（歸一化），職責為：解析

   先寫快取（redis）

   1. support（根據圖層命名），如：boua，該值作為處理器的判斷條件

5. 可提供轉換處理器，職責為：座標系轉換

   座標系轉換既可以在歸一化處理階段進行（保證不同空間參考的資料轉入同一空間參考），也可以在中間處理階段（統一進行資料轉換）或是釋放階段（根據釋放需求進行處理）

6. 提供資料合併處理器，職責為：跨圖幅的圖層合併

   處理邏輯，根據給定的關聯列進行空間物件合併（如：在行政境界中使用行政區劃編碼以及名稱進行空間合併）

   不同圖層可以具有不同的合併規則，由不同的例項進行提供對應的合併行為

7. 釋放邏輯也透過定義格式處理器完成處理，職責：格式處理

   1. 需指明格式與比例尺
   2. support，根據格式與比例尺，如：postgis:1000000，該值作為處理器的判斷條件

8. 不同比例尺建立不同的儲存名稱。如，儲存格式為PostGIS：

   表命名規則：圖層名稱_比例尺，小寫

   如：boua_1000000，意為1：1000000比例尺下boua資料層資料

💡 LayerDefine結構可作為所有策略選擇的入參，內部則根據自己依賴的條件進行決策。

Driver(register) → GDALDriverManager → DataSet → Layer → Feature → Geometry

Redis快取作為轉換過程中的資料儲存與交換容器

圖層定義拆分為兩個部分進行儲存，常規內容使用hash結構儲存，便於進行變更
prefix key : feature code : layer code : definition

• scale
• sourceSpatialRef
• sinkSpatialRef
• storage

屬性定義使用set結構儲存，便於序列化與反序列化

prefix key : feature code : layer code : definition : fields

• fieldDefinitions

圖層資料則使用list資料結構進行儲存，便於序列化與反序列化，同時便於資料追加
prefix key : feature code : layer code : features

• features

結構拆解圖

開發

環境準備

• 作業系統：Debian 11 bookworm（testing）
• IDE：IntelliJ IDEA 2022.1.4 (Community Edition)
• PostgreSQL：14（13+即可）
• JDK：openjdk version "11" 2018-09-25
• Maven：Apache Maven 3.8.6
• GDAL：GDAL 3.5.1, released 2022/06/30

GDAL環境準備

debian 11 testing，自行編譯gdal程式

/debian 11 testing，已經自行安裝了gdal程式，但是沒有java binding模組，所以還是需要自行編譯。

我今天用的是Fedora 38 workstation，可以直接 install gdal, gdal-java，而後將libgdalalljni.so檔案複製到/usr/java/packages/lib 目錄下（如果不存在，請自行建立）

sudo mkdir -p /usr/java/packages/lib & sudo ln -s /usr/lib/java/gdal/libgdalalljni.so /usr/java/packages/lib/libgdalalljni.so

@time: 2023.10.08

異常

下圖是經過pac,name合併後的資料（並去除了境外與海洋資料），如圖可見，仍然存在些許裂縫

經過排查發現，發現部分原始圖副之間便存在裂縫，且該裂縫的兩邊不一定是平行的。那麼這就會導致合併出現兩種情況

1. 待合併的多邊形之間存在裂縫，且多邊形互相不相交

   也就是隔海相望的那一種，此種情況是可以透過後續手段修復的，不過也是隻有在裂縫的定義存在於合併後的多邊形中（只不過可能被定義為洞）的情況下，如果不在的話，也不知道應該如何出來

2. 待合併的多邊形之間存在裂縫，且多邊形存在相交，並且交點在裂縫之上

   這種形式便不知道應該怎麼處理了。

💡 感覺本質上，是一個找到裂縫的過程，並將裂縫作為內容的一部分進行合併。但是如何確定裂縫的範圍也是不容易的，這可能和多邊形間相互位置關係有關係

成果

歸一化

即將分離的圖幅資料提取到統一位置存放，座標本來就是統一的，不需要額外操作

去除無效區域與資料

• 去除境外資料（pac < 100000）
• 去除海域面資料（pac = 250100）
• 去除空間無效資料（ST_IsValid(geometry)==false）

合併處理

即進行圖幅合併（使用pac與name）

💡 從此處可以發現，部分割槽域存在裂縫。
1、pac=632722
2、pac=140404
3、pac=140427
4、pac=520625
5、pac=520628

我同時對比了QGIS、ArcGIS中提供的Dissolve工具的融合效果（使用pac與name），就都存在裂縫情況上：
合併後的數量都是：2900

arcgis 優於 wukong 優於 qgis

ArcGIS

QGIS

WuKong

Source Code

wukong-github

小結

因做境界合併，開發了WuKong程式，但是其含義並不僅僅只是進行境界合併，他是一個基於國家開放的基礎地理資訊資料的（目前是基於1：100玩基礎地理資訊資料），可以進行歸一、轉換、釋放的通用程式。

雖然目前僅完成了基礎境界的合併（且存在瑕疵，但好歹也和QGIS的實現差不多不是），但這僅僅只是開始，還有更多的內容可以實現。

參考

全國地理資訊資源目錄服務系統（https://www.webmap.cn/）

資料來源於全國地理資訊資源目錄服務系統：www.webmap.cn

GB/T 13989-2012 國家基本比例尺地形圖分幅和編號

1：1000000地形圖的分幅

經度：72_{138（E），43}53（11）

緯度：0_{56（N），A}N（14）

GB/T 33183-2016 基礎地理資訊 1:50 000地形要素資料規範

1:50000地形要素屬性表結構

• 境界與政區

• 地名

1:50000地形要素資料要素內容與選取指標

走天涯徐小洋公眾號

2021版國家基礎地理資料庫批次下載與拼接方法

【資源分享】如何查詢靠譜的國標，全文免費看！全文免費看！全文免費看！

1：100萬基礎地理資料庫怎麼看？怎麼用？能幹嘛？

GDAL

vector-openfilegdb

Vector Data Model - GDAL documentation

ESRI File Geodatabase (OpenFileGDB) - GDAL documentation

Java bindings - GDAL documentation

Overview (GDAL/OGR 3.5.0 Java bindings API)

Other

https://github.com/modood/Administrative-divisions-of-China

ST_Centroid

Conda

Conda - Conda documentation

conda的安裝與使用 2.0版（2022-08-12更新）