記事作成時から、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

衛星データの1つの特徴である、複数の周波数帯(バンド)を一度に観測している、というもの。「課題に応じて変幻自在？ 衛星データをブレンドして見えるモノ・コト #マンガでわかる衛星データ」をご覧になった読者の中には、実際に自分でブレンドした衛星データを確認してみたい、と思われた方もいらっしゃるのではないでしょうか？
本記事では、Tellus OSと開発環境、それぞれを用いてブレンドする方法についてご紹介します。

1. 利用するデータ

今回は、複数バンドで観測しているLandsat-8のデータを用いた例を紹介します。

2. Tellus OSを用いたブレンドの方法

まず、見たい地域の画像を選択します。

画面右側にある「選択データリスト」をクリックすると、選択した画像情報を確認できます。
初期状態では、Rにバンド4、Gにバンド3、Bにバンド2が割り振られています。

Rに割り振るバンドを変更することで、ブレンドを変えることが可能です。

どのバンドを何にあてたらいいか分からない、という方もいるかと思います。その場合には、Presetをクリックすることで、Tellusで予め用意しているブレンドの組み合わせを確認できます。

例えば陸域/水域をクリックすると、水域と土が見えていそうなところがそれぞれ目立って描画されているのが分かります。他のブレンドがどのように描画されるのか、ぜひご自身で確かめてみてください。

3. 開発環境を用いたブレンドの方法

Tellus OSで見る場合と操作は基本的に同じです。RGBにそれぞれバンドを割り当てることで、好きなブレンドの画像を描画することが可能です。
衛星画像(タイル画像)を取得し、バンドを指定して描画するプログラムとしては下記のような形になります。
タイル画像の引っ張り方の解説については「【ゼロからのTellusの使い方】大きなタイル画像を取得するには？地図タイルを扱う上で知っておきたいTips4選」をご覧ください。

import requests
import json
import math
from skimage import io
from io import BytesIO
%matplotlib inline

TOKEN = '自身のトークンを貼り付けてください'

def search_landsat8_scenes(min_lon, min_lat, max_lon, max_lat):
    url = 'https://gisapi.tellusxdp.com/api/v1/landsat8/scene'
    headers = {
        'Authorization': 'Bearer ' + TOKEN
    }
    payload = {'min_lat': min_lat, 'min_lon': min_lon,
               'max_lat': max_lat, 'max_lon': max_lon}

    r = requests.get(url, params=payload, headers=headers)
    if not r.status_code == requests.codes.ok:
        r.raise_for_status()
    return r.json()


def get_Landsat8_blend_tile(z, x, y, option, params={'opacity':1, 'r':4, 'g':3, 'b':2, 
                            'rdepth':1, 'gdepth':1, 'bdepth':1, 'preset':None}):
    """
    Landsat-8のブレンドタイル画像を取得する
    Parameters
    ----------
    z : int 
        タイルのズーム率 
    x : int 
        タイルのX座標 
    y : int 
        タイルのY座標 
    option : dict 
        APIパスのオプション(path, row, productId)
    params : dict 
        クエリパラメータ
    Returns
    -------
    img: ndarray
        タイル画像
    """
    url = 'https://gisapi.tellusxdp.com/blend/landsat8/{}/{}/{}/{}/{}/{}.png'.format(
        option['path'], option['row'], option['productId'], z, x, y)
    headers = {
        'Authorization': 'Bearer ' + TOKEN
    }
    r = requests.get(url, params=params, headers=headers)
    if not r.status_code == requests.codes.ok:
        r.raise_for_status()
    return io.imread(BytesIO(r.content))
import numpy as np

