Updating table fields from a raster file

445
6
12-05-2018 03:08 AM
RenanBarroso
New Contributor II

Hi there! I have some raster files in grid format with signed integer pixels in which I want to add some fields in their table then compute values with some formulas in the new fields. I managed to add the fields but I'm struggling with the second part. I tried with arcpy.da.UpdateCursor but it returns a RuntimeError: 'in_table' is not a table or a feature class. It seems clear for me that I can not use this function with raster file. Is there a way to compute the values directly in the raster or do I have to try a different approach? How?

The code is:

import arcpy as ap
import os
import random_function as rf # change name after
from arcpy.sa import *

if ap.CheckExtension('Spatial') == 'Available':
     ap.CheckOutExtension('Spatial')
else: 
     print('Spatial Analyst extension is not available')

# Setting environments & paths
ap.env.workspace = os.path.join("D:\\", "master", "_thesis_data")
ap.env.scratchWorkspace = os.path.join(ap.env.workspace, "temp")
ap.env.overwriteOutput = True
inp_ws = ap.env.workspace
temp_ws = ap.env.scratchWorkspace 
lu_folder = os.path.join(inp_ws, "lu")

lu_scenarios = []
for dirpath, dirnames, filenames in ap.da.Walk(lu_folder, datatype='RasterDataset'):
     for f in filenames:
          lu_scenarios.append(os.path.join(dirpath,f))

# Setting input data 
climate = Raster("climate_c")
clim_soil = Raster("rc_climsoil_c")
tbl_uSOCr = os.path.join(inp_ws, "uSOCr.txt")
tbl_uSOCf = os.path.join(inp_ws, "uSOCf.txt")
tbl_uBCi = os.path.join(inp_ws, "uBCi.txt")

for (sc, lu_sc) in enumerate(lu_scenarios):
     lu = Raster(lu_sc)
     clim_lu = climate + lu
     soc_f1k = ReclassByTable(clim_lu, tbl_uSOCf, "id", "id", "soc_factor")
     soc_r1k = ReclassByTable(clim_soil, tbl_uSOCr, "id", "id", "soc_ref")
     soc_1k = (soc_f1k*soc_r1k)/1000 
     bc_1k = ReclassByTable(clim_lu, tbl_uBCi, "id", "id", "bcs")          
     tc_1k = (soc_1k + bc_1k)
     tc_1k_new_fields = ["tc", "tot_ha", "tot_cs"]
     for field in tc_1k_new_fields:
          ap.AddField_management(tc_1k, field, "FLOAT", 12, 4)
     with ap.da.UpdateCursor(tc_1k, "*") as cursor:
          tc_all_fields = [f.name for f in arcpy.ListFields(tc_1k)]               
          VALUE, COUNT, tc, tot_ha, tot_cs = row[1], row[2], row[3], row[4], row[5]
                 for row in cursor:
                    tc = VALUE/1000 
                    tot_ha = COUNT*2500  
                    tot_cs = tc*tot_ha  
                    cursor.updateRow(row)                                   
                    
     if ap.GetRasterProperties_management(tc_1k, "ALLNODATA"):
          tc_1k.save(temp_ws + "\\" + "0sc" + str(sc) + "_tc" + str(mc_counter+1)) # ==> TESTING DATA!!!
     print tc_1k‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍

Thanks again!

Tags (1)
0 Kudos
6 Replies
DanPatterson_Retired
MVP Esteemed Contributor

You could export the raster table to a standalone table or a table in a geodatabase

0 Kudos
RenanBarroso
New Contributor II

And how can I do that in arcpy? I was also wondering if there is a way to have the raster attribute tables in numpy array. I have tried the RasterToNumpyArray, but without success. The thing I want to do seems quite easy but I'm really struggling with it =/.

Thanks again!

0 Kudos
DanPatterson_Retired
MVP Esteemed Contributor

I use RasterToNumPyArray all the time.  What was your problem.

pth = r"C:\Temp\dem.tif"


pth = r"C:\Temp\myarray.tif"

a = RasterToNumPyArray(pth)

a.shape
Out[7]: (50, 50)

a
Out[8]: 
array([[68, 33, 78, ..., 53, 39, 89],
       [87, 88, 22, ..., 75,  7, 96],
       [69, 51, 77, ..., 84, 69, 83],
       ...,
       [73,  0,  1, ..., 29, 25, 88],
       [81, 60, 94, ..., 28, 23, 68],
       [47,  5, 55, ..., 19, 78,  9]])

a - a.mean()
 
array([[ 17.9472, -17.0528,  27.9472, ...,   2.9472, -11.0528,  38.9472],
       [ 36.9472,  37.9472, -28.0528, ...,  24.9472, -43.0528,  45.9472],
       [ 18.9472,   0.9472,  26.9472, ...,  33.9472,  18.9472,  32.9472],
       ...,
       [ 22.9472, -50.0528, -49.0528, ..., -21.0528, -25.0528,  37.9472],
       [ 30.9472,   9.9472,  43.9472, ..., -22.0528, -27.0528,  17.9472],
       [ -3.0528, -45.0528,   4.9472, ..., -31.0528,  27.9472, -41.0528]])
RenanBarroso
New Contributor II

Thank you! Now it is working and I've changed all my raster processing to NumpyArray. Much better! 

DanPatterson_Retired
MVP Esteemed Contributor

I have lots of NumPy examples on my blog post if you are interested in numpy functionality with arcpy and other modules.

XanderBakker
Esri Esteemed Contributor

The cursors in the da module do not provide access to raster attribute tables (RAT). The old cursors do. Although I have never used it to write to a RAT, I have used the old cursors to read from it...

https://community.esri.com/message/355845?commentID=355845#comment-633411