gp.calculate field expression syntax - PYTHON

[ThomasBurley]

Hi,
I have a script that just needs slight tweaking (I hope/think). It iterates through a workspace containing rasters, adds a new field called Proportion, then sums the raster value COUNT field in order to then determine the proportion of each raster value type. This is for land cover data - trying to calculate the proportion of land cover type per raster.

I think my syntax in calculate field is not quite right:

import sys
import string
import os

try:
  #for 9.2 and above
  import arcgisscripting
  gp = arcgisscripting.create(9.3)
  gp.AddMessage("\n" + "Using ArcMap 9.2 or above with arcgisscripting..." + "\n")
except:
  #for 9.0/9.1
  import win32com.client
  gp = win32com.client.Dispatch("esriGeoprocessing.GpDispatch.1")
  gp.AddMessage("\n" + "Using ArcMap 9.0/9.1 or above with win32com.client.Dispatch..." + "\n")
#
# Check for the necessary product

Licensed = gp.ProductInfo()
gp.addmessage("You are using the following product " + str(Licensed))

try: 
    gp.SetProduct(str(Licensed))
    gp.addmessage("Successfully set the product type to " + str(Licensed))
except:
    gp.addmessage("Could not find an available ArcGIS License")

# Load required toolboxes...
gp.toolbox = "management"

# Set the Workspace
gp.Workspace = sys.argv[1]
Workspace = str(gp.Workspace)
Workspace = string.replace(Workspace, "\\", "/")

# Set the Output Folder
OutFolder = sys.argv[2]
OutFolder = str(OutFolder)
OutFolder = string.replace(OutFolder, "\\", "/")


LogFile = open(OutFolder + "/" + "log.txt", "w")

try:
 #Generate a list of the grids in the grid folder
 rasters = gp.ListRasters("*", "ALL")
 #Iterate through each raster, add a new field called "Proportion", calculate the proportion of land cover type based on those two fields
 for raster in rasters:
  try: 
   #Add the new Proportion field
   gp.addfield (raster, "PROPORTION", "FLOAT", "2")
   #Use the searchcursor method to sum the values of the count field for calculating the proportions
   x = 0
   rows = gp.SearchCursor(raster)
   row = rows.Next()
   field = "COUNT"
   while row:
    x += row.getvalue(field)
    row = rows.next()
   
   print str(x) + " " + raster
   gp.calculatefield(raster,"PROPORTION","!COUNT! / float(x)","PYTHON")
   
   #rows = gp.SearchCursor(raster)
   #row = rows.Next()
   logfile.write("Just finished processing" + raster + "\n")
  except:
   gp.AddMessage(gp.GetMessages(2))
   ErrorMsg = gp.GetMessages(2)
   print ErrorMsg
   LogFile.write(ErrorMsg + "\n")
  
except:
 gp.AddMessage(gp.GetMessages(2))
 ErrorMsg = gp.GetMessages(2)
 logfile.write(ErrorMsg + "\n")
 

LogFile.close()

Specifically this part of the code:

    gp.calculatefield(raster,"PROPORTION","!COUNT! / float(x)","PYTHON")

I'm storing the sum of the COUNT field in the variable "x", and then want to use that in gp.calculatefield to divide the values in the COUNT field by the sum of the values stored in "x" in memory to determine the proportion of each land cover type


any help would be greatly appreciated - my programming experience is in VBA and VBscript, but I'm trying to break into Python

thanks,

Tom Burley
U.S. Geological Survey
Austin, TX

[KenHartling]

Hi Tom,

When you include the variable "x" within the quotes it will be treated as a string, 'x', and not as a variable.

Try this:

gp.calculatefield(raster,"PROPORTION","!COUNT! / " + x,"PYTHON")

Ken

[ThomasBurley]

Hi Ken,

I tried that but am still not getting success. I checked to verify that the COUNT field is numeric. The script seems to run OK, but I go and check the PROPORTION field and all the values are zero. I tried just assigning the value of the X variable (no calculation) and that worked, there's just something not working with determining the quotient of the COUNT field records based on the total stored in the variable X.


Here's my code again, and attached are the rasters I'm trying in a personal GDB:

import sys
import string
import os

try:
  #for 9.2 and above
  import arcgisscripting
  gp = arcgisscripting.create(9.3)
  gp.AddMessage("\n" + "Using ArcMap 9.2 or above with arcgisscripting..." + "\n")
except:
  #for 9.0/9.1
  import win32com.client
  gp = win32com.client.Dispatch("esriGeoprocessing.GpDispatch.1")
  gp.AddMessage("\n" + "Using ArcMap 9.0/9.1 or above with win32com.client.Dispatch..." + "\n")
#
# Check for the necessary product

