nav[aria-label="Primary Navigation"] { padding: 0; & ul { list-style: none; width: 100%; display: flex; flex-direction: row; justify-content: start; align-items: start; gap: 30px; padding: 0; & li { margin: 0; } & ul li { list-style: none; } } }

Multi batch spatial joins

2Quiker

I have over 20 layers that I need to spatially join, doing separate arcpy.SpatialJoin_analysis for each one would be time consuming. fc1 is the layer I want to multi batch the other layers in the database to.

I have tried the following but it only spatially joins the first layer in the database. Is it possible to do a multi batch spatial joins? I am going about it the wrong way?

import arcpy, os  
from arcpy import env
from datetime import datetime as d
startTime = d.now()

arcpy.env.workspace = r"C:\Temp\Default2.gdb"  
fc1 = r"C:\Temp\Default2.gdb\Parcels"

arcpy.env.overwriteOutput = True
fcs = arcpy.ListFeatureClasses() #gets the feature classes
for fc in fcs:
    print(fc)
    arcpy.SpatialJoin_analysis(fc1,fc, "SP_Test", "JOIN_ONE_TO_ONE", "KEEP_ALL", "", "INTERSECT", "", "")
print ('(Elapsed time: ' + str(d.now() - startTime)[:-3] + ')')

fc1

Find more posts tagged with

Accepted answers

All comments

Lyle-GIS-Analyst

I came to this thread (four years later woo) looking for a solution to do my own batch spatial join, where I could get the count of points within polygons with just the one polygon layer, but about 11 different point layers. This code currently makes a copy for each spatial join, but I think that can be easily worked around- I just don't want this to break right now and I spent 6 hours figuring this out so bear with me, I was also doing this within ArcGIS Pro's Python window, so I didn't need to set the workspace. I sincerely hope this helps others figure this out faster than I did:

import time #I used a function at the end to wait 3 seconds just because I find that
#when Python in Pro works a bit too fast it can mess up, especially with spatial joins. This isn't really necessary.
#first create a list with each layer you want to join to your main layer
pointfc = ["pointlayer1", #0 - I use notes after lists I iterate through to help me remember which is which, especially if I add some later that point to separate spots in either list. This example looks nicer than some real uses for it.
            "pointlayer2", #1
            "pointlayer3", #2
            "pointlayer4", #3
            "pointlayer5", #4
            "pointlayer6", #5
            "pointlayer7", #6
            "pointlayer8", #7
            "pointlayer9", #8
            "pointlayer10", #9
            "pointlayer11"] #10
#second, I made a list of field names, this is because I hate ESRI's field mapping function
fieldlist = ["fieldname1", #0 -> 0
             "fieldname2", #1 -> 1
             "fieldname3", #2 -> 2
             "fieldname4", #3 -> 3
             "fieldname5", #4 -> 4
             "fieldname6", #5 -> 5
             "fieldname7", #6 -> 6
             "fieldname8", #7 -> 7
             "fieldname9", #8 -> 8
             "fieldname10", #9 -> 9
             "fieldname11"] #10 ->10
#Next, I made a string list that I could alter later with the field mapping from the
#layer I am joining to, I left it as it was just to show how I made it but yours should not be the exact same
fieldmapping = ['id "ID" true true false 4 Long 0 0,First,#,polygonfc,id,-1,-1;',
   'sectorid "SectorID" true true false 50 Text 0 0,First,#,polygonfc,sectorid,0,49;',
   'publicstatus "PublicStatus" true true false 4 Long 0 0,First,#,polygonfc,publicstatus,-1,-1;',
   'market "market" true true false 10 Text 0 0,First,#,polygonfc,market,0,9;',
   'chassis "chassis" true true false 25 Text 0 0,First,#,polygonfc,chassis,0,24;']

#this is where we iterate through the points to join them one at a time
for i in range(len(pointfc)):
   print(fieldlist[i]) #not necessary, just used to help identify which step I was on
#here we append the next field mapping that is needed, this one is putting in the
#desired field name and getting a count of the "ticket" field because that's what I
#wanted to count, your field will look different but the idea is to get the right
#field name for the initial joining.
   fieldmapping.append(fieldlist[i]+' "'+fieldlist[i]+'"'+' true true false 12 Double 0 0,Count,#,i,ticket,0,11;')
