Error 010240 arcpy-numpy

2368
5
Jump to solution
04-16-2016 05:44 AM
IvanD_Ortenzio1
New Contributor III

I'm following this guide to extract pixels with certain rgb code (eg 255,0,0)  in a .tif raster using Anaconda 32bit:

RasterToNumPyArray—ArcPy Functions | ArcGIS for Desktop

But I receive this error:

os.chdir("C:\\Users\\ivn\\Desktop\\NumpArr")

filein = os.path.join(os.getcwd(),r'AR00199.tif')

fGDB = "fGDB.gdb"

if not os.path.exists(os.path.join(os.getcwd(),fGDB)):

arcpy.CreateFileGDB_management(os.getcwd(),fGDB)

fileoutgdb = os.path.join(os.getcwd(),r"fGDB.gdb\blocco")

arcpy.env.overwriteOutput = True

arcpy.env.outputCoordinateSystem = filein

arcpy.env.cellSize = filein

arcpy.env.workspace = os.path.join(os.getcwd(),r"fGDB.gdb")

arcpy.env.scratchWorkspace = os.path.join(os.getcwd(),r"fGDB.gdb")

inRas = arcpy.Raster(filein)

b1 = arcpy.sa.Raster(r'AR00199.tif\Band_1')

b2 = arcpy.sa.Raster(r'AR00199.tif\Band_2')

b3 = arcpy.sa.Raster(r'AR00199.tif\Band_3')

blocksize = 512

blockno = 0

for x in range(0,inRas.width,blocksize): 

    for y in range(0, inRas.height, blocksize):

        blockno += 1

        mx = inRas.extent.XMin + x * inRas.meanCellWidth

        my = inRas.extent.YMin + y * inRas.meanCellHeight

        lx = min([x + blocksize, inRas.width])

        ly = min([y + blocksize, inRas.height])

        b1nr = arcpy.RasterToNumPyArray(b1,arcpy.Point(mx, my),lx-x, ly-y)      

        b2nr = arcpy.RasterToNumPyArray(b2,arcpy.Point(mx, my),lx-x, ly-y)

        b3nr = arcpy.RasterToNumPyArray(b3,arcpy.Point(mx, my),lx-x, ly-y)

        rgb = np.array([b1nr, b2nr,b3nr], np.uint8)

   

         #Operations on blocks....

        ....

        ....

        ....

        myRasterBlock = arcpy.NumPyArrayToRaster(output_array,\

        arcpy.Point(mx, my),inRas.meanCellWidth,inRas.meanCellHeight)

        filetemp = '_'.join([fileoutgdb,str(blockno)])

        myRasterBlock.save(filetemp)

It save about 173 grids with blocksize = 512 and 73 whit blocksize = 1024 in gdb and then I receive the error... why?? Can someone help?

This is the result of the operation for a block:

Tags (2)
0 Kudos
1 Solution

Accepted Solutions
JoshuaBixby
MVP Esteemed Contributor

My guess is something is happening along the way with "#Operations on blocks" and arcpy.NumPyArrayToRaster() may not be returning a valid array.  I would do some debugging of myRasterBlock before saving it, e.g., insert some print statements to print out properties of the raster.  Also, print out filetemp .  If you could include the results of those print statements, it might help people provide more feedback.

View solution in original post

5 Replies
JoshuaBixby
MVP Esteemed Contributor

My guess is something is happening along the way with "#Operations on blocks" and arcpy.NumPyArrayToRaster() may not be returning a valid array.  I would do some debugging of myRasterBlock before saving it, e.g., insert some print statements to print out properties of the raster.  Also, print out filetemp .  If you could include the results of those print statements, it might help people provide more feedback.

View solution in original post

DanPatterson_Retired
MVP Esteemed Contributor

You use the term 'block' ... do you mean memory size block or a physical block as in the number of rows and columns.

Here is a 3D array and its construction

import numpy as np

a = np.ones((3,6,6))  # make a 3D array

a[0] = np.arange(36).reshape(6,6) # 36 sequential numbers

a[1] = 2  # assign all 2

a[2] = 3  # ditto but use 3

