请选择 进入手机版 | 继续访问电脑版
爱气象,爱气象家园! 

气象家园

 找回密码
 立即注册

QQ登录

只需一步,快速开始

新浪微博登陆

只需一步, 快速开始

搜索
查看: 4689|回复: 2

ncl做组图时出错

[复制链接]

新浪微博达人勋

发表于 2018-4-21 08:25:26 | 显示全部楼层 |阅读模式

登录后查看更多精彩内容~

您需要 登录 才可以下载或查看,没有帐号?立即注册 新浪微博登陆

x
想画一下不同时刻的风场,位温场以及dbz的组合图。出现两个问题,单个图里面没有叠加上位温,组合图里面只有dbz。。。。
从昨天到今天都没搞好,求大神指教一下做组合图到底是不是有啥规则

QQ图片20180421082430.jpg
QQ图片20180421082420.jpg
密码修改失败请联系微信:mofangbao

新浪微博达人勋

 楼主| 发表于 2018-4-21 08:46:14 | 显示全部楼层
问题解决了,我自己再调一下。好像是原始数据的问题。这个程序里面涉及了Wrf数据的提取,插值还有以台风为中心点的转换,有需要的自提
密码修改失败请联系微信:mofangbao
回复 支持 1 反对 0

使用道具 举报

新浪微博达人勋

 楼主| 发表于 2018-4-21 08:27:01 | 显示全部楼层
脚本如下
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"

