[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