GrADS-Note

'set gxout shaded'

という中身の ”sh.gs” を作成し、同じ場所でgradsを開いて

 open 〜
 sh
 d 変数名

と打ちこんでみると、シェードで変数が描かれる。
これはgrads上で ”sh” と打つと ”run sh.gs” が実行されるためである。※
→東北大学大学院理学研究科 流体地球物理学講座 公開情報/GrADS

これを使って、汎用性の高いコマンドをGrADSスクリプトに書き留めておき、指定の場所に置いておくことで、作業効率が格段に向上する。

※このため、コマンド名.gsというスクリプトを作ってしまうのはタブーである。例えばrun.gsなど。run hoge.gsと命令したつもりがrun.gs hoge.gsと命令されてしまう。

Web上で公開されている、種々のサイトのGrADSスクリプト（＋オリジナルも複数追加）をzipにまとめたので(2015年02月14日更新)、以下の手順で使用することができる。

https://researchmap.jp/multidatabases/multidatabas...

gradslib.zipをダウンロード、解凍。

unzip gradslib.zip

.bashrcに

export GASCRP=/・・・/gradslib　　#gradslibへの絶対パス

の一文を加える。

もしくは、既に設定されているGASCRPのディレクトリの中に入れる。

ついでに、basemap.gsで必要になるマスクデータも、gradslib/以下に置く。

wget ftp://cola.gmu.edu/grads/scripts/lpoly_lowres.asc
wget ftp://cola.gmu.edu/grads/scripts/lpoly_mres.asc
wget ftp://cola.gmu.edu/grads/scripts/lpoly_hires.asc
wget ftp://cola.gmu.edu/grads/scripts/opoly_lowres.asc
wget ftp://cola.gmu.edu/grads/scripts/opoly_mres.asc
wget ftp://cola.gmu.edu/grads/scripts/opoly_hires.asc
wget ftp://cola.gmu.edu/grads/scripts/lpoly_US.asc

x軸が11×8.5の座標の中でx=1.0〜10.0、y=2.0の位置にあるとき、
目盛りで6つに区切られた１区間内を5つに分ける補助目盛を表示
（色＝1、太さ＝4、補助区分＝5、目盛区分＝6、目盛りのサイズ＝0.05）

auxlinex 1.0 10.0 2.0 -c 1 -th 4 -d 5 -n 6 -s 0.05

目盛区分の半分のところから目盛りを書きたい場合
（軸ラベルを60度間隔で書いているが、開始が30度である場合など）

auxlinex 1.0 10.0 2.0 -c 1 -th 4 -d 1 -n 6 -s 0.07 -st 2

使用例

 'reinit'
 'set display color white'
 'c'
 'sdfopen precip.mon.ltm.nc'　　#CMAP降水量気候値
 'set mproj scaled'　　　　　　#x軸、y軸の位置がずれないように設定
 'set lon 0 360'
 'set lat -90 90'
 'set xlint 60'
 'set ylint 30'
 'set parea 1 10 2 7'
 'set xlopts 1 6 0.22'
 'set ylopts 1 6 0.22'
 'set gxout shaded'
 'color 0 10 0.2 -kind white->blue->lime'
 'set grads off'
 'd ave(precip,t=1,t=12)'
 'auxlinex 1.0 10.0 2.0 -c 1 -th 4 -d 5 -n 6 -s 0.05'　#南x軸に補助目盛りを表示
 'auxliney 1.0 2.0 7.0 -c 1 -th 4 -d 5 -n 6 -s 0.05'　#西y軸に補助目盛りを表示
 'auxlinex 1.0 10.0 7.0 -th 6 -d 1 -s -0.07'　#北x軸に目盛りを表示
 'auxlinex 1.0 10.0 7.0 -th 4 -d 5 -s -0.05'　#北x軸に補助目盛りを表示
 'auxliney 10.0 2.0 7.0 -th 6 -d 1 -s -0.08'　#東y軸に目盛りを表示
 'auxliney 10.0 2.0 7.0 -th 4 -d 5 -s -0.05'　#東y軸に補助目盛りを表示
 'xcbar 4.5 6.5 1.0 1.4 -fwidth 0.18 -fheight 0.19 -edge box -line off -fstep 25'

GrADS Script Library(日本語訳)

陸上/海上を塗りつぶす

 d var
 basemap L/O

陸上を色番号1(デフォルトは15)で塗りつぶし、境界の色を色番号2(デフォルトは0)にする

basemap L 色番号1 色番号2

