GIS奮闘記

現役GISエンジニアの技術紹介ブログ。主にPythonを使用。

Pythonで自分だけの世界地図を作ってみる ~オープンデータを自分の欲しい形にしてみよう~

本日は世界地図を作ってみようと思います。オープンデータ関連で検索すれば世界地図がダウンロードできるサイトがいくつか検索できます。ただし、すべて海外のサイトのため、表記が英語になっていて非常に見づらいです。なので、今回は日本人向けの世界地図を作ってみようと思います。

※今回使用したソースと作成したGDBGitHubに公開しています。
GitHub - sanvarie/MinnanoArcGIS

ダウンロード

このサイトからShapeファイルをダウンロードします。
http://www.naturalearthdata.com/

ダウンロードページ>Large scale data, 1:10m>Cultural からダウンロードします。
f:id:sanvarie:20160228162730p:plain

こんな感じの表記なので、うーん、といった感じでしょうか。ただ、グルジアがジョージアになっていたりと、データ自体は新しそうです。
f:id:sanvarie:20160228162824p:plain

実行環境

事前準備

  • 「Map.gdb」を作成します。

ポリゴンを格納するフィーチャクラスを作成する

ダウンロードした「ne_10m_admin_0_countries.shp」を「Map.gdb」にコピーします。

サンプルコード①

ダウンロードした「ne_10m_admin_0_countries.shp」を「Map.gdb」にコピーして、不要なカラムを削除するサンプルです。また、正式名称がないフィーチャに関しては名称で置き換えます(クック諸島とか)。

■copyWorldMap.py

# -*- coding: utf-8 -*-
import arcpy

#対象GDB
arcpy.env.workspace = "C:\ArcPySample\Map.gdb"

#コピーするフィーチャクラス名
inFeature = "WorldMap"

#ShapeをMap.gdbにコピー
arcpy.CopyFeatures_management(r"D:\python\WorldMap\ne_10m_admin_0_countries.shp", inFeature)

#コードブロック
codeblock = """
def UpdateColumn(FORMAL_EN,NAME_LONG):
    if FORMAL_EN == " ":
        return NAME_LONG
    else:
        return FORMAL_EN
"""

# 条件式を設定
expression = "UpdateColumn(!FORMAL_EN!,!NAME_LONG!)"

#フィールド演算(正式名称が半角スペースになっているものをNAME_LONGに置き換える)
arcpy.CalculateField_management(inFeature, "FORMAL_EN", expression, "PYTHON_9.3", codeblock)

#カラム追加
arcpy.AddField_management(inFeature, "NAME_EN", "Text",field_length = 60,field_alias="英語名")

#フィールド演算(FORMAL_ENをNAME_ENに置き換える)
arcpy.CalculateField_management(inFeature, "NAME_EN", "!FORMAL_EN!", "PYTHON_9.3")

field_list = []
fields = arcpy.ListFields(inFeature)
for field in fields:
    if field.type != "Geometry":
        if field.name != "NAME_EN" and field.name != "OBJECTID" and field.name != "Shape_Length" and field.name != "Shape_Area":
            field_list.append(field.name)

#不要なフィールドを削除
arcpy.DeleteField_management(inFeature,
                             field_list)

こんな感じになれば成功です。
f:id:sanvarie:20160229111608p:plain

世界地図に属性を付与する

属性カラムがすっきりしたところで、新たに属性を付与させようと思います。

対象サイト

外務省のサイトがよさげだったので今回はこのサイトから情報を取得しようと思います。世界の国々 | 外務省
ただし、世界と日本のデータを見る(世界の国の数,国連加盟国数,日本の大使館数など) | 外務省で確認すると、「世界の国の数=196か国。これは,現在,日本が承認している国の数である195か国に日本を加えた数です。」とあるようにフィーチャ数(255)と一致しません。残念ですが、一致するものだけに属性を付与したいと思います。

また、「独立年」がうまくとってこれないので、ここはカットしようと思います(泣)。Pandasでスクレイピングする限界ですかね・・・

サンプルコード②

必要なカラムを追加するサンプルです。

■addColumn.py

# -*- coding: utf-8 -*-
import arcpy

#対象GDB
arcpy.env.workspace = "C:\ArcPySample\Map.gdb"

#コピーするフィーチャクラス名
inFeature = "WorldMap"

#カラム追加
arcpy.AddField_management(inFeature, "NAME", "Text",field_length = 50,field_alias="国名")
arcpy.AddField_management(inFeature, "CAPITAL", "Text",field_length = 50,field_alias="首都")
arcpy.AddField_management(inFeature, "LANGUAGE", "Text",field_length = 50,field_alias="主要言語")
arcpy.AddField_management(inFeature, "AREASQUARE", "Double",field_length = 50,field_alias="面積(1,000平方キロ)")
arcpy.AddField_management(inFeature, "POPULATION", "Double",field_length = 50,field_alias="人口(100万人)")
arcpy.AddField_management(inFeature, "CURRENCY", "Text",field_length = 50,field_alias="通貨単位")
サンプルコード③

