Extremely slow processing time to map lat/lon data onto census shapefile with ArcPy?

Question

I am a newcomer to ArcGIS and Python.

What I am trying to do is very simple: map 400K business locations using lat/lon coordinates onto a shapefile of California census tracts, and then export the data from the spatial join into a text file so that I have a spreadsheet containing the census tract that each business falls into. The input file of lat/lon data is a csv with 400K records, 5 variables and is ~15MB.

I have successfully completed this task interactively in ArcMap 10.2.1, but I would like to have a record of my process for future replication, hence my foray into ArcPy. In ArcMap, this process takes a few minutes. However, my ArcPy script takes multiple hours to complete. It appears to be the CopyFeatures step that requires significant computer resources and so I am not including ArcMap details or ArcPy code beyond this step. Note also that I have run my full ArcPy script on a subset of my data (the first 20 records of my 15MB file) and it works fine and finishes in just a few seconds. Given that, I know the issue is likely with the size of my file, but I have not been able to determine or successfully implement a work around to decrease the memory usage of my code.

Are there any ways I can modify my script to decrease processing time for this simple task?

I'm hopeful there's a way to make this as efficient with ArcPy as it is with ArcMap.

ArcMap steps:

1. Import shapefile containing census tracts

   • Add data -> choose shapefile

2. Import and prepare lat/lon data for spatial join

   • File -> Add Data -> Add XY Data -> select my 15GB csv -> X Field (longitude) -> Y Field (latitude) -> Coordinate System (WGS 1984)

   • Right click on XY data -> Data -> Export data -> All Features -> Output feature class is a .shp saved in my GIS workspace -> upon completion of this step, a pop up displays asking whether want to add to map as layer and I choose Yes (don't think this is significant however)

ArcPy script:

# Import modules
import arcpy, csv
import os
from arcpy import env

# Set environment settings
env.workspace = "C:/Users/ccwright/Desktop/GIS_workspace"
arcpy.env.overwriteOutput = True

# Set local variables:
in_table = "infousa_latlon_small.csv"
x_coords = "longitude"
y_coords = "latitude"
out_layer = "lat_lon_layer"
feature_class = "lat_lon_feature_class.shp"
census_tracts_shp = "census_tracts_CA_only/US_tract_2000_CA_only.shp"
spatial_join = "spatial_join.shp"
spatial_join_output_txt = "spatial_join_data.txt"

# Set coordinate system to WGS-1984:
WKID = 4326 
sr = arcpy.SpatialReference()
sr.factoryCode = WKID
sr.create()
env.outputCoordinateSystem = sr

# Process: Make XY Event Layer
arcpy.MakeXYEventLayer_management(in_table, x_coords, y_coords, out_layer, sr, "")

# Process: Copy Features
arcpy.CopyFeatures_management(out_layer, feature_class, "", "0", "0", "0")

Answer 1

Try changing the CopyFeatures_management to FeatureClassToFeatureClass_conversion(). It also takes a feature layer and is essentially going to do the same thing you would when you right-click on your XY Event Layer in the TOC, go to Data and select Export.

Also, you could probably simplify setting your output spatial coordinate system by just directly referencing the spatial reference object. no need to create the object. The env object knows how to do all that stuff already.

# Set output coordinate system to WGS-1984:
env.outputCoordinateSystem = arcpy.SpatialReference(4326)

Answer 2

I suggest to change

arcpy.SpatialJoin_analysis(census_tracts_shp, feature_class, spatial_join, "JOIN_ONE_TO_MANY", "KEEP_ALL","", "INTERSECT", "", "")

to

arcpy.SpatialJoin_analysis(feature_class,census_tracts_shp, spatial_join, "JOIN_ONE_TO_ONE", "KEEP_ALL","", "INTERSECT", "", "")