006, 160, 50, 0, 0, W, 0, , 0, 0, INVEST, S,
WP, 14, 2005082906, , BEST, 0, 151N, 1532E, 25, 1002, TD, 0, , 0,
0, 0, 0, 1006, 160, 50, 0, 0, W, 0, , 0, 0, INVEST, S,
WP, 14, 2005082912, , BEST, 0, 151N, 1523E, 30, 1000, TD, 0, , 0,
0, 0, 0, 1005, 180, 45, 0, 0, W, 0, , 0, 0, NONAME, S,
WP, 14, 2005082918, , BEST, 0, 151N, 1512E, 35, 997, TS, 34, NEQ,
30, 40, 40, 30, 1004, 190, 20, 0, 0, W, 0, , 0, 0, NABI, M,
WP, 14, 2005083000, , BEST, 0, 150N, 1501E, 45, 991, TS, 34, NEQ,
80, 60, 50, 70, 1002, 190, 25, 0, 0, W, 0, , 0, 0, NABI, M,
WP, 14, 2005083006, , BEST, 0, 149N, 1490E, 55, 984, TS, 34, NEQ,
70, 70, 70, 70, 1002, 190, 15, 0, 0, W, 0, , 0, 0, NABI, M,
WP, 14, 2005083006, , BEST, 0, 149N, 1490E, 55, 984, TS, 50, NEQ,
20, 20, 20, 20, 1002, 190, 15, 0, 0, W, 0, , 0, 0, NABI, M,
WP, 14, 2005083012, , BEST, 0, 150N, 1480E, 65, 976, TY, 34, NEQ,
100, 85, 70, 75, 1000, 200, 10, 0, 0, W, 0, , 0, 0, NABI, M,
WP, 14, 2005083012, , BEST, 0, 150N, 1480E, 65, 976, TY, 50, NEQ,
20, 20, 20, 20, 1000, 200, 10, 0, 0, W, 0, , 0, 0, NABI, M,
WP, 14, 2005083018, , BEST, 0, 152N, 1471E, 75, 967, TY, 34, NEQ,
100, 110, 110, 90, 1003, 250, 10, 0, 0, W, 0, , 0, 0, NABI, D,
......
1,
2,
3,
4,
5,
6,
7,
8, 9,
10
先用NCL读取这个文本文件，提取出绘制台风路径最需要的信息。
由于这个文本文件各个变量间是由逗号隔开，用字符串形式将文本信息读入，然后用str_get_field函数读取数据可能是最方便的。
;********************************************************
; Load NCL scripts
load
'$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl'
load
'$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl'
load
'$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl'
;********************************************************
begin
fiTY = 'bwp142005.txt'
; 获取文本文件的行数，相应的还有numAsciiCol函数用于获取列数
nrow = numAsciiRow(fiTY)
YYYYMMDDHH = new(nrow, 'string')
lat = new(nrow, 'float')
lon = new(nrow, 'float')
vmax = new(nrow, 'integer')
mslp = new(nrow, 'integer')
data = asciiread(fiTY, -1, 'string')
YYYYMMDDHH = str_get_field(data, 3, ',')
lat = stringtofloat(str_get_field(data, 7, ',')) *0.1
lon = stringtofloat(str_get_field(data, 8, ',')) *0.1
vmax = stringtoint(str_get_field(data, 9, ','))
mslp = stringtoint(str_get_field(data, 10,
','))
DateChar = stringtochar(YYYYMMDDHH)
MM = chartostring(DateChar(:,5:6))
DD = charactertostring(DateChar(:,7:8))
HH = charactertostring(DateChar(:,9:10))
end
|
台风路径数据用一行存储一个时刻的所有信息，所以这里先获取文本文件的行数，然后根据文本行数预先定义变量，然后
asciiread一次将所有数据读入。
这里data得到的是一个字符串(string)类型的一维数组，而这种数据类型与FORTRAN的“字符串”是不同的，这是NCL与FORTRAN的差别之一：
比如“Hello world!”这个句子，NCL中，n='Hello
world!'，n为string类型，大小为8bytes (类似于一个指针)；FORTRAN中，f='Hello
world!'，f为character的类型的数组，大小由字符串的长度决定。
这种差别带来的最大的影响就是NCL中不能对string类型的变量做子串操作(即无法直接取出字符串中的某部分)，而FORTRAN可以，f(3:5)就是“llo”。因此，如果需要在NCL中做子串操作，需要
stringtochar，然后再
charactertostring，这就带来了一些麻烦。幸运的是，NCL有专门的字符串函数，使用这些函数就能完成很多复杂的操作。
str_get_field是NCL的字符串函数之一，是由用户自己指定一个分隔符，然后取出分隔后的某一段字符串。
比如这里，每一行的字符串用逗号分隔开后，第3段是时间，第7段是纬度，第8段经度，第9段风速，第10段为海平面最低气压。在分段读取数据后进行变量类型转换，便得到了需要的信息。最后再对时间字符串进行子串操作，获取了月日和时的字符串。
月、日、时的信息也可以用下面的方式从YYYYMMDDHH变量得到：
HHi = mod(stringtoint(YYYYMMDDHH), 100) DDi =
mod(stringtoint(YYYYMMDDHH)/100, 100)
MMi = mod(stringtoint(YYYYMMDDHH)/10000, 100)
|
注意这里时间都是整形变量。
有这些信息，就能画出台风路径了。
画图部分的代码如下：
先绘制底图，然后绘制线和点
;.....前略....
; plot
wks = gsn_open_wks('ps', 'BestTrack')
res = True
res@gsnDraw = False
res@gsnFrame = False
res@mpMinLatF = 10
res@mpMaxLatF = 50
res@mpMinLonF = 110
res@mpMaxLonF = 160
plot = gsn_csm_map_ce(wks,res)
resDot = True
resLine = True
dumDot= new(nrow, graphic)
dumLine = new(nrow, graphic)
; 绘制线
resLine@gsLineColor = 'black'
do i = 0, nrow-2
xx = (/ lon(i), lon(i+1)/)
yy = (/ lat(i), lat(i+1)/)
dumLine(i) = gsn_add_polyline(wks, plot,
xx, yy, resLine)
end do
; 绘制00时的点
resDot@gsMarkerColor = 'black'
resDot@gsMarkerIndex = 1
resDot@gsMarkerSizeF = 0.02
do i = 0, nrow-1
if (HH(i) .eq. '00') then
dumDot(i) =
gsn_add_polymarker(wks, plot, lon(i), lat(i), resDot)
end if
end do
draw(wks)
frame(wks)
|
完成后的结果就是这样：
![[转载]NCL绘制台风路径(1)](http://s11.sinaimg.cn/bmiddle/a1345723gc0c3c9b9d09a&690 "[转载]NCL绘制台风路径(1)")