basemap.gsでは以下のマスクファイルが必要。

wget ftp://cola.gmu.edu/grads/scripts/lpoly_lowres.asc
wget ftp://cola.gmu.edu/grads/scripts/lpoly_mres.asc
wget ftp://cola.gmu.edu/grads/scripts/lpoly_hires.asc
wget ftp://cola.gmu.edu/grads/scripts/opoly_lowres.asc
wget ftp://cola.gmu.edu/grads/scripts/opoly_mres.asc
wget ftp://cola.gmu.edu/grads/scripts/opoly_hires.asc
wget ftp://cola.gmu.edu/grads/scripts/lpoly_US.asc

set mpdsetで地図の解像度を変更しているときは、低解像度L、中解像度M、高解像度Hを表記する

basemap L/O 色番号1 色番号2 L/M/H

※オリジナルのスクリプトでは作業場所に海陸ascファイルが必要であったが、上記の方法でインストールしたbasemap.gsは73行目で海陸ascファイルのディレクトリを指定している。使用時はここを適宜書き直すこと。

※hiresのとき、北米大陸上にしか色が塗られないバグあり。未修正。

GrADS Script Library(日本語訳)

x, y = 3.0, 6.0　の位置に、色1、マーク0、スタイル1の凡例線と「line1」という文字を表示

cbar_line -x 3.0 -y 6.0 -c 1 -m 0 -l 1 -t "line1"

並べて表示

cbar_line -x 3.0 -y 6.0 -c 1 2 3 -m 0 -l 1 -t "line1" "line2" "line3"

x, y = 3.0, 6.0　の位置に、太さ6、スタイル1、色1の凡例線に、マーク1、サイズ0.13のマークをつけ、「line1」という太さ4、色1、サイズ0.17の文字を表示

cbar_line_set -x 3.0 -y 6.0 -n 6 -l 1 -c 1 -m 1 -z 0.13 -sn 4 -sc 1 -sz 0.17 -t "line1"

ヘルプ表示

cbar_line_set

マークに線を重ねたくないとき、-rオプションで線の両端を縮める

cbar_line_set -x 3.0 -y 6.0 -r 0.08 -n 6 -l 1 -c 1 -m 1 -z 0.13 -sn 4 -sc 1 -sz 0.17 -t "line1"

markplot.gsとtrackplot.gsを使うと、線・マークのサイズ・太さ・位置を細かく調整できて便利。

GrADS Reference Manual

シェードを描いたあとに、

cbarn

と打つと、図の下にカラーバーが表示される。
ただしデフォルトの設定だと縦横比の関係で数字が一部切れてしまうことがある。これを防ぐには

cbarn 1 0 5.5 0.18

等の微調整が必要になる。

cbarn scale num xmid ymid

　　　 scale：カラーバーの大きさ。defaultは1.0
　　　 num：カラーバーを縦にする場合は1、横にする場合は0とする。
　　　 xmid,ymid：カラーバーの位置の指定をする。数字はカラーバーの中心がvirtual pageのどこに位置するのかを指定している。

カラーバーを画面に収める調節の例

 set parea 0.8 10.2 1.0 8.0
 set lon 0 360
 set lat -90 90
 set xlint 30
 set ylint 30
 set xlopts 1 6 0.2
 set ylopts 1 6 0.2
 set gxout shaded
 color 〜
 set grads off
 d var
 cbarn 1.2 0 5.5 0.3

その他、長方形のcbar、隅に表示されるcbarcなどがある。

カラーバーの設定(xcbar.gs)

cc

中身

 'set display color white'
 'c'
 'xy'

Grads Script 集

nps/sps表示において東経90度線上に30度ごとに緯度を0.18の大きさの文字で表示

 set mproj nps/sps
 d var
 circlat 30 90 0.18

東経90度線上、30度ごと、0.18の大きさ、文字スタイル1

circlat 30 90 0.18 1

※-360〜720度の範囲に対応できるよう、オリジナルを修正

Grads Script 集

nps/sps表示において30度ごとに経度を表示

 set mproj nps/sps
 d var
 circlon 30

30度ごと、文字の大きさ0.18、文字のスタイル1

circlon 30 0.18 1

※数字の下を白く塗らず、90度の範囲に対応できるよう、オリジナルを修正

co

中身

 'set gxout contour'
 'set ccolor 1'

color.gsの説明（小玉さん）

最小値から最大値までを青→赤でシェード表示

 color 最小値 最大値 間隔
 d var