Licensed = gp.ProductInfo()
gp.addmessage("You are using the following product " + str(Licensed))

try: 
    gp.SetProduct(str(Licensed))
    gp.addmessage("Successfully set the product type to " + str(Licensed))
except:
    gp.addmessage("Could not find an available ArcGIS License")

# Load required toolboxes...
gp.toolbox = "management"

# Set the Workspace
gp.Workspace = sys.argv[1]
Workspace = str(gp.Workspace)
Workspace = string.replace(Workspace, "\\", "/")

# Set the Output Folder
OutFolder = sys.argv[2]
OutFolder = str(OutFolder)
OutFolder = string.replace(OutFolder, "\\", "/")


LogFile = open(OutFolder + "/" + "log.txt", "w")

try:
 #Generate a list of the grids in the grid folder
 rasters = gp.ListRasters("*", "ALL")
 #Iterate through each raster, add a new field called "Proportion", calculate the proportion of land cover type based on those two fields
 for raster in rasters:
  try: 
   #Add the new Proportion field
   gp.addfield (raster, "PROPORTION", "LONG", "2")
   #Use the searchcursor method to sum the values of the count field for calculating the proportions
   x = 0
   rows = gp.SearchCursor(raster)
   row = rows.Next()
   field = "COUNT"
   while row:
    x += row.getvalue(field)
    row = rows.next()
   
   print str(x) + " " + raster
   gp.calculatefield(raster,"PROPORTION","!COUNT! / " + long(x),"PYTHON")
   
   #rows = gp.SearchCursor(raster)
   #row = rows.Next()
   LogFile.write("Just finished processing" + raster + "\n")
  except:
   gp.AddMessage(gp.GetMessages(2))
   ErrorMsg = gp.GetMessages(2)
   print ErrorMsg
   LogFile.write(ErrorMsg + "\n")
  
except:
 gp.AddMessage(gp.GetMessages(2))
 ErrorMsg = gp.GetMessages(2)
 logfile.write(ErrorMsg + "\n")
 

LogFile.close()

thanks for the assistance with this....it's probably something small...hopefully

Tom Burley
U.S. Geological Survey
Austin, TX

[KenHartling]

Tom,

The zip file appears to be corrupt, but I've figured it out with my own data.

When I first took a look at your sample code you were adding an integer field to the raster attr table and then performing interger division ... but what you actually want to do is "true division".

Python, like C, performs "floor" division by default. No remainders are returned when using integers so for your workflow, with the code and data you were using, the field will always be calced to 0.

So, first change the field you create to something that will accept floats (I used DOUBLE).
And then you have two options:
1. Use '/' with floats.
2. Use a code block that switches python from "floor" division to "true" division.

Find attached the new script I put together using a code block.

Summary of changes:
Defined a variable, EXPRESSION, containing the code block:
EXPRESSION = '''
from __future__ import division
def realDiv(fld,val):
    return fld / val
'''

Add field as DOUBLE, not LONG:
#gp.addfield (raster, "PROPORTION", "LONG", "")
gp.addfield (raster, "PROPORTION", "DOUBLE", "")

Use the code block in calculatefield:
#gp.calculatefield(raster,"PROPORTION","!COUNT! / " + str(x),"PYTHON")
gp.calculatefield(raster,"PROPORTION","realDiv(!COUNT!,x)","PYTHON_9.3",EXPRESSION)

Please let me know how it goes.
Ken

[ThomasBurley]

Ken,
that worked, just tested it on my data. I would have never figured that out. Yeah I'm not familiar with C so that was completely new for me - this is really the first week I've messed with Python substantially.

I'm working through the online free Dive Into Python book and have been dinking around with some stuff I've wanted to automate with geoprocessing like this.

Do you have any online references for this? I don't quite understand what's going on with the Expression that's now being references via the codeblock parameter

I just did a quick google search and found this: http://www.python.org/dev/peps/pep-0238/

Does that cover the gist of what you changed and what I was missing?

I didn't see any references to the function you include after the from __future__ part with the realDiv function - what does that do? That's just handling the passed-in Count field variable as the numerator and the x sum variable as the denominator?
Is that a necessary item to include or would the future import division statement as described in the above url be sufficient?

thanks much for your help - this was driving me insane

Tom

[KenHartling]

Hi Tom,

You are not alone.  There has been such an debate going on in the python community over this that Python will be changing their default behaviour to "true" division in version 3, but we (and many others) will not be moving to Python 3 for quite a long time so this is the solution you should use until then.

Our doc with CalculateField examples, include using code blocks is here:
http://help.arcgis.com/en/arcgisdesktop/10.0/help/index.html#/Calculate_Field_examples/00170000004s0...

The link you sent is one of many discussing part of the issue.

Good luck.
Ken