OS Terrain50 to GeoTIFF

osterrain50.py

# https://gis.stackexchange.com/a/462686
# CC-BY-SA 4.0

import os
import zipfile
from osgeo import gdal

gdal.UseExceptions()

input = "/path/to/terr50_gagg_gb.zip"
output = input.replace(".zip", ".tif")

srs = "EPSG:27700"

zips = (z for z in zipfile.ZipFile(input).namelist() if z.endswith(".zip"))
vrts = []

for zip in zips:
    basename = os.path.basename(zip)[0:4].upper()
    asc = "/vsizip/{/vsizip/{%s}/%s}/%s.asc" % (input, zip, basename)
    vrt = f"/vsimem/{basename}.vrt"
    gdal.Translate(vrt, asc, outputType=gdal.GDT_Float32)
    vrts.append(vrt)

vrt = f"/vsimem/terr50_gagg_gb.vrt"
gdal.BuildVRT(f"/vsimem/terr50_gagg_gb.vrt", vrts)
gdal.Translate(output, vrt, outputSRS=srs, creationOptions=["TILED=YES", "COMPRESS=LZW", "PREDICTOR=2"])

osterrain50_qgsprocessing.py

# QGIS Processing script for OS Terrain 50 ASCII grid zip file
# Add to Processing->Scripts https://docs.qgis.org/latest/en/docs/user_manual/processing/scripts.html

import os
import zipfile
from osgeo import gdal

gdal.UseExceptions()
from qgis import processing
from qgis.processing import alg
from qgis.core import QgsRasterLayer, QgsProcessingOutputRasterLayer, QgsProcessingException, QgsProcessingContext, QgsProcessingUtils


@alg(
    name="osterrain", label=alg.tr("Convert OS Terrain data"), 
    group="rasterconversion", group_label=alg.tr("Raster Conversion")
)
@alg.input(
    type=alg.FILE, name="INPUT_ZIP", 
    label="Input OS Terrain 50 ASCII Grid zip file", extension="zip"
)
@alg.input(
    type=alg.RASTER_LAYER_DEST, name='OUTPUT_RASTER',
    label='Output OS Terrain 50 Raster'
)
@alg.input(
    type=alg.STRING, name='ADDITIONAL_CREATE_OPTIONS',
    label='Additional GDAL driver creation options',
    default="COMPRESS=DEFLATE, TILED=YES, PREDICTOR=2"
)
@alg.output(
    type=alg.RASTER_LAYER, name="OUTPUT_LAYER", label="Output Raster Layer"
)
def osterrain(instance, parameters, context, feedback, inputs):
    """Convert OS Terrain data"""  # Do not delete this docstring
    input = instance.parameterAsString(parameters, "INPUT_ZIP", context)
    output_layer =  instance.parameterAsOutputLayer(parameters, 'OUTPUT_RASTER', context)
    output_file = str(output_layer)
    opts = instance.parameterAsString(parameters, "ADDITIONAL_CREATE_OPTIONS", context)

    srs = "EPSG:27700"
    opts = [opt.strip() for opt in opts.split(",")]
    zips = (z for z in zipfile.ZipFile(input).namelist() if z.endswith(".zip"))
    vrts = []

    feedback.pushInfo("Building list of ASCII Grids")
    
    # Conversion code based on 
    # https://gis.stackexchange.com/a/462686
    # CC-BY-SA 4.0

    for zip in zips:
        basename = os.path.basename(zip)[0:4].upper()
        asc = "/vsizip/{/vsizip/{%s}/%s}/%s.asc" % (input, zip, basename)
        vrt = f"/vsimem/{basename}.vrt"
        gdal.Translate(vrt, asc, outputType=gdal.GDT_Float32)
        vrts.append(vrt)
    feedback.pushInfo(f"{len(vrts)} ASCII Grids found")

    vrt = f"/vsimem/terr50_gagg_gb.vrt"
    gdal.BuildVRT(f"/vsimem/terr50_gagg_gb.vrt", vrts)
    feedback.pushInfo(f"Creating {output_layer}")
    gdal.Translate(output_layer, vrt, outputSRS=srs, creationOptions=opts)

    return {"OUTPUT_LAYER": output_layer}