Advertisement
Guest User

Minimum Bounding Rectangle

a guest
May 1st, 2018
115
0
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
Julia 2.21 KB | None | 0 0
  1. module MBR
  2.  
  3. using Plots
  4.  
  5. """
  6.    rightturn(a, b, c)
  7.  
  8. Do `a`, `b`, and `c` constitute a right turn?
  9. """
  10. rightturn(a, b, c) = det([(b - a) (c - a)]) < 0
  11.  
  12. """
  13.    convexhull(points)
  14.  
  15. `points` is given as a vector of vectors.
  16. Returns an open convex polygon.
  17. """
  18. function convexhull(points)
  19.     n = length(points)
  20.     sorted = sort(points, lt = (p, q) -> p[1] < q[1] || (p[1] == q[1] && p[2] < q[2]))
  21.     hull = [sorted[1], sorted[2]]
  22.     for i = 3:n
  23.         while length(hull) >= 2 && !rightturn(hull[end-1], hull[end], sorted[i])
  24.             pop!(hull)
  25.         end
  26.         push!(hull, sorted[i])
  27.     end
  28.     lower = [sorted[end], sorted[end-1]]
  29.     for i = n-2:-1:1
  30.         while length(lower) >= 2 && !rightturn(lower[end-1], lower[end], sorted[i])
  31.             pop!(lower)
  32.         end
  33.         push!(lower, sorted[i])
  34.     end
  35.     append!(hull, lower[2:end-1])
  36. end
  37.  
  38. """
  39.    mbr(hull)
  40.  
  41. `hull` is a closed convex polygon, given as a vector of vectors.
  42. Returns a tuple of (area, corner, edgevector1, edgevector2).
  43. """
  44. function mbr(hull)
  45.     function rectangle(dx)
  46.         dy = [dx[2], -dx[1]]
  47.         xs = [dot(p, dx) for p in hull[2:end]]
  48.         ys = [dot(p, dy) for p in hull[2:end]]
  49.         min = [minimum(xs), minimum(ys)]
  50.         max = [maximum(xs), maximum(ys)]
  51.         dev = max - min
  52.         (abs(dev[1] * dev[2]), dx * min[1] + dy * min[2], dx * dev[1], dy * dev[2])
  53.     end
  54.     best = (Inf,)
  55.     for i = 2:length(hull)
  56.         rec = rectangle(normalize(hull[i] - hull[i-1]))
  57.         if rec[1] < best[1]
  58.             best = rec
  59.         end
  60.     end
  61.     best
  62. end
  63.  
  64. """
  65.    test(n)
  66.  
  67. Generate a random test with `n` points and show it on a plot.
  68. """
  69. function test(n)
  70.     points = [[rand(), rand()] for _ = 1:n]
  71.     scatter(map(p -> p[1], points), map(p -> p[2], points), label = "Data points",
  72.             xlims = (-0.5, 1.5), ylims = (-0.5, 1.5), aspect_ratio = :equal)
  73.     hull = convexhull(points)
  74.     push!(hull, hull[1])
  75.     plot!(map(p -> p[1], hull), map(p -> p[2], hull), label = "Convex hull")
  76.     rec = mbr(hull)
  77.     rec = [rec[2], rec[2] + rec[3], rec[2] + rec[3] + rec[4], rec[2] + rec[4], rec[2]]
  78.     plot!(map(p -> p[1], rec), map(p -> p[2], rec), label = "MBR")
  79. end
  80.  
  81. end
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement