[ncl-talk] help with plotting shapefile polygons

Dave Allured - NOAA Affiliate dave.allured at noaa.gov
Tue Apr 5 16:39:28 MDT 2022


You may have encountered a known bug whereby the X and Y coordinates from a
shapefile are swapped with each other.  We think this was caused by a
change in an underlying library.  Please see this report for details and a
possible easy workaround.  Let us know if that works.

    https://github.com/NCAR/ncl/issues/176

What version of NCL are you using, and what was the source of your NCL
installation?


On Tue, Apr 5, 2022 at 4:29 PM Vanessa Vincente via ncl-talk <
ncl-talk at mailman.ucar.edu> wrote:

> Hello,
>
> I am having trouble plotting polygon data from a shapefile of a burn scar
> perimeter over a map of Colorado with counties.  I am using the example
> shapefiles_15.ncl from the ncl page
> <https://www.ncl.ucar.edu/Applications/shapefiles.shtml>. I am able to
> successfully plot the state of Colorado and its counties, but not the burn
> scar perimeter.  Burn scar shapefile can be found here
> <https://ftp.wildfire.gov/public/incident_specific_data/rocky_mtn/2018/SpringCreek/IR/20180715/PerimeterFileProvidedByIncident/>.
>  Can you please help? Thanks in advance!
>
> Below is a copy of my ncl script:
>
> load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/shea_util.ncl"
> load "./shapefile_utils.ncl"
>
> ;----------------------------------------------------------------------
> ; This function creates a cylindrical equidistant map of Colorado.
> ;----------------------------------------------------------------------
>
> function create_colorado_map(wks,res,draw_ncl_outlines)
> local a, mpres
> begin
>   mpres                       = res
>
>   mpres at gsnMaximize           = True
>   mpres at gsnPaperOrientation   = "portrait"
>
>   mpres at gsnDraw               = False
>   mpres at gsnFrame              = False
>
>   mpres at mpFillOn              = False
>
> ;---Turn on fancier tickmark labels.
>   mpres at pmTickMarkDisplayMode = "Always"
>   mpres at tmXBLabelFontHeightF  = 0.008      ; smaller tickmark labels
>
> ;---Zoom in on area of interest
>   mpres at mpLimitMode           = "LatLon"
>   mpres at mpMinLatF             =   37
>   mpres at mpMaxLatF             =   41
>   mpres at mpMinLonF             = -109.05
>   mpres at mpMaxLonF             = -102.05
>   mpres at mpFillOn              = False
>
>   if(draw_ncl_outlines) then
>     mpres at mpOutlineOn           = True
>     mpres at mpOutlineBoundarySets    = "AllBoundaries"
>     mpres at mpDataBaseVersion        = "MediumRes"
>     mpres at mpDataSetName            = "Earth..4"      ; U.S. counties
>   else
>     mpres at mpOutlineOn              = False
>   end if
>
> ;---Create map.
>   map = gsn_csm_map(wks,mpres)
>
>   return(map)
> end
>
> ;--------------------------------------------------
> ; Main code
> ;--------------------------------------------------
> begin
>   wtype          = "png"              ; send graphics to PNG file
>   wtype at wkWidth  = 2000
>   wtype at wkHeight = 2000
>   wks = gsn_open_wks(wtype,"shapefiles")
>
>   ncl_version = get_ncl_version()
>
> ;---Create two maps of Colorado
>   res                   = True
>   res at tiMainFontHeightF = 0.015
>
>   res at tiMainString  = "Colorado counties - shapefile"
>   map_shp = create_colorado_map(wks,res,False)
>
>   res at tiMainString = "Colorado counties - NCL (version " + ncl_version +
> ")"
>   map_ncl  = create_colorado_map(wks,res,True)
>
> ;---Add shapefiles to one of the maps
>   lnres             = True
>   lnres at gsLineColor = "red"
>
>   plot_line_shp =
> gsn_add_shapefile_polylines(wks,map_shp,"./usa/gadm40_USA_2.shp",lnres)
>
> ;---Add a slightly transparent marker to both maps to show location of
> county updates
>   mkres                  = True
>   mkres at gsMarkerIndex    = 16        ; filled dot
>   mkres at gsMarkerOpacityF = 0.5       ; make the marker half transparent
>   mkres at gsMarkerColor    = "red"
>   mkres at gsMarkerSizeF    = 50.
>
>   counties_lat_center =   39.8
>   counties_lon_center = -104.9
>
>   plot_marker_ncl =
> gsn_add_polymarker(wks,map_ncl,counties_lon_center,counties_lat_center,mkres)
>   plot_marker_shp =
> gsn_add_polymarker(wks,map_shp,counties_lon_center,counties_lat_center,mkres)
>
> print("Adding polygons to ncl map and shp map...")
>
> filename = "Perimeter_20180715_SpringCreek_COHUX001313.shp"
>
>  ;---Set some options for the polygons
>    bscar = True
>    bscar at gsLineColor = "Blue"
>    bscar at gsLineThicknessF = 3.0  ; 3x thickness
>
> id1 = gsn_add_shapefile_polygons(wks,map_ncl,filename,bscar)
> id2 = gsn_add_shapefile_polygons(wks,map_shp,filename,bscar)
>
> ;---Panel both plots. Markers, lines, and polygons will be drawn.
>   pres             = True
>   pres at gsnMaximize = True
>   gsn_panel(wks,(/map_ncl,map_shp/),(/1,2/),pres)
>
> end
> --
> Vanessa Vincente
> Associate Scientist
> The COMET® Program
> University Corporation for Atmospheric Research
>
> *My working hours may differ from your own.  Please do not feel obligated
> to reply outside of standard business hours.*
> _______________________________________________
> ncl-talk mailing list
> ncl-talk at mailman.ucar.edu
> List instructions, subscriber options, unsubscribe:
> https://mailman.ucar.edu/mailman/listinfo/ncl-talk
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://mailman.ucar.edu/pipermail/ncl-talk/attachments/20220405/28e2bf3c/attachment.html>


More information about the ncl-talk mailing list