色のグラデーションの設定（降水量　茶色→白→緑）
偏差

 color 最小値 最大値 間隔 -kind saddlebrown->white->green
 d var

生の値

 color 最小値 最大値 間隔 -kind white->green
 d var

0近傍に色を塗りたくない

color ...

のあとに表示される、ccols= 16 17 18 19 20 21（6区分の場合）
をコピーし、中心近傍を0に変更した

set ccols 16 17 0 0 20 21

を入力したのち、dすればよい。

寒色系から暖色系

darkblue->blue->dodgerblue->skyblue->white->white->gold->darkorange->red->darkred

乾燥から湿潤

maroon->saddlebrown->darkgoldenrod->khaki->white->white->palegreen->lightgreen->limegreen->darkgreen

一定値以上とそれ以下を強調（降水量や水蒸気量の気候値など）

white->antiquewhite->mistyrose->lightpink->mediumvioletred->navy->darkblue->blue->dodgerblue->aqua

GrADS-aholic / Script / colormaps.gs

最小値265、最大値310、間隔1、カラーパレット"ocean"を使用

colormaps -l 265 310 1 -map ocean

色の並びを逆に

colormaps -flipped -l 265 310 1 -map ocean

使用できるカラーパレット一覧

color.gsの初期設定を saddlebrown->white->green に変えたもの

Colorbrewer website



6種類のカラーパレットを選べる。色の塗り分け数は、6、7、8、9に対応（デフォルトは7）。5以下の塗り分けは7とかで作ったccolsを抜き出せばよい。colorcbコマンドのあとに、set clevsでレベルを設定のこと（デフォルトは【-3 -2 -1 1 2 3】に設定）塗り分けタイプは、Kaye et al. (2012; GMD)のFig.4にもとづく。

使用方法

 colorcb 塗り分け数 タイプ
 set clevs -3 -2 … 
 d var

cタイプのカラーパレットを9段階（-4,-3,-2,-1,1,2,3,4）に設定

 colorcb 9 c
 set clevs -4 -3 -2 -1 1 2 3 4
 d var

カラーパレットを作成してくれるサイトOnline HCL Creatorを使うと、好みの色調のカラーパレットを
GrADSで使うときのスクリプトを作ってくれて便利。

経度a〜b°、緯度c〜d°に長方形を描く

drawbox a b c d

※11x8.5の画面格子で指定したい場合は、draw recf

Bin Guan script libraly

通常図の下・左にしか表示されない軸とラベルを上・右に表示

 d var
 faxis

スクリプト集（小玉さん）

3から5の値の範囲にハッチをかける

 d var
 hatch var 3 5

ななめ格子状のハッチ

 hatch var 3 5 -angle 45
 hatch var 3 5 -angle 135

粗さを0.1に上げ、ハッチの間隔を0.5に広げる

hatch var 3 5 -density 0.1 -int 0.5

※オリジナルのスクリプトを、ハッチ線のスタイル、色、太さを設定できるように書き換えた(ver.0.01r1_2)

ハッチ線のスタイル(デフォルト：1)、色(1)、太さ(3)を設定する

hatch var -stlye 1 -color 1 -thick 3

GrADS Script Library(日本語訳)

 map type
 d var

typeの種類：nps, sps, usa, usa2, n_amer, s_amer, africa, europe, euro2, asia, aust, robinson, c_pac, n_pac, n_atl, nps2, nps3

緯度140.0度、経度40.0度の地点にマーク3、サイズ0.2、色1のマークをプロット

markplot 140.0 40.0 -m 3 -s 0.2 -c 1

※set line 1 1 X でマークの線の太さを変更できる。

使用例（沖縄、大阪、東京、札幌の地点をプロット）

 'set lon 120 152'
 'set lat 20 48'
 'set mpdset hires'
 'set clevs -1e+30'
 'd var'
 'markplot 127.41 26.13 -c 1 -m 2 -s 0.25'
 'markplot 135.29 34.41 -c 2 -m 3 -s 0.25'
 'markplot 139.45 35.41 -c 3 -m 4 -s 0.25'
 'markplot 141.21 43.03 -c 4 -m 5 -s 0.25'

ここを参考にしました。

スクリプト集（小玉さん）

t=1から4までの4枚の図を並べる

 mul 2 2 1 1
 d var(t=1)
 mul 2 2 1 2
 d var(t=2)
 mul 2 2 2 1
 d var(t=3)
 mul 2 2 2 2
 d var(t=4)

