Re: Projection methods used in sas map datasets

Darrell_sas — Fri, 07 Nov 2014 19:22:07 GMT

You really should be using MAPSGFK instead of MAPSSAS. MAPSSAS is not changing. It will eventually stop shipping.

MAPSGFK is different than MAPSSAS in some ways. It is more consistant between maps.
The paper explaining the differences is: http://support.sas.com/rnd/papers/sasgf13/GfKMaps2013.pdf
(If this doesn't work, go to: http://support.sas.com/rnd/papers/sasgf13/GfKMaps2013.pdf and then look for the paper "Tips and Tricks: Using the new SAS Map data sets").

You can use use the MAPSGFK.PROJPARM data set to project your annotate data on a map that is already projected (See example 6C in the above documentation).
And you can chang your map with the insets like you had in MAPSSAS with the following program (It will soon be on MapsOnLine):

/*******************************************************************

| SAS (R) is a registered trademark of SAS Institute Inc.

| SAS Institute does not assume responsibility for the accuracy of

| any material presented in this file. |

\*******************************************************************/

/* Re-projecting US, USCITY and USCENTER Datasets */

/*******************************************************************

| This program is a revised version to be used with V8 datasets.

| The code was used to reproject US, USCITY and USCENTER data

| sets from the STATES data set and the LONG/LAT variables in

| the USCITY and USCENTERS data sets.

| It was necessary to re-project these datasets with the purchase

| of new city data.

|*******************IMPORTANT NOTE *******************************

| PLEASE BE AWARE THAT THIS PROGRAM CANNOT BE USED ASIS. ANY NEW

| INFORMATION SHOULD FIRST BE PROCESSED OR UNNECESSARY STEPS

| DELETED. THE NEW USCITY DATA IS NOT INCLUDED.

| |

| This program uses datasets from SAS MAPS library and new city

| information from USCITY dataset.In the end, this job produces |

| three work data sets: US,USCITY, USCENTER and graphics output.

| In the mean time the following naming convention is used for

| work ds: |

| |

| Prefix - Indicates geographic area |

| US - All 50 states |

| CONT - 48 states |

| AK - Alaska |

| HI - Hawaii |

| Middle - Projection / reduction |

| U - Unprojected |

| P - Projected |

| R - Projected and reduced |

| blank - Projected low density |

| Suffix - Indicates data's type(s) |

| MAP - Only map points |

| CTY - Cities |

| CTR - Centers |

| ALL - Map, cities, and centers |

| |

*******************************************************************/

data USUMAP;

set MAPSGFK.US_STATES(drop=x y where=(state ne 72));

FLAG = 1;

X = LONG;

Y = LAT;

if STATE = 02 and

SEGMENT = 1 then output;

else if ((STATE = 25 | STATE = 36 | STATE = 53) and (

SEGMENT = 1)) then output;

else output;

run;

/*********** Process new city information **************************/

data USUCTY;

set mapsGFK.uscity(drop=x y);

X = LONG;

Y = LAT;

SEGMENT = 99;

FLAG = 2;

run;

data USUCTR;

set MAPSGFK.USCENTER(drop=x y);

X = LONG;

Y = LAT;

SEGMENT = 98;

FLAG = 3;

run;

data USUALL;

set USUMAP USUCTY USUCTR;

run;

/*------------------ Seperate Alaska and Hawaii ---------------------*/

data CONTUALL AKUALL HIUALL;

set USUALL;

if STATE = 15 then output HIUALL;

else if STATE = 2 then output AKUALL;

else output CONTUALL;

run;

/*-------------------------- Project -------------------------------*/

proc gproject data=CONTUALL out=CONTPALL degrees eastlong;

id STATE;

run;

proc gproject data=AKUALL out=AKPALL degrees eastlong;

id STATE;

run;

proc gproject data=HIUALL out=HIPALL degrees eastlong;

id STATE;

run;

/*------------------ Determine Scaling Factors ---------------------*/

title 'Range of Alaska';

proc means noprint MIN MAX data=AKPALL;/* find range of proj. Alaska */

var X Y;

output out=AKMEAN min=XMIN YMIN max=XMAX YMAX;

run;

title 'Range of Hawaii';

proc means noprint MIN MAX data=HIPALL; /*find range of proj. Hawaii */

var X Y;

output out=HIMEAN min=XMIN YMIN max=XMAX YMAX;

run;

title 'Range of Continental US';

proc means noprint MIN MAX data=CONTPALL;

var X Y;

output out=CONTMEAN min=CXMIN CYMIN max=CXMAX CYMAX;

run;

/*------------------ Move Alaska and Hawaii -------------------------*/

data AKPALL; /* rescale Alaska to fit under California */

/* --------------------------------------------------------

Scale X to (CXMIN, CXMIN + 0.200 * (CXMAX - CXMIN))

Scale Y to (CYMIN, CYMIN + 0.285 * (CYMAX - CYMIN))

------------------------------------------------------- */

set AKPALL;

if _n_ = 1 then do;

set CONTMEAN;

FACT1 = CXMIN + 0.200 * (CXMAX - CXMIN);

FACT2 = CYMIN + 0.285 * (CYMAX - CYMIN);

retain FACT1 FACT2;

set AKMEAN;

AX = (FACT1 - CXMIN) / (XMAX - XMIN);

AY = (FACT2 - CYMIN) / (YMAX - YMIN);

BX = CXMIN - AX * XMIN;

BY = CYMIN - AY * YMIN;

retain AX BX AY BY;

end;

X = X * AX + BX;

Y = Y * AY + BY;

drop FACT1-FACT2 BX BY AX AY XMAX XMIN YMAX YMIN;

drop CXMAX CXMIN CYMAX CYMIN;

run;

data HIPALL; /* rescale Hawaii to fit under New Mexico */

/* ---------------------------------------------

Scale X to (CXMIN + 0.200 * (CXMAX - CXMIN),

CXMIN + 0.300 * (CXMAX - CXMIN))

Scale Y to (CYMIN, CYMIN + 0.143 * (CYMAX - CYMIN))

--------------------------------------------- */

set HIPALL;

if _n_ = 1 then do;

set CONTMEAN;

FACT1 = CXMIN + 0.200 * (CXMAX - CXMIN);

FACT2 = CXMIN + 0.300 * (CXMAX - CXMIN);

FACT3 = CYMIN + 0.143 * (CYMAX - CYMIN) ;

retain FACT1 FACT2 FACT3;

set HIMEAN;

AX = (FACT2 - FACT1) / (XMAX - XMIN);

AY = (FACT3 - CYMIN) / (YMAX - YMIN);

BX = FACT1 - AX * XMIN;

BY = CYMIN - AY * YMIN;

retain AX BX AY BY;

end;

X = X * AX + BX;

Y = Y * AY + BY;

drop FACT1-FACT3 BX BY AX AY XMAX XMIN YMAX YMIN;

drop CXMAX CXMIN CYMAX CYMIN;

run;

/*-------------------- Separate Data Sets --------------------------*/

data AKPCTY AKPCTR AKPMAP;

set AKPALL;

if FLAG = 1 then output AKPMAP;

else if FLAG = 2 then output AKPCTY;

else output AKPCTR;

run;

data HIPCTY HIPCTR HIPMAP;

set HIPALL;

if FLAG = 1 then output HIPMAP;

else if FLAG = 2 then output HIPCTY;

else output HIPCTR;

run;

data CONTPCTY CONTPCTR CONTPMAP;

set CONTPALL;

if FLAG = 1 then output CONTPMAP;

else if FLAG = 2 then output CONTPCTY;

else output CONTPCTR;

run;

/*-------------------- Assign Densities ----------------------------*/

title 'Reducing Maps';

proc greduce data=HIPMAP out=HIRMAP n1=114;

id STATE;

run;

proc greduce data=AKPMAP out=AKRMAP n1=109;

id STATE;

run;

proc greduce data=CONTPMAP out=CONTRMAP n1=84;

id STATE;

run;

data USRMAP;

set CONTRMAP AKRMAP HIRMAP;

run;

/*-------------------- Strip Maps to Low Density -------------------*/

data HIMAP;

set HIRMAP;

if DENSITY <= 1;

run;

data AKMAP;

set AKRMAP;

if DENSITY <= 1;

run;

data CONTMAP;

set CONTRMAP;

if DENSITY <= 1;

run;

/*-------------------- Recombine States ----------------------------*/

data USMAP; /* add repostioned Alaska and Hawaii to US */

set CONTMAP AKMAP HIMAP;

run;

data USCTY; /* add repostioned Alaska and Hawaii to US */

set CONTPCTY AKPCTY HIPCTY;

run;

data USCTR;/* add repostioned Alaska and Hawaii to US */

set CONTPCTR AKPCTR HIPCTR;

run;

/*-------------------- Sort Data Sets by STATE ---------------------*/

proc sort data=USMAP out=USMAP;

by STATE;

run;

proc sort data=USCTY out=USCTY;

by STATE;

run;

proc sort data=USCTR out=USCTR;

by STATE descending OCEAN;

run;

/*----------------------- Look It Over -----------------------------*/

set USMAP;

label X = 'Projected Longitude in Radians';

label Y = 'Projected Latitude in Radians';

label SEGMENT = 'State Segement Number';

label STATE = 'State FIPS Code';

keep X Y SEGMENT STATE;

run;

data USCTY; /* add temporary variables */

length COLOR $ 8;

set USCTY;

XSYS = '2';

YSYS = '2';

TEXT = '*';

COLOR = 'RED';

run;

set USCTY;

keep X Y LONG LAT CITY CNTYFIPS CITYFIPS STATE CAPITAL

ALT POP PLACE FEATYPE;

run;

set USCTR;

label X = 'Projected Longitude from Radians';

label Y = 'Projected Latitude from Radians';

label LONG = 'Unprojected Longitude in Degrees';

label LAT = 'Unprojected Latitude in Degrees';

label OCEAN = 'Y/N';

label STATE = 'State FIPS Code';

keep X Y LONG LAT ID OCEAN STATE STATECODE;

run;

title 'Contents of US Map';

proc contents data=US;

run;

title 'Contents of USCITY';

proc contents data=USCITY;

run;

title 'Contents of USCENTER';

proc contents data=USCENTER;

run;

/*---------------------- Test Data Sets --------------*/

DATA CTCENTERS(KEEP=state city flag

XSYS YSYS WHEN FUNCTION TEXT X Y COLOR SIZE);

LENGTH FUNCTION COLOR $8 TEXT $25;

RETAIN

WHEN 'A' XSYS YSYS '2' ;

set uscity;

by state;

FLAG='C';

FUNCTION='LABEL';

SIZE=.65;

TEXT=city;

color='RED';

* TEXT='*';

RUN;

DATA SYMBCTR(KEEP=state city flag

XSYS YSYS WHEN FUNCTION TEXT X Y COLOR SIZE);

LENGTH FUNCTION COLOR $8 TEXT $25;

RETAIN

WHEN 'A' XSYS YSYS '2' ;

set uscity;

by state;

FLAG='C';

FUNCTION='LABEL';

SIZE=.65;

TEXT='*';

color='RED';

RUN;

DATA CENTERS(KEEP=state

XSYS YSYS WHEN FUNCTION TEXT X Y COLOR SIZE);

LENGTH FUNCTION COLOR $8 TEXT $25;

RETAIN WHEN 'A' XSYS YSYS '2' FLAG 0;

set uscenter;

FUNCTION='LABEL';

SIZE=.65;

TEXT=statecode;

color='BLACK';

if ocean='Y' then do;

position='6'; output;

function='move';

flag=1;

end;

/* If the labeling coordinates are outside the state (OCEAN='Y'), Annotate */

/* adds the label and prepares to draw the leader line. Note: OCEAN is a */

/* character variable and is therefore case sensitive. OCEAN='Y' must specify */

/* an uppercase Y. */

/* When external labeling is in effect, Annotate */

/* draws the leader line and resets the flag. */

else if flag=1 then do;

function='draw'; size=.5;

flag=0;

end;

output;

RUN;

filename GSASFILE 'us-project9.gif';

goptions reset=pattern dev=GIF

cback=white

noprompt

rotate=LANDSCAPE

gaccess=gsasfile

gsfmode=append

colors=(green);

footnote "As of &SYSDATE";

ods listing close;

ods html file="us-project9.html" (title="USA State Map")

nogtitle nogfootnote;

%macro show2(dsn,note,st,anno);

PROC GMAP

%if &st ne %str( ) %then %do;

DATA=&dsn(where=(state in(&st) ) )

map=&dsn(where=(state in(&st ) ))

ANNO=&anno;

%end;

%else %do;

DATA=&dsn

map=&dsn

ANNO=&anno;

%end;

ID state ;

CHORO state / DISCRETE NOLEGEND COUTLINE=BLACK;

PATTERN V=s r=100;

title; footnote "&note";

run;

%mend;

%show2(%str(work.us),%str(work.us reduced using text),17,

%str(ctcenters(where=(state in(17) and city='Chicago'))));

%show2(%str(work.us),%str(work.us reduced using symbols),17,

%str(symbctr(where=(state in(17) and city='Chicago'))));

%show2(%str(work.usrmap),%str(work.us unreduced using text),17,

%str(ctcenters(where=(state in(17)))));

%show2(%str(work.usrmap),%str(work.us unreduced using symbols),17,

%str(symbctr(where=(state in(17)))));

%show2(%str(work.us),%str(work.us using text),%str(29),

%str(ctcenters(where=(state in(29)))));

%show2(%str(work.us),%str(work.us using symbols),%str(29),

%str(symbctr(where=(state in(29)))));

%show2(%str(work.us),%str(work.us using text),%str(29),

%str(ctcenters(where=(state in(29)

and city in:('Maryvil' 'St.Lou' 'Kansas') ))));

%show2(%str(work.us),%str(work.us using symbols),%str(29),

%str(symbctr(where=(state in(29)

and city in:('Maryvil' 'St.Lou' 'Kansas') ))));

%show2(%str(work.us),%str(work.us using text),,centers);

quit;

ods html close;

ods listing;

topic Re: Projection methods used in sas map datasets in Graphics Programming

Projection methods used in sas map datasets

Re: Projection methods used in sas map datasets

Re: Projection methods used in sas map datasets

Re: Projection methods used in sas map datasets

Re: Projection methods used in sas map datasets

Re: Projection methods used in sas map datasets