NCL画图问题

sss123

尽头的尽头

只给错误，不给代码，恐怕很难帮到你

sss123

尽头的尽头 发表于 2019-9-26 08:42
只给错误，不给代码，恐怕很难帮到你

load "$NCARG_ROOT/lib/ncarg/nclscripts/esmf/ESMF_regridding.ncl"
begin
;-------open WRF_out files-----------------------
dir = "/state/partition1/home/zhugf/WRF_out_6_1.10/"
f = addfile( dir + "wrfout_d01_2018-06-08_18^%00^%00.nc", "r")
f1 = addfile("real_obs_T_UV.nc","r")
;-----read varibles from wrf files---------
  it        = 0     ; first time step
  precip       = wrf_user_getvar(f,"RAINNC",it)    ; accumulated total precipitation
  precip@lat2d = wrf_user_getvar(f,"XLAT",it)   ; latitude/longitude
  precip@lon2d = wrf_user_getvar(f,"XLONG",it)  ; required for plotting
  
  T       = wrf_user_getvar(f,"T2",it)    ; accumulated total precipitation
  T@lat2d = wrf_user_getvar(f,"XLAT",it)   ; latitude/longitude
  T@lon2d = wrf_user_getvar(f,"XLONG",it)  ; required for plotting
  printVarSummary(precip)
  print(min(precip))
  print(max(precip))
  ;-------------read varibles from nc files----------------------
  time = f1->time
  ob_T2 = short2flt(f1->t2m)
  ;ob_p  = f1->tp
  printVarSummary(ob_T2)

   utc_date= cd_calendar(time, 0)
   year   = tointeger(utc_date(:,0))   ; Convert to integer for
   month  = tointeger(utc_date(:,1))    ; use sprinti
   day    = tointeger(utc_date(:,2))    ; use sprinti
   hour   = tointeger(utc_date(:,3))    ; use sprinti
   ;print(day)
   ;print(hour)
   
   iYYYY = ind(day.eq.8.and.hour.eq.18)
   ob_T2_6_08 = ob_T2(iYYYY,:,:)
   copy_VarMeta(ob_T2(0,:,:), ob_T2_6_08)
   printVarSummary(ob_T2_6_08)
;;;get the source lat/lon grid-----------------------------------
  lat2d = wrf_user_getvar(f,"XLAT",it)   ; latitude/longitude
  lon2d = wrf_user_getvar(f,"XLONG",it)  ; required for plotting
  dims = dimsizes(lat2d)
  nlat = dims(0)
  nlon = dims(1)
  
  minlon  = min(T@lon2d)-1
  minlat  = min(T@lat2d)-1
  maxlon  = max(T@lon2d)+1
  maxlat  = max(T@lat2d)+1

  ;--------get the destination lat/lon grid------------------------
  lat  = f1->latitude
  lon  = f1->longitude

   Opt = True
   Opt@SrcFileName     = "WRF_SCRIP.nc"      ; Name of source and
   Opt@DstFileName     = "Rectilinear.nc"    ; destination files
   Opt@Overwrite            = True
   Opt@DstGridType          = "rectilinear"
   Opt@DstGridLat           = lat
   Opt@DstGridLon           = lon

   Opt@InterpMethod         = "bilinear"
   Opt@SrcRegional          = True
   Opt@DstRegional          = True

    Opt@SkipSrcGrid   = True
    Opt@SkipDstGrid   = True
    Opt@SkipWgtGen    = True

   T_regrid = ESMF_regrid(T,Opt)     ; Do the regridding for TMP
;---Reset 0 values to missing values.
    T_regrid@_FillValue = default_fillvalue(typeof(T_regrid))

    T_regrid            = where(T_regrid.eq.0.0,T_regrid@_FillValue,\
                                  T_regrid)


     printVarSummary(T_regrid)

  ;--------------------------difference-----------------------------------------
   T_diff = T_regrid - ob_T2_6_08
   copy_VarMeta(ob_T2_6_08, T_diff)
   printVarSummary(T_diff)
  ;------------------------creat plot------")
  wks = gsn_open_wks("pdf", "wrf_ERA_T_diff")
  res                      = True

  res@gsnMaximize          = True
  res@gsnDraw              = False
  res@gsnFrame             = False
  res@gsnAddCyclic        = True

  res@mpMinLatF     = 15.;min(hgt@lat2d)
  res@mpMaxLatF     = 40.;max(hgt@lat2d)
  res@mpMinLonF     = 90.;min(hgt@lon2d)
  res@mpMaxLonF     = 120;max(hgt@lon2d)
  
  res@mpOutlineOn           = True
  res@mpFillOn              = False
  res@mpDataBaseVersion     = "MediumRes"       ; better and more map outlines
  ;res@mpDataSetName         = "Earth..4"
  res@mpFillBoundarySets     = "National"
  ;res@mpOutlineBoundarySets = "AllBoundaries"
  ;res@mpOutlineSpecifiers   =(/"China:states"/)
  res@mpGeophysicalLineColor ="black"
  res@gsnLeftString         ="2018_6_8_T2_diff"
  res@cnFillOn             = True
  res@cnFillPalette = "BlueDarkRed18";NCV_blue_red";WhViBlGrYeOrRe"  
  res@cnLinesOn            = False
  ;res@cnLevelSpacingF      = 5.           ; NCL chose 4.0
  res@lbOrientation        = "Vertical"
  res@tiMainString         = "wrfout_d01_2018-06-08_18:00:00"

  ;---Change contour levels to better match the color map being used
   res@cnLevelSelectionMode =  "ManualLevels"   
   res@cnMinLevelValF       = -20;-1.3;-1.5
   res@cnMaxLevelValF       = 20;1.3; 1.5
   res@cnLevelSpacingF      =  1
   res@cnSmoothingOn = True
   res@cnSmoothingTensionF = 0.002
  plot = gsn_csm_contour_map(wks,ob_T2_6_08,res)
  draw(plot)
  frame(wks)
end

以上是代码，我是算wrf模拟的数据和观测数据的差值，观测数据，插值之后的数据都可以正常画图，但是同样的代码作差之后的数据就会出现这个报错，不知道什么原因，还请指教！

尽头的尽头

sss123 发表于 2019-9-26 10:02
load "$NCARG_ROOT/lib/ncarg/nclscripts/esmf/ESMF_regridding.ncl"
begin
;-------open WRF_out fil ...

看了一下，一个是你的数据应该不是全球数据，gsnAddCyclic = True应该设为false吧，还有差值的话注意等值线间隔，应该很小才对

sss123

sss123 发表于 2019-9-26 10:02
load "$NCARG_ROOT/lib/ncarg/nclscripts/esmf/ESMF_regridding.ncl"
begin
;-------open WRF_out fil ...

好的，已经解决了，非常感谢