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'm trying to run a script that will create cost paths for a large number of randomly generated points then assign the cost to those points. But my script is beginning to create cost path rasters that don't contain a path. I'm pretty new at python scripting, and the code does need to be cleaned up a bit and there are a few print statements left over from where I was trying to find bugs, but my main concern is getting the cost paths created, then getting the total cost into the points. Every once in awhile the script will work completetly without error, but it takes 3-4 hours to run on my computer, so running it 5 times to get it to work isn't an option.

# Import system modules for coupling with ArcGIS, create geoprocessing link,
# acquire licenses, and load the necessary toolboxes
# copied from lab 4 script by Dr. Aldrich

import sys, string, os, arcgisscripting, arcpy
arcpy.CheckOutExtension("Spatial")


#set the work location for the script

arcpy.env.workspace = "C:\Users\Ryan\Desktop\GIS\python\project_test\project_test_gdb.gdb"

#script will overwrite existing data

arcpy.env.overwriteOutput = 1

#classes
##add to this as necessary

class randompoint:
    pass





#variables
##add to this as necessary
dem1 = "grdn36w089_13"
dest1 = "destination_one"
pointarea = "random_point_polygon"
point_names = []
costpath_names = []


#definitions
##add to this as necessary




#creates 100 random points, names each point sequencially, and adds that name to the list "point_names"
def create_points ():
    pointcount = 0
    while pointcount < 100:
        pointcount = pointcount + 1
        print "making point " + str(pointcount)
        pointname = "random_point" + str(pointcount)
        point_names.append(pointname)
        arcpy.CreateRandomPoints_management ("", pointname, pointarea, "", 1, "", "")

#creates cost paths for points, names them sequencially, and adds that name to the list "costpath_names"        
def create_costpaths ():
    for point in point_names:
        print "creating costpath for point " + point
        cpname = "costpath_for_" + point
        costpath = arcpy.sa.CostPath (point, "cost_for_first_destination", "back_link_for_first_destination", "", "")
        costpath.save ("costpath_for_" + point)
        costpath_names.append (cpname)

        #failed ideas
        #costpath = None 

#assigns the values from costpaths to their points using lists
#eventually replace the lists with classes that contain the name of the point and costpath
def assignvalues ():

    checklists ()

    counter = 0

    print "assigning values"
    while counter <= len(point_names) - 1:

        print counter

        pointname = point_names[counter]

        costpathname = costpath_names[counter]

        print pointname
        print costpathname

        arcpy.AddField_management (pointname, "cost1", "DOUBLE", "", "", "", "", "", "", "")

        findcursor = arcpy.SearchCursor (costpathname, "", "", "PATHCOST", "")

        replacecursor = arcpy.UpdateCursor (pointname, "", "", "cost1", "")


        for item in findcursor:
            newvalue = item.getValue ("PATHCOST")
            if newvalue != 0:
                print newvalue
                for row in replacecursor:
                    print newvalue
                    print row.getValue ("cost1")
                    row.setValue ("cost1", newvalue)
                    print row.getValue ("cost1")
                    replacecursor.updateRow(row) 

        counter += 1

#checks the lists to make sure they are the same length
def checklists ():
    print "checking lists"
    if len(point_names) == len(costpath_names):
        print "lists ok"
    else:
        raw_input ("Error: lists do not match, check script.  Press any key to continue")


#script

#generate slope raster

print "generating slope raster"
slope = arcpy.sa.Slope (dem1, "PERCENT_RISE", "")  

slope.save ("slope_for_cost")  


#generate cost distance raster

print "generating cost distance raster"
cost1 = arcpy.sa.CostDistance (dest1, slope, "", "")

cost1.save ("cost_for_first_destination") 

#generate back link raster

print "generating backlink raster"
backlink1 = arcpy.sa.CostBackLink (dest1, cost1, "", "")  

backlink1.save ("back_link_for_first_destination") 


create_points()

create_costpaths ()

assignvalues ()


raw_input('press enter to close')
share|improve this question
1  
I doubt that using a slope raster in a CostDistance calculation will do what you expect, because it does not account for the direction in which each cell is crossed: the "impedance" will be the same regardless of whether the cell is crossed going uphill, downhill, or sideways. It is unusual for that to be an accurate model of anything. –  whuber Jun 25 '13 at 18:59
    
It's more of a test of the process. I needed something to generate a cost distance raster from, and slope was simple enough. Eventually it will use a number of things to determine cost. –  ryan Jun 25 '13 at 19:20
add comment

Your Answer

 
discard

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

Browse other questions tagged or ask your own question.