Merge with custom method in rasterio

Question

I would like to merge a few rasters by using a custom method, and not the default one (reverse painting). For this purpose I am trying to figure out how to pass the custom callable.

For example, I would like to merge raster cells by using the cell with the maximum value. Therefore here is my callback function:

def custom_merge(old_data, new_data, old_nodata, new_nodata, index=None, roff=None, coff=None):
    new_data = np.maximum(old_data, new_data)

And here is how I am calling it:

mosaic, out_trans = merge(src_files_to_mosaic_opened, method=custom_merge)

I am not getting an error, however I am gettin a raster with all zeros.

How do I make this work?

Answer 1

Setting new_data data to a new variable (what you did) or returning a result (as suggested in a comment) will not work.

The documentation could be clearer, but you need to update (in place) the old data array old_data:

def function(old_data, etc...):

old_data array_like

array to update with new_data

Here's a worked example:

import numpy as np
import rasterio as rio
from rasterio.merge import merge
def custom_merge_works(old_data, new_data, old_nodata, new_nodata, index=None, roff=None, coff=None):
    old_data[:] = np.maximum(old_data, new_data)  # <== NOTE old_data[:] updates the old data array in place
def custom_merge_doesnt_work(old_data, new_data, old_nodata, new_nodata, index=None, roff=None, coff=None):
    return np.maximum(old_data, new_data)
with rio.open('test1.tif') as test1, rio.open('test2.tif') as test2:
    arr1, arr2 = test1.read(), test2.read()
mosaic, out_trans = merge([test1, test2], method=custom_merge_doesnt_work)
print(&quot;All zeros\n&quot;, mosaic)

mosaic, out_trans = merge([test1, test2], method=custom_merge_works)
print(&quot;Expected result\n&quot;, mosaic)

Output:

All zeros
 [[[0. 0. 0. ... 0. 0. 0.]
  [0. 0. 0. ... 0. 0. 0.]
  [0. 0. 0. ... 0. 0. 0.]
  ...
  [0. 0. 0. ... 0. 0. 0.]
  [0. 0. 0. ... 0. 0. 0.]
  [0. 0. 0. ... 0. 0. 0.]]]
Expected result
 [[[ 5.  7.  5. ...  9.  6.  5.]
  [ 5.  5. 10. ...  5. 10.  7.]
  [ 5.  8.  9. ...  7.  6.  5.]
  ...
  [ 5.  5.  6. ...  5.  9.  5.]
  [10.  7. 10. ...  8.  8.  5.]
  [ 7.  6.  5. ...  6.  6.  9.]]]