#I recognize the "i,ticket" part there maybe is not being used correctly, but this still
#worked right so I really don't know what to say
#then we make the holy grail of field mapping without using ESRI's function
#fmj is Field Mapping Join, and this just puts the list back together into
#that one ugly long string for the join code to look at normally
   fmj = str(''.join(fieldmapping))
   print(fmj) #I used this to make sure I was looking at things right, it isn't necessary
   if i == 0: #because I was going through each one separately and making a new layer each time, the first time had to be slightly different.
        arcpy.analysis.SpatialJoin(
            target_features="polygonfc",
            join_features=pointfc[i],
            out_feature_class=r"C:/blah/blah/path/blah/polygonfc"+str(i), #used str(i) here just to differentiate the step and notate which point set I was on
            join_operation="JOIN_ONE_TO_ONE",
            join_type="KEEP_ALL",
            field_mapping=fmj,#here goes that field mapping variable, looks so much nicer this way
            match_option="INTERSECT",
            search_radius=None,
            distance_field_name="",
            match_fields=None
            )
   else: #this is the actual, normal operation, it targets the previous join to go off of.
        arcpy.analysis.SpatialJoin(
            target_features="polygonfc"+str(i-1), #grabs previous join to build off of
            join_features=pointfc[i],
            out_feature_class=r"C:/blah/blah/path/blah/polygonfc"+str(i),
            join_operation="JOIN_ONE_TO_ONE",
            join_type="KEEP_ALL",
            field_mapping=fmj,
            match_option="INTERSECT",
            search_radius=None,
            distance_field_name="",
            match_fields=None
            )
   fieldmapping[i+5] = (fieldlist[i]+' "'+fieldlist[i]+'"'+' true true false 12 Double 0 0,First,#,polygonfc,'+fieldlist[i]+',0,11;') #this part is very important
#the i+5 is based on the number of attributes present in the feature, which grows by
#one each time this iterates through since we are keeping the count of each previous join.
#We mainly have to make sure the end part changes from "Count,#,i,ticket,0,11;" to "First,#,polygonfc,'+fieldlist[i]+',0,11;"
#so that we don't accidentally set all fields in the joins to the last join done.
#This part tripped me up quite a bit until I figured that out.
   time.sleep(3) #here is that 3 second sleep, you may not need it, I may not need it, but I am more comfy with it there
#after that it will loop back up and grab the next join until there are no more joins.
#this code isn't perfect, I know there is likely a better way to not make a number of features equal to the joins done
#but it works and it is a start, and I wanted to provide that to others who may be able to work with it better.

2Quiker

I have tried to implement some field mapping to change some field names and it works on the first one but then errors. I am assuming that it is because it's is trying to change the field again? How can I by pass this on the second go around?

Error line 33.

FieldMappings: Error in getting field map from field mapping for GetFieldMap

import arcpy, os  
from arcpy import env
from datetime import datetime as d
startTime = d.now()

###Working 7/2/2021
arcpy.env.overwriteOutput = True

arcpy.env.workspace = r"C:\Temp\Default.gdb"  
fc1 = r"C:\Temp\fc1.shp"


fcs = arcpy.ListFeatureClasses() #gets the feature class from database
counter = 1

for fc in fcs:
###Field Mapping
    fieldmappings = arcpy.FieldMappings()

    # Add all fields from inputs.
    fieldmappings.addTable(fc1)
    fieldmappings.addTable(fc)

    # Name fields you want. Could get these names programmatically too.
    keepers = ["Field1", "Field2", "Field3", "Field4","Field5", "Field6"] # etc.

    # Remove all output fields you don't want.
    for field in fieldmappings.fields:
        if field.name not in keepers:
            fieldmappings.removeFieldMap(fieldmappings.findFieldMapIndex(field.name))
    
    # Rename fields
    fieldmap = fieldmappings.getFieldMap(fieldmappings.findFieldMapIndex("DXF_LAYER"))
    field.name = "Field1A"
    field.type = 'Text'
    #field.aliasName = "mean_city_pop"
    
    fieldmap.outputField = field
    #fieldmap.mergeRule = "Count"
    fieldmappings.replaceFieldMap(fieldmappings.findFieldMapIndex("Field1"),fieldmap)
    ####
    ####
    
    #print(fc)
    out_fc = "SP_Test" + str(counter)

    spjoin= arcpy.SpatialJoin_analysis(fc1,fc, out_fc, "JOIN_ONE_TO_ONE", "KEEP_ALL", fieldmappings, "INTERSECT", "", "")
    counter += 1
    fc1 = out_fc
print ('(Elapsed time: ' + str(d.now() - startTime)[:-3] + ')')

Anonymous User

Yes, that could help.

altfields is a helper method that takes a list of tuples and runs them through arcpy.AlterField_management() to rename them:

