記事作成時から、Tellusからデータを検索・取得するAPIが変更になっております。該当箇所のコードについては、以下のリンクをご参照ください。
https://www.tellusxdp.com/ja/howtouse/access/traveler_api_20220310_
firstpart.html

2022年8月31日以降、Tellus OSでのデータの閲覧方法など使い方が一部変更になっております。新しいTellus OSの基本操作は以下のリンクをご参照ください。
https://www.tellusxdp.com/ja/howtouse/tellus_os/start_tellus_os.html

Tellus上で扱える衛星データの中には、特定のバンド幅を組み合わせることで植物の植生や地表の降雪などの様子をよりわかりやすく可視化できます。

今回は、琵琶湖の湖上に繁茂する水草を衛星画像が検出しようとした論文「LANDSATの衛星画像を利用した滋賀県琵琶湖南湖での水草の検出」を参考に、Tellus上でもLandsat-8の衛星画像を使い、バンド幅と呼ばれる光の波長の組み合わせを変えた画像を表示し、琵琶湖上の水草をわかりやすく可視化できるか試してみようと思います。

(1)琵琶湖の水草とは

琵琶湖では2009年ごろから、オオバナミズキンバイという外来植物が繁茂していることが問題となっています。このオオバナミズキンバイは葉や茎の断片からでも根を伸ばせるほど生命力が強く、急速に分布を拡大しています。

繁殖した水草が排水設備に絡まると水害を引き起こすなどの危険性があり、滋賀県ではこの水草を取り除くためにおよそ３億円の費用をかけるなど、行政上の大きな問題となっています。

琵琶湖が抱える環境問題については、下記の動画もよろしければご覧ください。

(2)水草の検出方法と使用するデータについて

今回、湖上の水草を検出するにあたり、参考にした論文と同じLandsatの衛星画像を利用することとしました。Landsat-8は我々が日常で目にする「可視光線」以外にも近赤外線や熱赤外など、様々な波長帯の光を撮影しています。近赤外域の観測結果には、特に植物による反射が反映されることで知られています。そこため、近赤外の波長帯を画像として扱うことで、植物の分布をよりくっきりと可視化できます。

今回使用するLandsat-8の画像を検索するために、まずはTellus OS上で琵琶湖周辺の地図に移動し、左上のメニューから衛星画像の「シーン」を検索してみます。

「シーン」とは衛星画像が地球の周囲を周回している際に撮影した画像の一つひとつのことをさします。ここでは、なるべく雲によって遮られていない画像を選択するために、「被雲率」（撮影した画像内の雲の割合）の上限を1%に設定して検索しています。

(3)実際にTellusのデータを利用して水草を検出してみた

それでは実際にLandsat-8の衛星画像を用いて、琵琶湖での水草の繁茂位置を確認していきましょう。

まずは必要なライブラリ群をインストールしておきます。

import os, json, requests, math
import numpy as np
import cv2
from skimage import io
from io import BytesIO
import matplotlib.pyplot as plt
import json

%matplotlib inline

まず先にシーンIDから「バンド別幅API」を検索するためのパスを取得する関数を定義します。シーンIDとは、先ほどのシーン検索結果にあった「LC811XXXXXXXXXX」のような文字列のことを指します。

TOKEN = "ここにAPIトークンを設定してください"

def get_landsat_tile_path_from_scene_id(id):
    url = "https://gisapi.tellusxdp.com/api/v1/landsat8/get_scene/{}".format(id)
    headers = {
        "content-type": "application/json",
        "Authorization": "Bearer " + TOKEN
    }
    
    r = requests.get(url, headers = headers)
    data = json.loads(r.content.decode('utf-8'))
    return data["tile_path"]

バンド別画像合成APIからLandsat-8の画像を取得する関数を定義しています。
APIの仕様については下記のページも参考にしてください。
https://www.tellusxdp.com/ja/dev/api

それでは、実際に取得した画像をみてみましょう。

def get_landsat_image(scene_id, zoom, xtile, ytile, opacity=1, r=3, g=2, b=1, rdepth=1, gdepth=1, bdepth=1, preset=None):
    host = "https://gisapi.tellusxdp.com"
    path = "/blend/" + get_landsat_tile_path_from_scene_id(scene_id).format(z=zoom, x=xtile, y=ytile)
    queries = "?opacity={}&r={}&g={}&b={}&rdepth={}&gdepth={}&bdepth={}".format(opacity, r, g, b, rdepth, gdepth, bdepth)
    if preset is not None:
        queries += "&preset=" + preset
    url = host + path + queries
    
    headers = {
        "content-type": "application/json",
        "Authorization": "Bearer " + TOKEN
    }
    
    r = requests.get(url, headers=headers)

    return io.imread(BytesIO(r.content))

