surprising result for RasterizeLayer with MERGE_ALG="ADD", "ALL_TOUCHED" for polygon

[mathause]

Expected behavior and actual behavior.

When rasterizing a single polygon with options=["MERGE_ALG=ADD", "ALL_TOUCHED"] I expect the maximum value to be the burn in value. However, the borders get added in addition to the core - see the example below.

• Of course this is not relevant for a single polygon (I could just remove the "MERGE_ALG=ADD"), but it also happens for several polygons. I use ADD to separate overlapping regions.
• I already reported in rasterio (surprising result of rasterize with all_touched=True and merge_alg=add rasterio/rasterio#2925) but I think it's more relevant here.
• Similar to gdal_rasterize with -add introduces additional values in the output #1307 but for polygons

Steps to reproduce the problem.

from osgeo import gdal, ogr, osr

# Setup working spatial reference
sr_wkt = 'LOCAL_CS["arbitrary"]'
sr = osr.SpatialReference(sr_wkt)

# Create a memory raster to rasterize into.

target_ds = gdal.GetDriverByName("MEM").Create("", 15, 15, 3, gdal.GDT_Byte)
target_ds.SetGeoTransform((15, -1, 0, 15, 0, -1))

target_ds.SetProjection(sr_wkt)

# Create a memory layer to rasterize from.

rast_ogr_ds = ogr.GetDriverByName("Memory").CreateDataSource("wrk")
rast_mem_lyr = rast_ogr_ds.CreateLayer("poly", srs=sr)

# Add a polygon.

wkt_geom = "POLYGON((12.5 2.5, 12.5 12.5, 2.5 12.5, 2.5 2.5, 12.5 2.5))"

feat = ogr.Feature(rast_mem_lyr.GetLayerDefn())
feat.SetGeometryDirectly(ogr.Geometry(wkt=wkt_geom))

rast_mem_lyr.CreateFeature(feat)

# Run the algorithm.

err = gdal.RasterizeLayer(
    target_ds,
    [2, 1],
    rast_mem_lyr,
    burn_values=[80, -1],
    options=["MERGE_ALG=ADD", "ALL_TOUCHED"],
)

assert err == 0, "got non-zero result code from RasterizeLayer"

gdal.GetDriverByName("GTiff").CreateCopy("rasterize_1.tif", target_ds)

_, maxval = target_ds.GetRasterBand(2).ComputeRasterMinMax()
assert maxval == 80, maxval

Operating system

Linux 6.2.0-33-generic #33~22.04.1-Ubuntu SMP PREEMPT_DYNAMIC x86_64 GNU/Linux

GDAL version and provenance

version 3.7.0 installed from conda-forge

[mathause]

Thanks a lot!