CAPE calculation

Hello

   I am trying to calculate CAPE and other variables from my WRF data fiel using the code below. However I get the following message and I don't know why.

 capecalc3d: Outside of lookup table bounds. prs,thte= 950.832343750000 NaN

Is it because I have missing values? If so, how can I fix it? Or what other explanation there is for the error?

Thank you,
Cecille

; Example script to produce standard plots for a WRF quarter_ss run

load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
;load "$NCARG_ROOT/lib/ncarg/nclscripts/wrf/WRFUserARW.ncl"
load "./WRFUserARW.ncl"

begin

dir = "Morrison"
;scheme = "MOR"
grid = "1km"
;city = "Jasper"
perturbation = "delta_4"
period = "past"

path="/project/hypocenter/i/cloudmp/Cecille/WRF_output/"+dir+"/"+grid+"/"+perturbation+"/"

a = addfile(path+"wrfout_d01_0001-01-01_00:00:00_"+period+".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_cape")
  gsn_define_colormap(wks,"WhViBlGrYeOrReWh")

; Set some basic resources
  res = True
  res@MainTitle = "REAL-TIME WRF"

  pltres = True
  mpres = True
  mpres@mpGeophysicalLineColor = "Black"
  mpres@mpNationalLineColor = "Black"
  mpres@mpUSStateLineColor = "Black"
  mpres@mpGridLineColor = "Black"
  mpres@mpLimbLineColor = "Black"
  mpres@mpPerimLineColor = "Black"
  mpres@mpGridSpacingF = 10

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
; Which times and how many time steps are in the data set?

  ;times = wrf_user_getvar(a,"times",-1) ; get all times in the file
  times = wrf_user_list_times(a)
  ntimes = dimsizes(times) ; number of times in the file

 it = ntimes-1 ; only interested in last time
  res@TimeLabel = times(it) ; set Valid time to use on plot

  cape2d = wrf_user_getvar(a,"cape_2d",it)
  cape3d = wrf_user_getvar(a,"cape_3d",it)

  print ("cape3d")

; extract the data from the multi-dimensional fields
  cape = cape3d(0,:,:,:)
  cin = cape3d(1,:,:,:)
  mcape = cape2d(0,:,:)
  mcin = cape2d(1,:,:)
  lcl = cape2d(2,:,:)
  lfc = cape2d(3,:,:)

print ("Variables Calculated")

  opts = res
  opts@cnFillOn = True
  opts@cnLineLabelsOn = False

  opts@FieldTitle = "mcape"
  opts@ContourParameters = (/ 500., 3000., 250./)
  contour = wrf_contour(a,wks,mcape,opts)
  plot = wrf_map_overlays(a,wks,(/contour/),pltres,mpres)

  opts@FieldTitle = "mcin"
  opts@ContourParameters = (/ 0., 125., 25./)
  contour = wrf_contour(a,wks,mcin,opts)
  plot = wrf_map_overlays(a,wks,(/contour/),pltres,mpres)

  opts@FieldTitle = "lfc"
  opts@ContourParameters = (/ 200., 3800., 400./)
  contour = wrf_contour(a,wks,lfc,opts)
  plot = wrf_map_overlays(a,wks,(/contour/),pltres,mpres)

  opts@FieldTitle = "lcl"
  opts@ContourParameters = (/ 200., 3800., 400./)
  contour = wrf_contour(a,wks,lcl,opts)
  plot = wrf_map_overlays(a,wks,(/contour/),pltres,mpres)

  opts@FieldTitle = "cape"
  opts@ContourParameters = (/ 250., 3250., 250./)
  do il = 0,9 ; only interested in the first 10 levels
    ilev = il + 1
    opts@PlotLevelID = "level " + ilev
    contour = wrf_contour(a,wks,cape(il,:,:),opts)
    plot = wrf_map_overlays(a,wks,(/contour/),pltres,mpres)
  end do

  opts@FieldTitle = "cin"
  opts@ContourParameters = (/ 0., 300., 25./)
  do il = 0,9 ; only interested in the first 10 levels
    ilev = il + 1
    opts@PlotLevelID = "level " + ilev
    contour = wrf_contour(a,wks,cin(il,:,:),opts)
    plot = wrf_map_overlays(a,wks,(/contour/),pltres,mpres)
  end do

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

end

> ********************************************************
> Cecille M. Villanueva Birriel
> Ph.D. Candidate
> Cloud Microphysics Research Group
> Purdue University
> Department of Earth, Atmospheric, & Planetary Sciences
> 550 Stadium Mall Drive, West Lafayette, IN 47907-2051
> email: cvillanu@purdue.edu
> ********************************************************

_______________________________________________
ncl-talk mailing list
List instructions, subscriber options, unsubscribe:
http://mailman.ucar.edu/mailman/listinfo/ncl-talk
Received on Tue Sep 4 15:01:16 2012