TIGER1.sml

# TIGER1.SML
# Create city boundary vector object from TIGER data
# INPUT: vector object imported from U.S. Census Bureau TIGER Line files
# OUTPUT: city boundary vector
# extracts only those TIGER lines that are on city boundaries and assembles
# them into polygons with the city name as attribute. 
### Declarations ###
vector Input, Output;
class DATABASE db;		# handle for polygon database in output vector
class DBTABLEINFO tinfo; # handle for new table in polygon database
class DATABASE linedb;	# handle for line database in input vector
class DBTABLEINFO linetable;	# handle for required table in TIGER line database 
numeric numLines;		# number of lines in Input
numeric numpolys;		# counter for number of polygons created in output vector
numeric left, right;		# variables for city code from left and right side of line
numeric numPoints;			# number of vertices in line to be copied
Array numeric xarray[1];	# arrays to hold coordinates of line vertices
Array numeric yarray[1];
Array numeric records[1];			# array to hold record(s) attached to extracted line
numeric numberofthem;	# number of records attached to input line; should be 1
string name$;				# string variable to hold name of city
numeric recordnumber;	# number of record added to polygon city table
### Process ###
clear();
GetInputVector( Input );
### get output vector and initialize vector toolkit so elements can be added
GetOutputVector( Output, "VectorToolkit,Polygonal" );
CopySubobjects( Input, Output, "GEOREF" );			# copy georeference to output
### open polygon database and create a table to hold city names;
### set attachment type and add field to table
db = OpenVectorPolyDatabase( Output );
tinfo = TableCreate( db, "CITY_NAMES", "Created by SML script" );
tinfo.OneRecordPerElement = 1;		
TableAddFieldString( Output.poly.CITY_NAMES, "City_Names", 40 );
### get handle for line database in Input vector
linedb = OpenVectorLineDatabase( Input );
### check if required line table exits and continue processing lines if it does
if (	!TableExists(linedb, "Geo_Names_P"	) > 0	)  {
		PopupMessage( "Geo_Names_P table missing; abort processing." );
		return;
		}	
	else {
	### get handle for required line table
	linetable = DatabaseGetTableInfo( linedb, "Geo_Names_P" );
	### loop through input lines to identify & copy city boundary lines
	numeric j;			# loop counter = input line element number being processed
	numpolys = 0;		# initialize polygon counter
	numLines = NumVectorLines( Input );
	for j = 1 to numLines {
		SetStatusBar( j, numLines );
 
		### get city codes for left and right side of input line
		left = ( Input.line[j].Basic_Data.FIPS_Pub55Pla_L );
		right = ( Input.line[j].Basic_Data.FIPS_Pub55Pla_R );
		### if city codes don't match, get line vertex list and use to add line to
		### output vector; if added line closes a polygon, a new polygon is created automatically.
		if (left != right) {
			numPoints = GetVectorLinePointList( Input, xarray, yarray, j );
			VectorAddLine( Output, numPoints, xarray, yarray );
			}
		### check if new line closed a polygon by comparing number of polygons to polygon counter;
		### if so, get city name from line database and add to polygon database
		if  ( NumVectorPolys( Output ) > numpolys ) {
			numpolys = numpolys + 1;	# increment polygon counter
			### get record(s) attached to the input line and store city name in name$
			numberofthem = TableReadAttachment( linetable, j, records );
			name$ = TableReadFieldStr( linetable, "Geographic_Name", records[1] );
			if  ( ( name$ == "" )  &&  ( numberofthem > 1 ) ) #hack for bad data in Nebraska otoe county
					name$ = TableReadFieldStr( linetable, "Geographic_Name", records[2] );
			### create new record in polygon table with city name
			recordnumber = TableNewRecord( tinfo, name$ );
			records[1] = recordnumber;
			### attach this new record to the last polygon 
			TableWriteAttachment( tinfo, numpolys, records, 1);
			}
		}
	}
### process output vector to remove island polygons inside city polygons
array numeric islands[1];  # array to hold list of island element numbers for current polygon
numeric size;
size = 100;
array numeric deleteisland[size]; # array to hold cumulative list of island element numbers
numeric numofislands;
numeric k;		# counter for cumulative number of islands (k -1) and index for deleteisland array
k = 1;
for each poly in Output {
	### identify island polygons, delete their attached records, and record
	### element numbers in array deleteisland
	numofislands =  GetVectorPolyIslandList( Output, islands );
	if ( numofislands > 0 ) {
		for  j = 1 to numofislands {
				TableReadAttachment( tinfo, islands[j], records );
				RecordDelete( tinfo, records[1] );
				deleteisland[k] = islands[j];
				k += 1;
				if ( k >= size ) {
					ResizeArrayPreserve( deleteisland, size + 100 );
					size += 100;
					}
			}
		}
	}
if ( k > 1 )
	VectorDeletePolys( Output, deleteisland, k-1);
print( "The number of island polygons deleted was", k-1 );
SetStatusBar( 0, 10 );
VectorValidate( Output );   # validate vector topology
CloseVector( Output );