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 over 100 features in a shapefile and I want to create individual viewsheds for each and then run processes on each of those resulting viewsheds in ArcGIS using Python.

  • Is it possible to write a script using a cursor to run through each "row" and use that individual observer point to run the viewshed? Then update the cursor to the next row and rerun it until all features have been exhausted?

Or would I have to first Split the shapefile into over 100 files and run each individually?

  • Also, is it possible to designate the filenames for the resulting viewsheds (or shapefiles) using a specific field?
share|improve this question
2  
Instead of running 100 operations you could reduce this to about seven using the Observer Points Tool, which will perform the viewshed calculation for up to 16 observers at a time. There's a chance it could be a little more efficient. –  whuber Apr 22 '13 at 17:18

3 Answers 3

up vote 5 down vote accepted

I've never used the viewshed tool so I can't speak to the specifics of using that tool, however in regards to your suggestion in this question I will add an answer here to build on KHibma's solution.

The first thing we need to address is the format of your SQL query: "[OID] = count". As it stands, this query will fail since "count", as represented here, is a string and strings need to be surrounded by single quotes: "[OID] = 'count'". Now, the second reason this query will fail is that OID is an integer field and as such there will be no OID value equal to "count". This is where string formatting comes into play: '"OID" = {}'.format(i) or '"OID" = {0}'.format(i) for those using Python 2.6.x (Use quotation marks to enclose the fieldname when using a shapefile. See ArcGIS help.).

To avoid certain pitfalls of making queries on the OBJECTID field, we can nest the loop into a cursor to make sure that we are only querying OID values that exist:

arcpy.MakeFeatureLayer_management ("C:/data/pts.shp", "pts")

with arcpy.da.SearchCursor('pts',['OID']) as cursor:
    for row in cursor:
        oid = row.OID # or whatever your OBJECTID field is called
        arcpy.SelectLayerByAttribute_management ("pts", "NEW_SELECTION", '"OID" = {}'.format(oid))
        outViewshed = Viewshed("elevation","pts",2,"CURVED_EARTH",0.15)
        outViewshed.save("C:/sapyexamples/output/outvwshd"+oid)
share|improve this answer
    
I am now getting the error: AttributeError: 'da.SearchCursor' object has no attribute 'OID'. It breaks at "oid = cursor.OID... Thoughts? Thanks for everything btw! –  Spatial Pariah Apr 22 '13 at 21:36
    
You'll need to determine the actual name of you object id field and substitute that for ‛OID‛. –  Jason Apr 23 '13 at 2:27

OK, borrowing from the answer from @Jason, I was able to tweak it a bit to finally come up with a script that works great. I also switched it to use FID instead of OID. Either way, it will work.

arcpy.MakeFeatureLayer_management (inPoints, "pts")

with arcpy.da.SearchCursor('pts',['FID']) as cursor:
    for row in cursor:
        fid = str(row[0]) ### Must be string
        print fid
        arcpy.SelectLayerByAttribute_management ("pts", "NEW_SELECTION", '"FID" = {}'.format(fid))
        outViewshed = Viewshed(inDEM,"pts",1,"CURVED_EARTH",0.15)
        outViewshed.save("C:/Python/output/viewsheds/"+fid)

Perhaps this will help someone in the future. Thanks everyone!

share|improve this answer
    
Good one. I spaced out on the row iterator! Glad you got it worked out. –  Jason Apr 23 '13 at 2:29

You could do a cursor, but its more straight forward to use the select by attributes tool. Simply put it in a loop ... Here's some untested code (you'll probably have to tinker with the expression in the select: (note I just add the count to the output name, if you want to get the specific value from a field then you'll have to open a cursor)

arcpy.MakeFeatureLayer_management ("C:/data/pts.shp", "pts")
count = int(arcpy.GetCount_management("pts").getOutput(0)) 
i=0
while i < count:
   arcpy.SelectLayerByAttribute_management ("pts", "NEW_SELECTION", "[OID] = i")
   outViewshed = Viewshed("elevation","pts",2,"CURVED_EARTH",0.15)
   outViewshed.save("C:/sapyexamples/output/outvwshd"+count)
   i+=1

Or you could do it all in a model. The "advanced" model at the bottom of this topic explains how to do it. Just cut out the tabulate area part, and do just the viewsheds with the iterator. Using the model only tool "get field value" you can use inline variables to update the final output name. This would actually be my preferred way to tackle this problem (model over script)

share|improve this answer
    
I am betting this will work just fine. Thanks! –  Spatial Pariah Apr 22 '13 at 17:41
    
I am getting a syntax error at "for i < count:" and the "<" is highlighted? I have closed all other arguments correctly, so I am not sure what's going on here. –  Spatial Pariah Apr 22 '13 at 17:55
    
doh me - changed "for" to "while" –  KHibma Apr 22 '13 at 17:56
    
That helped, but now I get an invalid expression. I changed [OID] to 'OID' since it seems you use quotes to call up featurelayers. I had it print count immediately following the initial "while" statement. I get a number of 57? So should 'OID' actually equal "i" instead? –  Spatial Pariah Apr 22 '13 at 18:12

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.