Create a polygon and test whether point(s) lie within?

[genemroz]

I need to create a polygon and test whether a number of points lie within.  I've created the polygon  (from lat/long coordinates transformed to a Cartesian xy grid)  using the POLYGON function available in the  SG Annotation Function Dictionary.  But I don't know how to test whether a x,y point lies within or without.  Will PROC GINSIDE accept a polygon generated by the Polygon Function?  It seems dubious given that the annotation function is about drawing and not about "mapping" as seems to be a requirement for PROC GINSIDE.  Looking for some guidance on how best to proceed.

Thank you,

Gene

[Reeza]

So you're trying to now convert your polygon to a spatial file format that GINSIDE will be able to use with your point data?

[genemroz]

I think that's what needs to happen but I don't know how to execute it.

Gene

[ballardw · Accepted solution]

GINSIDE wants a map data set.

A map dataset comprises ordered pairs of coordinates with an identifier for the polygon. Best is to provide a closed polygon where the last pair of coordinates is the same as the beginning.

data examplemap;
   input id x y;
datalines;
1   0 0
1   0 10
1   5 15
1   8 10
1   8 0
1   0 0
2   12  12
2   12  20
2   18  20
2   18  12
2   12  12
;


data exampledata;
   input otherid x y;
datalines;
1  0  .5
2  1  4
3  0  20
4  4  12
5  10 10
6  13 14
7  20 1
;

proc ginside map =examplemap
             data=exampledata
             out =result
             includeborder
;
   id id;
run;

      

The output data set has the coordinates of the Exampledata, and indicator 1/0 if the value is on the border, and the OTHERID of the data is matched with the ID from the Map data when the point is inside the area for the ID.

So the Map data set doesn't have to be too complex as long as the Map and Data sets use the same coordinate systems.

[genemroz]

Thanks, ballardw, for your prompt reply.  I think I can make this solution work for my application and will mark it as accepted.  But if I run in to trouble, well, I'll be back...with more questions.

Thanks again,

Gene

[genemroz]

Well, shoot!  I marked the wrong solution as "accepted".  How do I fix this?

Gene

[GraphGuy]

Proc Ginside requires x/y variables, but there's nothing saying they have to be anything 'special' (such as projected x/y values, etc). Therefore you could simply rename or duplicate your long/lat variables as x/y. 🙂

Here's a simple example, where I just duplicate long/lat as x/y:

data my_coords;
input lat long;
x=long; y=lat;
datalines;
35.8252831 -78.7590949
38.8975995 -77.0366965
;
run;

data my_map; set mapsgfk.us_states (drop=x y);
x=long; y=lat;
run;

proc ginside data=my_coords map=my_map out=my_coords;
id statecode;
run;

proc print data=my_coords;
var lat long statecode;
run;