[ncl-talk] Help With Cross Section Script

Kerwyn Texeira ktish86 at gmail.com
Fri Aug 26 11:52:20 MDT 2016


Hi ncl-talk,

I'm having issues trying to plot ageostrophic winds cross section across a
terrain.  I'm using a 12km wrf data.  When I use 36km no terrain with the
same script, I did not get any errors and I got very nice plots.

The 12km wrf data, to me, is very good because I was able to get other
vertical cross sections across the terrain using Matlab.

Whatever help anyone can give me would be greatly apreciated

Thanks in advance!

Error Message:

fatal:VectorPlotDraw: VVECTR - VECTOR NDC LENGTH TOO GREAT

fatal:VectorPlotDraw: error drawing vectors

fatal:VectorPlotDraw: draw error

fatal:PlotManagerDraw: error in plot draw

fatal:_NhlPlotManagerDraw: Draw error

Varaible Outputs:

Variable: ua
Type: float
Total Size: 2690400 bytes
            672600 values
Number of Dimensions: 3
Dimensions and sizes:    [bottom_top | 50] x [south_north | 118] x
[west_east | 114]
Coordinates:
Number Of Attributes: 6
  FieldType :    104
  MemoryOrder :    XYZ
  description :    x-wind component
  units :    m s-1
  stagger :
  coordinates :    XLONG XLAT

Variable: va
Type: float
Total Size: 2690400 bytes
            672600 values
Number of Dimensions: 3
Dimensions and sizes:    [bottom_top | 50] x [south_north | 118] x
[west_east | 114]
Coordinates:
Number Of Attributes: 6
  FieldType :    104
  MemoryOrder :    XYZ
  description :    y-wind component
  units :    m s-1
  stagger :
  coordinates :    XLONG XLAT

Variable: z
Type: float
Total Size: 2690400 bytes
            672600 values
Number of Dimensions: 3
Dimensions and sizes:    [bottom_top | 50] x [south_north | 118] x
[west_east | 114]
Coordinates:
Number Of Attributes: 6
  FieldType :    104
  MemoryOrder :    XYZ
  description :    Height
  units :    m
  stagger :
  coordinates :    XLONG XLAT

Variable: p
Type: float
Total Size: 2690400 bytes
            672600 values
Number of Dimensions: 3
Dimensions and sizes:    [bottom_top | 50] x [south_north | 118] x
[west_east | 114]
Coordinates:
Number Of Attributes: 6
  coordinates :    XLONG XLAT
  stagger :
  units :    hPa
  description :    Pressure
  MemoryOrder :    XYZ
  FieldType :    104

Variable: ur
Type: float
Total Size: 2690400 bytes
            672600 values
Number of Dimensions: 3
Dimensions and sizes:    [bottom_top | 50] x [lat | 118] x [lon | 114]
Coordinates:
            lat: [31.90399932861328..44.84910202026367]
            lon: [-128.5808258056641..-111.4191741943359]
Number Of Attributes: 8
  _FillValue :    9.96921e+36
  FieldType :    104
  MemoryOrder :    XYZ
  description :    x-wind component
  units :    m s-1
  stagger :
  remap :    remapped via ESMF_regrid_with_weights: Bilinear remapping
  missing_value :    9.96921e+36

Variable: vr
Type: float
Total Size: 2690400 bytes
            672600 values
Number of Dimensions: 3
Dimensions and sizes:    [bottom_top | 50] x [lat | 118] x [lon | 114]
Coordinates:
            lat: [31.90399932861328..44.84910202026367]
            lon: [-128.5808258056641..-111.4191741943359]
Number Of Attributes: 8
  _FillValue :    9.96921e+36
  FieldType :    104
  MemoryOrder :    XYZ
  description :    y-wind component
  units :    m s-1
  stagger :
  remap :    remapped via ESMF_regrid_with_weights: Bilinear remapping
  missing_value :    9.96921e+36

Variable: zr
Type: float
Total Size: 2690400 bytes
            672600 values
Number of Dimensions: 3
Dimensions and sizes:    [bottom_top | 50] x [lat | 118] x [lon | 114]
Coordinates:
            lat: [31.90399932861328..44.84910202026367]
            lon: [-128.5808258056641..-111.4191741943359]
Number Of Attributes: 8
  _FillValue :    9.96921e+36
  FieldType :    104
  MemoryOrder :    XYZ
  description :    Height
  units :    m
  stagger :
  remap :    remapped via ESMF_regrid_with_weights: Bilinear remapping
  missing_value :    9.96921e+36

Variable: press
Type: float
Total Size: 2690400 bytes
            672600 values
