気象衛星ひまわり#

衛星データは,放射輝度等をセンサで測定しています。バイアス等の修正を行ったあと,物理量に変換するとともに,座標変換をして格子点に内挿し使いやすい形に変換して提供されています。データ処理の程度をレベルといいます。処理前のセンサの解像度そのままのデータをレベル0,補正や時空間の情報が付加されたものをレベル1,物理量に変換したものをレベル2,一様な時空間に内挿されたものをレベル3,モデルの出力や複数の測定を組み合わせたものをレベル4といいます。ここでは,千葉大学環境リモートセンシング研究センター(CEReS)から提供されているレベル3の気象衛星ひまわり8号のデータを使います。リリースノートを参考に読み解いていきましょう。

ひまわり8号は可視,近赤外,赤外の16バンド(波長帯)を観測しています。CEReSのレベル3赤外データは6,000×6,000点で水平解像度は約2 kmです。ファイルはFTPサーバから提供されています。Python標準ライブラリのftplibを使うこともできますが,1行で済むcurlコマンドを使うことにします。シェルのコマンドを呼び出すには,最初に!を付けます。

server = "ftp://hmwr829gr.cr.chiba-u.ac.jp/gridded/FD/V20190123"
yyyymm = "202107"
ddhhnn = "221200"
band = "tir"
ch = "01"
ext = "fld.geoss.bz2"
fname = f"{yyyymm}{ddhhnn}.{band}.{ch}.{ext}"
url = f"{server}/{yyyymm}/{band.upper()}/{fname}"
!curl -O $url

データにはセンサが測定したカウント値(デジタル値)が格納されています。放射輝度への変換テーブルをダウンロードし展開します。

import numpy as np


cfname = "count2tbb_v102.tgz"
url = f"http://www.cr.chiba-u.jp/databases/GEO/H8_9/FD/{cfname}"
!curl -O $url
!tar zxvf $cfname
_, tbb = np.loadtxt(f"count2tbb_v102/{band}.{ch}", unpack=True)

tarはファイルやディレクトリをまとめたアーカイブファイルを扱うコマンドです。ここではgzip圧縮を解凍するとともに,ディレクトリやファイルに展開しています。 numpy.loadtxt()メソッドを使って変換テーブルを読みます。最初の_は返り値を無視することを表しています。

座標系を定義します。

import numpy as np
 
lon_min, lon_max = 85, 205
lat_min, lat_max = -60, 60
resolution = 0.02
lon = np.arange(lon_min, lon_max, resolution)
lat = np.arange(lat_max, lat_min, -resolution)

データファイルを読んで等価黒体温度に変換します。

import bz2
 
n = 6000
with bz2.BZ2File(fname) as bz2file:
   dataDN = np.frombuffer(bz2file.read(), dtype=">u2").reshape(n, n)
   dataTBB = tbb[dataDN]

ファイルはbzip2圧縮されているので,標準ライブラリbz2を使って開きます。read()メソッドでバイナリデータをメモリに読み,numpy.frombuffer()メソッドで数値に変換しています。バイナリデータは,バイト数とバイトの並び順序(エンディアン)を指定する必要があります。dtype=">u2"で型を指定していますが,>はビッグエンディアン,uは符号なし,2はバイト数を表しています。なお,現在のほとんどの計算機に使われているプロセッサのエンディアンはリトルエンディアンですが,ビッグエンディアンはかつての大型計算機で多く使われていました。ビッグエンディアンは,ネットワークのプロトコルで用いられているため,ネットワークオーダーと呼ばれることもあります。

まず,ヒストグラムを描いてみましょう。

import matplotlib.pyplot as plt
import xarray as xr

xr_tbb = xr.DataArray(np.float32(dataTBB), name="tbb",
                      coords = {
                          'lat':('lat', lat, {'units': 'degrees_north'}),
                          'lon':('lon', lon, {'units': 'degrees_east'})},
                      dims = ['lat', 'lon'])
xr_tbb.plot.hist()
plt.show()
_images/88ce3f159790bac96452f519ec26f6d8ac2956355f8ce21326fc540800a2a305.png

\(T_{\mathrm{BB}}\)の範囲は\(200\)から\(300\)で良さそうです。

xr_tbb.loc[40:0, 110:160].plot.imshow(vmin=200, vmax=300,
                  interpolation="None", cmap="gray_r", figsize=[9,7])
plt.show()
_images/b7150a378345b51beaf8967bdaae466c1855d25aca8a1db4769735a14b6977a9.png