error de muestreo utilizando gdalwarp en un VRT subconjunto?

Question

Actualmente estoy intentando preparar un subconjunto de SRTM que incluya la conversión de altura de EGM96 a WGS84 utilizando gdalwarp. Me he encontrado con algunos desplazamientos extraños en parte de la imagen.
He aquí un ejemplo reproducible:
Prepare el archivo necesario para la conversión de altura:

cd /usr/share/proj
sudo wget https://download.osgeo.org/proj/vdatum/egm96_15/egm96_15.gtx
sudo chmod 644 egm96_15.gtx

Descargue dos mosaicos SRTM 3Sec y combínelos en un VRT:

wget http://srtm.csi.cgiar.org/wp-content/uploads/files/srtm_5x5/TIFF/srtm_39_03.zip
wget http://srtm.csi.cgiar.org/wp-content/uploads/files/srtm_5x5/TIFF/srtm_38_03.zip
gdalbuildvrt dem.vrt /vsizip/srtm_39_03.zip/srtm_39_03.tif \
/vsizip/srtm_38_03.zip/srtm_38_03.tif \
-te 9.497235 45.432812 13.183811 47.336189

Utilice gdalwarp para crear dos nuevos archivos, uno que es una simple conversión del contenido del archivo original y otro que contiene las alturas convertidas al elipsoide WGS84:

gdalwarp -r bilinear -of GTiff dem.vrt dem.tif
gdalwarp -r bilinear -of GTiff -s_srs EPSG:4326+5773 -t_srs EPSG:4326 \
dem.vrt dem_ellp.tif

Calcule su diferencia como GTiff y JPEG para obtener la diferencia EGM96-WGS84:

gdal_calc.py -A dem.tif -B dem_ellp.tif --outfile=dem_ellp_diff.tif --calc="A-B"
gdal_translate -of JPEG -scale -outsize 50% 50% \
dem_ellp_diff.tif dem_ellp_diff.jpg

El resultado es el siguiente:

Al oeste y al este, la imagen parece la esperada con simples líneas de degradado. Sin embargo, en el centro parece un archivo de aspecto DEM. Si se examinan más detenidamente los dos archivos, se puede ver un desplazamiento de píxeles hacia el Este introducido por gdalwarp en el segundo archivo (altura del elipsoide).
¿Alguien conoce este problema? No aparece cuando el MDE no es subconjunto durante la creación del archivo vrt, es decir, llamando a gdalbuildvrt sin el te argumento:

Además, las imágenes tienen el mismo aspecto cuando se utiliza un CRS espacial diferente, por ejemplo, UTM a través de
-t_srs EPSG:32632

El problema tampoco se presenta cuando se establece alternativamente la extensión a través de gdalwarp en lugar de gdalbuildvrt :

gdalbuildvrt dem2.vrt /vsizip/srtm_39_03.zip/srtm_39_03.tif \
/vsizip/srtm_38_03.zip/srtm_38_03.tif
gdalwarp -r bilinear -of GTiff -te 9.497235 45.432812 13.183811 47.336189 \
dem2.vrt dem2.tif
gdalwarp -r bilinear -of GTiff -s_srs EPSG:4326+5773 -t_srs EPSG:4326 \
-te 9.497235 45.432812 13.183811 47.336189 dem2.vrt dem2_ellp.tif
gdal_calc.py -A dem2.tif -B dem2_ellp.tif --outfile=dem2_ellp_diff.tif \
--calc="A-B"
gdal_translate -of JPEG -scale -outsize 50% 50% \ 
dem2_ellp_diff.tif dem2_ellp_diff.jpg

GDAL v2.2.3, Ubuntu 18.04, libproj v4.9.3-2

Answer 1

Acabo de probarlo de nuevo con la versión 2.4.0 de GDAL de ubuntugis unstable.
Aquí funciona perfectamente y al parecer se arregló en algún momento entre la versión 2.2.3 y la 2.4.0.