I've written a Python script that reclassifies a raster using a set of conditional statements into 4 classes:
rec_ras = Con(costmin <= 12, 12, Con((costmin > 12) & (costmin <= 24), 24, Con((costmin > 24) & (costmin <= 30), 30, Con((costmin > 30) & (costmin <= 60), 60))))
I then needed to change the classes into 7 classes:
rec_ras = Con(costmin <= 5, 5, Con((costmin > 5) & (costmin <= 10), 10, Con((costmin > 10) & (costmin <= 15), 15, Con((costmin > 15) & (costmin <= 20), 20, Con((costmin > 20) & (costmin <= 25), 25, Con((costmin > 25) & (costmin <= 30), 30, Con((costmin > 30) & (costmin <= 60), 60)))))))
Is ther a better way of creating a Python function that can create my conditional statement as and if when I need to change the classes before running it. The following is part of a larger Python Class.
Steve... I venture is was because of the speed and pure beauty of numpy
Hi Dan
I found that the new script that I wrote using Numpy has a further advantage in that it supports reclassifying an integer raster to float raster.
''' Created on Mar 19, 2016 @author: PeterW ''' # import site-packages and modules import os import numpy as np import arcpy from arcpy.sa import * # set arguments input_raster = arcpy.Raster(r"E:\Projects\2016\G112030\SAA\Rasters\euc_road_clip") output_folder = r"E:\Projects\2016\G112030\SAA\Rasters" # set environment settings arcpy.env.overwriteOutput = True arcpy.env.outputCoordinateSystem = input_raster # check out extensions arcpy.CheckOutExtension("Spatial") def reclass_raster(input_raster, output_folder): lower_left = arcpy.Point(input_raster.extent.XMin, input_raster.extent.YMin) cell_size = input_raster.meanCellWidth saa_array = arcpy.RasterToNumPyArray(input_raster) reclass_array = np.zeros_like(saa_array) bins = [0, 10, 25, 50, 100, 13069.5234375] new_bins = [1, 1.25, 1.5, 2, 5] new_classes = zip(bins[:-1], bins[1:], new_bins) print(new_classes) for rc in new_classes: q1 = (saa_array >= rc[0]) q2 = (saa_array < rc[1]) z = np.where(q1 & q2, rc[2], 0) reclass_array = reclass_array + z reclass_raster = arcpy.NumPyArrayToRaster(reclass_array, lower_left, x_cell_size=cell_size, value_to_nodata=0) reclass_raster.save(os.path.join(output_folder, "cost_road")) reclass_raster(input_raster, output_folder) # check in extensions arcpy.CheckInExtension("Spatial")
Nice... further reducing steps and showing how arcpy and numpy get along well.
Hi Steve
I was originally using map algebra Conditional Statement to reclassify my output raster and was looking for a way to build a function for the new classes that I wanted the raster to be reclassified to. Darren and Dan suggested that Numpy would allow me to generate a function to define the classes more cleanly. I must admit I had not considered the reclassify tool, RemapRange would have worked just just as well.
You numpy fanboys. I mean really.
Another thing about spatial analyst tools is they are somewhat more likely to work when you throw a 5GB raster at them. As I have been discovering recently processing big chunks of 10m elevation data!
RemapRange—Spatial Analyst module | ArcGIS for Desktop
Four lines, my friends! (OK, Dan, six lines with the imports)
>>> import arcpy >>> from arcpy.sa import * >>> bins = [0, 10, 25, 50, 100, 13069.5234375] >>> new_bins = [1, 1.25, 1.5, 2, 5] >>> remap = [[bins, bins[k+1], int(new_bins >>> for kk in remap: print(kk) ... [0, 10, 100] [10, 25, 125] [25, 50, 150] [50, 100, 200] [100, 13069.5234375, 500] >>> out_raster = Reclassify("in_raster", "VALUE", RemapRange(remap)) * 0.01* 100)] for k in range(len(bins) - 1)]
missing the import lines of course
rows and columns of raster?
assuming it is float of course...
just curious
Okay, I added the import lines for you. This is Python map algebra, I don't need to worry about the raster dimensions!
I really liked Peter's nifty construction of the remap table below, but I needed to scale the data for the float output he was looking for, so I had to go with my somewhat clunkier (but still so pythonesque) list comprehension.
new_classes = zip(bins[:-1], bins[1:], new_bins) # cool
Oh... it just changes the 'view' of the raster, it doesn't actually produce a new reclassed one (if that is correct)
yes, one of several ways to generate it using LCs, but be cautious of zip
When people start using python 3.x + make sure to use list(zip(....)) or you will get a zip object... which is useful, or not depending upon what you need to do.
bins = [0, 10, 25, 50, 100, 13079]
new_bins = [1, 1.25, 1.5, 2, 5]
nc = [ (bins[i-1], bins, new_bins[i-1]) for i in range(1,len(bins)) ] # safe in python 2.x and 3.x
or
zip(bins[:-1],bins[1:],new_bins) # python 2.x only
or
list(zip(blah)) # in python 3.x
yields
[(0, 10, 1), (10, 25, 1.25), (25, 50, 1.5), (50, 100, 2), (100, 13079, 5)]
Sounds like yet another thing to do slightly differently to make your code work in both.
list(zip(blah)) it is.