/*******************************************************************\
| Copyright (C) 2005 by SAS Institute Inc., Cary, NC, USA. |
| |
| 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. |
\*******************************************************************/
/********************************************************************
* Plot Geocoded Points with PROC GMAP *
* Example 5 from SUGI 30 Presentation: *
* "Cheap Geocoding: SAS/GIS and Free TIGER Data" *
********************************************************************/
/********************************************************************
* Set location and libref for data and html map page. *
********************************************************************/
%let PATH=C:\Temp;
libname CENTROID "&PATH";
/********************************************************************
* Create a customer address data set. *
********************************************************************/
data CENTROID.CUSTOMERS (label = 'Customer Addresses');
input address $ 1-22 /* House address (not used) */
name $ 24-48 /* Customer name */
city $ 49-63 /* City name (not used) */
statecode $ 64-65 /* State abbreviation (not used) */
zip /* Customer ZIP Code (numeric) */
revenue; /* Total revenue for customer */
format revenue dollar9.0; /* Show revenue in dollars */
label revenue = 'Revenue'; /* Label to appear on map legend */
cust_ID = _n_; /* Assign customer ID number */
cards;
123 E. Wine Street J. Cheever Loophole Mullins SC 29574 985.33
611 Jackson Street Cuthbert J. Twillie Vidalia GA 30474 2533.25
300 North Odum Street Prof. Quincy Wagstaff Pembroke NC 28372 1899.45
5575 Peachtree Pkwy. Kaspar Gutman Norcross GA 30022 33.44
1525 N. Shenandoah Ave Dr. Hugo Z. Hackenbush Front Royal VA 22630 987.02
2310 Cannon Bridge Rd Charlie Allnut Cordova SC 29039 1459.65
305 Mill Creek Drive Larson E. Whipsnade Fuquay Varina NC 27526 4566.98
1911 W. Main St. Guillermo Ugarte Richmond VA 23220 17.45
2313 Dale Ave. Dewey, Cheatem & Howe Roanoke VA 24013 989.20
7825 Old Middlesex Rd Capt. Geoffrey Spaulding Bailey NC 27807 1798.33
1935 E. Frederick St Joel Cairo Blacksburg SC 29702 2566.45
900 Stony Landing Rd Charles Blutoski Easley SC 29641 783.25
26 Ancho Street Brigid O'Shaughnessy Burlingame CA 94000 543.99
1373 London Bridge Rd Fred C. Dobbs Virginia Beach VA 23451 898.25
626 Mohican Trail Rufus T. Firefly Wilmington NC 28409 123.45
216 12th Street Otis B. Driftwood Lynchburg VA 24504 98.56
35 W. Chatham Str Elwood Blues Pittsboro NC 27312 661.13
100 W. Church Street Duke Mantee Cherryville NC 28021 4531.87
;
run;
/********************************************************************
* Reorder customer data by ZIP for merging. *
********************************************************************/
proc sort data = CENTROID.CUSTOMERS;
by ZIP;
run;
/********************************************************************
* Assign ZIP centroid coordinates to customer addresses. *
********************************************************************/
data CENTROID.CUSTOMERS (label = 'Customers Located by ZIP Centroid Geocoding')
CENTROID.MISSING (label = 'Customers Not Found by ZIP');
merge CENTROID.CUSTOMERS (in=a) /* Customer addresses */
SASHELP.ZIPCODE; /* ZIP and centroid X/Y data */
by ZIP; /* ZIP is the linkage variable */
if a; /* Keep only customer obs */
/*****************************************************************
* If this customer has a bad ZIP Code, set this obs aside. *
*****************************************************************/
if X = . or Y = . then output CENTROID.MISSING;
/*****************************************************************
* Otherwise this customer's ZIP Code was found and has valid *
* X/Y values. *
*****************************************************************/
else output CENTROID.CUSTOMERS;
run;
/********************************************************************
* Determine which states contain customers. The base map data *
* will be subset to display only these states. *
********************************************************************/
proc freq data = CENTROID.CUSTOMERS noprint;
tables state / out = _states_ (label = 'States containing customers');
run;
/********************************************************************
* Generate a where clause to subset the map data for only the *
* states containing customers. The clause is placed in the macro *
* variable &CoverWhere for use below. *
********************************************************************/
data _NULL_;
set _STATES_ (rename = (state = curState)) end = last;
length stateList $ 256; /* States with customers */
retain stateList "";
/*****************************************************************
* If the first state, initialize the where clause string. *
*****************************************************************/
if _n_ = 1 then stateList = trim(left(curState));
/*****************************************************************
* If not the first state, append it to the string. *
*****************************************************************/
else stateList = compress(stateList || "," || curState);
/*****************************************************************
* If final state, put where clause string into macro variable. *
*****************************************************************/
if last then do;
call symput( 'StateList', trim(stateList) );
end;
run;
/********************************************************************
* Use where clause to subset US map for states with customers. *
********************************************************************/
data CENTROID.MAP;
set MAPS.STATES; /* Base map of US states */
where state in (&stateList) & /* List of states wanted */
density < 3; /* Thin the point data */
run;
/********************************************************************
* Create an ANNOTATE data set to draw stars at customer sites. *
********************************************************************/
data CENTROID.ANNOTATE; /* Annotate data set */
set CENTROID.CUSTOMERS; /* Customer data with X/Y */
length html $ 500 /* Chart-tip text string */
function style color $ 8;
xsys = '2';
ysys = '2';
hsys = '3';
when = 'A';
/*****************************************************************
* Convert decimal degrees to radians to match map data. *
*****************************************************************/
x = -x * (atan(1)/45);
y = y * (atan(1)/45);
/*****************************************************************
* Draw customer point marker with a star. *
*****************************************************************/
function = 'label';
style = 'special';
text = 'M'; /* Use star symbol */
color = 'cx2e19dd'; /* Star color */
size = 5; /* Star height */
position = '5'; /* Center star at X/Y */
/*****************************************************************
* Add chart-tip text for this point. This information will be *
* displayed when you hover over the customer point. *
*****************************************************************/
html = 'alt='|| quote(
'Name: '||trim(left(name))||'0D'x||
'City: '||trim(left(city))||', '||trim(left(statecode))||'0D'x||
'Revenue: '||trim(left(put(revenue, dollar8.2)))||' ');
run;
/********************************************************************
* Combine map and annotate data sets. *
* They must be projected together. *
********************************************************************/
data _COMBINED_;
set CENTROID.MAP CENTROID.ANNOTATE;
run;
/********************************************************************
* Project the map and annotate data. *
********************************************************************/
proc gproject data = _COMBINED_
out = _COMBINED_;
id state;
run;
/********************************************************************
* Separate projected map and annotate into two data sets. *
********************************************************************/
data CENTROID.MAP CENTROID.ANNOTATE; /* Separate data sets */
set _COMBINED_; /* Combined projected data */
if density = . then
output CENTROID.ANNOTATE; /* Points */
else
output CENTROID.MAP; /* Base map */
run;
/********************************************************************
* Use SAS/GRAPH GMAP Procedure and ODS to create web page. *
********************************************************************/
ods results off;
ods html close; /* Start new ODS output */
ods html path="&PATH" /* Location of output file */
body="Customers.html"; /* Output file name */
goptions reset=all /* Clear graphics settings */
device=gif /* Output image type */
cback=cxc3bdb6; /* Map background color */
pattern1 v=solid c=cxeeeee6 r=50; /* State fill color */
title1 h=2 f=centb "Customer Locations"; /* Map titles */
title2 h=1.5 f=centbi "Geocoded by ZIP Centroid with SAS";
footnote1 h=1.3 j=l f=centb " Map by PROC GMAP";
footnote2 h=1.3 j=l f=centb " Hover to see point values";
/********************************************************************
* Generate the map. *
********************************************************************/
proc gmap data = CENTROID.MAP /* Response data set needed */
map = CENTROID.MAP /* Projected map data */
anno = CENTROID.ANNOTATE; /* Projected point data */
id state; /* Map variable */
choro state / /* Draw choropleth map */
coutline = cx739436 /* State outline color */
des = "Customers" /* grseg entry label */
nolegend; /* Suppress legend */
run;
quit;
ods listing close; /* End ODS output */
/********************************************************************
* Open the map in your default browser. Hover over the address *
* points to see the customer information. *
********************************************************************/
dm 'wbrowse "&PATH\Customers.html"';
|
