Skip to content

readOGR and PostGIS: "Error getting field 0"

8 messages · Mathieu Basille, Roger Bivand

#
Dear R users,

I have a working PostGIS database that I routinely access through R, using 
readOGR or readGDAL. I have this point layer, however, that I can't read in R:

 > ctrdcoast <- readOGR("PG:dbname=florida host=localhost user=pguser", 
layer = "gis.world_ctrdcoast")
OGR data source with driver: PostgreSQL
Source: "PG:dbname=florida host=localhost user=pguser", layer: 
"gis.world_ctrdcoast"
with 35398 features and 2 fields
Error in readOGR("PG:dbname=florida host=localhost user=pguser", layer = 
"gis.world_ctrdcoast") :
   Error getting field 0

The layer is basically a table with just 2 fields: geom (a point feature) 
and gid (a serial primary key). I can export this layer without any problem 
using ogr2ogr, and the layer looks fine to me (it displays fine in QGIS):

ogr2ogr -f "ESRI Shapefile" ctrd.shp PG:"dbname=florida host=localhost 
user=pguser" "gis.world_ctrdcoast"

Any idea where this is going wrong? Thanks in advance for any suggestion.
Mathieu.


 > sessionInfo()
R version 2.15.1 (2012-06-22)
Platform: x86_64-pc-linux-gnu (64-bit)

locale:
  [1] LC_CTYPE=fr_FR.UTF-8       LC_NUMERIC=C 
LC_TIME=fr_FR.UTF-8
  [4] LC_COLLATE=fr_FR.UTF-8     LC_MONETARY=fr_FR.UTF-8 
LC_MESSAGES=fr_FR.UTF-8
  [7] LC_PAPER=C                 LC_NAME=C                  LC_ADDRESS=C 

[10] LC_TELEPHONE=C             LC_MEASUREMENT=fr_FR.UTF-8 
LC_IDENTIFICATION=C

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base

other attached packages:
[1] raster_2.0-31 rgdal_0.7-22  sp_1.0-2

loaded via a namespace (and not attached):
[1] compiler_2.15.1 grid_2.15.1     lattice_0.20-6  tools_2.15.1


# SELECT PostGIS_Full_Version();
 
     postgis_full_version 

-------------------------------------------------------------------------
  POSTGIS="2.1.0SVN r10597" GEOS="3.3.3-CAPI-1.7.4" PROJ="Rel. 4.7.1, 23 
September 2009" GDAL="GDAL 2.0dev, released 2011/12/29" LIBXML="2.8.0" 
LIBJSON="UNKNOWN" TOPOLOGY RASTER
#
On Fri, 30 Nov 2012, Mathieu Basille wrote:

            
The error message is from ogrReadColumn(), line 221 in src/ogrsource.cpp, 
for example here:

https://r-forge.r-project.org/scm/viewvc.php/pkg/src/ogrsource.cpp?view=markup&root=rgdal

What does ogrInfo() report? Does it make a difference if you add an extra 
field? I'm wondering whether in this case the GID is not being seen as a 
field as such, but the FID of the point. Possibly ogr2ogr steps round 
this. I don't have a PostGIS system running, so may need help debugging. 
If you debug(readOGR), what happens here:

         if (drop_unsupported_fields) {
              iflds <- as.integer((1:ogr_info$nitems)-1)
              iflds <- iflds[keep]
              fldnms <- ogr_info$iteminfo$name[keep]
              if (any(!keep)) warning(paste("Fields dropped:",
                  paste(ogr_info$iteminfo$name[!keep], collapse=" ")))
         } else {
              if (any(!keep)) stop(paste("Unsupported field type:",
                  paste(ogr_info$iteminfo$typeName[!keep], collapse=" ")))
              iflds <- as.integer((1:ogr_info$nitems)-1)
              fldnms <- ogr_info$iteminfo$name
         }

I think ogr_info$nitems is 1, so we get 0, but if ogr_info$nitems is 0, 
(0, -1).

Roger

  
    
#
Dear Roger,

Thanks a lot for your very fast answer! Here are the requested outputs:

Le 30/11/2012 10:39, Roger Bivand a ?crit :
$ ogrinfo "PG:dbname=florida host=localhost user=pguser" 
gis.world_ctrdcoast -summary
INFO: Open of `PG:dbname=florida host=localhost user=pguser'
       using driver `PostgreSQL' successful.

