关于台风中心点绘图求助……

陪着你的路人

唉~折腾了好几天了……换了好几次思路，实在是搞不懂……直接上程序吧……
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"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/shea_util.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/wrf/WRFUserARW.ncl"
undef("create_hurricane_symbol")
function create_hurricane_symbol(wks)
begin
  mstring = "p"
  fontnum = 37
  xoffset = 0.0
  yoffset = 0.0
  ratio   = 1.0
  size    = 1.0
  angle   = 0.0
  new_index = NhlNewMarker(wks, mstring, fontnum, xoffset, yoffset, \
                           ratio, size, angle)
  return(new_index)
end
;************************************************

begin

  a = addfile("./wrfout_d01_2011-07-25_00:00:00.nc","r")
  type = "png"
  res=True
  res@MainTitle="trackplot"

  times=wrf_user_list_times(a)
  ntimes=dimsizes(times)

  wks=gsn_open_wks(type,"trackplot")
  gsn_define_colormap(wks,"rainbow+white+gray")

  plot=gsn_csm_map(wks,res)

alon    =wrf_user_getvar(a,"lon",0)
lon    =new(dimsizes(alon(0,:)),float)
lat    =new(dimsizes(alon(:,0)),float)

  mres            = True
  mres@gsMarkerIndex    = create_hurricane_symbol(wks)
  mres@gsMarkerSizeF    = 8
  mres@gsMarkerColor    = "green"
do  it=0,ntimes-1
  slp = wrf_user_getvar(a,"slp",it)    ; slp
    wrf_smooth_2d( slp, 6 )            ; smooth slp
  locmin = local_min(slp,False,0.)
    x=lon(locmin@xi)                ; get lat/lon points of minima
    y=lat(locmin@yi)
    gsn_polymarker(wks,plot,x,y,mres)

end do
end
理论上说应该很简单的啊……可是为啥还是这么复杂……报错……
fatal:Subscript out of range, error in subscript #0
fatal:An error occurred reading xf
fatal:Execute: Error occurred at or near line 4219 in file $NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl

fatal:Execute: Error occurred at or near line 4460 in file $NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl

fatal:Execute: Error occurred at or near line 51 in file tracksplot1.ncl
我就是把最小值的程序和绘制台风图标的程序放一块儿了，因为想要画出路径来，所以用了一个循环，然后就开始报错了……我调试的时候发现是第一个时次就错了……还不是循环的问题……不懂了……http://www.ncl.ucar.edu/Applications/Scripts/minmax_1.ncl里面的例子应该差不多呀……
求教各位大虾了……

mofangbao

我还以为是grads的。。。ncl还没用过，昨天倒是用grads实现了一个类似的

浪迹天涯

看了一下，Subscript out of range, error in subscript ，Error occurred at or near line 51 in file tracksplot1.ncl
第51行有问题，下标超出了~要是兰溪在就好了，他做过这方面的，可惜他已经走了~

陪着你的路人

浪迹天涯 发表于 2011-12-7 22:43

                               
登录/注册后可看大图

看了一下，Subscript out of range, error in subscript ，Error occurred at or near line 51 in file tra ...

我也知道，但是就是改不对……

浪迹天涯

陪着你的路人 发表于 2011-12-7 22:57

                               
登录/注册后可看大图

我也知道，但是就是改不对……

你可以用print出来看看的嘛，唉，我也走了~

陪着你的路人

浪迹天涯 发表于 2011-12-7 23:07

                               
登录/注册后可看大图

你可以用print出来看看的嘛，唉，我也走了~

嗯嗯~这下改对了，是维数没有弄对……

随缘

陪着你的路人 发表于 2011-12-8 11:55

                               
登录/注册后可看大图

嗯嗯~这下改对了，是维数没有弄对……

楼主能否把改对的脚本放上来，我们学习学习啊，谢谢啊！具体是怎么改的，我也是遇到维数不对，不知道怎么改

陪着你的路人

随缘 发表于 2012-9-19 21:50

                               
登录/注册后可看大图

楼主能否把改对的脚本放上来，我们学习学习啊，谢谢啊！具体是怎么改的，我也是遇到维数不对，不知道怎 ...

不好意思一直没有上，这是我自己写的最后程序，跟最开始的算法有改动，你看看合不合用吧
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"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/shea_util.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/wrf/WRFUserARW.ncl"

;********************Symbol Settings****************************
undef("create_hurricane_symbol")
function create_hurricane_symbol(wks)
begin
  mstring = "p"
  fontnum = 37
  xoffset = 0.0
  yoffset = 0.0
  ratio   = 1.0
  size    = 1.0
  angle   = 0.0
  new_index = NhlNewMarker(wks, mstring, fontnum, xoffset, yoffset, \
                           ratio, size, angle)
  return(new_index)
end
;*******************Symbol Settings End*****************************

begin

  a = addfile("./wrfout_d01_2011-07-25_00:00:00.nc","r")
  b = addfile("./geo_em.d01.nc","r")   ; Open a file

  type = "png"

  pltres=True
  mpres=True

  times=wrf_user_list_times(a)
  ntimes=dimsizes(times)

  wks=gsn_open_wks(type,"tracksplot")
  gsn_define_colormap(wks,"rainbow+white+gray")


