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.

I have been trying for some time now to create a script to extract the polygons within an extent, but am unable to get the same result as when I use the user interface manual method.

Maybe someone can see the fault of my script, I have stripped it down to the minimum for simplicity.

import arcpy

layer = arcpy.GetParameterAsText(0)
top = arcpy.GetParameterAsText(1)
left = arcpy.GetParameterAsText(2)
bottom = arcpy.GetParameterAsText(3)
right = arcpy.GetParameterAsText(4)
outputFeature = arcpy.GetParameterAsText(5)

# Determine if Layer is a Layer
desc = arcpy.Describe(layer)
if hasattr(desc, "layer"):
     lyr = layer
     arcpy.AddMessage('Parameter was layer: {0}'.format(lyr))
else:
     lyr = arcpy.MakeFeatureLayer_management(layer,"in_memory\\source_layer")
     arcpy.AddMessage('Parameter made to layer: {0}'.format(lyr))

sr = arcpy.CreateSpatialReference_management("",lyr,"","","","")

lowerLeft = arcpy.Point(left,bottom)
upperLeft = arcpy.Point(left,top)
upperRight = arcpy.Point(right,top)
lowerRight = arcpy.Point(right,bottom)

a = arcpy.Array([lowerLeft,lowerRight,upperRight,upperLeft,lowerLeft])
thepoly = arcpy.Polygon(a,sr)

arcpy.AddMessage('Selecting form lyr by polygon: {0}'.format(thepoly))
selection = arcpy.SelectLayerByLocation_management(lyr, "WITHIN", thepoly, "", "NEW_SELECTION")

matchcount = int(arcpy.GetCount_management(selection).getOutput(0)) 
if matchcount == 0:
    arcpy.AddMessage('no features matched spatial and attribute criteria')
else:
    if matchcount < 1000:
        arcpy.CopyFeatures_management(selection, outputFeature)
        arcpy.AddMessage('{0} objects that matched criteria written to {1}'.format(matchcount, outputFeature))

Thanks for your help.

Follow on after some responses

