Projecting Map of US and Canada from GFK Map Datasets

[fastb]

I am creating a map of the US and Canada (excluding lakes and some states and provinces).

I have successfully created the maps of each country, but they are not projected to the right position (see example).

Can someone help me with the PROC GPROJECT statement that I need to use to get these in the right position?

I'm not very familiar with that PROC so I tried several options and they made it worse instead of better.

Thanks!

DM 'LOG;CLEAR;OUT;CLEAR;'; OPTIONS Missing = '' LS=128 PAGESIZE=36 NONUMBER; 
goptions reset=all border cback=white htitle=12pt htext=10pt;

proc gremove data=mapsgfk.canada out=canada_cleaned;
by id1; id id; run;

proc contents data=canada_cleaned; quit;

data map_us; length country $8.;set mapsgfk.us/*states*/;
if statecode = 'AK' or statecode = 'PR' or statecode = 'HI' then delete;
country = 'US';
run;

data map_can; length country $8.; set canada_cleaned;
if ID1 = 'CA-60' or ID1 = 'CA-61' or ID1 = 'CA-62' then delete;
if lake ne 0 then delete;
country = 'CAN';
run;

data map2(drop=statecode ID1);
set map_can map_us;
if country = 'CAN' then state_province_code = ID1;
if country ne 'CAN' then state_province_code=statecode;
run;

FILENAME plots "C:\Map Files";
goptions reset=all 
DEVICE=gif 
xpixels=1200
ypixels=1000
transparency 
GSFNAME=plots; 
ods _all_ close;
ods listing; 
pattern1 v=s color=white;

proc format;
value $state_province_code other = ' ';
run;

proc gmap data=map2 map=map2 all;
id state_province_code; choro state_province_code / cdefault = white nolegend uniform;
format state_province_code $state_province_code. ;
run; quit;
quit;

[ballardw]

Same hint applies from the previous question: Check the max and min of the x y variables by country.

Since the mapsgfk.us does not contain lat and long it is hard to project to the same coordinate system:

data map_us; length country $8.;set mapsgfk.us_states;
if statecode = 'AK' or statecode = 'PR' or statecode = 'HI' then delete;
country = 'US';
run;

data map_can; length country $8.; set canada_cleaned;
if ID1 = 'CA-60' or ID1 = 'CA-61' or ID1 = 'CA-62' then delete;
if lake ne 0 then delete;
country = 'CAN';
run;

data map2(drop=statecode ID1);
set  map_can map_us ;
if country = 'CAN' then state_province_code = ID1;
if country ne 'CAN' then state_province_code=statecode;
run;

proc gproject data= map2 latlong deg eastlong out=mapproj
;
id  state_province_code;
run;


proc format;
value $state_province_code other = ' ';
run;

proc gmap data=mapproj map=mapproj all;
id state_province_code; choro state_province_code / cdefault = white nolegend uniform;
format state_province_code $state_province_code. ;
run; quit;
quit;

The Gproject LATLONG says to use coordinates the Lat Long variables, Deg that they are in degrees, EASTLONG because the default values will have east/west swapped otherwise.

[fastb]

I was initially using the maps.states and maps.canada4 datasets.

However, I switched to mapsgfk.us and mapsgfk.canada because the mapsgfk.canada4 dataset includes the LAKE variable so I can exclude the lakes from the Canada map.

Is there any way to prevent the lakes from displaying using the maps.canada4 dataset?  

[ballardw]

Which lakes? When I ran the modified code I show with the mapsgfk.us_states the only lakes I see are the Great Lakes. And if you remove them you are likely to cause international incidents 😀 Without the Great Lakes I think anyone is going to question the map appearances.

I'm not sure if there is an easy to find map with the political boundary that runs through most of those lakes. You can spend some time adjusting coordinates of polygons to estimate such a line and then remove the appropriate current lake shores.

[GraphGuy · Accepted solution]

See my sample code in https://communities.sas.com/t5/Graphics-Programming/Random-Lines-Across-Map-from-PROC-GMAP/m-p/80536...

[fastb]

Thanks for your help with this GraphGuy!

Made great progress using your code suggestions.  Now I am making the two country borders thicker and darker than the state/province borders.  For some reason the US/Canada border is a double line above Montana.  Can you identify what is causing that?

Thank you!

DM 'LOG;CLEAR;OUT;CLEAR;'; OPTIONS Missing = '' LS=128 PAGESIZE=36 NONUMBER; 
goptions reset=all border cback=white htitle=12pt htext=10pt;

data map_us; set mapsgfk.us_states (where=(statecode not in ('AK' 'HI') and density<=4));
my_id=statecode; country='US'; run;

data map_canada; set mapsgfk.canada (where=(id1 not in ('CA-60' 'CA-61' 'CA-62') and density<=1)); run;

proc gremove data=map_canada out=map_canada;
by id1; id id; run;

data map_canada; set map_canada;
my_id=id1; country='CA'; run;


data mapall; length my_id $8;
set map_us map_canada;
run;

proc sql noprint;
create table mapall_attr as
select unique country, my_id, segment, 1 as fakedata
from mapall;
quit; run;

proc gproject data=mapall out=mapall_proj latlong eastlong degrees project=gnomon;
id country my_id;
run;

proc gremove data=mapall_proj out=anno_outline;
by country notsorted;
id my_id;
run;

data anno_outline; set anno_outline;
by country segment notsorted; length function $8 color $8;
color='gray33'; style='mempty'; when='a'; xsys='2'; ysys='2';
if first.segment then function='poly'; else function='polycont';
run; 


FILENAME plots "C\Output Files";
goptions reset=all 
DEVICE=gif 
xpixels=1200
ypixels=1000
transparency 
GSFNAME=plots; 
ods _all_ close;
ods listing; 
pattern1 v=solid color=white;

proc gmap map=mapall_proj data=mapall_attr anno=anno_outline all;
id country my_id;
choro fakedata / levels=1 nolegend coutline=graydd;
run; quit;

[GraphGuy]

I suspect that one of the two maps doesn't have enough points along the border to follow the desired projection's curve. To try to fix it, you could either use a less-curvy projection (such as cylindri), or you could leave more points in the map (by eliminating the where clause's "density<=something").

The longer explanation: 

The reason it doesn't have enough points is (probably) that when the maps were processed by SAS (before  being shipped to the customers), Proc Greduce was probably run against the unprojected coordinates ... and the coordinates along that part of the US/Canada border pretty much have the same latitude, and therefore Proc Greduce thinks it doesn't need to keep very many data points along that border. But once you start using a projection, and those borders are shown as 'curves' then you need more coordinates to follow the curve.

One of the final projects I inherited and worked on before retiring was the processing of the maps. I changed the algorithm to project the map first, before running Proc Greduce, to try to reduce/eliminate this problem. I believe that the new maps created using my algorithm are the ones currently available for download (and are likely newer than the ones you have in your mapsgfk library). You might try downloading them, and see if it gives you better resolution along that border: https://support.sas.com/en/knowledge-base/maps-geocoding/maps-online.html