[ncl-talk] Draw only outer boundaries of a shapefiles polygon

[Rick Brownrigg]

Hi,

Unfortunately, shapefiles do not contain topological information, and in
particular, adjacency info; the lines that make up a border are repeated
for the two areas it separates. Thus there's not an easy way to readily
determine which borders are interior to your outer regions, and which are
borders between two outer regions. Since you are drawing the outer regions
in distinctly different colors, perhaps just not drawing any edges will
suffice?

Rick

On Thu, May 30, 2019 at 1:33 PM Tabish Ansari <tabishumaransari at gmail.com>
wrote:

> Hi
>
> I'm trying to map three distinct regions, two of which are a collection of
> several districts, and create boundaries for them on a map. I have
> shapefiles for each district but drawing polylines using all the member
> districts would also draw internal boundaries which I don't want. I only
> want to draw outer boundaries for each region. Please see the attached
> plot.
>
> My script is pasted here:
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
> *begin  wks  = gsn_open_wks("png","pulse-districts") ; send graphics to
> PNG file  res                     = True  res at gsnDraw             = False
>     ; don't draw yet  res at gsnFrame            = False       ; don't advance
> frame yet  res at gsnMaximize         = True       ; maximize plot in frame;
>  res at gsnPaperOrientation = "portrait"  res at vpKeepAspect        = True
> res at mpProjection        = "LambertConformal"    ; choose projection
> res at mpLambertParallel1F = 20          ; first parallel
> res at mpLambertParallel2F = 50          ; second parallel
> res at mpLambertMeridianF  = 110         ; meridian  res at mpLimitMode         =
> "Corners"   ; corner method of zoom  res at mpLeftCornerLatF  = 35
> res at mpLeftCornerLonF  = 109.35  res at mpRightCornerLatF = 44.41
> res at mpRightCornerLonF = 124.94  res at pmTickMarkDisplayMode = "Always"  ;
> turn on tickmarks  res at tmXBMajorOutwardLengthF = 0.0
> res at tmXTMajorOutwardLengthF = 0.0  res at tmYLMajorOutwardLengthF = 0.0
> res at tmYRMajorOutwardLengthF = 0.0   res at tiMainFontHeightF  = 0.025
> res at tiYAxisFontHeightF = 0.01  res at tiXAxisFontHeightF = 0.01
> res at tmXBLabelFontHeightF = 0.01  plot = gsn_csm_map(wks,res)   ; Create
> map, but don't draw it
> yet.;*************************************************; Section to add
> polygons to map.;*************************************************  f =
> addfile("/home/tabish/Shapefiles/China/CHN_adm2.shp", "r")   ; Open
> shapefile  ;; Read data off shapefile;  segments = f->segments  geometry =
> f->geometry  segsDims = dimsizes(segments)  geomDims = dimsizes(geometry);;
> Read global attributes  ;  geom_segIndex = f at geom_segIndex  geom_numSegs  =
> f at geom_numSegs  segs_xyzIndex = f at segs_xyzIndex  segs_numPnts  =
> f at segs_numPnts  lines       = new(segsDims(0),graphic)   ; Array to hold
> polygons  numFeatures = geomDims(0)  name =  f->NAME_2  lon   = f->x  lat
> = f->y  plres             = True       ; resources for polylines
> plres at gsEdgesOn   = True       ; draw border around polygons
> plres at gsEdgeColor = "grey"    ;  plres at gsFillColor      = "Navy";
>  plres at gsFillOpacityF    = 0.5    segNum = 0  do i=0, numFeatures-1       ;
> color assignment (probably a better way to do this?)     if
> (name(i).eq."Beijing") then        plres at gsFillColor = "tomato"
>  else if
> (name(i).eq."Langfang").or.(name(i).eq."Baoding").or.(name(i).eq."Shijiazhuang").or.(name(i).eq."Xingtai").or.(name(i).eq."Handan").or.(name(i).eq."Tangshan").or.(name(i).eq."Tianjin").or.(name(i).eq."Cangzhou").or.(name(i).eq."Hengshui").or.(name(i).eq."Dezhou").or.(name(i).eq."Binzhou").or.(name(i).eq."Dongying").or.(name(i).eq."Liaocheng").or.(name(i).eq."Jinan").or.(name(i).eq."Zibo")then
>       plres at gsFillColor = "limegreen"          else
> plres at gsFillColor = "dodgerblue"     end if     end if     startSegment =
> geometry(i, geom_segIndex)     numSegments  = geometry(i, geom_numSegs)
>  do seg=startSegment, startSegment+numSegments-1        startPT =
> segments(seg, segs_xyzIndex)        endPT = startPT + segments(seg,
> segs_numPnts) - 1        lines(segNum) = gsn_add_polygon(wks, plot,
> lon(startPT:endPT),  \
>  lat(startPT:endPT), plres)        segNum = segNum + 1     end do
> plres at gsFillColor      = "Grey"  end doprint("done making
> polygons")delete(lat)delete(lon);---Attach some dummy text strings  txres
>             = True  txres at txFontHeightF = 0.020  txres at txJust        =
> "CenterCenter";STARTING TO MARK CITIES WHERE APEC CONTROLS WERE
> IMPLEMENTEDlat=40.15lon=116.4074text = "Beijing"text_id1 =
> gsn_add_text(wks,plot,text,lon,lat,txres)lat=38lon=116.2text = "Near
> Neighbourhood"text_id8 =
> gsn_add_text(wks,plot,text,lon,lat,txres)lat=41lon=112.2text = "Far
> Neighbourhood"text_id8 = gsn_add_text(wks,plot,text,lon,lat,txres);
> Maximize output in frame. This will draw everything and advance; the
> frame.;  maximize_output(wks,False)end*
>
> Thanks very much.
>
> Cheers,
>
> Tabish
>
> Tabish U Ansari
> PhD student, Lancaster Environment Center
> Lancaster Univeristy
> Bailrigg, Lancaster,
> LA1 4YW, United Kingdom
> _______________________________________________
> ncl-talk mailing list
> ncl-talk at ucar.edu
> List instructions, subscriber options, unsubscribe:
> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20190530/48b07297/attachment.html>