Python Google Earth Engine（GEE）を使用した地図上に、指定した区域を描く

　Google Earth Engine（GEE）は、地球の様子を衛星画像によって可視化、分析できる地理空間分析のプラットフォームです。クラウドベースでデータを取得します。
　本記事では、指定した区域を描く雛形コードを載せました。下図例では滋賀県の市町村の境界を色分けして地図上に重ね合わせたものです。JupyterLabを使用してインタラクティブに操作していきます。

■1. ライブラリのインストール

anaconda環境下の場合

conda install -c conda-forge earthengine-api
conda install -c conda-forge geemap

pipの場合

pip install earthengine-api
pip install geemap

■2. gcloud CLI のインストール

　Google Cloud CLI は、Google Cloud リソースを操作するツールです。コマンドラインやスクリプト等の一般的なプラットフォームでタスクを実行できます。

cloud.google.com

上記でインストール後に、JupyterLabにて下図のような手順で使用できるようになります。

import ee
ee.Authenticate()
ee.Initialize()

上記の処理は初回のみです。

■3. 地図を区切るための区域データgeojsonの取得

手順は下記リンクを参照下さい。この緯度と経度の位置座標データから、所望の場所を可視化したり、衛星データを取得することができるようになります。

hk29.hatenablog.jp

■本プログラム

上記1, 2, 3を実施後、JupyterLabを用いて以降のコードをインタラクティブに実行して進めます。

import ee
import geemap
import geopandas as gpd
import pandas as pd
import matplotlib.pyplot as plt
import matplotlib.colors as mcolors
css_colors = mcolors.CSS4_COLORS
import japanize_matplotlib

# GeoJSON形式のファイルを読み込む
geojson_path = 'N03-20210101_25_GML/N03-21_25_210101.geojson' # 滋賀県
gdf = gpd.read_file(geojson_path)
gdf

# 必要な市町村だけ抽出
# target_list = ['大阪市']
# gdf2 = gdf.query('N03_003 in @target_list')
gdf2 = gdf.copy()
gdf2

# 抽出した区域を図示
#gdf2.plot(column='N03_004', legend=False, cmap='OrRd', figsize=(6, 6))
ax = gdf2.boundary.plot()
ax.set_axis_off()

# 抽出したデータの市町村名一覧
my_group = gdf2.groupby("N03_004").groups.keys()
print(len(my_group))
my_group

# 抽出した地図領域データをgeojson形式で保存
data_file = 'shiga_ken.geojson'
gdf2.to_file(data_file, driver='GeoJSON') # geojson形式

# 市町村データからポリゴンデータを抽出する
import json
with open(data_file, encoding='utf-8') as f:
    data = json.load(f)

# 市町村名と座標の辞書を作成する
city_dict = {}
for feature in data['features']:
    geometry = feature['geometry']
    city_name = feature['properties']['N03_004'] # 市町村名のキーはデータによって異なる場合がある
    if city_name not in city_dict:
        city_dict[city_name] = [] # 市町村名が辞書になければ空のリストを追加する
    if geometry['type'] == 'Polygon':
        city_dict[city_name].append(geometry['coordinates'])
    elif geometry['type'] == 'MultiPolygon':
        city_dict[city_name].extend(geometry['coordinates'])

# 市町村名と座標の辞書を表示する
#print(city_dict)

# GEE初期化
geemap.ee_initialize()

# データセットの選択
dataset = ee.ImageCollection("USGS/NLCD")

# データセットから一つの画像を取得
image = dataset.first()

# Mapオブジェクトの作成
Map = geemap.Map()

# 市町村別の色で領域を枠で囲む（画像のレイヤーの前に追加する）
for i, (key, val) in enumerate(city_dict.items(), start=1):
    print(i, key)
    coords = val[0][0]
    Map.addLayer(ee.Geometry.Polygon(coords), 
                 {'color': list(css_colors.keys())[i*2-1]}, # 市町村名に対応する色を選択する
                 key,
                )

# データの範囲を指定して画像をクリップしてマップに追加する
for key, val in city_dict.items():
    #print(key)
    coords = val[0][0]
    clipped_image = image.clip(ee.Geometry.Polygon(coords))
    Map.addLayer(clipped_image, {'min': 0, 'max': 100, 'opacity': 0.1}, 'NLCD Land Cover', True) # opacityは小さいほど透明度が高い

# マップを表示
Map