Take the 2-minute tour ×
Geographic Information Systems Stack Exchange is a question and answer site for cartographers, geographers and GIS professionals. It's 100% free, no registration required.

Basically it's to perform a one-to-many join in a python script, from a dbf to a shapefile. I saw (Greg Krakow's method) using QGIS, but I can't use it (I'm working on it, but QGIS is not yet recognize as an "official and usable" option within our organisation).

My dbf contains many records with the same ID and my source shapefile contains the geometry for each ID.

So, my approach is to create from scratch a shapefile ("outFile" in the script) by copying the structure of my dbf (adding field by field), using arcpy.FieldList()

dbfFieldList = arcpy.ListFields(dbf)

for field in dbfFieldList:
    fieldTypeShp = validFieldType(field.type)
    arcpy.AddField_management (outFile, field.name, fieldTypeShp, "", "", field.length)

validFieldType(field.type) is a function "converting" the field type of my dbf to a valid field type for arcpy. For example, I have a "smallint" which came back as "SHORT" so I could use AddField_management...

So I figured that if I could read my dbf line by line, I could, field by field (since both my dbf and shp have the exact same structure), udpate the value of my shapefile.

#-- Seach cursor for my .dbf
searchCur = arcpy.SearchCursor(dbf)

for sCur in searchCur:

    #-- Insert cursor for my shapefile
    cur = arcpy.InsertCursor(outFile)
    row = cur.newRow()

    #-- Reading for the current row, the fields value, one by one.
    for field in dbfFieldList:
        fieldName = field.name
        val = sCur.getValue(fieldName)
        row.setValue(fieldName,val)

    #-- Inserting the new row into my shapefile
    cur.InsertRow(row) 

