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

TelluSARでできること～干渉SARってなに？～

Tellusマーケットの公式ツールとリリースした「TelluSAR（てるーさー）」
このツールは、2時期の衛星画像を使って地盤の沈下や隆起が分かるツールです。

Tellusマーケットをご覧いただくと2種類のTelluSARがあることがお分かりいただけると思います。

アドインと書いてある方が、Tellus OSの中でボタンをクリックしていくことで解析を行うことができるツールで、APIとついている方はコマンドラインベースで解析を行うツールです。それぞれの用途に合わせて、適切なものを選択していただければと思います。

アドインはTellus OS上での操作をしていただければ解析が進められると思いますので、本記事では特にAPIの方のツールの使い方について解説して行きたいと思います。

ここで、簡単にTelluSARで地盤の沈下や隆起が分かる原理を説明して起きたいと思います。

TelluSARで用いているのは「干渉SAR」という技術です。

SAR（さー）というのは、衛星に搭載されるセンサの一つで、衛星から電波を発射して地表面で跳ね返って来る電波を観測することで、地表面の様子を知ることができるセンサです。

電波は文字通り波なので、2時期のSARデータを比べると波の微妙なズレを検知することができます（上図）。これが、地盤の沈下や隆起が分かる仕組みとなります。

詳しく知りたいという方は、宙畑の別の記事で解説していますので、ぜひご覧ください。

TelluSARのAPIを使ってみる

それでは、TelluSARのAPIを利用して、衛星画像から地盤の沈下や隆起を計測してみましょう。

今回利用するTelluSARのAPIリファレンスはこちらです。
TelluSARのAPIリファレンス

事前準備

まずTellusマーケットで「TelluSAR（API）」を購入します（無料で購入できます）。

マーケットトークン発行APIを実行するには、APIアクセストークンが必要です。
APIアクセストークンはダッシュボードの「開発環境」から発行することができます。

ここまで準備ができたら、マーケットで購入したAPIを利用するために、マーケットトークン発行APIを実行しましょう。

import requests, datetime, json

TOKEN = '（ご自身で発行したトークンを入力してください）'
provider_id = 'tellus-product'
tool_label = 'tellusar-api-f'

# マーケットトークン有効期限（30分に設定）
# デフォルトでは5分、最長で60分まで設定可能
expires_at = (datetime.datetime.now(datetime.timezone(datetime.timedelta(hours=+9))) + datetime.timedelta(minutes=30)).isoformat()
def get_market_token(payload={}):
    url = 'https://sdk.tellusxdp.com/api/manager/v1/auth/api_access_token/token'
    headers = {
        'Authorization': 'Bearer ' + TOKEN,
        'Content-type': 'application/json'
    }
    r = requests.post(url, headers=headers, data=json.dumps(payload))
    if r.status_code is not 200:
        raise ValueError('status error({}).'.format(r.status_code))
    return json.loads(r.content) 

# マーケットトークンを発行する
ret = get_market_token({'provider_id': provider_id, 'tool_label': tool_label, 'expires_at':expires_at})
market_token = ret['token']
print(market_token)

これで事前準備は完了です。

ここから、APIの利用方法を紹介していきます。
シーン検索・干渉解析処理・処理結果の取得の3ステップを順番に実行していきましょう。

TelluSARのAPI (1)シーン検索

まずは利用するシーンを決めましょう。
利用できるシーン一覧は、参照APIで確認することができます。

# 利用できるシーン一覧を取得する

def get_free_scene(market_token, payload={}):
    url = 'https://tellusar.tellusxdp.com/api/v1/search'
    headers = {
        'Authorization': 'Bearer ' + market_token
    }
    r = requests.get(url, headers=headers, params=payload)
    if r.status_code is not 200:
        print(r.status_code)
    return json.loads(r.content)

get_free_scene_result = get_free_scene(market_token, {})
print(get_free_scene_result['data']['scenes'])

