GrADS繪圖學習技巧與實例

1、nanjing university of information science& technology 阿木 以下技巧總結(jié)都是筆者從學習實踐過程中總結(jié)出來的，基本的問題。不求全面，希望對讀者學習有用，如果有問題，敬請留言指正，以促進交流學習！1、 軟件綜述：grads軟件是一款繪圖軟件除了繪制圖形，還可以提取數(shù)據(jù)，主要應(yīng)用是在大氣科學中，當然只要是數(shù)據(jù)處理成grads能夠讀取的數(shù)據(jù)文件就可以進行相關(guān)繪圖。軟件版本問題，軟件本身不是很大，我接觸到1.8、1.9、2.0版本的，1.8版本的安裝很多情況還要修改環(huán)境變量、1.9版本的不識別sdfopen命令，最穩(wěn)定的版本是2.0版本，所以筆者推薦

2、學習者安裝2.0版本，選擇默認安裝路徑就可以。2、 文件類型簡述 ：grads處理的是網(wǎng)格數(shù)據(jù)，可以處理的數(shù)據(jù)類型有：grd、grib、nc（海洋常用的數(shù)據(jù)），cdf（雷達衛(wèi)星數(shù)據(jù)），其中nc、cdf數(shù)據(jù)都是自帶描述文件，不需要ctl，grib數(shù)據(jù)要通過命令生成ctl、index數(shù)據(jù)才可以調(diào)用，常用的是grd數(shù)據(jù)，需要ctl。3、 數(shù)據(jù)文件轉(zhuǎn)換：grads軟件識別的數(shù)據(jù)是二進制無格式數(shù)據(jù)，文件類型是binary，寫入和生成時是不需要格式的如 read(20) sst(i,j,iz,it),20為文件號，通常是十進制數(shù)據(jù)與grd數(shù)據(jù)間轉(zhuǎn)換，這里給一個grd轉(zhuǎn)換成txt數(shù)據(jù)的fortran程序：p

