¿Crear una cuadrícula Esri y rellenar los valores de las celdas utilizando Python/GDAL o ArcObjects?

Question

En términos de velocidad de procesamiento, ¿cuál es la forma más rápida de crear y luego añadir valores a una trama?

Los valores provendrán originalmente de archivos ASCII que pueden tener un tamaño de unos 500MB - idealmente me gustaría crear un GRID de Esri pero no es esencial.

Tengo Python (GDAL) y/o ArcOjects (9.3.1) en mi arsenal

No puedo utilizar la herramienta ASCII to Raster debido al formato de los valores en el archivo ASCII . Estoy buscando una solución para poder escribir los valores en el raster directamente.

Como ejemplo puedo tener un XYZ asccii como:

    @File_Version: 3
@Coordinate_Type_is: 3
@Export_Type_is: 1
@Number_of_Projects 1
@Project_Type_Name: 3 cstone25 
@Project_Unit_is: meters
#File_Version____________-> 3
#Project_Name____________-> cstone25
#Project_Type____________-> 3
#Project_Unit____________-> meters
#OpenWorks_Project_______-> EXPLORATION_UTM31
#Master_Project_______-> 
#Coordinate_type_________-> 3
#Number_of_points_in_hzd_-> 1
zpxc1_ufa_nra_env_sele+20+100_Max_Peak_Amp
#Horizon_internal_id_____->_120
#Horizon_extremes_are____->   7245798912.00000
#Horizon_onset_is__Minimum_____-> 1
#Horizon_type_is_Time______-> 1
#Horizon_color_is________-> 0 255 0
 6148.0 17392.0     529940.17623    6274994.18251 394123616.00000
 6148.0 17396.0     529940.83286    6275019.17388 326897280.00000
 6148.0 17400.0     529941.48949    6275044.16526 567459712.00000
 6148.0 17404.0     529942.14613    6275069.15663 429482368.00000
 6148.0 17408.0     529942.80276    6275094.14801 133182280.00000
 6148.0 17412.0     529943.45939    6275119.13938 32126824.00000
 6148.0 17416.0     529944.11603    6275144.13076 17940496.00000
 6148.0 17420.0     529944.77266    6275169.12213 320429248.00000
 6148.0 17424.0     529945.42929    6275194.11351 453848192.00000
 6148.0 17428.0     529946.08592    6275219.10488 336978496.00000
 6148.0 17432.0     529946.74256    6275244.09626 290245120.00000
 6148.0 17436.0     529947.39919    6275269.08763 323590976.00000
 6148.0 17440.0     529948.05582    6275294.07901 249979616.00000
 6148.0 17444.0     529948.71246    6275319.07038 221433360.00000
 6148.0 17448.0     529949.36909    6275344.06176 244277504.00000

o podría tener una cabecera que los valores:

!
!     FILE NAME :  Rob_Grid                
!     FORMATTED FILE CREATION DATE:  APR 25 2011
!     FORMATTED FILE CREATION TIME:  17:11
!
@Rob_Grid HEADER                         ,      GRID,         5
      15,  0.1000000E+31,            ,      7,      1
     245,    232,   262410.0    ,   273960.0    ,   6426990.    ,   6439190.    
    200.0000    ,   0.000000    ,   0.000000    
@
  0.1000000E+31  0.1000000E+31  0.1000000E+31  0.1000000E+31  0.1000000E+31
  0.1000000E+31  0.1000000E+31  0.1000000E+31  0.1000000E+31  0.1000000E+31
  0.1000000E+31  0.1000000E+31  0.1000000E+31   1023.936       1021.771    
   1020.476       1019.986       1019.994       1020.656       1022.158    
   1024.154       1026.193       1027.803       1028.656       1028.808    
   1028.548       1028.080       1027.858       1028.684       1030.780    

NOTA: esta pregunta se refiere a cómo crear y rellenar los valores de la trama lo más rápido posible

Answer 1

Construí un Python script con GDAL/NumPy que convierte algún formato ASCII personalizado (similar al tuyo) en un GeoTIFF para ArcGIS u otro software GIS. Todo lo que necesitas saber para hacer esto se proporciona en este tutorial de GDAL .

Básicamente, el script (que no puedo publicar) lee los datos en estructuras de datos NumPy, y hace lo siguiente:

1. Determina la tupla GeoTransformación a partir de los orígenes X e Y y las separaciones de la malla. Utilizo algo como (xdim.min() - dx/2, dx, 0.0, ydim.min() - dy/2, 0.0, dy) (ambos dx y dy son positivos, por lo que utilizo min() )
2. Estructurar las alturas del raster en un array NumPy con forma (len(ydim), len(xdim)) , a continuación, utilice WriteArray para cada banda
3. Si tiene coberturas rasterizadas incompletas, asegúrese de establecer SetNoDataValue tanto en la trama como en la matriz NumPy

Por último, como nota especial para ArcGIS: cuando vea sus datos, éstos pueden aparecer en negro, por lo que debe Calcular estadísticas . He añadido este paso en mi Python script también.

Answer 2

No estoy seguro de si es el el más rápido , pero podría considerar el uso de IRasterEdit .

La muestra aquí podría modificarse para escribir en un formato diferente (en lugar de TIFF), proporcionando el formato adecuado a IRasterWorkspace2.CreateRasterDataSet . Creo que también podría querer hacer un bucle en lugar de escribir todo el raster con una llamada a IRasterEdit.Write, como se hace en el ejemplo.

Creo que esto también se puede hacer utilizando IPixelBlockCursor pero los comentarios me dan la impresión de que se ha dejado de utilizar: "Una interfaz heredada en ArcGIS 9.2".