Lire une table à partir d'une géodatabase fichier ESRI (.gdb) à l'aide de R

21

J'essaie de lire une table directement d'une géodatabase fichier ESRI dans R. Un exemple de fichier de données peut être téléchargé ici . La base de données contient une classe d'entités ponctuelles (Zone9_2014_01_Broadcast) et deux tables liées (Zone9_2014_01_Vessel et Zone9_2014_01_Voyage). Vous pouvez lire le fichier de formes dans R à l'aide readOGRdu rgeospackage:

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")

Les deux tableaux liés indiquent également quand vous utilisez ogrListLayersou ogrInfo. Cependant, ogrInfodonne un avertissement:

Message d'avertissement: Dans ogrInfo ("Zone9_2014_01.gdb", layer = "Zone9_2014_01_Vessel"): ogrInfo: toutes les fonctionnalités NULL

Et si vous essayez d'utiliser readOGRsur les tables, vous obtenez une erreur:

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

Erreur dans readOGR (dsn = "Zone9_2014_01.gdb", layer = "Zone9_2014_01_Vessel"): aucune fonctionnalité trouvée En outre: Message d'avertissement: dans ogrInfo (dsn = dsn, layer = layer, encoding = encoding, use_iconv = use_iconv,: ogrInfo: toutes les fonctionnalités NULL

Ainsi, il apparaît que seules les caractéristiques géographiques peuvent être lues par readOGR. Existe-t-il un moyen d'importer les tables directement dans R ou est-ce la seule solution pour les exporter d'abord depuis ArcGIS en tant que fichiers * .dbf (ou * .txt) comme dans cette réponse?

En outre, si quelqu'un peut fournir des appels de R à un script python qui automatise l'exportation de fichiers * csv (de préférence) ou * .dbf, ce serait une solution de contournement acceptable. La solution doit simplement être évolutive et automatisée.