上述したサイトから各国の情報を取得して「WorldMap」フィーチャクラスに属性を付与するサンプルです。ものすごい力技になってしまいましたが。。。

■setCountryData.py

# -*- coding: utf-8 -*-
import arcpy
import pandas as pd

def updateAttribute(row,cur,conTable):
    row.setValue("NAME", conTable[0])
    row.setValue("CAPITAL", conTable[2])
    row.setValue("LANGUAGE", conTable[4])

    AreaSquere = conTable[5]
    #変な文字列がまじっているので
    if isinstance(AreaSquere, float):
        row.setValue("AREASQUARE", AreaSquere)
    elif isinstance(AreaSquere, long):
        row.setValue("AREASQUARE", AreaSquere)
    #なぜかdatetimeとして認識されるものがあるのでそれは無視
    elif isinstance(AreaSquere, datetime.date):
        pass
    else:
        if AreaSquere.find(u"平方キロ") > -1:
            row.setValue("AREASQUARE", AreaSquere[0:AreaSquere.find(u"平方キロ")])

    population = conTable[6]
    #変な文字列がまじっているので
    if isinstance(population, float):
        row.setValue("POPULATION", population)
    else:
        if population.find(u"約") > -1:
            if population.find(u"人") > -1:
                population = population[population.find(u"約")+1:]
                row.setValue("POPULATION", population[0:population.find(u"人")])
            else:
                row.setValue("POPULATION", population[population.find(u"約")+1:])
        elif population.find(u"人") > -1:
            row.setValue("POPULATION", population[0:population.find(u"人")])

    currency = conTable[7]
    if currency.find(u"(") > 1:
        row.setValue("CURRENCY", currency[0:currency.find(u"(")])
    else:
        row.setValue("CURRENCY", currency)

    cur.updateRow(row)

tableList = []
countryList = []

#対象のURL
countryList.append("http://www.mofa.go.jp/mofaj/kids/ichiran/i_europe.html")
countryList.append("http://www.mofa.go.jp/mofaj/kids/ichiran/i_africa.html")
countryList.append("http://www.mofa.go.jp/mofaj/kids/ichiran/i_chuto.html")
countryList.append("http://www.mofa.go.jp/mofaj/kids/ichiran/i_asia.html")
countryList.append("http://www.mofa.go.jp/mofaj/kids/ichiran/i_oceania.html")
countryList.append("http://www.mofa.go.jp/mofaj/kids/ichiran/i_n_america.html")
countryList.append("http://www.mofa.go.jp/mofaj/kids/ichiran/i_latinamerica.html")

#対象のフィーチャクラス
arcpy.env.workspace = "C:\ArcPySample\Map.gdb"
spatial_reference=arcpy.SpatialReference(4326)

for l in countryList:
    #HTMLを読込
    df = pd.io.html.read_html(l)
    tableList.append(df[0])

#データフレームを結合
conTable = pd.concat(tableList, ignore_index=True)
conTable.columns = ["NAME","ENGLISH_NAME","CAPITAL","INDE_YEAR","LANGUAGE","AREASQUARE","POPULATION","CURRENCY"]