But I have an error executing the function "setValue". If I remove that line and the cur.InsertRow (or put a "print val" it works perfectly. I have created a shapefile with the same number of row as my dbf, but I can't add the field values.

After that part, I would use the ID field to retrieve the geometry from my source shapefile with the same ID, recreating my 1-to-many join.

If only I could make that part work...

Any ideas?

Thanks!

share|improve this question
1  
What is the field name, field type, and value of the line/field that's bombing your script? What is the error that's being thrown? Just to help efficiency of the script I would create the insert cursor outside of your search cursor loop; then at the end you'll want a "del cur,searchCur" to delete the insert and search cursor objects. –  Jason Oct 17 '12 at 19:23
    
Hard to know if it's the value/field or it's a scripting error? The name is OID, type = OID and value = 0... –  fgcartographix Oct 17 '12 at 19:59
    
Can you post a sample of your data? I have some ideas that may work but I'll need to test them when I get into the office tomorrow. –  sgrieve Oct 17 '12 at 20:07
4  
I think the OID field in a shapefile/feature class is an autonumber field maintained by ArcGIS...I don't think it will accept a value from the insert cursor. If you want to keep the OID from your source dbf table for reference you would need to create a new text/integer/float field to hold that piece of data. –  Jason Oct 17 '12 at 20:21
    
I've added a «if field.type <> "OID":» and it didn't crash, but it didn't add the values to the shapefile. Oh! I deleted my insertRow... brb ;) –  fgcartographix Oct 18 '12 at 11:41
show 4 more comments

3 Answers

up vote 3 down vote accepted

You should be able to overcome the problem by adding in some error handling. Without the exact error message, I can only assume that you are running into issues with attempting to manually fill restricted fields such as OID, as @Jason pointed out. Try the following and see if this works out for you:

dbfFieldList = arcpy.ListFields(dbf)

for field in dbfFieldList:
    fieldTypeShp = validFieldType(field.type)
    try:
        arcpy.AddField_management (outFile, field.name, fieldTypeShp, "", "", field.length)
    except:
        print "Field, '" + field.name + "' cannot be added."

searchCur = arcpy.SearchCursor(dbf)
cur = arcpy.InsertCursor(outFile)

for sCur in searchCur:


    #-- Create a new row in the outfile
    row = cur.newRow()

    #-- Reading for the current row, the fields value, one by one.
    for field in dbfFieldList:
        val = sCur.getValue(field.name)
        try:
            row.setValue(field.name,val)
        except:
            print "Cannot populate field, '" + field.name + "'."

    #-- Inserting the new row into my shapefile
    cur.insertRow(row) 
share|improve this answer
    
That pretty much solve it! –  fgcartographix Oct 18 '12 at 13:08
add comment

Assuming you want to create duplicate geometries for those features in the shapefile for which there exist multiple rows having the same ID in the DBF table (one to many relationship), you could do this using Make Query Table, followed by Copy Features to make a permanent feature class or shapefile.

See the example 2 of Examples of queries with the Make Query Table tool for a relevant example:

Example 2

You could implement this pretty easily in ModelBuilder.

You could also do this once interactively using ArcToolbox and then right-click each of the results in the Results window and Copy As Python Snippet to get the syntax and paste it into a new Python script (you'll need to add import arcpy at the top and add arcpy. in front of the geoprocessing function calls).

Note: You will need to first copy the shapefile and DBF file into an intermediate geodatabase as Make Query Table only operates on data from an ArcSDE geodatabase, a file geodatabase, a personal geodatabase, or an OLE DB connection, and all input tables must be in the same workspace.

You can accomplish this intermediate step using additional geoprocessing tools as part of your script or model. You might even be able to use the in-memory workspace for this for better performance and to avoid an intermediate file geodatabase, though I wouldn't suggest this if your shapefile or table are very large.

share|improve this answer
    
I'll try that one. I didn't want to create a geodatabase because that script is only a single part of a greater process, but I guess I could create one, use the tool, export to shp and the delete it or keep it as I may need it later. –  fgcartographix Oct 18 '12 at 12:31
    
I don't see how that would be a problem. The script/model itself could (should) manage the lifetime of the intermediate data. –  blah238 Oct 18 '12 at 17:44
add comment

Finally, it worked! Only, it's **ing slow. I have about a record/second (adding all the field values + find the geometry and adding it to my record), but 5807 record to to. So it took about 2h to run. But with your help, I've done a one-to-many join in arcpy so... Here's the code (sorry for the comments in French. I'll come back to edit my answer to translate them when I'll have the chance :

# -*- coding: utf-8 -*-
###############################################################################
from sys import argv
import arcpy
import os, os.path
###############################################################################

#-- Le dbf contient des types de champs qui ne sont pas compatibles avec les shapefiles.
#-- Pour s'en assurer, on prend le type de champ du dbf et on le convertis en qqchose de compatible

def validFieldType(fieldType):
    if fieldType == "OID":
        return "LONG"
    elif fieldType == "String":
        return "TEXT"
    elif fieldType == "Integer":
        return "LONG"
    elif fieldType == "Double":
        return "DOUBLE"
    elif fieldType == "SmallInteger":
        return "SHORT"
    else:
        return "TEXT"   # On peut entrer presque n'importe quoi en texte alors dans le doute...

#-- Reçoit un nom de fichier, un champ dans lequel rechercher et un ID à trouver.
#-- Retourne la géométrie

#script, dbf, shp = argv

dbf = "D:/#PROJETS_AMT/#RABATTEMENT/2012/Donnees/CITLA_2012au/__Competitivite.dbf"
dbfIdFld = "ID"
shp = "D:/#PROJETS_AMT/#RABATTEMENT/2012/Donnees/CITLA_2012au/__CITLA_2012au_lignes.dbf"
shpIdFld = "To_ID"
arcpy.env.workspace = "D:/#PROJETS_AMT/#RABATTEMENT/2012/Donnees/CITLA_2012au/"

#-- On permet l'overwrite de l'output pour que si ça plante, on ait juste à repartir
#-- le script et que ça réécrire par-dessus les fichiers déjà créés.
arcpy.env.overwriteOutput


#-- On fait la liste des champs du dbf qui va nous servir à recréer un shapefile identique.
dbfFieldList = arcpy.ListFields(dbf)

#-- On va chercher la géométrie du shapefile pour savoir quel type de géométrie créer
descSHP = arcpy.Describe(shp)
descDBF = arcpy.Describe(dbf)
dbfName = descDBF.name[0:-4]
shpName = descSHP.name[0:-4]



#-- J'aime bien mettre des "__" au début de mes outputs pour les retrouver facilement
#-- dans mes dossiers, mais en mergeant des noms de fichiers, parfois ça fait trop long
#-- C'est juste pour éviter d'avoir des "___" dans mes noms de fichiers.
outFile = dbfName + "_" + shpName + ".shp"
outFile = outFile.replace("___","__")

#-- On permet l'overwrite de l'output pour que si ça plante, on ait juste à repartir
#-- le script et que ça réécrire par-dessus les fichiers déjà créés.
listOutput = arcpy.ListFiles("*.shp")

for output in listOutput:
    if output == outFile:
        print "Fichier existant... suppression!"
        arcpy.Delete_management(outFile)


#-- Création du shapefile de sortie
arcpy.CreateFeatureclass_management(arcpy.env.workspace, outFile, descSHP.shapeType, "", "", "", descSHP.spatialReference, "", "", "", "")

#-- On rajoute les mêmes fichiers que le dbf
for field in dbfFieldList:
    if field.name == "ID":
        arcpy.DeleteField_management(outFile, "ID") #-- Un shp "vierge" contient déjà un champ ID. On veut le remplacer par le nôtre.
    fieldTypeShp = validFieldType(field.type)
    arcpy.AddField_management (outFile, field.name, fieldTypeShp, "", "", field.length, "", "", "", "")


#-- On ajoute dans le fichier le même
searchCur = arcpy.SearchCursor(dbf)

#-- Pour suivre la progression
total = arcpy.GetCount_management(dbf)
i = 1

for sCur in searchCur:

    print "Traitement de la ligne %i/%i" %(i, total)

    #-- On ajoute une nouvelle ligne dans le shapefile de sortie
    cur = arcpy.InsertCursor(outFile)
    row = cur.newRow()

    #-- On passe champ par champ pour copier les valeurs

    for field in dbfFieldList:
            fieldName = field.name
            val = sCur.getValue(fieldName)
            try:
                row.setValue(field.name,val)
            except:
                #print "  Cannot populate field %s. Moving on!" %(str(field.name))

    #-- Maintenant on cherche dans le fichier à joindre la géométrie correspondant à l'ID.
    srcCur = arcpy.SearchCursor(shp)

    for src in srcCur:
        if src.getValue(shpIdFld) == sCur.getValue(dbfIdFld):
            row.Shape = src.Shape

    cur.insertRow(row)    
    del cur, row, src, srcCur

    i+=1        
del sCur, searchCur

There's still some hardcoded value, but I'm sure you'll forgive me.

Now, I guess I could have done it with a LEFT JOIN in PostGIS? Would have been easier, no?

Thanks!!

share|improve this answer
1  
That doesn't really surprise me. Cursors have been notoriously slow, which is why at 10.1 they made faster cursors in the arcpy.da module. Still, I suggest using geoprocessing tools where they can do the same job, probably do it faster, and do it with a lot less code. –  blah238 Oct 18 '12 at 17:48
    
Try creating your insert cursor, once, outside of your search cursor loop - you will still call row = cur.newRow() for each record, but this way you should be able to reduce your run time. –  Jason Oct 18 '12 at 18:08
    
I tried to put «cur = arcpy.InsertCursor(outFile)» before «for sCur in searchCur:» but I got an error : «cur is not defined» ?? –  fgcartographix Oct 18 '12 at 19:11
1  
For the record, I imported my dbf and shp into PostGIS and with a single query "CREATE TABLE test AS SELECT compet.*, lignes.geom AS geom FROM compet LEFT JOIN lignes ON compet.id = lignes.to_id;", got my result in 181 ms... –  fgcartographix Oct 18 '12 at 19:13
1  
@Jason Yes. Besides, my error message pointed me to the line where I was calling cur for my new row. Didn't search long after that. My computer must had something yesterday... Recreate the query in PostGIS in 170 ms this morning ;) I know I said in my intro that QGIS is not an option, but I said to myself "scr*w them", that what's working better ;) –  fgcartographix Oct 19 '12 at 13:45
show 1 more comment

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.