Layer name: gis.world_ctrdcoast
Geometry: Unknown (any)
Feature Count: 35398
Extent: (-179.916500, -59.083000) - (179.968500, 83.535000)
Layer SRS WKT:
(unknown)
FID Column = gid
Geometry Column = geom

I noticed there was a call to ogrInfo in the readOGR function. After this 
line, here is what's reported:

Browse[2]> ogr_info
Source: "PG:dbname=florida host=localhost user=pguser", layer: 
"gis.world_ctrdcoast"
Driver: PostgreSQL number of rows 35398
Feature type: wkbPoint with 2 dimensions
NA
Number of fields: 0
[1] name     type     length   typeName
<0 rows> (or 0-length row.names)
As a matter of fact, it does! If I add a serial field named 'bla', it works 
perfectly:

 > ctrdcoast <- readOGR("PG:dbname=florida host=localhost user=pguser", 
layer = "gis.world_ctrdcoast")
OGR data source with driver: PostgreSQL
Source: "PG:dbname=florida host=localhost user=pguser", layer: 
"gis.world_ctrdcoast"
with 35398 features and 1 fields
Feature type: wkbPoint with 2 dimensions

 > summary(ctrdcoast)
Object of class SpatialPointsDataFrame
Coordinates:
                 min      max
coords.x1 -179.9165 179.9685
coords.x2  -59.0830  83.5350
Is projected: NA
proj4string : [NA]
Number of points: 35398
Data attributes:
    Min. 1st Qu.  Median    Mean 3rd Qu.    Max.
       1    8850   17700   17700   26550   35400

 > head(ctrdcoast at data)
   bla
1   1
2   2
3   3
4   4
5   5
6   6
It looks like a correct guess, given the output of ogrinfo. In addition, 
the exported shapefile from ogr2ogr present one field (named FID), while if 
I load the table into QGIS, I can see the field with its correct name 
(gid). Admittedly, while I'm pretty sure that QGIS uses GDAL to load 
tables, I have no idea how it deals with it.
Browse[2]> ogr_info$nitems
[1] 0
Browse[2]> iflds
[1]  0 -1
Browse[2]> fldnms
character(0)
Does that mean that we necessarily need a field in the table (in addition 
to geom + gid)? It's not hard to 'fix' by adding a nonsensical field, but I 
find it quite counterintuitive.

Thanks again for the help!
Mathieu.

  
    
#
On Fri, 30 Nov 2012, Mathieu Basille wrote:

            
Thanks very much, this is the information needed, that there are no 
fields. It appears that at least in this driver, readOGR() needs to be 
able to return a geometry-only object, or to create a GID column in the 
data.frame. I'll let you know when I've committed to R-forge, and will be 
very grateful if you can check out the source package, install, and test 
whether the changes have helped. It may take a little time to get there, 
but I'll give it priority.

Roger

  
    
#
Le 30/11/2012 11:27, Roger Bivand a ?crit :
I'd be happy to test when it's ready! However, no need to hurry: in the 
meantime, it is still possible (and easy) to add a dummy field to make 
readOGR work.

Thank you very much for you work!
Mathieu.

  
    
#
On Fri, 30 Nov 2012, Mathieu Basille wrote:

            
OK, thanks. SVN revision 389 committed to R-forge project rgdal, checkout 
from:

svn checkout svn://scm.r-forge.r-project.org/svnroot/rgdal/pkg

and rename pkg -> rgdal (to avoid getting multiple branches), then try
R CMD check rgdal, if OK INSTALL, and try your problem.

Roger

  
    
#
Roger,

It works beautifully! Checkout and compilation went successfully without 
any problem, and now readOGR can read layers with no field. As a matter of 
fact, it doesn't even need a primary key, just the geometry field. If there 
is a primary key, it is used as 'FID', if there is no primary key, a 'FID' 
field is created starting at 0.

This is perfect, thanks a lot!

Mathieu.


Le 30/11/2012 11:56, Roger Bivand a ?crit :

  
    
#
On Fri, 30 Nov 2012, Mathieu Basille wrote:

            
Thanks for checking, I'll put a release submission together.

Roger