Thank you for your answer, and I would have expected you to be correct. However what ever I try (and I've tried lots of different combinations) I don't get the desired results. Unless I use the GUI interface to SelectByLocation

Screen shot of GUI, just to clarify

The selected features produced by the GUI interface is exactly what I would expect. I've change the code a bit to output thepoly as feature class (BoundingAOI) so I could check that input parameters where correctly used.

import arcpy

try:

    layer = arcpy.GetParameterAsText(0)
    top = arcpy.GetParameter(1)
    left = arcpy.GetParameter(2)
    bottom = arcpy.GetParameter(3)
    right = arcpy.GetParameter(4)
    outputFeature = arcpy.GetParameterAsText(5)

    # Determine if Layer is a Layer
    desc = arcpy.Describe(layer)
    if hasattr(desc, "layer"):
         lyr = layer
         arcpy.AddMessage('Parameter was layer: {0}'.format(lyr))
    else:
         lyr = arcpy.MakeFeatureLayer_management(layer,"in_memory\\source_layer")
         arcpy.AddMessage('Parameter made to layer: {0}'.format(lyr))

    sr = arcpy.CreateSpatialReference_management("",lyr,"","","","")

    lowerLeft = arcpy.Point(left,bottom)
    upperLeft = arcpy.Point(left,top)
    upperRight = arcpy.Point(right,top)
    lowerRight = arcpy.Point(right,bottom)

    a = arcpy.Array([lowerLeft,lowerRight,upperRight,upperLeft,lowerLeft])
    thepoly = arcpy.Polygon(a)

    outFeatureClass = arcpy.CreateFeatureclass_management(arcpy.env.workspace, "BoundingAOI", "Polygon", "#","#","#", sr)
    arcpy.AddMessage('Created {0}'.format(outFeatureClass))
    cur = arcpy.da.InsertCursor(outFeatureClass, ["SHAPE@"])
    cur.insertRow([thepoly])

    arcpy.AddMessage('Selecting form lyr by polygon: {0}'.format(thepoly))
    selection = arcpy.SelectLayerByLocation_management(lyr, "COMPLETELY_WITHIN", outFeatureClass, "", "NEW_SELECTION")

    matchcount = int(arcpy.GetCount_management(selection).getOutput(0)) 
    if matchcount == 0:
        arcpy.AddMessage('no features matched spatial and attribute criteria')
    else:
        if matchcount < 1000:
            arcpy.CopyFeatures_management(selection, outputFeature)
            arcpy.AddMessage('{0} objects that matched criteria written to {1}'.format(matchcount, outputFeature))

except Exception, e:
    import traceback
    map(arcpy.AddError, traceback.format_exc().split("\n"))
    arcpy.AddError(str(e))

I have to admit that I'm new to ArcGIS and Python (background in Fortran, C, perl and C#), so thank you for the exception handling hint.

Is there any way of step-by-step debuging arcpy scripts?

Thanks again for your responce. Edit:

Thanks for all the comments and suggestions, so far:-

  • Have applied service pack 1 to ArcGIS 10.1, so am now running 10.1sp1
  • Got PyScripter running and able to debug from ArcMap (did need to copy the contents of the rpyc.zip file from the PyScripter install to a rpyc directory under ArcGIS installation of Python C:\Python27\ArcGIS10.1\lib\site-packages\rpyc on my system), this maybe of use to someone
  • Cast the double input parameters as double as suggested
  • No errors from the script
  • Below is a screen shot of the script GUI as requested

result pane and script GUI

My input parameters are (on the comand line)

"D:\Current\PaddockGRASP VS2012\PaddockGRASP Rev 1\PaddockGRASP.gdb\QLD_CADASTRE_DCDB_20" -18.933367 145.278049 -19.735903 146.119615 "D:\Current\PaddockGRASP VS2012\PaddockGRASP Rev 1\PaddockGRASP.gdb\OutputFeatures"

Edit:

It appears after all this that the input dataset had some sort of corruption (ArcSDE can't find anything wrong, not sure why or what is actually wrong), found another older version of the dataset and the problem went away. Frustrating, educational, and hopefully useful to someone else too. Thanks to everyone for the help.

share|improve this question
2  
what do you mean by being unable to get the same result? does the script crash or selections are different? If it crashes, what error do you get? –  jareks Feb 8 '13 at 13:01
    
I'm refering to the fact that I can do the selection using the GUI interface of SelectByLocation in ArcMap but not using the arcpy script. –  David Ahrens Feb 10 '13 at 23:57
    
@David, it is important that you describe what went wrong, especially if there are any error messages. –  blah238 Feb 11 '13 at 1:18
1  
It would be helpful to see what parameters you are specifying in the "Add Script" GUI. –  Aaron Feb 11 '13 at 1:37
    
There are no errors, just wrong results, I'll put the script GUI in the orginal message –  David Ahrens Feb 11 '13 at 4:43

2 Answers 2

The problem is that these variables (left, top, bottom, right):

lowerLeft = arcpy.Point(left,bottom)
upperLeft = arcpy.Point(left,top)
upperRight = arcpy.Point(right,top)
lowerRight = arcpy.Point(right,bottom)

need to be numeric not text. Surround your code with a try/except error handeling will help you figure out potential coding errors e.g.

try:

  # your code in here

except Exception, e:
  import traceback
  map(arcpy.AddError, traceback.format_exc().split("\n"))
  arcpy.AddError(str(e))
share|improve this answer
1  
I disagree with the huge try/except, more often than not it just ends up obscuring errors. Let them bubble up normally. –  blah238 Feb 11 '13 at 0:19
    
You can pass off a numeric value to an object using, for example: right = arcpy.GetParameterAsText(0). I don't think that is the issue. Rather, I suspect the problem is with setting the parameters within the "Add Script" GUI. –  Aaron Feb 11 '13 at 1:35

I think @artwork21's suggestion was to try casting your left/top/bottom/right parameters to numbers, e.g.:

top = float(arcpy.GetParameterAsText(1))
left = float(arcpy.GetParameterAsText(2))
bottom = float(arcpy.GetParameterAsText(3))
right = float(arcpy.GetParameterAsText(4))

However, this is not necessary according to my testing (see below).

I would suggest stepping through your code in a Python debugger, such as PyScripter's. You can supply your arguments as command line parameters or hardcode them for testing purposes:

Command Line Parameters menu item Command Line Parameters dialog

More on debugging tools here: How can I effectively debug arcpy Python scripts?

Edit: I tested your original script as a script tool, and it works fine for me. How have you set up your script tool parameters? This works for me (all but the last are Inputs, the last is an Output parameter). I am on 10.1 SP1, what version are you at?

Script tool parameters dialog

Also tested hardcoded parameter values as a standalone script in PyScripter:

layer = r"C:\GISData\12345.shp"
top = "394605.709522"
left = "13411007.111189"
bottom = "392946.646919"
right = "13413137.265642"
outputFeature = r"C:\GISData\12345_copyFeaturesByRectangle.shp"
share|improve this answer
    
Thank you for the links to PyScripter, which I've now got working. Just checked my ArcGIS version and found it at 10.1, in the process on appling the service pack, so we're on the same version. Get back to you with the results –  David Ahrens Feb 11 '13 at 4:38
    
Thank you for testing the script. –  David Ahrens Feb 11 '13 at 6:30
    
Thank you for testing the script. Can't for the life of me figure out what is going on. Just tried changing my input and output layers to shapefiles (as you have), same result. I figure cause the AOIBoundary layer which I'm currently outputing can be inported and used in ArcMap the problem is not there, and it works for you. One very frustrated indivdual here... Could ArcMap installation be broken (doesn't make any sence in my books, cause the GUI method does work). The debugger is a bit difficult cause it only presents Result objects, but the step thru does as expected. Thanks for your time. –  David Ahrens Feb 11 '13 at 6:39

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.