まず、気象衛星が何を見ているのかの概略について、教材ページ [気象衛星「ひまわり」] を見てほしい。
「ひまわり」は、日本の気象庁が運用してきた気象衛星のシリーズ名である。 現在、ひまわり 8号と 9号が運用されている。 その運用は民間会社が担当しているが、発注もとは気象庁である。
「ひまわり」についての基礎的な情報は、気象庁 気象衛星センターのサイトにある。
研究者むけのデータ公開は、気象庁から直接ではなく、4つの機関からされている。
そのうち、千葉大学 環境リモートセンシング研究センターでは、 画像データを物理量の観測値としてつかう人を想定し、 緯度経度格子に変換したプロダクトを作成している。 位置合わせの精度が (自動処理のうちでは) よい。 位置合わせの論文の書誌情報は、下記「リリースノート」中にある。
CEReSの「衛星データ」のページの 「フル公開データ / 気象庁ひまわりシリーズ」の部に「ひまわり8/9号」の節に「データの各種仕様に関するリリースノート」というリンクがある。 そのリンクさきはつぎの文書で、文書名は「リリースノート」ではないが、 ここでは「リリースノート」とよぶことにする。
このリリースノートに記述されたバンド名「TIR 01」と 「EXT 01」のデータ (期間は2015年7月から2021年3月まで) を CEReS から提供していただいた。
画像ファイルに収録された「カウント値」は、 赤外のばあいは輝度温度、可視のばあいは反射 [注] を仮にあらわすものである。 カウント値からそれぞれの物理量の変換表が用意されている。
【[注] 可視センサーがしめす物理量は、千葉大のリリースノートでは「分光反射率」と 書かれている。他のドキュメントでは「アルベド」と書かれていることもある。 本来、反射率あるいは「アルベド」とは、反射放射の強さを入射放射の強さでわったものである。 ところが、衛星観測の可視センサーが観測した量をさすときは、 反射放射の強さを、太陽が真上にあるばあいの入射放射の強さ (一定値) でわったものを そう呼んでしまうことが多い。 ([アルベド (albedo) (2014-07-02)] の記事で「衛星観測独特の用法」と呼んだもの)。 わたしはこれを「アルベド」(albedo) とか「反射率」(reflectance) と呼びたくない。 ここではくるしまぎれに、単に「反射」 (reflection) と呼んでおく。 「規格化された反射」(normalized reflection) というべきかもしれない。】
このページで紹介するのは、1時間に1回 (毎正時) の画像を処理するばあいである。 10分ごとの画像を処理した経験はあるが、そのための変更はふくんでいない。 (別のところで紹介するかもしれない。)
また、たまたま正時の画像が千葉大のアーカイブにない (観測がうまくいかなかったと思われる) ばあいには、 10分後の画像で代用したことがあるが、そのための kiridasi のプログラムの変更点 (何年何月何日何時のばあいは「分」の部分が '10' であるファイルを読め、 というような個別の指示を追加した) は省略している。 (切り出し・物理量への変換をすませた画像のファイル名では「分」は省略しているので、 map の処理は変えなくてすむが、時刻のちがいはファイル名からはわからなくなる。)
赤外画像 (千葉大CEResのいう「TIR 01」、気象庁のいう「Band 13」) について、 緯度経度範囲によって部分を切り出し、数値を輝度温度に変換し、 空間分布を地図上に作図するところまでやってみた。
実際には、GrADS と Python を併用して試行錯誤したのだが、 このページには、Python でとおして処理する手順を説明する。
作業は Linux を動かしているパソコンでおこなった。 Windows でもできるはずだが、切り出しと輝度温度への変換は、Linux のほうが適していると思う。 その結果を転送して Windows で作図することはできる。 その場合のプログラム変更は、 入力データ (輝度温度の数値)、補助入力データ (県境)、と 出力データ (画像) を置く ディレクトリ名を指定するところだけでよいはずである。
切り出す緯度経度範囲は、北緯 20 ～ 50 度、東経 120 ～ 150 度とした。 ここでは便宜上「拡大日本領域」 (greater Japan region) とよび、「gjapan」という符丁で示す。
カウント値から物理量 (赤外は輝度温度、可視は反射) への変換テーブルが、 http://www.cr.chiba-u.jp/databases/GEO/H8_9/FD/count2tbb_v102.tgz にふくまれた 「tir.01」の変換テーブルを適用し、輝度温度の格子点データを書きだした。
【ひまわり9号の運用にともない、千葉大の変換テーブルは改訂され、 http://www.cr.chiba-u.jp/databases/GEO/H8_9/FD/count2tbb_v103.tgz に「tir.01.H08」(8号用) と「tir.01.H09」(9号用) が含まれており、 さらに「tir.01」があるが、これは9号用のようである。 それぞれの日時に、ひまわり8号、9号のどちらが運用されているかによって使いわけるべきものである。】
輝度温度のデータの形式は、4バイト浮動小数点数の単純バイナリデータとした。 (Big endian か little endian かは実行した計算機の CPU によってきまる。 わたしのばあいは Intel CPU なので little endian になる。) ファイル名は「tbb_gjapan2015070700.r4l」のようにつけることにした。 数字列は、年(4けた)、月(2けた)、日(2けた)、UTCの時(2けた) である。観測時刻の「分」は無視している。 「r4l」 (最後の字はエルの小文字) は、4バイト浮動小数点数 little endian であることをしめす、わたしの個人的符丁である。
期間が長くなるとファイル数が多くなるので、1か月ごとにディレクトリをわけ、 その名まえは「201507」のように年(4けた)と月(2けた)の数字にした。
輝度温度の有効な値はかならず正なので、 データが欠損であるところを埋める値としては「-1.0」を採用した。 材料の画像データがなかった日時には、全部が「-1.0」である画像データを作成することにした。
千葉大から受け取ったファイルは、bzip2 圧縮がかかっている。 Python では bz2 パッケージの BZ2File 関数でファイルを開き、 NumPy の frombuffer 関数をつかうことによって、圧縮をといた内容を NumPy の配列に読みこむことができる。
多数の時点について、切り出しと輝度温度への変換をつづけておこなう Python のプログラム をつくった。
これは 2015年 7月7日から 2021年 3月31日までの 1時間ごと (毎正時) のデータを続けて処理する形になっているが、 わたしが 2023年8月現在つかっている Linux マシンでは、 1日ぶんの24シーンを処理するのに約30分かかるので、その日数倍の時間がかかる。 期間を変更するには、yearstart, monthstart, daystart, hourstart, yearend, monthend, dayend, hourend に代入する値を変更する。
この作業でえられた輝度温度のデータは、つぎのような ctl ファイルをつくれば GrADS で読むこともできる。 (この例の ctlファイルは、2015年7月から2021年3月までの期間のデータに対応するものである。)
Python で、cartopy の地図投影と contourf によるぬりわけをつかって、 輝度温度の分布を作図するプログラムをつくった。
region という変数の値としてあたえる文字列によって指定できるようにした。
また、プログラムの「# 作図する緯度経度」 につづく部分で、 あらたな region を追加して、 その端の経度・緯度の値lonmin, lonmax, latmin, latmax を指定することもできる。 ただし、プログラムで region の値を参照しているところのふるまいの設定も必要になる。
[注] ここに示した画像の例は、旧版の作図プログラムによるもので、 現行版の結果は細部がちがいます。(作図時のバージョンでは set_over, set_under が有効でした。)
台風 2019年15号 が関東に上陸していたときである。
熊谷で地上の最高気温の記録 41.1℃ が出たときである。 熊谷付近は雲がなく輝度温度が 300 K 以上になっている。 名古屋付近に輝度温度が 300 K 以上のところがもっと広くある。 (GrADS での作図では名古屋付近に 310 K 以上のところが見られたが、この図では見られない。 ちがいの理由はこれから検討したい。)
浜松で地上の最高気温の記録 41.1℃ が出たときである。 浜松付近に輝度温度が 300 K 以上のところがあるが、あまり広くない。 むしろ岐阜付近に 300 K 以上のところがひろがっている。
冬型の一例。本州の太平洋側は晴れているが、日本海側には雲がある。 (もっとわかりやすい例にさしかえるかもしれない。)
雲がない太平洋岸で、陸面のほうが海面より温度が高いことがわかる。
高分解能 可視画像 (千葉大CEReSのいう「EXT 01」、気象庁のいう「Band 03」) について、 緯度経度範囲によって部分を切り出し、数値を「反射」[仮称] に変換し、 空間分布を地図上に作図するところまで、手順を同様に整理した。
空間分解能が高いのはくわしい情報を得るためにはよいことだが、 同じ空間のひろがりに対するデータ量が多くなる。 ここであつかう EXT 画像の画素の辺の長さは赤外画像の 4分の1 なので、 同じひろがりの画像のデータ量は、赤外画像の 16倍になる。 計算機のディスク容量にもかぎりがあるから、一度にあつかえる時刻の数がすくなくなる。
【ひまわり 8号の可視画像を大まかに見るには、この EXT 画像データから出発するよりも、 千葉大 CEReS の gridded のうちでも「4KM」のデータ (画素間隔を 約 4 km までまびいて物理量に変換したもの) をつかうか、 NICT、JAXA、高知大などで画像化ずみのものを見るのがよさそうだ。】
多数の時点について、切り出しと「反射」への変換をつづけておこなう Python のプログラム をつくった。
緯度経度範囲は、上記の sjapan (限定日本領域) 北緯 30 ～ 50 度、東経 125 ～ 145 度 とした。 (赤外のばあいの gjapan よりも狭い範囲にした。)
日時の範囲は、yearstart, monthstart, daystart, hourstart, yearend, monthend, dayend, hourend で指定する。 利用可能なデータは、2015年 7月 7日 2時 (UTC) から 2021年 3月 31日 23時 (UTC) まであるが、 プログラム例の現状では、2020年 7月 7日の 0時 (UTC) からその日の 23時 (UTC) までをあつかうようになっている。 なお、夜は暗いので、ここでは処理対象を UTC の 0時から 7時および 23時 (JST の 8時から16時まで) の画像にかぎった。
物理量への変換は、count2tbb_v102.tgz にふくまれた 「ext.01」によった。 2023年8月の画像処理から、count2tbb_v103.tgz にふくまれた「ext.01.H08」「ext.01.H09」によることにした。 【なお、可視のばあい、センサーの感度の変化 (劣化) が無視できない。 ひまわり 8号 については、それを考慮した変換表もつくられていて、 CEReSのリリースノートにも紹介されている。 しかしここではそれを採用せず、一定値の変換表によることにした。】
変換結果の物理量は「反射」 (仮称) だが、百分率になっていて、 原則として 0 と 100 のあいだの値をとる。 【同じ千葉大 CEReS 作成のデータでも、他の衛星、 たとえば [GOES 9 データのばあい] は 百分率でなく単純な比で、0と 1の間の値である。 このちがいは、物理量への変換の際には統一せず、作図の際に対処することにした。】
この作業でえられた「反射」のデータは、つぎのような ctl ファイルをつくれば GrADS で読むこともできる。 (この例の ctlファイルは、2020年 7月7日の0時 (UTC) から同じ日の 7時 (UTC)021年3月までの データに対応するものである。 コメントとして、2015年7月7日から2021年3月までの期間全体の日時をいれているが、その形で動作させたことはない。)
Python で、cartopy の地図投影と contourf によるぬりわけをつかって、 「反射」の分布を作図するプログラムをつくった。
日時の範囲の指定の方法は、上の切り出しと「反射」への変換のプログラムと同様である。
色の指定は、増田の習慣で、数値が大きいほど濃い赤、小さいほど濃い青としている。 可視の反射の強さのばあいは、厚い雲で赤く、海面が見えているところで青くなる。 同じ雲が可視と赤外でちがう色になることが多い。
気象衛星画像の習慣と同様に、反射が強いところを白く、弱いところを黒くする グレースケールの colorf も用意し、コメントにしてある。
[注] ここに示した画像の例は、旧版の作図プログラムによるもので、現行版の結果は細部がちがいます。 (作図時のバージョンでは、set_over, set_under が有効でした。)
「反射」の分布の例として、輝度温度の例と同じ日時のものをあげておく。
