StackOverflow2013

Note that there are some explanatory texts on larger screens.

plurals

PO
text
Body
copied!<p>If you want to restrict the second set of points to one of the tiles of the tessellation, you can use <code>tile.list</code> to have a description of each tile, and then check which points are in this tile (there are many functions to do so: in the following example, I use <code>secr::pointsInPolygon</code>).</p> <pre><code># Sample data x <- matrix( rnorm(20), nc = 2 ) y <- matrix( rnorm(1000), nc=2 ) # Tessellation library(deldir) d <- deldir(x[,1], x[,2]) plot(d, wlines="tess") # Pick a cell at random cell <- sample( tile.list(d), 1 )[[1]] points( cell$pt[1], cell$pt[2], pch=16 ) polygon( cell$x, cell$y, lwd=3 ) # Select the points inside that cell library(secr) i <- pointsInPolygon( y, cbind( c(cell$x,cell$x[1]), c(cell$y,cell$y[1]) ) ) points(y[!i,], pch=".") points(y[i,], pch="+") # Compute a tessellation of those points dd <- deldir(y[i,1], y[i,2]) plot(dd, wlines="tess", add=TRUE) </code></pre> <p><img src="https://i.stack.imgur.com/yT8Tg.png" alt="Tessellation inside a cell of another tessellation"></p> <p>If, instead, you want to translate and rescale the points to fit them into the tile, that is trickier.</p> <p>We need to somehow estimate how far away from the tile the points are: to this end, let us define a few auxilliary functions to compute, first the distance from a point to a segment, then the distance from a point to a polygon.</p> <pre><code>distance_to_segment <- function(M, A, B) { norm <- function(u) sqrt(sum(u^2)) lambda <- sum( (B-A) * (M-A) ) / norm(B-A)^2 if( lambda <= 0 ) { norm(M-A) } else if( lambda >= 1 ) { norm(M-B) } else { N <- A + lambda * (B-A) norm(M-N) } } A <- c(-.5,0) B <- c(.5,.5) x <- seq(-1,1,length=100) y <- seq(-1,1,length=100) z <- apply( expand.grid(x,y), 1, function(u) distance_to_segment( u, A, B ) ) par(las=1) image(x, y, matrix(z,nr=length(x))) box() segments(A[1],A[2],B[1],B[2],lwd=3) library(secr) distance_to_polygon <- function(x, poly) { closed_polygon <- rbind(poly, poly[1,]) if( pointsInPolygon( t(x), closed_polygon ) ) return(0) d <- rep(Inf, nrow(poly)) for(i in 1:nrow(poly)) { A <- closed_polygon[i,] B <- closed_polygon[i+1,] d[i] <- distance_to_segment(x,A,B) } min(d) } x <- matrix(rnorm(20),nc=2) poly <- x[chull(x),] x <- seq(-5,5,length=100) y <- seq(-5,5,length=100) z <- apply( expand.grid(x,y), 1, function(u) distance_to_polygon( u, poly ) ) par(las=1) image(x, y, matrix(z,nr=length(x))) box() polygon(poly, lwd=3) </code></pre> <p>We can now look for a transformation of the form </p> <pre><code>x --> lambda * x + a y --> lambda * y + b </code></pre> <p>that minimizes the (sum of the squared) distances to the polygon. That is actually not sufficient: we are likely to end up with scaling factor lambda equal to (or close to) zero. To avoid this, we can add a penalty if lambda is small.</p> <pre><code># Sample data x <- matrix(rnorm(20),nc=2) x <- x[chull(x),] y <- matrix( c(1,2) + 5*rnorm(20), nc=2 ) plot(y, axes=FALSE, xlab="", ylab="") polygon(x) # Function to minimize: # either the sum of the squares of the distances to the polygon, # if at least one point is outside, # or minus the square of the scaling factor. # It is not continuous, but (surprisingly) that does not seem to be a problem. f <- function( p ) { lambda <- log( 1 + exp(p[1]) ) a <- p[2:3] y0 <- colMeans(y) transformed_points <- t( lambda * (t(y)-y0) + a ) distances <- apply( transformed_points, 1, function(u) distance_to_polygon(u, x) ) if( all(distances == 0) ) - lambda^2 else sum( distances^2 ) } # Minimize this function p <- optim(c(1,0,0), f)$par # Compute the optimal parameters lambda <- log( 1 + exp(p[1]) ) a <- p[2:3] y0 <- colMeans(y) # Compute the new coordinates transformed_points <- t( lambda * (t(y)-y0) + a ) # Plot them segments( y[,1], y[,2], transformed_points[,1], transformed_points[,2], lty=3 ) points( transformed_points, pch=3 ) library(deldir) plot( deldir( transformed_points[,1], transformed_points[,2] ), wlines="tess", add=TRUE ) </code></pre> <p><img src="https://i.stack.imgur.com/Ym5VG.png" alt="Shifting and rescaling a set of points to put them inside a polygon"></p>

Querying!

Guidance

An individual column

Larger individual text columns get their own page to allow for proper reading.

SQuiL has stopped working due to an internal error.

If you are curious you may find further information in the browser console, which is accessible through the devtools (F12).

Reload