[ncl-talk] Overlaying Contour on WRF Domain

[Sean Egan]

Hi,

I have an HDF file that contains volcanic ash column densities in a
900x3712 array. Most of the elements are fill values. In addition, there
are two arrays of latitudes and longitudes, both 900x3712.

I am trying to use gsn_contour or gsn_csm_contour to plot the ash
concentrations over a WRF model domain. Currently, I'm unable to get NCL to
reference the ash array with the latitude and longitude arrays. I have
posted my script below and can upload the HDF file if needed.

Any tips as to why I cannot get NCL to reference the lat and lon arrays
correctly would be greatly appreciated. I believe it may have something to
do with the sfXArray and sfYArray resources....

Thanks,
Sean


Script Start

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/wrf/WRFUserARW.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"

begin

a = addfile("file1.hdf","r")
b = addfile("wrfout_d01_2010-04-16_00:00:00","r")
wks=gsn_open_wks("x11","testing")

latS = a->pixel_latitude
latN = latS * 0.002746665852839747
lonS = a->pixel_longitude
lonN = lonS * 0.005493331705679495
ash  = a->volcld_ret_14_15_16_ash_mass_loading

lat = b->XLAT(i0,:,:) ;Not a nested domain
lon = b->XLONG(0,:,:)
dims = dimsizes(lat)
nlat = dims(0)
nlon = dims(1)

gsn_define_colormap(wks,"WhBlGrYeRe")
res                         = True
res at gsnMaximize             = False
res at gsnDraw                 = False        ; Turn off for now.
res at gsnFrame                = False        ; Will draw later
res at tiMainOn                = False

mpres                              = res
mpres at mpProjection                 = "LambertConformal"
mpres at mpLimitMode                  = "Corners"
mpres at mpLeftCornerLatF             = lat(0,0)
mpres at mpLeftCornerLonF             = lon(0,0)
mpres at mpRightCornerLatF            = lat(nlat-1,nlon-1)
mpres at mpRightCornerLonF            = lon(nlat-1,nlon-1)
mpres at mpLambertParallel1F          = b at TRUELAT1
mpres at mpLambertParallel2F          = b at TRUELAT2
mpres at mpLambertMeridianF           = b at STAND_LON
mpres at mpDataBaseVersion            = "MediumRes"
mpres at mpOutlineBoundarySets        = "GeophysicalAndUSStates"
mpres at mpFillOn                     = False
mpres at mpOutlineOn                  = True
mpres at mpCountyLineThicknessF       = 2.5
mpres at mpGeophysicalLineThicknessF  = 2.5
mpres at mpNationalLineThicknessF     = 2.5
mpres at mpProvincialLineThicknessF   = 2.5
mpres at mpUSStateLineThicknessF      = 2.5
mpres at mpCountyLineColor            = "Black"
mpres at mpGeophysicalLineColor       = "Black"
mpres at mpNationalLineColor          = "Black"
mpres at mpUSStateLineColor           = "Black"

cnres                       = res
cnres at tfDoNDCOverlay        = False ; Testing
cnres at cnFillOn              = True
cnres at lbLabelBarOn          = True
cnres at cnLinesOn             = False
cnres at lbOrientation         = "Vertical"
cnres at sfDataArray           = ash
cnres at sfXArray              = lonN
cnres at sfYArray              = latN
cnres at sfMissingValueV       = ash at _FillValue
cnres at cnLevelSelectionMode  = "ManualLevels"
cnres at cnFillMode            = "RasterFill"
cnres at cnRasterSmoothingOn   = False
cnres at cnFillOpacityF        = 0.75
cnres at cnMinLevelValF        = 0.
cnres at cnMaxLevelValF        = 2.0
cnres at cnLevelSpacingF       = 0.05

map         = gsn_csm_map(wks,mpres)
ash_contour = gsn_csm_contour(wks,ashNN,cnres)

overlay(map,ash_contour)
draw(map)
frame(wks)

end
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20160908/d6063c9b/attachment.html