Re: question about country outline

Shinn,

Its a bit hard to say why. If I recall your previous plot correctly,
the shapefile contained boundaries of internal regions, in addition to
the country boundary. This likely confuses the in/out testing. If you
are trying to do something like example shapefile_5.ncl, then what is
needed is a coherent list of the coordinates of just the country
boundary. In shapefile_5.ncl, this was obtained through a shapefile
that contained *only* that boundary.

It would be easiest if you can find such a shapefile for your area of
interest. If not however, you might be able to extract just the
required polygon by subsetting the shapefile based upon a field in the
database. Trying doing a "printVarSummary" on the shapefile file-
variable. Your looking for a field that might encode boundary types,
and most importantly, you want to make sure that the geometry_type
attribute is "polygon". Example shapefile_1.ncl shows how to query
the shapefile's database.

If that's not clear and you want to send me your shapefiles, I'll take
a look and see what might be possible.

Good luck!
Rick

On Dec 9, 2010, at 10:50 PM, shinn wong wrote:

> Dear Rick,
> Thank you very much for your help. Now I can draw the country
> outline, but I've got another question. I want to check whether a
> point is within the country. So I use the function
> gc_inout(lat,lon,mrb_lat,mrb_lon). But there are points outside the
> country that show positive results. Does anyone know why? Thanks in
> advance for any given help.
> Shinn
> 2010/12/10 Rick Brownrigg <brownrig@ucar.edu>
> Hi Shinn,
>
> The trouble stems from the single call to:
>
> line_data = gsn_add_polyline(wks, map_data, mrb_lon, mrb_lat, lnres)
>
> which will treat all the coordinates in mrb_lot/mrb_lat as a
> continuous line. You want to draw the line segments that make up the
> countries individually, as in shapefile_3.ncl. In particular, you
> want to replace the above call with something like:
> segments = f->segments
> geometry = f->geometry
> segsDims = dimsizes(segments)
> geomDims = dimsizes(geometry)
> geom_segIndex = f@geom_segIndex
> geom_numSegs = f@geom_numSegs
> segs_xyzIndex = f@segs_xyzIndex
> segs_numPnts = f@segs_numPnts
> numFeatures = geomDims(0)
> lines = new(segsDims(0),graphic) ; array to hold polylines
> segNum = 0 ; Counter for adding polylines
> do i=0, numFeatures-1
> 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_polyline(wks, map_data,
> mrb_lon(startPT:endPT), \
>
> mrb_lat(startPT:endPT), lnres)
> segNum = segNum + 1
> end do
> end do
>
> Check out the discussion at the top of the shapefiles page if the
> above code doesn't make sense. Hope that helps...
> Rick
>
>
>
>
>
> On Dec 9, 2010, at 2:00 AM, shinn wong wrote:
>
>> Dear all,
>> I've tried to plot the outline boundary of some countries using
>> data obtained from DIVA-GIS. The code and the output file can be
>> seen in the attachment. However, the output is quite different from
>> that of example 4 in the shapefile example page. Does anyone know
>> how can i just plot the outline of the country like those in
>> example 4 instead of all the data points? How can i get the
>> coordinate of the outline of the country? Many thanks.
>> Shinn
>> <
>> aa1
>> .ncl
>> >
>> <
>> shapefiles1.000001
>> .png>_______________________________________________
>> ncl-talk mailing list
>> List instructions, subscriber options, unsubscribe:
>> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>
>

_______________________________________________
ncl-talk mailing list
List instructions, subscriber options, unsubscribe:
http://mailman.ucar.edu/mailman/listinfo/ncl-talk
Received on Fri Dec 10 09:09:18 2010