TransWikia.com

Raster equivalent of "Recalculate Feature Class Extent"

Geographic Information Systems Asked on June 21, 2021

Using Arc Desktop 10.7, I create an output raster using the Extract by Attributes tool. However, the output raster retains the same x/y spatial extent as the input raster, even though its pixels occupy a different extent.

How do I recalculate the output raster’s new spatial extent?

3 Answers

To the contrary of feature dataset which can have holes between two features, a raster is a matrix of predefined size (number of rows and columns) where each cell has a value (this value can be NoData). So even if you see only the data value of your raster, its extent is also taking into account the NoData pixels.

Changing the extent property of a raster will therefore reduce its size (and you'll loose information, but in your case this information is not useful as far as I understand, so it is OK).

To most straightforward tool to change the extent of a raster is the "clip raster" tool, but the challenge is to determine the extent that you want. This can be done by converting the values that you consider as useful into polygons (raster to polygon) then use the extent of the feature class. Make sure that you use "snap raster" to your original raster when you use this extent, otherwise your pixels are likely to be slightly shifted (and thus resampled)

Notes:

  1. if you know "a priori" the extent of your output raster of interest, you can define it in the tool environment which will be much more efficient than running the tool on the full raster and then clipping.

  2. raster extents are always rectangles, so you will probably still have some NoData pixels.

  3. my links are pointing to ArcGIS pro, but in this case it is the same in ArcGIS 10.7.

  4. converting raster to polygon can be very heavy, especialy if you have many isolated pixels of different values. You could first use the raster calculator to set all your "valid" pixel values = 1, and the other as NoData (e.g. Con(raster>5,1) )

  5. A quick and dirty way to achieve this (if you don't care about slightly shifted pixels and do not need automation) is to zoom on the area of interest, right click on the raster, export data and select the extent that is "same as dataframe"

Answered by radouxju on June 21, 2021

The new raster's pixels don't "occupy a different extent", rather more are filed with a NO DATA value, so aren't obviously there. While running the Extract by Attributes tool (or any tool), you can set the processing extent in the environment settings.

  1. With the Extract by Attributes GP tool open, click the "Environments" button.
  2. Expand "Processing Extent"
  3. Pick from the options (same as display, same as a specific layer..., etc)
  4. Click OK, and run the tool.

Per the documentation, "The extent of the output dataset will typically be larger than the Extent setting to account for features or cells that pass through the extent rectangle". Nevertheless, if the processing extent is much smaller than the original raster, don't worry, your output raster will have a smaller extent.

Answered by JMers on June 21, 2021

Here's a version of @whubers method for ArcGIS 10.1+ as a python toolbox (.pyt).

import arcpy

class Toolbox(object):
    def __init__(self):
        """Define the toolbox (the name of the toolbox is the name of the
        .pyt file)."""
        self.label = "Raster Toolbox"
        self.alias = ""

        # List of tool classes associated with this toolbox
        self.tools = [ClipNoData]


class ClipNoData(object):
    def __init__(self):
        """Clip raster extent to the data that have values"""
        self.label = "Clip NoData"
        self.description = "Clip raster extent to the data that have values. "
        self.description += "Method by Bill Huber - https://gis.stackexchange.com/a/55150/2856"

        self.canRunInBackground = False

    def getParameterInfo(self):
        """Define parameter definitions"""
        params = []

        # First parameter
        params+=[arcpy.Parameter(
            displayName="Input Raster",
            name="in_raster",
            datatype='GPRasterLayer',
            parameterType="Required",
            direction="Input")
        ]

        # Second parameter
        params+=[arcpy.Parameter(
            displayName="Output Raster",
            name="out_raster",
            datatype="DERasterDataset",
            parameterType="Required",
            direction="Output")
        ]

        return params

    def isLicensed(self):
        """Set whether tool is licensed to execute."""
        return arcpy.CheckExtension('spatial')==u'Available'

    def execute(self, parameters, messages):
        """See https://gis.stackexchange.com/a/55150/2856
           ##Code comments paraphrased from @whubers GIS StackExchange answer
        """
        try:
            #Setup
            arcpy.CheckOutExtension('spatial')
            from arcpy.sa import *
            in_raster = parameters[0].valueAsText
            out_raster = parameters[1].valueAsText

            dsc=arcpy.Describe(in_raster)
            xmin=dsc.extent.XMin
            ymin=dsc.extent.YMin
            mx=dsc.meanCellWidth
            my=dsc.meanCellHeight
            arcpy.env.extent=in_raster
            arcpy.env.cellSize=in_raster
            arcpy.AddMessage(out_raster)

            ## 1. Convert the grid into a single zone by equating it with itself
            arcpy.AddMessage(r'1. Convert the grid into a single zone by equating it with itself...')
            zones = Raster(in_raster) == Raster(in_raster)

            ## 2. Create a column index grid by flow-accumulating a constant grid
            ##    with value 1. (The indexes will start with 0.) Multiply this by
            ##    the cellsize and add the x-coordinate of the origin to obtain
            ##    an x-coordinate grid.
            arcpy.AddMessage(r'Create a constant grid...')
            const = CreateConstantRaster(1)

            arcpy.AddMessage(r'2. Create an x-coordinate grid...')
            xmap = (FlowAccumulation(const)) * mx + xmin

            ## 3. Similarly, create a y-coordinate grid by flow-accumulating a
            ##    constant grid with value 64.
            arcpy.AddMessage(r'3. Create a y-coordinate grid...')
            ymap = (FlowAccumulation(const * 64)) * my + ymin

            ## 4. Use the zone grid from step (1) to compute the zonal min and
            ##    max of "X" and "Y"
            arcpy.AddMessage(r'4. Use the zone grid from step (1) to compute the zonal min and max of "X" and "Y"...')

            xminmax=ZonalStatisticsAsTable(zones, "value", xmap,r"IN_MEMORYxrange", "NODATA", "MIN_MAX")
            xmin,xmax=arcpy.da.SearchCursor(r"IN_MEMORYxrange", ["MIN","MAX"]).next()

            yminmax=ZonalStatisticsAsTable(zones, "value", ymap,r"IN_MEMORYyrange", "NODATA", "MIN_MAX")
            ymin,ymax=arcpy.da.SearchCursor(r"IN_MEMORYyrange", ["MIN","MAX"]).next()

            arcpy.Delete_management(r"IN_MEMORYxrange")
            arcpy.Delete_management(r"IN_MEMORYyrange")

            # Write out the reduced raster
            arcpy.env.extent = arcpy.Extent(xmin,ymin,xmax+mx,ymax+my)
            output = Raster(in_raster) * 1
            output.save(out_raster)

        except:raise
        finally:arcpy.CheckInExtension('spatial')

Answered by user2856 on June 21, 2021

Add your own answers!

Ask a Question

Get help from others!

© 2024 TransWikia.com. All rights reserved. Sites we Love: PCI Database, UKBizDB, Menu Kuliner, Sharing RPP