Stand alone script fails using spatial analyst tools

ElizabethLewis · ‎03-05-2012

Hi,

I'm trying to create a stand alone script to automate watershed analysis. So far it works up until it gets to the spatial analyst
tools but then fails. Can anyone show me where I have gone wrong please? Here is my code:

import arcpy
import os
import time
import sys
from arcpy import env
from arcpy.sa import *

env.workspace = "H:/PhD/PythonLessons/Scripts/testIt"
catchmentName = "CatchmentName"
catchmentRaster = catchmentName + ".tif"
catchmentRasterProj = catchmentName + "proj.tif"
catchment100 = catchmentName + "100.tif"
catchment100Fill = catchmentName + "100Fill.tif"
catchment100Flow = catchmentName + "100Flow.tif"
catchment100Acc = catchmentName + "100Acc.tif"
catchment100PP = catchmentName + "100PP.tif"
catchment100WS = catchmentName + "100WS.tif"
folder = "H:/PhD/PythonLessons/Scripts/testIt"
folderAdd = "H:/PhD/PythonLessons/Scripts/testIt/"
outlet = "outlet.shp"
projection = "C:/Program Files (x86)/ArcGIS/Desktop10.0/Coordinate Systems/Projected Coordinate Systems/National Grids/Europe/British National Grid.prj"
# identify guage
#create pourpoint shapefile
# identify catchment outline
# select appropriate testIt

# have a folder with all of the little ascii testIt for a catchment. This folder will be called *catchment name*
print "checking extensions"

try:
 arcpy.CheckExtension("Spatial") == "Available"
 arcpy.CheckOutExtension("Spatial")
 print "done"
except:
 print "piddle"
 

print "converting the ascii files to raster files"

asciiList = os.listdir(folder)    

for x in range(len(asciiList)-1, -1, -1):
 if not asciiList.endswith(".asc"):
   asciiList.pop(x)


for file in asciiList:  
 print folderAdd + file
 print folderAdd + file[:-4]
 
 try:
  arcpy.ASCIIToRaster_conversion(folderAdd+ file, folderAdd+ file[:-4] + ".tif", "INTEGER")
 except:
  print "Processing " + file + "FAILED"
  
 time.sleep(3)
 
tifList = os.listdir(folder)

print tifList

for x in range(len(tifList) -1, -1, -1):
    if not tifList.endswith(".tif"):
        tifList.pop(x)

print tifList

tifString = ""
for tif in tifList:
 tifString += tif + ";"
 
tifString = tifString[:-1]

print tifString

print "mosaic to new raster"

try:
 arcpy.MosaicToNewRaster_management(tifString, folderAdd, catchmentRaster, "#", "16_BIT_SIGNED", "#", "1" , "MEAN")
except:
 print "Error with mosaic to new raster"
 
time.sleep(3)


print "Giving it the right coordinates"
try:
 arcpy.DefineProjection_management(folderAdd + catchmentRaster, projection)
 arcpy.ProjectRaster_management(folderAdd + catchmentRaster, folderAdd + catchmentRasterProj, projection, "BILINEAR", "#", "#","#",projection)
except: 
 print "Error giving it the right coordinates"

time.sleep(3)
 
print "Changing resolution"

try:
 outAggreg = Aggregate(catchmentRasterProj, 10, "MEDIAN", "TRUNCATE", "DATA")
 outAggreg.save(folderAdd + catchment100)
except:
 print "Error changing resolution"
 
time.sleep(3)
 
print "Filling sinks"

try:
 outFill = Fill(catchment100)
 outFill.save(folderAdd + catchment100Fill)
except:
 print "error filling sinks"
 
print "Calculating fow direction"

time.sleep(3)

try:
 outFlowDirection = FlowDirection(catchment100Fill, "NORMAL")
 outFlowDirection.save(folderAdd + catchment100Flow)
except:
 print "Error calculating flow direction"
 
 
print "Calculating flow accumulation"

time.sleep(3)

try:
 outFlowAccumulation = FlowAccumulation(catchment100Flow)
 outFlowAccumulation.save(folderAdd + catchment100Acc)
except:
 print "Error calculating flow accumulation"

print "Creating pour point snap"

time.sleep(3)

try:
 outSnapPour = SnapPourPoint(outlet, catchment100Acc, 5, "FID")
 outSnapPour.save(folderAdd + catchment100PP)