このコードを実行すると、利用できるシーンの情報が大量に表示されます。
今回は例として、地盤の沈下や隆起がわかりやすい、新燃岳の噴火口のシーンを利用していきます。

シーンは、BeforeとAfterで2枚必要です。

この2枚は撮影した場所だけでなく、撮影した際の衛星の軌道や複数ある観測モードが同一である必要があります。同じ条件で撮影されていないと、わずかな電波のずれを観測することができないためです。条件の合う2枚のペアのことを「可干渉ペア」と言います。
TelluSARにはこの「可干渉ペア」を見つけるためのAPIが用意されています。

まずは、シーンの一覧から1枚目のシーンを選びます。今回は以下のシーンとします。

# 新燃岳（2018年3月6日）のシーンID
scene_id = 'ALOS2204710630-180308'

続いて、取得したシーンIDに対応する「可干渉ペア」のシーンを探します。

# ペアの候補となるシーン一覧を取得する
# after_scene_id

def get_scene_after(scene_id, market_token, payload={}):
    url = 'https://tellusar.tellusxdp.com/api/v1/search/{}/afters'.format(scene_id)
    headers = {
        'Authorization': 'Bearer ' + market_token
    }
    r = requests.get(url, headers=headers, params=payload)
    if r.status_code is not 200:
        print(r.status_code)
    return json.loads(r.content)

get_scene_after_result = get_scene_after(scene_id, market_token, {})
print(get_scene_after_result)

レスポンス

{'data': {'scenes': [{'scene_id': 'ALOS2210920630-180419', 'observation_datetime': '2018-04-19T15:11:28+00:00', 'polarisations': ['HH'], 'left_bottom_lon': 130.537, 'left_bottom_lat': 31.62, 'right_top_lon': 131.26, 'right_top_lat': 32.325}, {'scene_id': 'ALOS2212990630-180503', 'observation_datetime': '2018-05-03T15:11:27+00:00', 'polarisations': ['HH'], 'left_bottom_lon': 130.538, 'left_bottom_lat': 31.62, 'right_top_lon': 131.261, 'right_top_lat': 32.325}, {'scene_id': 'ALOS2235760630-181004', 'observation_datetime': '2018-10-04T15:11:28+00:00', 'polarisations': ['HH', 'HV'], 'left_bottom_lon': 130.536, 'left_bottom_lat': 31.619, 'right_top_lon': 131.262, 'right_top_lat': 32.325}, {'scene_id': 'ALOS2254390630-190207', 'observation_datetime': '2019-02-07T15:11:30+00:00', 'polarisations': ['HH'], 'left_bottom_lon': 130.535, 'left_bottom_lat': 31.62, 'right_top_lon': 131.26, 'right_top_lat': 32.325}, {'scene_id': 'ALOS2268880630-190516', 'observation_datetime': '2019-05-16T15:11:27+00:00', 'polarisations': ['HH'], 'left_bottom_lon': 130.545, 'left_bottom_lat': 31.633, 'right_top_lon': 131.258, 'right_top_lat': 32.336}], 'next': None}}

2枚目のscene_idの値が取得できています。

このとき上手く結果を出すコツとしては、あまり時期の離れていない（まずは3か月とか）シーンを選ぶということがあります。1年以上時期が離れてしまうと、地盤の変化だけでなく地表面の状態が変わってくる可能性が大きくなり上手く干渉しないケースがあるためです。

# 2枚目のscene_id
# 新燃岳（2018年4月19日）のシーンID
# ALOS2210920630-180419

これで、処理に必要な2つのシーン情報が検索できました。

TelluSARのAPI (2)処理

続いて、先ほどまでの手順で取得した2つのシーンを使って、処理をしてみましょう。

# 処理を実行する
# work_idが返却される

def request_work(market_token, payload={}):
    url = 'https://tellusar.tellusxdp.com/api/v1/works'
    headers = {
        'Authorization': 'Bearer ' + market_token
    }
    r = requests.post(url, headers=headers, params=payload)
    print(r.url)
    if r.status_code is not 200:
        print(r.status_code)
    return json.loads(r.content)

