台风涡度垂直剖面图

香猫猫

辐合在低层，都到不了500hPa,不知道哪里出问题啦，求指点

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/wrf/WRFUserARW.ncl"
begin
        a = addfile("wrfout_d01_2017-08-21_18:00:00","r")

        type = "png"
        wks = gsn_open_wks(type,"prevor")
        gsn_define_colormap(wks,"WhiteBlueGreenYellowRed")
        times  = wrf_user_getvar(a,"times",-1) ; get times in the file
    ntimes = dimsizes(times)         ; number of times in the file
       
        do it = 14,16,2
        print("Working on time: " + times(it) )
        print(it)
        u  = wrf_user_getvar(a,"ua",it)
        v  = wrf_user_getvar(a,"va",it)
        p  = wrf_user_getvar(a, "pressure",it)
       
        latm2d=a->XLAT(0,:,:)
lonm2d=a->XLONG(0,:,:)
fstart=min(latm2d)
fend=max(latm2d)
dims=dimsizes(latm2d)
lat = fspan(fstart, fend, dims(0))
lat@units = "degrees_north"
fstart=min(lonm2d)
fend=max(lonm2d)
lon = fspan(fstart, fend, dims(1))
lon@units = "degrees_east"
printVarSummary(lon)
        ;printVarSummary(p)
       
        vor = uv2dv_cfd(u, v, lat, lon, 1)
        vor@long_name = "Vorticity"
        vor@units = "10 -5 s -1"
        ;printVarSummary(vor)
       
        plane=(/111.7,21.7/)
        ;plane1=(/120,21.7/)
        ;rh_plane2 = wrf_user_intrp3d(rh,p,"v",plane,90.0,False)
       
         vor_plane = wrf_user_intrp3d(vor,p,"v",plane,90.0,False)
         ;printVarSummary(vor_plane)
     p_plane = wrf_user_intrp3d(p, p,"v",plane, 90.,False)
    lon_plane = wrf_user_intrp2d(lonm2d,plane,90.,False)
        ;print(lon_plane)
  ;X_plane = wrf_user_intrp2d(lon,plane,0,False)
;xmin=min(X_plane)
;xmax=max(X_plane)

    vor_plane!0="level"
    vor_plane!1="lon"
    vor_plane&level = p_plane(:, 0)
    vor_plane&lon = lon_plane
        ;vor_plane@_FillValue=-999
        ;printVarSummary(vor_plane)
        ;print(vor_plane&level)
    ; MAKE PLOTS
       
        scale = 1.e05
         vor_plane=vor_plane*scale
         
        opts = True           
      opts@cnLinesOn = False              
      opts@cnFillOn = True  
      opts@cnLineLabelsOn = False
      opts@cnLineLabelFontHeightF = 0.01
      opts@cnLineLabelInterval = 1      ; Dengzhixian Biaozhu JianGe
      opts@cnLineLabelPerimOn  = False  ; Dengzhixian Biaozhu Biankuang de Xianyin Kaiguan
      opts@gsnSpreadColorEnd    = -3    ; End third from the last color in color map
      opts@cnLineLabelPlacementMode = "Computed"
      opts@cnLineLabelDensityF = 0.0
      opts@lbTitleOn            = False
      opts@pmLabelBarOrthogonalPosF = -0.05
      opts@cnLevelSelectionMode = "ExplicitLevels"
      opts@cnLevels             = (/ -20,-16,-12, -8, -4, 0, 4,  8,  12,  16,  20,  24, 28, 32, 36/)
      opts@cnFillColors         = (/  80 ,64, 48, 32, 16, 10, 140,160,180, 190, 200,220,230,240,255/)
          ;opts@SubFieldTitle =  times(it)
          ;opts@SubFieldTitle =  times(it)
          opts@tiMainString = times(it)
          
          opts@tmXBMode = "Manual"
          opts@tmXBTickStartF = 124
          opts@tmXBTickEndF = 112
          opts@tmXBTickSpacingF = 2
          ;opts@tmXBValues =
          ;opts@tmXBMinorOn=True
     plot  = gsn_csm_pres_hgt(wks,vor_plane(1:71,{112:124}),opts )
         
         x=117.4
         y=1000
         polyres               = True          ; poly marker mods desired
     polyres@gsMarkerIndex = 16            ; choose circle as polymarker
     polyres@gsMarkerSizeF = 40.0           ; select size to avoid streaking
     polyres@gsMarkerColor = (/"GREEN"/)   ; choose color
     gsn_polymarker(wks,plot,x,y,polyres)
        end do
end