1. はじめに
皆さんは、テクスチャ付きの3D地形データを簡単に、しかも世界中で取得したいと思ったことはありませんか？実際、地理院地図は日本国内向けに3Dモデル（OBJやSTL形式）とPNG画像のペアでテクスチャ付きの3D地形データを提供していますが、データが日本国内に限定されています。一方で、世界中のデータを収集しようとすると、それぞれの衛星データの解像度や座標参照系が合わなかったり、手動でデータを集める必要があったりと、さまざまな問題があります（例: USGS EarthExplorer）。そこでGoogle Earth Engine (GEE)を利用すれば、効率的に世界中のDEMと衛星画像を整合させることができるのです。今回は、GEEを活用して、デジタル標高モデル（DEM）とSentinel-2衛星画像を統合的に取得・処理・可視化・エクスポートする方法をご紹介します。GEEを使うことで、以下のようなメリットがあります：
- 自動化されたデータ収集
- 高精度なデータ整合
- 広範囲なデータカバレッジ
- クラウドベースの処理
2. 実装
2-1 データセット
今回使用したデータセットは以下になります：
-
デジタル標高モデル (DEM):
- USGS/SRTMGL1_003
- 解像度: 約30メートル
- NASAのシャトルレーダートポグラフィミッション（SRTM）によって取得された全球の標高データ
- USGS/SRTMGL1_003
-
衛星画像:
- COPERNICUS/S2_SR
- 解像度: 10メートル（特定のバンド）
- 欧州宇宙機関（ESA）のSentinel-2衛星による表面反射率データ（Level-2A）
- COPERNICUS/S2_SR
2-2 コードの実装
今回使用するスクリプトは以下のGitHubリポジトリから確認・利用できます：
2-3 スクリプトの説明
GEE Python APIを使ってDEMとSentinel-2の衛星画像を統合・処理・可視化・エクスポートする方法を詳しく解説します。実装もGEE APIが便利な関数をたくさん用意してくれてるおかげで簡単にできます。
2-3-1 衛星画像データの取得
まず、Sentinel-2の画像コレクションを取得し、必要な条件でフィルタリングします。
import ee
import geemap
# 認証と初期化
ee.Authenticate()
ee.Initialize(project='google_cloud_project_name')
# エリアの定義（座標は適宜変更してください）
geometry = ee.Geometry.Rectangle([-122.5, 37.0, -121.5, 38.0])
# Sentinel-2データの取得とフィルタリング
sentinel2 = ee.ImageCollection('COPERNICUS/S2_SR') \
.filterBounds(geometry) \
.filterDate('2020-01-01', '2022-12-31') \
.filter(ee.Filter.lt('CLOUDY_PIXEL_PERCENTAGE', 20)) \
.select(['B2', 'B3', 'B4', 'B8', 'SCL']) # クラウドマスキング用にSCLバンドを含む
-
ee.Authenticate()でGEEに認証し、ee.Initialize()でプロジェクトを初期化 -
ee.Geometry.Rectangleで解析対象の地域を矩形で指定 -
ee.ImageCollection()で衛星画像データを取得-
'COPERNICUS/S2_SR'：Sentinel-2のSurface Reflectanceデータを指定 -
.filterBounds(geometry)：定義したエリアに基づいて画像をフィルタリング -
.filterDate('2020-01-01', '2022-12-31')：指定した期間内の画像を取得 -
.filter(ee.Filter.lt('CLOUDY_PIXEL_PERCENTAGE', 20))：雲カバー率が20%未満の画像に絞り込み - .select(['B2', 'B3', 'B4', 'B8', 'SCL'])：必要なバンドを選択('SCL'バンドはクラウドマスキング用)
-
2-3-2 クラウドマスキング
# クラウドマスキング関数の定義
def maskS2clouds_L2A(image):
scl = image.select('SCL')
# クラウドとクラウドシャドウをマスク
mask = scl.neq(3).And(scl.neq(9)).And(scl.neq(8)).And(scl.neq(7))
return image.updateMask(mask)
# クラウドマスキングの適用
sentinel2 = sentinel2.map(maskS2clouds_L2A)
-
scl = image.select('SCL')：シーン分類レイヤー（SCL）を選択 - マスクの作成
-
scl.neq(3)：クラウドシャドウ（値3）でないピクセル -
scl.neq(9)：高確率の雲（値9）でないピクセル -
scl.neq(8)：中確率の雲（値8）でないピクセル -
scl.neq(7)：未分類（値7）でないピクセル - これらをAnd演算子で結合し、すべての条件を満たすピクセルのみを保持
-
-
image.updateMask(mask)でマスクを画像に適用し、雲や影のあるピクセルを透明にします -
.map(maskS2clouds_L2A)で画像コレクション全体に関数を適用
2-3-3 統合衛星画像の作成
クラウドマスクされた画像コレクションから、メディアン（中央値）を計算して一つの合成画像を作成します。
# 画像コレクションのメディアン合成
s2_median = sentinel2.median().clip(geometry)
-
.median()で各ピクセル位置において、画像コレクション内の値の中央値を計算 -
.clip(geometry)で成画像をエリアの範囲内にクリッピング
2-3-4 DEMのリサンプリングと整合
DEMデータをSentinel-2の解像度と座標系に合わせて再サンプリングします。
# DEMの再サンプリング
dem_resampled = dem.reproject(
crs=s2_median.projection(),
scale=10
)
-
.reproject()で画像の座標参照系（CRS）と解像度（スケール）を変更-
crs=s2_median.projection()：Sentinel-2画像と同じCRSを設定。これにより、DEMとSentinel-2画像が正確に重ね合わせられます。 -
scale=10：解像度を10メートル(Sentinel-2の高解像度バンドと同じ)に設定。
-
2-3-4 データのエクスポート
DEMと衛星画像それぞれをtiff(タグ付き画像)ファイルとしてGoogle Driveにエクスポートします。
# DEMのエクスポート
dem_export = ee.batch.Export.image.toDrive(
image=dem_resampled,
description='DEM_export',
folder='EarthEngineExports',
fileNamePrefix='DEM_export',
scale=10,
region=geometry.getInfo()['coordinates'],
crs=s2_median.projection(),
maxPixels=1e13
)
dem_export.start()
# Sentinel-2画像のエクスポート
s2_export = ee.batch.Export.image.toDrive(
image=s2_median.select(['B2', 'B3', 'B4']), # 必要なバンドのみを選択
description='Sentinel2_export',
folder='EarthEngineExports',
fileNamePrefix='Sentinel2_export',
scale=10,
region=geometry.getInfo()['coordinates'],
crs=s2_median.projection(),
maxPixels=1e13
)
s2_export.start()
-
.batch.Export.image.toDrive()でEarth Engineの画像をGoogleドライブにエクスポートする関数-
image：エクスポートする画像オブジェクト。DEMとSentinel-2それぞれで異なります -
description：エクスポートタスクの説明 -
folder：エクスポート先のGoogleドライブ内のフォルダ名 -
fileNamePrefix：エクスポートされるファイル名のプレフィックス -
scale：エクスポートする画像の解像度(ここでは10メートルに設定) -
region：エクスポートする領域の座標(geometry.getInfo()['coordinates']でAOIの座標情報を取得) -
crs：エクスポートする画像の座標参照系をSentinel-2画像と合わせる -
maxPixels：エクスポート可能な最大ピクセル数を設定
-
- エクスポートの開始：
-
.start()：エクスポートタスクを開始(タスクの進行状況はEarth EngineのTaskタブで確認できます)
-
2-3-5 データの可視化
geemapライブラリを使ってデータを地図上で可視化します。
# geemap地図の作成
Map = geemap.Map(center=[37.5, -122.0], zoom=10)
# DEMの可視化パラメータ
dem_vis_params = {
'min': 0,
'max': 3000,
'palette': ['blue', 'green', 'yellow', 'orange', 'red']
}
Map.addLayer(dem_resampled, dem_vis_params, 'DEM')
# Sentinel-2の可視化パラメータ
s2_vis_params = {
'bands': ['B4', 'B3', 'B2'], # 真カラー合成
'min': 0,
'max': 3000
}
Map.addLayer(s2_median, s2_vis_params, 'Sentinel-2')
# マップをHTMLファイルとして保存
Map.addLayerControl()
Map.to_html('map.html')
-
geemap.Mapでgeemap地図オブジェクトを作成します-
center：地図の中心点を設定（今回は緯度37.5、経度-122.0（カリフォルニア州中央部付近）に設定） -
zoom：ズームレベルを10に設定
-
- DEMパラメーターの設定：
-
dem_vis_params：-
'min'：表示する最小値（0メートル） -
'max'：表示する最大値（3000メートル） -
'palette'：色のグラデーションを指定
-
-
- Map.addLayer()でDEMを地図に追加
- ラベル：'DEM'
- 衛星画像パラメーターの設定：
-
s2_vis_params：-
'bands'：真カラー合成のために赤（B4）、緑（B3）、青（B2）のバンドを指定 -
'min'：反射率の最小値（0） -
'max'：反射率の最大値（3000）
-
-
-
Map.addLayer()でメディアン合成されたSentinel-2画像を地図に追加- ラベル：
'Sentinel-2'
- ラベル：
- マップの保存：
-
Map.addLayerControl()で地図上にレイヤーコントロール機能（レイヤーの表示/非表示を切り替えるパネル）を追加 -
Map.to_html('map.html')で作成した地図をインタラクティブなHTMLファイルとして保存
-
3. 実行結果
スクリプトを実行すると、map.htmlというファイルが生成されます。このファイルをウェブブラウザで開くと、以下のような地図が表示されます。DEMと衛星画像が正確に重なっていることが確認できます！
4. まとめ
今回紹介したスクリプトを使えば、Google Earth Engine を活用して、世界中のDEMとSentinel-2衛星画像を効率的に取得・処理・可視化・エクスポートすることができます。ありそうで意外と出回っていなかったものなので皆さんにもぜひ活用していただきたいです！
参考文献
https://earthengine.google.com/
https://qiita.com/ksasao/items/7718234c9cbc19abd3bd