// 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 (
	"fmt"
	"log"
	"math"
	"os"

	svg "github.com/ajstarks/svgo"
)

type Triangulation struct {
	Points     []Vertex
	ConvexHull []Vertex
	Triangles  []int32
	Halfedges  []int32
}

// Triangulate returns a Delaunay triangulation of the provided points.
func Triangulate(points []Vertex) (*Triangulation, error) {
	t := newTriangulator(points)
	err := t.triangulate()
	return &Triangulation{points, t.convexHull(), t.triangles, t.halfedges}, err
}

func (t *Triangulation) ConcaveHull(tooLong float64) []int32 {

	tooLong *= tooLong

	var candidates []int32

	for i, num := 0, len(t.Triangles)/3; i < num; i++ {
		idx := i * 3
		var max float64
		vs := t.Triangles[idx : idx+3]
		for j, vj := range vs {
			k := (j + 1) % 3
			if l := t.Points[vj].SquaredDistance2D(t.Points[vs[k]]); l > max {
				max = l
			}
		}
		if max > tooLong {
			candidates = append(candidates, int32(i))
		}
	}

	removed := map[int32]bool{}

	isBorder := func(n int32) bool {
		n *= 3
		for i := int32(0); i < 3; i++ {
			e := n + i
			if o := t.Halfedges[e]; o < 0 || removed[o/3] {
				return true
			}
		}
		return false
	}

	var newCandidates []int32

	for len(candidates) > 0 {
		log.Printf("info: candidates: %d\n", len(candidates))

		oldRemoved := len(removed)

		for _, i := range candidates {

			if isBorder(i) {
				removed[i] = true
			} else {
				newCandidates = append(newCandidates, i)
			}
		}

		if oldRemoved == len(removed) {
			break
		}

		candidates = newCandidates
		newCandidates = newCandidates[:0]
	}

	log.Printf("info: candidates left: %d\n", len(candidates))
	log.Printf("info: triangles: %d\n", len(t.Triangles)/3)
	log.Printf("info: triangles to remove: %d\n", len(removed))

	var edges [][2]int32

	for i, num := int32(0), int32(len(t.Triangles)/3); i < num; i++ {
		if removed[i] {
			continue
		}
		n := i * 3
		for j := int32(0); j < 3; j++ {
			e := n + j
			if t.Halfedges[e] < 0 || removed[e/3] {
				edges = append(edges, [2]int32{
					t.Triangles[e],
					t.Triangles[n+(j+1)%3],
				})
			}
		}
	}

	for i := range edges {
		fmt.Printf("%d - %d\n", edges[i][0], edges[i][1])
	}

	log.Printf("num of border triangles: %d\n", len(edges))
	log.Printf("len convex hull: %d\n", len(t.ConvexHull))

	var (
		minX float64 = math.MaxFloat64
		minY float64 = math.MaxFloat64
		maxX float64 = -math.MaxFloat64
		maxY float64 = -math.MaxFloat64
	)

	for i := range edges {
		p0 := t.Points[edges[i][0]]
		p1 := t.Points[edges[i][1]]

		minX = math.Min(p0.X, minX)
		maxX = math.Max(p0.X, maxX)
		minY = math.Min(p0.Y, minY)
		maxY = math.Max(p0.Y, maxY)

		minX = math.Min(p1.X, minX)
		maxX = math.Max(p1.X, maxX)
		minY = math.Min(p1.Y, minY)
		maxY = math.Max(p1.Y, maxY)
	}
	linear := func(x1, y1, x2, y2 float64) func(float64) float64 {
		// y1 = x1*m + b
		// y2 = x2*m + b
		// b = y1 - x1*m
		// y1-y2 = (x1-x2)*m
		// m = (y1-y2)/(x1-x2) for x1 != x2

		if x1 == x2 {
			return func(float64) float64 { return (y1 + y2) * 0.5 }
		}

		m := (y1 - y2) / (x1 - x2)
		b := y1 - x1*m

		return func(x float64) float64 { return m*x + b }
	}
	xf := linear(minX, 0, maxX, 1499)
	yf := linear(minY, 0, maxY, 1499)

	f, err := os.Create("/tmp/hull.svg")
	if err == nil {

		canvas := svg.New(f)
		canvas.Start(1500, 1500)

		for i := range edges {
			p0 := t.Points[edges[i][0]]
			p1 := t.Points[edges[i][1]]
			x1 := int(math.Floor(xf(p0.X)))
			y1 := int(math.Floor(yf(p0.Y)))
			x2 := int(math.Floor(xf(p1.X)))
			y2 := int(math.Floor(yf(p1.Y)))
			canvas.Line(x1, y1, x2, y2, "stroke:black")
		}

		canvas.End()
		f.Close()
	}
	return nil
}

func (t *Triangulation) TriangleSlices() [][]int32 {
	sl := make([][]int32, len(t.Triangles)/3)
	var j int
	for i := range sl {
		sl[i] = t.Triangles[j : j+3]
		j += 3
	}
	return sl
}

func (t *Triangulation) Tin() *Tin {

	min := Vertex{math.MaxFloat64, math.MaxFloat64, math.MaxFloat64}
	max := Vertex{-math.MaxFloat64, -math.MaxFloat64, -math.MaxFloat64}

	vertices := t.Points

	for _, v := range vertices {
		min.Minimize(v)
		max.Maximize(v)
	}

	return &Tin{
		Vertices:  vertices,
		Triangles: t.TriangleSlices(),
		Min:       min,
		Max:       max,
	}
}

func (t *Triangulation) area() float64 {
	var result float64
	points := t.Points
	ts := t.Triangles
	for i := 0; i < len(ts); i += 3 {
		p0 := points[ts[i+0]]
		p1 := points[ts[i+1]]
		p2 := points[ts[i+2]]
		result += area(p0, p1, p2)
	}
	return result / 2
}

// Validate performs several sanity checks on the Triangulation to check for
// potential errors. Returns nil if no issues were found. You normally
// shouldn't need to call this function but it can be useful for debugging.
func (t *Triangulation) Validate() error {
	// verify halfedges
	for i1, i2 := range t.Halfedges {
		if i1 != -1 && t.Halfedges[i1] != i2 {
			return fmt.Errorf("invalid halfedge connection")
		}
		if i2 != -1 && t.Halfedges[i2] != int32(i1) {
			return fmt.Errorf("invalid halfedge connection")
		}
	}

	// verify convex hull area vs sum of triangle areas
	hull1 := t.ConvexHull
	hull2 := ConvexHull(t.Points)
	area1 := polygonArea(hull1)
	area2 := polygonArea(hull2)
	area3 := t.area()
	if math.Abs(area1-area2) > 1e-9 || math.Abs(area1-area3) > 1e-9 {
		return fmt.Errorf("hull areas disagree: %f, %f, %f", area1, area2, area3)
	}

	// verify convex hull perimeter
	perimeter1 := polygonPerimeter(hull1)
	perimeter2 := polygonPerimeter(hull2)
	if math.Abs(perimeter1-perimeter2) > 1e-9 {
		return fmt.Errorf("hull perimeters disagree: %f, %f", perimeter1, perimeter2)
	}

	return nil
}