海老ノート

Google Earth Engine 苦闘の記録

石川チャレンジ2:なかなか始まらない

そう簡単にはいかない

いや、わかってましたよ。そんなサクサクはいかないだろうなーと。でも思った以上に考えることやらケースが多くて、勉強・経験不足が露呈します。どこらへんがなのか、というのは追々として、これまでのあらすじを書いてみましょうか。

作業手順

  1. QGISで植生図のカテゴリをリクラス
  2. QGISでランダムにポイント発生させて、ポイントに植生クラスを付与、シェープファイルに書き出し
  3. シェープファイルをGEEに読み込む
  4. GEEで特徴量として使う衛星データを整理して、読み込んだシェープファイルにそのデータをくっつける
  5. シェープファイルの各点のデータをcsvで書き出し。ダウンロードしてRでランダムフォレストの設定を調整
  6. 調整した数値などを使って、GEE上でランダムフォレスト。完成。

というわけで、結構手数があるわけです。で、5まで行ったと思ったら違ってたから1に戻るとかですね、そういうのがたびたび起きるわけで、あちゃーとなるわけです。モチベーションは下がるけど、一度通った道は早く感じるというのか、スクリプトはあるし思ったより楽なんだけど、モチベーションは下がる。慣れていれば落とし穴やらひっかけ問題にかからず、サクッと行くんでしょうが…。

GEEでの衛星データの整理

QGIS周りの作業は大胆に割愛。シェープファイルのアップロード法はこちらの説明の通り。

Sentinel2の前処理

GEEのサイトを見たらpythonの雲除去コードがあったので、こちらを借用。これをjavascriptに書き換えます(っつってもdefをvarに変えるとかそれくらい)。最後に雲を除去する関数を実行して、そこからNDVIとNDWI、各バンドのmedian、10分位、90分位を計算するようにしています。日付を微妙に設定しているのは、冬は雪で真っ白かもしれないから、そうじゃない時期だけにしたいから↓

https://code.earthengine.google.com/eccc6a6b863ab6ea8c3522cdc085f4dc

地形と夜間光も使ってみたい

地形データはJAXAALOSから、最近流行ってるっぽい夜間光データはNOAAのものを使いました。JAXAALOSは使い方にひと手間いるようですが、説明があるのでそのとおりやればできた。

// dem
var dsm = ee.ImageCollection("JAXA/ALOS/AW3D30/V3_2").select('DSM').filterBounds(aoi);
var elevation = dsm.mosaic() // 標高

var proj = dsm.first().select(0).projection();
var slope = ee.Terrain.slope(dsm.mosaic()
                             .setDefaultProjection(proj))
                             .rename('slope'); // 傾斜

var aspect = ee.Terrain.aspect(dsm.mosaic()
                             .setDefaultProjection(proj))
                             .rename('aspect'); // 斜面方位

var mtpi = ee.Image('CSP/ERGo/1_0/Global/ALOS_mTPI').select('AVE').rename('mtpi');

var dem = elevation
             .addBands(slope)
             .addBands(aspect)
             .addBands(mtpi);
Sentinel1はよくわからないのでそのまま

Sentinel1の地形補正ですが、こちらはSRTMのDEMを使って入るらしい、というところまでわかりましたが、あとはちょっと…タイルごとに書き出す部分を削除したのを除いてVenter先生のコードをほぼそのまま使ってみます。出てくるのはmasterstackというマルチバンドのラスタで、ここに「vh、vv、dpol」「メディアン、標準偏差」、「軌道上向き、下向き」の3×2×2=12通りのデータが入っています。

https://code.earthengine.google.com/8b7f12488967850eedbcb2930c439a82

見た感じ大体OKでしたが、白山の周りではちょっとしんどそうな谷間がありました。あと、海上の船?みたいな点が出てくるのも気になります。Sentinel1を弄り回すのはそのうちやってみたい。

sampleRegionsでデータを取得

公式分類ガイドと似たような感じでデータをとって、それをエクスポート。

// Overlay the points on the imagery to get training.
var sampledPoints = imageSentinel.select(bands).sampleRegions({  // imageSentinelというラスタのbandsというバンドから取得
  collection: points, // 点データのシェープファイルの指定
  scale: 10
});

// Google Driveに書き出し
Export.table.toDrive({
  collection: sampledPoints,
  description: 'sentinel2', // 書き出すファイル名
  fileFormat: 'CSV'
});

ということで、ここまでで手順5の半分まできた。こうしてようやく分析に入れるわけです。(続く)