0

I have two raster layers and a function sand.tif, clay.tif and def soiltexturalclass(sand,clay):. The function is a conditional function which can be viewed here.

The goal is to create an output raster from the result of the function. I was able to write a function that successfully creates a new array from my function. When trying to run the below function though, I receive a blank raster.

IE outras = soiltexturalclass(sand.tif,clay.tif)

import numpy as np

class SoilTexture():

    def __init__(self):
        self.name = "NRCS Soil Texture Classification"
        self.description = "Computes soil texture given clay and sand content."
        
    def getParameterInfo(self):
        return [
            { 
                'name': 'Clay',
                'dataType': 'raster',
                'value': None,
                'required': True,
                'displayName': "Clay %",
                'description': "Single-band raster where pixel values represent clay content in mass/mass."
            },
            {
                'name': 'Sand',
                'dataType': 'raster',
                'value': None,
                'required': True,
                'displayName': "Sand %",
                'description': "Single-band raster where pixel values represent sand content in mass/mass."
            },
        ]
    
    def getConfiguration(self, **scalars):
        return {
            'inheritProperties': 2 | 4 | 8,     # inherit all but the pixel type from the input raster
            'invalidateProperties': 2 | 4 | 8,  # reset any statistics and histogram that might be held by
                                                #   the parent dataset (because this function modifies pixel values).
            'inputMask': True                   # We need the input raster mask in .updatePixels().
        }

    def updateRasterInfo(self, **kwargs):
        kwargs['output_info']['bandCount'] = 1
        kwargs['output_info']['pixelType'] = 'u1'
        kwargs['output_info']['statistics'] = ()
        kwargs['output_info']['noData'] = np.array([0], 'u1')
        return kwargs

    def updatePixels(self, tlc, shape, props, **pixelBlocks):

        sand = np.array(pixelBlocks['sand_pixels'], dtype='f4', copy=False)
        clay = np.array(pixelBlocks['clay_pixels'], dtype='f4', copy=False)

        sand = sand/10
        clay = clay/10
        silt = (100 - sand - clay).astype('f4', copy=False)
        def TextureClass(sand, clay, silt):
    
        #if sand + clay > 100 or sand < 0 or clay < 0:
        #    raise Exception('Inputs adds over 100% or are negative')

            if silt + 1.5*clay < 15:
                textural_class = 'sand'

            elif silt + 1.5*clay >= 15 and silt + 2*clay < 30:
                textural_class = 'loamy sand'

            elif (clay >= 7 and clay < 20 and sand > 52 and silt + 2*clay >= 30) or (clay < 7 and silt < 50 and silt + 2*clay >= 30):
                textural_class = 'sandy loam'

            elif clay >= 7 and clay < 27 and silt >= 28 and silt < 50 and sand <= 52:
                textural_class = 'loam'

            elif (silt >= 50 and clay >= 12 and clay < 27) or (silt >= 50 and silt < 80 and clay < 12):
                textural_class = 'silt loam'

            elif silt >= 80 and clay < 12:
                textural_class = 'silt'

            elif clay >= 20 and clay < 35 and silt < 28 and sand > 45:
                textural_class = 'sandy clay loam'

            elif clay >= 27 and clay < 40 and sand > 20 and sand <= 45:
                textural_class = 'clay loam'

            elif clay >= 27 and clay < 40 and sand <= 20:
                textural_class = 'silty clay loam'

            elif clay >= 35 and sand > 45:
                textural_class = 'sandy clay'

            elif clay >= 40 and silt >= 40:
                textural_class = 'silty clay'

            elif clay >= 40 and sand <= 45 and silt < 40:
                textural_class = 'clay'

            else:
                textural_class = 'na'

            return textural_class

        vect_function = np.vectorize(TextureClass)

        pixelBlocks['output_pixels'] = vect_function(sand, clay, silt).astype(props['pixelType'], copy=False)
        return pixelBlocks
1
  • 1
    I don't know of a way to define functions like that inside the Raster Calculator, it's not set up like the Field Calculator. But, since you already have your function in python, you could use RasterToNumPyArray, run the calculations in python, and save a new raster from there. pro.arcgis.com/en/pro-app/latest/arcpy/functions/…
    – Brennan
    Aug 5 at 17:50

0

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.