IT pass HikiWiki - [Exp2023]NetCDF データの可視化 基礎 Diff

((<"スケジュール表・各回資料 (08/09)"|[Exp2023]スケジュール表・各回資料#08-2F09>))

{{toc}}

((<データと可視化ソフトウェアの準備|[Exp2023]データと可視化ソフトウェアの準備>))でダウンロードした NetCDF データを可視化してみましょう.

#============================================================================
= GPhys 概説

NetCDF 形式のファイルにはメタデータを含めることができるため, その内容を上手く解釈すれば, そのファイルだけで「簡単に」図を描くことができます.
そのようなことができるソフトウェアは様々な業界で複数開発されていますが, 本実習では((<地球流体電脳倶楽部|URL:http://www.gfd-dennou.org/>))の有志が開発してきた ((<GPhys|URL:http://www.gfd-dennou.org/library/ruby/products/gphys/>)) を用います.
GPhys は, 地球流体電脳倶楽部有志が 1990 年代から開発してきた Fortran ライブラリである((<地球流体電脳ライブラリ (DCL)|URL:http://www.gfd-dennou.org/library/dcl/>)) をオブジェクト指向スクリプト言語 ((<Ruby|URL:https://www.ruby-lang.org/ja/>)) から呼び出し (((<Ruby-DCL|URL:http://ruby.gfd-dennou.org/products/ruby-dcl/>))), 手軽に可視化できるように工夫して開発されたソフトウェアです.
日本各地の大気海洋研究者に広く使われているソフトウェアの一つです.


#============================================================================
= 最低限可視化

GPhys/DCL を使うと凝った図を描くこともできますが, 必要最低限だけを含めた簡単なスクリプトから説明していきます.


下のスクリプトは((<データと可視化ソフトウェアの準備|[Exp2023]データと可視化ソフトウェアの準備>))でダウンロードした NetCDF ファイル, air.2019.nc, air.2023.nc, を最低限の手順で可視化するスクリプトです.
このスクリプトの内容は下のようになります.

((<tone_xy.rb|URL:/exp/fy2023/230809/practice_gphys/tone_xy.rb>))

    01  #!/usr/bin/ruby
    02  # 使用するライブラリの読み込み. (以下 2 行は「決まり文句」.)
    03  require "numru/ggraph"
    04  include NumRu
    05
    06  # NetCDF ファイル "air.2019.nc" "air.2023.nc" から変数 "air" を読み, GPhys オブジェクト gp に
    07  # 格納
    08  gp = GPhys::IO.open( "air.2019.nc", "air.2023.nc", "air" )
    09
    10  # 画面を開く (open)
    11  #   引数の 1 は画面への描画を表す
    12  #          2 はファイルへの出力を表す
    13  #            (デフォルトでは出力は pdf 形式でファイル名は dcl.pdf)
    14  DCL.gropn(1)
    15
    16  # 描画画面を準備
    17  #   itr の 1 は横軸, 縦軸ともに線形を表す
    18  #          2 は横軸が線形軸, 縦軸が対数軸を表す
    19  #          3 は横軸が対数軸, 縦軸が線形軸を表す
    20  #          4 は横軸, 縦軸ともに対数軸を表す
    21  GGraph.set_fig( 'itr'=> 1 )
    22
    23  # トーン (色付け) で描画
    24  #   第一引数は描画するデータの GPhys オブジェクト
    25  #   第二引数の true はこの時点では「決まり文句」
    26  GGraph.tone( gp, true )
    27  # カラーバーを描画
    28  GGraph.color_bar
    29
    30  # 画面を閉じる (close)
    31  DCL.grcls

最低限の手順としては下のようになります.
* ruby で実行する宣言 (決まり文句; shebang, シバン) (1 行目)
* 使用するライブラリの読み込み (決まり文句) (3-4 行目)
* 変数の読み込み (8 行目)
* 画面を開く (10 行目)
* 描画画面の準備 (16 行目)
* 描画 (26, 28 行目)
* 画面を閉じる (31 行目)

なお, スクリプト中の "#" はコメントを記述するための記号です.
"#" 以降はコメントとして扱われて実行されません.
スクリプトを作成する際には処理内容を理解しやすいようにコメントを書くと良いでしょう.


このスクリプトを使って可視化してみましょう.
air.2019.nc air.2023.nc と同じディレクトリにスクリプトを準備して下のように実行しましょう.

  $ ruby tone_xy.rb

下のような図が表示されるはずです.

((<URL:/exp/fy2023/230809/practice_gphys/tone_xy.png>))

air.2019.nc air.2023.nc に保存されている変数 air は大気温度のデータで, この図は 2019 2023 年 1 月 1 日の 1000 hPa (mbar) 面における大気温度の経度-緯度分布を表しています.

= 可視化される次元の選択

((<NetCDF データの内容確認|[Exp2022]NetCDF データの内容確認>))で確認したように, air.2019.nc air.2023.nc に保存されている変数 air は, 経度(lon), 緯度(lat), 鉛直層(圧力, level), 時間(time) の 4 次元変数です.
つまり, 上のスクリプトでは, 4 次元変数の air から何の指定もないにも関わらず経度-緯度分布が描かれています.
特に指定がない場合には, 最初の 2 次元 (経度と緯度) の分布が可視化されます.
このとき, 3 次元目以降 (鉛直層(level) と時間(time)) の軸は最初の値が採用されます.

== (一つの) 軸の指定

では, 4 次元のうちの一つの次元の値を陽に指定してみましょう.

下のスクリプトは先のスクリプトに少しだけ手を入れて, 経度 135 度のデータを切り出して可視化するスクリプトです.

((<tone_yp_1.rb|URL:/exp/fy2023/230809/practice_gphys/tone_yp_1.rb>))

    01  #!/usr/bin/ruby
    02  # 使用するライブラリの読み込み. (以下 2 行は「決まり文句」.)
    03  require "numru/ggraph"
    04  include NumRu
    05
    06  # NetCDF ファイル "air.2019.nc" "air.2023.nc" から変数 "air" を読み, GPhys オブジェクト gp に
    07  # 格納
    08  gp = GPhys::IO.open( "air.2019.nc", "air.2023.nc", "air" )
    09
    10  # GPhys オブジェクト gp 内の変数 "air" から lon=135 のデータを切り出す (cut)
    11  gp = gp.cut('lon'=>135)
    12
    13  # 画面を開く (open)
    14  #   引数の 1 は画面への描画を表す
    15  #          2 はファイルへの出力を表す
    16  #            (デフォルトでは出力は pdf 形式でファイル名は dcl.pdf)
    17  DCL.gropn(1)
    18
    19  # 描画画面を準備
    20  #   itr の 1 は横軸, 縦軸ともに線形を表す
    21  #          2 は横軸が線形軸, 縦軸が対数軸を表す
    22  #          3 は横軸が対数軸, 縦軸が線形軸を表す
    23  #          4 は横軸, 縦軸ともに対数軸を表す
    24  GGraph.set_fig( 'itr'=> 1 )
    25
    26  # トーン (色付け) で描画
    27  #   第一引数は描画するデータの GPhys オブジェクト
    28  #   第二引数の true はこの時点では「決まり文句」
    29  GGraph.tone( gp, true )
    30  # カラーバーを描画
    31  GGraph.color_bar
    32
    33  # 画面を閉じる (close)
    34  DCL.grcls

先のスクリプトと比べると, 11 行目に切り出しの処理が入っています.
GPhys では, GPhys オブジェクトに作用できるメソッド cut を用いることで特定の経度を切り出すことができます.

このスクリプトを実行してみましょう.

  $ ruby tone_yp_1.rb

下のような図が表示されるはずです.

((<URL:/exp/fy2023/230809/practice_gphys/tone_yp_1.png>))

この図は 2019 2023 年 1 月 1 日の経度 135 度におけるの大気温度の緯度-圧力分布を表しています.

元々, 変数 air は, 経度, 緯度, 鉛直層(圧力), 時間の 4 次元変数でしたが, 上のスクリプトでは経度を指定したため, GPhys オブジェクト gp は 11 行目で緯度, 鉛直層(圧力), 時間の 3 次元となりました.
そのオブジェクトを用いて可視化すると, 緯度(lat)と鉛直層(圧力,level)の分布が可視化されます.
このとき 3 次元目の時間軸は最初の値が採用されます.

=== 練習問題
* 赤道における経度-圧力分布の図を描いてみましょう.


== 時間軸の指定と二つの軸の指定

さらに軸の値を指定してみましょう.

ここでは, 先に指定したように経度 135 度を指定し, さらに日時を指定し, 2019 2023 年 8 月 1 日の分布を可視化してみましょう.

日時の指定にも cut メソッドを使うことができますが, 少し工夫が必要です.
それは, ((<データと可視化の実習準備|[Exp2023]データと可視化の実習準備>))で確認したように, air.2019.nc air.2023.nc に保存されている時間軸の単位が

  time:units = "hours since 1800-01-01 00:00:0.0"

であり, 軸の値が

time = 1919712, 1919736, 1919760, 1919784, 1919808, 1919832, 1919856,
    1919880, 1919904, 1919928, 1919952, 1919976, 1920000, 1920024, 1920048,
    1920072, 1920096, 1920120, 1920144, 1920168, 1920192, 1920216, 1920240,
    ...

という値になっているためです.
ほとんどの人は 2019 2023 年 8 月 1 日が 1800 年 1 月 1 日から何時間後であるのかすぐに計算することはできません.

しかし, GPhys ではこの日時の変換を実現する方法が用意されており,
cut メソッドに

    'time'=>DateTime.parse("2019-08-01 'time'=>DateTime.parse("2023-08-01 00:00:0.0")

のように指定します.

日時の指定に用いているのは DateTime.parse メソッド (正確には DateTime クラスの parse メソッド) です.
このメソッドに引数として "2019-08-01 "2023-08-01 00:00:0.0" を与えたことで, GPhys オブジェクトが持っている日時の単位との関係が自動的に判定され, 日時を正しく指定することができます.

したがって, 先のスクリプトファイルの 12 行目を

    gp = gp.cut('lon'=>135)
    gp = gp.cut('time'=>DateTime.parse("2019-08-01 gp.cut('time'=>DateTime.parse("2023-08-01 00:00:0.0"))

のように cut を続けることで, 経度を切り出し, それに続けて時間を切り出すことができます.

あるいは, cut メソッドは 1 行に複数を続けることもできるため,

    gp = gp.cut('lon'=>135, 'time'=>DateTime.parse("2019-08-01 'time'=>DateTime.parse("2023-08-01 00:00:0.0"))

と書くこともできます.
または, cut メソッドは 1 行に複数を続けることもできるため,

    gp = gp.cut('lon'=>135).cut('time'=>DateTime.parse("2019-08-01 gp.cut('lon'=>135).cut('time'=>DateTime.parse("2023-08-01 00:00:0.0"))

と書くこともできます.


上のようにスクリプトを書き換えて実行してみましょう
(((<書き換えたスクリプトの例 (tone_yp_2.rb)|URL:/exp/fy2023/230809/practice_gphys/tone_yp_2.rb>))).
下のような図が表示されるはずです.

((<URL:/exp/fy2023/230809/practice_gphys/tone_yp_2.png>))

=== 問
* 先のスクリプトで可視化された結果 (2019 (2023 年 1 月 1 日の分布) とはどこがどのように異なるでしょうか?

=== 練習問題
* 春夏秋冬の緯度-圧力分布を描いてみましょう. そしてその違いの原因を考えてみましょう.
* 赤道の 1000 hPa 圧力面における経度-時間分布の図を描いてみましょう.


=== 補足: 日時の値の確認
上に書いたように,

  time:units = "hours since 1800-01-01 00:00:0.0"

の単位を持つ time の値

time = 1919712, 1919736, ...

が何年, 何月, 何日に対応するのかは, 俄かにはわからないでしょう.

これは, ruby スクリプトでは下のようにすることで確認することができます.

  p DateTime.parse("1800-01-01 00:00:0.0") + 1919712/24.0

上のようにすると, 端末 (ターミナル) に下のように表示されるはずです.

  #<DateTime: 2019-01-01T00:00:00+00:00 2023-01-01T00:00:00+00:00 ((2458485j,0s,0n),+0s,2299161j)>

これで, time = 1919712 hours since 1800-01-01 00:00:0.0 が
2019 2023 年 1 月 1 日 0 時 0 分 0 秒に対応することがわかります.

上の命令文は下のような意味になります.

  DateTime.parse("1800-01-01 00:00:0.0")

で time の単位にある基準日時の文字列を解釈してオブジェクトに
しています.

  + 1919712/24.0

で, その基準日時に (time の値の) 日数を足しています.
(time 変数の単位は hours ですが, ここでオブジェクトに足す数値は日数なので 24 で割っています.)

* 参考ページ:
  * ((<GPhys, GGraphチュートリアル (その4)|URL:http://ruby.gfd-dennou.org/products/gphys/tutorial2/ruby_basic.htm>))


== 三つの軸の指定と線グラフ

さらに軸の値を指定してみましょう.

ここでは, 先に指定したように経度 135 度, 2019 2023 年 8 月 1 日を指定し, さらに緯度 0 度を指定して分布を可視化してみましょう.

軸の指定のところを下のようにすると良いでしょう.

    gp = gp.cut('lon'=>135, 'time'=>DateTime.parse("2019-08-01 'time'=>DateTime.parse("2023-08-01 00:00:0.0")).cut('lat'=>0)

なお, cut を作用させる順番は入れ替えることができます.

しかし, この変更だけでスクリプトを実行するとエラーとなります.
それは, 元々 4 次元だった変数・GPhys オブジェクトに対して, 経度, 緯度, 時間の軸を指定したことで変数・GPhys オブジェクトが 1 次元になったためです.
1 次元のデータに対してトーンを使った 2 次元の色分けの図を描くことはできませんので, 可視化するメソッドを変更しなければなりません.

具体的には,

    GGraph.tone( gp, true )
    GGraph.color_bar

を

    GGraph.line( gp, true )

に変えてみましょう.
line は折れ線グラフを描くメソッドです.

上のようにスクリプトを書き換えて実行してみましょう
(((<書き換えたスクリプトの例 (line_p.rb)|URL:/exp/fy2023/230809/practice_gphys/line_p.rb>))).
下のような図が表示されるはずです.

((<URL:/exp/fy2023/230809/practice_gphys/line_p.png>))

これは経度 135 度, 緯度 0 度, 2019 2023 年 8 月 1 日の大気温度の圧力分布を表しています.

しかし, 圧力は鉛直方向の軸ですから作図するときには縦軸として描きたいものです.
そのような場合のために, 縦軸と横軸を入れ替える方法が用意されています.

    GGraph.line( gp, true )

に下のようにして縦軸と横軸を入れ替える指定を追加しましょう.

    GGraph.line( gp, true, "exchange"=>true )

書き換えたスクリプトを実行してみましょう
(((<書き換えたスクリプトの例 (line_p-ex.rb)|URL:/exp/fy2023/230809/practice_gphys/line_p-ex.rb>))).
下のような図になるでしょう.

((<URL:/exp/fy2023/230809/practice_gphys/line_p-ex.png>))

さらに, この図では縦軸である圧力軸は線形ですが対数軸にしたいかもしれません.
(対数軸にすることの意味は何でしょうか?)
その場合には,

    GGraph.set_fig( 'itr'=> 1 )

を下のように書き換えて itr の値を変更しましょう.

    GGraph.set_fig( 'itr'=> 2 )

書き換えたスクリプトを実行してみましょう
(((<書き換えたスクリプトの例 (line_logp-ex.rb)|URL:/exp/fy2023/230809/practice_gphys/line_logp-ex.rb>))).

下のような図になるでしょう.

((<URL:/exp/fy2023/230809/practice_gphys/line_logp-ex.png>))


#============================================================================
= 複数の折れ線の重ね描き

上の 2019 2023 年 8 月 1 日における経度 135 度, 緯度 0 度の温度の圧力分布のスクリプト
(((<書き換えたスクリプトの例 (line_logp-ex.rb)|URL:/exp/fy2023/230809/practice_gphys/line_logp-ex.rb>)))
を基にして複数の折れ線を描いてみましょう.

下のスクリプトを見てみましょう.

((<line_logp-ex_multi.rb|URL:/exp/fy2023/230809/practice_gphys/line_logp-ex_multi.rb>))

    01  #!/usr/bin/ruby
    02  # 使用するライブラリの読み込み. (以下 2 行は「決まり文句」.)
    03  require "numru/ggraph"
    04  include NumRu
    05
    06  # NetCDF ファイル "air.2019.nc" "air.2023.nc" から変数 "air" を読み, GPhys オブジェクト gp に
    07  # 格納
    08  gp = GPhys::IO.open( "air.2019.nc", "air.2023.nc", "air" )
    09
    10  # GPhys オブジェクト gp 内の変数 "air" から
    11  # lon=135, lat=35, および 2019 2023 年 8 月 1 日, 12 月 1 日のデータを切り出 す (cut)
    12  gp1 = gp.cut('lon'=>135, 'lat'=>35, 'time'=>DateTime.parse("2019-08-01 'time'=>DateTime.parse("2023-08-01 00:00:0.0"))
    13  gp2 = gp.cut('lon'=>135, 'lat'=>35, 'time'=>DateTime.parse("2019-12-01 'time'=>DateTime.parse("2023-12-01 00:00:0.0"))
    14
    15  # 画面を開く (open)
    16  #   引数の 1 は画面への描画を表す
    17  #          2 はファイルへの出力を表す
    18  #            (デフォルトでは出力は pdf 形式でファイル名は dcl.pdf)
    19  DCL.gropn(1)
    20
    21  # 描画画面を準備
    22  #   itr の 1 は横軸, 縦軸ともに線形を表す
    23  #          2 は横軸が線形軸, 縦軸が対数軸を表す
    24  #          3 は横軸が対数軸, 縦軸が線形軸を表す
    25  #          4 は横軸, 縦軸ともに対数軸を表す
    26  GGraph.set_fig( 'itr'=> 2 )
    27
    28  # 折れ線グラフを描画
    29  #   第一引数は描画するデータの GPhys オブジェクト
    30  #   第二引数の true はこの時点では「決まり文句」
    31  GGraph.line( gp1, true , "exchange"=>true )
    32  #   第二引数は重ね書きするため false
    33  GGraph.line( gp2, false, "exchange"=>true )
    34
    35  # 画面を閉じる (close)
    36  DCL.grcls

このスクリプトでは, 経度 135 度, 緯度 35 度における 2019 2023 年 8 月 1 日と 12 月 1 日の温度分布を重ね描きします.

スクリプトを変更して実行してみると下のような図になるでしょう.

((<URL:/exp/fy2023/230809/practice_gphys/line_logp-ex_multi.png>))

二つの日時における温度分布を重ね描きすることができました.

2 本の線を描くために 2 本の線のためのデータ (GPhys オブジェクト) が必要なので, 12, 13 行目で元データ (GPhys オブジェクト gp) からそれぞれのオブジェクト (gp1, gp2) を作成しています.
そしてそのオブジェクトは 31, 33 行目に line メソッドに渡すことで二本の線を重ね描きします.
ただし, 最初の折れ線を描く際 (31 行目) には line メソッドの第二引数を true にしていましたが, 二本目の折れ線を描く際 (33 行目) には false にします.
これは, line メソッドの第二引数が true の場合には「図の最初の線」と解釈され, 第二引数が false の場合は「図に追記する線」と解釈されるためです.
もし二本目の折れ線を描く際 (33 行目) の第二引数を true とすると, その線は「図の最初の線」になりますから, 次の画面 (ウィンドウ) に線が描かれます. (一枚目の図が書き終わった後でウィンドウをクリックすると二枚目の図が表示され, 新しい図が表示されます.)

なお, 上のスクリプトに同じように第二引数を false にして line メソッドを使う行を追記すれば, 三本目, 四本目, ... と, さらに多くの折れ線を追加することができます.

しかし, 上の例では二本の折れ線グラフが両方とも黒実線で区別がつきにくいです.
折れ線の色を変えるには

  GGraph.line( gp2, false, "exchange"=>true )

を下のようにすると良いでしょう.

  GGraph.line( gp2, false, "exchange"=>true, 'index'=>20 )

index のオプションの 20 は赤線を表します.
色は index オプションに指定する数字の十, 百の位によって変化します.
(1 の位は線の太さを表します.)

* 参考ページ
  * ((<ラインインデクス|URL:http://www.gfd-dennou.org/library/dcl/dcl-f77doc/Japanese/f77/grph1/node19.html>))

=== 練習問題

上の複数の折れ線を重ね描きする例において, index オプションに与える数字を変えることによってどのように色が変わるか試してみましょう.


#============================================================================
= 可視化するデータ範囲の指定

ここまでの実習では, cut メソッドを使って次元(軸)を指定してきた.
しかし, cut メソッドは, 次元(軸)の範囲を指定するためにも使うことができる.

最初に用いた((<このスクリプト|URL:/exp/fy2023/230809/practice_gphys/tone_xy.rb>))では全球の温度分布を可視化しましたが,
日本付近のみを切り出すことを考える.

そのためには, スクリプトを下のように変更します.

((<tone_xy_reg.rb|URL:/exp/fy2023/230809/practice_gphys/tone_xy_reg.rb>))

    01  #!/usr/bin/ruby
    02  # 使用するライブラリの読み込み. (以下 2 行は「決まり文句」.)
    03  require "numru/ggraph"
    04  include NumRu
    05
    06  # NetCDF ファイル "air.2019.nc" "air.2023.nc" から変数 "air" を読み, GPhys オブジェクト gp に
    07  # 格納
    08  gp = GPhys::IO.open( "air.2019.nc", "air.2023.nc", "air" )
    09
    10  # cut メソッドを使って経度, 緯度の範囲を切り出し.
    11  gp = gp.cut('lon'=>90..180, 'lat'=>0..90)
    12
    13  # 画面を開く (open)
    14  #   引数の 1 は画面への描画を表す
    15  #          2 はファイルへの出力を表す
    16  #            (デフォルトでは出力は pdf 形式でファイル名は dcl.pdf)
    17  DCL.gropn(1)
    18
    19  # 描画画面を準備
    20  #   itr の 1 は横軸, 縦軸ともに線形を表す
    21  #          2 は横軸が線形軸, 縦軸が対数軸を表す
    22  #          3 は横軸が対数軸, 縦軸が線形軸を表す
    23  #          4 は横軸, 縦軸ともに対数軸を表す
    24  GGraph.set_fig( 'itr'=> 1 )
    25
    26  # トーン (色付け) で描画
    27  #   第一引数は描画するデータの GPhys オブジェクト
    28  #   第二引数の true はこの時点では「決まり文句」
    29  GGraph.tone( gp, true )
    30  # カラーバーを描画
    31  GGraph.color_bar
    32
    33  # 画面を閉じる (close)
    34  DCL.grcls

最初のスクリプトとは, 11 行目のみ異なる.
11 行目では, 経度と緯度の範囲を切り出している.
cut メソッドでは,

  '次元名'=>開始値..終了値

の形式で範囲を指定できます.

上のスクリプトを実行してみましょう.
下のような図になるでしょう.

((<URL:/exp/fy2023/230809/practice_gphys/tone_xy_reg.png>))

しかし, この図では切り出された領域がよくわかりません.
このスクリプトに少し手を入れて
((<このスクリプト (tone_xy_reg_coast.rb)|URL:/exp/fy2023/230809/practice_gphys/tone_xy_reg_coast.rb>))のように変更すると, 海岸線が描かれ, 下のような図になるでしょう.

((<URL:/exp/fy2023/230809/practice_gphys/tone_xy_reg_coast.png>))

このように海岸線を描く設定の詳細は後で説明します.


#============================================================================
= 等値線の描画

上ではトーンを用いた可視化を行ってきましたが二次元分布であれば色分けによる表示だけでなく等値線を用いて表示したいこともあるでしょう.
トーンによる色付けは視覚的に大変わかりやすいですが, 各所の値を確認するにはカラーバーを参照しなければなりません.
一方で, 等値線は値を確認しやすい利点があります.
GPhys では等値線を描くメソッドとして contour が用意されています.

最初に用いた((<このスクリプト|URL:/exp/fy2023/230809/practice_gphys/tone_xy.rb>))を
そのためには, スクリプトを下のように変更します.

((<contour_xy.rb|URL:/exp/fy2023/230809/practice_gphys/contour_xy.rb>))

    01  #!/usr/bin/ruby
    02  # 使用するライブラリの読み込み. (以下 2 行は「決まり文句」.)
    03  require "numru/ggraph"
    04  include NumRu
    05
    06  # NetCDF ファイル "air.2019.nc" "air.2023.nc" から変数 "air" を読み, GPhys オブジェク ト gp に
    07  # 格納
    08  gp = GPhys::IO.open( "air.2019.nc", "air.2023.nc", "air" )
    09
    10  # 画面を開く (open)
    11  #   引数の 1 は画面への描画を表す
    12  #          2 はファイルへの出力を表す
    13  #            (デフォルトでは出力は pdf 形式でファイル名は dcl.pdf)
    14  DCL.gropn(1)
    15
    16  # 描画画面を準備
    17  #   itr の 1 は横軸, 縦軸ともに線形を表す
    18  #          2 は横軸が線形軸, 縦軸が対数軸を表す
    19  #          3 は横軸が対数軸, 縦軸が線形軸を表す
    20  #          4 は横軸, 縦軸ともに対数軸を表す
    21  GGraph.set_fig( 'itr'=> 1 )
    22
    23  # 等値線で描画
    24  #   第一引数は描画するデータの GPhys オブジェクト
    25  #   第二引数の true はこの時点では「決まり文句」
    26  GGraph.contour( gp, true )
    27
    28  # 画面を閉じる (close)
    29  DCL.grcls

変更されている箇所は 26 行目です.
このスクリプトでは, tone の代わりに contour を用いています.
また, 色付けしませんのでカラーバーを描画しません.

上のスクリプトを実行してみましょう.
下のような図になるでしょう.

((<URL:/exp/fy2023/230809/practice_gphys/contour_xy.png>))


ここまでやってくるとトーンと等値線を同時に描きたいと思うかもしれません.
そのためにはトーンに等値線を上描きします.

((<ここ|URL:#複数の折れ線の重ね描き>))では, 一つの図に複数の折れ線を重ね描きしました.
そのためには, 最初の折れ線とその後の折れ線で line メソッドの第二引数を変える (最初は true, 二本目以降は false) のでした.
同じようにすることでトーンに等値線を上描きすることができます.

つまり, 上のスクリプトの 26 行目,

    GGraph.contour( gp, true )

を,

    GGraph.tone( gp, true )
    GGraph.color_bar
    GGraph.contour( gp, false )

と変更すると良いでしょう.

ここで,

    GGraph.contour( gp, false )

の第二引数が true ではなく false であることに注意してください.
これは, 等値線を上描きすることを意味しています.
(既に説明していますが, この部分を true にすると「次の図」と解釈されて, 二枚目の図として等値線図が描かれます.)
なお, トーンと等値線を描く順番は重要です. 等値線を先に描いてトーンを後に描くと, トーンによって等値線が上書きされてしまうため等値線が見えなくなります.

スクリプトを変更して実行してみましょう
(((<変更したスクリプト|URL:/exp/fy2023/230809/practice_gphys/tone-contour_xy.rb>))).
下のような図になるでしょう.

((<URL:/exp/fy2023/230809/practice_gphys/tone-contour_xy.png>))

トーンの上に等値線が描かれて見やすくわかりやすい図になりました.


なお, トーンに等値線を上書きしたいときには, 上の方法を使う他に tone_and_contour メソッドを使うこともできますので好みの方法を
使うと良いでしょう.

  GGraph.tone_and_contour( gp, true )


#============================================================================
= 風速ベクトルの描画

次に, 風速ベクトルを描いてみましょう.
GPhys でベクトル図を描くには vector メソッドを使うことができます.
最初に用いた((<このスクリプト|URL:/exp/fy2023/230809/practice_gphys/tone_xy.rb>))を
下のように変更します.

((<vector_xy.rb|URL:/exp/fy2023/230809/practice_gphys/vector_xy.rb>))

    01  #!/usr/bin/ruby
    02  # 使用するライブラリの読み込み. (以下 2 行は「決まり文句」.)
    03  require "numru/ggraph"
    04  include NumRu
    05
    06  # NetCDF ファイルから東西風と南北風の変数を読み, GPhys オブジェクトに格 納
    07  gpu = GPhys::IO.open( "uwnd.2019.nc", "uwnd.2023.nc", "uwnd" )
    08  gpv = GPhys::IO.open( "vwnd.2019.nc", "vwnd.2023.nc", "vwnd" )
    09
    10  # 画面を開く (open)
    11  #   引数の 1 は画面への描画を表す
    12  #          2 はファイルへの出力を表す
    13  #            (デフォルトでは出力は pdf 形式でファイル名は dcl.pdf)
    14  DCL.gropn(1)
    15
    16  # 描画画面を準備
    17  #   itr の 1 は横軸, 縦軸ともに線形を表す
    18  #          2 は横軸が線形軸, 縦軸が対数軸を表す
    19  #          3 は横軸が対数軸, 縦軸が線形軸を表す
    20  #          4 は横軸, 縦軸ともに対数軸を表す
    21  GGraph.set_fig( 'itr'=> 1 )
    22
    23  # ベクトルを描画
    24  #   第一引数は描画するベクトルの x 成分の GPhys オブジェクト
    25  #   第二引数は描画するベクトルの y 成分の GPhys オブジェクト
    26  #   第三引数の true は新しい図を表す
    27  GGraph.vector( gpu, gpv, true, "xintv"=>2, "yintv"=>2, "factor"=>2, "unit_vect"=>true )
    28
    29  # 画面を閉じる (close)
    30  DCL.grcls

注目する個所は二か所です.

一つは変数の読み込み部分で, 二次元ベクトルを描くためには二成分のデータが必要ですから, 7, 8 行目で二つの NetCDF ファイルから二つの変数を読んでそれぞれ別の GPhys オブジェクトに格納しています.

もう一か所はベクトルの描画部分です.
27 行目では vector にベクトルの x 成分と y 成分の GPhys オブジェクトを渡してベクトルを描画しています.

さらに, この例では描画される図を見やすくするためにいくつかのオプションを指定しています.
* xintv, yintv はベクトルを描く際のそれぞれ x, y 方向のデータの間隔です.
  上の例では, 指定なしではベクトルがたくさん書かれるために見難くなってしまいます.
* factor はベクトルの長さを調節する係数です.
* unit_vect は図に描かれているベクトルの長さを示す「標本」を描く指定です.

スクリプトを変更して実行してみると下のような図になるでしょう.

((<URL:/exp/fy2023/230809/practice_gphys/vector_xy.png>))

=== 練習問題
* 温度の分布と風速ベクトルの分布を重ね描きしてみましょう.