work_result = request_work(market_token, {
    'before_scene_id': 'ALOS2204710630-180308',
    'after_scene_id': 'ALOS2210920630-180419',
    'polarisation': 'HH',
    'nlook_rg': 5,
    'nlook_az': 7,
    'filter': 0,
})
print(work_result)

すでに別のユーザーが処理したことのあるペアであれば、結果はすぐに表示されます(‘exist_flag’: True)。

一方でまだ解析が行われていないペアの場合は、処理に時間がかかる（30分～数時間程度）ため、ここではすぐ結果は表示されず、work_idが返却されます。
処理が終わった後で、work_idを利用して処理結果を参照しましょう。

TelluSARのAPI (3)処理結果の取得

最後に、処理結果を参照しましょう。
利用するwork_idには、先ほど処理APIから返却されたwork_idを設定しています。

# work_idを指定して処理結果を取得する

def get_work(work_id, market_token, payload={}):
    url = 'https://tellusar.tellusxdp.com/api/v1/works/{}'.format(work_id)
    headers = {
        'Authorization': 'Bearer ' + market_token
    }
    r = requests.get(url, headers=headers, params=payload)
    print(r.url)
    if r.status_code is not 200:
        print(r.status_code)
    return json.loads(r.content) 

work_id = work_result['data']['work_id']
get_work_result = get_work(work_id, market_token, {})
print(get_work_result)

レスポンス

https://tellusar.tellusxdp.com/api/v1/works/224
{'data': {'work_id': 224, 'before_scene_id': 'ALOS2204710630-180308', 'after_scene_id': 'ALOS2210920630-180419', 'polarisation': 'HH', 'nlook_rg': 5, 'nlook_az': 7, 'filter': 0, 'progress_status': 2, 'request_date': '2020-10-28T12:25:48+00:00', 'complete_date': '2020-10-28T12:25:48+00:00'}}

また、以下のように、処理結果の一覧を参照するAPIを利用することも可能です。

# 処理結果の一覧を取得する

def get_work_list(market_token, payload={}):
    url = 'https://tellusar.tellusxdp.com/api/v1/works'
    headers = {
        'Authorization': 'Bearer ' + market_token
    }
    r = requests.get(url, headers=headers, params=payload)
    print(r.url)
    if r.status_code is not 200:
        print(r.status_code)
    return json.loads(r.content) 

get_work_result_list = get_work_list(market_token, {})
print(get_work_result_list)

最後に、処理結果をpng形式の画像で表示してみましょう。
APIのリクエストパラメータに指定するz,x,yの値を取得する方法は、以下の記事を参考にしてください。

# 処理結果の差分干渉画像（png形式）を表示する

def get_finge_diff(work_id, z, x, y, market_token, payload={}):
    url = 'https://tellusar.tellusxdp.com/api/v1/works/{}/pngs/fringe_diffs/{}/{}/{}.png'.format(work_id, z, x, y)
    headers = {
        'Authorization': 'Bearer ' + market_token
    }
    r = requests.get(url, headers=headers, params=payload)
    print(r.url)
    if r.status_code is not 200:
        print(r.status_code)
    return io.imread(BytesIO(r.content))

# 新燃岳の座標を指定
z = 10
x = 884
y = 416

# 差分干渉画像を取得して表示
res = get_finge_diff(work_id, z, x, y, market_token, {})
io.imshow(res)

計算結果を眺めてみる

無事カラフルな画像を出力することができましたが、このままではどこが隆起しているのか、沈下しているのかよく分かりません。

このカラフルな画像の読み解き方を解説していきます。

ここでは、簡単に画像の読み方について説明しますが、詳しくはTelluSARのユーザーガイドをぜひご覧下さい。

STEP0：画像の観測条件を確認する

先ほど選んだ画像がどのような条件で観測されているのかを確認します。

例えば、以前の記事を使って、観測日時とフレーム番号で該当の画像の詳細情報を出力します。

import requests
from pprint import pprint