Coton.Rockwood
la source
2
Avez-vous vu la nouvelle intégration de R et ArcGIS? r-arcgis.github.io peut-être quelque chose d'utile pour votre travail.
Alex Tereshenkov
Merci pour la suggestion ... J'en avais vu la mention à un moment donné, mais je ne l'ai jamais examinée plus à fond. Ce serait peut-être le bon moment pour le faire!
Cotton.Rockwood
@AlexTereshenkov, si vous voulez rédiger une réponse courte pour cette solution, je l'accepterai car c'est ce que je cherchais.
Cotton.Rockwood
1
Il semble que le pont R-ArcGIS mentionné par @AlexTereshenkov ait la fonctionnalité de lire les tableaux directement dans R. L'intégration nécessite ArcGIS Desktop> 10.3.1 (ou ArcGIS Pro) et R> 3.2. La version 64 bits R ne peut être utilisée qu'avec le géotraitement en arrière-plan 64 bits (et autorise uniquement l'utilisation depuis ArcGIS, pas depuis R) ou ArcGIS Pro. Une fois les liaisons installées, vous pouvez utiliser le package arcgisbbindingen R. La fonction arc.open()ouvrira la table en tant que arc.dataset-class object. Pour ouvrir directement en tant que data.table, utilisez la fonction arc.select.
Cotton.Rockwood
bon à savoir. J'ai ajouté une réponse juste pour fermer le fil, mais vous avez tout compris par vous-même, alors acceptez mais ne votez pas: D
Alex Tereshenkov

Réponses:

18

Je suis un peu en retard à la fête, mais ça peut maintenant être lu par sf, avec

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

Il renvoie un avertissement (pas de géométrie d'entité présente) mais aussi un data.frame avec la table. Voir le fil de discussion qui a commencé ici: https://stat.ethz.ch/pipermail/r-sig-geo/2018-F February / 026344.html

Edzer Pebesma
la source
excellent! Merci Edzer ... content de voir ça et l'évolution de sf !!
Cotton.Rockwood
étrange, je n'ai pas pu faire ça sur 3 machines: j'obtiens une erreur, pas un avertissement?
Matifou
1
vous devrez installer la version de développement à partir de github, à partir de la source, ou attendre la sortie de la
version 0.6-1 le
Mieux vaut tard pour la fête que jamais! Je suis venu à cette fête il y a environ 2 ans et j'ai mis en œuvre l'une des solutions précédentes. Je suis juste allé chercher une sfsolution et Google m'a ramené avec joie à cette même soirée avec une solution super utile (j'ai donc ajouté mon vote positif à cette question).
D. Woods
9

J'utilise GDAL 2.0.2 qui est "livré" avec le support FDGB et sans tiers un pilote FGDB pour enquêter sur ce genre de choses. L'environnement de test est Debian Jessie 64 bits.

En bref, il semble que la "couche" Zone9_2014_01_Vesselcontient des données d'attribut pures et la couche Zone9_2014_01_Broadcastcontient des données de position. Vous pouvez utiliser une solution de contournement dans R via un appel système et la conversation de la GDB vers un conteneur de fichiers de formes (dernier script à la fin de la réponse).

Voici les étapes de l'enquête:

$ 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

Comme vous le voyez, le champ Geometryest défini sur None. Vous pouvez convertir les données dans un fichier de forme en utilisant ogr2ogret obtenir uniquement un fichier d'attribut dbase:

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

Zone9_2014_01_Vessel.dbf

Des géométries (positions) se trouvent dans la couche 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

Navire et voyage ne contenant aucune donnée de position selon le protocole des messages AIS .

Voici la solution de contournement complète dans R en utilisant un appel système pour que la GDB forme la conversation et le package foreignpour lire les 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
la source
merci @huckfinn! Il s'agit d'une bonne solution de contournement. J'ai pas mal de fichiers (dont beaucoup sont beaucoup plus gros que l'exemple), donc je vais essayer et voir comment la conversion en fichier de formes affecte les temps de traitement. J'espère aussi qu'il y a une solution directe dans R, mais si personne ne répond par une, je choisirai la vôtre comme réponse.
Cotton.Rockwood
3

Je ne sais pas si vous pouvez le faire avec readOGR mais essayez

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

Si cela ne fonctionne pas, essayez ogr2ogrdirectement, ce qui peut exporter des non-géométries vers la table. (Peut-être essayez le package R gdalUtilspour l'exécuter, une fois que vous avez terminé votre processus.)

mdsumner
la source
1
Malheureusement, readOGRn'a pas la capacité de lire les tables gdb.
Aaron
1
C'est probablement le cas maintenant.
mdsumner
Toujours pas au rgdal 1.2-8, gdal 2.0.1
gregmacfarlane
Cela s'appelle OpenFileGDB dans ogrDrivers () $ name, peut-être que vous essayez de lire un raster? Cela est toujours en cours d'implémentation, de toute façon si vous voulez savoir que vous pouvez poster une question avec des détails sur votre système et ce que vous avez essayé.
mdsumner
3

Il existe une intégration récemment publiée entre le R et ArcGIS d'Esri, appelée R ArcGIS Tools . Il fournit une intégration entre R et ArcGIS permettant d'accéder de manière interchangeable aux outils R et aux ressources ArcGIS. Vous devriez pouvoir accéder aux classes / tables d'entités de géodatabase avec cette intégration.

Des exemples d'outils R sont disponibles ici et des exemples d'outils illustrant l'utilisation de R dans les scripts de géotraitement sont ici .

Alex Tereshenkov
la source
1

Cette fonction personnalisée suit essentiellement le chemin tracé par @huckfinn mais utilise la gdalUtilsbibliothèque suggérée par @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)
}

Exécutez-le comme ceci:

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 vous ne l'avez pas déjà gdalinstallé, vous devrez l'installer pour y accéder gdalUtils. Vous pouvez trouver des binaires et des instructions pour l'installation de 'gdal' ici .

D. Woods
la source