module MBR

using LinearAlgebra
using Plots

"""
    rightturn(a, b, c)

Do `a`, `b`, and `c` constitute a right turn?
"""
rightturn(a, b, c) = det([(b - a) (c - a)]) < 0

"""
    convexhull(points)

`points` is given as a vector of vectors.
Returns an open convex polygon.
"""
function convexhull(points)
    n = length(points)
    sorted = sort(points, lt = (p, q) -> p[1] < q[1] || (p[1] == q[1] && p[2] < q[2]))
    hull = [sorted[1], sorted[2]]
    for i = 3:n
        while length(hull) >= 2 && !rightturn(hull[end-1], hull[end], sorted[i])
            pop!(hull)
        end
        push!(hull, sorted[i])
    end
    lower = [sorted[end], sorted[end-1]]
    for i = n-2:-1:1
        while length(lower) >= 2 && !rightturn(lower[end-1], lower[end], sorted[i])
            pop!(lower)
        end
        push!(lower, sorted[i])
    end
    append!(hull, lower[2:end-1])
end

"""
    mbr(hull)

`hull` is a closed convex polygon, given as a vector of vectors.
Returns a tuple of (area, corner, edgevector1, edgevector2).
"""
function mbr(hull)
    function rectangle(dx)
        dy = [dx[2], -dx[1]]
        xs = [dot(p, dx) for p in hull[2:end]]
        ys = [dot(p, dy) for p in hull[2:end]]
        min = [minimum(xs), minimum(ys)]
        max = [maximum(xs), maximum(ys)]
        dev = max - min
        (abs(dev[1] * dev[2]), dx * min[1] + dy * min[2], dx * dev[1], dy * dev[2])
    end
    best = (Inf,)
    for i = 2:length(hull)
        rec = rectangle(normalize(hull[i] - hull[i-1]))
        if rec[1] < best[1]
            best = rec
        end
    end
    best
end

"""
    test(n)

Generate a random test with `n` points and show it on a plot.
"""
function test(n)
    points = [[rand(), rand()] for _ = 1:n]
    scatter(map(p -> p[1], points), map(p -> p[2], points), label = "Data points",
            xlims = (-0.5, 1.5), ylims = (-0.5, 1.5), aspect_ratio = :equal)
    hull = convexhull(points)
    push!(hull, hull[1])
    plot!(map(p -> p[1], hull), map(p -> p[2], hull), label = "Convex hull")
    rec = mbr(hull)
    rec = [rec[2], rec[2] + rec[3], rec[2] + rec[3] + rec[4], rec[2] + rec[4], rec[2]]
    plot!(map(p -> p[1], rec), map(p -> p[2], rec), label = "MBR")
end

end