Implementing merge() in ArcPy map algebra

curtvprice — Mon, 01 Sep 2014 16:17:12 GMT

The merge function (which we've used for years in GRID, ArcView, and ArcGIS 8x, 9x) simply doesn't exist in ArcPy map algebra, unless I'm missing something...

Here's my hack. Note this is set up so it should work inside a map algebra expression. It also returns an integer raster if all the inputs are integer, otherwise float.

Just asking the community (and especially Esri SA team folks - I know you're out there) whether this is a viable approach.

import os
import arcpy


def MergeRasters(raster_list, out_raster=None):  
    """ArcPy Map Algebra implementation of the old merge() function. 

     Curtis Price, cprice@usgs.gov 
    """  
    
    # output raster name
    if not out_raster:
        out_raster = arcpy.CreateScratchName(
            "mrg", "", "raster", arcpy.env.workspace)  
        
    # determine raster type for output
    ptype = "32_BIT_SIGNED"  
    for ras in raster_list:  
        if not arcpy.Describe(ras).isInteger:  
            ptype = "32_BIT_FLOAT"  
            break  
        
    # create temporary raster
    tmpRaster = arcpy.CreateScratchName("", "mrg", "raster",  
                                          arcpy.env.scratchFolder)  
    arcpy.CreateRasterDataset_management(os.path.dirname(tmpRaster),  
                              os.path.basename(tmpRaster),  
                              pixel_type=ptype)  
    arcpy.Mosaic_management(raster_list, tmpRaster, "FIRST")  
    
    # Copy the raster, enforcing environment (extent, snap, mask, etc)  
    mRaster = ApplyEnvironment(tmpRaster)
    mRaster.save(out_raster)
    
    # Clean up  
    arcpy.Delete_management(tmpRaster)  
    return mRaster

Note, this may work sort of like the merge() function, but is not true arcpy map algebra local operator, so we lose out on the benefits of using it with complex statements like

Int(Merge(raster1, SetNull(raster2, 2), raster3 + 2))

Please vote up my idea: Bring back Spatial Analyst Merge tool

UPDATES

curtvprice - script bug fix

curtvprice - Used ApplyEnvironment function instead of just multiplying by one

curtvprice - added note about complex statements and link to ideas.esri.com

Re: Implementing merge() in ArcPy map algebra

curtvprice — Sun, 12 Dec 2021 07:08:59 GMT

Here's a new implementation. I think this is a much better replacement for merge().

(I did get some help on this from one of the Spatial Analyst developer luminaries.)

import arcpy
from arcpy.sa import *
arcpy.CheckOutExtension("spatial")

def Merge(raster_list): 
    """Merge a list of rasters. The first non-NoData value found (cell by cell) is kept.
       Merge([raster1, raster2, raster3, ...])

    This is another implementation of ArcGIS 9x Map Algebra merge() function"""
    if type(raster_list) == str:
        raster_list = raster_list.split(";")
    if len(raster_list) < 2:
        raise Exception("More than one raster required")
    rlist = [Raster(arcpy.Describe(ras).catalogPath) for ras in raster_list]
    merge_ras = rlist[0]
    for ras in rlist[1:]:
        merge_ras = Con(IsNull(merge_ras), ras, merge_ras)
    return merge_ras

topic Implementing merge() in ArcPy map algebra in ArcGIS Spatial Analyst Questions

Implementing merge() in ArcPy map algebra

Re: Implementing merge() in ArcPy map algebra