Arcpy Reclassify by Standard Deviation

I would need some help with the Reclassify tool. I have written a script that produces rasters which I need to reclassify. My problem is, I want the values to be reclassified by standard deviation which I don't know how to do. I would like to make two new classes:

values < mean + 2 * standard deviation

values >= mean + 2 * standard deviation

I need help with putting this into a script.

Try something like this:

import arcpy


Set the workspace to the folder where you have your data

arcpy.env.workspace="Path:\to\your\data\"
for raster in arcpy.ListRasters():


Create an object storing the st. dev. for the raster file:

    stdev_object=arcpy.GetRasterProperties_management(raster, "STD")


Store the numeric value of st. dev. in the variable stdev :

    stdev=stdev_object.getOutput[0]


Same thing for mean value:

    mean_object=arcpy.GetRasterProperties_management(raster, "MEAN")
mean=mean.getOutput[0]
threshold=mean + 2* stdev


You can use Raster Calculator for this instead of Reclassify:

    Output_raster="some\path\and\name\for\output\raster"
expression="Con(" + raster + "<=" + str(threshold) + ",1,2"


The string resulting from the above concatenation is an expression that tells Raster Calculator: if a raster cell value is <= than threshold change that value to 1, else change it to 2 (you can of course any values you like).

Now we actually reclassify the raster:

    arcpy.gp.RasterCalculator_sa(expression, Output_raster)


Done!

This works if you have your rasters in a ArcGIS native format, otherwise you need to convert them first. More info here: Get Raster Properties tool

rast = arcpy.sa.Raster(pathtorasterfile)


Find the mean, standard deviations, minimum and maximum:

meanValue = rast.mean
std = rast.standardDeviation
minR = rast.minimum - .1 #just a little buffer room
maxR = rast.maximum + .1 #just a little buffer room


Setup the bounds:

target = meanValue + (2*std)


Create a remaprange object:

remapRangeValues = arcpy.sa.RemapRange([[minR,target-.00001,0],[target,maxR,1]])


Reclassify and save the new raster:

outReclass = arcpy.sa.Reclassify(rast, "VALUE", remapRangeValues)
outReclass.save(reclassifiedname)


I didn't test this, and I wasn't sure what you wanted to reclassify the values as. You may change the 0 or 1 in the RemapRange object to get the values you want.

• Do you mean minR = rast.minimum - .1? – nmpeterson Apr 14 '16 at 19:11
• Thank you! The problem with this is that RemapRange does not seem to take in variables, only actual values which I cannot hard code because I have a set of rasters and they all have different values. – danielaandrea Apr 14 '16 at 19:40
• That's strange. I use it with variables in tools where I don't know what the values will be. Could you edit your question and post your existing code for troubleshooting? – dslamb Apr 14 '16 at 19:43

You can convert your raster to a numpy Array, then use numpy.where and convert back to raster:

import arcpy
import numpy as np

rasterfile = r"C:\GIS\data\DEM50m\grid50m\nh_63_4.tif"
new_raster = r"C:\GIS\data\DEM50m\grid50m\nh_63_4_reclassed.tif"

arr = arcpy.RasterToNumPyArray(rasterfile)
arr = np.where(arr<arr.mean()+2*arr.std(), 1, 2)

newarr = arcpy.NumPyArrayToRaster(in_array=arr, lower_left_corner=r.extent.lowerLeft, x_cell_size=r.meanCellWidth, y_cell_size=r.meanCellHeight)
newarr.save(new_raster)