julien
julien

Reputation: 113

Create numerous lines in Simple Features from list of coordinates in R

I am trying to detect whether pairs of objects (trees) are separated by roads or lie on the same side of them. I have downloaded my road network and think I more or less understand how to use st_intersects. So all I am missing are line segments between the pairs of trees I am considering in order to test intersections with the roads..

However, I cannot seem to figure out how to create lines between my objects. I have a large number of pairs (300K+), so must be able to do this programmatically, whereas all the examples I am finding seem to be "hand coded".

Suppose the following two matrices, containing the coordinates of the "origin" and "destination" of each pair.

    orig = matrix(runif(20),ncol=2)
    dest = matrix(runif(20),ncol=2)

In this example, I need to create 10 lines: one between orig[1,] and dest[1,], another (distinct) one between orig[2,] and dest[2,], etc. My understanding is that I should be using st_multilinestring, but I cannot figure out how to formulate the call. Typically, I either end up with "XYZM" points, or with a multi-segment line starting at orig[1,] and terminating at dest[10,] after going through all other coordinates. And when it is not one of these outcomes, it is a whole host of errors.

Is st_multilinestring what I should be using and if so, how does one do this? Thanks!!

Upvotes: 1

Views: 636

Answers (2)

SymbolixAU
SymbolixAU

Reputation: 26258

Here's a way to construct the sfc / sf object using library(sfheaders)

library(sf)
library(sfheaders)

## If you add a pseudo-id column
orig <- cbind( orig, 1:nrow( orig ) )
dest <- cbind( dest, 1:nrow( dest ) )

## you can bind these matrices together
m <- rbind( orig, dest )

## set the order by the 'id' column
m <- m[ order( m[,3] ), ]

## then use `sfheaders` to create your linestrings

sfc <- sfheaders::sfc_linestring(
  obj = m
  , linestring_id = 3 ## 3rd column
)

sfc

# Geometry set for 10 features 
# geometry type:  LINESTRING
# dimension:      XY
# bbox:           xmin: 0.01952919 ymin: 0.04603703 xmax: 0.9172785 ymax: 0.9516615
# epsg (SRID):    NA
# proj4string:    NA
# First 5 geometries:
# LINESTRING (0.7636528 0.2465392, 0.05899529 0.7...
# LINESTRING (0.6435893 0.9158161, 0.01952919 0.1...
# LINESTRING (0.05632407 0.3106372, 0.03306822 0....
# LINESTRING (0.1978259 0.07432209, 0.2907429 0.0...
# LINESTRING (0.1658199 0.6436758, 0.1407145 0.75...

Upvotes: 1

Spacedman
Spacedman

Reputation: 94192

Loop over rows of your origin and destination matrices using lapply and create a vector of LINESTRING objects:

> lines = do.call(st_sfc,
              lapply(
                1:nrow(orig),
                function(i){
                   st_linestring(
                      matrix(
                        c(orig[i,],dest[i,]), ncol=2,byrow=TRUE)
                   )
                }
               )
              )

This gives you this:

> lines
Geometry set for 10 features 
geometry type:  LINESTRING
dimension:      XY
bbox:           xmin: 0.06157865 ymin: 0.007712881 xmax: 0.967166 ymax: 0.9864812
epsg (SRID):    NA
proj4string:    NA
First 5 geometries:
LINESTRING (0.6646269 0.1545195, 0.8333102 0.40...
LINESTRING (0.5588124 0.5166538, 0.3213998 0.08...
LINESTRING (0.06157865 0.6138778, 0.06212246 0....
LINESTRING (0.202455 0.4883115, 0.5569435 0.986...
LINESTRING (0.3120373 0.8189916, 0.8499419 0.73...

Let's check we got all that the right way round. Where's the fourth line come from and going to?

> orig[4,]
[1] 0.2024550 0.4883115
> dest[4,]
[1] 0.5569435 0.9864812

which looks like the coordinates in the fourth LINESTRING output.

You can then st_intersects this with another set of features and see which of these cross them.

(You might also need to add the coordinate system to them...)

Upvotes: 0

Related Questions