except:
 print "Error snapping pour point"
 
time.sleep(3)

print "Creating watershed"

try:
 outWatershed = Watershed(catchment100Flow, catchment100PP)
 outWatershed.save(folderAdd + catchment100WS)
except:
 print "Error creating watershed"
 
time.sleep(3)
 
print "Creating files for SHETRAN prepare"

try:
 arcpy.RasterToASCII_conversion(catchment100, folderAdd + catchmentName + "100DEM.txt")
 arcpy.RasterToASCII_conversion(catchment100WS, folderAdd + catchmentName + "100Mask.txt")
except:
 print "error creating files for SHETRAN prepare"
 
time.sleep(3)
 
print "Done!"

raw_input("close to exit")

Thanks very much

ElizabethLewis · ‎03-07-2012

Ok, I set

env.scratchWorkspace = env.workspace

and tried this set up

try:
        outAggreg = Aggregate(catchmentRasterProj, 10, "MEDIAN", "TRUNCATE", "DATA")
 outAggreg.save(catchment100)
except Exception as e:
    print e.message

for my spatial analyst tools.

It now gives the error 'failed to open raster dataset'. I'm not sure why it won't open it. I checked that it would open in ArcMap 10 and it does that just fine. I've also checked that the file name is right. I checked the license and that I'd imported the right modules.

Also, when I run the, as well as an error it creates a new folder within the workspace called 'aggrega_catc1' but I thought that setting the scratch workspace to the env workspace would stop this?

Thanks again

curtvprice · ‎03-07-2012

It now gives the error 'failed to open raster dataset'. I'm not sure why it won't open it.

To see if it's a path issue, make sure it exists. (If you don't specify a full path, it will look in env.workspace for it.)

if not arcpy.exists(catchmentRasterProj): 
    raise Exception, "%s does not exist" % catchmentRasterProj

In the same manner, the output is also written to the current workspace, so instead of

outFill.save(folderAdd + catchment100Fill)

you can just specify

outFill.save(catchment100Fill)

as well as an error it creates a new folder within the workspace called 'aggrega_catc1' but I thought that setting the scratch workspace to the env workspace would stop this?

This is a temporary raster; all tools write a temporary raster to the scratch workspace when raster tools run. As I said the point of setting the workspaces the same is so when the tool successfully runs this temp raster is renamed (instead of copied) to your output raster when you execute outraster.save.

ElizabethLewis · ‎03-07-2012

I tried arcpy.Exists and it can find the file so it's not my file path that's wrong. Aggregate still won't open the file raising the error 'failed to open raster dataset'.

Thanks for explaining the temporary files.

DorothyMortenson · ‎03-15-2012

I had been trying to do something similar with the SnapPointPoint. I am converting from the gp to arcpy.
What I found was I was trying to create a grid and it didn't like it.
So I created a filegeodatabase for my results to go into. This seems to have solved my issue.
All of my original data is in GRID from arcgis 9.3 I plan to convert it, but I wanted to get the python scripts working first.

curtvprice · ‎03-15-2012

I had been trying to do something similar with the SnapPointPoint. I am converting from the gp to arcpy.
What I found was I was trying to create a grid and it didn't like it.
So I created a filegeodatabase for my results to go into. This seems to have solved my issue.
All of my original data is in GRID from arcgis 9.3 I plan to convert it, but I wanted to get the python scripts working first.

The grid data format is very path sensitive - you need to be very careful to avoid long grid names or names that do not follow the very restrictive grid/coverage naming rules. On the plus side, the grid data format is usually faster to process (though this depends on the tool).

ElizabethLewis · ‎03-22-2012

Hi,

I've finally resolved this issue. It turns out that the arc modules wouldn't run on the H drive that I was using so I saved everything to the C drive instead and it all ran fine.

Here's the code for future reference:

import arcpy
import os
import time
import sys
from arcpy import env
from arcpy.sa import *

