calling gdal_merge into python script

Question

I am new to python and breaking my head on how to use gdal_merge in a (spyder) python script. I use Windows 10, Python 3.6 and have the gdal tool installed from osgeo4w. I realize many other posts describe this problem but none could help me resolving this issue.

When i call the gdal module from cmd, it works like a charm:

 python "C:\OSGeo4W64\bin\gdal_merge.py" -o merged.tif input1.tif input2.tif

However, I cannot get it working properly in a python script (spyder).

A first approach produces an output but not with the right name: it produces an 'out.tif' file and not the 'merged.tif' file as I request:

import sys
sys.path.append('C:\\OSGeo4W64\\bin')
import gdal_merge as gm
gm.main(['-o', 'merged.tif', 'N00W078.tif', 'N00W079.tif'])

A second approach simply produces no output:

import subprocess
merge_command = ["python", "C:\OSGeo4W64\bin\gdal_merge.py", "-o", "merged.tif", "input1.tif", "input2.tif"]
subprocess.call(merge_command,shell=True)

Any thoughts on how to resolve this issue?

Answer 1

Add an empty argument in the first approach (because gdal_merge.py parses arguments starting from 1 and not 0):

import sys
sys.path.append('C:\\OSGeo4W64\\bin')
import gdal_merge as gm
gm.main(['', '-o', 'merged.tif', 'N00W078.tif', 'N00W079.tif'])

Join the path of gdal_merge.py in the second approach:

import os, subprocess
gm = os.path.join('C:\\','OSGeo4W64','bin','gdal_merge.py')
merge_command = ["python", gm, "-o", "merged.tif", "input1.tif", "input2.tif"]
subprocess.call(merge_command,shell=True)

Answer 2

Windows 10, Python 3.9, QGIS 3.20 (installed with OSGeo4W). This code works fine:

import osgeo_utils.gdal_merge
output_file_path = r'insert path here'
input_files_path = r'insert path here'
parameters = ['', '-o', output_file_path] + input_files_path + ['-separate', '-co', 'COMPRESS=LZW']
osgeo_utils.gdal_merge.main(parameters)

the -separate flag puts each input in a different band in the raster, if you need the normal merge where the output raster has one band don't use it

parameters = ['', '-o', output_file_path] + input_files_path + ['-co', 'COMPRESS=LZW']
osgeo_utils.gdal_merge.main(parameters)