var1に対するvar2の散布図上に回帰直線をひく(色1、線種1、太さ5)

 set gxout scatter
 set t 1 30
 d var1;var2
 regline var1 var2 -c 1 -l 1 -t 5

GrADSスクリプトでよく使う文まとめ→散布図のオプション

sh

中身

'set gxout shaded'

エマグラム（Skew-Tダイアグラム）の描画
Bob Hart/Software/plotskew.gs

Bin Guan's GrADS Script Library

まず、taylor.gsとともにparsestr.gsfとqdims.gsfを用意。

'open hoge.ctl'
'set parea 1 5 1 5'
'set line 1 1 5'
'set strsiz 0.15 0.16'
'taylor -s 0.8 1.3 1.7 -r 0.7 0.95 0.9 -t A B C -z 0.25 -c 1 2 4 -m 2 6 8'

[微修正]

35行目、tmpdirを'./tmp'に。オリジナルでは/tmp/ディレクトリに一時ファイルを格納する仕組みになっていたが、トラブルを避けるため。

set string 1 l -1を、全てset string 1 l 4に変更。文字が細くて見づらい。

58行目、tick_length=0.05を0.09に変更。目盛りを少し長めに。

緯度120.0度、経度20.0度から緯度140.0度、経度40.0度に色1、線種1、太さ6の線をプロット

trackplot 120.0 20.0 140.0 40.0 -c 1 -l 1 -t 6

使用例（沖縄、大阪、東京、札幌の地点間の線をプロット）

 #markplotと同じ工程#
 'd var'
 'markplot 127.41 26.13 -c 1 -m 2 -s 0.25'
 'markplot 135.29 34.41 -c 2 -m 3 -s 0.25'
 'trackplot 127.41 26.13 135.29 34.41 -c 2 -l 1 -t 6'
 'markplot 139.45 35.41 -c 3 -m 4 -s 0.25'
 'trackplot 135.29 34.41 139.45 35.41 -c 2 -l 1 -t 6'
 'markplot 141.21 43.03 -c 4 -m 5 -s 0.25'
 'trackplot 139.45 35.41 141.21 43.03 -c 2 -l 1 -t 6'

ここを参考にしました。

http://www.jamstec.go.jp/frsgc/research/iprc/nona/...

位置（3, 0.5）に、スケール1, 4の太さ5のベクトル凡例を描く

 set line 1 1 5
 set cthick 5
 vec var1;var2 -SCL 1 4 -P 3 0.5

IT memo

アニメーションをクリックで操作　左クリック：進める　右：戻る　中：終了

 set t 最小値 最大値
 xamin var -pause

1秒ごとに表示

 set t 最小値 最大値
 xamin var -sec 1

〜間隔で表示（12ならt=1,t=13,t=15・・・）

 set t 最小値 最大値
 xamin var -sec 間隔

スクリプト集（小玉さん）

x座標2.0〜9.0、y座標0.6〜1.0の範囲に、ラベルの幅0.18、高さ0.19(デフォルトの1.5倍)、端の形三角、ラベル表示間隔2、境界線ありでカラーバーを表示

xcbar 2.0 9.0 0.6 1.0 -fwidth 0.24 -fheight 0.26 -edge triangle -fstep 2 -line on

境界線を白抜き

-line on -lc 0

• line on　のとき線の太さを変える

 set line 1 1 6　　　# 線の太さを6に
 xcbar 2.0 9.0 0.6 1.0 -line on

Bin Guan script libraly

マスクアウトしない年を指定

 set ccolor 1
 d var
 yrmask 1981 1982 1983 1985 1986 1988 1990 1994 1995 1996 1999 2000 2003 2005 2009 2010 2012 a
 # 時系列データaとして、指定した年は1、それ以外は欠損値のデータをdefine
 set ccolor 2
 d var*a
 # 指定した年のみ描く

＃このgsを使うよりも、全ての値が0の線を引く方法のほうが汎用性がある。 ＃

 set cthick 3
 set ccolor 1
 set cmark 0
 d const(lev,0)
 set cthick 5

＃など。＃




通常表示されない折れ線グラフの0線を表示

 set vrange -3 3
 zeroline

※事前にvrangeを設定すること

色番号、マークの種類、線種、太さ、線を引く値を変更したい場合
デフォルトではzeroline 1 0 1 6 0.0

zeroline color mark linetype thickness value

ヘルプ表示

zeroline -help

Bin Guan script libraly