option = {}
scene_id = "LC81100352016082LGN02"
zoom = 12
x = 3595
y  = 1620
img = get_landsat_image(scene_id, zoom, x, y)

io.imshow(img)

まずはオプションを指定せずデフォルトの可視光域の画像を取得してみました。

琵琶湖ぞいの陸地から、湖上に白いもやのような箇所が確認できます。こちらが琵琶湖沿いに繁茂する水草の分布でしょうか。

次に、APIであらかじめ用意されているバンド幅の組み合わせであるpresetを指定して画像を取得してみましょう。観測範囲の植物をより目立たせると言われているNatural Colorという合成方法（青：Band1、緑：Band2、赤：Band3）をpresetオプションに指定してみます。

# use composition preset 
img_natural = get_landsat_image(scene_id, zoom, x, y, preset="natural")

io.imshow(img_natural)

こちらはpresetにnaturalカラー画像を指定して取得した画像になります。たしかに、陸地にある植物の領域が緑色に浮き上がっていることがわかりますね。

しかし湖上の水草の部分も多少判別できるものの、陸地の箇所ほど緑色になっているわけではなく、そこまで見やすくなったとは言い難いですね。

次に、植生の活性度を広 域に分析しようとする際によく用いられるNDVI（Normalized Difference Vegetation Index:正規化植生指標）をpresetに指定した画像を見てみましょう。

NDVIは、「NDVI＝（IR－Red）/（IR＋Red）」という計算式で植生の状況を把握することを目的として考案された指標です。緑植生の茂り具合をより目立たせられます。

# use composition preset nvdi
img_ndvi = get_landsat_image(scene_id, zoom, x, y, preset="ndvi")

io.imshow(img_ndvi)

こちらの画像も先ほどとnatural color合成と同じく、陸地上の緑はくっきりと認識できるのですが、湖上の水草については、あまり目立たせることには成功していないようです。

次に、false color合成（青：Band2、緑：Band3、赤：Band4）を見てみましょう。

img_false = get_landsat_image(scene_id, zoom, x, y, preset="false")

io.imshow(img_false)

湖上にある水草の範囲が青白く浮き上がって見えて、これまで利用したpresetの中では比較的水草の分布位置が見やすいように思います。

(4)おまけ：preset変更以外の画像処理「RGB値反転」

せっかくなのでpresetを変える以外の画像処理もしてみたいと思います。

ここでは画像のRGB値を反転させた画像を生成してみます。実際のコードと生成した画像は下記になります。元の画像として上記のfalse color合成で取得した画像を使用しています。

# nega
from PIL import Image
import numpy as np

width = 255
height = 255
img2 = Image.new('RGB', (width, height))

# 色を反転する
reverse_color_pixels = 255 - img_false
for y in range(height):
  for x in range(width):
    # 反転した色の画像を作成する
    r,g,b,o = reverse_color_pixels[y][x]

    img2.putpixel((x,y), (r,g,b))

plt.imshow(np.array(img2))

ややおどろおどろしい画像になったようにも思いますが、赤みがかったことで湖上の水草の分布がより浮かび上がったように思います。

(5)まとめ

今回はLandsat-8の衛星画像から、バンドの組み合わせを変えることで、水草が繁茂していそうな場所を見つけることにチャレンジしました。

論文や報告書から、実際にその場所に水草が生えていそうか確認するに本記事では留めましたが、どのような条件下であればどのようなモノが映っているのかを検証するためには、実際に現場に行って正解データを一度は確認することが大切になります。

確認した結果と衛星データの関係性、つまりはどのような状態の時にどのように映るのか、ということを把握することが、衛星データの理解への第一歩。
今回の例の場合には、バンド組合せごとに湖や水草がどのように見えるか、ということにチャレンジしていますが、
・水草の種類
・水草の生えている水深
・湖の水質(濁り具合)
・気象条件
など、さまざまな条件の違いにより、画像での見え方は変わってきます。

どのような条件下でどのような見え方か、地道にモデルを組み立てることで、最終的に衛星画像から現場の状況を類推できるようになるのです。

また、Landsat-8以外にもバンド幅を組み合わせられる衛星としてAVNIR-2やASNARO-1といった衛星があります。これら別の衛星の画像では水草の様子がどのように違って見えるのか、また別の機会に検証してみたいと思います。

光の波長の組み合わせ次第で、通常の可視光では見えないものが可視化できるのは衛星データの面白いところ。ぜひ皆様の気になる場所で実践されてみてください。