begin
  f   = addfile("chanchu_2km.nc","r")
   
   times = wrf_user_getvar(f,"times",-1)  ; get all times in the file
   ntimes = dimsizes(times)  
   ti=(/180,276,372,392/);14日03时,15日03时 , ;16日03时,16日18时

     km=41
     lq=243
     lp=243
     nx=201
     ny=201
     nr=100                        
     nt=360                     
     vmis=1.0e30
     dbz1=new((/4,41,243,243/), float)
     tc1 =new((/4,41,243,243/), float)
     u1 =new((/4,41,243,243/), float)
     v1 =new((/4,41,243,243/), float)
    ; mslp=new((/2/), "float")
     owa1=new((/4,ny,nx/),"float") ;dbz
     owa2=new((/4,ny,nx/),"float")  ;tc
     owa3=new((/4,ny,nx/),"float")  ;v
     owa4=new((/4,ny,nx/),"float")  ;v

     do i =0,3
      
      dbz           =wrf_user_getvar(f, (/"dbz","1","1"/),ti(i))
      slp           =wrf_user_getvar(f, "slp", ti(i))
      u             =wrf_user_getvar(f, "ua",ti(i))
      v             =wrf_user_getvar(f, "va",ti(i))
      tc            =wrf_user_getvar(f, "tc", ti(i))
      
      H=(/0.1,0.2,0.3,0.4,0.5,0.75,1.0,1.5,2.0,2.5,3.0,3.5,4.0,4.5,5.0,5.5,6.0,6.5,7.0,7.5,8.0,8.5,9.0,9.5,10.0,10.5,11.0,11.5,12.0,12.5,13.0,13.5,14.0,14.5,15.0,15.5,16.0,16.5,17.0,17.5,18.0/)
      dims=dimsizes(slp)
      opts          =True
      opts@time     =i
      dbz1(i,:,:,:) =wrf_user_vert_interp(f,dbz, "ght_agl", H, opts)
      tc1(i,:,:,:) =wrf_user_vert_interp(f,tc, "ght_agl", H, opts)
      u1(i,:,:,:) =wrf_user_vert_interp(f,u, "ght_agl", H, opts)
      v1(i,:,:,:) =wrf_user_vert_interp(f,v, "ght_agl", H, opts)
      
      slp1D         =ndtooned(slp)
      mslp          =min(slp1D)
      ii            =ind_resolve(ind(slp1D.eq.mslp),dims)
      jj            =ndtooned(ii)
      tyc           =jj(0)
      txc           =jj(1)

      do j=0,ny-1
            do k=0,nx-1
           owa1(i,j,k)=dbz1(i,6,tyc-nr+j,txc-nr+k);dbz
           owa2(i,j,k)=tc1(i,6,tyc-nr+j,txc-nr+k);tc
           owa3(i,j,k)=u1(i,6,tyc-nr+j,txc-nr+k);u
           owa4(i,j,k)=v1(i,6,tyc-nr+j,txc-nr+k);v

           end do
        end do
     end do
   
    ;绘图设置
      wks  =gsn_open_wks("png", "dbzzutu")
      gsn_define_colormap(wks,"WhViBlGrYeOrReWh")
      res                          =True
      pltres                       =True
      pltres@PanelPlot             =True
      res@NoHeaderFooter           =True
      res@pmLabelBarOrthogonalPosF =-0.05
      res@lbTitleOn                =False

      plots      =new(4, graphic)
     
      res@tiMainOn        =False
      res@tiXAxisOn       =False
      res@tiYAxisOn       =False
      res@vpHeightF       =0.5
      res@vpWidthF        =0.5
      res@cnInfoLabelOn   =False
      res@tmXBMode        ="Explicit"
      res@tmXBValues      =(/0,25,50,75,100,125,150,175,200/)
      res@tmXBLabels      =(/"-200","-150","-100","-50","0","50","100","150","200"/)
      res@tmYLMode        ="Explicit"
      res@tmYLValues      =(/0,25,50,75,100,125,150,175,200/)
      res@tmYLLabels      =(/"-200","-150","-100","-50","0","50","100","150","200"/)

     do i = 0, 3
         ;设置dbz
      resd                      =res
      resd@cnFillOn             = True
      resd@cnLinesOn            = False
      resd@cnLevelSelectionMode ="ManualLevels"
      ;resd@cnMinLevelValF       =5.
      ;resd@cnLevelSpacingF      =5.
     ; resd@cnMaxLevelValF       =75.
      resd@ContourParameters = (/ 5.,65.,5./)
      resd@lbOrientation        = "Horizontal"
      counter_dbz           = wrf_contour(f,wks,owa1(i,:,:),resd)    ; plot only lowest level
      delete(resd)

      ;设置风场
      resd                       =res
      resd@vcGlyphStyle          ="FillArrow"
      resd@vcFillArrowFillColor  ="black"
      resd@vcFillArrowWidthF     =0.1
      resd@vcLineArrowThicknessF =4
      resd@vcMinMagnitudeF       =2.0
      resd@vcMinDistanceF        =0.04
      ;参考箭头
      resd@vcRefAnnoOn           =True
      resd@vcFillArrowsOn       = True
      resd@vcRefMagnitudeF       =40
      resd@vcFillArrowEdgeColor = "black"
      resd@vcRefLengthF          =0.035
      resd@vcRefAnnoParallelPosF =0.95
      resd@vcRefAnnoFontHeightF  =0.015
      resd@vcRefAnnoSide         ="Right"
      resd@vcRefAnnoString1On    =False
      resd@vcRefAnnoString2      ="40 m/s"
      vector                =wrf_vector(f, wks,owa3(i,:,:),owa4(i,:,:),resd)
      delete(resd)
       ;设置tc
       resd                      =res
       resd@cnFillOn             = False
       resd@cnLinesOn            = True
       resd@cnLevelSelectionMode ="ManualLevels"
       resd@cnLineThicknesses    =4
       resd@cnLineDashPattern    =0
       resd@cnLineColors          ="red"
       ;resd@cnMinLevelValF      =5.
       resd@cnLevelSpacingF      =5.
       counter_tc             = wrf_contour(f,wks,owa2(i,:,:),resd)   
       delete(resd)
   
      plots(i) =wrf_overlays(f, wks,(/counter_dbz,vector,counter_tc/),pltres)
     ; plots(1) =wrf_overlays(f, wks,(/counter_dbz(1),counter_tc(1),vector(1)/),pltres)
    ;  plots(2) =wrf_overlays(f, wks,(/counter_dbz(2),counter_tc(2),vector(2)/),pltres)
     ; plots(3) =wrf_overlays(f, wks,(/counter_dbz(3),counter_tc(3),vector(3)/),pltres)
      delete(counter_dbz)
      delete(counter_tc)
      delete(vector)
      end do

     pntlre                 =True
     pntlre@gsnPanelYWhiteSpacePercent  =3
     pntlre@gsnPanelScalePlotIndex      =1;指定某幅图来确定整个图的比例
     pntlre@gsnPanelFigureStrings       =(/"(a) 14/03","(b) 15/03","(c) 16/03","(d) 16/18"/)
     pntlre@amJust                      ="TopLeft"
     pntlre@gsnPanelFigureStringsFontHeightF =0.015
     pntlre@gsnPanelFigureStringsPerimOn    =False

     gsn_panel(wks,(/plots/),(/2,2/),pntlre)

  end

                        
   
     
   

     
   
   


   
密码修改失败请联系微信:mofangbao
您需要登录后才可以回帖 登录 | 立即注册 新浪微博登陆

本版积分规则

Copyright ©2011-2014 bbs.06climate.com All Rights Reserved.  Powered by Discuz! (京ICP-10201084)

本站信息均由会员发表,不代表气象家园立场,禁止在本站发表与国家法律相抵触言论

快速回复 返回顶部 返回列表