¿Cómo solucionar el proceso de mosaico grande que está fallando?

Question

Tengo que hacer un mosaico de unos 550Gb de imágenes tif y el software que he probado sigue fallando. El área se ha dividido en zonas, de modo que la más pequeña tiene aproximadamente 200 mosaicos.

He utilizado las últimas versiones de ERDAS (Imagine y Mapper), ArcINFO y Global Mapper en un Intel Xeon E31245 de 3,30 gigahercios, DELL, 16GB RAM, Win 7 Professional de 64 bits. Mullti-core (4 en total), Hyper-threaded (8 en total) máquina. Mi C tiene 700GB libres y el D tiene 1,5TB.

Estoy estudiando la posibilidad de utilizar Grass (nunca lo he hecho antes) pero i.image.mosaic sólo parece manejar 4 archivos... algunos de los míos tienen 600 mosaicos. Cualquier otra opción o software de código abierto para probar?

Siento tener que añadir que no podemos utilizar un conjunto de datos en mosaico (o su equivalente en otro software), ya que necesitamos crear zonas con áreas sin datos definidas como ecw para que puedan abrirse en cualquier software SIG y combinarse con datos de menor resolución/antiguos cuando no existan nuevos datos sin problemas.

Un ejemplo de cómo se ven algunos archivos en mosaico en diferentes softwares. Global Mapper/ERDAS están bien pero no es correcto en arcgis.

--- INFORMACIÓN MÁS ANTIGUA---

Perdón por el dibujo tan tosco. Así que tener las áreas coloreadas como 5 zonas minimizará las áreas sin datos en la AOI más grande.

En arcgis el código es el siguiente (esto se ejecuta como un modelo y no en python ya que no puedo conseguir que tome la entrada tifList).

arcpy.MosaicToNewRaster_management(tifList+";" +mask,RootOutput,"Tile1.tif","PROJCS['GDA_1994_MGA_Zone_55',GEOGCS['GCS_GDA_1994',DATUM['D_GDA_1994',SPHEROID['GRS_1980',6378137.0,298.257222101]],PRIMEM['Greenwich',0.0],UNIT['Degree',0.0174532925199433]],PROJECTION['Transverse_Mercator'],PARAMETER['False_Easting',500000.0],PARAMETER['False_Northing',10000000.0],PARAMETER['Central_Meridian',147.0],PARAMETER['Scale_Factor',0.9996],PARAMETER['Latitude_Of_Origin',0.0],UNIT['Meter',1.0]]","16_BIT_UNSIGNED","0.5","3","MAXIMUM","#")

# Replace a layer/table view name with a path to a dataset (which can be a layer file) or create the layer/table view within the script
# The following inputs are layers or table views: "test2"

arcpy.CopyRaster_management(OutputFile,RootOutput+"Tile1b.tif","#","256","256","NONE","NONE","16_BIT_UNSIGNED")

donde tifList debería ser leído desde un archivo csv pero esto no funcionó en python así que estoy ejecutando lo anterior en un modelo en su lugar...

Tengo más de 1,5TB de espacio libre en mi disco, pero el proceso se bloquea con un error 9999.

¿Procesarían incluso 100 baldosas? -Es decir, ¿deberíamos considerar la posibilidad de dividir más las zonas?

Answer 1

Tendré que secundar las sugerencias de @blah238 de utilizar algún otro método de acceso a los datos que no sea la creación de una sola imagen en mosaico. Una simple suposición diría que no hay un escritorio por ahí que pueda manejar la cantidad de datos que tendría que procesar con el fin de mosaico todos esos azulejos.
Para desglosarlo, es probable que haya dos lugares en los que se está quedando sin espacio.

1. El primero va a estar en su memoria RAM. Con el fin de mosaico datos, usted está fusionando todo en un solo archivo, lo que significa idealmente, que todo el archivo se mantendrá en la memoria. Usted no tiene 550GB de RAM, lo que significa que habrá lectura/escritura desde la memoria virtual. Eso es suficiente para que el proceso se bloquee ahí mismo.
2. El otro problema es probablemente el espacio en el disco duro. Muchas de las operaciones de rasterización crean archivos temporales en el directorio "espacio de trabajo" que son bastante grandes. que son bastante grandes. Puede haber 2 o 3 o más de estos antes de que el conjunto de datos final de datos final, por lo que es posible utilizar todo el espacio de disco durante el procesamiento. espacio en el disco durante el procesamiento.

Ahora, otras soluciones. Como se mencionó en los comentarios anteriores, existe la opción de crear un Conjunto de datos de mosaico . Este conjunto de datos le permitirá no sólo tratar todos los mosaicos individuales como una única imagen sin fisuras, sino que también mantiene los metadatos sobre los mosaicos individuales que contiene. También le permite realizar operaciones de rasterización como Hillshade .

