Cómo usar gdal2tiles en una imagen tiff personalizada recibida de proveedores para generar mosaicos

10

He estado luchando un poco para generar mosaicos para una imagen de alta resolución que tenemos. La imagen actual que tenemos es una imagen muy grande (+ 20GB), guardada como un archivo GeoTiff.Gran imagen de GTiff

Me gustaría generar los mosaicos utilizando la utilidad de línea de comando gdal2tiles y luego abrirlos y verlos en Cesium, utilizando el proveedor de imágenes TMS para proporcionar los mosaicos. Usando gdalinfo, aquí están algunos de los detalles de la imagen:

Driver: GTiff/GeoTIFF
Files: image.tif
Size is 52250, 56119
Coordinate System is:
PROJCS["WGS 84 / UTM zone 35S",
    GEOGCS["WGS 84",
        DATUM["WGS_1984",
            SPHEROID["WGS 84",6378137,298.257223563,
                AUTHORITY["EPSG","7030"]],
            AUTHORITY["EPSG","6326"]],
        PRIMEM["Greenwich",0],
        UNIT["degree",0.0174532925199433],
        AUTHORITY["EPSG","4326"]],
    PROJECTION["Transverse_Mercator"],
    PARAMETER["latitude_of_origin",0],
    PARAMETER["central_meridian",27],
    PARAMETER["scale_factor",0.9996],
    PARAMETER["false_easting",500000],
    PARAMETER["false_northing",10000000],
    UNIT["meters",1],
    AUTHORITY["EPSG","32735"]]
Origin = (606276.000000000000000,7197873.000000000000000)
Pixel Size = (0.500000000000000,-0.500000000000000)
Metadata:
  AREA_OR_POINT=Area
  TIFFTAG_MAXSAMPLEVALUE=13165
  TIFFTAG_MINSAMPLEVALUE=1
  TIFFTAG_RESOLUTIONUNIT=2 (pixels/inch)
  TIFFTAG_SOFTWARE=ERDAS IMAGINE
  TIFFTAG_XRESOLUTION=1
  TIFFTAG_YRESOLUTION=1
Image Structure Metadata:
  INTERLEAVE=PIXEL