def altfields(inLayer, flds):
    """
    Support function to alter fields
    """
    for f in flds:
        arcpy.AlterField_management(inLayer, f[0], f[1], f[1])

its use is something like:

altfields(yourlayer, [('DMSLat', 'latitude'), ('panel', 'firm_panel')])

'updateAddress' is the target featureclass that you want to join the other featureclasses to.

'subdivision' is the feature class that you want to join to the target featureclass.

Both of these are just variables pointing to featureclasses, which you'll change to match yours.

I see ESRI's example has changed and has improved some over the one tree example, but there is a lot of additional things happening in the example that is distracting from being able to clearly see how to make a field map.

If you want to rename the fields in the fieldmapping, you can do so by doing this:

fldMap = arcpy.FieldMappings()
... # map like fields and exclude fields you dont want 

#-----------------------
# create specific mapped fields:
specificfldmap = arcpy.FieldMap() # create FieldMap object

# add the target featureclass and its field name that you want to map.
specificfldmap.addInputField(targetfc, 'targetfcfldname') 

# add the joining featureclass and its field name that you want to map.
specificfldmap.addInputField(joinfc, 'joinfcfldname') 

# get the outputfield object
specfldobj = specificfldmap.outputField 

# set the desired outputfield name
specfldobj.name = 'newoutputfieldname'

# replace the outputfield object with the modifed version
specificfldmap.outputField = specfldmapobj 

# add your fieldmap to your fieldMappings object
fldMap.addFieldMap(specificfldmap)

# proceed with using the fldMap

You can iterate over a list of tuples for this part if you know the fieldnames and what you want to be mapped to in the target featureclass.

2Quiker

Where are you defining 'altfields', 'updateAddress' and 'subdivisions'?

2Quiker

To have one feature class with certain fields and rename some fields. I can keep certain fields with the following but I need to change some fields names and just the whole fieldmappings gives me a headache.

fieldmappings = arcpy.FieldMappings()

    # Add all fields from inputs.
    fieldmappings.addTable(fc1)
    fieldmappings.addTable(fc)

    # Name fields you want. Could get these names programmatically too.
    keepers = ["Field1", "Field2", "Field3", "Field4","Field5", "Field6"] # etc.

    # Remove all output fields you don't want.
    for field in fieldmappings.fields:
        if field.name not in keepers:
            fieldmappings.removeFieldMap(fieldmappings.findFieldMapIndex(field.name))

I need to be able to rename "Field1", "Field2", "Field3", "Field4","Field5", "Field6" to something else, "AFieldA", "BFieldB", "CFieldC", "DFieldD","EFieldE", "FFieldF". No specifically these names but to something else. I also need to use a merge rule and join delimiter to one layer in the database.

blake.terhune

Is your goal to have one feature class with everything smashed into it (as many spatial joins) or to have separate feature classes, each spatial joined to only Parcels?

2Quiker

I guess I am confused about the who process. I was thinking./hoping I could do a spatial join all at once and spatial join would add the new fields from every feature class to the end of fc1. But it sounds like it will create a new feature class for every feature class in the database, maybe I am better off just doing individual spatial joins so I can do the field mapping.

import arcpy, os  
from arcpy import env
import traceback

from datetime import datetime as d
startTime = d.now()

arcpy.env.workspace = r"C:\Temp\Default2.gdb"  
fc1 = r"C:\Temp\Taxparcels1.shp"

arcpy.env.overwriteOutput = True
fcs = arcpy.ListFeatureClasses() #gets the feature class names 

# Creating a for loop to iterate over featureclasses
for i in range(0,len(fcs)):
    fcs_s = fcs[i]
    print (fcs_s)
    arcpy.SpatialJoin_analysis(fc1,fcs_s, "SP_Test", "JOIN_ONE_TO_ONE", "KEEP_ALL", "", "INTERSECT", "", "")

DavidPike

Yes sorry, merge actually makes no sense, I'm not sure why I thought it would work.

The below code should keep writing the outputs into your workspace, with the last feature class being the one you need e.g. "SP_Test20". You may want to adapt it to process these in memory and then return only the final feature.

Also as Blake say's - remove fc1 from your workspace and put it into a separate geodatabase, as it will be iterating over itself (or put some logic i.e. if fc !=... : , but I would say the best thing to do is just take it out of there).

This is very simple and won't account for removing shape_length, shape_area etc. or any field mappings and join rules, it simply returns the first matching feature.

import arcpy, os  
from arcpy import env
from datetime import datetime as d
startTime = d.now()

