Getting country using longitude and latitude

Question

First of all I need to mention that I do not have much experience using SAS in geo data. I have a dataset with longitude and latitude and I need to see the distribution by countries. I tried to Google and the only solution I found is this: https://communities.sas.com/t5/Base-SAS-Programming/Is-there-a-way-to-get-County-using-latitude-and-longitude/td-p/145233 It works pretty well for USA data, but the ratio for the rest of the world is very poor (below 10%). I use this code:

proc ginside data=gps map=maps.world out=gpscounties;
    id id cont;
run;

Where GPS is my data with longitude and latitude. Is this 10% ratio is normal or I miss something? Or may be there are other better ways to solve my problem?

I'm on SAS 9.4 TS1M2.

Here is an example:

data gps;
  input x y;
datalines;
143.21  -33.494
119.306 26.0614
119.306 26.0614
143.21  -33.494
113.25  23.1167
139.751 35.685
113.25  23.1167
133.935 34.6617
133.935 34.6617
133.051 35.4722
;
run;

Can you provide a small dataset of examples that don't work (or, that sporadically work) ? Feel free to muck them up so they're not as precise to protect PII/privacy/etc., or even get generic ones by looking around google maps for a while, but a runnable example in the question of exactly what you're doing would help. Use DATALINES in the question to input them, please, rather than linking to a document or pictures. — Joe
Also, please include your version of SAS (9.2, 9.3, 9.4, and TS1M3 or whatever the release is - should show up in your intialization). — Joe

Joe Joe · Accepted Answer · 2016-10-13T17:15:47

At least part of your problem is you have lat/long but you need polar radians.

data gps;
  input x y;
   x=x*arcos(-1)/180;
  x=x*(-1);
  y=y*arcos(-1)/180;
  datalines;
143.21  -33.494
119.306 26.0614
119.306 26.0614
143.21  -33.494
113.25  23.1167
139.751 35.685
113.25  23.1167
133.935 34.6617
133.935 34.6617
133.051 35.4722
;
run;

However, that's not all of your problem as it still doesn't work for me. I suspect you need to use GPROJECT to get the cartesian coordinates into x/y coordinates.

Since you're using 9.4, though, you have the GFK maps, which should make this much easier as they save their projection information!

See the following:

data gps;
  input x y;
  rownum = _n_;
  datalines;
143.21  -33.494
119.306 26.0614
119.306 26.0614
143.21  -33.494
113.25  23.1167
139.751 35.685
113.25  23.1167
133.935 34.6617
133.935 34.6617
133.051 35.4722
;
run;

proc gproject data=gps out=projected
  parmin=mapsgfk.projparm parmentry=world;
  id rownum;
run;


proc ginside data=projected map=mapsgfk.world out=countries;
  id idname id;
run;

That seems to work for me.

Getting country using longitude and latitude

1 Answers