I have a large .tif raster with integer values. I wish to generate Jenks classification breaks using code and then use those breaks to reclassify my .tif file.
I came across a block of code that I am using to define my Jenks function. Its source is: https://gist.github.com/drewda/1299198
So far I have imported my raster, converted it to a (233,679) 2D array and then converts to a (158207,) 1D array.
Then I attempt to run the function on the array but get the error:
Traceback (most recent call last):
File "C:\Users\agates44\Desktop\FinalProject\script1.py", line 80, in <module>
new = getJenksBreaks(oneDarr,9)
File "C:\Users\agates44\Desktop\FinalProject\script1.py", line 40, in getJenksBreaks
mat2.append(temp)
MemoryError
My code is below
# Define the folder path
folderPath = r'C:\Users\agates44\Desktop\FinalProject'
# Import modules
import arcpy
import numpy as np
# Set workspace
from arcpy import env
env.workspace = folderPath
# Allow overwrite
env.overwriteOutput = True
# Get Input Raster properties
inRas = arcpy.Raster('smallorange.tif')
lowerLeft = arcpy.Point(inRas.extent.XMin,inRas.extent.YMin)
cellSize = inRas.meanCellWidth
# Convert Raster to numpy array
arr = arcpy.RasterToNumPyArray(inRas,nodata_to_value=0)
print arr.shape
print arr.dtype
oneDarr = arr.ravel()
print oneDarr.shape
print oneDarr.dtype
def getJenksBreaks(dataList, numClass):
dataList.sort()
mat1 = []
for i in range(0,len(dataList)+1):
temp = []
for j in range(0,numClass+1):
temp.append(0)
mat1.append(temp)
mat2 = []
for i in range(0,len(dataList)+1):
temp = []
for j in range(0,numClass+1):
temp.append(0)
mat2.append(temp)
for i in range(1,numClass+1):
mat1[1][i] = 1
mat2[1][i] = 0
for j in range(2,len(dataList)+1):
mat2[j][i] = float('inf')
v = 0.0
for l in range(2,len(dataList)+1):
s1 = 0.0
s2 = 0.0
w = 0.0
for m in range(1,l+1):
i3 = l - m + 1
val = float(dataList[i3-1])
s2 += val * val
s1 += val
w += 1
v = s2 - (s1 * s1) / w
i4 = i3 - 1
if i4 != 0:
for j in range(2,numClass+1):
if mat2[l][j] >= (v + mat2[i4][j - 1]):
mat1[l][j] = i3
mat2[l][j] = v + mat2[i4][j - 1]
mat1[l][1] = 1
mat2[l][1] = v
k = len(dataList)
kclass = []
for i in range(0,numClass+1):
kclass.append(0)
kclass[numClass] = float(dataList[len(dataList) - 1])
countNum = numClass
while countNum >= 2:#print "rank = " + str(mat1[k][countNum])
id = int((mat1[k][countNum]) - 2)
#print "val = " + str(dataList[id])
kclass[countNum - 1] = dataList[id]
k = int((mat1[k][countNum] - 1))
countNum -= 1
return kclass
new = getJenksBreaks(oneDarr,9)
print new
### used after running getJenksBreaks()
##def classify(value, breaks):
## for i in range(1, len(breaks)):
## if value < breaks[i]:
## return i
## return len(breaks) - 1
How can I figure out which line of my code is using all the memory and what is a workaround I can use?
****** EDIT The function ran (although it took forever) on the raster size listed above. But that was only a subset of my much larger dataset.
(2400, 2337) 2D array (5608800,) 1D array
This is when I get the memory error.
