Convex hull: Difference between revisions
(→=={{header|Racket}}==: stub added) |
|||
Line 130: | Line 130: | ||
=={{header|Racket}}== |
=={{header|Racket}}== |
||
Also an implementation of https://en.wikibooks.org/wiki/Algorithm_Implementation/Geometry/Convex_hull/Monotone_chain (therefore kinda {{trans|Go}} |
|||
<lang racket> |
<lang racket>#lang typed/racket |
||
(define-type Point (Pair Real Real)) |
|||
</lang> |
|||
(define-type Points (Listof Point)) |
|||
(: ⊗ (Point Point Point -> Real)) |
|||
(define/match (⊗ o a b) |
|||
[((cons o.x o.y) (cons a.x a.y) (cons b.x b.y)) |
|||
(- (* (- a.x o.x) (- b.y o.y)) (* (- a.y o.y) (- b.x o.x)))]) |
|||
(: Point<? (Point Point -> Boolean)) |
|||
(define (Point<? a b) |
|||
(cond [(< (car a) (car b)) #t] [(> (car a) (car b)) #f] [else (< (cdr a) (cdr b))])) |
|||
;; Input: a list P of points in the plane. |
|||
(define (convex-hull [P:unsorted : Points]) |
|||
;; Sort the points of P by x-coordinate (in case of a tie, sort by y-coordinate). |
|||
(define P (sort P:unsorted Point<?)) |
|||
;; for i = 1, 2, ..., n: |
|||
;; while L contains at least two points and the sequence of last two points |
|||
;; of L and the point P[i] does not make a counter-clockwise turn: |
|||
;; remove the last point from L |
|||
;; append P[i] to L |
|||
;; TB: U is identical with (reverse P) |
|||
(define (upper/lower-hull [P : Points]) |
|||
(reverse |
|||
(for/fold ((rev : Points null)) |
|||
((P.i (in-list P))) |
|||
(let u/l : Points ((rev rev)) |
|||
(match rev |
|||
[(list p-2 p-1 ps ...) #:when (not (positive? (⊗ p-2 P.i p-1))) (u/l (list* p-1 ps))] |
|||
[(list ps ...) (cons P.i ps)]))))) |
|||
;; Initialize U and L as empty lists. |
|||
;; The lists will hold the vertices of upper and lower hulls respectively. |
|||
(let ((U (upper/lower-hull (reverse P))) |
|||
(L (upper/lower-hull P))) |
|||
;; Remove the last point of each list (it's the same as the first point of the other list). |
|||
;; Concatenate L and U to obtain the convex hull of P. |
|||
(append (drop-right L 1) (drop-right U 1)))) ; Points in the result will be listed in CCW order.) |
|||
(module+ test |
|||
(require typed/rackunit) |
|||
(check-equal? |
|||
(convex-hull |
|||
(list '(16 . 3) '(12 . 17) '(0 . 6) '(-4 . -6) '(16 . 6) '(16 . -7) '(16 . -3) '(17 . -4) |
|||
'(5 . 19) '(19 . -8) '(3 . 16) '(12 . 13) '(3 . -4) '(17 . 5) '(-3 . 15) '(-3 . -9) |
|||
'(0 . 11) '(-9 . -3) '(-4 . -2) '(12 . 10))) |
|||
(list '(-9 . -3) '(-3 . -9) '(19 . -8) '(17 . 5) '(12 . 17) '(5 . 19) '(-3 . 15))))</lang> |
|||
{{out}} |
{{out}} |
Revision as of 22:56, 11 January 2017
Find the points which form a convex hull from a set of arbitrary two dimensional points.
For example, given the points (16,3), (12,17), (0,6), (-4,-6), (16,6), (16,-7), (16,-3), (17,-4), (5,19), (19,-8), (3,16), (12,13), (3,-4), (17,5), (-3,15), (-3,-9), (0,11), (-9,-3), (-4,-2) and (12,10) the convex hull would be (-9,-3), (-3,-9), (19,-8), (17,5), (12,17), (5,19) and (-3,15).
- See also
Go
<lang go>package main
import ( "fmt" "image" "sort" )
// ConvexHull returns the set of points that define the
// convex hull of p in CCW order starting from the left most.
func (p points) ConvexHull() points {
// From https://en.wikibooks.org/wiki/Algorithm_Implementation/Geometry/Convex_hull/Monotone_chain
// with only minor deviations.
sort.Sort(p)
var h points
// Lower hull for _, pt := range p { for len(h) >= 2 && !ccw(h[len(h)-2], h[len(h)-1], pt) { h = h[:len(h)-1] } h = append(h, pt) }
// Upper hull for i, t := len(p)-2, len(h)+1; i >= 0; i-- { pt := p[i] for len(h) >= t && !ccw(h[len(h)-2], h[len(h)-1], pt) { h = h[:len(h)-1] } h = append(h, pt) }
return h[:len(h)-1] }
// ccw returns true if the three points make a counter-clockwise turn func ccw(a, b, c image.Point) bool { return ((b.X - a.X) * (c.Y - a.Y)) > ((b.Y - a.Y) * (c.X - a.X)) }
type points []image.Point
func (p points) Len() int { return len(p) } func (p points) Swap(i, j int) { p[i], p[j] = p[j], p[i] } func (p points) Less(i, j int) bool { if p[i].X == p[j].X { return p[i].Y < p[i].Y } return p[i].X < p[j].X }
func main() { pts := points{ {16, 3}, {12, 17}, {0, 6}, {-4, -6}, {16, 6}, {16, -7}, {16, -3}, {17, -4}, {5, 19}, {19, -8}, {3, 16}, {12, 13}, {3, -4}, {17, 5}, {-3, 15}, {-3, -9}, {0, 11}, {-9, -3}, {-4, -2}, {12, 10}, } hull := pts.ConvexHull() fmt.Println("Convex Hull:", hull) }</lang>
- Output:
Convex Hull: [(-9,-3) (-3,-9) (19,-8) (17,5) (12,17) (5,19) (-3,15)]
Haskell
<lang Haskell> import Data.Ord;import Data.List;(x,y)=((!!0),(!!1)); compareFrom o l r=compare((x l-x o)*(y r-y o))((y l-y o)*(x r-x o)); distanceFrom from to = ((x to-x from)**2+(y to-y from)**2)**(1/2); convexHull points =
let{o = minimum points; presorted=sortBy(compareFrom o)(filter(/=o)points); collinears=groupBy(((EQ==).).compareFrom o)presorted; outmost=map(maximumBy(comparing(distanceFrom o)))collinears; }in dropConcavities[o]outmost;
dropConcavities (left:lefter) (right:righter:rightest) =
case compareFrom left right righter of{ LT -> dropConcavities (right:left:lefter) (righter:rightest); EQ -> dropConcavities (left:lefter) (righter:rightest); GT -> dropConcavities lefter (left:righter:rightest); };dropConcavities output lastInput = lastInput++output;
example=convexHull[[16,3],[12,17],[0,6],[-4,-6],[16,6],[16,-7],[16,-3],
[17,-4],[5,19],[19,-8],[3,16],[12,13],[3,-4],[17,5],[-3,15],[-3,-9], [0,11],[-9,-3],[-4,-2],[12,10]]
</lang>
J
Restated from the implementation at http://kukuruku.co/hub/funcprog/introduction-to-j-programming-language-2004 which in turn is a translation of http://dr-klm.livejournal.com/42312.html
<lang J>counterclockwise =: ({. , }. /: 12 o. }. - {.) @ /:~ crossproduct =: 11"_ o. [: (* +)/ }. - {. removeinner =: #~ 1, 0 > 3 crossproduct\ ], 1: hull =: [: removeinner^:_ counterclockwise</lang>
Example use:
<lang J> hull 16j3 12j17 0j6 _4j_6 16j6 16j_7 16j_3 17j_4 5j19 19j_8 3j16 12j13 3j_4 17j5 _3j15 _3j_9 0j11 _9j_3 _4j_2 12j10 _9j_3 _3j_9 19j_8 17j5 12j17 5j19 _3j15</lang>
Python
An approach that uses the shapely library:
<lang python>from __future__ import print_function from shapely.geometry import MultiPoint
if __name__=="__main__": pts = MultiPoint([(16,3), (12,17), (0,6), (-4,-6), (16,6), (16,-7), (16,-3), (17,-4), (5,19), (19,-8), (3,16), (12,13), (3,-4), (17,5), (-3,15), (-3,-9), (0,11), (-9,-3), (-4,-2), (12,10)]) print (pts.convex_hull)</lang>
- Output:
POLYGON ((-3 -9, -9 -3, -3 15, 5 19, 12 17, 17 5, 19 -8, -3 -9))
Racket
Also an implementation of https://en.wikibooks.org/wiki/Algorithm_Implementation/Geometry/Convex_hull/Monotone_chain (therefore kinda
<lang racket>#lang typed/racket (define-type Point (Pair Real Real)) (define-type Points (Listof Point))
(: ⊗ (Point Point Point -> Real)) (define/match (⊗ o a b)
[((cons o.x o.y) (cons a.x a.y) (cons b.x b.y)) (- (* (- a.x o.x) (- b.y o.y)) (* (- a.y o.y) (- b.x o.x)))])
(: Point<? (Point Point -> Boolean)) (define (Point<? a b)
(cond [(< (car a) (car b)) #t] [(> (car a) (car b)) #f] [else (< (cdr a) (cdr b))]))
- Input
- a list P of points in the plane.
(define (convex-hull [P:unsorted : Points])
;; Sort the points of P by x-coordinate (in case of a tie, sort by y-coordinate). (define P (sort P:unsorted Point<?)) ;; for i = 1, 2, ..., n: ;; while L contains at least two points and the sequence of last two points ;; of L and the point P[i] does not make a counter-clockwise turn: ;; remove the last point from L ;; append P[i] to L ;; TB: U is identical with (reverse P) (define (upper/lower-hull [P : Points]) (reverse (for/fold ((rev : Points null)) ((P.i (in-list P))) (let u/l : Points ((rev rev)) (match rev [(list p-2 p-1 ps ...) #:when (not (positive? (⊗ p-2 P.i p-1))) (u/l (list* p-1 ps))] [(list ps ...) (cons P.i ps)])))))
;; Initialize U and L as empty lists. ;; The lists will hold the vertices of upper and lower hulls respectively. (let ((U (upper/lower-hull (reverse P))) (L (upper/lower-hull P))) ;; Remove the last point of each list (it's the same as the first point of the other list). ;; Concatenate L and U to obtain the convex hull of P. (append (drop-right L 1) (drop-right U 1)))) ; Points in the result will be listed in CCW order.)
(module+ test
(require typed/rackunit) (check-equal? (convex-hull (list '(16 . 3) '(12 . 17) '(0 . 6) '(-4 . -6) '(16 . 6) '(16 . -7) '(16 . -3) '(17 . -4) '(5 . 19) '(19 . -8) '(3 . 16) '(12 . 13) '(3 . -4) '(17 . 5) '(-3 . 15) '(-3 . -9) '(0 . 11) '(-9 . -3) '(-4 . -2) '(12 . 10))) (list '(-9 . -3) '(-3 . -9) '(19 . -8) '(17 . 5) '(12 . 17) '(5 . 19) '(-3 . 15))))</lang>
- Output:
silence implies tests pass (and output is as expected)
zkl
<lang zkl>// Use Graham Scan to sort points into a convex hull // https://en.wikipedia.org/wiki/Graham_scan, O(n log n) // http://www.geeksforgeeks.org/convex-hull-set-2-graham-scan/ // http://geomalgorithms.com/a10-_hull-1.html fcn grahamScan(points){
N:=points.len(); # find the point with the lowest y-coordinate, x is tie breaker p0:=points.reduce(fcn([(a,b)]ab,[(x,y)]xy){
if(b<y)ab else if(b==y and a<x)ab else xy });
#sort points by polar angle with p0, ie ccw from p0 points.sort('wrap(p1,p2){ ccw(p0,p1,p2)>0 });
# We want points[0] to be a sentinel point that will stop the loop. points.insert(0,points[-1]); M:=1; # M will denote the number of points on the convex hull. foreach i in ([2..N]){ # Find next valid point on convex hull. while(ccw(points[M-1], points[M], points[i])<=0){
if(M>1) M-=1; else if(i==N) break; # All points are collinear else i+=1;
} points.swap(M+=1,i); # Update M and swap points[i] to the correct place. } points[0,M]
}
- Three points are a counter-clockwise turn if ccw > 0, clockwise if
- ccw < 0, and collinear if ccw = 0 because ccw is a determinant that
- gives twice the signed area of the triangle formed by p1, p2 and p3.
fcn ccw(a,b,c){ // a,b,c are points: (x,y)
((b[0] - a[0])*(c[1] - a[1])) - ((b[1] - a[1])*(c[0] - a[0]))
}</lang> <lang zkl>pts:=List( T(16,3), T(12,17), T(0,6), T(-4,-6), T(16,6), T(16, -7), T(16,-3),T(17,-4), T(5,19), T(19,-8), T(3,16), T(12,13), T(3,-4), T(17,5), T(-3,15), T(-3,-9), T(0,11), T(-9,-3), T(-4,-2), T(12,10), ) .apply(fcn(xy){ xy.apply("toFloat") }).copy(); hull:=grahamScan(pts); println("Convex Hull (%d points): %s".fmt(hull.len(),hull.toString(*)));</lang>
- Output:
Convex Hull (7 points): L(L(-3,-9),L(19,-8),L(17,5),L(12,17),L(5,19),L(-3,15),L(-9,-3))