view pkg/octree/hull.go @ 2483:620038ade708 octree-diff

Incorporated fogleman's fast Delaunay triangulation adjuted to our vertex model. License: MIT Home: https://github.com/fogleman/delaunay
author Sascha L. Teichmann <sascha.teichmann@intevation.de>
date Fri, 01 Mar 2019 15:33:27 +0100
parents
children
line wrap: on
line source

// Copyright (C) 2018 Michael Fogleman
//
// Permission is hereby granted, free of charge, to any person obtaining
// a copy of this software and associated documentation files (the "Software"),
// to deal in the Software without restriction, including without limitation
// the rights to use, copy, modify, merge, publish, distribute, sublicense,
// and/or sell copies of the Software, and to permit persons to whom the
// Software is furnished to do so, subject to the following conditions:
//
// The above copyright notice and this permission notice shall be included
// in all copies or substantial portions of the Software.
//
// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
// OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
// FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
// THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
// LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
// OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.

package octree

import "sort"

func cross2D(p, a, b Vertex) float64 {
	return (a.X-p.X)*(b.Y-p.Y) - (a.Y-p.Y)*(b.X-p.X)
}

// ConvexHull returns the convex hull of the provided points.
func ConvexHull(points []Vertex) []Vertex {
	// copy points
	pointsCopy := make([]Vertex, len(points))
	copy(pointsCopy, points)
	points = pointsCopy

	// sort points
	sort.Slice(points, func(i, j int) bool {
		a := points[i]
		b := points[j]
		if a.X != b.X {
			return a.X < b.X
		}
		return a.Y < b.Y
	})

	// filter nearly-duplicate points
	distinctPoints := points[:0]
	for i, p := range points {
		if i > 0 && p.SquaredDistance2D(points[i-1]) < eps {
			continue
		}
		distinctPoints = append(distinctPoints, p)
	}
	points = distinctPoints

	// find upper and lower portions
	var U, L []Vertex
	for _, p := range points {
		for len(U) > 1 && cross2D(U[len(U)-2], U[len(U)-1], p) > 0 {
			U = U[:len(U)-1]
		}
		for len(L) > 1 && cross2D(L[len(L)-2], L[len(L)-1], p) < 0 {
			L = L[:len(L)-1]
		}
		U = append(U, p)
		L = append(L, p)
	}

	// reverse upper portion
	for i, j := 0, len(U)-1; i < j; i, j = i+1, j-1 {
		U[i], U[j] = U[j], U[i]
	}

	// construct complete hull
	if len(U) > 0 {
		U = U[:len(U)-1]
	}
	if len(L) > 0 {
		L = L[:len(L)-1]
	}
	return append(L, U...)
}