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
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