- 积分
- 6920
- 贡献
-
- 精华
- 在线时间
- 小时
- 注册时间
- 2013-1-1
- 最后登录
- 1970-1-1
|
登录后查看更多精彩内容~
您需要 登录 才可以下载或查看,没有帐号?立即注册
x
本帖最后由 lysx 于 2014-8-29 19:29 编辑
今天看到清风大侠的转帖——“GrADS中站点分布图的画法【转】http://bbs.06climate.com/forum.php?mod=viewthread&tid=7386&fromuid=21651”,心动不已,然而在站点附近只能用汉语拼音来标记站点名称,不能不说是一个缺憾。幸好最近研究了grads下的汉字输出问题,于是自己试着也写了一段程序以实现同样的功能。
以下代码,由本人免费送上
PROGRAM station
USE DFPORT
IMPLICIT NONE
REAL :: lat, lon, temp
CHARACTER ( LEN = 20 ) :: sta
INTEGER :: k, status
OPEN ( 99 , FILE = 'cn_station.gs', STATUS = 'UNKNOWN', ACTION = 'WRITE', IOSTAT = status )
IF ( status .EQ. 0 ) THEN
WRITE ( 99, '(A)' ) "'reinit'"
WRITE ( 99, '(A)' ) "'enable print cn_station.gmf'"
WRITE ( 99, '(A)' ) "'set parea 0.3 10.7 0.5 8.0'"
WRITE ( 99, '(A)' ) "'open fnl_20130301_00_00.ctl'"
WRITE ( 99, '(A)' ) "'set grads off'"
WRITE ( 99, '(A)' ) "'set grid off'"
WRITE ( 99, '(A)' ) "'set xlopts 1 5 0.12'"
WRITE ( 99, '(A)' ) "'set ylopts 1 5 0.12'"
WRITE ( 99, '(A)' ) "'set xlab off'"
WRITE ( 99, '(A)' ) "'set ylab off'"
WRITE ( 99, '(A)' ) "'set font 1'"
WRITE ( 99, '(A)' ) "'set mpdset cnmap cnriver'"
WRITE ( 99, '(A)' ) "'set map 15 1 1'"
WRITE ( 99, '(A)' ) "'set lon 73 136'"
WRITE ( 99, '(A)' ) "'set lat 15 55'"
WRITE ( 99, '(A)' ) "'set cmin 10000000.0'"
WRITE ( 99, '(A)' ) "'d hgtprs'"
WRITE ( 99, '(A)' ) "'q gxinfo'"
WRITE ( 99, '(A)' ) "xx=sublin(result,3)"
WRITE ( 99, '(A)' ) "yy=sublin(result,4)"
WRITE ( 99, '(A)' ) "x1=subwrd(xx,4)"
WRITE ( 99, '(A)' ) "x2=subwrd(xx,6)"
WRITE ( 99, '(A)' ) "x=x1+(x2-x1)/2"
WRITE ( 99, '(A)' ) "y1=subwrd(yy,4)"
WRITE ( 99, '(A)' ) "y2=subwrd(yy,6)"
WRITE ( 99, '(A)' ) "y=y1-0.5"
WRITE ( 99, '(A)' ) "'axis -type b -position i -start 75 -end 135 -interval 10 -sinterval 5 -suffix `3.'"
WRITE ( 99, '(A)' ) "'axis -type l -position i -start 20 -end 50 -interval 5 -sinterval 5 -suffix `3.'"
WRITE ( 99, '(A)' ) "'axis -type t -position i -start 75 -end 135 -interval 10 -sinterval 5 -suffix `3.'"
WRITE ( 99, '(A)' ) "'axis -type r -position i -start 20 -end 50 -interval 5 -sinterval 5 -suffix `3.'"
OPEN ( 19, FILE = '160_station.txt', STATUS = 'OLD', ACTION = 'READ', IOSTAT = status )
IF ( status .EQ. 0 ) THEN
DO
READ ( 19, *, IOSTAT = status ) temp, temp, lat, lon, sta
IF ( status .NE. 0 ) EXIT
WRITE ( 99, '(A8, 2F10.4, A10 )' ) "'q w2xy ", lon, lat, "'"
WRITE ( 99, '(A)' ) "x=subwrd(result,3)"
WRITE ( 99, '(A)' ) "y=subwrd(result,6)"
k = LEN_TRIM ( sta )
WRITE ( 99, '(A18, A<K>, A15 )' ) " 'writehz 'x' 'y' ", TRIM ( sta ), " 2 1 2'"
WRITE ( 99, '(A)' ) "'draw mark 2 'x' 'y' 0.05'"
END DO
ELSE
WRITE ( *, '(/, A, /)' ) 'WARNING:: the file 160_station.txt is not opened!!!'
END IF
CLOSE ( 19 )
WRITE ( 99, '(A)' ) " 'southsea_last 0 11 0 8.5' "
WRITE ( 99 , '(A)' ) "'printim cn_station.png white X1366 Y768'"
WRITE ( 99, '(A)' ) "'print'"
WRITE ( 99, '(A)' ) "'disable print'"
WRITE ( 99, '(A)' ) "'quit'"
ELSE
WRITE ( *, '(/, A, /)' ) 'WARNING:: the file cn_station.gs is not opened!!!'
END IF
CLOSE ( 99 )
CALL SYSTEM ( 'grads -lbc cn_station.gs' )
END PROGRAM station
PS:红色代码为FORTRAN调用grads的方式
看看图吧:
|
-
评分
-
查看全部评分
|