Re: Plotting HGT in 2D Cross Section

Lauren,

What kind of difficulty are you having? Do you get an error message, or a bad plot, or ?

Can you provide the wrfout file you're using?

--Mary

On Feb 26, 2014, at 12:46 PM, Lauren Wheeler <lauren.bronwyn.wheeler@gmail.com> wrote:

> Hi all,
> I'm trying to plot the output from an idealized WRF simulation, we have modified the 2d_hill example so that the simulation is actually 3D and we can input idealized topography.
>
> What I'm trying to do is use the example script provided by WRF/NCL to plot a cross section through the center of the ridge but I'm having difficulty adding the terrain height in as a contour. Below is my script which works before I try and add in the topo, any suggestions?
>
> Thanks,
> Lauren Wheeler
>
>
> ;Plots from 3D Hill Simulation
>
> load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
> load "$NCARG_ROOT/lib/ncarg/nclscripts/wrf/WRFUserARW.ncl"
> ;load "./WRFUserARW.ncl"
>
> begin
> ;
> ; The WRF ARW input file.
> ; This needs to have a ".nc" appended, so just do it.
>
> a = addfile("wrfout_d01_0001-01-03_00:00:00"+".nc","r")
>
> ; We generate plots, but what kind do we prefer?
>
> ; type = "x11"
> type = "pdf"
> ; type = "ps"
> ; type = "ncgm"
>
> wks = gsn_open_wks(type,"plt_Hill3D_06”)
>
>
> ; Set some Basic Plot options
> res = True
> res@MainTitle = "WRF Hill 3D"
> res@InitTime = False
> res@Footer = False
>
> pltres = True
>
>
> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
>
> ; What times and how many time steps are in the data set?
>
> times = wrf_user_list_times(a) ; get times in the file
> ntimes = dimsizes(times) ; number of times in the file
>
> ; The specific plane we want to plot data on
>
> plane = (/ 100., 0./) ; (x,y) point for vertical plane
> angle = 90.0
> pii = 3.14159
> aspect_ratio = .7
>
> ; This is the big loop over all of the time periods to process.
>
> do it = 2,ntimes-1,2
>
> time = it
> res@TimeLabel = times(it)
> res@AspectRatio = aspect_ratio
>
>
> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
> ; First get the variables we will need
>
> u_in = wrf_user_getvar(a,"ua",time) ; ua is u averaged to mass points
> v_in = wrf_user_getvar(a,"va",time) ; va is v averaged to mass points
> w_in = wrf_user_getvar(a,"wa",time) ; vertical velocity
> th_in = wrf_user_getvar(a,"th",time) ; get temperature (C)
> z_in = wrf_user_getvar(a, "z",time) ; grid point height
> ter = wrf_user_getvar(a,"HGT",time) ; need terrain height sometimes
>
> u = u_in(0:34,:,:)
> v = v_in(0:34,:,:)
> w = w_in(0:34,:,:)
> th = th_in(0:34,:,:)
> z = z_in(0:34,:,:)
>
> res@tmYLMode = "Explicit"
> res@tmYLValues = fspan(0.,100.,5)
> res@tmYLLabels = sprintf("%.1f",fspan(0.,30.,5))
>
> u_plane = wrf_user_intrp3d( u,z,"v",plane,angle,False)
> v_plane = wrf_user_intrp3d( v,z,"v",plane,angle,False)
> w_plane = wrf_user_intrp3d( w,z,"v",plane,angle,False)
> th_plane = wrf_user_intrp3d(th,z,"v",plane,angle,False)
> ter_plane = wrf_user_intrp3d(u,z,”ter”,plane,angle,False)
>
> vel_normal = u_plane*cos(2.*pii*angle/360.) - v_plane*sin(2.*pii*angle/360.)
> vel_tangent = u_plane*sin(2.*pii*angle/360.) + v_plane*cos(2.*pii*angle/360.)
> vel_tangent = vel_tangent - 10.
>
> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
>
> ; Horizontal Wind Perturbation Contour Fill
> opts_th = res
> opts_th@FieldTitle = u_in@description
> opts_th@cnFillOn = True
> opts_th@gsnSpreadColorEnd = -10
> opts_th@PlotOrientation = u_plane@Orientation
> contour_th = wrf_contour(a,wks,u_plane,opts_th)
>
>
> ; Theta
> ; opts_th = res
> ; opts_th@FieldTitle = th_in@description
> ; opts_th@cnFillOn = True
> ; opts_th@gsnSpreadColorEnd = -10
> ; opts_th@PlotOrientation = th_plane@Orientation
> ; contour_th = wrf_contour(a,wks,th_plane,opts_th)
>
> ; Vertical Velocity
> opts_w = res
> opts_w@FieldTitle = w_in@description
> contour_w = wrf_contour(a,wks, w_plane,opts_w)
>
> ; Vel Tangent
> opts_vt = res
> opts_vt@FieldTitle = "Perturbation u"
> opts_vt@UnitLabel = "m/s"
> contour_vt = wrf_contour(a,wks,vel_tangent,opts_vt)
>
> ;Topography
> opts_ter = res
> contour_ter = wrf_contour(a,weks,ter_plane,opts_ter)
>
>
> ; plot = wrf_overlays(a,wks,(/contour_th, contour_vt/),pltres)
> plot = wrf_overlays(a,wks,(/contour_th, contour_w/),pltres)
> plot = wrf_overlays(a,wks,(/contour_th, contour_w, contour_ter/),pltres)
> ; ************************************************************
>
> end do ; end of the time loop
>
> end
>
> _______________________________________________
> ncl-talk mailing list
> List instructions, subscriber options, unsubscribe:
> http://mailman.ucar.edu/mailman/listinfo/ncl-talk

_______________________________________________
ncl-talk mailing list
List instructions, subscriber options, unsubscribe:
http://mailman.ucar.edu/mailman/listinfo/ncl-talk
Received on Thu Feb 27 10:35:06 2014