Leer una tabla de una geodatabase de archivos ESRI (.gdb) usando R

21

Estoy tratando de leer una tabla directamente desde una geodatabase de archivos ESRI en R. Aquí se puede descargar un archivo de datos de ejemplo . La base de datos contiene una clase de entidad de puntos (Zone9_2014_01_Broadcast) y dos tablas vinculadas (Zone9_2014_01_Vessel y Zone9_2014_01_Voyage). Puede leer el archivo de forma en R usando readOGRel rgeospaquete:

library(rgeos)
library(downloader)

download("https://coast.noaa.gov/htdata/CMSP/AISDataHandler/2014/01/Zone9_2014_01.zip", dest="Zone9_2014_01.zip", mode="wb")
unzip("Zone9_2014_01.zip", exdir = ".")

#  Not Run (loads large point file)
#  broadcast <- readOGR(dsn = "Zone9_2014_01.gdb", layer = "Zone9_2014_01_Broadcast")

Las dos tablas vinculadas también muestran cuándo usa ogrListLayerso ogrInfo. Sin embargo, ogrInfoda una advertencia:

Mensaje de advertencia: En ogrInfo ("Zone9_2014_01.gdb", layer = "Zone9_2014_01_Vessel"): ogrInfo: todas las funciones NULL

Y si intentas usarlo readOGRen las tablas, obtienes un error:

vessel <- readOGR(dsn = "Zone9_2014_01.gdb", layer = "Zone9_2014_01_Vessel")