Corner Coordinates:
Upper Left  (  606276.000, 7197873.000) ( 28d 3'21.59"E, 25d19'55.12"S)
Lower Left  (  606276.000, 7169813.500) ( 28d 3'29.55"E, 25d35' 7.17"S)
Upper Right (  632401.000, 7197873.000) ( 28d18'55.92"E, 25d19'47.60"S)
Lower Right (  632401.000, 7169813.500) ( 28d19' 5.85"E, 25d34'59.57"S)
Center      (  619338.500, 7183843.250) ( 28d11'13.23"E, 25d27'27.58"S)
Band 1 Block=512x512 Type=UInt16, ColorInterp=Gray
Band 2 Block=512x512 Type=UInt16, ColorInterp=Undefined
Band 3 Block=512x512 Type=UInt16, ColorInterp=Undefined
Band 4 Block=512x512 Type=UInt16, ColorInterp=Undefined

Mi primer intento fue usar gdal_translate para georreferenciar la imagen, y luego usar gdalwarp para cambiar la proyección a EPSG: 3857, como lo requiere Cesium (consulte la referencia de la API)

gdal_translate -of VRT -a_srs EPSG:4326 -gcp 606275 7197875 28.055987 -25.331974 -gcp 606275 7169814 28.058200 -25.585326 -gcp 632400.5 7197875 28.31553 -25.329876 -gcp 632400.5 7169814 28.318286 -25.583209 image.tif newImage1.vrt
gdalwarp -of VRT -t_srs EPSG:3857 newImage1.vrt newImage2.vrt

Sin embargo, obtengo muchos de los siguientes errores:

ERROR 1: latitud o longitud excedieron los límites

Otro método que probé fue usar gdal2tiles directamente y generar los mosaicos:

gdal2tiles.py image.tif

Esto creó una carpeta con una subcarpeta (etiquetada 18) como el único nivel de zoom en el que se crearon los mosaicos. Sin embargo, las imágenes que obtengo aquí están completamente "equivocadas" y "borrosas".

Un ejemplo de uno de los mosaicos:

ingrese la descripción de la imagen aquí

¿Alguna sugerencia para generar mosaicos para esta imagen de imagen grande de un área específica usando gdal2tiles para que pueda cargarla y verla en cesio?

Actualizar

Entonces, después de probar la sugerencia de @ iant, utilicé los siguientes comandos:

gdalwarp -co TILED=YES -co COMPRESS=DEFLATE -co BIGTIFF=YES -t_srs EPSG:3857 image.tif newImage.tif

Este comando funcionó perfectamente bien hasta el final donde obtuve el siguiente error:

ERROR 1: TIFFFillTile: error de lectura en la fila 43520, col 47104; obtuvo 35788250 bytes, esperado 37421449

No estoy seguro de qué significaba este error, lo dejé por el momento y todavía obtuve una imagen final "newImage.tif", producida por el paso gdalwarp. Usando esto llamé a gdal2tiles.py

gdal2tiles.py newImage.tif

Esto produjo una carpeta con subcarpetas 10-18 (y no solo un nivel de zoom 18 como obtuve anteriormente). También se lee perfectamente en Cesium, sin ningún error de consola, pero la imagen todavía se ve "incorrecta":

Imagen cargada en cesio

Estoy considerando que mi problema puede ser como @ user30184 ha sugerido "... los datos de origen no son adecuados para gdal2tiles". Sin embargo, hasta que nuestro proveedor pueda proporcionarnos algo para usar con gdal, esto es todo lo que tengo.

Estaba considerando eliminar una de las bandas para evitar que Gdal interfiriera en la última banda como un canal alfa. ¿Alguna sugerencia?

esfuerzo
fuente
¿Por qué quieres georreferenciar la imagen? Ya tiene toda la información de CRS dentro.
AndreJ

Respuestas:

7

Creo que todo lo que necesitas hacer es reproyectarlo usando:

gdalwarp -co TILED=YES -co COMPRESS=DEFLATE -t_srs EPSG:3857 newImage.tif image.tif

y luego en mosaico:

gdal2tiles.py newImage.tif

Si su archivo es muy grande, demorará un poco.

Ian Turton
fuente
Puede controlar los niveles de zoom con el parámetro -z, ¿ya lo intentó? Y tenga en cuenta que su imagen tiene 4 bandas que pueden generar resultados inesperados, especialmente porque maneja datos de 16 bits. Es posible que primero se necesite algo de preprocesamiento
usuario 30184
Gracias por tu respuesta @iant. Voy a intentar esto y ver qué pasa. ¿Podrías explicar un poco más en tu respuesta qué significan las opciones que has seleccionado? Según la documentación, la opción -co "pasa una opción de creación al controlador de formato de salida". Entonces, ¿estás agregando más propiedades al archivo tiff?
esfuerzo
@ user30184 No intenté eso todavía no. ¿Cómo sabría cuáles deberían ser los niveles de zoom? ¿O puedo especificarlos como quiera? Pensé que al dejar fuera esta opción, dejaría que el script determinara el nivel de zoom en función del área seleccionada.
esfuerzo
@ user30184 también mencionó que podría obtener resultados inesperados con una imagen de 4 bandas y datos de 16 bits. ¿Por qué exactamente? ¿No está en el formato correcto para el procesamiento de mosaico gdal? Si es así, la imagen se obtuvo directamente del proveedor, ¿cuáles serían los pasos para obtener el archivo en el formato correcto? A saber, un archivo tiff GDAL? (Si puedo decir algo así)
esfuerzo
Tengo los siguientes problemas al ejecutar gdalwarp como sugirió @iant: >>>>>>>>>>>> ERROR 1: TIFFFillTile: Error de lectura en la fila 43520, col 47104; obtuvo 35788250 bytes, esperado 37421449 ERROR 1: TIFFReadEncodedTile () falló. ERROR 1: pleiades_merge05m_2015-06-19.tif, banda 1: IReadBlock falló en X offset 86, Y offset 109 ERROR 1: GetBlockRef falló en X block offset 86, Y block offset 109 >>>>>>>>>>> > ¿Alguna sugerencia para arreglar esto?
esfuerzo
2

Supongo que su imagen es uno de los productos de 4 bandas de Airbus DS:

http://www.intelligence-airbusds.com/en/4951-which-spectral-mode-do-i-choose

Gdal2tiles está hecho para dividir imágenes visuales comunes en mosaicos png. Dichas imágenes usan 8 bits por banda y tienen una banda (escala de grises), 3 bandas (rojo-verde-azul) de 4 bandas (reg-verde-azul + alfa).

Diría que su pregunta es irrelevante en gran medida porque sus datos de origen no son adecuados para gdal2tiles. Puede solucionar los problemas inmediatos que tiene en este momento, pero el resultado final aún no será bueno si no vuelve a procesar sus datos.

La razón del mosaico no atractivo que adjuntó a su pregunta puede ser que la cuarta banda de datos se interpreta como un canal alfa.

usuario30184
fuente
Gracias @ user30184, he estado leyendo sobre algunos recursos y he llegado a un pensamiento similar. Creo que sería mejor pedirles a nuestros proveedores que nos proporcionen archivos tiff "compatibles con GDAL", pero hasta que nos respondan, esto es todo lo que tenemos. Estaba considerando eliminar una de las bandas para evitar que Gdal interfiriera en la última banda como un canal alfa. ¿Alguna sugerencia?
esfuerzo
Usar gdal_translate para cortar un pequeño subconjunto de la imagen gdal_translate -srcwin 20000 20000 1000 1000 original.tif sample.tifdebería hacerlo (grandes compensaciones para evitar las áreas de nodata). Abra esta pequeña imagen con QGIS y podrá jugar con la configuración de visualización rápidamente. Mi suposición sobre el canal alfa probablemente sea incorrecta, de lo contrario el resultado debería verse colorido, no gris.
user30184
Gracias @ user30184, hice lo que sugirió y lo abrí con éxito en QGIS. Consulte este enlace: drive.google.com/open?id=0B97NtaPJrVz-anRYQmxjZFludk0 ¿Cómo puedo ahora "depurar" mi problema? ¿Usando QGIS para hacer gdalwarp y gdal2tiles?
esfuerzo