for i in conTable.index:
    dfEnglishName = conTable.ix[i][1]
    cursor = arcpy.UpdateCursor("WorldMap")
    for row in cursor:
        if dfEnglishName == row.NAME_EN:
            updateAttribute(row,cursor,conTable.ix[i])
        #色々力技
        elif dfEnglishName.find("of") > -1:
            if dfEnglishName.find("Australia") > -1 and row.NAME_EN.find("Australia") > -1:
                updateAttribute(row,cursor,conTable.ix[i])
            elif dfEnglishName.find("Vatican") > -1 and row.NAME_EN.find("Vatican") > -1:
                updateAttribute(row,cursor,conTable.ix[i])
            elif dfEnglishName.find("Hungary") > -1 and row.NAME_EN.find("Hungary") > -1:
                updateAttribute(row,cursor,conTable.ix[i])
            elif dfEnglishName.find("Verde") > -1 and row.NAME_EN.find("Verde") > -1:
                updateAttribute(row,cursor,conTable.ix[i])
            elif dfEnglishName.find("Gambia") > -1 and row.NAME_EN.find("Gambia") > -1:
                updateAttribute(row,cursor,conTable.ix[i])
            elif dfEnglishName.find("Comoros") > -1 and row.NAME_EN.find("Comoros") > -1:
                updateAttribute(row,cursor,conTable.ix[i])
            elif dfEnglishName.find("Principe") > -1 and row.NAME_EN.find("Principe") > -1:
                updateAttribute(row,cursor,conTable.ix[i])
            elif dfEnglishName.find("Nepal") > -1 and row.NAME_EN.find("Nepal") > -1:
                updateAttribute(row,cursor,conTable.ix[i])
            elif dfEnglishName.find("Viet") > -1 and row.NAME_EN.find("Viet") > -1:
                updateAttribute(row,cursor,conTable.ix[i])
            elif dfEnglishName.find("Bahamas") > -1 and row.NAME_EN.find("Bahamas") > -1:
                updateAttribute(row,cursor,conTable.ix[i])
            elif dfEnglishName == "Republic of Guinea" and row.NAME_EN == "Republic of Guinea":
                updateAttribute(row,cursor,conTable.ix[i])
            elif dfEnglishName == "Independent State of Papua New Guinea" and row.NAME_EN == "Independent State of Papua New Guinea":
                updateAttribute(row,cursor,conTable.ix[i])
            elif dfEnglishName.find("Cote") > -1  and row.NAME_EN == "Republic of Ivory Coast":
                updateAttribute(row,cursor,conTable.ix[i])
            elif dfEnglishName == "Republic of Congo" and row.NAME_EN == "Republic of Congo":
                updateAttribute(row,cursor,conTable.ix[i])
            elif dfEnglishName == "Democratic Republic of the Congo" and row.NAME_EN == "Democratic Republic of the Congo":
                updateAttribute(row,cursor,conTable.ix[i])
            elif dfEnglishName.find("Togo") > -1 and row.NAME_EN == "Togolese Republic":
                updateAttribute(row,cursor,conTable.ix[i])
            elif dfEnglishName.find("Lebanon") > -1 and row.NAME_EN == "Lebanese Republic":
                updateAttribute(row,cursor,conTable.ix[i])
            elif dfEnglishName == "Republic of Korea" and row.NAME_EN == "Republic of Korea":
                updateAttribute(row,cursor,conTable.ix[i])
            elif dfEnglishName.find(row.NAME_EN[row.NAME_EN.find("of")+3:]) > -1 and row.NAME_EN.find("of") > -1 \
                 and dfEnglishName != "Republic of Korea" and dfEnglishName != "Republic of Congo" and dfEnglishName != "Independent State of Papua New Guinea" \
                 and dfEnglishName != "Republic of Equatorial Guinea" and dfEnglishName != "Republic of Guinea-Bissau" and dfEnglishName != "Democratic Republic of the Congo":
                updateAttribute(row,cursor,conTable.ix[i])
        elif dfEnglishName.find("Spain") > -1 and row.NAME_EN.find("Spain") > -1:
                updateAttribute(row,cursor,conTable.ix[i])
        elif dfEnglishName.find("Brunei") > -1 and row.NAME_EN.find("Brunei") > -1:
            updateAttribute(row,cursor,conTable.ix[i])
        elif dfEnglishName.find("Nevis") > -1 and row.NAME_EN.find("Nevis") > -1:
            updateAttribute(row,cursor,conTable.ix[i])
        elif row.NAME_EN.find(dfEnglishName) > -1 and row.NAME_EN.find("Kingdom") == -1  and dfEnglishName.find("of") == -1 \
             and row.NAME_EN != "South Georgia and South Sandwich Islands" and row.NAME_EN != "Indian Ocean Territories" and row.NAME_EN != "British Indian Ocean Territory":
            updateAttribute(row,cursor,conTable.ix[i])

こんな感じになりました。面積や人口がうまく入らないところもあったのですが、とりあえずよしとします。。。
f:id:sanvarie:20160302100327p:plain

サンプルコード④

Shape_Length、Shape_Areaの位置が気持ち悪いのでフィーチャクラスをコピーします。

■copyFeatures.py

# -*- coding: utf-8 -*-
import arcpy

#対象GDB
arcpy.env.workspace = "C:\ArcPySample\Map.gdb"

#ShapeをMap.gdbにコピー
arcpy.CopyFeatures_management("WorldMap", "WorldMap2")

OKですね。後は「WolrMap」フィーチャクラスを削除して「WolrMap2」を「WolrMap」にリネームしてもOKだと思います。
f:id:sanvarie:20160302101018p:plain

課題

  1. 「独立年」がとってこれなかった
  2. 「面積」「人口」が一部とってこれなかった
  3. 承認していない国の情報がない
  4. プログラムが力技すぎ

課題は多いです。ですが、とりあえず、日本が承認している国の日本語名を入れることができたので最低限のことはできたかと思います。台湾がないのは中国との政治的なあれですかね。今後、この世界地図を使って色々情報を追加していきたいですね。「サッカー日本代表と各国との対戦成績を分析」みたいに。

以上、今回は「Pythonで自分だけの世界地図を作ってみる ~オープンデータを自分の欲しい形にしてみよう~」でした。