arcpy.env.workspace = r"C:\ArcGISHome\multipleSpatialJoin\b.gdb" 
fc1 = r"C:\ArcGISHome\multipleSpatialJoin\a.gdb\a"

arcpy.env.overwriteOutput = True
fcs = arcpy.ListFeatureClasses() #gets the feature classes
counter = 1
for fc in fcs:
    print(fc)
    out_fc = "SP_Test" + str(counter)
    arcpy.SpatialJoin_analysis(fc1,fc, out_fc, "JOIN_ONE_TO_ONE", "KEEP_ALL", "", "INTERSECT", "", "")
    counter += 1
    fc1 = out_fc
print ('(Elapsed time: ' + str(d.now() - startTime)[:-3] + ')')

2Quiker

ok, I see what the "Counter" is doing (SP_Test1, SP_Test2, SP_Test3). The Fc1 layer has 95K features after merging all there is 361,320k features in the merged layer. I must be missing something?

Anonymous User

Here is the code that handles the spatial joining process:

def fieldmapping(targetFc, joinFc, keepFields):
    """
    Support function create a field map between two featureclasses
    """
    fldMap = arcpy.FieldMappings()
    # Creating field maps for the two files
    fldMap.addTable(targetFc)
    fldMap.addTable(joinFc)

    # Adds the fields from the targetFC to the list of fields that you want to keep
    keepFields.extend(i.name for i in arcpy.ListFields(targetFc))

    # Removing unwanted fields
    for field in fldMap.fields:
        if field.name not in keepFields or field.name in ['TARGET_FID', 'Join_Count']:
            fldMap.removeFieldMap(fldMap.findFieldMapIndex(field.name))

    return fldMap
        
# ----------------- 
def joinsubdivision():
    # List of fields that you want to keep from the input featureclass 
    subd = ['sub_code', 'blk_tr_no', 'lot_no', 'sub_name', 'plat_year', 'plat_month']
    # helper function to map the target featureclass to the input featureclasses fields that you want to keep.
    subdMap = fieldmapping(updatedAddress, subdivision, subd)
    lyrInMem = arcpy.SpatialJoin_analysis(inLayer, subdivision, os.path.join(in_memory, 'subdiv'), 'JOIN_ONE_TO_ONE', 'KEEP_ALL', subdMap, 'INTERSECT')
     # If needed, alter field names (alternative to having to create a specific field map for the field)
    altfields(lyrInMem, [('blk_tr_no', 'blk_tr')])

Hope this gives you an idea to get started and if you have any questions please feel free to ask.

2Quiker

When you say 'all joins' do you meaning making 20+ separate spatial joins or

arcpy.SpatialJoin_analysis(fc1,fc, "memory\SP_Test" + str(counter), "JOIN_ONE_TO_ONE", "KEEP_ALL", "", "INTERSECT", "", "")

2Quiker

If you don't mind sharing that would be great. I hate and think the way ESRI's traditional field mapping is so counterintuitive.

blake.terhune

Spatial join is an inherently intensive operation and there's not much you can do to improve that. Indexes, multithreading, maybe using memory workspace. The looping is how I would go about doing it. You could consider using field mappings when you copy the data to memory and keep only the necessary fields but that's a lot more work because you have to have a list of keep fields for every feature class.

Also, looking at your code, it looks like you'll end up spatial joining fc1 to its self when it gets to that in the fcs list.

DavidPike

it looks like SP_test is being overwritten each time. I would just produce all the joins in memory then merge.

import arcpy, os  
from arcpy import env
from datetime import datetime as d
startTime = d.now()

arcpy.env.workspace = r"C:\Temp\Default2.gdb"  
fc1 = r"C:\Temp\Default2.gdb\Parcels"

arcpy.env.overwriteOutput = True
fcs = arcpy.ListFeatureClasses() #gets the feature classes
counter = 1
for fc in fcs:
    print(fc)
    arcpy.SpatialJoin_analysis(fc1,fc, "SP_Test" + str(counter), "JOIN_ONE_TO_ONE", "KEEP_ALL", "", "INTERSECT", "", "")
    counter += 1
print ('(Elapsed time: ' + str(d.now() - startTime)[:-3] + ')')

Anonymous User

If you just need certain fields from the other featureclasses, I have a few lines of code that I use in a process to field map the desired fields and skip the rest to keep the target fc slim throughout 14 joins in a loop. If you want to see it as an idea I can post it.

Join the Conversation

Signed in members can post, follow updates, and more. New here? Register a free account.

Quick Links

Community Help
Find useful guides, FAQs, and documents to help you navigate and make the most of Esri Community.