>>> a  # like an image
array([[[  0.,   1.,   2.,   3.,   4.,   5.],
        [  6.,   7.,   8.,   9.,  10.,  11.],
        [ 12.,  13.,  14.,  15.,  16.,  17.],
        [ 18.,  19.,  20.,  21.,  22.,  23.],
        [ 24.,  25.,  26.,  27.,  28.,  29.],
        [ 30.,  31.,  32.,  33.,  34.,  35.]],

       [[  2.,   2.,   2.,   2.,   2.,   2.],
        [  2.,   2.,   2.,   2.,   2.,   2.],
        [  2.,   2.,   2.,   2.,   2.,   2.],
        [  2.,   2.,   2.,   2.,   2.,   2.],
        [  2.,   2.,   2.,   2.,   2.,   2.],
        [  2.,   2.,   2.,   2.,   2.,   2.]],

       [[  3.,   3.,   3.,   3.,   3.,   3.],
        [  3.,   3.,   3.,   3.,   3.,   3.],
        [  3.,   3.,   3.,   3.,   3.,   3.],
        [  3.,   3.,   3.,   3.,   3.,   3.],
        [  3.,   3.,   3.,   3.,   3.,   3.],
        [  3.,   3.,   3.,   3.,   3.,   3.]]])

Now you can block it in 3 or 2 dimensions etc.  This is the top left 3x3 block for all 3 'bands'

>>> a[:,:3,:3]
array([[[  0.,   1.,   2.],
        [  6.,   7.,   8.],
        [ 12.,  13.,  14.]],

       [[  2.,   2.,   2.],
        [  2.,   2.,   2.],
        [  2.,   2.,   2.]],

       [[  3.,   3.,   3.],
        [  3.,   3.,   3.],
        [  3.,   3.,   3.]]])

So far so good?  Here is the average for each of those 3x3 blocks for the first band (just scrunched up the answer)

>>> np.mean(a[0])       # 17.5   overall for band 1

>>> np.mean(a[0,:3,:3])  #  7.0   top left 3x3 block band 0

>>> np.mean(a[1,:3,:3])  #  2.0   top left 3x3 block band 1

>>> np.mean(a[2,:3,:3])  #  3.0   top left 3x3 block band 2

Now you normally work with blocks based on the number of rows and columns rather than 'stride size'.  But you can block an array based on the band and the number of rows and columns rather than that arbitrary unexplained 512 block.

here is some code to work with blocks to get you thinks about how to tackle to problem

import numpy as np
from numpy.lib.stride_tricks import as_strided

def block_array(a, block=(3,3)): 
    """Provide a 2D block view of a 2D array or a slice of a 3D array.
    No error checking made. Columns
    and rows outside of the block are truncated.
    """ 
    r, c = block  # 3x3 block default 
    shape = (a.shape[0]/r, a.shape[1]/c) + block 
    strides = (r*a.strides[0], c*a.strides[1]) + a.strides 
    b_a = as_strided(a, shape=shape, strides=strides) 
    return b_a 

# assign all 3a = np.ones((3,6,6))   # make a 3D array
a[0] = np.arange(36).reshape(6,6) # 36 sequential numbers
a[1] = 2  # assign all 2
a[2] = 3

block = (3,3)
b_0 = block_array(a[0],block)
b_1 = block_array(a[1],block)
b_2 = block_array(a[2],block)

DanPatterson_Retired
MVP Esteemed Contributor

ps...

I forgot to add, that the simple average through all 3 bands is simply

>>> np.set_printoptions(precision=2,threshold=10,edgeitems=3,linewidth=75,suppress=True)

>>> np.mean(a,axis=0,keepdims=False)

array([[  1.67,   2.  ,   2.33,   2.67,   3.  ,   3.33],

       [  3.67,   4.  ,   4.33,   4.67,   5.  ,   5.33],

       [  5.67,   6.  ,   6.33,   6.67,   7.  ,   7.33],

       [  7.67,   8.  ,   8.33,   8.67,   9.  ,   9.33],

       [  9.67,  10.  ,  10.33,  10.67,  11.  ,  11.33],

       [ 11.67,  12.  ,  12.33,  12.67,  13.  ,  13.33]])

>>>

IvanD_Ortenzio1
New Contributor III

I solved the problem by changing some lines of code... here is the result:

Dan Patterson thanks for the suggestions. I will keep in mind for future elaborations.

curtvprice
MVP Esteemed Contributor

Please share your fixed code, so others can benefit. Thanks!

0 Kudos