Plotting HGT in 2D Cross Section

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
Received on Wed Feb 26 12:46:37 2014