def get_combined_image(get_image_func, z, topleft_x, topleft_y, 
                       size_x=1, size_y=1, option={}, params={}):
    """
    結合されたタイル画像を取得する
    Parameters
    ----------
    get_image_func : function 
        タイル画像取得メソッド
        引数はz, x, y, option, params
    z : int 
        タイルのズーム率 
    topleft_x : int 
        左上のタイルのX座標 
    topleft_y : int 
        左上のタイルのY座標 
    size_x : int 
         タイルを経度方向につなぎわせる枚数 
    size_y : int 
        タイルを緯度方向につなぎわせる枚数
    option : dict 
        APIパスのオプション(z,x,y除く)
    params : dict 
        クエリパラメータ
    Returns
    -------
    combined_image: ndarray
        結合されたタイル画像
    """
    rows = []
    blank = np.zeros((256, 256, 4), dtype=np.uint8)
    
    for y in range(size_y):
        row = []
        for x in range(size_x):
            try:
                img = get_image_func(z, topleft_x + x, topleft_y + y, option, params)
            except Exception as e:
                img = blank
                
            row.append(img)
        rows.append(np.hstack(row))
    return  np.vstack(rows)

z = 12
(x, y) = (3641, 1614)
option = {'path':107, 'row':35, 'productId':'LC08_L1TP_107035_20151009_20170403_01_T1'}
params={'r':4, 'g':3, 'b':2}
combined = get_combined_image(get_Landsat8_blend_tile, z, x, y, 2, 1, option, params)

io.imshow(combined)

上記の例ではまずは普段目にするTrue Colorの画像を表示しました。

以下の各パラメータを変更することで、取得する画像とブレンドを変更できます。
params={‘r’:7, ‘g’:6, ‘b’:5}
のrgbに割り振るバンド(番号)を変更してみましょう。

例えば、以下のバンドの組み合わせからは、熱源が確認しやすくなる、と言われています。赤く目立っているところが熱源と思しきところです。

z = 12
(x, y) = (3641, 1614)
option = {'path':107, 'row':35, 'productId':'LC08_L1TP_107035_20151009_20170403_01_T1'}
params={'r':7, 'g':6, 'b':5}
combined = get_combined_image(get_Landsat8_blend_tile, z, x, y, 2, 1, option, params)

io.imshow(combined)

他にも、例えば植物を目立たせたい場合には以下のようにバンドを割り当てます。

z = 12
(x, y) = (3641, 1614)
option = {'path':107, 'row':35, 'productId':'LC08_L1TP_107035_20151009_20170403_01_T1'}
params={'r':5, 'g':4, 'b':3}
combined = get_combined_image(get_Landsat8_blend_tile, z, x, y, 2, 1, option, params)

io.imshow(combined)

水や土を目立たせたい場合は以下のブレンドを試してみてください。

z = 12
(x, y) = (3641, 1614)
option = {'path':107, 'row':35, 'productId':'LC08_L1TP_107035_20151009_20170403_01_T1'}
params={'r':5, 'g':6, 'b':4}
combined = get_combined_image(get_Landsat8_blend_tile, z, x, y, 2, 1, option, params)

io.imshow(combined)

4. まとめ

同じ場所の画像であっても、選択するバンドの組み合わせを変えることで、異なる見え方が可能になります。物事の本質を捉えるときには、多角的に見てみる、ということが重要だと言われています。複数の見方ができる衛星データは、多角的に見てみることにとても向いています。

Tellusで実装されているpresetはほんの一部、もっと他にも組み合わせはあります。各波長の特性(どの波長帯は何を強く反射するのか)を理解して試してみる、ということが重要です。ただ、いきなり波長の特性を理解しながらブレンドを探すのは難しいと思いますので、いろいろと試しに変えてみて、何が目立って見えているのか、ということから適したバンドの組み合わせを見つけるのも良いかもしれません。
皆さまもご自身のブレンドを見つけるべく、割り当てるバンドの組み合わせを変えながら衛星データを見てみてはいかがでしょうか。

衛星データプラットフォーム「Tellus」の登録はこちらから。
※TellusでJupyter Notebookを使う方法は、「【ゼロからのTellusの使い方】Tellusの開発環境でAPI引っ張ってみた」をご覧ください。

参考記事