Error en readOGR (dsn = "Zone9_2014_01.gdb", layer = "Zone9_2014_01_Vessel"): no se encontraron características Además: Mensaje de advertencia: En ogrInfo (dsn = dsn, layer = layer, codificación = codificación, use_iconv = use_iconv,: ogrInfo: todas las características NULL

Por lo tanto, parece que solo readOGR puede leer las características geográficas. ¿Hay alguna forma de importar las tablas directamente en R o es la única solución para exportarlas primero desde ArcGIS como archivos * .dbf (o * .txt) como en esta respuesta?

Además, si alguien puede proporcionar llamadas de R a un script de Python que automatice la exportación de archivos * csv (preferiblemente) o * .dbf, sería una solución aceptable. La solución solo necesita ser escalable y automatizada.

Algodón Rockwood
fuente
2
¿Has visto la nueva integración de R y ArcGIS? r-arcgis.github.io quizás algo útil para su trabajo.
Alex Tereshenkov
Gracias por la sugerencia ... Había visto mención de ello en un momento, pero nunca lo examiné más a fondo. ¡Quizás ahora sea un buen momento para hacerlo!
Cotton.Rockwood
@AlexTereshenkov, si desea escribir una respuesta breve para esta solución, la aceptaré, ya que es lo que estaba buscando.
Cotton.Rockwood
1
Parece que el puente R-ArcGIS que mencionó @AlexTereshenkov tiene la funcionalidad de leer tablas directamente en R. La integración requiere ArcGIS Desktop> 10.3.1 (o ArcGIS Pro) y R> 3.2. R de 64 bits solo se puede usar con el geoprocesamiento de fondo de 64 bits (y solo permite el uso desde ArcGIS, no desde R) o ArcGIS Pro. Una vez que se instalan los enlaces, puede usar el paquete arcgisbbindingen R. La función arc.open()abrirá la tabla como un arc.dataset-class object. Para abrir directamente como a data.table, use la función arc.select.
Cotton.Rockwood
bueno saber. Agregué una respuesta solo para cerrar el hilo, pero has resuelto todo por tu cuenta, así que acepta pero no votes: D
Alex Tereshenkov

Respuestas:

18

Llego un poco tarde a la fiesta, pero ahora puedo leer esto sf, con

vessel <- sf::st_read(dsn = "Zone9_2014_01.gdb", layer = "Zone9_2014_01_Vessel")

Devuelve una advertencia (no hay geometrías de entidades presentes) sino también un data.frame con la tabla. Vea el hilo que comenzó aquí: https://stat.ethz.ch/pipermail/r-sig-geo/2018-February/026344.html

Edzer Pebesma
fuente
¡excelente! Gracias Edzer ... me alegro de ver esto y la evolución de sf !!
Cotton.Rockwood
extraño, no pude ejecutar esto en 3 máquinas: ¿aparece un error, no una advertencia?
Matifou
1
deberás instalar la versión de desarrollo desde github, desde la fuente, o esperar hasta la versión 0.6-1 el próximo mes
Edzer Pebesma
¡Más vale tarde a la fiesta que nunca! Vine a esta fiesta hace ~ 2 años y tuve una implementación de una de las soluciones anteriores. Acabo de buscar una sfsolución y Google felizmente me trajo de vuelta a esta misma fiesta con una solución muy útil (así que felizmente agregué mi voto a esta pregunta).
D. Woods el
9

Utilizo GDAL 2.0.2 que se "envía" con soporte FDGB y sin un tercero un controlador FGDB para investigar esas cosas. El entorno de prueba es Debian Jessie de 64 bits.

En resumen, parece que la "capa" Zone9_2014_01_Vesselcontiene datos de atributos puros y la capa Zone9_2014_01_Broadcastcontiene datos de posición. Puede usar una solución alternativa dentro de R a través de una llamada al sistema y la conversación del GDB a un contenedor de archivos de forma (último script al final de la respuesta).

Aquí están los pasos de investigación:

$ mkdir ~/dev.d/gis-se.d/gdb 
$ cd ~/dev.d/gis-se.d/gdb
$ wget https://coast.noaa.gov/htdata/CMSP/AISDataHandler/2014/01/Zone9_2014_01.zip
$ unzip Zone9_2014_01.zip
$ ogrinfo Zone9_2014_01.gdb Zone9_2014_01_Vessel | head -20
Had to open data source read-only.
INFO: Open of `Zone9_2014_01.gdb'
      using driver `OpenFileGDB' successful.

Layer name: Zone9_2014_01_Vessel
Geometry: None <---------------------------- HERE 
Feature Count: 1282
Layer SRS WKT:
(unknown)
FID Column = OID
MMSI: Integer (0.0)
IMO: Integer (0.0)
CallSign: String (255.0)
Name: String (255.0)
VesselType: Integer (0.0)
Length: Integer (0.0)
Width: Integer (0.0)
DimensionComponents: String (255.0)
OGRFeature(Zone9_2014_01_Vessel):1
  MMSI (Integer) = 367603345

Como puede ver, el campo Geometryestá establecido en None. Puede convertir los datos a un archivo de forma usando ogr2ogry obtener también solo un archivo de atributo dbase:

$ ogr2ogr -f 'ESRI SHAPEFILE' test Zone9_2014_01.gdb Zone9_2014_01_Vessel
$ ls test

Zone9_2014_01_Vessel.dbf

Las geometrías (posiciones) se pueden encontrar en la capa Zone9_2014_01_Broadcast.

$ ogr2ogr -f 'ESRI SHAPEFILE' test Zone9_2014_01.gdb
$ ls test

Zone9_2014_01_Broadcast.dbf  
Zone9_2014_01_Broadcast.shp  
Zone9_2014_01_Broadcast.prj  
Zone9_2014_01_Broadcast.shx  
Zone9_2014_01_Vessel.dbf
Zone9_2014_01_Voyage.dbf

Buque y viaje que no contienen datos de posición según el protocolo de mensajes AIS .

Aquí, la solución completa en R usando una llamada al sistema para que el GDB modele la conversación y el paquete foreignlea los dbf:

# Load module to get readOGR
require('rgdal');

# Load module to get read.dbf
require('foreign');

# goto the directory with the GDB stuff
setwd('~/dev.d/gis-se.d/gdb')

# Conversation to a shapefile container 
system("ogr2ogr -f 'ESRI SHAPEFILE' test Zone9_2014_01.gdb") 

# read the vessels
vessel <- read.dbf('test/Zone9_2014_01_Vessel.dbf');

# read hte voyage data
voyage <- read.dbf('test/Zone9_2014_01_Voyage.dbf');

# read the geometries in broad cast
broadcast <- readOGR('test/Zone9_2014_01_Broadcast.shp','Zone9_2014_01_Broadcast')

OGR data source with driver: ESRI Shapefile
Source: "test/Zone9_2014_01_Broadcast.shp", layer: "Zone9_2014_01_Broadcast"
with 1639274 features
It has 10 fields

# is vessel OK?    
head(vessel)

MMSI IMO CallSign Name VesselType Length Width   DimensionC
1 367603345  NA     <NA> <NA>         50     20     6     7,13,3,3
2 563000574  NA     <NA> <NA>         70    276    40 188,88,20,20
3 367449580  NA     <NA> <NA>         31     28    10     9,19,5,5
4 367302503  NA     <NA> <NA>         31     20     8     8,12,4,4
5 304003909  NA     <NA> <NA>         71    222    32 174,48,21,11
6 210080027  NA     <NA> <NA>         71    294    32 222,72,22,10

# is voyage OK?
head(voyage)

VoyageID           Destinatio Cargo Draught        ETA  StartTime    EndTime      MMSI
1       12                 KAKE    50      20       <NA> 2014-01-01       <NA> 367603345
2       23             YOKOHAMA    70     125 2014-01-11 2014-01-01 2014-01-30 563000574
3       38         KETCHIKAN AK    31      40 2014-11-12 2014-01-01       <NA> 367449580
4       52 CLARENCE STRAIT LOGS    31      30 2014-09-12 2014-01-01       <NA> 367302503
5       62               JP TYO    71      90 2014-01-13 2014-01-01 2014-01-31 304003909
6       47           VOSTOCHNYY    71     106 2014-01-13 2014-01-01       <NA> 210080027
huckfinn
fuente
gracias @huckfinn! Esta es una buena solución. Tengo bastantes archivos (muchos de los cuales son mucho más grandes que el ejemplo), así que lo probaré y veré cómo la conversión a un shapefile afecta los tiempos de procesamiento. También tengo esperanzas, hay una solución directa en R, pero si nadie responde con una, seleccionaré la suya como respuesta.
Cotton.Rockwood
3

No estoy seguro si puede hacer esto con readOGR pero intente

vessel <- readOGR(dsn = "Zone9_2014_01.gdb", layer = "Zone9_2014_01_Vessel", dropNULLGeometries = FALSE)

Si eso no funciona, intente ogr2ogrdirectamente, que puede exportar no geometrías a la tabla. (Tal vez intente con el paquete R gdalUtilspara ejecutar eso, una vez que haya finalizado su proceso).

mdsumner
fuente
1
Desafortunadamente, readOGRno tiene la capacidad de leer tablas gdb.
Aaron
1
Probablemente sí, ahora.
mdsumner
Todavía no a partir de rgdal 1.2-8, gdal 2.0.1
gregmacfarlane
Se llama OpenFileGDB en ogrDrivers () $ name, ¿tal vez estás intentando leer un ráster? Eso todavía se está implementando, de cualquier manera, si desea averiguar, puede publicar una pregunta con detalles sobre su sistema y lo que ha intentado.
mdsumner
3

Existe una integración lanzada recientemente entre R y ArcGIS de Esri, llamada R ArcGIS Tools . Proporciona integración entre R y ArcGIS, lo que permite acceder indistintamente a las herramientas R y a los recursos de ArcGIS. Debería poder acceder a las clases / tablas de entidades de geodatabase con esta integración.

Las herramientas de muestra R están disponibles aquí y las herramientas de muestra que ilustran el uso de R en scripts de geoprocesamiento están aquí .

Alex Tereshenkov
fuente
1

Esta función personalizada básicamente sigue la ruta descrita por @huckfinn pero usa la gdalUtilsbiblioteca sugerida por @mdsumner.

read_GDB_Layer <- function(dsn, layerName, overwrite = T) {
   conversionDir <- tempdir()

   gdalUtils::ogr2ogr(src_datasource_name = dsn, 
                      dst_datasource_name = conversionDir, 
                      f = "ESRI Shapefile", layer = layerName, 
                      verbose = T, overwrite = overwrite)

   df <- foreign::read.dbf(file.path(conversionDir, paste0(layerName, ".dbf")))

   return(df)
}

Ejecútelo así:

vsl <- read_GDB_Layer(dsn = "Zone9_2014_01.gdb", layerName = "Zone9_2014_01_Vessel")
vyg <- read_GDB_Layer(dsn = "Zone9_2014_01.gdb", layerName = "Zone9_2014_01_Voyage")

Si aún no lo ha gdalinstalado, deberá instalarlo para poder acceder gdalUtils. Puede encontrar binarios e instrucciones para la instalación de 'gdal' aquí .

D. Woods
fuente