t=1からt=100までのvar1とvar2の変動の相関係数

 set t 1 100
 fcorr var1 var2
 set t 1
 d fcorrout

なぜか365_day_calendarを設定してるとエラーになる。

Bin Guan script libraly

t=1からt=100までの線形回帰の値、傾き、標準偏差

 set t 1 100
 ltrend var a b c
 d a　　　# 線形回帰直線の表示(1次元の場合のみ)
 set t 1
 d b　　　# 傾きの表示
 d c　　　# 標準偏差の表示

set t 1 34
ltrend hoge a b c
criteria=b*57.2*sqrt(34-2)/sqrt(34)/c

ltrendではRMSE（y_i - Y_iの2乗の平均の平方根）を得ることができるので、これを利用する。
57.2の部分には、説明変数の偏差平方和（x_i - x_barの2乗のΣ）の平方根を入れる。criteriaが有意水準に対応するtの値よりも大きければ有意。

varのt=aからt=bまでの期間の各月気候値climを計算

 monclim var clim 1 240
 set t 1 12
 d clim

中身

 function clim(args)
 var=subwrd(args,1)
 varclim=subwrd(args,2)
 tst=subwrd(args,3)
 ten=subwrd(args,4)
 tst2=tst+12
 
 'set t 'tst' 'tst2''
 'define 'varclim'=ave('var',t+0,t='ten',1yr)'
 'modify 'varclim' seasonal'
 
 'set t 'tst''

pinterp.gs

GrADSスクリプト関数

gsf(GrADS Script functions)

公開されているpinterp.gsは拡張子がgsですが、pinterp.gsfという名前にする必要があります。pinterp.gsしかない場合

cp (GASCRPで指定しているディレクトリ)/pinterp.gs (〃)/pinterp.gsf

として下さい。

800hPa面の値を表示

 rc=gsfallow("on")　　　# GrADSスクリプト関数を使用
 'open var.ctl'
 'set x 1 ???'
 'set y 1 ???'
 'set lev 1000 700'　　　　# 800hPaを挟む二つのレベルを含むように範囲を設定
 'define a='pinterp(var, lev, 800)　　# 800hPa面に内挿したvarをaと定義
 'set lev 800'　　　　# 適当に高さを一つに設定
 'd a'

※必ずlevが減っていく座標でないと使用できません。

GrADS用「私」作ユーティリティー

set gxout stat

 d var         # 変数を表示
 pstat var     # 画面上部とターミナルウィンドウに統計量（最小値、最大値、
                 平均値、標準偏差）を表示

Bin Guan script libraly

t-6からt+6までの13ステップ移動平均

 rmean var -6 6 a
 d a

Bin Guan script libraly

t=1からt=100までの期間の標準偏差

 set t 1 100
 rms var a
 set t 1
 d a

IT memo

変数varの結果をnum列に並べてアスキー形式のファイルhoge.txtに出力

textout num hoge.txt var

あらかじめ、grd2txt.fの実行ファイルgrd2txtを

$GADDIR

で設定されているディレクトリに置いておく必要がある。

http://cookbooks.opengrads.org/index.php?title=Rec...

zinterp.gs

GrADSスクリプト関数

gsf(GrADS Script functions)

公開されているzinterp.gsは拡張子がgsですが、zinterp.gsfという名前にする必要があります。zinterp.gsしかない場合、

cp (GASCRPで指定しているディレクトリ)/zinterp.gs (〃)/zinterp.gsf

として下さい。

高度1,000mの値を表示

 rc=gsfallow("on")　　　# GrADSスクリプト関数を使用
 'open var.ctl'
 'set x 1 ???'
 'set y 1 ???'
 'set lev 0 2000'　　　　# 1,000mを挟む二つのレベルを含むように範囲を設定
 'define a='zinterp(var, lev, 1000)　　# 1,000mの高度に内挿したvarをaと定義
 'set lev 1000'　　　　# 適当に高さを一つに設定
 'd a'

※必ずlevが増えていく座標でないと使用できません。
※海洋のようにlevの値が減っていく座標のときは、ctlファイルに書かれたlevの順番を逆にし、options zrevを表記すると計算できます。

 am　# アジアモンスーン域
 eu　# ユーラシア大陸・大西洋
 pi　# 熱帯太平洋・インド洋
 pj　# PJパターン・南アジア
 po　# 熱帯太平洋
 wnp　# 西部北太平洋
 zen　# 全球

必要に応じてgsファイルを自作して置いてみましょう。Fortranのサブルーチンのような使い方ができます。