[R-sig-Geo] minimum area bounding rectangles

Barry Rowlingson b.rowlingson at lancaster.ac.uk
Wed Mar 16 00:47:31 CET 2011


Shazam:

minbb <- function(pts){

  ch = chull(pts)
  pts=pts[ch,]
  np = nrow(pts)
  pts=rbind(pts,pts[1,])
  minbba = Inf ; bbth = NA; rotmin = NA
  for(i in 1:np){
    th = pi-atan2(pts[i+1,2]-pts[i,2],pts[i+1,1]-pts[i,1])
    prot = rotxy(pts,th)
    bba = diff(range(prot[,1])) * diff(range(prot[,2]))
    if(bba < minbba){
      xyb=cbind(
        c(min(prot[,1]),max(prot[,1]),max(prot[,1]),min(prot[,1])),
        c(min(prot[,2]),min(prot[,2]),max(prot[,2]),max(prot[,2]))
        )
      xyb=rbind(xyb,xyb[1,])
      xyb= rotxy(xyb,-th)

      rotmin=prot
      minbba = bba
      bbth = th
    }
  }
  return(list(minbba=minbba,theta=th,pts=rotmin,box=xyb))

}

rotxy = function(pts,angle){
  co = cos(angle)
  si = sin(angle)
  cbind(co * pts[,1] - si * pts[,2], si * pts[,1] + co * pts[,2])
}

And usage:

 > pts=cbind(runif(60),runif(60))
 > bb=minbb(pts)
 > plot(bb$box,type="l")
 > points(pts)
 > pts2=rotxy(pts,pi/1.6)
 > bb2=minbb(pts2)
 > plot(bb2$box,type="l")
 > points(pts2)


> Thanks Barry for the great tips.  The MATLAB code will be a huge help -
> should be relatively straightforward to implement in R.

 Wrote this myself with a helpful rotate function nabbed from
spatstat. Could be optimised a bit, I reckon.

Barry



More information about the R-sig-Geo mailing list