La otra opción que yo recomendaría, basándome en tu comentario de querer tener las zonas separadas sería crear un Catálogo de tramas . Un Catálogo Raster es esencialmente una capa de grupo. Se pueden añadir múltiples conjuntos de datos ráster a ella. Pueden ser gestionados en una geodatabase, e importar los rásteres, o simplemente crear un conjunto de datos no gestionado en el que el Catálogo de Rásteres mantiene las rutas a los conjuntos de datos ráster originales. Cuando se carga esta capa en ArcMap, se pueden establecer las propiedades de visualización para que sólo se cargue un determinado número de mosaicos ráster a la vez, o establecer la escala de visualización y la resolución.
Actualmente estoy utilizando un catálogo rasterizado para muestrear un conjunto de más de 100 GB de fotos aéreas. El rendimiento es muy bueno. Si usted está buscando un tipo diferente de almacenamiento de datos simplemente para los medios de gestión de un gran número de azulejos, entonces yo realmente lo recomiendo.

Este es el código que puede utilizar para crear un Catálogo Raster y luego importar un espacio de trabajo de azulejos en él:

import arcpy
import os

newdir = r"c:\data"
dbname = "Aerialphotos.gdb"
gdbdir = os.path.join(newdir, dbname)
rcat = "AerialCatalog"

arcpy.CreateRasterCatalog_management(gdbdir, rcat,
                                     "NAD 1983 StatePlane California VI FIPS 0406 (US Feet).prj", 
                                     "NAD 1983 StatePlane California VI FIPS 0406 (US Feet).prj",
                                     "MAX_FILE_SIZE_4GB", "1000", "3000", "9000",
                                     "UNMANAGED", "")

#Load all raster datasets in workspace to Raster Catalog
rcatdir = os.path.join(gdbdir, rcat)
rastertiledir = os.path.join(newdir, "tiles")

arcpy.WorkspaceToRasterCatalog_management(rastertiledir, rcatdir,
                                          "INCLUDE_SUBDIRECTORIES",
                                          "PROJECT_ONFLY")

Espero que esto ayude.

------------- Editar

Aquí hay un gráfico de los azulejos manejados por mi catálogo de trama. Tenga en cuenta que puede elegir que se muestren los wireframes o los datos rasterizados. El catálogo de tramas tiene una tabla de atributos a la que puede añadir campos, por ejemplo, si quiere añadir designaciones de zonas como en su gráfico. Entonces, podría elegir mostrar sólo los rásteres de una zona en particular.
Cuando se imprime en un gráfico desde la vista de diseño, se utiliza la resolución completa de los rásteres, por lo que no hay pérdida de calidad en la impresión.

Aquí está el mismo gráfico, pero mostrando algunos de los datos rasterizados, junto con algunos wireframes.

Answer 2

Sé que llego tarde a la fiesta. Pero aquí está mi sugerencia.

1) tamaño de la imagen
Si tus originales de 550 GB están sin comprimir, deberías convertirlos en archivos tiff comprimidos en jpeg. Guárdalos individualmente (no fusionados). Puedes comprimir usando arcgis, gdal, lo que quieras. La compresión te llevará a unos 23GB. No cree pirámides/vistas generales todavía. Para comprimir puedes usar cualquier programa de gis que te guste, pero a mi me gusta usar gdal así que el comando es básicamente este:

gdalwarp -multi -wm 512 --config  GDAL_CACHEMAX 512 -co compress=jpeg -co tiled=yes -co jpeg_quality=70 -co PHOTOMETRIC=YCBCR -co INTERLEAVE=PIXEL uncompressed.tif compressed.tif

Puedes hacer fácilmente un archivo bat que recorra todos tus tiffs descomprimidos. Me gusta usar gdalwarp para comprimir mis imágenes en lugar del habitual gdal_translate, porque es más rápido (usando la opción multi para los multinúcleos, y -wm para tener mucha memoria).

2) manipulación como una sola imagen
Puedes crear un mosaico "virtual" utilizando el formato gdal vrt. Esto es compatible con arcgis, qgis, mapserver, etc. No estoy seguro sobre el mapeador global y el mapinfo. El formato .vrt es sólo un único archivo xml con la lista de sus imágenes. Es un solo comando para crear:

gdalbuildvrt nameofmosaic.vrt compressed_tif_folder\*.tif

Este archivo tiene un tamaño de unos pocos kb.

3) visualización del exceso de velocidad
Hay que construir pirámides/vistas. Para ello, utilice su software preferido. Manteniendo las herramientas de gdal puedes:

gdaladdo -ro -r average --config COMPRESS_OVERVIEW JPEG --config JPEG_QUALITY_OVERVIEW 70 nameofmosaic.vrt 2 4 8 16 32 64 128

Esto llevará mucho tiempo. Prepárate para esperar de 2 a 3 días de procesamiento ininterrumpido.

4) utilizar el mosaico
Cargue el mosaico virtual en su programa GIS. Será rápido porque está leyendo los resúmenes que están en un solo archivo como un ecw. Cuando se acerca a la resolución real de las imágenes, entonces sólo se leerán las pocas visibles de las imágenes comprimidas, y eso también es muy rápido.

5) manejo de las áreas sin datos que se muestran en negro
Tienes 3 soluciones para esto: i) usar un formato de archivo que maneje nodatos, lo cual va a ser complicado; o ii) usar una banda alfa o iii) un archivo de máscara. Puedes crear una banda alfa automáticamente en el paso 2 diciéndole a GDAL que quieres que las áreas de nodatos estén en la banda alfa - sólo tienes que añadir la opción -addalpha:

