Main Page/Research/Smart Maps/dwelling units

Create a smart map raster that represents the count or density of dwelling units

1. add a field to the table for dwelling units per km².

   greater precision than per ft²

2. calculate the field, e.g.

   du_per_km2 = [num_units]*10763910.42/[shape_area]

3. convert to raster using the field

   make a note of cell size, e.g., 10 m cell size = 32.808 ft

4. focalsum of the raster

   note the radius, e.g., 10 cells = 100 m

   gives sum of du/km² for the circle

5. divide by area of the circle to give density

   note unit conversions

The following source code can be used. Requirements:

1. Existence of drive W for output. Use subst in Windows for this.
2. Modify script in ArcToolbox to include arguments/parameters, matching the sys.argv[] lines.

Note: this assumes you want the density to be in units/ha.

<source lang="python">

1. ---------------------------------------------------------------------------
2. focal_du_per_ha.py
3. Created on: Tue May 18 2010 02:48:36 PM
4. (generated by ArcGIS/ModelBuilder)
5. Usage: focal_du_per_ha <parcels> <Value_field> <Cellsize> <Neighborhood> <du_per_ha>
6. ---------------------------------------------------------------------------

1. Import system modules

import sys, string, os, arcgisscripting

1. Create the Geoprocessor object

gp = arcgisscripting.create()

1. Check out any necessary licenses

gp.CheckOutExtension("spatial")

1. Load required toolboxes...

gp.AddToolbox("C:/Program Files (x86)/ArcGIS/ArcToolbox/Toolboxes/Spatial Analyst Tools.tbx") gp.AddToolbox("C:/Program Files (x86)/ArcGIS/ArcToolbox/Toolboxes/Conversion Tools.tbx")

1. output dir

outdir = "w:/" if os.path.exists(outdir) == 0: gp.AddMessage("Path does not exist.") exit()

1. Script arguments...

parcels = sys.argv[1] if parcels == '#':

parcels = "PIN_ADDRESS" # provide a default value if unspecified

Value_field = sys.argv[2] if Value_field == '#':

Value_field = "du_per_km2" # provide a default value if unspecified

Cellsize = sys.argv[3] if Cellsize == '#':

Cellsize = "32.808" # provide a default value if unspecified

StrCellAreaMeters = str((float(Cellsize) / 3.2808)**2)

Neighborhood = sys.argv[4] if Neighborhood == '#':

Neighborhood = "Circle 10 CELL" # provide a default value if unspecified

du_per_ha = sys.argv[5] if du_per_ha == '#':

du_per_ha = outdir + "du_per_ha.img" # provide a default value if unspecified

1. Local variables...

duperkm2 = outdir + "duperkm2.img" FocalSumR = outdir + "FocalSumR.img"

1. Process: Polygon to Raster...

gp.AddMessage("Polygon to raster\n") gp.PolygonToRaster_conversion(parcels, Value_field, duperkm2, "CELL_CENTER", "NONE", Cellsize)

1. Process: Focal Sum...

gp.AddMessage("FocalSum\n") gp.FocalStatistics_sa(duperkm2, FocalSumR, Neighborhood, "SUM", "DATA")

1. Process: DU per hectare...

gp.AddMessage("Du per ha\n") gp.SingleOutputMapAlgebra_sa(FocalSumR + " / ( (Pi) * sqr( " + StrCellAreaMeters + " ))", du_per_ha) </source>