當您在分析資料時，清晰的圖表，是表現成果最有效的方式，雖然市面上已有許多相當方便的科學繪圖軟體，但對大氣科學上而言，其仍有所不足之處，例如等值線很難看、風標及流線場無法畫、通用性及價錢等考量。目前氣象上最通用的就屬NCAR GRAPHIC，其雖可滿足我們的需求，但每次都要寫程式，又有一堆堆的設定，對一般的"懶人"及要求時效的人來說，是相當麻煩的。而Grads 是一套專門為氣象方面所用的繪圖軟體，優點為:
- 簡單，只要給一定的BINARY DATA，不須寫程式，只需線上下指令即可。
- 方便，目前在各種平台工作站如 DEC, IBM, HP, Linux 等都可支援(要X-WINDOWS)，PC下DOS亦可，目前亦準備有Windows的版本。
- 易學，您也許不相信，看這demo，把您的常用資料放進來試試看，一小時以內即可畫出來。
- 經濟，目前是freeware。
雖有以上的便利，可能是過於便利，因此缺少了一些圖形的細節控制，例如許多都已是default，要改很難，筆者建議您若是做研究分析，時間就是一種競爭，用Grads是很好的工具來協助您，但是若成果要發表成論文時時，NCAR GRAPHIC還是較好的選擇，除非你不大計較這些東西。
目前各工作站及PC幾乎皆可支援，詳情可見
http://grads.iges.org/grads
或中央大學大氣系軟體資料庫
ftp://www.atm.ncu.edu.tw/grads/
內有最新訊息及hypertext,ps,text,doc檔的說明，約150頁。(為節省資源，用WWW看就好)
其實根本不需要印，只要看看 grads reference card version就可了，
ftp://www.atm.ncu.edu.tw/grads/manual/
ftp: grads.iges.org
您要最新的版本可在此發現。
開始學 GRADS 可分為六個東西
- 資料處理程式: 將資料格式轉成 BINARY檔 。
- OPEN(8,FILE='samp.dat',FORM='UNFORMATTED',ACCESS='DIRECT',
RECL=100)
WRITE (8,REC=1) X
- Binary 檔: 可為 *.bin。
- 描述 data 格式的檔: *.ctl。
- 繪圖控制檔: *.gs。
- 轉圖成為 postscript。(gxps, gxpsc, gxpscw)
- 印 postscript圖。(用postscript印表機或用ghostview)
若您在PC下請在autoexec.bat內加入兩個路徑
path .........;c:\grads
append ......;c:\grads
若您使用 PC 版本，妳可利用網路芳鄰連至 ATMOS 網域內，選擇 SOFTWARE-FREE當成妳的 H 碟，此時在 H:\下鍵入 RUN 此時包括 FORTRAN 設定及 GRADS 設定皆已完成。
執行grads就可以進入 Grads繪圖模式內，出現Grads的提示符號鍵入
ga>run maps.gs
你就可以看見一些DEMO的圖，您要看更多可打test即可(在PC版才有)
您在IBM上以X-windows內打grads即可。
GRADS 繪圖樣式有兩種選擇，landscape(-A4) 及 portait(A4)在打入GRADS 後即要妳選擇繪圖樣式，妳也可加入一些指令簡化流程，
-b run batch mode 此時不會秀出任何圖形
-l run in landscape mode
-p run in portrait mode
-c 執行後者所提供的 command 指令，例如 *.gs
ex:
grads -c "run batch.gs"
grads -blc "run batch.gs"
原始範例
我們看看 maps.gs 內容，'reinit'，reinit是對系統的設定值重新定義，而 ’ ’內則表示這是在Grads繪圖狀態下的指令，相對於若你不以run maps.gs的方式進行，直接在 ga>reinit 一行一行鍵入也是可以。接著在maps.gs內開了三個檔案
'open sa.ctl'
'open ob.ctl'
'open wx.ctl'
在*.ctl 檔內包含了二個很重要的部分，一是定義圖形的格式，另一個是定義所要繪出的變數，我們來看 sa.clt 內容:
DSET sa.dat
TITLE 0.5 degree analyzed grids
UNDEF -9.9E33
XDEF 81 LINEAR -110 0.5
YDEF 51 LINEAR 25 0.5
ZDEF 1 linear 1 1
TDEF 1 LINEAR 16z16feb1993 1hr
VARS 5
slp 0 99 sea level pressure
ts 0 99 surface temp
ds 0 99 surface dwpts
us 0 99 surface u winds
vs 0 99 surface v winds
ENDVARS
一開始的 DSET sa.dat 是定義你所要用的資料檔檔名,
VARS 5
slp 0 99 sea level pressure
ts 0 99 surface temp
ds 0 99 surface dwpts
us 0 99 surface u winds
vs 0 99 surface v winds
ENDVARS
表示有五個變數，內含的就是定義所有變數的名稱和格式，例slp 0 99 sea level pressure 其中0和99表示時間和空間都只有一層，而 sea level pressure是變數的說明（若時間和空間不止一層則可用相關loop 指令做成動畫，在此不做介紹）。實際上應寫成
Time 1,Level ? ,Variable slp
即第一個是代表時間，第二個是代表層數（編者不知99的意義，Grads.doc內沒有說明詳細，大概是預留的未知個數)。
TITLE這僅是對此檔案的說明，幫助了解這一個檔案的作用。
UNDEF無效值或可忽略的資料值，在此為-9.9E33。
XDEF 表示經度上Ｘ軸的點數，在此有81個點，起始值由-110開始，間值是0.5，增加的方式是線性(Linear)。若不是線性，則用 LEVELS -110,-100,80,......，逐一列出。
YDEF 緯度上......和上面相同。
ZDEF 高度上......和上面相同。
TDEF 時間上......和上面相同。但時間表示的格式是 hh:mmZddmmmyyyy
hh:hour
mm:minutes
dd:day
mmm:month
yyyy:year
例如：12Z1JAN1990,14:20Z22JAN1987等。
做完*.ctl後，如何把資料變數放進去，把變數場繪出來。底下先說明一維陣列的方式
REAL X(100)
DO 10 I=1,100
X(I)=I
10 CONTINUE
OPEN(8,FILE='samp.dat',FORM='UNFORMATTED',ACCESS='DIRECT', RECL=100)
WRITE (8,REC=1) X
STOP
END
必須注意的是RECL，在PC中實數是4個byes因此必須乘上4，才是所須要的大小。底下則是一個三維矩陣
三維矩陣的排列方式：
REAL Z(72,46,16)
OPEN (8,FILE='grads.dat',FORM='UNFORMATTED', ACCESS='DIRECT',RECL=72*46)
IREC=1
DO 10 I=1,16
WRITE(8,REC=IREC) ((Z(J,K,I),J=1,72),K=1,46)
IREC=IREC+1
10 CONTINUE
ob.ctl，在Grads內對於測站圖形的處理相當的美觀，因此如何將測站的資料融入圖形內是很重要，對此Grads自有一套格式來解決。
dtype station
stnmap ob.map
Year Month Stid Lat Lon Rainfall
1980 1 QQQ 34.3 -85.5 123.3
1980 1 RRR 44.2 -84.5 87.1
1980 1 SSS 22.4 -83.5 412.8
1980 1 TTT 33.4 -82.5 23.3
1980 2 QQQ 34.3 -85.5 145.1
1980 2 RRR 44.2 -84.5 871.4
1980 2 SSS 22.4 -83.5 223.1
1980 2 TTT 33.4 -82.5 45.5
.
.
.
A FORTRAN program in DEC FORTRAN to write this data set
in GrADS format might be:
CHARACTER*8 STID
C
OPEN (8,NAME='rain.ch')
OPEN (10,NAME='rain.dat',FORM='UNFORMATTED',
& RECORDTYPE='STREAM')
C
IFLAG = 0
C
C Read and Write
C
10 READ (8,9000,END=90) IYEAR,IMONTH,STID,RLAT,RLON,RVAL
9000 FORMAT (I4,3X,I2,2X,A8,3F8.1)
IF (IFLAG.EQ.0) THEN
IFLAG = 1
IYROLD = IYEAR
IMNOLD = IMONTH
ENDIF
C
C If new time group, write time group terminator.
C Assuming no empty time groups.
C
IF (IYROLD.NE.IYEAR.OR.IMNOLD.NE.IMONTH) THEN
NLEV = 0
WRITE (10) STID,RLAT,RLON,TIM,NLEV,NFLAG
ENDIF
IYROLD = IYEAR
IMNOLD = IMONTH
C
C Write this report
C
TIM = 0.0
NLEV = 1
NFLAG = 1
WRITE (10) STID,RLAT,RLON,TIM,NLEV,NFLAG
WRITE (10) RVAL
GO TO 10
C
C On end of file write last time group terminator.
C
90 CONTINUE
NLEV = 0
WRITE (10) STID,RLAT,RLON,TIM,NLEV,NFLAG
STOP
END
另外一個檔案wx.ctl 則依此類推。
常用功能摘要(請見附表)
Graphic Type
set gxout bar BAR圖
set gxout winbarb winbarb
set gxout contour 畫等值線
set gxout shaded 陰影
set gxout vector 風向量圖
set gxout sream 氣流線
set gxout line 直條圖
set gxout fgrid Shaded grid boxes
Graphic options
set csmooth on|off 是否用cubic interpolation內差
set cmin 0 最小值
set cmax 20 最大值
set cint 5 距隔值
set clevs lev1 lev2 .... 等值線值
Page control 圖形大小控制
Map projections/drawing 地圖投影
'set mproj latlon' 定義值，經緯度直角座標。
scaled 緯度投影。
nps north ploar stereographic．
sps south polar sterographic．
off 不做地圖投影。
'set mpvals -105 -75 25 50' 設定繪出座標的範圍，在此表示105W~75W,25N~50N。
'set mpdset lowers' 內定值，
mres 表示有州界
hires 表示含國界
nmaps
'set poli on|off 內定值是on ，表示上一項有否要作用。
筆者範例 ftp://www.atm.ncu.edu.tw/grads/example/1
以上是照原文翻過來的方法，在此筆者用一個自己嘗試的範例來說明，筆者的資料是環保署桃園站一個月的資料，目標是看其一個月的每日變化，筆者用上述 z(i,j,k)的寫法將原始環保署data成了一個 hs9417.bin檔因此描述此檔的file可寫為 epa.ctl
DSET hs9417.bin
TITLE model data
UNDEF -999.0
XDEF 24 LINEAR 0 24 -----> x軸為24小時，24點
YDEF 31 LINEAR 0 31 -----> y軸為31日，31點
ZDEF 1 linear 1 1 -----> 沒有z故以1 表示
TDEF 1 LINEAR 12z4feb1994 12hr
VARS 16 -----> 16個變數
so2 0 99 precipitable water ----->依順序為16個變數，內容懶得寫故都一樣
co 0 99 lifted index
o3 0 99 sea level pressure
pm 0 99 precip
nox 0 99 z1000
no 0 99 z850
no2 0 99 t850
thc 0 99 u850
nmhc 0 99 v850
wspeed 0 99 rh850
wdir 0 99 z500
U 0 99 z500 ----> 畫風標其只認得分量
V 0 99 z500
temp 0 99 u500
rainfall 0 99 rain
ch4 0 99 ch4
ENDVARS
接下來可畫圖了，進入grads後
ga>open plot.ctl
您可看到 lon .....,lat.....等，因為其會將其當成經緯度來看，我們先不管，座標軸示等一下說明。可以先用query file來顯示plot.ctl內宣告的內容，若此時要畫 O3的月內日變化圖可直接 ga>d o3，此即可繪出，但如前所述，其會將其當成經緯度來看，故座標會怪怪的，因此可ga> set xaxis 1 24 1，表x軸圖示為1至24間隔1，同時還有世界地圖在內，故須 set mpdraw off，再將網格去掉 set grid off，clear一下，再d 03就可畫出一張不錯的圖，若要看每日中午12時的 O3日變化，只需set x 12，再 d O3 即可看到一個剖面圖。當要畫好幾個圖時，這樣太麻煩了，可以寫成批次檔如下所示
epa.gs
'reinit'
'clear'
'open plot.ctl'
'enable print test' ---->將圖形初出至test檔
'set xaxis 0 24 1'
'set yaxis 0 31 1'
'set mpdraw off'
'set csmooth on' -----> smooth 點
'set grid off'
'draw title EPA 17 O3'
'draw xlab Hour'
'draw ylab Day(May)'
'set t 1'
'set z 1'
'd o3'
'print' ----> 印在test上
pull dummy ----> 有點像call wait
'clear'
-----> 下一張
'set xaxis 0 24 1'
'set yaxis 0 31 1'
'set mpdraw off'
'set csmooth on'
'set grid off'
'draw title EPA 17 O3'
'draw xlab Hour'
'draw ylab Day(May)'
'set t 1'
'set z 1'
'd u;v' -----> 風的 vector
'print'
pull dummy
'clear'
'disable print' -----> 結束輸出圖形
故只要 ga> run plot.gs即可做完一堆工作，此時會產生test，將其用gxps轉成postscript檔，在pc上用ghostview，工作站用xv或 ghostview即可看出。
ftp://www.atm.ncu.edu.tw/grads/example/2
這是幫我老闆投到 JGR 的PAPER所弄的，妳可去最近的1997 PAPER上找到，特殊的是為非等間距，無效值為陰影區及三張圖疊在一起，不過要加一些特殊符號還是很麻煩的，所以GRADS還是有其侷限性!
ftp://www.atm.ncu.edu.tw/software/grads/example/3
由於我不是做大尺度的，這是謝榮傑所提供的，主要是在 ibm與 dec機器在存取binary 時之差異，妳可傳至 IBM試試
在 IBM, HP, SUN 的位元數值，是採取 4BYTES=32BIT的存取位元資料，是IEEE的標準，其是屬於 big endian，若是在 ms-dos下開資料時用 *4則可適用於 big endian.
在 DEC, PC linux則是 little endian，與big endian 不同，但在dec 上，也可讀 big endian code，只要在 open 這個敘述中加入 convert='big_endian' 即可
open(10,file='10.bin',status='old',form='unformatted',convert='big_endian')
ftp://www.atm.ncu.edu.tw/grads/example/4
這是東亞地區硫化物排放之模擬。
ftp://www.atm.ncu.edu.tw/grads/example/skew
斜溫圖，是氣象很獨特的一種圖形，一般而言，目前除了 NCAR Graphic及PC上的繪圖軟體，一般是相當難繪製的。但皇天不負苦心人，賓州大學氣象系研究生Bob Hart (Graduate Student, Penn State University Meteorology, http://www.ems.psu.edu/~hart/skew.html )終於將 Grads 版斜溫圖寫出來，並公開予大眾使用。
目前經測試，其函數及各項數值均為正確，但目前只能在工作站上 run，pc部份仍有些問題待解決。
關於斜溫圖的程式範例詳見 ftp://www.atm.ncu.edu.tw/grads/example/skew
須注意的是風速可以設一個很小值，絕不能為0，這是許多錯誤嘗試的經驗！
