5

Edited: I have added login details to a demo SQL Server 2017.

I am querying an SQL Server from R using the FreeTDS driver (preferred over the Microsoft odbc driver as it supports Windows authentication in Linux).

I'm struggling to interpret the geometry column with st_read.

The connection is initiated as such:

library(odbc)
library(DBI)
library(sf)

username <- 'SO-user'
sql_server_ip <- '35.214.169.110'
password <- 'SQLaskingfortrouble?!'

con <- dbConnect(odbc(),
                 Driver = "FreeTDS",
                 Server = sql_server_ip,
                 Database = "stackoverflow-example",
                 UID = username,
                 PWD = password,
                 Port = 1433)
# I have also tried Driver = "ODBC Driver 17 for SQL Server"

Here's a query result:

data <- dbGetQuery(con, 'SELECT TOP 2
                         objectid, geom,
                         geom.STGeometryType() geom_type,
                         geom.STSrid STSrid
                         FROM spatialdata')

objectid    geom          geom_type     stsRID
864         blob[416 B]   LineString    25832
865         blob[416 B]   LineString    25832

class(data$geom)
"blob" "vctrs_list_of" "vctrs_vctr"

The raw/blob geometry looks like hex:

rawToHex(data$geom[1])
e86400000104180000004003780b38331f416b4df3ea9ecf5741a04c150c38331f4159ca32f09ecf574160ba...

Trying to use st_read, so I can get a useful sf dataframe out of it, gives me a range of error messages, depending on the options I pass to the function:

st_read(con, geometry_column = "geom", 
        query = "SELECT TOP 10 objectid, geom FROM \"spatialdata\"")

Error: nanodbc/nanodbc.cpp:1617: 00000: [FreeTDS][SQL Server]Invalid object name 'spatial_ref_sys'.  [FreeTDS][SQL Server]Statement(s) could not be prepared. 
<SQL> 'select srtext from spatial_ref_sys where srid = 70647808'
In addition: Warning message:
In CPL_crs_from_input(x) :
  GDAL Error 6: EPSG PCS/GCS code 70647808 not found in EPSG support files.  Is this a valid EPSG coordinate system?

It looks like the EPSG (25832) is not being recognised correctly.

What is the correct construction of st_read in my case?

Bonus question: How can a bounding box be efficiently appended to st_read?

fifthace
  • 506
  • 1
  • 10
  • 33
  • possible duplicate? https://gis.stackexchange.com/questions/57935/reading-ogc-geometry-data-from-sql-server-using-r – arcee123 Mar 30 '20 at 18:26
  • Thanks for this. However, this answer doesn't produce an sf object. It's also incredibly slow, because readWKT isn't vectorized. – fifthace Mar 30 '20 at 21:05
  • 1
    Could you share the data or some sort of [reproducible example](https://github.com/tidyverse/reprex) ? – agila Apr 02 '20 at 12:27

2 Answers2

3

Have you tried the WKB representation?

mytest <- st_read(con, geometry_column = "geom", 
                         query = "
select geometry::STGeomFromText('POLYGON ((0 0, 1.5 0, 1.5 1.5, 0 1.5, 0 0))', 0).STAsBinary() as geom
union all
select geometry::Parse('CIRCULARSTRING(0 0, 1 1.10, 2 2.3246, 0 7, -3 2.3246, -1 2.1082, 0 0)').STAsBinary() as geom;
")

plot(mytest)
lptr
  • 1
  • 2
  • 6
  • 16
0

The issue here is that you want to translate the format in which MSSQL stores geometry data to something that the sf package can understand. STAsBinary() should work.

data <- dbGetQuery(con, 'SELECT TOP 2
                     objectid, 
                     geom.STAsBinary(),
                     geom.geom_type,
                     geom.STSrid STSrid
                     FROM spatialdata')

Then st_read(data) should work

M. Leon
  • 11
  • 3