def search_file(params={}, next_url=''):
    if len(next_url) > 0:
        url = next_url
    else:
        url = 'https://file.tellusxdp.com/api/v1/origin/search/palsar2-l11'
        
    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 r.json()

ret = search_file({'after': '2018-03-08T00:00:00','before': '2018-03-08T23:59:59','frame_number': '0630'})

pprint(ret['count'])
pprint(ret['items'])

結果として、以下のような情報を得られます。

[{'bbox': [130.536, 31.62, 131.261, 32.325],
  'beam_number': 'U2-6',
  'begin_datetime': '2018-03-08T15:11:25.036000+00:00',
  'center_datetime': '2018-03-08T15:11:30.036000+00:00',
  'code': 'UBS',
  'coordinates': [[130.674, 31.62],
                  [131.261, 31.714],
                  [131.128, 32.325],
                  [130.536, 32.231]],
  'dataset_id': 'ALOS2204710630-180308',
  'date_added': '2019-12-18T09:34:25.186298',
  'end_datetime': '2018-03-08T15:11:35.036000+00:00',
  'files': ['1', '2', '3', '4', '5', '6'],
  'frame_number': 630,
  'mode': 'SM1',
  'observation_datetime': '2018-03-08T15:11:30.036000+00:00',
  'observation_direction': 'R',
  'off_nadir_angle': 29.5,
  'orbit_direction': 'A',
  'path_number': 130,
  'polarisations': ['HH'],
  'publish_link': 'https://file.tellusxdp.com/api/v1/origin/publish/palsar2-l11/ALOS2204710630-180308'}]

STEP1：衛星がどちら側を飛んでいるか

まず、確認するのは衛星がどこを飛んでいるか、です。
SAR衛星は、その観測方法から必ず斜めに撮影をしています。

右と左のどちらから撮影しているのかは、出力した詳細情報の中の「observation_direction」から確認できます。
「R」が右方向で、「L」が左方向です。

また、「orbit_direction」が「D」であれば「降交軌道」、「A」であれば「昇交軌道」となります。

今回の場合は、「昇交軌道」で「右方向」ということになります。
つまり、以下のような配置となります。

衛星は新燃岳の西側を南から北へ通過しています。

色はBefore/Afterの変化の度合いを表します。
2つの波のずれを表現するためーπ～πという値で表現されます。

2つの波がずれなく全く同じであればゼロ、半分ずれるとπということです。もう半分ずれるとまたぴったり重なるのでゼロになってしまいます。

したがって、先ほど調べた衛星から飛んでいる側から何本の縞があるかを数えて、最終的な変位を考えることになります。

今回の場合、広い範囲で色が変化しているので、これは擾乱の影響と思われます。

その中でも、新燃岳の山頂付近はピンク色になっており、これは凡例に照らすと3.14=π の半分程度、すなわち0.5π程度位相が変化しているように見えます。

また、その向きは先ほど確認した通り、西側（画面では左側）から読んでいき、青からピンクへ変化しているため、正の向きであり、衛星から離れていく方向、すなわち沈下していると考えられます。

もっとも、今回のケースでは水蒸気による変化とも考えられ（このあたりが判読の難しいところですね、、、）、地殻変動が生じた場合、これとよく似た位相変化が見られます。

今回のケースでは、2018年3月6日午後2時27分頃、九州の霧島連山新燃岳で爆発的噴火が発生したことが確認されており、その直後の画像となっています。

霧島連山新燃岳の爆発的噴火（2018年3月6日）
https://www.jma-net.go.jp/sat/himawari/obsimg/image_volc.html#obs_j20180306

まとめ

今回は、Tellus初の公式解析ツールであるTelluSARについてご紹介しました。

これまで宙畑でも干渉SARについてはご紹介してきましたが、「自分でやるのは難しそうだな」と思われていた方でも、このツールを用いることで、簡単に解析を行うことができます。
ぜひ、チャレンジしてみてください！
https://www.tellusxdp.com/market/tool_detail/tellus-product/288