Trazando un gran número de objetos en un gráfico de uso de R y formato WKT

Question

Yo soy la recuperación de miles de geometría de líneas en formato WKT de un Servidor Microsoft SQL server base de datos. Me gustaría trazar todas las líneas en un mapa, pero estoy teniendo algunas dificultades con el código. ¿Cómo debo modificar el código siguiente para llegar a mi meta?

library(RODBC)

png(file="examplex.png", width=600, height=480)
setwd("C:/ArcR")

con <- odbcDriverConnect('driver={SQL Server};server=SERVERNAME;database=DBNAME;trusted_connection=true')

objects_1 <- sqlQuery(con, 'SELECT TOP (1) Shape.STAsText() as ShapeWKT FROM TABLENAME ;')

things <- vector("list", 1)

z = 0
for(line in objects_1$ShapeWKT)
{
  {
    things[[z+1]]<-readWKT(line)
  }
  z = z + 1
}

plot(things[[1]]) 

dev.off()

Sí, el código hace referencia a sólo un registro devuelto, pero yo quería empezar la pregunta con el código en la orden de trabajo, para aquellos que no saben la respuesta, pero todavía son capaces de pruebas de posibilidades.

Answer 1

Desde readWKT no leer más de un objeto a un tiempo, usted todavía tiene que recorrer cada fila de la dataframe regresó de sqlQuery, pero usted puede hacerlo todo desde dentro de R en lugar de las llamadas de sistema de ogr2ogr. He aquí un ejemplo de uso de SpatialPoints:

library(RODBC)
library(rgeos)

con <- odbcDriverConnect('driver={SQL Server}; server=SERVERNAME; database=DBNAME; trusted_connection=true')
objects_1 <- sqlQuery(con, 'SELECT Shape.STAsText() AS ShapeWKT, OBJECTID FROM TABLE;')

Establecer un SpatialPointsDataFrame:

point.sp <- SpatialPointsDataFrame(readWKT(objects_1$ShapeWKT[1]), data=data.frame(OBJECTID=objects_1$OBJECTID[1]))

Iterar a través de las filas restantes.

for (n in 2:length(objects_1$OBJECTID)) {
  point.sp <- rbind(point.sp, 
                    SpatialPointsDataFrame(readWKT(objects_1$ShapeWKT[n]), 
                    data.frame(OBJECTID=objects_1$OBJECTID[n])))
}

Carecen de su conjunto de datos, me encontré con algo de código de ejemplo:

obj_wkt <- c("POINT (30 10)", "POINT (30 11)", "POINT (31 11)", "POINT (31 10)")
obj_id <- 1:4
objects_1 <- data.frame(ShapeWKT=obj_wkt, OBJECTID=obj_id)

plot(point.sp, axes = 1)

Usted debe configurar el inicial SpatialPointsDataFrame a ser la longitud de la objects_1 (con ceros) o esto ralentizará la iteración aún más.

Answer 2

Esta respuesta requiere que los bombardeos a cabo ogr2ogr. Estoy seguro de que tiene que haber una mejor manera, pero por ahora, esto funciona.

library(RODBC)

xmin <- xxxxxxx  #insert your number here
xmax <- xxxxxxx  #insert your number here
ymin <- xxxxxx   #insert your number here
ymax <- xxxxxx   #insert your number here

png(file="examplex.png", width=600, height=480)
setwd("C:/ArcR")

con <- odbcDriverConnect('driver={SQL     Server};server=SERVERNAME;database=DBNAME;trusted_connection=true')

objects_1 <- sqlQuery(con, 'SELECT Shape.STAsText() AS ShapeWKT, OBJECTID FROM TABLE;')

cat("WKT;OBJECTID",file="outfile.csv",sep="\n")

for(i in 1:nrow(objects_1)) {
  row <- objects_1[i,]

  cat('"',file="outfile.csv",append=TRUE)
  cat(toString(row$ShapeWKT),file="outfile.csv",append=TRUE)
  cat('";',file="outfile.csv",append=TRUE)
  cat(row$OBJECTID,file="outfile.csv",sep="\n",append=TRUE)
}

shell('ogr2ogr -f "ESRI Shapefile" mydata.shp "C:\\ArcR\\Outfile.csv"')

Inspections.mp <- readShapeLines("mydata")

#I put a shapefile called 'boundary' in the working directory.  Its not important, but its how I made my results.
border.mp <- readShapePoly("Boundary")

plot.new()
par(mar=c(0,0,0,0))
plot(border.mp, col="#000000FF", xlim = c(xmin, xmax), ylim = c(ymin, ymax), border="#000000FF")
lines(Inspections.mp, col="blue",lwd=0.25)

dev.off()

Estoy bastante seguro de que estoy haciendo mal uso de 'gato' aquí, pero eh... eh.