Number of Dimensions: 3
Dimensions and sizes:    [bottom_top | 50] x [lat | 118] x [lon | 114]
Coordinates:
            lat: [31.90399932861328..44.84910202026367]
            lon: [-128.5808258056641..-111.4191741943359]
Number Of Attributes: 8
  _FillValue :    9.96921e+36
  stagger :
  units :    hPa
  description :    Pressure
  MemoryOrder :    XYZ
  FieldType :    104
  remap :    remapped via ESMF_regrid_with_weights: Bilinear remapping
  missing_value :    9.96921e+36

Variable: uv
Type: float
Total Size: 5380800 bytes
            1345200 values
Number of Dimensions: 4
Dimensions and sizes:    [2] x [50] x [118] x [114]
Coordinates:
Number Of Attributes: 1
  _FillValue :    9.96921e+36
(0)    uv: min=-6821.45   max=7810.25

Variable: u
Type: float
Total Size: 2690400 bytes
            672600 values
Number of Dimensions: 3
Dimensions and sizes:    [50] x [118] x [114]
Coordinates:
Number Of Attributes: 1
  _FillValue :    9.96921e+36
(0)    u: min=-5291.31   max=5693.54

Variable: v
Type: float
Total Size: 2690400 bytes
            672600 values
Number of Dimensions: 3
Dimensions and sizes:    [50] x [118] x [114]
Coordinates:
Number Of Attributes: 1
  _FillValue :    9.96921e+36
(0)    v:min=-6821.45    max=7810.25

Variable: uageo
Type: float
Total Size: 2690400 bytes
            672600 values
Number of Dimensions: 3
Dimensions and sizes:    [50] x [118] x [114]
Coordinates:
Number Of Attributes: 1
  _FillValue :    9.96921e+36
(0)    uageo:min=-5693.41  max=5288.43

Variable: vageo
Type: float
Total Size: 2690400 bytes
            672600 values
Number of Dimensions: 3
Dimensions and sizes:    [50] x [118] x [114]
Coordinates:
Number Of Attributes: 1
  _FillValue :    9.96921e+36
(0)    vageo:min=-7812.72  max=6823.08


My Script:

;======================================================================
; ESMF_regrid.ncl

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"
load "$NCARG_ROOT/lib/ncarg/nclscripts/esmf/ESMF_regridding.ncl"

begin

;---Specify method to be used
    InterpMethod= "patch"

;---Input file
    srcDirName  = "./"
    srcFileName = "wrfout_d01_2014-01-11_13_00_00.nc"
    srcFilePath =  srcDirName + srcFileName

   pltres = True

;---Wgt File: WRF to Rectilinear
    wgtDirName  = "./"
    wgtFileName = "WRF_to_Rect.WgtFile_"+InterpMethod+".nc"
    wgtFilePath =  wgtDirName + wgtFileName


;---Retrieve either one level, or all levels. Use '-1' for all.
    sfile = addfile(srcFilePath,"r")

    FirstTime = True
    FirstTimeMap = True
    times = wrf_user_getvar(sfile, "times", -1)
    ntimes = dimsizes(times)

    mdims = getfilevardimsizes(sfile, "P")
    nd = dimsizes(mdims)


    ua   = wrf_user_getvar(sfile,"ua",0)    ; On mass grid
    va   = wrf_user_getvar(sfile,"va",0)
    z    = wrf_user_getvar(sfile,"z", 0)
    p    = wrf_user_getvar(sfile,"pressure", 0)


    printVarSummary(ua)                      ;
(Time,bottom_top,south_north,west_east)
    printVarSummary(va)                      ;
(Time,bottom_top,south_north,west_east)
    printVarSummary(z)
    printVarSummary(p)



;---Regrid the wind components to a rectilinear grid
    ur   = ESMF_regrid_with_weights(ua,wgtFilePath,False)
    vr   = ESMF_regrid_with_weights(va,wgtFilePath,False)
    zr   = ESMF_regrid_with_weights(z,wgtFilePath, False)
    press = ESMF_regrid_with_weights(p,wgtFilePath,False)


    printVarSummary(ur)
    printVarSummary(vr)
    printVarSummary(zr)
    printVarSummary(press)

;---Compute the geostrophic winds  on the rectilinear grid
    lat2d = sfile->XLAT(0,:,:)               ; (south_north,west_east)
    lon2d = sfile->XLONG(0,:,:)

;---Generate the same rectilinear grid used to generate the weight file
    dims  = dimsizes(lat2d)
    nlat  = dims(0)
    nlon  = dims(1)

    lat = fspan(min(lat2d), max(lat2d) ,nlat)
    lon = fspan(min(lon2d), max(lon2d) ,nlon)

;---Calculate the geostrophic winds on a rectilinear grid
    uv = z2geouv(zr, lat, lon, 1)
    printVarSummary(uv)
    print("uv: min="+min(uv)+"   max="+max(uv))