env.workspace = "C:/temp/testIt"
env.scratchWorkspace = env.workspace
catchmentName = "CatchmentName"
catchmentRaster = catchmentName + ".tif"
catchmentRasterProj = catchmentName + "proj.tif"
catchment100 = catchmentName + "100.tif"
catchment100Fill = catchmentName + "100Fill.tif"
catchment100Flow = catchmentName + "100Flow.tif"
catchment100Acc = catchmentName + "100Acc.tif"
catchment100PP = catchmentName + "100PP.tif"
catchment100WS = catchmentName + "100WS.tif"
folder = "C:/temp/testIt"
folderAdd = "C:/temp/testIt/"
outlet = "outlet.shp"
projection = "C:/Program Files (x86)/ArcGIS/Desktop10.0/Coordinate Systems/Projected Coordinate Systems/National Grids/Europe/British National Grid.prj"
# identify guage
#create pourpoint shapefile
# identify catchment outline
# select appropriate testIt

# have a folder with all of the little ascii testIt for a catchment. This folder will be called *catchment name*
print "checking extensions"

try:
 arcpy.CheckExtension("Spatial") == "Available"
 arcpy.CheckOutExtension("Spatial")
 print "done"
except:
 print "Didn't work!"
 

print "converting the ascii files to raster files"

asciiList = os.listdir(folder)    

for x in range(len(asciiList)-1, -1, -1):
 if not asciiList.endswith(".asc"):
   asciiList.pop(x)


for file in asciiList:  
 print folderAdd + file
 print folderAdd + file[:-4]
 
 try:
  arcpy.ASCIIToRaster_conversion(folderAdd+ file, folderAdd+ file[:-4] + ".tif", "INTEGER")
 except:
  print "Processing " + file + "FAILED"
  
 time.sleep(3)
 
tifList = os.listdir(folder)

print tifList

for x in range(len(tifList) -1, -1, -1):
    if not tifList.endswith(".tif"):
        tifList.pop(x)

print tifList

tifString = ""
for tif in tifList:
 tifString += tif + ";"
 
tifString = tifString[:-1]

print tifString

print "mosaic to new raster"

try:
 arcpy.MosaicToNewRaster_management(tifString, folderAdd, catchmentRaster, "#", "16_BIT_SIGNED", "#", "1" , "MEAN")
except:
 print "Error with mosaic to new raster"
 
time.sleep(3)


print "Giving it the right coordinates"
try:
 arcpy.DefineProjection_management(catchmentRaster, projection)
 arcpy.ProjectRaster_management(catchmentRaster, catchmentRasterProj, projection, "BILINEAR", "#", "#","#",projection)
except Exception as e:
    print e.message


time.sleep(3)
 
print "Changing resolution"

try:
 outAggreg = Aggregate(folderAdd + catchmentRasterProj, 10, "MEDIAN", "TRUNCATE", "DATA")
 outAggreg.save(catchment100)
except Exception as e:
    print e.message
 
time.sleep(3)
 
print "Filling sinks"

try:
 outFill = Fill(catchment100)
 outFill.save(catchment100Fill)
except Exception as e:
    print e.message
 
print "Calculating fow direction"

time.sleep(3)

try:
 outFlowDirection = FlowDirection(catchment100Fill, "NORMAL")
 outFlowDirection.save(catchment100Flow)
except Exception as e:
    print e.message
 
 
print "Calculating flow accumulation"

time.sleep(3)

try:
 outFlowAccumulation = FlowAccumulation(catchment100Flow)
 outFlowAccumulation.save(catchment100Acc)
except Exception as e:
    print e.message

print "Creating pour point snap"

time.sleep(3)

try:
 outSnapPour = SnapPourPoint(outlet, catchment100Acc, 5, "FID")
 outSnapPour.save(catchment100PP)
except Exception as e:
    print e.message
 
time.sleep(3)

print "Creating watershed"

try:
 outWatershed = Watershed(catchment100Flow,catchment100PP)
 outWatershed.save(catchment100WS)
except Exception as e:
    print e.message
 
time.sleep(3)
 
print "Creating files for SHETRAN prepare"

try:
 arcpy.RasterToASCII_conversion(catchment100, catchmentName + "100DEM.txt")
 arcpy.RasterToASCII_conversion(catchment100WS, catchmentName + "100Mask.txt")
except Exception as e:
    print e.message
 
time.sleep(3)
 
print "Done!"

raw_input("close to exit")

Thanks for the help!

MichaelGleason · ‎06-21-2013

I'm having a very similar problem over a year later -- arcpy.sa modules fail with 999999 errors on the S: drive, but run without issue on C:. I'm experiencing this at ArcGIS 10.0 SP 3 and ArcGIS 10.1 SP 1.

This seems like it should really be a bug report.