Re: Re:

Moti,

I ran your script with my sample wrfout data, and it worked well.
The warning message you've got seems to be related to the data values in
your wrfout file.
Could you check your wrf run first? I would check u,v winds in wrfout.nc
using "ncview".
Again, the script itself has no problem.

So-Young

Motilal Mittal wrote:
> Hi Mary,
>
> When I execute the command
> ncl < wrf_real.ncl
>
> I get the error
> in compute_uvmet 143 143 30 144 144
> computing velocities
> warning:ContourPlotInitialize: scalar field is constant; ContourPlot
> not possible:[errno=1102]
> warning:ContourPlotSetValues: Data values out of range of levels set
> by EXPLICITLEVELS mode
>
> I am attaching the script wrf_real.ncl for your reveiw.
>
> I am still in learning mode and want to understand how to resolve this
> error.
>
> Thank you.
>
> Moti
> Sr. Scientist, Ohio Supercomputer Center
> Graduate Faculty, Graduate Program in Environmental Science (OSU)
> Columbus, OH 43212-1163
> www.osc.edu/pcrm/emissions
> voice:(614)292-2552 (O); (614)789-9384 (H) fax:(614)292-7168
> e-mail: moti_at_osc.edu, mittal.1_at_osu.edu
> ------------------------------------------------------------------------
>
> external WRF "./wrf_user_fortran_util_0.so"
> ; Script to produce standard plots for a WRF real-data run,
> ; with the ARW coordinate dynamics option.
>
> load "WRFOptions.ncl" ; set basic plot options here
> ;load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
> load "gsn_code.ncl"
> ;load "/local/LINUX_PC/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
> load "WRFPlot.ncl"
> load "WRFUserARW.ncl"
> load "SkewTFunc.ncl"
>
> begin
> ;
> ; The WRF ARW input file.
> ; This needs to have a ".nc" appended, so just do it.
>
> a = addfile("../EXAMPLES/wrfout.nc","r")
>
> ; We generate plots, but what kind do we prefer?
>
> type = "x11"
> ; type = "pdf"
> ; type = "ps"
> ; type = "ncgm"
>
> wks = gsn_open_wks(type,"wrf_real")
>
> ; Debug information.
>
> debug = False
> ; debug = True
> icount = 0
>
> ; Basic Plot Information
>
> res_at_MainTitle = "REAL-TIME WRF"
>
> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
>
> ; What times and how many time steps are in the data set?
>
> FirstTime = True
> times = wrf_user_list_times(a) ; get times in the file
> ntimes = dimsizes(times) ; number of times in the file
>
> ; The specific pressure levels that we want the data
> ; interpolated to.
>
> pressure_levels = (/ 850., 700., 500., 300./) ; pressure levels to plot
> nlevels = dimsizes(pressure_levels) ; number of pressure levels
>
> ; This is the big loop over all of the time periods to process.
>
> do it = 0,ntimes-1,2
> ; do it = 0,2,2
>
> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
>
> ;Save some Time Information and create the Map Background
>
> time = it
> if (FirstTime) then
> times_sav = times(it)
> end if
> res_at_TimeLabel = times(it)
>
> map = wrf_map(wks,a,res)
>
> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
> ; First get the variables we will need
>
> slvl = wrf_user_getvar(a,"slvl",time) ; psl
> wrf_user_filter2d(slvl, 3) ; smooth psl
> qv = wrf_user_getvar(a,"QVAPOR",time) ; Qv
> tc = wrf_user_getvar(a,"tc",time) ; T in C
> td = wrf_user_getvar(a,"td",time) ; Td in C
> if ( it .eq. 0 ) then
> tc2 = tc(0,:,:) ; Use lowest T at time zero
> else
> tc2 = wrf_user_getvar(a,"T2",time) ; T2 in Kelvin
> tc2 = tc2-273.16 ; T2 in C
> end if
> if (it .eq. 0) then
> td2 = td(0,:,:) ; Use lowest Td at time zero
> else
> td2 = wrf_user_getvar(a,"td2",time) ; Td2 in C
> end if
> u = wrf_user_getvar(a,"ua",time) ; u averaged to mass points
> v = wrf_user_getvar(a,"va",time) ; v averaged to mass points
> uvm = wrf_user_getvar(a,"umeta",time) ; umet,vmet averaged to mass points
> if (it .eq. 0) then
> u10 = u(0,:,:) ; Use lowest level at time 0
> v10 = v(0,:,:)
> else
> u10 = wrf_user_getvar(a,"U10",time) ; u at 10 m, mass point
> v10 = wrf_user_getvar(a,"V10",time) ; v at 10 m, mass point
> end if
> w = wrf_user_getvar(a,"wa",time) ; vertical velocity, averaged to half levels
> w = 100.*w ; vert vel in cm/s
> p = wrf_user_getvar(a, "p",time) ; pressure is our vertical coordinate
> z = wrf_user_getvar(a, "Z",time) ; grid point height
> rh = wrf_user_getvar(a,"rh",time) ; relative humidity
> ter = wrf_user_getvar(a,"HGT",time) ; terrain height
>
> ; Get non-convective, convective and total precipitation
> ; Calculate tendency values
> rain_exp = wrf_user_getvar(a,"RAINNC",time)
> rain_con = wrf_user_getvar(a,"RAINC",time)
> rain_tot = rain_exp + rain_con
> if( FirstTime ) then
> rain_exp_save = rain_exp
> rain_con_save = rain_con
> rain_tot_save = rain_tot
> end if
> rain_exp_tend = rain_exp - rain_exp_save
> rain_con_tend = rain_con - rain_con_save
> rain_tot_tend = rain_tot - rain_tot_save
> ; Bookkeeping, just to allow the tendency at the next time step
> rain_exp_save = rain_exp
> rain_con_save = rain_con
> rain_tot_save = rain_tot
>
> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
>
> ; ************************************************************
>
> ; 1) Surface Temp + SLP + wind barbs
> ; Surface Td + wind barbs
>
>
> ; Plotting options for T
>
> tc_plane = 1.8*tc2+32. ; Turn temperature into Fahrenheit - Typical of a surface map
> opts_tc = res
> opts_tc_at_FieldTitle = "Surface Temperature"
> opts_tc_at_UnitLabel = "F"
> opts_tc_at_ContourParameters = (/ -20., 90., 5./)
> opts_tc_at_cnFillOn = True
> opts_tc_at_gsnSpreadColorEnd = -3 ; End third from the last color in color map
> contour_tc = wrf_contour(a,wks,tc_plane,opts_tc)
> print_opts("opts_tc", opts_tc, debug)
>
>
> ; Plotting options for Td
>
> td_plane = 1.8*td2+32. ; Turn temperature into Fahrenheit - Typical of a surface map
> opts_td = res
> opts_td_at_FieldTitle = "Surface Dew Point Temp"
> opts_td_at_UnitLabel = "F"
> opts_td_at_ContourParameters = (/ -20., 90., 5./)
> opts_td_at_cnFillOn = True
> opts_td_at_cnLinesOn = True
> opts_td_at_cnLineLabelsOn = True
> opts_td_at_cnLineLabelBackgroundColor = -1
> opts_td_at_gsnSpreadColorEnd = -3 ; End third from the last color in color map
> contour_td = wrf_contour(a,wks,td_plane,opts_td)
> print_opts("opts_td", opts_td, debug)
>
> ; Plotting options for SLP
>
> opts_psl = res
> ;opts_psl_at_FieldTitle = "Sea Level Pressure"
> ;opts_psl_at_UnitLabel = "mb"
> opts_psl_at_ContourParameters = (/ 900., 1100., 4. /)
> opts_psl_at_cnLineColor = "NavyBlue"
> opts_psl_at_cnHighLabelsOn = True
> opts_psl_at_cnLowLabelsOn = True
> opts_psl_at_cnLineLabelBackgroundColor = -1
> opts_psl_at_gsnContourLineThicknessesScale = 2.0
> contour_psl = wrf_contour(a,wks,slvl,opts_psl)
> print_opts("opts_psl", opts_psl, debug)
>
> ; Plotting options for Wind Vectors
>
> u_plane = u10*1.94386 ; Turn wind into knots - Typical of a surface map
> v_plane = v10*1.94386
> opts_vct = res
> opts_vct_at_FieldTitle = "Winds"
> opts_vct_at_UnitLabel = "kts"
> opts_vct_at_NumVectors = 47
> vector = wrf_vector(a,wks,u_plane, v_plane,opts_vct)
> print_opts("opts_vct", opts_vct, debug)
>
> ; MAKE PLOTS
>
> wrf_map_overlay(wks,map,(/contour_tc,contour_psl,vector/),False)
> print_header(icount,debug)
>
> wrf_map_overlay(wks,map,(/contour_td,vector/),False)
> print_header(icount,debug)
>
> ; Delete options and fields, so we don't have carry over
> delete(opts_tc)
> delete(opts_td)
> delete(opts_vct)
> delete(opts_psl)
> delete(tc_plane)
> delete(td_plane)
> delete(u_plane)
> delete(v_plane)
>
> ; ************************************************************
>
> ; 2) 2d precip plots
>
> ; Plotting options for Sea Level Pressure
>
> opts_psl = res
> opts_psl_at_ContourParameters = (/ 900., 1100., 2. /)
> opts_psl_at_cnLineColor = "Blue"
> opts_psl_at_cnInfoLabelOn = False
> opts_psl_at_cnHighLabelsOn = True
> opts_psl_at_cnLowLabelsOn = True
> opts_psl_at_gsnContourLineThicknessesScale = 1.0
> contour_psl = wrf_contour(a,wks,slvl,opts_psl)
> print_opts("opts_psl",opts_psl, debug)
>
>
> ; Plotting options for Precipitation
>
> opts_r = res
> opts_r_at_UnitLabel = "mm"
> opts_r_at_cnLevelSelectionMode = "ExplicitLevels"
> opts_r_at_cnLevels = (/ .1, .2, .4, .8, 1.6, 3.2, 6.4, \
> 12.8, 25.6, 51.2, 102.4/)
> opts_r_at_cnFillColors = (/"White","White","DarkOliveGreen1", \
> "DarkOliveGreen3","Chartreuse", \
> "Chartreuse3","Green","ForestGreen", \
> "Yellow","Orange","Red","Violet"/)
> opts_r_at_cnInfoLabelOn = False
> opts_r_at_cnConstFLabelOn = False
> print_opts("opts_r", opts_r, debug)
>
>
> ; MAKE PLOTS
> opts_r_at_cnFillOn = True
>
> ; Total Precipitation (color fill)
> opts_r_at_FieldTitle = "Total Precipitation"
> contour_tot = wrf_contour(a,wks, rain_tot, opts_r)
>
> ; Precipitation Tendency (color fill)
> opts_r_at_FieldTitle = "Precipitation Tendency"
> opts_r_at_SubFieldTitle = "from " + times_sav + " to " + times(it)
> contour_tend = wrf_contour(a,wks, rain_tot_tend,opts_r)
>
> ; Non-Convective Precipitation Tendency (color fill)
> opts_r_at_FieldTitle = "Explicit Precipitation Tendency"
> opts_r_at_SubFieldTitle = "from " + times_sav + " to " + times(it)
> contour_res = wrf_contour(a,wks,rain_exp_tend,opts_r)
>
> ; Convective Precipitation Tendency (contour lines)
> opts_r_at_FieldTitle = "Param Precipitation Tendency"
> opts_r_at_SubFieldTitle = "from " + times_sav + " to " + times(it)
> opts_r_at_cnFillOn = False
> opts_r_at_cnLineColor = "Red4"
> contour_prm = wrf_contour(a,wks,rain_con_tend,opts_r)
>
>
> ; Total Precipitation Tendency + SLP
>
> wrf_map_overlay(wks,map,(/contour_tend,contour_psl/),False)
> print_header(icount,debug)
>
> ; Total Precipitation
>
> wrf_map_overlay(wks,map,contour_tot,False)
> print_header(icount,debug)
>
> ; Non-Convective and Convective Precipiation Tendencies
>
> wrf_map_overlay(wks,map,(/contour_res,contour_prm/),False)
> print_header(icount,debug)
>
>
> ; Delete options, so we don't have carry over
> delete(opts_r)
> delete(opts_psl)
>
> ; ************************************************************
>
> ; 3) 3d plots: height, rh, vert vel, temp
>
> ; Loop over the selected pressure levels to plot
>
> do level = 0,nlevels-1
>
> pressure = pressure_levels(level)
>
> tc_plane = wrf_user_intrp3d(tc,p,ter,"h",pressure,0.)
> z_plane = wrf_user_intrp3d( z,p,ter,"h",pressure,0.)
> rh_plane = wrf_user_intrp3d(rh,p,ter,"h",pressure,0.)
> u_plane = wrf_user_intrp3d( u,p,ter,"h",pressure,0.)
> v_plane = wrf_user_intrp3d( v,p,ter,"h",pressure,0.)
> spd = (u_plane*u_plane + v_plane*v_plane)^(0.5) ; m/sec
> u_plane = u_plane*1.94386 ; kts
> v_plane = v_plane*1.94386 ; kts
>
>
> ; Plotting options for T
>
> opts_tc = res
> opts_tc_at_cnLineColor = "Red"
> opts_tc_at_ContourParameters = (/ 5.0 /)
> opts_tc_at_cnInfoLabelOrthogonalPosF = 0.06
> opts_tc_at_gsnContourLineThicknessesScale = 2.0
>
> ; Plotting options for Geopotential Heigh
>
> opts_z = res
> opts_z_at_cnLineColor = "NavyBlue"
> opts_z_at_gsnContourLineThicknessesScale = 3.0
>
> ; Plotting options for RH
>
> opts_rh = res
> opts_rh_at_ContourParameters = (/ 10., 90., 10./)
> opts_rh_at_cnFillOn = True
> opts_rh_at_cnFillColors = (/"White","White","White", \
> "White","Chartreuse","Green",\
> "Green3","Green4", \
> "ForestGreen","PaleGreen4"/)
>
>
>
> ; Plotting options for Wind Speed
>
> opts_spd = res
> opts_spd_at_FieldTitle = "Wind Speed"
> opts_spd_at_UnitLabel = "m/s"
> opts_spd_at_cnLineColor = "MediumSeaGreen"
> opts_spd_at_ContourParameters = (/ 10. /)
> opts_spd_at_cnInfoLabelOrthogonalPosF = 0.06
> opts_spd_at_gsnContourLineThicknessesScale = 3.0
>
> ; Plotting options for Wind Vectors
>
> opts_vct = res
> opts_vct_at_FieldTitle = "Winds"
> opts_vct_at_UnitLabel = "kts"
> opts_vct_at_NumVectors = 47
>
> ; MAKE PLOTS
>
> if ( pressure .eq. 850 ) then
> ; temp, rh, height, barbs
>
> opts_z_at_ContourParameters = (/ 20.0 /)
>
> contour_rh = wrf_contour(a,wks,rh_plane,opts_rh)
> print_opts("opts_rh", opts_rh, debug)
> contour_tc = wrf_contour(a,wks,tc_plane,opts_tc)
> print_opts("opts_tc", opts_tc, debug)
> contour_height = wrf_contour(a,wks, z_plane,opts_z)
> print_opts("opts_z", opts_z, debug)
> vector = wrf_vector(a,wks,u_plane, v_plane,opts_vct)
> print_opts("opts_vct", opts_vct, debug)
>
> wrf_map_overlay(wks,map,(/contour_rh,contour_tc,contour_height, \
> vector/),False)
> print_header(icount,debug)
>
> end if
>
> if ( pressure .eq. 700 ) then
> ; temp, height, barbs
>
> opts_z_at_ContourParameters = (/ 30.0 /)
>
> contour_tc = wrf_contour(a,wks,tc_plane,opts_tc)
> print_opts("opts_tc", opts_tc, debug)
> contour_height = wrf_contour(a,wks, z_plane,opts_z)
> print_opts("opts_z", opts_z, debug)
> vector = wrf_vector(a,wks,u_plane, v_plane,opts_vct)
> print_opts("opts_vct", opts_vct, debug)
>
> wrf_map_overlay(wks,map,(/contour_tc,contour_height, \
> vector/),False)
> print_header(icount,debug)
>
> end if
>
>
> if ( pressure .eq. 500 ) then
> ; temp, height, barbs
>
> opts_z_at_ContourParameters = (/ 60.0 /)
>
> contour_tc = wrf_contour(a,wks,tc_plane,opts_tc)
> print_opts("opts_tc", opts_tc, debug)
> contour_height = wrf_contour(a,wks, z_plane,opts_z)
> print_opts("opts_z", opts_z, debug)
> vector = wrf_vector(a,wks,u_plane, v_plane,opts_vct)
> print_opts("opts_vct", opts_vct, debug)
>
> wrf_map_overlay(wks,map,(/contour_tc,contour_height, \
> vector/),False)
> print_header(icount,debug)
>
> end if
>
> if ( pressure .eq. 300 ) then
> ; windspeed, height, barbs
>
> opts_z_at_ContourParameters = (/ 60.0 /)
>
> contour_spd = wrf_contour(a,wks,spd,opts_spd)
> print_opts("opts_spd", opts_spd, debug)
> contour_height = wrf_contour(a,wks, z_plane,opts_z)
> print_opts("opts_z", opts_z, debug)
> vector = wrf_vector(a,wks,u_plane, v_plane,opts_vct)
> print_opts("opts_vct", opts_vct, debug)
>
> wrf_map_overlay(wks,map,(/contour_spd,contour_height, \
> vector/),False)
> print_header(icount,debug)
>
> end if
>
> end do ; END OF 3D LOOP
>
> ; Delete options and fields, so we don't have carry over
> delete(opts_tc)
> delete(opts_z)
> delete(opts_rh)
> delete(opts_spd)
> delete(opts_vct)
> delete(z_plane)
> delete(tc_plane)
> delete(rh_plane)
> delete(u_plane)
> delete(v_plane)
> delete(spd)
>
> ; ************************************************************
>
> ; 4) xsection plots
>
> do ip = 1, 2
> ; |
> ; |
> ; Two different xsections, angle=90 is |, angle=0 is ------
> ; |
> ; |
> if(ip .eq. 1) then
> angle = 90.
> Orientation = "Cross-Section Orientation : West-East"
> else
> angle = 0.
> Orientation = "Cross-Section Orientation : South-North"
> end if
>
> ; Options for XY Plots
> opts_xy = res
> opts_xy_at_tiXAxisString = "Number of Grid Points"
> opts_xy_at_tiYAxisString = "Height (km)"
> opts_xy_at_PlotType = "XY"
> opts_xy_at_AspectRatio = 0.75
> opts_xy_at_PlotOrientation = Orientation
> opts_xy_at_cnMissingValPerimOn = True
> opts_xy_at_cnMissingValFillColor = 0
> opts_xy_at_cnMissingValFillPattern = 11
> opts_xy_at_tmYLMode = "Explicit"
> opts_xy_at_tmYLValues = ispan(0,100,10)
> opts_xy_at_tmYLLabels = ispan(0,20,2)
> opts_xy_at_tiXAxisFontHeightF = 0.020
> opts_xy_at_tiYAxisFontHeightF = 0.020
> opts_xy_at_tmXBMajorLengthF = 0.02
> opts_xy_at_tmYLMajorLengthF = 0.02
>
>
> ; Build planes of data for temp and rh
>
> dimsrh = dimsizes(rh)
> plane = new(2,float)
> plane = (/ dimsrh(2)/2, dimsrh(1)/2 /)
>
> rh_plane = wrf_user_intrp3d(rh,z,ter,"v",plane,angle)
> tc_plane = wrf_user_intrp3d(tc,z,ter,"v",plane,angle)
>
> ; Plotting options for RH
>
> opts_rh = opts_xy
> opts_rh_at_ContourParameters = (/ 10., 90., 10. /)
> opts_rh_at_cnFillOn = True
> opts_rh_at_cnFillColors = (/"White","White","White", \
> "White","Chartreuse","Green", \
> "Green3","Green4", \
> "ForestGreen","PaleGreen4"/)
>
> ; Plotting options for Temperature
>
> opts_tc = opts_xy
> opts_tc_at_ContourParameters = (/ 5. /)
>
> ; Get the contour info for the rh and temp
>
> contour_tc = wrf_contour(a,wks,tc_plane,opts_tc)
> print_opts("opts_tc", opts_tc, debug)
> contour_rh = wrf_contour(a,wks,rh_plane,opts_rh)
> print_opts("opts_rh", opts_rh, debug)
>
> ; Overlay the fields
>
> wrf_overlay(wks,(/contour_rh,contour_tc/),False)
> print_header(icount,debug)
>
> ; Delete options and fields, so we don't have carry over
> delete(opts_tc)
> delete(opts_rh)
> delete(tc_plane)
> delete(rh_plane)
> end do ; end of xsection plot rotation
>
> ; ************************************************************
>
> ; 5) skew-T plots
>
> u = uvm(0,:,:,:)*1.94386
> v = uvm(1,:,:,:)*1.94386
>
> do ip = 12, 13
>
> ; Define a few skew-T plotting options
>
> skewtOpts = True
> dataOpts = True
> dataOpts_at_Parcel = 0
> dataOpts_at_WspdWdir = False ; wind speed and dir [else: u,v]
> dataOpts_at_PlotWindH = False ; plot wind barbs at h lvls [pibal; special]
> dataOpts_at_HspdHdir = True ; wind speed and dir [else: u,v]
>
> ; Choose the lat/lon of the skew-T plots
>
> ip_lats = (/ 39.77, 37.77, 41.12, 41.30, 39.07, 35.22, \
> 37.23, 35.22, 31.95, 32.50, 43.56, 38.94, \
> 40.0 /)
> ip_lons = (/-104.87, -99.97, -100.67, -95.90, -95.62, -97.45, \
> -93.38, -101.72, -102.20, -97.30, -116.24, -77.44, \
> -82.90 /)
> ip_locs = (/"Denver (DEN)", "Dodge City (DDC)", "North Platte (LBF)", \
> "Omaha (OMA)", "Topeka (TOP)", "Norman (OUN)", \
> "Springfield (SGF)", "Amarillo (AMA)", "Midland (MAF)", \
> "Fort Worth (FWD)", "Boise (BOI)", "Washington DC (IAD)", \
> "Columbus (CMH)"/)
>
> ip_locs = ip_locs + ", Valid at " + times(it)
> locr = wrf_user_find_ij_lat_long(a, ip_lats(ip-1), ip_lons(ip-1))
> loc = floattointeger(locr)
>
> ; Get the skew-T background
>
> skewtOpts_at_tiMainString = ip_locs(ip-1)
> skewt_bkgd = skewT_BackGround (wks, skewtOpts)
>
> ; Draw the skew-T plot
>
> draw (skewt_bkgd)
> skewT_data = skewT_PlotData(wks, skewt_bkgd, p(:,loc(0), loc(1)), \
> tc(:,loc(0), loc(1)), \
> td(:,loc(0), loc(1)), \
> z(:,loc(0), loc(1)), \
> -u(:,loc(0), loc(1)), \
> -v(:,loc(0), loc(1)), \
> dataOpts)
>
> ; Close the frame
>
> frame(wks)
> print_header(icount,debug)
>
> ; Delete the temporary arrays
>
> delete(skewtOpts)
> delete(dataOpts)
> delete(skewT_data)
> delete(skewt_bkgd)
>
> end do
>
> ; ************************************************************
>
> times_sav = times(it)
> FirstTime = False
> end do ; end of the time loop
>
> end
>
> ------------------------------------------------------------------------
>
> _______________________________________________
> ncl-talk mailing list
> ncl-talk_at_ucar.edu
> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>
Received on Tue Nov 07 2006 - 13:55:18 MST