まず、気象衛星が何を見ているのかの概略について、教材ページ [気象衛星「ひまわり」] を見てほしい。
「ひまわり」は、日本の気象庁が運用してきた気象衛星のシリーズ名である。 現在、ひまわり 8号と 9号が運用されている。 その運用は民間会社が担当しているが、発注もとは気象庁である。
「ひまわり」についての基礎的な情報は、気象庁 気象衛星センターのサイトにある。
研究者むけのデータ公開は、気象庁から直接ではなく、4つの機関からされている。
そのうち、千葉大学 環境リモートセンシング研究センターでは、 画像データを物理量の観測値としてつかう人を想定し、 緯度経度格子に変換したプロダクトを作成している。 位置合わせの精度が (自動処理のうちでは) よい。 位置合わせの論文の書誌情報は、下記「リリースノート」中にある。
CEReSの「衛星データ」のページの 「フル公開データ / 気象庁ひまわりシリーズ」の部の のところに「ひまわり8/9号」の節に「データの各種仕様に関するリリースノート」というリンクがある。 そのリンクさきはつぎの文書で、文書名は「リリースノート」ではないが、 ここでは「リリースノート」とよぶことにする。
このリリースノートに記述された
のデータ (期間は2015年7月から2021年3月まで) を CEReS から提供していただいた。
リリースノートの情報と見くらべながら、 GrADS の ctl ファイルを作成して、読むことができた。
データは、2バイト符号なし整数、big endian で収録されている。
データが 7月7日の 00 UTC からはじまっているので、 TDEF には、7月7日から31日まで1日24回に対応する 600 を指定した。 31日の月の全部の日時のデータがあるならば 744、 30日ならば 720、29日ならば 696, 28日ならば 672になる。
他の月に適用するばあい、ctl ファイルの変更は、TDEFの行だけでよい。 TDEF に続く数値が、上にのべた日時の個数であり、 そのあと LINEAR のあとに、最初にあらわれる日時を、GrADSの形式で指定する。
実際には画像のない時刻があるが、 対象期間の途中に画像ファイルがなくても、データ欠損にはなるが、処理は継続できる。
なお、DSET の行のファイルのパス名指定は、 CEReSからのディスクの内容のコピーを tir_dist の下に置き、 その横ならび (同じ親ディレクトリの下) のディレクトリで作業することを想定した。
EXT データは、つぎのような ctl ファイルをつくれば GrADS で読めるようだ。 しかし、ファイルが大きいので、動作が不確実である。 地図上にカウント値の分布を表示するところまで行く場合と、 欠損でない数値があるのに「値が全部欠損である」ようなメッセージが出て 図が出ないばあいがあった。
EXT 01 のデータは 2015年 7月 7 日 02 UTC からはじまっているので、 ctl ファイルのはじめもその時刻にした。
画像ファイルに収録された「カウント値」は、 赤外のばあいは輝度温度、可視のばあいは反射 [注] を仮にあらわすものである。 カウント値からそれぞれの物理量の変換表が用意されている。
【[注] 可視センサーがしめす物理量は、千葉大のリリースノートでは「分光反射率」と 書かれている。他のドキュメントでは「アルベド」と書かれていることもある。 本来、反射率あるいは「アルベド」とは、反射放射の強さを入射放射の強さでわったものである。 ところが、衛星観測の可視センサーが観測した量をさすときは、 反射放射の強さを、太陽が真上にあるばあいの入射放射の強さ (一定値) でわったものを そう呼んでしまうことが多い。 ([アルベド (albedo) (2014-07-02)] の記事で「衛星観測独特の用法」と呼んだもの)。 わたしはこれを「アルベド」(albedo) とか「反射率」(reflectance) と呼びたくない。 ここではくるしまぎれに、単に「反射」 (reflection) と呼んでおく。 「規格化された反射」(normalized reflection) というべきかもしれない。】
赤外画像 (千葉大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 によるぬりわけをつかって、 輝度温度の分布を作図するプログラムをつくった。
対象期間は yearstart, monthstart, daystart, hourstart, yearend, monthend, dayend, hourend で指定する。 (例文では、2015年 7月 7日の 0時 (UTC) から、その日の 23時 (UTC) までになっている。)
対象領域として、(gjapan のままでもよいのだが、わざと変えて) つぎの2とおりのものをつくってみた。
【[注] いま作業につかっている Linux マシンにはいっているバージョンの matplotlib では、 contourf で set_under, set_over が無効になってしまう。 その対策として、levels[] と colorf[] の要素をふやしている。 set_under, set_over が有効なばあいの手続きはコメント行にしてある。】
実行するマシンに応じて変更が必要なのは、 datapath, imagefilepath (gkanto のばあいはさらに shapefilepath) だけでよいはずである。
(画素数が多いばあいは、contourf では能率がわるくなり、 pcolormesh をつかったほうがよさそうだが、わたしはまだその方法を修得していない。)
ひとまず、およそ中央になる 250-260 K を薄い黄色、 そこから輝度温度が高いほうを赤系統、低いほうを青系統の段彩にしてみた。
赤は地表面 (陸面・海面) または背の低い雲、 青は背の高い雲に、だいたい対応する。青が濃いほど雲頂が高い。
[注 (2023-08-02)] ここに示した画像の例は、旧版の作図プログラムによるもので、現行版の結果は細部がちがうかもしれません。
台風 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」によった。 (今後は、count2tbb_v103.tgz にふくまれた「ext.01.H08」「ext.01.H09」によるべきである。) 【可視のばあい、センサーの感度の変化 (劣化) が無視できない。 ひまわり 8号 については、それを考慮した変換表もつくられていて、 CEReSのリリースノートにも紹介されている。 しかしここではそれを採用せず、一定値の変換表によることにした。】
変換結果の物理量は「反射」 (仮称) だが、百分率になっていて、 原則として 0 と 100 のあいだの値をとる。 (同じ千葉大 CEReS 作成のデータでも、他の衛星、 たとえば [GOES 9 データのばあい] は 百分率でなく単純な比で、0と 1の間の値である。 このちがいには、物理量への変換の際には統一せず、作図の際に対処することにした。)
Python で、cartopy の地図投影と contourf によるぬりわけをつかって、 「反射」の分布を作図するプログラムをつくった。
日時の範囲の指定の方法は、上の切り出しと「反射」への変換のプログラムと同様である。
対象領域はつぎの2とおりのものをつくってみた。
【[注] いま作業につかっている Linux マシンにはいっているバージョンの matplotlib では、 contourf で set_under, set_over が無効になってしまう。 その対策として、levels[] と colorf[] の要素をふやしている。 set_under, set_over が有効なばあいの手続きはコメント行にしてある。】
[注 (2023-08-02)] ここに示した画像の例は、旧版の作図プログラムによるもので、現行版の結果は細部がちがうかもしれません。
「反射」の分布の例として、輝度温度の例と同じ日時のものをあげておく。
