由經緯度查找鄉鎮市區名稱

 許多應用程式可以顯示使用者所在地的各種資訊，例如筆者住在台北市內湖區，開啟應用程式就顯示內湖區的天氣狀況，讓人倍感親切。其原理是以行動裝置取得使用者所在的經緯度，再由經緯度取得鄉鎮市區行政區域名稱，就可做各種應用了！

由國土測繪圖資取得鄉鎮市區邊界資料

要判斷經緯度屬於哪一個鄉鎮市區，必須先知道鄉鎮市區的邊界資料才能判斷，國土測繪圖資公開資料有提供鄉鎮市區邊界資料。

開啟 https://whgis.nlsc.gov.tw/Opendata/Files.aspx 網頁，按「鄉鎮市區界線(TWD97經緯度)1070516」項目、格式為「SHP」右方下載鈕。

解壓縮後將 <TOWN_MOI_1070516.shp>、<TOWN_MOI_1070516.shx>、<TOWN_MOI_1070516.dbf> 三個檔複製到本專題的 Python 檔案資料夾中。

安裝套件

本專題需要使用 Fiona 及Shapely 套件，但試了許多方法都無法安裝成功。後來找到資料提及這兩個套件是 geopandas 的相關套件，結果安裝了 geopandas 後，Fiona 及Shapely 套件也都安裝完成了！

pip install geopandas

安裝完後最好使用 pip list 檢查 Fiona 及Shapely 套件是否安裝完成。

pip list

邊界資料的資料結構

首先讀取邊界資料：(邊界資料檔與 Python 程式在同一資料夾)

import fiona

from shapely.geometry import shape, Point

import os

module_dir = os.path.dirname(__file__)  #取得目前目錄

collection = fiona.open(os.path.join(module_dir, 'TOWN_MOI_1070516.shp'))

讀取的資料在 collection 變數中。

撰寫程式碼顯示邊界資料：

for f in collection:

    print(f)

下面是其中一筆資料：

{

'type': 'Feature', 

'id': '367', 

'geometry': 

   {'type': 'Polygon', 

    'coordinates': 

      [[(120.30864905600004, 22.65505427000005), 

        (120.308820749, 22.65661072200004), 

        ………

        (120.30864905600004, 22.65505427000005)]]}, 

'properties': OrderedDict([

    ('TOWNID', 'E03'), 

    ('TOWNCODE', '64000030'), 

    ('COUNTYNAME', '高雄市'), 

    ('TOWNNAME', '左營區'), 

    ('TOWNENG', 'Zuoying District'), 

    ('COUNTYID', 'E'), 

    ('COUNTYCODE', '64000')])

}

較重要的資料有：「geometry」項目表示這是一個多邊形，其中有幾百個經緯度畫出多邊形邊界，此多邊形就是鄉鎮市區邊界；「properties」項目的「COUNTYNAME」為縣市名稱，「TOWNNAME」為鄉鎮市區名稱。

完整程式碼

1 def search(x, y):  #尋找鄉鎮

 2     global shapes, townnames

 3     return next((townnames[town_id]  #如果鄉鎮區域包含傳入的經緯度就傳回townnames[town_id]

 4                  for town_id in shapes  #逐一尋找各鄉鎮

 5                  if shapes[town_id].contains(Point(x, y))), None)

 6 

 7 import fiona

 8 from shapely.geometry import shape, Point

 9 import os

10 

11 lng = float(input('輸入經度：'))

12 lat = float(input('輸入緯度：'))

13 module_dir = os.path.dirname(__file__)  #取得目前目錄

14 collection = fiona.open(os.path.join(module_dir, 'TOWN_MOI_1070516.shp'))

15 shapes = {}

16 townnames = {}

17 

18 for f in collection:

19     town_id = f['properties']['TOWNCODE']  #鄉鎮代碼

20     shapes[town_id] = shape(f['geometry'])  #鄉鎮界限經緯度

21     townnames[town_id] = f['properties']['COUNTYNAME'] + ',' + f['properties']['TOWNNAME'] #search函式傳回值

22 

23 print(search(float(lng), float(lat)))

第 7-9 列含入套件。

第 11-12 列讓使用者輸入經度和緯度。

第 13-14 列讀取鄉鎮市區邊界資料。

第 15 列 shapes 串列儲存邊界資料，16 列 townnames 儲存縣市及鄉鎮市區名稱。

第 18-21 列以迴圈建立各種串列資料：19 列建立鄉鎮市區代碼，20 列建立邊界資料，21 列建立縣市及鄉鎮市區名稱。

第 1-5 列 search 函式可由傳入的經緯度尋找鄉鎮市區名稱：4-5 列逐一尋找經緯度是否在鄉鎮市區範圍內，直到找到為止，找到後由 第 3 列傳回縣市及鄉鎮市區名稱。

第 23 列印出傳回的縣市及鄉鎮市區名稱。

執行結果

下載 news5-2.zip