gdalbuildvrt -addalpha nameofmosaic.vrt compressed_tif_folder\*.tif

El problema de las bandas alfa es que se comprimen mal. Así que tus panorámicas van a ser más grandes. Si eso te parece bien, entonces ya está hecho.

Si quiere crear un archivo de máscara, entonces es un poco más complicado. Y me parece que esto no encaja en la presente pregunta.

Así que espero que esto ayude. Puedes encontrar información sobre las herramientas gdal buscando en Google. Hay muchas cosas interesantes.

Answer 3

550gb de datos TIF de entrada son fácilmente manejados por un solo archivo ECW. Tenemos muchos clientes que comprimen conjuntos de datos mucho más grandes que esto, así que no piense que el formato no es capaz en esta área.

Su estrategia de dividir el proyecto en pequeños mosaicos para minimizar el área nula es también un buen enfoque para tomar con la versión de formato actual, ya que reducirá el tiempo de compresión

Su ejemplo incluye una referencia a datos de entrada de 16 bits sin signo. Yo recomendaría reescalar a 8bit si es posible (dependiendo de sus requisitos)

Por favor, explique por qué no ha podido procesar su proyecto con IMAGINE o ERMapper, ya que sin esta información no puedo ayudarle. O mejor aún, póngase en contacto con el equipo de soporte local

Tenga en cuenta que al utilizar el formato ESRI Mosaic Dataset, lo que las respuestas anteriores no mencionan es el requisito de generar la capa de pirámide/vista general. Sin ella, el rendimiento se verá considerablemente afectado. Es probable que pueda crear los archivos equivalentes de ECW en el mismo tiempo, pero con una calidad de imagen mejorada y unos requisitos de almacenamiento de salida significativamente menores.

Answer 4

Aunque está claro que es mejor utilizar una de las otras opciones mencionadas, podrías probar lo siguiente:

gdalbuildvrt index.vrt *.tif
gdal_translate -of "GTiff" -co "COMPRESS=LZW" -co "TILED=YES" -co "BIGTIFF=YES" index.vrt out.tif

Esto construye un formato virtual GDAL y luego convierte a un único GeoTiff.

Answer 5

Eso me suena bastante, también producimos grandes archivos ECW individuales de 500 a 1TB de archivos TIF. Pero yo no duraría en ArcGIS (ArcObjects y el motor de geoprocesamiento), ya que no es capaz de mosaico de esta cantidad de una manera fiable. Si usted quiere permanecer en el mundo de ESRI recomendaría a mosaico trozos de unos 50 GB o incluso más pequeños a la vez a un conjunto de datos raster almacenados en una base de datos de archivos. La herramienta de mosaico tiende a bloquearse después de un tiempo, así que es una buena idea dejar que ArcGIS libere memoria después de algunos GigaBytes.

Otra posibilidad es utilizar una geodatabase SDE de empresa o de grupo de trabajo. Con SDE se obtienen las antiguas herramientas de línea de comandos de SDE, que se basan en una robusta arquitectura de C++ y no en el poco fiable material de ArcObjects. Con el comando "sderaster -o mosaic...", puede hacer mosaicos a un RasterDataset hasta que su almacén de Base de Datos esté lleno. También hay comandos para construir pirámides y estadísticas para el RasterDataset, de lo contrario no es muy útil, porque la mayoría de los clientes no pueden mantener las imágenes en la memoria al leerlas, como blah238 mencionó anteriormente. Pero las pirámides (de hecho la indexación espacial) deberían resolver este problema.

Pero estas soluciones no le ayudan con MapInfo con seguridad. Usted mencionó que ya ha probado ERDAS Mapper. Esa es también la herramienta que yo preferiría. Ya hemos hecho un mosaico de 16000 archivos TIF, cada uno de 50 MB de tamaño en conjunto, que son 800GB.Luego lo comprimimos a un solo ECW con una relación de compresión de 1:20 que resultó en un archivo ECW de 30GB. Me pregunto si esto no funciona para usted...

Al menos todo el proceso se ejecutaba en un Pentium 4 de un solo núcleo a 1,6 GHz con 2 GB de RAM, por lo que el hardware no debería ser el problema. Estamos utilizando Windows Server 2003 (u otro sistema operativo de servidor) porque utiliza mejor los recursos de hardware. Tenga en cuenta que todo el proceso de compresión necesita mucho tiempo. Nuestra máquina estuvo trabajando unas 5 semanas en ese único archivo, y como a veces también se colapsaba tuvimos que hacerlo varias veces, pero al final conseguimos nuestro archivo ECW.

No conozco otro sistema o mecanismo para almacenar grandes cantidades de rásteres de una manera mayormente neutral para el vendedor. Todas las formas mencionadas anteriormente son muy específicas de ESRI. Por lo menos con Oracle RASTER y una implementación bastante similar en PostGIS, hay dos variantes de base de datos que son, bueno también no proveedor neutral, pero abierto a través de la interfaz SQL/MM.

Espero que esto ayude un poco.