trovare il punto di intersezione in R

Ho 2 vettori:

set.seed(1) x1 = rnorm(100,0,1) x2 = rnorm(100,1,1) 

Voglio tracciare questi come linee e quindi trovare i punti di intersezione delle linee, anche se ci sono più punti di intersezione, quindi voglio localizzarli ognuno.

inserisci la descrizione dell'immagine qui

Mi sono imbattuto in una domanda simile e ho cercato di risolvere questo problema usando spatstat , ma non ero in grado di convertire la mia cornice dati combinata contenente entrambi i valori vettoriali in psp object .

Se hai letteralmente solo due vettori di numeri casuali, puoi usare una tecnica piuttosto semplice per ottenere l’intersezione di entrambi. Trova tutti i punti in cui x1 è sopra x2 e quindi sotto di esso nel punto successivo o viceversa. Questi sono i punti di intersezione. Quindi usa le rispettive pendenze per trovare l’intercetta per quel segmento.

 set.seed(1) x1=rnorm(100,0,1) x2=rnorm(100,1,1) # Find points where x1 is above x2. above<-x1>x2 # Points always intersect when above=TRUE, then FALSE or reverse intersect.points<-which(diff(above)!=0) # Find the slopes for each line segment. x1.slopes<-x1[intersect.points+1]-x1[intersect.points] x2.slopes<-x2[intersect.points+1]-x2[intersect.points] # Find the intersection for each segment. x.points<-intersect.points + ((x2[intersect.points] - x1[intersect.points]) / (x1.slopes-x2.slopes)) y.points<-x1[intersect.points] + (x1.slopes*(x.points-intersect.points)) # Plot. plot(x1,type='l') lines(x2,type='l',col='red') points(x.points,y.points,col='blue') 

inserisci la descrizione dell'immagine qui

Ecco un codice di intersezione segmento segmento alternativo,

 # segment-segment intersection code # http://paulbourke.net/geometry/pointlineplane/ ssi <- function(x1, x2, x3, x4, y1, y2, y3, y4){ denom <- ((y4 - y3)*(x2 - x1) - (x4 - x3)*(y2 - y1)) denom[abs(denom) < 1e-10] <- NA # parallel lines ua <- ((x4 - x3)*(y1 - y3) - (y4 - y3)*(x1 - x3)) / denom ub <- ((x2 - x1)*(y1 - y3) - (y2 - y1)*(x1 - x3)) / denom x <- x1 + ua * (x2 - x1) y <- y1 + ua * (y2 - y1) inside <- (ua >= 0) & (ua <= 1) & (ub >= 0) & (ub <= 1) data.frame(x = ifelse(inside, x, NA), y = ifelse(inside, y, NA)) } # do it with two polylines (xy dataframes) ssi_polyline <- function(l1, l2){ n1 <- nrow(l1) n2 <- nrow(l2) stopifnot(n1==n2) x1 <- l1[-n1,1] ; y1 <- l1[-n1,2] x2 <- l1[-1L,1] ; y2 <- l1[-1L,2] x3 <- l2[-n2,1] ; y3 <- l2[-n2,2] x4 <- l2[-1L,1] ; y4 <- l2[-1L,2] ssi(x1, x2, x3, x4, y1, y2, y3, y4) } # do it with all columns of a matrix ssi_matrix <- function(x, m){ # pairwise combinations cn <- combn(ncol(m), 2) test_pair <- function(i){ l1 <- cbind(x, m[,cn[1,i]]) l2 <- cbind(x, m[,cn[2,i]]) pts <- ssi_polyline(l1, l2) pts[complete.cases(pts),] } ints <- lapply(seq_len(ncol(cn)), test_pair) do.call(rbind, ints) } # testing the above y1 = rnorm(100,0,1) y2 = rnorm(100,1,1) m = cbind(y1, y2) x = 1:100 matplot(x, m, t="l", lty=1) points(ssi_matrix(x, m)) 

inserisci la descrizione dell'immagine qui

Risposta tardiva, ma qui c’è un metodo “spaziale” che utilizza il pacchetto SP e RGEOS. Ciò richiede che sia x che y siano numerici (o che possano essere convertiti in numerici). La proiezione è arbitraria, ma l’epsg: 4269 sembrava funzionare bene:

 library(sp) library(rgeos) # dummy x data x1 = rnorm(100,0,1) x2 = rnorm(100,1,1) #dummy y data y1 <- seq(1, 100, 1) y2 <- seq(1, 100, 1) # convert to a sp object (spatial lines) l1 <- Line(matrix(c(x1, y1), nc = 2, byrow = F)) l2 <- Line(matrix(c(x2, y2), nc = 2, byrow = F)) ll1 <- Lines(list(l1), ID = "1") ll2 <- Lines(list(l2), ID = "1") sl1 <- SpatialLines(list(ll1), proj4string = CRS("+init=epsg:4269")) sl2 <- SpatialLines(list(ll2), proj4string = CRS("+init=epsg:4269")) # Calculate locations where spatial lines intersect int.pts <- gIntersection(sl1, sl2, byid = TRUE) int.coords <- [email protected] # Plot line data and points of intersection plot(x1, y1, type = "l") lines(x2, y2, type = "l", col = "red") points(int.coords[,1], int.coords[,2], pch = 20, col = "blue")