[ncl-talk] Overlaying Contour on WRF Domain
Sean Egan
sdegan at alaska.edu
Thu Sep 8 10:54:38 MDT 2016
Hi Mary,
Thanks for the tips. Sorry about the typo - the NN shouldn't proceed the
ash variable. I've attached (and copy/pasted) the script as well as the
variable summaries to this e-mail. If you would like to look at the HDF
file, you can download it by connecting via FTP to emc2.iarc.uaf.edu with
username set to anonymous. No password is required. The file is named
file1.hdf.
Thanks,
Sean
---------------------------------
VARIABLE SUMMARIES
---------------------------------
ncl 61> printVarSummary(ash)
Variable: ash
Type: float
Total Size: 13363200 bytes
3340800 values
Number of Dimensions: 2
Dimensions and sizes: [lines | 900] x [elements | 3712]
Coordinates:
Number Of Attributes: 10
algorithm_name : VOLCAT volcanic ash retrieval
algorithm_version : none
algorithm_index : 10
reference : Pavolonis and Heidinger
units : none
scaling_method : 0
scale_factor : 1
add_offset : 0
_FillValue : -999
hdf_name : volcld_ret_14_15_16_ash_mass_loading
ncl 62> printVarSummary(latN)
Variable: latN
Type: float
Total Size: 13363200 bytes
3340800 values
Number of Dimensions: 2
Dimensions and sizes: [900] x [3712]
Coordinates:
Number Of Attributes: 1
_FillValue : -32768
ncl 63> printVarSummary(lonN)
Variable: lonN
Type: float
Total Size: 13363200 bytes
3340800 values
Number of Dimensions: 2
Dimensions and sizes: [900] x [3712]
Coordinates:
Number Of Attributes: 1
_FillValue : -32768
---------------------------------
MAX AND MIN VALUES
---------------------------------
ncl 64> printMinMax(ash,0)
(0) min=0 max=4.58412
ncl 65> printMinMax(latN,0)
(0) min=27.5848 max=81.2629
ncl 66> printMinMax(lonN,0)
(0) min=-79.6368 max=79.6368
---------------------------
START OF SCRIPT
---------------------------
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"
a = addfile("file1.hdf","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
b = addfile("wrfout_d01_2010-04-16_00:00:00","r")
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 = wrf_map_resources(b,mpres)
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 sfXArray = lonN
cnres at sfYArray = latN
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,ash,cnres)
overlay(map,ash_contour)
draw(map)
frame(wks)
end
-------------------------
SCRIPT ERRORS
-------------------------
pacman13 % ncl testing.ncl
Copyright (C) 1995-2013 - All Rights Reserved
University Corporation for Atmospheric Research
NCAR Command Language Version 6.1.2
The use of this software is governed by a License Agreement.
See http://www.ncl.ucar.edu/ for more details.
warning:mpNestTime is not a valid resource in map at this time
warning:start_lat is not a valid resource in map at this time
warning:start_lon is not a valid resource in map at this time
warning:end_lat is not a valid resource in map at this time
warning:end_lon is not a valid resource in map at this time
On Thu, Sep 8, 2016 at 8:19 AM, Mary Haley <haley at ucar.edu> wrote:
> Sean,
>
> First, I recommend using "wrf_map_resources" to set the necessary WRF map
> resources for you. This will allow you to get rid of all this code:
>
> 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"
>
> and replace it with one line:
>
> mpres = wrf_map_resources(b,mpres)
>
> The other resources you are setting:
>
> 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"
>
> are also set by "wrf_map_resources", but since you are using different
> values than what wrf_map_resources uses, I recommend leaving these alone.
> Just make sure you set all these resources *after* you call
> wrf_map_resources. You can always do a "print(mpres)" after you call
> wrf_map_resources to see exactly what values it is setting for you.
>
> As for the "ash" array you are trying to plot: please look at your latN
> and lonN arrays to make sure they look correct, and that they in "degrees"
> units. That is, the lat values should range from -90 to 90, and the
> longitude values from -180 to 180 or 0 to 360:
>
> printMinMax(latN,0)
> printMinMax(lonN,0)
>
> You have the sfXArray and sfYArray arrays set correctly, but you do NOT
> want to set this:
>
> cnres at sfDataArray = ash
>
> You also do not need this, because the gsn_csm_contour function will set
> it for you:
>
> cnres at sfMissingValueV = ash at _FillValue
>
> I'm not sure if this is a typo:
>
> ash_contour = gsn_csm_contour(wks,ashNN,cnres)
>
> but shouldn't "asnNN" be "ash"?
>
> Finally, you want to be sure that "ash" doesn't already have some kind of
> coordinate arrays attached to it.
>
> If you continue to have problems, it would really help if we could get the
> latest version of your script and your data, or, if you can't provide the
> data, then insert a bunch of "printVarSummary" and "printMinMax" statements
> so we can see what your data looks like:
>
> printVarSummary(ash)
> printVarSummary(latN)
> printVarSummary(lonN)
> printMinMax(ash,0)
> printMinMax(latN,0)
> printMinMax(lonN,0)
>
> On Thu, Sep 8, 2016 at 9:57 AM, Genito Amos Maure <genito.maure at uem.mz>
> wrote:
>
>> Hi Sean.
>> Is that "i" before the zero in the expression
>> lat = b->XLAT(i0,:,:) ;Not a nested domain
>> a typo?
>>
>>
>>
>> On Thu, 8 Sep 2016 07:48:52 -0800, Sean Egan wrote
>> > 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
>>
>>
>>
>> _______________________________________________
>> ncl-talk mailing list
>> ncl-talk at ucar.edu
>> List instructions, subscriber options, unsubscribe:
>> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>>
>
>
--
PO1 Sean D. Egan, US Navy
PhD Candidate
Geophysical Institute
University of Alaska Fairbanks
903 Koyukuk Drive
Fairbanks, AK 99775
sdegan at alaska.edu | 907 474-5483
IARC Office 334 | www.gi.alaska.edu/~sdegan
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20160908/5fc2d019/attachment.html
-------------- next part --------------
A non-text attachment was scrubbed...
Name: testing.ncl
Type: application/octet-stream
Size: 2166 bytes
Desc: not available
Url : http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20160908/5fc2d019/attachment.obj
More information about the ncl-talk
mailing list