/*******************************************************************\
| 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. |
\*******************************************************************/
/********************************************************************
* Street Level Geocoding *
* Example 2 from SUGI 30 Presentation: *
* "Cheap Geocoding: SAS/GIS and Free TIGER Data" *
* *
* Notes: *
* 1) You will have to download the Delaware TIGER files from *
* the Census Bureau prior to running this program. They *
* should be placed in the folder referenced by the *
* &TIGER_DATA macro variable. The .RT1 and .RT2 files *
* should be unzipped. *
* 2) If you run this code a second time, you must first delete *
* the old files from the MAP_DATA location. *
********************************************************************/
/********************************************************************
* Set paths for your map and data. *
********************************************************************/
%let TIGER_DATA=C:\Temp; /* Downloaded TIGER files */
%let MAP_DATA=C:\Temp; /* Map spatial data sets and catalog */
%let LOOKUP_DATA=C:\Temp; /* Geocoding lookup data sets */
%let ADDRESS_DATA=C:\Temp; /* Address data set to geocode */
/********************************************************************
* Set librefs. *
********************************************************************/
libname STREET "&MAP_DATA";
libname LOOKUP "&LOOKUP_DATA";
libname ADDRESS "&ADDRESS_DATA";
/********************************************************************
* Create a customer address data set. *
********************************************************************/
data ADDRESS.CUSTOMERS;
input address $ 1-23 /* House number and street name */
customer $ 24-48 /* Customer name */
zip /* Customer ZIP Code (numeric) */
city $ 55-70 /* City name */
state $ 71-72 /* State abbreviation */
revenue; /* Total revenue for customer */
format revenue dollar9.0; /* Show revenue in dollars */
label revenue = 'Revenue'; /* Label appears on map legend */
cust_ID = _n_; /* Assign customer ID number */
cards;
115 E. Water St J. Cheever Loophole 19901 Dover 985.33
760 Moose Lodge Road Cuthbert J. Twillie 19934 Camden 2533.25
200 S. Madison Str Prof. Quincy Wagstaff 19801 Wilmington DE 1899.45
4701 Limestone Road Kaspar Gutman 19808 Wilmington 33.44
2117 N 4th St Dr. Hugo Z. Hackenbush 19363 Oxford PA 987.02
133 Silver Lake Dr Charlie Allnut 19971 Rehoboth Beach DE 1459.65
11 SE Front Street Larson E. Whipsnade 19963 Milford DE 4566.98
402 Nylon Boulevard Guillermo Ugarte . Seaford DE 17.45
363 E Commerce St Dewey, Cheatem & Howe . Smyrna DE 989.20
5595 Polly Branch Rd Capt. Geoffrey Spaulding 19975 Selbyville DE 1798.33
1209 Coastal Highway Joel Cairo 19944 Fenwick Island DE 2566.45
2899 Arthursville Rd Charles Blutoski 19953 Hartly DE 783.25
41 Bramhall St Fred C. Dobbs . DE 898.25
9320 Old Racetrack Rd Rufus T. Firefly . Delmar DE 123.45
281 W Commerce Str Otis B. Driftwood 19955 Kenton 98.56
211 Blue Ball Road Gordon Miller 21921 Elkton MD 661.13
3893 Turkey Point Rd Duke Mantee 19980 Woodside DE 4531.87
;
run;
/********************************************************************
* Determine the states that the addresses are in. *
********************************************************************/
proc freq data = ADDRESS.CUSTOMERS;
tables state / nocum;
run;
/********************************************************************
* Look in the Output Window and note that most of the addresses *
* are in Delaware. Go to the Census Bureau web site and download *
* the TIGER files for Delaware's three counties: *
* http://www.census.gov/geo/www/tiger/ *
* Put the files in the directory referenced by the TIGER_DATA *
* macro variable and unzip the .RT1 and .RT2 files so they can be *
* imported into SAS/GIS. *
********************************************************************/
/********************************************************************
* Set batch import parameters common to all counties. *
********************************************************************/
%let imp_type = TIGER; /* Import source file type */
%let maplib = STREET; /* Libref for imported map */
%let mapcat = DELAWARE; /* Catalog for imported entries */
%let mapname = CUSTOMERS; /* Name of imported map entry */
%let spalib = STREET; /* Libref to write spatial data sets */
%let spaname = DELAWARE; /* Base name for spatial data sets */
/********************************************************************
* Import first Delaware county and create initial map data. *
********************************************************************/
%let cathow = CREATE; /* Catalog not existing */
%let spahow = CREATE; /* Initialize data sets */
filename tiger1 "&TIGER_DATA\TGR10001.RT1"; /* County TIGER file 1 */
filename tiger2 "&TIGER_DATA\TGR10001.RT2"; /* County TIGER file 2 */
dm 'af c = sashelp.gisimp.batch.scl;'; /* Invoke import */
/********************************************************************
* Import second Delaware county and add to initial map data. *
********************************************************************/
%let cathow = UPDATE; /* Catalog exists */
%let spahow = APPEND; /* Data sets exist */
filename tiger1 "&TIGER_DATA\TGR10003.RT1"; /* County TIGER file 1 */
filename tiger2 "&TIGER_DATA\TGR10003.RT2"; /* County TIGER file 2 */
dm 'af c = sashelp.gisimp.batch.scl;'; /* Invoke import */
/********************************************************************
* Import third Delaware county add to map data. *
********************************************************************/
filename tiger1 "&TIGER_DATA\TGR10005.RT1"; /* County TIGER file 1 */
filename tiger2 "&TIGER_DATA\TGR10005.RT2"; /* County TIGER file 2 */
dm 'af c = sashelp.gisimp.batch.scl;'; /* Invoke import */
/********************************************************************
* Modify the imported chains data set to remove character values *
* from the address variables on the street chains. *
********************************************************************/
data STREET.DELAWAREC;
set STREET.DELAWAREC; /* Chains data set */
retain chars 'ABCDEFGHIJKLMNOPQRSTUVWXYZ'; /* Unwanted characters */
if _n_ > 1 then do; /* Skip header record */
FrAddL = compress(upcase(FrAddL), chars);
FrAddR = compress(upcase(FrAddR), chars);
ToAddL = compress(upcase(ToAddL), chars);
ToAddR = compress(upcase(ToAddR), chars);
end;
run;
/********************************************************************
* Delaware base map is complete. *
* Set up input parameters for the batch geocoder. *
********************************************************************/
%gcbatch(
glib = LOOKUP, /* Street data sets libref */
newdata = YES, /* (Re)create lookup data? */
mname = STREET.DELAWARE.CUSTOMERS, /* Base map of streets */
geod = ADDRESS.CUSTOMERS, /* Address data to geocode */
nv = CUSTOMER, /* Who's at the address */
av = ADDRESS, /* House number/street name */
cv = CITY, /* City name */
sv = STATE, /* State abbreviation */
zv = ZIP, /* ZIP Code (5-digit) */
pv = TRACT); /* AREA var wanted from map */
/********************************************************************
* Invoke the batch geocoding program. *
********************************************************************/
dm 'AF C=SASHELP.GIS.GEOCODEB.SCL';
/********************************************************************
* Show geocoding results on a bar chart. *
********************************************************************/
ods results off;
ods html close; /* Start new ODS output */
ods html path="&ADDRESS_DATA" /* Location of output file */
body="Geocoding Results.html"; /* Output file name */
goptions reset=all /* Clear graphics settings */
device=gif /* Output image type */
ftext=swiss; /* Chart font */
pattern1 color=pay value=solid; /* Bar color */
axis1 label=(height=1.3 /* x-axis label size */
'Address Status'); /* x-axis label text */
axis2 label=(angle=-90 rotate=90 /* y-axis label is vertical */
height=1.3 /* y-axis label size */
'Percent'); /* y-axis label text */
title1 h=1.6 "Delaware Customers"; /* Chart titles */
title2 h=1.3 "Geocoding Results";
/********************************************************************
* Generate the bar chart. *
********************************************************************/
proc gchart data = ADDRESS.CUSTOMERS; /* Address data set */
vbar _status_ / /* Midpoint (x-axis) var */
descending /* Order of results */
type = pct /* Response (y-axis) var */
outside = pct /* Label on top of bars */
inside = freq /* Label inside of bars */
coutline = black /* Bar outline color */
maxis = axis1 /* x-axis */
raxis = axis2; /* y-axis */
run;
quit;
ods html close; /* End ODS output */
/********************************************************************
* Open html chart in default browser. *
********************************************************************/
dm 'wbrowse "&ADDRESS_DATA\Geocoding Results.html"';
/********************************************************************
* Set up input parameters for the batch importer to add the *
* geocoded addresses to the Delaware map as points. *
********************************************************************/
%let imp_type = GENPOINT; /* Import customers as points */
%let maplib = STREET; /* Map libref */
%let mapcat = DELAWARE; /* Map catalog name */
%let mapname = CUSTOMERS; /* Map entry name */
%let spalib = STREET; /* Spatial data set libref */
%let spaname = DELAWARE; /* Spatial entry name */
%let cathow = UPDATE; /* Use existing entries */
%let spahow = APPEND; /* Add to spatial data sets */
%let nidvars = 0; /* Put points in one layer */
%let infile = ADDRESS.CUSTOMERS; /* Geocoded address data set */
/********************************************************************
* Invoke the batch import program. *
********************************************************************/
dm 'af c=SASHELP.GISIMP.BATCH.SCL';
/********************************************************************
* Use the GIS Procedure to modify the Delaware map layers and *
* add a theme to the address points using customer revenue. *
********************************************************************/
proc gis; /* Invoke Procedure */
catalog STREET.DELAWARE; /* Map location */
spatial DELAWARE; /* Spatial entry name */
/*****************************************************************
* Modify visual parameters of customer layer. *
*****************************************************************/
layer update CUSTOMERS / /* Layer entry name */
where = 'node = 1' /* Layer definition */
des = 'Geocoded addresses' /* Layer entry label */
default = (point = (color = cxff0000 /* Symbol color */
font = marker /* Symbol font */
character = 'V' /* Point symbol (star) */
size = 10)); /* Symbol height */
/*****************************************************************
* Add a theme using customer revenue to the layer. *
*****************************************************************/
/*--- Add a theme using customer revenue to the layer. ----------*/
layer update CUSTOMERS / /* Layer entry name */
theme = (create /* Make a new theme */
link = cust_rev /* Data set link name */
dataset = ADDRESS.CUSTOMERS /* Theme data set */
datavar = cust_ID /* Linking variable */
themevar = Revenue /* Thematic variable */
nlevels = 3 /* Number of ranges */
point = ((level = first /* First theme range: */
color = cxff0000 /* Symbol color */
size = 6) /* Symbol height */
(level = last /* Third theme range: */
color = cx00a900 /* Symbol color */
size = 16) /* Symbol height */
blendcolor /* Middle range color */
blendsize)); /* Middle point size */
/*****************************************************************
* Modify visual parameters of tract layer. *
*****************************************************************/
layer update TRACT / /* Layer entry name */
default = (area = (color=cxdddddd) /* Fill color */
outline = (color=cx739436));/* Outline color */
/*****************************************************************
* Modify visual parameters of county layer. *
*****************************************************************/
layer update COUNTY / /* Layer entry name */
des = 'Delaware Counties' /* Label for entry */
default = (line = (width=2)); /* Boundary line size */
run;
/*****************************************************************
* Modify the map visually. *
*****************************************************************/
map update CUSTOMERS / /* Map entry name */
layers = (COUNTY TRACT CUSTOMERS /* Set layer order */
STREET WATER RAIL BG )
layerson = (COUNTY TRACT CUSTOMERS) /* Layers to display */
select = (CUSTOMERS) /* Selectable layers */
cback = cxcae1e9 /* Background color */
des = 'Geocoded Customer Map' /* Label for entry */
action = (create /* Make action to open */
type = spatial /* Spatial Info Window */
name = Spatial /* when a map feature */
when = immediate); /* is selected */
run;
/*****************************************************************
* Add labels to the map window. *
*****************************************************************/
maplabel create /
force /* Replace existing DS */
dataset = STREET.LABELS /* label data set name */
text = 'Customer Locations' /* Label text */
map = STREET.DELAWARE.CUSTOMERS /* Map entry */
attach_to = window /* Pin label to window */
position = (top right) /* Label location */
color = cxff0000 /* Text color */
offscale = (5 km/cm on) /* Turn off on zoom in */
notransparent /* Show in text block */
font = 'Swis721 Blk BT-12pt-Roman-Normal'; /* Text font */
maplabel create /
text = 'Geocoded by SAS/GIS' /* Label text */
map = STREET.DELAWARE.CUSTOMERS /* Map entry */
attach_to = window /* Pin label to window */
position = (bottom left) /* Label location */
color = cxff0000 /* Text color */
offscale = (5 km/cm on) /* Turn off on zoom in */
notransparent /* Show in text block */
font = 'Swis721 Blk BT-11pt-Roman-Normal'; /* Text font */
/*****************************************************************
* Label each CUSTOMER point with the customer name. *
*****************************************************************/
layerlabel create /
layer = CUSTOMERS /* Layer to label */
map = STREET.DELAWARE.CUSTOMERS /* Map entry */
sas_var = customer /* Values for text */
color = cxff0000 /* Text color */
onscale = (2.9 km/cm on) /* Scale to turn on at */
notransparent; /* Show in text block */
/*****************************************************************
* Set scale for layer to turn on as we zoom. *
*****************************************************************/
layer update STREET / /* Layer entry name */
onscale = (2.9 km/cm on); /* Scale to turn on at */
run;
quit;
/********************************************************************
* Open the map in SAS/GIS. See the SUGI presentation paper for *
* information on using the map. *
********************************************************************/
dm 'gis map = STREET.DELAWARE.CUSTOMERS;';
|