3、arameter(nx=56,ny=41,nz=1,nt=360)dimension sst(nx,ny,nz,nt)real sstopen(15,file=sst.grd,form=binary) ！固定的用form=binary就是二進制數(shù)據(jù)open(16,file=sst.txt) ！新建txt文件do it=1,nt do iz=1,nz read(15) (sst(i,j,iz,it),i=1,nx),j=1,ny) ！read后只有文件號，數(shù)據(jù)是無格式的 enddoenddodo it=1,nt do iz=1,nz write(16,*) (sst(i,j,iz,it),i=1

4、,nx),j=1,ny) ！輸出時是txt文件可直接看的數(shù)據(jù)，有格式輸出，有 * enddoenddoclose(15)close(16)end 寫程序時：注意格點數(shù)要與數(shù)據(jù)對應(yīng)，如：上程序?qū)?yīng)的數(shù)據(jù)是 經(jīng)度90200，緯度-2060，時間：1971.012000.12共360個月的海面溫度數(shù)據(jù)，數(shù)據(jù)格點精度 2*2 ，nx=(200-90)/2+1,ny=(60-(-20)/2+1,nt=360，nz=1，大氣的數(shù)據(jù)要根據(jù)數(shù)據(jù)的層次確定幾層。4、 grd 、ctl、gs、nc詳述grd文件：grd數(shù)據(jù)不可直接看，為二進制無格式數(shù)據(jù)，簡單的說只有1和0，而且數(shù)據(jù)間沒有間隔，grads識別grd

5、是根據(jù)ctl進行劃分的，根據(jù)ctl中的經(jīng)度、緯度、層次、時間，精度進行數(shù)據(jù)分塊。ctl實例：dset c:datasst.grdundef -9.99e+33title sea surface tempturexdef 56 linear 90 2ydef 41 linear -20 2tdef 360 linear jan1971 1mozdef 1 levels 0vars 1sst 0 99 surface sea tempture endvars上例數(shù)據(jù)的數(shù)據(jù)順序是（以下是數(shù)據(jù)對應(yīng)的經(jīng)緯度）緯度 經(jīng)度（t=1）-20 90200的 緯度20s的從90e200的56個數(shù)據(jù)-18 9020

6、0的 緯度18s的從90e200的56個數(shù)據(jù).60 90200的 緯度60n的從90e200的56個數(shù)據(jù)以上為一層的數(shù)據(jù)，接下來是t=2，t=3t=360的數(shù)據(jù)，每個時間點的每一層是如上格式，編程時讀取和寫如的數(shù)據(jù)循環(huán)順序依次是：時間、層次、緯度、經(jīng)度，讀者參照3中的fortran程序加以理解。ctl文件：具體其他指導書上都有，我這里強調(diào)的是 sst 后面的0表示一層，如果是兩層以上則是2,3，1層是固定用法，sst后面的99是默認設(shè)置；undef -9.99e+33 此處的數(shù)值決定了軟件將文件中的那些值認定為不繪制的數(shù)據(jù)，所以這個值一定要與數(shù)據(jù)對應(yīng)。gs 文件，批命令文件，與ctl一樣是用記

7、事本編寫，另存為.gs文件，文件都是命令，方便大段的命令編寫、修改，很常用，尤其時繪制的圖要求比較多時必須用，免得在命令窗口重復(fù)輸入命令浪費時間，如下例子：reinitopen c:datauv.ctlset t 7set lev 850set grid offset vpage 0 8.5 0 3.6set parea 0.3 8.5 0.3 3.5set gxout vectorset grads offd u;vset vpage 0 8.5 3.6 7.2set grads offset gxout barbd u;vset vpage 0 8.5 7.2 11set grads of

8、fset gxout gridd u;vprintim c:imagesgxout3uv850.png white; gs文件以 ；結(jié)尾，分號后面不能有空格，除了循環(huán)命令外都需要 單引號 將命令引起來。nc文件：此類數(shù)據(jù)不需要ctl可直接用sdfopen命令直接打開，往往需要知道數(shù)據(jù)文件中的各個緯度特征，可以用：q ctlinfo 命令查詢該文件的內(nèi)置ctl，這樣一切都會很清楚，尤其是數(shù)據(jù)精度：截圖如下（下例即2*2的網(wǎng)格經(jīng)度）：以上為基礎(chǔ)知識介紹，以下為技巧命令1、 繪圖時，圖的時間下標可用 set grads off 命令關(guān)掉，網(wǎng)格用 set grid off 命令關(guān)掉，需注意的是，網(wǎng)格關(guān)

9、一次就一直有效，而下標則是每次繪完圖就自動開啟，所以建議讀者在每次的 繪圖命令之前加set grads off 。2、3、 reinit命令是讓窗口恢復(fù)到剛打開時的界面，會恢復(fù)所有的set，同時關(guān)掉所有的打開文件，為為防止前面打開文件的干擾，建議讀者在每個gs文件的第一條就加上這條命令（參考gs描述的例子）。4、5、 數(shù)據(jù)維度設(shè)定： set lon 90 set x 1 set lat 80 set y 20 set lev 1000 set z 1 set t 1 12 set time jan1970 dec1970 上面兩種設(shè)定等價，左邊是實際維度設(shè)定法，右邊是給點設(shè)定法，讀者需要確定具

10、體的格點數(shù)。需要注意的是，在設(shè)定全球尺度時，經(jīng)度0和360是同一個格點，所以 set 0 360會出錯，這時設(shè)定格點的方法比較好：set x 1 180。6、 vpage和parea的區(qū)別：vapge是對整個繪圖區(qū)分塊，需要幾張圖就劃分成幾塊，給英寸時不需要給標注、標題預(yù)留空間，如：要橫著繪制兩張圖 那么就是set vpage 0 5.5 0 8.5（11*8.5） 繪制第一張圖 set vpage 5.5 11 0 8.5 繪制第二張圖parea 是描述的繪圖時圖形的四根邊線的大小，而且是虛頁的尺寸如：畫一張圖 set parea 1.0 10.2 0.8 7.8 四個值的范圍 取決于vap

11、ge的長度01.010.211，如果是上例中的兩張圖，那么：set parea 0.5 4.8 0.8 7.9，00.54.85.5 ，數(shù)值范圍取決于水平、垂直的長度，與起點無關(guān)如 set vpage 0 11 3.5 7 set parea 0.8 10 0.3 3.2 y上滿足00.33.2(7-3.5) 即可，實際y長度決定范圍，與起點無關(guān)。 5、set annot 5 8set xlopts 3 5 0.18set ylopts 3 5 0.18上面這兩條命令是設(shè)定的x、y軸下標數(shù)字的顏色、粗細、字號set annot 5 8是設(shè)定坐標軸線、標題的顏色 粗細，會重置xlopts的部分設(shè)

12、定，希望讀者注意，命令間的互相干擾6、cbar cbarn命令cabr cbarn 這里只講解cbarn 命令，cbarn sf vert xmid ymid 其中sf 為標尺,1為全尺寸0.5為半尺寸；vert為放置位置，0為水平，1為垂直；xmind 、ymid為色標的中心位置如：cbarn 0.5 0 2.5 1.5 色標 半尺度長，水平放置，色標的中心英寸坐標（2.5,1.5）。如果只是cabrn或cbar命令，會按照默認的全長，水平或數(shù)值取決于圖的哪邊空位大繪制。強調(diào)一點，加色標和加標題都是在繪制出圖形后才能加，因為色標是根據(jù)陰影圖確定對應(yīng)色值，標題根據(jù)圖形大小確定標題位置。7、 c

13、min、cmax是命令是用于繪制大于或小于某數(shù)值的線或區(qū)域圖形，在每次繪圖之后會重設(shè)，所以如果多次使用一定要每次d之前加上，這點與set grads off相似。8、求12個月每個月的海溫距平值set t 1 12asst=ave(sst,t+1404,t=1764,12)modify asst seasonalset time jan1971d asstset t 1405 1764nasst=sst-asstset time jan1998d nasstset t 1 12 是設(shè)定asst變量有12個時間序列，每個時間格點放一個平面的平均值，asst=ave(sst,t+1404,t=17

14、64,12) 定義變量asst放每個月的平均值，隨著t從1變化到12，一次求的每個月的平均值，起始時間是t=1405，終止時間是1764，t=1時，相當于1+13+25+37，即每年的1月份的值求平均（這里t=1+1404=1405是1971年1月，因此是求19712000年每年1月份的平均值），t=2,3.4,，12 與1同理。如此將12個月的平均值都放在了asst里。modify asst seasonal 本來12個月的平均值知識放在了112的時間序列里，（這里調(diào)用的文件數(shù)據(jù)開始時間是1854年1月），所以asst的值只是在設(shè)定112以內(nèi)才能畫出，假了這條命令，可以將asst的時間序列擴

15、展到所有時間里，使得每一年的每個月對應(yīng)都是該月30年的平均值，方便后面求距平。需要注意的，如果set t 1 4 ，那么這條命令的作用是每年的asst的14月值是一樣的，是該月的平均值。由于grads本身軟件有一定的問題，當你不是設(shè)定一段時間而是一個時間點如：set t 1，步、不用modify命令，所有的時間序列都會有該值set t 1405 1764nasst=sst-asst這兩句是求30年（1971 012000 12）每個月的海溫距平值，共360各月的距平值。這個不難理解，不做解釋。以下為實例：1、利用所提供的數(shù)據(jù)文件，繪制出2003年7月60-150e、0-40n區(qū)域內(nèi)700hpa

16、流線圖，且地形高度場超過2000米以上用黑色陰影顯示，并且給出相應(yīng)標題（請包含姓名拼音與學號），最終將圖形保存。所有命令編寫于.gs文件中。reinitopen c:datadxgd.ctlopen c:datauv.ctlset grid offset grads offset map 3 1 5set xlopts 3 5 0.18set ylopts 3 5 0.18set parea 0.7 10.4 0.6 7.8set t 1set lon 60 150set lat 0 40 set gxout shadedset rbcols 1 1 1 1 1 1 1set cmin 200

17、0d hset gxout streamset time jul2003set grads offset lev 700d u.2;v.2set annot 5 8draw title dingxiaoli 20081331001printim c:imagesuv200307850.png white2、利用所提供的數(shù)據(jù)文件，畫出2002年1-12月120 e、 0-40n 200hpa緯向風的緯度-時間剖面圖。圖型要求：（1）緯向風為西風時填色，東風繪制等值線，給出色標，0值線加粗；（2）x軸標注為“time”,y軸標注為“l(fā)at”,標題標注為“u 100-120e”。（3）最終將圖形保存

18、為gmf格式。所有命令編寫于.gs文件中reinitopen c:datauv.ctlset grid offset grads offset xlopts 3 5 0.16set ylopts 3 5 0.16set parea 1 10.4 1.2 7.9set lon 120set lat 0 40set lev 200set t 1 12set annot 5 8enable print c:imagesuv200.gmfset xyrev onset gxout shadedset cmin 0d ucbarn 1 0 5.5 0.3set grads offset gxout co

19、ntourset cmax 0set grads offd uset clevs 0set ccolor 2set cthick 10set grads offd udraw xlab timedraw ylab latdraw title u 100-120eprintdisable print;3、十二個月的nc數(shù)據(jù)整合到一個文件中（fortran）這個程序可用于整合數(shù)據(jù)parameter(nx=360,ny=181,nz=26,nt=12)dimension temp(nx,ny,nz,nt)open(1,file=argo_200901.grd,form=binary)open(2,f

20、ile=argo_200902.grd,form=binary)open(3,file=argo_200903.grd,form=binary)open(4,file=argo_200904.grd,form=binary)open(5,file=argo_200905.grd,form=binary)open(6,file=argo_200906.grd,form=binary)open(7,file=argo_200907.grd,form=binary)open(8,file=argo_200908.grd,form=binary)open(9,file=argo_200909.grd,

21、form=binary)open(10,file=argo_200910.grd,form=binary)open(11,file=argo_200911.grd,form=binary)open(12,file=argo_200912.grd,form=binary)open(37,file=argo_2009.grd,form=binary)do iz=1,nz read(1) (temp(i,j,iz,1),i=1,nx),j=1,ny)read(2) (temp(i,j,iz,2),i=1,nx),j=1,ny)read(3) (temp(i,j,iz,3),i=1,nx),j=1,n

22、y)read(4) (temp(i,j,iz,4),i=1,nx),j=1,ny)read(5) (temp(i,j,iz,5),i=1,nx),j=1,ny)read(6) (temp(i,j,iz,6),i=1,nx),j=1,ny) read(7) (temp(i,j,iz,7),i=1,nx),j=1,ny) read(8) (temp(i,j,iz,8),i=1,nx),j=1,ny)read(9) (temp(i,j,iz,9),i=1,nx),j=1,ny)read(10) (temp(i,j,iz,10),i=1,nx),j=1,ny)read(11) (temp(i,j,iz,11),i=1,nx),j=1,ny)read(12) (temp(i,j,iz,12),i=1,nx),j=1,ny)enddodo it=1,nt do iz=1,nz write(37) (temp(i,j