;----Calculate the geostrophic u wind component on a rectilinear grid
    u = uv(0,:,:,:)
    printVarSummary(u)
    print("u: min="+min(u)+"   max="+max(u))

;----Calculate the geostrophic v wind component on a rectinlinear grid
   v = uv(1,:,:,:)
   printVarSummary(v)
   print("v:min="+min(v)+"    max="+max(v))

;---Calculate the uageo component
 uageo = ur - u
   printVarSummary(uageo)
   print("uageo:min="+min(uageo)+"  max="+max(uageo))

;---Calculate the vageo component
  vageo = vr - v
  printVarSummary(vageo)
  print("vageo:min="+min(vageo)+"  max="+max(vageo))
;----------------------------------------------------------------------
  opt = False
  plane = new(4, float)
  plane = (/  63,37  , 63,72   /)
  uageo_plane = wrf_user_intrp3d(uageo, p, "v", plane, 0.0, opt)
  vageo_plane = wrf_user_intrp3d(vageo, p, "v", plane, 0.0, opt)
  x_plane  = wrf_user_intrp2d(lat2d, plane, 0.0, opt)

;--Let's create nice labels - only have to do this one
   if ( FirstTime ) then
          zmax = 200.     ;  Place top at model top or close to zmax
          zz = wrf_user_intrp3d(p,p,"v",plane,0.,opt)
          z_ind = ind(.not.ismissing(zz(:,0)))
          zmin = zz(z_ind(0),0)
          delete(z_ind)
          nice_levs = floor((zmin-zmax)/50)*50
          zmax = zmin - nice_levs
          dims = dimsizes(zz)
          zmax_pos = dims(0)-1
          do imax = 1,dims(0)-1
            if ( .not.ismissing(zz(imax,0)) .and. zz(imax,0) .ge. zmax )
then
              zmax_pos = imax
            end if
          end do
          zspan = zmax_pos
          zmax = zz(zmax_pos,0)
          nz = floattoint((zmin-zmax)/50+1)
          FirstTime = False
        end if


;-- x-axis labeles
   dimsX = dimsizes(x_plane)
      xmin  = x_plane(0)
      xmax  = x_plane(dimsX(0)-1)
      xspan = dimsX(0)-1
      nx    = floattoint( (xmax-xmin)/0.5 + 1)


; create plots: Note some defaults chaanged in NCL v6.1.0
;************************************************
    wks_type = "png"
    wks_type at wkWidth = 2500
    wks_type at wkHeight = 2500
   wks = gsn_open_wks("png","ESMF")

  gsn_define_colormap(wks,"matlab_jet")     ; this is the default v6.1.0
onward

   res = True
   res at tiXAxisString = "Latitude"
   res at tiYAxisString = "Pressure (mb)"
   res at tmXTOn        = False
   res at tmYROn        = False
   res at tmXBMode      = "Explicit"
   res at tmXBValues    = fspan(0, xspan, nx)
   res at tmXBLabels    = sprintf("%.1f",fspan(xmin,xmax,nx))
   res at tmXBTickSpacingF = 1.0
   res at tmXBLabelFontHeightF  = 0.015
   res at tmYLMode      = "Explicit"
   res at tmYLValues    = fspan(0, zspan, nz)
   res at tmYLLabels    = sprintf("%.1f",fspan(zmin, zmax, nz))
   res at tiXAxisFontHeightF = 0.015
   res at tiYAxisFontHeightF = 0.015
   res at tmXBMajorLengthF   = 0.02
   res at tmYLMajorLengthF   = 0.02
   res at tmYLLabelFontHeightF = 0.015


    res1                      = res
    res1 at gsnDraw              = False           ; don't draw
    res1 at gsnFrame             = False           ; don't advance frame
   ; res1 at gsnAddCyclic         = False           ; regional data
    res1 at vcLineArrowThicknessF = 3.0
    res1 at vcRefMagnitudeF       = 10
    res1 at vcRefLengthF          = 0.018
   ; res1 at gsnLeftString         = ""
   ; res1 at gsnRightString         = ""
    res1 at vcGlyphStyle          = "LineArrow"
    res1 at vcMinDistanceF        = 0.05
    res1 at vcRefAnnoOn           = True
    res1 at gsnMaximize           = True
    res1 at tiMainString          = "Ageostrophic Winds Cross Section on Jan
11 at 13UTC"
    res1 at tiMainFontHeightF     = 0.018



ageo = gsn_csm_vector(wks, uageo_plane(0:zmax_pos,:),
vageo_plane(0:zmax_pos,:), res1)


;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;

draw(ageo)
frame(wks)

end
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20160826/30caf6b7/attachment.html 


More information about the ncl-talk mailing list