Embedded Files
雲はいらぬ
雲はいらぬ
前回、Sentinel-2のmedianコンポジット画像をつくってみましたが、かなりガスった画になりました。これは雲が写り込んだピクセルも含めて統計値をとったためです。事前に雲を除去して、cloud-freeなピクセルだけで改めてmedianをとれば、より綺麗な画が得られると予想されますのでやってみましょう。
Sentinel-2 level 1Cプロダクトは品質保証（quality assurance）データが付帯しており、その中のdense cloud, cirrus cloudのマスクが利用できます[1]。GEEでは当該プロダクトのImageCollectionの中の、QA60というバンドがそれです。もう一度データサマリーを開き、下の方までスクロールしていってみましょう。
衛星データのQAではよくありますが、ビットマスクです。2進数で、10番目のビットがdense(opaque) cloudフラグ、11番目のフラグがcirrus cloudフラグです。これらがともに0なら、そのピクセルはcloud-freeだ、と判定できます。
この判定処理を、まずひとつのImageについて行う関数を実装しましょう。こういう場合、ビット演算を利用するのが便利です。少し複雑に見えるかもしれませんが、以下のような関数で実現できます。
function構文自体はJavaScriptの仕様です。imageという仮引数（Image型を想定）をとるcloudMaskingという関数を定義しています。
2行目でimageからQA60のバンドを抽出し、qaというImageを作っています。
3, 4行目はビットシフト演算で、それぞれ10ビット目、11ビット目を1にした整数を作っています。
5行目で雲マスクを作っています。Image型のメソッドが連鎖しており、.bitwiseAndはビットAND演算、.eq(x)は「ピクセル値がxでないときに0を割り当てる」論理演算、.andは「前後に指定したImageのピクセル値がともに0以外であればTrue(1)、そうでなければFalse(0)を返す」論理演算です。ややこしいですが要するに、２つの雲フラグがともに0であればTrue(1),そうでなければFalse(0)となるImageが変数maskに格納されます。
6行目は返り値で、「マスクが更新されたimage」を返します。Image型にはデータ以外の付帯情報として「マスク」をつけることができます。マスクは1枚の画像で、この画像のピクセル値が0以下であれば、Imageの当該ピクセルは後の処理から除外されます（正であれば処理される）。
というわけでまとめますと、cloudMaskingにImageを与えると、雲ピクセルのみマスクされた（同じ構造をもつ）Imageが帰ってきます。
さて、あとはこの関数を、medianコンポジットをとる前に、ImageCollectionがもつ全てのImageに対して適用すればいいのですが、どうしましょうか。
for文的ななにかで回す？
いいえ、ImageCollectionには.map() という神メソッドがあります:
.map() は引数に与えた関数をImageCollectionのもつすべてのImageに作用させ、結果をImageCollection型で返します。ここでは一括雲マスクをかけたことになります。
というわけで、前回のコードを以下のように変更して実行してみましょう。
だいぶ綺麗になりました。まだ全体が青っぽくボヤッとしてない？と思う方、これをどうにかしたければ大気補正が必要です。Level 1CデータはTOA反射率、つまりレイリー散乱等によるパスラジアンスや大気中での減衰の影響がのったままのデータになりますので注意してください。とはいえこのデータでも伐採の様子は十分見えると思います。
森を見たい
森を見たい
森林とそれ以外を際立たせるために、植生指数NDVIを計算してみましょう。次のような関数を作ります:
.expressionメソッドはImage中の複数のバンドを指定してラスタ演算を行います。書き方は感じ取ってください。.select() メソッドで、Image内の特定のバンドだけを抜き出すことができます。
では描画しましょう。
森林は緑っぽく、水域は青、それ以外は黄色っぽくなるようなカラースケールにしてあります。右上のウィンドウでInspectorタブを選択し、マップをクリックすると、カーソルのある場所のNDVIや緯度経度情報を示してくれますので、実際にいろいろな土地被覆のNDVIがどんな値になっているかチェックしてみましょう。
ちなみにトゥルーカラー画像も重なって表示されており、描画画面のLayersからNDVIのチェックを外せば、先ほどのトゥルーカラー画像が見えるはずです。次回はいよいよ時系列で森林の変化を見てみましょう。
[1] 2019年1月現在では、cloud shadow（雲影）のマスクまでは提供されていないようです。雲影も意外と解析に影響するのでやっかいです。どうしても気になる方は、独自に処理アルゴリズムを作りmap演算を追加しましょう。
Google Sites
Report abuse