;################Finding Typhoon Eyes Settings##############
lon2d        = a->XLONG(0,:,:)
lat2d        = a->XLAT(0,:,:)
fx        = new(ntimes,integer)
fy        = new(ntimes,integer)
temp        = new(ntimes,integer)

  mres                                = True
  mres@gsMarkerIndex                = create_hurricane_symbol(wks)
  mres@gsMarkerSizeF                = 15
  mres@gsMarkerColor                = "green"

  pres                         = True               ; polyline resources
  pres@gsLineThicknessF        = 2.0                ; line thickness
  pres@gsLineColor                = "Blue"

  mpres@gsnDraw                 = False
  mpres@gsnFrame                = False
  mpres@gsnMaximize                = True
  mpres@tiMainString         = "Typhoon track"

  mpres@pmLabelBarWidthF         = 0.6
  mpres@pmLabelBarHeightF        = 0.1
  mpres@mpMinLatF                = min(lat2d)
  mpres@mpMaxLatF                = max(lat2d)
  mpres@mpMinLonF                = min(lon2d)
  mpres@mpMaxLonF                = max(lon2d)

;################Finding Typhoon Eyes Settings End##############

;#################Finding Typhoon Eyes##########################

do it=0,ntimes-1
  slp2d = wrf_user_getvar(a,"slp",it)    ; slp
    wrf_smooth_2d( slp2d, 6 )            ; smooth slp
  dims                = dimsizes(slp2d)
  slp1d        = ndtooned(slp2d)
  if (it .eq. 0) then
        slp        = slp2d(14:(dims(0)-50),100:(dims(1)-5))
  else
        if (fx(it-1) .gt. (dims(0)-16)) then
                  slp        = slp2d((fx(it-1)-15):fx(it-1),(fy(it-1)-15):(fy(it-1)+15))
        else
                if (fx(it-1) .lt. 14) then
                          slp        = slp2d(fx(it-1):(fx(it-1)+15),(fy(it-1)-15):(fy(it-1)+15))
                else
                        if (fy(it-1) .gt. (dims(1)-16)) then
                                  slp        = slp2d((fx(it-1)-15):(fx(it-1)+15),(fy(it-1)-15):fy(it-1))
                        else
                                if (fy(it-1) .lt. 14) then
                                          slp        = slp2d((fx(it-1)-15):(fx(it-1)+15),fy(it-1):(fy(it-1)+15))
                                else
                                          slp        = slp2d((fx(it-1)-15):(fx(it-1)+15),(fy(it-1)-15):(fy(it-1)+15))
                                end if
                        end if
                end if
        end if
  end if
  minij        = ind_resolve(ind(slp1d.eq.min(slp)),dims)
  fx(it)        = minij(0,0)
  fy(it)        = minij(0,1)
  temp(it)        = get1Dindex(slp1d,slp2d(fx(it),fy(it)))
  delete(slp)
end do
;#################Finding Typhoon Eyes End#######################

;###########################Drawing##############################
  ter        = b->HGT_M(0,:,:)                      ; Read the variable to memory
  ter@lon2d        = b->XLONG_M(0,:,:)
  ter@lat2d        = b->XLAT_M(0,:,:)
  terrain = gsn_csm_contour_map(wks,ter,mpres)
  dot  = new(ntimes,graphic)    ; Make sure each gsn_add_polyxxx call
  line = new(ntimes-1,graphic)    ; is assigned to a unique variable.

do it=0,ntimes-2
  xx         = (/lon2d(fx(it),fy(it)),lon2d(fx(it+1),fy(it+1))/)
  yy         = (/lat2d(fx(it),fy(it)),lat2d(fx(it+1),fy(it+1))/)
  line(it)        = gsn_add_polyline(wks,terrain,xx,yy,pres)
end do
lon1d        = ndtooned(lon2d)
lat1d        = ndtooned(lat2d)
do it=0,ntimes-1,10
  dot(it)=gsn_add_polymarker(wks,terrain,lon1d(temp(it)),lat1d(temp(it)),mres)
end do

;##################draw the track##########
do  it=1,ntimes-1
;  gsn_polyline(wks,terrain,(/lon2d(fx(it-1),fy(it-1)),lon2d(fx(it),fy(it))/),(/lat2d(fx(it-1),fy(it-1)),lat2d(fx(it),fy(it))/),pres)
end do
do  it=0,ntimes-1,10
;  gsn_polymarker(wks,terrain,lon2d(fx(it),fy(it)),lat2d(fx(it),fy(it)),mres)
end do
draw(terrain)
frame(wks)
end

若当块

新手，过来学习一下。

sunlight

陪着你的路人 发表于 2012-10-23 11:00

                               
登录/注册后可看大图

不好意思一直没有上，这是我自己写的最后程序，跟最开始的算法有改动，你看看合不合用吧
load "$NCARG_R ...

很赞！
弱弱的想问一下
if (it .eq. 0) then
        slp        = slp2d(14:(dims(0)-50),100:(dims(1)-5)) 。。。。。。。
这些判断语句的目的是什么？